Free cookie consent management tool by TermsFeed Policy Generator

source: branches/2994-AutoDiffForIntervals/HeuristicLab.Problems.DataAnalysis.Symbolic/3.4/Interpreter/SymbolicDataAnalysisExpressionTreeILEmittingInterpreter.cs @ 17217

Last change on this file since 17217 was 17209, checked in by gkronber, 5 years ago

#2994: merged r17132:17198 from trunk to branch

File size: 39.1 KB
RevLine 
[6732]1#region License Information
2/* HeuristicLab
[17209]3 * Copyright (C) Heuristic and Evolutionary Algorithms Laboratory (HEAL)
[6732]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;
[7615]23using System.Collections.Generic;
[6741]24using System.Linq;
[6732]25using System.Reflection;
26using System.Reflection.Emit;
27using HeuristicLab.Common;
28using HeuristicLab.Core;
[6740]29using HeuristicLab.Data;
[6732]30using HeuristicLab.Encodings.SymbolicExpressionTreeEncoding;
[6740]31using HeuristicLab.Parameters;
[16565]32using HEAL.Attic;
[6732]33
34namespace HeuristicLab.Problems.DataAnalysis.Symbolic {
[16565]35  [StorableType("426718E3-2A57-4CA4-98A1-65EDD0B0BDBF")]
[6732]36  [Item("SymbolicDataAnalysisExpressionTreeILEmittingInterpreter", "Interpreter for symbolic expression trees.")]
[8436]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);
[7842]40
[8436]41    #region method infos
42    private static MethodInfo listGetValue = typeof(IList<double>).GetProperty("Item", new Type[] { typeof(int) }).GetGetMethod();
43
[6732]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) });
[16668]47    private static MethodInfo tanh = typeof(Math).GetMethod("Tanh", new Type[] { typeof(double) });
[6732]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) });
[6755]51    private static MethodInfo round = typeof(Math).GetMethod("Round", new Type[] { typeof(double) });
[7842]52    private static MethodInfo sqrt = typeof(Math).GetMethod("Sqrt", new Type[] { typeof(double) });
[16356]53    private static MethodInfo abs = typeof(Math).GetMethod("Abs", new Type[] { typeof(double) });
[6732]54
[7842]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) });
[14826]70    private static MethodInfo string_eq = typeof(string).GetMethod("Equals", new Type[] { typeof(string) });
[8436]71    #endregion
[7842]72
[6732]73    private const string CheckExpressionsWithIntervalArithmeticParameterName = "CheckExpressionsWithIntervalArithmetic";
[13248]74    private const string CheckExpressionsWithIntervalArithmeticParameterDescription = "Switch that determines if the interpreter checks the validity of expressions with interval arithmetic before evaluating the expression.";
[7615]75    private const string EvaluatedSolutionsParameterName = "EvaluatedSolutions";
[7842]76
[6732]77    public override bool CanChangeName {
78      get { return false; }
79    }
[7842]80
[6732]81    public override bool CanChangeDescription {
82      get { return false; }
83    }
84
85    #region parameter properties
[13248]86    public IFixedValueParameter<BoolValue> CheckExpressionsWithIntervalArithmeticParameter {
87      get { return (IFixedValueParameter<BoolValue>)Parameters[CheckExpressionsWithIntervalArithmeticParameterName]; }
[6732]88    }
[7615]89
[13248]90    public IFixedValueParameter<IntValue> EvaluatedSolutionsParameter {
91      get { return (IFixedValueParameter<IntValue>)Parameters[EvaluatedSolutionsParameterName]; }
[7615]92    }
[6732]93    #endregion
94
95    #region properties
[13248]96    public bool CheckExpressionsWithIntervalArithmetic {
97      get { return CheckExpressionsWithIntervalArithmeticParameter.Value.Value; }
98      set { CheckExpressionsWithIntervalArithmeticParameter.Value.Value = value; }
[6732]99    }
[13248]100    public int EvaluatedSolutions {
101      get { return EvaluatedSolutionsParameter.Value.Value; }
102      set { EvaluatedSolutionsParameter.Value.Value = value; }
[7615]103    }
[6732]104    #endregion
105
106    [StorableConstructor]
[16565]107    private SymbolicDataAnalysisExpressionTreeILEmittingInterpreter(StorableConstructorFlag _) : base(_) { }
[7842]108
[8436]109    private SymbolicDataAnalysisExpressionTreeILEmittingInterpreter(SymbolicDataAnalysisExpressionTreeILEmittingInterpreter original, Cloner cloner) : base(original, cloner) { }
[6732]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.") {
[13248]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)));
[6732]119    }
120
[13248]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
[7615]128    [StorableHook(HookType.AfterDeserialization)]
129    private void AfterDeserialization() {
[13248]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));
[7615]144    }
145
146    #region IStatefulItem
147    public void InitializeState() {
[13248]148      EvaluatedSolutions = 0;
[7615]149    }
150
151    public void ClearState() {
152    }
153    #endregion
154
[13251]155    private readonly object syncRoot = new object();
[12509]156    public IEnumerable<double> GetSymbolicExpressionTreeValues(ISymbolicExpressionTree tree, IDataset dataset, IEnumerable<int> rows) {
[13248]157      if (CheckExpressionsWithIntervalArithmetic)
[8436]158        throw new NotSupportedException("Interval arithmetic is not yet supported in the symbolic data analysis interpreter.");
159
[13251]160      lock (syncRoot) {
161        EvaluatedSolutions++; // increment the evaluated solutions counter
162      }
[8436]163      var state = PrepareInterpreterState(tree, dataset);
[6741]164
[8436]165      Type[] methodArgs = { typeof(int), typeof(IList<double>[]) };
166      DynamicMethod testFun = new DynamicMethod("TestFun", typeof(double), methodArgs, typeof(SymbolicDataAnalysisExpressionTreeILEmittingInterpreter).Module);
[6741]167
[8436]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
[12509]181    private InterpreterState PrepareInterpreterState(ISymbolicExpressionTree tree, IDataset dataset) {
[8436]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) {
[6732]186        if (instr.opCode == OpCodes.Variable) {
[8436]187          var variableTreeNode = (VariableTreeNode)instr.dynamicNode;
[9828]188          instr.data = doubleVariableNames[variableTreeNode.VariableName];
[6732]189        } else if (instr.opCode == OpCodes.LagVariable) {
[8436]190          var laggedVariableTreeNode = (LaggedVariableTreeNode)instr.dynamicNode;
[9828]191          instr.data = doubleVariableNames[laggedVariableTreeNode.VariableName];
[6732]192        } else if (instr.opCode == OpCodes.VariableCondition) {
[8436]193          var variableConditionTreeNode = (VariableConditionTreeNode)instr.dynamicNode;
[9828]194          instr.data = doubleVariableNames[variableConditionTreeNode.VariableName];
[6732]195        } else if (instr.opCode == OpCodes.Call) {
196          necessaryArgStackSize += instr.nArguments + 1;
197        }
198      }
[8436]199      return new InterpreterState(code, necessaryArgStackSize);
[6732]200    }
201
[12509]202    private void CompileInstructions(ILGenerator il, InterpreterState state, IDataset ds) {
[6732]203      Instruction currentInstr = state.NextInstruction();
204      int nArgs = currentInstr.nArguments;
205
206      switch (currentInstr.opCode) {
207        case OpCodes.Add: {
208            if (nArgs > 0) {
[6809]209              CompileInstructions(il, state, ds);
[6732]210            }
211            for (int i = 1; i < nArgs; i++) {
[6809]212              CompileInstructions(il, state, ds);
[6732]213              il.Emit(System.Reflection.Emit.OpCodes.Add);
214            }
215            return;
216          }
217        case OpCodes.Sub: {
218            if (nArgs == 1) {
[6809]219              CompileInstructions(il, state, ds);
[6732]220              il.Emit(System.Reflection.Emit.OpCodes.Neg);
221              return;
222            }
223            if (nArgs > 0) {
[6809]224              CompileInstructions(il, state, ds);
[6732]225            }
226            for (int i = 1; i < nArgs; i++) {
[6809]227              CompileInstructions(il, state, ds);
[6732]228              il.Emit(System.Reflection.Emit.OpCodes.Sub);
229            }
230            return;
231          }
232        case OpCodes.Mul: {
233            if (nArgs > 0) {
[6809]234              CompileInstructions(il, state, ds);
[6732]235            }
236            for (int i = 1; i < nArgs; i++) {
[6809]237              CompileInstructions(il, state, ds);
[6732]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);
[6809]245              CompileInstructions(il, state, ds);
[6732]246              il.Emit(System.Reflection.Emit.OpCodes.Div);
247              return;
248            }
249            if (nArgs > 0) {
[6809]250              CompileInstructions(il, state, ds);
[6732]251            }
252            for (int i = 1; i < nArgs; i++) {
[6809]253              CompileInstructions(il, state, ds);
[6732]254              il.Emit(System.Reflection.Emit.OpCodes.Div);
255            }
256            return;
257          }
258        case OpCodes.Average: {
[6809]259            CompileInstructions(il, state, ds);
[6732]260            for (int i = 1; i < nArgs; i++) {
[6809]261              CompileInstructions(il, state, ds);
[6732]262              il.Emit(System.Reflection.Emit.OpCodes.Add);
263            }
[14809]264            il.Emit(System.Reflection.Emit.OpCodes.Ldc_R8, (double)nArgs);
[6732]265            il.Emit(System.Reflection.Emit.OpCodes.Div);
266            return;
267          }
[16356]268        case OpCodes.Absolute: {
269            CompileInstructions(il, state, ds);
270            il.Emit(System.Reflection.Emit.OpCodes.Call, abs);
271            return;
272          }
[6732]273        case OpCodes.Cos: {
[6809]274            CompileInstructions(il, state, ds);
[6732]275            il.Emit(System.Reflection.Emit.OpCodes.Call, cos);
276            return;
277          }
278        case OpCodes.Sin: {
[6809]279            CompileInstructions(il, state, ds);
[6732]280            il.Emit(System.Reflection.Emit.OpCodes.Call, sin);
281            return;
282          }
283        case OpCodes.Tan: {
[6809]284            CompileInstructions(il, state, ds);
[6732]285            il.Emit(System.Reflection.Emit.OpCodes.Call, tan);
286            return;
287          }
[16670]288        case OpCodes.Tanh: {
[16668]289            CompileInstructions(il, state, ds);
290            il.Emit(System.Reflection.Emit.OpCodes.Call, tanh);
291            return;
292          }
[6732]293        case OpCodes.Power: {
[6809]294            CompileInstructions(il, state, ds);
295            CompileInstructions(il, state, ds);
[6755]296            il.Emit(System.Reflection.Emit.OpCodes.Call, round);
[6732]297            il.Emit(System.Reflection.Emit.OpCodes.Call, power);
298            return;
299          }
300        case OpCodes.Root: {
[6809]301            CompileInstructions(il, state, ds);
[6755]302            il.Emit(System.Reflection.Emit.OpCodes.Ldc_R8, 1.0); // 1 / round(...)
[6809]303            CompileInstructions(il, state, ds);
[6755]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;
[6732]308          }
309        case OpCodes.Exp: {
[6809]310            CompileInstructions(il, state, ds);
[6732]311            il.Emit(System.Reflection.Emit.OpCodes.Call, exp);
312            return;
313          }
314        case OpCodes.Log: {
[6809]315            CompileInstructions(il, state, ds);
[6732]316            il.Emit(System.Reflection.Emit.OpCodes.Call, log);
317            return;
318          }
[7842]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          }
[16356]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          }
[7842]331        case OpCodes.SquareRoot: {
332            CompileInstructions(il, state, ds);
333            il.Emit(System.Reflection.Emit.OpCodes.Call, sqrt);
334            return;
335          }
[16356]336        case OpCodes.CubeRoot: {
337            CompileInstructions(il, state, ds);
[16911]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
[16356]346            il.Emit(System.Reflection.Emit.OpCodes.Ldc_R8, 1.0 / 3.0);
347            il.Emit(System.Reflection.Emit.OpCodes.Call, power);
[16911]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);
[16356]354            return;
355          }
[7842]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          }
[16360]431        case OpCodes.AnalyticQuotient: {
[16356]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);
[17120]438            il.Emit(System.Reflection.Emit.OpCodes.Add); // 1+x2*x2
[16356]439            il.Emit(System.Reflection.Emit.OpCodes.Call, sqrt);
440            il.Emit(System.Reflection.Emit.OpCodes.Div);
441            return;
442          }
[6732]443        case OpCodes.IfThenElse: {
444            Label end = il.DefineLabel();
445            Label c1 = il.DefineLabel();
[6809]446            CompileInstructions(il, state, ds);
[16585]447            il.Emit(System.Reflection.Emit.OpCodes.Ldc_R8, 0.0); // > 0
[6732]448            il.Emit(System.Reflection.Emit.OpCodes.Cgt);
449            il.Emit(System.Reflection.Emit.OpCodes.Brfalse, c1);
[6809]450            CompileInstructions(il, state, ds);
[6732]451            il.Emit(System.Reflection.Emit.OpCodes.Br, end);
452            il.MarkLabel(c1);
[6809]453            CompileInstructions(il, state, ds);
[6732]454            il.MarkLabel(end);
455            return;
456          }
457        case OpCodes.AND: {
458            Label falseBranch = il.DefineLabel();
459            Label end = il.DefineLabel();
[6809]460            CompileInstructions(il, state, ds);
[6732]461            for (int i = 1; i < nArgs; i++) {
[16585]462              il.Emit(System.Reflection.Emit.OpCodes.Ldc_R8, 0.0); // > 0
[6732]463              il.Emit(System.Reflection.Emit.OpCodes.Cgt);
464              il.Emit(System.Reflection.Emit.OpCodes.Brfalse, falseBranch);
[6809]465              CompileInstructions(il, state, ds);
[6732]466            }
[16585]467            il.Emit(System.Reflection.Emit.OpCodes.Ldc_R8, 0.0); // > 0
[6732]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();
[6809]482            CompileInstructions(il, state, ds);
[6732]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);
[16585]487              il.Emit(System.Reflection.Emit.OpCodes.Ldc_R8, 0.0); // <= 0       
[6732]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);
[6809]492              CompileInstructions(il, state, ds);
[6732]493            }
494            il.MarkLabel(resultBranch);
[16585]495            il.Emit(System.Reflection.Emit.OpCodes.Ldc_R8, 0.0); // > 0
[6732]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: {
[6809]507            CompileInstructions(il, state, ds);
[16585]508            il.Emit(System.Reflection.Emit.OpCodes.Ldc_R8, 0.0); // > 0
[6732]509            il.Emit(System.Reflection.Emit.OpCodes.Cgt);
[14809]510            il.Emit(System.Reflection.Emit.OpCodes.Conv_R8); // convert to float64
[6732]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          }
[10791]518        case OpCodes.XOR: {
519            CompileInstructions(il, state, ds);
[16585]520            il.Emit(System.Reflection.Emit.OpCodes.Ldc_R8, 0.0);
[10791]521            il.Emit(System.Reflection.Emit.OpCodes.Cgt);// > 0
522
523            for (int i = 1; i < nArgs; i++) {
524              CompileInstructions(il, state, ds);
[16585]525              il.Emit(System.Reflection.Emit.OpCodes.Ldc_R8, 0.0);
[10791]526              il.Emit(System.Reflection.Emit.OpCodes.Cgt);// > 0
527              il.Emit(System.Reflection.Emit.OpCodes.Xor);
528            }
[14809]529            il.Emit(System.Reflection.Emit.OpCodes.Conv_R8); // convert to float64
530
[10791]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          }
[6732]537        case OpCodes.GT: {
[6809]538            CompileInstructions(il, state, ds);
539            CompileInstructions(il, state, ds);
540
[6732]541            il.Emit(System.Reflection.Emit.OpCodes.Cgt); // 1 (>) / 0 (otherwise)
[14809]542            il.Emit(System.Reflection.Emit.OpCodes.Conv_R8); // convert to float64
[6732]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: {
[6809]550            CompileInstructions(il, state, ds);
551            CompileInstructions(il, state, ds);
[6732]552            il.Emit(System.Reflection.Emit.OpCodes.Clt);
[14809]553            il.Emit(System.Reflection.Emit.OpCodes.Conv_R8); // convert to float64
[6732]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: {
[6809]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);
[6849]566            var prevLaggedContext = state.InLaggedContext;
567            state.InLaggedContext = true;
[6809]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);
[6849]573            state.InLaggedContext = prevLaggedContext;
[6809]574            return;
[6732]575          }
576        case OpCodes.Integral: {
[6809]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);
[6849]583            var prevLaggedContext = state.InLaggedContext;
584            state.InLaggedContext = true;
[6809]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            }
[6849]595            state.InLaggedContext = prevLaggedContext;
[6809]596            return;
[6732]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: {
[6809]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;
[6849]611            var prevLaggedContext = state.InLaggedContext;
612            state.InLaggedContext = true;
[6809]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);
[6849]642            state.InLaggedContext = prevLaggedContext;
[6809]643            return;
[6732]644          }
645        case OpCodes.Call: {
[7842]646            throw new NotSupportedException(
647              "Automatically defined functions are not supported by the SymbolicDataAnalysisTreeILEmittingInterpreter. Either turn of ADFs or change the interpeter.");
[6732]648          }
649        case OpCodes.Arg: {
[7842]650            throw new NotSupportedException(
651              "Automatically defined functions are not supported by the SymbolicDataAnalysisTreeILEmittingInterpreter. Either turn of ADFs or change the interpeter.");
[6732]652          }
653        case OpCodes.Variable: {
[6741]654            VariableTreeNode varNode = (VariableTreeNode)currentInstr.dynamicNode;
[8798]655            il.Emit(System.Reflection.Emit.OpCodes.Ldarg_1); // load columns array
[9828]656            il.Emit(System.Reflection.Emit.OpCodes.Ldc_I4, (int)currentInstr.data);
[8798]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
[6849]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
[8798]679              il.Emit(System.Reflection.Emit.OpCodes.Pop); // column reference
[6849]680              il.Emit(System.Reflection.Emit.OpCodes.Ldc_R8, double.NaN);
681              il.MarkLabel(normalResult);
682            }
[6732]683            return;
684          }
685        case OpCodes.LagVariable: {
[6809]686            var nanResult = il.DefineLabel();
687            var normalResult = il.DefineLabel();
[6770]688            LaggedVariableTreeNode varNode = (LaggedVariableTreeNode)currentInstr.dynamicNode;
[8798]689            il.Emit(System.Reflection.Emit.OpCodes.Ldarg_1); // load columns array
[9828]690            il.Emit(System.Reflection.Emit.OpCodes.Ldc_I4, (int)currentInstr.data);
[8798]691            // load correct column of the current variable
692            il.Emit(System.Reflection.Emit.OpCodes.Ldelem_Ref);
[6770]693            il.Emit(System.Reflection.Emit.OpCodes.Ldc_I4, varNode.Lag); // lag
[6849]694            il.Emit(System.Reflection.Emit.OpCodes.Ldarg_0); // rowIndex
695            il.Emit(System.Reflection.Emit.OpCodes.Add); // actualRowIndex = rowIndex + sampleOffset
[6809]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);
[6770]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);
[6809]705            il.Emit(System.Reflection.Emit.OpCodes.Br, normalResult);
706            il.MarkLabel(nanResult);
[8798]707            il.Emit(System.Reflection.Emit.OpCodes.Pop); // sample index
708            il.Emit(System.Reflection.Emit.OpCodes.Pop); // column reference
[6809]709            il.Emit(System.Reflection.Emit.OpCodes.Ldc_R8, double.NaN);
710            il.MarkLabel(normalResult);
[6770]711            return;
[6732]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: {
[7842]722            throw new NotSupportedException("Interpretation of symbol " + currentInstr.dynamicNode.Symbol.Name +
723                                            " is not supported by the SymbolicDataAnalysisTreeILEmittingInterpreter");
[6732]724          }
[7842]725        default:
726          throw new NotSupportedException("Interpretation of symbol " + currentInstr.dynamicNode.Symbol.Name +
727                                          " is not supported by the SymbolicDataAnalysisTreeILEmittingInterpreter");
[6732]728      }
729    }
730
[7842]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;
[8430]756      else if (x <= 0 && (Math.Floor(x) - x).IsAlmost(0)) return double.NaN;
[7842]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);
[8430]783      return shi;
[7842]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    }
[6732]821  }
822}
Note: See TracBrowser for help on using the repository browser.