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

Last change on this file since 14811 was 14811, checked in by mkommend, 2 years ago

#2442: Merged r12807, r13039, r13139, r13140, r13141, r13222, r13247, r13248, r13251, r13254, r13255, r13256, r1326, r13288, r13313, r13314, r13315, r13318, r14282, r14809, r14810 into stable.

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