Free cookie consent management tool by TermsFeed Policy Generator

source: stable/HeuristicLab.Problems.DataAnalysis.Symbolic/3.4/Interpreter/SymbolicDataAnalysisExpressionTreeILEmittingInterpreter.cs @ 17108

Last change on this file since 17108 was 17105, checked in by mkommend, 5 years ago

#2520: Merged 16584, 16585,16594,16595, 16625, 16658, 16659, 16672, 16707, 16729, 16792, 16796, 16797, 16799, 16819, 16906, 16907, 16908, 16933, 16945, 16992, 16994, 16995, 16996, 16997, 17014, 17015, 17017, 17020, 17021, 17022, 17023, 17024, 17029, 17086, 17087, 17088, 17089 into stable.
File size: 39.1 KB
Line 
1#region License Information
2/* HeuristicLab
3 * Copyright (C) 2002-2019 Heuristic and Evolutionary Algorithms Laboratory (HEAL)
4 *
5 * This file is part of HeuristicLab.
6 *
7 * HeuristicLab is free software: you can redistribute it and/or modify
8 * it under the terms of the GNU General Public License as published by
9 * the Free Software Foundation, either version 3 of the License, or
10 * (at your option) any later version.
11 *
12 * HeuristicLab is distributed in the hope that it will be useful,
13 * but WITHOUT ANY WARRANTY; without even the implied warranty of
14 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
15 * GNU General Public License for more details.
16 *
17 * You should have received a copy of the GNU General Public License
18 * along with HeuristicLab. If not, see <http://www.gnu.org/licenses/>.
19 */
20#endregion
21
22using System;
23using System.Collections.Generic;
24using System.Linq;
25using System.Reflection;
26using System.Reflection.Emit;
27using HeuristicLab.Common;
28using HeuristicLab.Core;
29using HeuristicLab.Data;
30using HeuristicLab.Encodings.SymbolicExpressionTreeEncoding;
31using HeuristicLab.Parameters;
32using HEAL.Attic;
33
34namespace HeuristicLab.Problems.DataAnalysis.Symbolic {
35  [StorableType("426718E3-2A57-4CA4-98A1-65EDD0B0BDBF")]
36  [Item("SymbolicDataAnalysisExpressionTreeILEmittingInterpreter", "Interpreter for symbolic expression trees.")]
37  public sealed class SymbolicDataAnalysisExpressionTreeILEmittingInterpreter : ParameterizedNamedItem, ISymbolicDataAnalysisExpressionTreeInterpreter {
38    private static readonly Type thisType = typeof(SymbolicDataAnalysisExpressionTreeILEmittingInterpreter);
39    internal delegate double CompiledFunction(int sampleIndex, IList<double>[] columns);
40
41    #region method infos
42    private static MethodInfo listGetValue = typeof(IList<double>).GetProperty("Item", new Type[] { typeof(int) }).GetGetMethod();
43
44    private static MethodInfo cos = typeof(Math).GetMethod("Cos", new Type[] { typeof(double) });
45    private static MethodInfo sin = typeof(Math).GetMethod("Sin", new Type[] { typeof(double) });
46    private static MethodInfo tan = typeof(Math).GetMethod("Tan", new Type[] { typeof(double) });
47    private static MethodInfo tanh = typeof(Math).GetMethod("Tanh", new Type[] { typeof(double) });
48    private static MethodInfo exp = typeof(Math).GetMethod("Exp", new Type[] { typeof(double) });
49    private static MethodInfo log = typeof(Math).GetMethod("Log", new Type[] { typeof(double) });
50    private static MethodInfo power = typeof(Math).GetMethod("Pow", new Type[] { typeof(double), typeof(double) });
51    private static MethodInfo round = typeof(Math).GetMethod("Round", new Type[] { typeof(double) });
52    private static MethodInfo sqrt = typeof(Math).GetMethod("Sqrt", new Type[] { typeof(double) });
53    private static MethodInfo abs = typeof(Math).GetMethod("Abs", new Type[] { typeof(double) });
54
55    private static MethodInfo airyA = thisType.GetMethod("AiryA", new Type[] { typeof(double) });
56    private static MethodInfo airyB = thisType.GetMethod("AiryB", new Type[] { typeof(double) });
57    private static MethodInfo gamma = thisType.GetMethod("Gamma", new Type[] { typeof(double) });
58    private static MethodInfo psi = thisType.GetMethod("Psi", new Type[] { typeof(double) });
59    private static MethodInfo dawson = thisType.GetMethod("Dawson", new Type[] { typeof(double) });
60    private static MethodInfo expIntegralEi = thisType.GetMethod("ExpIntegralEi", new Type[] { typeof(double) });
61    private static MethodInfo sinIntegral = thisType.GetMethod("SinIntegral", new Type[] { typeof(double) });
62    private static MethodInfo cosIntegral = thisType.GetMethod("CosIntegral", new Type[] { typeof(double) });
63    private static MethodInfo hypSinIntegral = thisType.GetMethod("HypSinIntegral", new Type[] { typeof(double) });
64    private static MethodInfo hypCosIntegral = thisType.GetMethod("HypCosIntegral", new Type[] { typeof(double) });
65    private static MethodInfo fresnelCosIntegral = thisType.GetMethod("FresnelCosIntegral", new Type[] { typeof(double) });
66    private static MethodInfo fresnelSinIntegral = thisType.GetMethod("FresnelSinIntegral", new Type[] { typeof(double) });
67    private static MethodInfo norm = thisType.GetMethod("Norm", new Type[] { typeof(double) });
68    private static MethodInfo erf = thisType.GetMethod("Erf", new Type[] { typeof(double) });
69    private static MethodInfo bessel = thisType.GetMethod("Bessel", new Type[] { typeof(double) });
70    private static MethodInfo string_eq = typeof(string).GetMethod("Equals", new Type[] { typeof(string) });
71    #endregion
72
73    private const string CheckExpressionsWithIntervalArithmeticParameterName = "CheckExpressionsWithIntervalArithmetic";
74    private const string CheckExpressionsWithIntervalArithmeticParameterDescription = "Switch that determines if the interpreter checks the validity of expressions with interval arithmetic before evaluating the expression.";
75    private const string EvaluatedSolutionsParameterName = "EvaluatedSolutions";
76
77    public override bool CanChangeName {
78      get { return false; }
79    }
80
81    public override bool CanChangeDescription {
82      get { return false; }
83    }
84
85    #region parameter properties
86    public IFixedValueParameter<BoolValue> CheckExpressionsWithIntervalArithmeticParameter {
87      get { return (IFixedValueParameter<BoolValue>)Parameters[CheckExpressionsWithIntervalArithmeticParameterName]; }
88    }
89
90    public IFixedValueParameter<IntValue> EvaluatedSolutionsParameter {
91      get { return (IFixedValueParameter<IntValue>)Parameters[EvaluatedSolutionsParameterName]; }
92    }
93    #endregion
94
95    #region properties
96    public bool CheckExpressionsWithIntervalArithmetic {
97      get { return CheckExpressionsWithIntervalArithmeticParameter.Value.Value; }
98      set { CheckExpressionsWithIntervalArithmeticParameter.Value.Value = value; }
99    }
100    public int EvaluatedSolutions {
101      get { return EvaluatedSolutionsParameter.Value.Value; }
102      set { EvaluatedSolutionsParameter.Value.Value = value; }
103    }
104    #endregion
105
106    [StorableConstructor]
107    private SymbolicDataAnalysisExpressionTreeILEmittingInterpreter(StorableConstructorFlag _) : base(_) { }
108
109    private SymbolicDataAnalysisExpressionTreeILEmittingInterpreter(SymbolicDataAnalysisExpressionTreeILEmittingInterpreter original, Cloner cloner) : base(original, cloner) { }
110    public override IDeepCloneable Clone(Cloner cloner) {
111      return new SymbolicDataAnalysisExpressionTreeILEmittingInterpreter(this, cloner);
112    }
113
114    public SymbolicDataAnalysisExpressionTreeILEmittingInterpreter()
115      : base("SymbolicDataAnalysisExpressionTreeILEmittingInterpreter", "Interpreter for symbolic expression trees.") {
116      Parameters.Add(new FixedValueParameter<BoolValue>(CheckExpressionsWithIntervalArithmeticParameterName,
117        "Switch that determines if the interpreter checks the validity of expressions with interval arithmetic before evaluating the expression.", new BoolValue(false)));
118      Parameters.Add(new FixedValueParameter<IntValue>(EvaluatedSolutionsParameterName, "A counter for the total number of solutions the interpreter has evaluated", new IntValue(0)));
119    }
120
121    public SymbolicDataAnalysisExpressionTreeILEmittingInterpreter(string name, string description)
122      : base(name, description) {
123      Parameters.Add(new FixedValueParameter<BoolValue>(CheckExpressionsWithIntervalArithmeticParameterName,
124        "Switch that determines if the interpreter checks the validity of expressions with interval arithmetic before evaluating the expression.", new BoolValue(false)));
125      Parameters.Add(new FixedValueParameter<IntValue>(EvaluatedSolutionsParameterName, "A counter for the total number of solutions the interpreter has evaluated", new IntValue(0)));
126    }
127
128    [StorableHook(HookType.AfterDeserialization)]
129    private void AfterDeserialization() {
130      var evaluatedSolutions = new IntValue(0);
131      var checkExpressionsWithIntervalArithmetic = new BoolValue(false);
132      if (Parameters.ContainsKey(EvaluatedSolutionsParameterName)) {
133        var evaluatedSolutionsParameter = (IValueParameter<IntValue>)Parameters[EvaluatedSolutionsParameterName];
134        evaluatedSolutions = evaluatedSolutionsParameter.Value;
135        Parameters.Remove(EvaluatedSolutionsParameterName);
136      }
137      Parameters.Add(new FixedValueParameter<IntValue>(EvaluatedSolutionsParameterName, "A counter for the total number of solutions the interpreter has evaluated", evaluatedSolutions));
138      if (Parameters.ContainsKey(CheckExpressionsWithIntervalArithmeticParameterName)) {
139        var checkExpressionsWithIntervalArithmeticParameter = (IValueParameter<BoolValue>)Parameters[CheckExpressionsWithIntervalArithmeticParameterName];
140        Parameters.Remove(CheckExpressionsWithIntervalArithmeticParameterName);
141        checkExpressionsWithIntervalArithmetic = checkExpressionsWithIntervalArithmeticParameter.Value;
142      }
143      Parameters.Add(new FixedValueParameter<BoolValue>(CheckExpressionsWithIntervalArithmeticParameterName, CheckExpressionsWithIntervalArithmeticParameterDescription, checkExpressionsWithIntervalArithmetic));
144    }
145
146    #region IStatefulItem
147    public void InitializeState() {
148      EvaluatedSolutions = 0;
149    }
150
151    public void ClearState() {
152    }
153    #endregion
154
155    private readonly object syncRoot = new object();
156    public IEnumerable<double> GetSymbolicExpressionTreeValues(ISymbolicExpressionTree tree, IDataset dataset, IEnumerable<int> rows) {
157      if (CheckExpressionsWithIntervalArithmetic)
158        throw new NotSupportedException("Interval arithmetic is not yet supported in the symbolic data analysis interpreter.");
159
160      lock (syncRoot) {
161        EvaluatedSolutions++; // increment the evaluated solutions counter
162      }
163      var state = PrepareInterpreterState(tree, dataset);
164
165      Type[] methodArgs = { typeof(int), typeof(IList<double>[]) };
166      DynamicMethod testFun = new DynamicMethod("TestFun", typeof(double), methodArgs, typeof(SymbolicDataAnalysisExpressionTreeILEmittingInterpreter).Module);
167
168      ILGenerator il = testFun.GetILGenerator();
169      CompileInstructions(il, state, dataset);
170      il.Emit(System.Reflection.Emit.OpCodes.Conv_R8);
171      il.Emit(System.Reflection.Emit.OpCodes.Ret);
172      var function = (CompiledFunction)testFun.CreateDelegate(typeof(CompiledFunction));
173
174      IList<double>[] columns = dataset.DoubleVariables.Select(v => dataset.GetReadOnlyDoubleValues(v)).ToArray();
175
176      foreach (var row in rows) {
177        yield return function(row, columns);
178      }
179    }
180
181    private InterpreterState PrepareInterpreterState(ISymbolicExpressionTree tree, IDataset dataset) {
182      Instruction[] code = SymbolicExpressionTreeCompiler.Compile(tree, OpCodes.MapSymbolToOpCode);
183      Dictionary<string, int> doubleVariableNames = dataset.DoubleVariables.Select((x, i) => new { x, i }).ToDictionary(e => e.x, e => e.i);
184      int necessaryArgStackSize = 0;
185      foreach (Instruction instr in code) {
186        if (instr.opCode == OpCodes.Variable) {
187          var variableTreeNode = (VariableTreeNode)instr.dynamicNode;
188          instr.data = doubleVariableNames[variableTreeNode.VariableName];
189        } else if (instr.opCode == OpCodes.LagVariable) {
190          var laggedVariableTreeNode = (LaggedVariableTreeNode)instr.dynamicNode;
191          instr.data = doubleVariableNames[laggedVariableTreeNode.VariableName];
192        } else if (instr.opCode == OpCodes.VariableCondition) {
193          var variableConditionTreeNode = (VariableConditionTreeNode)instr.dynamicNode;
194          instr.data = doubleVariableNames[variableConditionTreeNode.VariableName];
195        } else if (instr.opCode == OpCodes.Call) {
196          necessaryArgStackSize += instr.nArguments + 1;
197        }
198      }
199      return new InterpreterState(code, necessaryArgStackSize);
200    }
201
202    private void CompileInstructions(ILGenerator il, InterpreterState state, IDataset ds) {
203      Instruction currentInstr = state.NextInstruction();
204      int nArgs = currentInstr.nArguments;
205
206      switch (currentInstr.opCode) {
207        case OpCodes.Add: {
208            if (nArgs > 0) {
209              CompileInstructions(il, state, ds);
210            }
211            for (int i = 1; i < nArgs; i++) {
212              CompileInstructions(il, state, ds);
213              il.Emit(System.Reflection.Emit.OpCodes.Add);
214            }
215            return;
216          }
217        case OpCodes.Sub: {
218            if (nArgs == 1) {
219              CompileInstructions(il, state, ds);
220              il.Emit(System.Reflection.Emit.OpCodes.Neg);
221              return;
222            }
223            if (nArgs > 0) {
224              CompileInstructions(il, state, ds);
225            }
226            for (int i = 1; i < nArgs; i++) {
227              CompileInstructions(il, state, ds);
228              il.Emit(System.Reflection.Emit.OpCodes.Sub);
229            }
230            return;
231          }
232        case OpCodes.Mul: {
233            if (nArgs > 0) {
234              CompileInstructions(il, state, ds);
235            }
236            for (int i = 1; i < nArgs; i++) {
237              CompileInstructions(il, state, ds);
238              il.Emit(System.Reflection.Emit.OpCodes.Mul);
239            }
240            return;
241          }
242        case OpCodes.Div: {
243            if (nArgs == 1) {
244              il.Emit(System.Reflection.Emit.OpCodes.Ldc_R8, 1.0);
245              CompileInstructions(il, state, ds);
246              il.Emit(System.Reflection.Emit.OpCodes.Div);
247              return;
248            }
249            if (nArgs > 0) {
250              CompileInstructions(il, state, ds);
251            }
252            for (int i = 1; i < nArgs; i++) {
253              CompileInstructions(il, state, ds);
254              il.Emit(System.Reflection.Emit.OpCodes.Div);
255            }
256            return;
257          }
258        case OpCodes.Average: {
259            CompileInstructions(il, state, ds);
260            for (int i = 1; i < nArgs; i++) {
261              CompileInstructions(il, state, ds);
262              il.Emit(System.Reflection.Emit.OpCodes.Add);
263            }
264            il.Emit(System.Reflection.Emit.OpCodes.Ldc_R8, (double)nArgs);
265            il.Emit(System.Reflection.Emit.OpCodes.Div);
266            return;
267          }
268        case OpCodes.Absolute: {
269            CompileInstructions(il, state, ds);
270            il.Emit(System.Reflection.Emit.OpCodes.Call, abs);
271            return;
272          }
273        case OpCodes.Cos: {
274            CompileInstructions(il, state, ds);
275            il.Emit(System.Reflection.Emit.OpCodes.Call, cos);
276            return;
277          }
278        case OpCodes.Sin: {
279            CompileInstructions(il, state, ds);
280            il.Emit(System.Reflection.Emit.OpCodes.Call, sin);
281            return;
282          }
283        case OpCodes.Tan: {
284            CompileInstructions(il, state, ds);
285            il.Emit(System.Reflection.Emit.OpCodes.Call, tan);
286            return;
287          }
288        case OpCodes.Tanh: {
289            CompileInstructions(il, state, ds);
290            il.Emit(System.Reflection.Emit.OpCodes.Call, tanh);
291            return;
292          }
293        case OpCodes.Power: {
294            CompileInstructions(il, state, ds);
295            CompileInstructions(il, state, ds);
296            il.Emit(System.Reflection.Emit.OpCodes.Call, round);
297            il.Emit(System.Reflection.Emit.OpCodes.Call, power);
298            return;
299          }
300        case OpCodes.Root: {
301            CompileInstructions(il, state, ds);
302            il.Emit(System.Reflection.Emit.OpCodes.Ldc_R8, 1.0); // 1 / round(...)
303            CompileInstructions(il, state, ds);
304            il.Emit(System.Reflection.Emit.OpCodes.Call, round);
305            il.Emit(System.Reflection.Emit.OpCodes.Div);
306            il.Emit(System.Reflection.Emit.OpCodes.Call, power);
307            return;
308          }
309        case OpCodes.Exp: {
310            CompileInstructions(il, state, ds);
311            il.Emit(System.Reflection.Emit.OpCodes.Call, exp);
312            return;
313          }
314        case OpCodes.Log: {
315            CompileInstructions(il, state, ds);
316            il.Emit(System.Reflection.Emit.OpCodes.Call, log);
317            return;
318          }
319        case OpCodes.Square: {
320            CompileInstructions(il, state, ds);
321            il.Emit(System.Reflection.Emit.OpCodes.Ldc_R8, 2.0);
322            il.Emit(System.Reflection.Emit.OpCodes.Call, power);
323            return;
324          }
325        case OpCodes.Cube: {
326            CompileInstructions(il, state, ds);
327            il.Emit(System.Reflection.Emit.OpCodes.Ldc_R8, 3.0);
328            il.Emit(System.Reflection.Emit.OpCodes.Call, power);
329            return;
330          }
331        case OpCodes.SquareRoot: {
332            CompileInstructions(il, state, ds);
333            il.Emit(System.Reflection.Emit.OpCodes.Call, sqrt);
334            return;
335          }
336        case OpCodes.CubeRoot: {
337            CompileInstructions(il, state, ds);
338            var c1 = il.DefineLabel();
339            var end = il.DefineLabel();
340
341            il.Emit(System.Reflection.Emit.OpCodes.Dup); // x
342            il.Emit(System.Reflection.Emit.OpCodes.Ldc_R8, 0.0);
343            il.Emit(System.Reflection.Emit.OpCodes.Clt); // x < 0?
344            il.Emit(System.Reflection.Emit.OpCodes.Brfalse, c1);
345            il.Emit(System.Reflection.Emit.OpCodes.Neg); // x = -x
346            il.Emit(System.Reflection.Emit.OpCodes.Ldc_R8, 1.0 / 3.0);
347            il.Emit(System.Reflection.Emit.OpCodes.Call, power);
348            il.Emit(System.Reflection.Emit.OpCodes.Neg); // -Math.Pow(-x, 1/3)
349            il.Emit(System.Reflection.Emit.OpCodes.Br, end);
350            il.MarkLabel(c1);
351            il.Emit(System.Reflection.Emit.OpCodes.Ldc_R8, 1.0 / 3.0);
352            il.Emit(System.Reflection.Emit.OpCodes.Call, power);
353            il.MarkLabel(end);
354            return;
355          }
356        case OpCodes.AiryA: {
357            CompileInstructions(il, state, ds);
358            il.Emit(System.Reflection.Emit.OpCodes.Call, airyA);
359            return;
360          }
361        case OpCodes.AiryB: {
362            CompileInstructions(il, state, ds);
363            il.Emit(System.Reflection.Emit.OpCodes.Call, airyB);
364            return;
365          }
366        case OpCodes.Bessel: {
367            CompileInstructions(il, state, ds);
368            il.Emit(System.Reflection.Emit.OpCodes.Call, bessel);
369            return;
370          }
371        case OpCodes.CosineIntegral: {
372            CompileInstructions(il, state, ds);
373            il.Emit(System.Reflection.Emit.OpCodes.Call, cosIntegral);
374            return;
375          }
376        case OpCodes.Dawson: {
377            CompileInstructions(il, state, ds);
378            il.Emit(System.Reflection.Emit.OpCodes.Call, dawson);
379            return;
380          }
381        case OpCodes.Erf: {
382            CompileInstructions(il, state, ds);
383            il.Emit(System.Reflection.Emit.OpCodes.Call, erf);
384            return;
385          }
386        case OpCodes.ExponentialIntegralEi: {
387            CompileInstructions(il, state, ds);
388            il.Emit(System.Reflection.Emit.OpCodes.Call, expIntegralEi);
389            return;
390          }
391        case OpCodes.FresnelCosineIntegral: {
392            CompileInstructions(il, state, ds);
393            il.Emit(System.Reflection.Emit.OpCodes.Call, fresnelCosIntegral);
394            return;
395          }
396        case OpCodes.FresnelSineIntegral: {
397            CompileInstructions(il, state, ds);
398            il.Emit(System.Reflection.Emit.OpCodes.Call, fresnelSinIntegral);
399            return;
400          }
401        case OpCodes.Gamma: {
402            CompileInstructions(il, state, ds);
403            il.Emit(System.Reflection.Emit.OpCodes.Call, gamma);
404            return;
405          }
406        case OpCodes.HyperbolicCosineIntegral: {
407            CompileInstructions(il, state, ds);
408            il.Emit(System.Reflection.Emit.OpCodes.Call, hypCosIntegral);
409            return;
410          }
411        case OpCodes.HyperbolicSineIntegral: {
412            CompileInstructions(il, state, ds);
413            il.Emit(System.Reflection.Emit.OpCodes.Call, hypSinIntegral);
414            return;
415          }
416        case OpCodes.Norm: {
417            CompileInstructions(il, state, ds);
418            il.Emit(System.Reflection.Emit.OpCodes.Call, norm);
419            return;
420          }
421        case OpCodes.Psi: {
422            CompileInstructions(il, state, ds);
423            il.Emit(System.Reflection.Emit.OpCodes.Call, psi);
424            return;
425          }
426        case OpCodes.SineIntegral: {
427            CompileInstructions(il, state, ds);
428            il.Emit(System.Reflection.Emit.OpCodes.Call, sinIntegral);
429            return;
430          }
431        case OpCodes.AnalyticQuotient: {
432            CompileInstructions(il, state, ds); // x1
433            CompileInstructions(il, state, ds); // x2
434
435            il.Emit(System.Reflection.Emit.OpCodes.Dup);
436            il.Emit(System.Reflection.Emit.OpCodes.Mul); // x2*x2
437            il.Emit(System.Reflection.Emit.OpCodes.Ldc_R8, 1.0);
438            il.Emit(System.Reflection.Emit.OpCodes.Add); // 1+x2*x2
439            il.Emit(System.Reflection.Emit.OpCodes.Call, sqrt);
440            il.Emit(System.Reflection.Emit.OpCodes.Div);
441            return;
442          }
443        case OpCodes.IfThenElse: {
444            Label end = il.DefineLabel();
445            Label c1 = il.DefineLabel();
446            CompileInstructions(il, state, ds);
447            il.Emit(System.Reflection.Emit.OpCodes.Ldc_R8, 0.0); // > 0
448            il.Emit(System.Reflection.Emit.OpCodes.Cgt);
449            il.Emit(System.Reflection.Emit.OpCodes.Brfalse, c1);
450            CompileInstructions(il, state, ds);
451            il.Emit(System.Reflection.Emit.OpCodes.Br, end);
452            il.MarkLabel(c1);
453            CompileInstructions(il, state, ds);
454            il.MarkLabel(end);
455            return;
456          }
457        case OpCodes.AND: {
458            Label falseBranch = il.DefineLabel();
459            Label end = il.DefineLabel();
460            CompileInstructions(il, state, ds);
461            for (int i = 1; i < nArgs; i++) {
462              il.Emit(System.Reflection.Emit.OpCodes.Ldc_R8, 0.0); // > 0
463              il.Emit(System.Reflection.Emit.OpCodes.Cgt);
464              il.Emit(System.Reflection.Emit.OpCodes.Brfalse, falseBranch);
465              CompileInstructions(il, state, ds);
466            }
467            il.Emit(System.Reflection.Emit.OpCodes.Ldc_R8, 0.0); // > 0
468            il.Emit(System.Reflection.Emit.OpCodes.Cgt);
469            il.Emit(System.Reflection.Emit.OpCodes.Brfalse, falseBranch);
470            il.Emit(System.Reflection.Emit.OpCodes.Ldc_R8, 1.0); // 1
471            il.Emit(System.Reflection.Emit.OpCodes.Br, end);
472            il.MarkLabel(falseBranch);
473            il.Emit(System.Reflection.Emit.OpCodes.Ldc_R8, 1.0); // -1
474            il.Emit(System.Reflection.Emit.OpCodes.Neg);
475            il.MarkLabel(end);
476            return;
477          }
478        case OpCodes.OR: {
479            Label trueBranch = il.DefineLabel();
480            Label end = il.DefineLabel();
481            Label resultBranch = il.DefineLabel();
482            CompileInstructions(il, state, ds);
483            for (int i = 1; i < nArgs; i++) {
484              Label nextArgBranch = il.DefineLabel();
485              // complex definition because of special properties of NaN 
486              il.Emit(System.Reflection.Emit.OpCodes.Dup);
487              il.Emit(System.Reflection.Emit.OpCodes.Ldc_R8, 0.0); // <= 0       
488              il.Emit(System.Reflection.Emit.OpCodes.Ble, nextArgBranch);
489              il.Emit(System.Reflection.Emit.OpCodes.Br, resultBranch);
490              il.MarkLabel(nextArgBranch);
491              il.Emit(System.Reflection.Emit.OpCodes.Pop);
492              CompileInstructions(il, state, ds);
493            }
494            il.MarkLabel(resultBranch);
495            il.Emit(System.Reflection.Emit.OpCodes.Ldc_R8, 0.0); // > 0
496            il.Emit(System.Reflection.Emit.OpCodes.Cgt);
497            il.Emit(System.Reflection.Emit.OpCodes.Brtrue, trueBranch);
498            il.Emit(System.Reflection.Emit.OpCodes.Ldc_R8, 1.0); // -1
499            il.Emit(System.Reflection.Emit.OpCodes.Neg);
500            il.Emit(System.Reflection.Emit.OpCodes.Br, end);
501            il.MarkLabel(trueBranch);
502            il.Emit(System.Reflection.Emit.OpCodes.Ldc_R8, 1.0); // 1
503            il.MarkLabel(end);
504            return;
505          }
506        case OpCodes.NOT: {
507            CompileInstructions(il, state, ds);
508            il.Emit(System.Reflection.Emit.OpCodes.Ldc_R8, 0.0); // > 0
509            il.Emit(System.Reflection.Emit.OpCodes.Cgt);
510            il.Emit(System.Reflection.Emit.OpCodes.Conv_R8); // convert to float64
511            il.Emit(System.Reflection.Emit.OpCodes.Ldc_R8, 2.0); // * 2
512            il.Emit(System.Reflection.Emit.OpCodes.Mul);
513            il.Emit(System.Reflection.Emit.OpCodes.Ldc_R8, 1.0); // - 1
514            il.Emit(System.Reflection.Emit.OpCodes.Sub);
515            il.Emit(System.Reflection.Emit.OpCodes.Neg); // * -1
516            return;
517          }
518        case OpCodes.XOR: {
519            CompileInstructions(il, state, ds);
520            il.Emit(System.Reflection.Emit.OpCodes.Ldc_R8, 0.0);
521            il.Emit(System.Reflection.Emit.OpCodes.Cgt);// > 0
522
523            for (int i = 1; i < nArgs; i++) {
524              CompileInstructions(il, state, ds);
525              il.Emit(System.Reflection.Emit.OpCodes.Ldc_R8, 0.0);
526              il.Emit(System.Reflection.Emit.OpCodes.Cgt);// > 0
527              il.Emit(System.Reflection.Emit.OpCodes.Xor);
528            }
529            il.Emit(System.Reflection.Emit.OpCodes.Conv_R8); // convert to float64
530
531            il.Emit(System.Reflection.Emit.OpCodes.Ldc_R8, 2.0); // * 2
532            il.Emit(System.Reflection.Emit.OpCodes.Mul);
533            il.Emit(System.Reflection.Emit.OpCodes.Ldc_R8, 1.0); // - 1
534            il.Emit(System.Reflection.Emit.OpCodes.Sub);
535            return;
536          }
537        case OpCodes.GT: {
538            CompileInstructions(il, state, ds);
539            CompileInstructions(il, state, ds);
540
541            il.Emit(System.Reflection.Emit.OpCodes.Cgt); // 1 (>) / 0 (otherwise)
542            il.Emit(System.Reflection.Emit.OpCodes.Conv_R8); // convert to float64
543            il.Emit(System.Reflection.Emit.OpCodes.Ldc_R8, 2.0); // * 2
544            il.Emit(System.Reflection.Emit.OpCodes.Mul);
545            il.Emit(System.Reflection.Emit.OpCodes.Ldc_R8, 1.0); // - 1
546            il.Emit(System.Reflection.Emit.OpCodes.Sub);
547            return;
548          }
549        case OpCodes.LT: {
550            CompileInstructions(il, state, ds);
551            CompileInstructions(il, state, ds);
552            il.Emit(System.Reflection.Emit.OpCodes.Clt);
553            il.Emit(System.Reflection.Emit.OpCodes.Conv_R8); // convert to float64
554            il.Emit(System.Reflection.Emit.OpCodes.Ldc_R8, 2.0); // * 2
555            il.Emit(System.Reflection.Emit.OpCodes.Mul);
556            il.Emit(System.Reflection.Emit.OpCodes.Ldc_R8, 1.0); // - 1
557            il.Emit(System.Reflection.Emit.OpCodes.Sub);
558            return;
559          }
560        case OpCodes.TimeLag: {
561            LaggedTreeNode laggedTreeNode = (LaggedTreeNode)currentInstr.dynamicNode;
562            il.Emit(System.Reflection.Emit.OpCodes.Ldarg_0); // row -= lag
563            il.Emit(System.Reflection.Emit.OpCodes.Ldc_I4, laggedTreeNode.Lag);
564            il.Emit(System.Reflection.Emit.OpCodes.Add);
565            il.Emit(System.Reflection.Emit.OpCodes.Starg, 0);
566            var prevLaggedContext = state.InLaggedContext;
567            state.InLaggedContext = true;
568            CompileInstructions(il, state, ds);
569            il.Emit(System.Reflection.Emit.OpCodes.Ldarg_0); // row += lag
570            il.Emit(System.Reflection.Emit.OpCodes.Ldc_I4, laggedTreeNode.Lag);
571            il.Emit(System.Reflection.Emit.OpCodes.Sub);
572            il.Emit(System.Reflection.Emit.OpCodes.Starg, 0);
573            state.InLaggedContext = prevLaggedContext;
574            return;
575          }
576        case OpCodes.Integral: {
577            int savedPc = state.ProgramCounter;
578            LaggedTreeNode laggedTreeNode = (LaggedTreeNode)currentInstr.dynamicNode;
579            il.Emit(System.Reflection.Emit.OpCodes.Ldarg_0); // row -= lag
580            il.Emit(System.Reflection.Emit.OpCodes.Ldc_I4, laggedTreeNode.Lag);
581            il.Emit(System.Reflection.Emit.OpCodes.Add);
582            il.Emit(System.Reflection.Emit.OpCodes.Starg, 0);
583            var prevLaggedContext = state.InLaggedContext;
584            state.InLaggedContext = true;
585            CompileInstructions(il, state, ds);
586            for (int l = laggedTreeNode.Lag; l < 0; l++) {
587              il.Emit(System.Reflection.Emit.OpCodes.Ldarg_0); // row += lag
588              il.Emit(System.Reflection.Emit.OpCodes.Ldc_I4_1);
589              il.Emit(System.Reflection.Emit.OpCodes.Add);
590              il.Emit(System.Reflection.Emit.OpCodes.Starg, 0);
591              state.ProgramCounter = savedPc;
592              CompileInstructions(il, state, ds);
593              il.Emit(System.Reflection.Emit.OpCodes.Add);
594            }
595            state.InLaggedContext = prevLaggedContext;
596            return;
597          }
598
599        //mkommend: derivate calculation taken from:
600        //http://www.holoborodko.com/pavel/numerical-methods/numerical-derivative/smooth-low-noise-differentiators/
601        //one sided smooth differentiatior, N = 4
602        // y' = 1/8h (f_i + 2f_i-1, -2 f_i-3 - f_i-4)
603        case OpCodes.Derivative: {
604            int savedPc = state.ProgramCounter;
605            CompileInstructions(il, state, ds);
606            il.Emit(System.Reflection.Emit.OpCodes.Ldarg_0); // row --
607            il.Emit(System.Reflection.Emit.OpCodes.Ldc_I4_M1);
608            il.Emit(System.Reflection.Emit.OpCodes.Add);
609            il.Emit(System.Reflection.Emit.OpCodes.Starg, 0);
610            state.ProgramCounter = savedPc;
611            var prevLaggedContext = state.InLaggedContext;
612            state.InLaggedContext = true;
613            CompileInstructions(il, state, ds);
614            il.Emit(System.Reflection.Emit.OpCodes.Ldc_R8, 2.0); // f_0 + 2 * f_1
615            il.Emit(System.Reflection.Emit.OpCodes.Mul);
616            il.Emit(System.Reflection.Emit.OpCodes.Add);
617
618            il.Emit(System.Reflection.Emit.OpCodes.Ldarg_0); // row -=2
619            il.Emit(System.Reflection.Emit.OpCodes.Ldc_I4_2);
620            il.Emit(System.Reflection.Emit.OpCodes.Sub);
621            il.Emit(System.Reflection.Emit.OpCodes.Starg, 0);
622            state.ProgramCounter = savedPc;
623            CompileInstructions(il, state, ds);
624            il.Emit(System.Reflection.Emit.OpCodes.Ldc_R8, 2.0); // f_0 + 2 * f_1 - 2 * f_3
625            il.Emit(System.Reflection.Emit.OpCodes.Mul);
626            il.Emit(System.Reflection.Emit.OpCodes.Sub);
627
628            il.Emit(System.Reflection.Emit.OpCodes.Ldarg_0); // row --
629            il.Emit(System.Reflection.Emit.OpCodes.Ldc_I4_M1);
630            il.Emit(System.Reflection.Emit.OpCodes.Add);
631            il.Emit(System.Reflection.Emit.OpCodes.Starg, 0);
632            state.ProgramCounter = savedPc;
633            CompileInstructions(il, state, ds);
634            il.Emit(System.Reflection.Emit.OpCodes.Sub); // f_0 + 2 * f_1 - 2 * f_3 - f_4
635            il.Emit(System.Reflection.Emit.OpCodes.Ldc_R8, 8.0); // / 8
636            il.Emit(System.Reflection.Emit.OpCodes.Div);
637
638            il.Emit(System.Reflection.Emit.OpCodes.Ldarg_0); // row +=4
639            il.Emit(System.Reflection.Emit.OpCodes.Ldc_I4_4);
640            il.Emit(System.Reflection.Emit.OpCodes.Add);
641            il.Emit(System.Reflection.Emit.OpCodes.Starg, 0);
642            state.InLaggedContext = prevLaggedContext;
643            return;
644          }
645        case OpCodes.Call: {
646            throw new NotSupportedException(
647              "Automatically defined functions are not supported by the SymbolicDataAnalysisTreeILEmittingInterpreter. Either turn of ADFs or change the interpeter.");
648          }
649        case OpCodes.Arg: {
650            throw new NotSupportedException(
651              "Automatically defined functions are not supported by the SymbolicDataAnalysisTreeILEmittingInterpreter. Either turn of ADFs or change the interpeter.");
652          }
653        case OpCodes.Variable: {
654            VariableTreeNode varNode = (VariableTreeNode)currentInstr.dynamicNode;
655            il.Emit(System.Reflection.Emit.OpCodes.Ldarg_1); // load columns array
656            il.Emit(System.Reflection.Emit.OpCodes.Ldc_I4, (int)currentInstr.data);
657            // load correct column of the current variable
658            il.Emit(System.Reflection.Emit.OpCodes.Ldelem_Ref);
659            il.Emit(System.Reflection.Emit.OpCodes.Ldarg_0); // rowIndex
660            if (!state.InLaggedContext) {
661              il.Emit(System.Reflection.Emit.OpCodes.Call, listGetValue);
662              il.Emit(System.Reflection.Emit.OpCodes.Ldc_R8, varNode.Weight); // load weight
663              il.Emit(System.Reflection.Emit.OpCodes.Mul);
664            } else {
665              var nanResult = il.DefineLabel();
666              var normalResult = il.DefineLabel();
667              il.Emit(System.Reflection.Emit.OpCodes.Dup);
668              il.Emit(System.Reflection.Emit.OpCodes.Ldc_I4_0);
669              il.Emit(System.Reflection.Emit.OpCodes.Blt, nanResult);
670              il.Emit(System.Reflection.Emit.OpCodes.Dup);
671              il.Emit(System.Reflection.Emit.OpCodes.Ldc_I4, ds.Rows);
672              il.Emit(System.Reflection.Emit.OpCodes.Bge, nanResult);
673              il.Emit(System.Reflection.Emit.OpCodes.Call, listGetValue);
674              il.Emit(System.Reflection.Emit.OpCodes.Ldc_R8, varNode.Weight); // load weight
675              il.Emit(System.Reflection.Emit.OpCodes.Mul);
676              il.Emit(System.Reflection.Emit.OpCodes.Br, normalResult);
677              il.MarkLabel(nanResult);
678              il.Emit(System.Reflection.Emit.OpCodes.Pop); // rowIndex
679              il.Emit(System.Reflection.Emit.OpCodes.Pop); // column reference
680              il.Emit(System.Reflection.Emit.OpCodes.Ldc_R8, double.NaN);
681              il.MarkLabel(normalResult);
682            }
683            return;
684          }
685        case OpCodes.LagVariable: {
686            var nanResult = il.DefineLabel();
687            var normalResult = il.DefineLabel();
688            LaggedVariableTreeNode varNode = (LaggedVariableTreeNode)currentInstr.dynamicNode;
689            il.Emit(System.Reflection.Emit.OpCodes.Ldarg_1); // load columns array
690            il.Emit(System.Reflection.Emit.OpCodes.Ldc_I4, (int)currentInstr.data);
691            // load correct column of the current variable
692            il.Emit(System.Reflection.Emit.OpCodes.Ldelem_Ref);
693            il.Emit(System.Reflection.Emit.OpCodes.Ldc_I4, varNode.Lag); // lag
694            il.Emit(System.Reflection.Emit.OpCodes.Ldarg_0); // rowIndex
695            il.Emit(System.Reflection.Emit.OpCodes.Add); // actualRowIndex = rowIndex + sampleOffset
696            il.Emit(System.Reflection.Emit.OpCodes.Dup);
697            il.Emit(System.Reflection.Emit.OpCodes.Ldc_I4_0);
698            il.Emit(System.Reflection.Emit.OpCodes.Blt, nanResult);
699            il.Emit(System.Reflection.Emit.OpCodes.Dup);
700            il.Emit(System.Reflection.Emit.OpCodes.Ldc_I4, ds.Rows);
701            il.Emit(System.Reflection.Emit.OpCodes.Bge, nanResult);
702            il.Emit(System.Reflection.Emit.OpCodes.Call, listGetValue);
703            il.Emit(System.Reflection.Emit.OpCodes.Ldc_R8, varNode.Weight); // load weight
704            il.Emit(System.Reflection.Emit.OpCodes.Mul);
705            il.Emit(System.Reflection.Emit.OpCodes.Br, normalResult);
706            il.MarkLabel(nanResult);
707            il.Emit(System.Reflection.Emit.OpCodes.Pop); // sample index
708            il.Emit(System.Reflection.Emit.OpCodes.Pop); // column reference
709            il.Emit(System.Reflection.Emit.OpCodes.Ldc_R8, double.NaN);
710            il.MarkLabel(normalResult);
711            return;
712          }
713        case OpCodes.Constant: {
714            ConstantTreeNode constNode = (ConstantTreeNode)currentInstr.dynamicNode;
715            il.Emit(System.Reflection.Emit.OpCodes.Ldc_R8, constNode.Value);
716            return;
717          }
718
719        //mkommend: this symbol uses the logistic function f(x) = 1 / (1 + e^(-alpha * x) )
720        //to determine the relative amounts of the true and false branch see http://en.wikipedia.org/wiki/Logistic_function
721        case OpCodes.VariableCondition: {
722            throw new NotSupportedException("Interpretation of symbol " + currentInstr.dynamicNode.Symbol.Name +
723                                            " is not supported by the SymbolicDataAnalysisTreeILEmittingInterpreter");
724          }
725        default:
726          throw new NotSupportedException("Interpretation of symbol " + currentInstr.dynamicNode.Symbol.Name +
727                                          " is not supported by the SymbolicDataAnalysisTreeILEmittingInterpreter");
728      }
729    }
730
731    public static double AiryA(double x) {
732      if (double.IsNaN(x)) return double.NaN;
733      double ai, aip, bi, bip;
734      alglib.airy(x, out ai, out aip, out bi, out bip);
735      return ai;
736    }
737
738    public static double AiryB(double x) {
739      if (double.IsNaN(x)) return double.NaN;
740      double ai, aip, bi, bip;
741      alglib.airy(x, out ai, out aip, out bi, out bip);
742      return bi;
743    }
744    public static double Dawson(double x) {
745      if (double.IsNaN(x)) return double.NaN;
746      return alglib.dawsonintegral(x);
747    }
748
749    public static double Gamma(double x) {
750      if (double.IsNaN(x)) return double.NaN;
751      return alglib.gammafunction(x);
752    }
753
754    public static double Psi(double x) {
755      if (double.IsNaN(x)) return double.NaN;
756      else if (x <= 0 && (Math.Floor(x) - x).IsAlmost(0)) return double.NaN;
757      return alglib.psi(x);
758    }
759
760    public static double ExpIntegralEi(double x) {
761      if (double.IsNaN(x)) return double.NaN;
762      return alglib.exponentialintegralei(x);
763    }
764
765    public static double SinIntegral(double x) {
766      if (double.IsNaN(x)) return double.NaN;
767      double si, ci;
768      alglib.sinecosineintegrals(x, out si, out ci);
769      return si;
770    }
771
772    public static double CosIntegral(double x) {
773      if (double.IsNaN(x)) return double.NaN;
774      double si, ci;
775      alglib.sinecosineintegrals(x, out si, out ci);
776      return ci;
777    }
778
779    public static double HypSinIntegral(double x) {
780      if (double.IsNaN(x)) return double.NaN;
781      double shi, chi;
782      alglib.hyperbolicsinecosineintegrals(x, out shi, out chi);
783      return shi;
784    }
785
786    public static double HypCosIntegral(double x) {
787      if (double.IsNaN(x)) return double.NaN;
788      double shi, chi;
789      alglib.hyperbolicsinecosineintegrals(x, out shi, out chi);
790      return chi;
791    }
792
793    public static double FresnelCosIntegral(double x) {
794      if (double.IsNaN(x)) return double.NaN;
795      double c = 0, s = 0;
796      alglib.fresnelintegral(x, ref c, ref s);
797      return c;
798    }
799
800    public static double FresnelSinIntegral(double x) {
801      if (double.IsNaN(x)) return double.NaN;
802      double c = 0, s = 0;
803      alglib.fresnelintegral(x, ref c, ref s);
804      return s;
805    }
806
807    public static double Norm(double x) {
808      if (double.IsNaN(x)) return double.NaN;
809      return alglib.normaldistribution(x);
810    }
811
812    public static double Erf(double x) {
813      if (double.IsNaN(x)) return double.NaN;
814      return alglib.errorfunction(x);
815    }
816
817    public static double Bessel(double x) {
818      if (double.IsNaN(x)) return double.NaN;
819      return alglib.besseli0(x);
820    }
821  }
822}
Note: See TracBrowser for help on using the repository browser.