source: branches/HeuristicLab.LinqExpressionTreeInterpreter/HeuristicLab.Problems.DataAnalysis.Symbolic/3.4/Interpreter/SymbolicDataAnalysisExpressionTreeILEmittingInterpreter.cs @ 13141

Last change on this file since 13141 was 13141, checked in by bburlacu, 6 years ago

#2442: Merged files from trunk and updated project file. Implemented missing operations in the CompiledTreeInterpreter: Integral, Derivative, Lag, TimeLag. Adapted lambda signature to accept an array of List<double> in order to make it easier to work with compiled trees. Changed value parameters to fixed value parameters and adjusted interpreter constructors and after serialization hooks. Removed function symbol.

From the performance point of view, compiling the tree into a lambda accepting a double[][] parameter (an array of arrays for the values of each double variable), accessed with Expression.ArrayIndex is the fastest, but it can be cumbersome to provide the data as a double[][]. Therefore the variant with List<double>[] was chosen. Internally, for each variable node the List's underlying double array is used, result in an overall decent speed compromise.

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