Index: branches/2925_AutoDiffForDynamicalModels/HeuristicLab.Problems.DynamicalSystemsModelling/3.3/OdeParameterIdentification.cs
===================================================================
--- branches/2925_AutoDiffForDynamicalModels/HeuristicLab.Problems.DynamicalSystemsModelling/3.3/OdeParameterIdentification.cs	(revision 16601)
+++ branches/2925_AutoDiffForDynamicalModels/HeuristicLab.Problems.DynamicalSystemsModelling/3.3/OdeParameterIdentification.cs	(revision 16602)
@@ -177,5 +177,5 @@
     public void CreateSolution(Problem problem, string[] modelStructure, int maxIterations, IRandom rand) {
       var parser = new InfixExpressionParser();
-      var trees = modelStructure.Select(expr => Convert(parser.Parse(expr))).ToArray();
+      var trees = modelStructure.Select(expr => parser.Parse(expr)).ToArray();
       var names = problem.Encoding.Encodings.Select(enc => enc.Name).ToArray();
       if (trees.Length != names.Length) throw new ArgumentException("The number of expressions must match the number of target variables exactly");
@@ -190,68 +190,68 @@
     }
 
-    private ISymbolicExpressionTree Convert(ISymbolicExpressionTree tree) {
-      return new SymbolicExpressionTree(Convert(tree.Root));
-    }
+    // private ISymbolicExpressionTree Convert(ISymbolicExpressionTree tree) {
+    //   return new SymbolicExpressionTree(Convert(tree.Root));
+    // }
 
 
     // for translation from symbolic expressions to simple symbols
-    private static Dictionary<Type, string> sym2str = new Dictionary<Type, string>() {
-      {typeof(Addition), "+" },
-      {typeof(Subtraction), "-" },
-      {typeof(Multiplication), "*" },
-      {typeof(Sine), "sin" },
-      {typeof(Cosine), "cos" },
-      {typeof(Square), "sqr" },
-    };
-
-    private ISymbolicExpressionTreeNode Convert(ISymbolicExpressionTreeNode node) {
-      if (sym2str.ContainsKey(node.Symbol.GetType())) {
-        var children = node.Subtrees.Select(st => Convert(st)).ToArray();
-        return Make(sym2str[node.Symbol.GetType()], children);
-      } else if (node.Symbol is ProgramRootSymbol) {
-        var child = Convert(node.GetSubtree(0));
-        node.RemoveSubtree(0);
-        node.AddSubtree(child);
-        return node;
-      } else if (node.Symbol is StartSymbol) {
-        var child = Convert(node.GetSubtree(0));
-        node.RemoveSubtree(0);
-        node.AddSubtree(child);
-        return node;
-      } else if (node.Symbol is Division) {
-        var children = node.Subtrees.Select(st => Convert(st)).ToArray();
-        if (children.Length == 1) {
-          return Make("%", new[] { new SimpleSymbol("θ", 0).CreateTreeNode(), children[0] });
-        } else if (children.Length != 2) throw new ArgumentException("Division is not supported for multiple arguments");
-        else return Make("%", children);
-      } else if (node.Symbol is Constant) {
-        return new SimpleSymbol("θ", 0).CreateTreeNode();
-      } else if (node.Symbol is DataAnalysis.Symbolic.Variable) {
-        var varNode = node as VariableTreeNode;
-        if (!varNode.Weight.IsAlmost(1.0)) throw new ArgumentException("Variable weights are not supported");
-        return new SimpleSymbol(varNode.VariableName, 0).CreateTreeNode();
-      } else throw new ArgumentException("Unsupported symbol: " + node.Symbol.Name);
-    }
-
-    private ISymbolicExpressionTreeNode Make(string op, ISymbolicExpressionTreeNode[] children) {
-      if (children.Length == 1) {
-        var s = new SimpleSymbol(op, 1).CreateTreeNode();
-        s.AddSubtree(children.First());
-        return s;
-      } else {
-        var s = new SimpleSymbol(op, 2).CreateTreeNode();
-        var c0 = children[0];
-        var c1 = children[1];
-        s.AddSubtree(c0);
-        s.AddSubtree(c1);
-        for (int i = 2; i < children.Length; i++) {
-          var sn = new SimpleSymbol(op, 2).CreateTreeNode();
-          sn.AddSubtree(s);
-          sn.AddSubtree(children[i]);
-          s = sn;
-        }
-        return s;
-      }
-    }
+    // private static Dictionary<Type, string> sym2str = new Dictionary<Type, string>() {
+    //   {typeof(Addition), "+" },
+    //   {typeof(Subtraction), "-" },
+    //   {typeof(Multiplication), "*" },
+    //   {typeof(Sine), "sin" },
+    //   {typeof(Cosine), "cos" },
+    //   {typeof(Square), "sqr" },
+    // };
+
+    // private ISymbolicExpressionTreeNode Convert(ISymbolicExpressionTreeNode node) {
+    //   if (sym2str.ContainsKey(node.Symbol.GetType())) {
+    //     var children = node.Subtrees.Select(st => Convert(st)).ToArray();
+    //     return Make(sym2str[node.Symbol.GetType()], children);
+    //   } else if (node.Symbol is ProgramRootSymbol) {
+    //     var child = Convert(node.GetSubtree(0));
+    //     node.RemoveSubtree(0);
+    //     node.AddSubtree(child);
+    //     return node;
+    //   } else if (node.Symbol is StartSymbol) {
+    //     var child = Convert(node.GetSubtree(0));
+    //     node.RemoveSubtree(0);
+    //     node.AddSubtree(child);
+    //     return node;
+    //   } else if (node.Symbol is Division) {
+    //     var children = node.Subtrees.Select(st => Convert(st)).ToArray();
+    //     if (children.Length == 1) {
+    //       return Make("%", new[] { new SimpleSymbol("θ", 0).CreateTreeNode(), children[0] });
+    //     } else if (children.Length != 2) throw new ArgumentException("Division is not supported for multiple arguments");
+    //     else return Make("%", children);
+    //   } else if (node.Symbol is Constant) {
+    //     return new SimpleSymbol("θ", 0).CreateTreeNode();
+    //   } else if (node.Symbol is DataAnalysis.Symbolic.Variable) {
+    //     var varNode = node as VariableTreeNode;
+    //     if (!varNode.Weight.IsAlmost(1.0)) throw new ArgumentException("Variable weights are not supported");
+    //     return new SimpleSymbol(varNode.VariableName, 0).CreateTreeNode();
+    //   } else throw new ArgumentException("Unsupported symbol: " + node.Symbol.Name);
+    // }
+
+    // private ISymbolicExpressionTreeNode Make(string op, ISymbolicExpressionTreeNode[] children) {
+    //   if (children.Length == 1) {
+    //     var s = new SimpleSymbol(op, 1).CreateTreeNode();
+    //     s.AddSubtree(children.First());
+    //     return s;
+    //   } else {
+    //     var s = new SimpleSymbol(op, 2).CreateTreeNode();
+    //     var c0 = children[0];
+    //     var c1 = children[1];
+    //     s.AddSubtree(c0);
+    //     s.AddSubtree(c1);
+    //     for (int i = 2; i < children.Length; i++) {
+    //       var sn = new SimpleSymbol(op, 2).CreateTreeNode();
+    //       sn.AddSubtree(s);
+    //       sn.AddSubtree(children[i]);
+    //       s = sn;
+    //     }
+    //     return s;
+    //   }
+    // }
     #endregion
   }
Index: branches/2925_AutoDiffForDynamicalModels/HeuristicLab.Problems.DynamicalSystemsModelling/3.3/Problem.cs
===================================================================
--- branches/2925_AutoDiffForDynamicalModels/HeuristicLab.Problems.DynamicalSystemsModelling/3.3/Problem.cs	(revision 16601)
+++ branches/2925_AutoDiffForDynamicalModels/HeuristicLab.Problems.DynamicalSystemsModelling/3.3/Problem.cs	(revision 16602)
@@ -71,4 +71,5 @@
       get { return (IFixedValueParameter<IntValue>)Parameters[MaximumLengthParameterName]; }
     }
+
     public IFixedValueParameter<IntValue> MaximumParameterOptimizationIterationsParameter {
       get { return (IFixedValueParameter<IntValue>)Parameters[MaximumParameterOptimizationIterationsParameterName]; }
@@ -196,6 +197,4 @@
     public override double Evaluate(Individual individual, IRandom random) {
       var trees = individual.Values.Select(v => v.Value).OfType<ISymbolicExpressionTree>().ToArray(); // extract all trees from individual
-                                                                                                      // write back optimized parameters to tree nodes instead of the separate OptTheta variable
-                                                                                                      // retreive optimized parameters from nodes?
 
       var problemData = ProblemData;
@@ -208,8 +207,7 @@
         int totalSize = 0;
         foreach (var episode in TrainingEpisodes) {
-          double[] optTheta;
-          double nmse;
-          OptimizeForEpisodes(trees, problemData, targetVars, latentVariables, random, new[] { episode }, MaximumParameterOptimizationIterations, NumericIntegrationSteps, OdeSolver, out optTheta, out nmse);
-          individual["OptTheta_" + eIdx] = new DoubleArray(optTheta); // write back optimized parameters so that we can use them in the Analysis method
+          // double[] optTheta;
+          double nmse = OptimizeForEpisodes(trees, problemData, targetVars, latentVariables, random, new[] { episode }, MaximumParameterOptimizationIterations, NumericIntegrationSteps, OdeSolver);
+          // individual["OptTheta_" + eIdx] = new DoubleArray(optTheta); // write back optimized parameters so that we can use them in the Analysis method
           eIdx++;
           totalNMSE += nmse * episode.Size;
@@ -218,13 +216,12 @@
         return totalNMSE / totalSize;
       } else {
-        double[] optTheta;
-        double nmse;
-        OptimizeForEpisodes(trees, problemData, targetVars, latentVariables, random, TrainingEpisodes, MaximumParameterOptimizationIterations, NumericIntegrationSteps, OdeSolver, out optTheta, out nmse);
-        individual["OptTheta"] = new DoubleArray(optTheta); // write back optimized parameters so that we can use them in the Analysis method
+        // double[] optTheta;
+        double nmse = OptimizeForEpisodes(trees, problemData, targetVars, latentVariables, random, TrainingEpisodes, MaximumParameterOptimizationIterations, NumericIntegrationSteps, OdeSolver);
+        // individual["OptTheta"] = new DoubleArray(optTheta); // write back optimized parameters so that we can use them in the Analysis method
         return nmse;
       }
     }
 
-    public static void OptimizeForEpisodes(
+    public static double OptimizeForEpisodes(
       ISymbolicExpressionTree[] trees,
       IRegressionProblemData problemData,
@@ -235,16 +232,29 @@
       int maxParameterOptIterations,
       int numericIntegrationSteps,
-      string odeSolver,
-      out double[] optTheta,
-      out double nmse) {
-
+      string odeSolver) {
+
+      var constantNodes = trees.Select(t => t.IterateNodesPrefix().OfType<ConstantTreeNode>().ToArray()).ToArray();
+      var initialTheta = constantNodes.Select(nodes => nodes.Select(n => n.Value).ToArray()).ToArray();
 
       // optimize parameters by fitting f(x,y) to calculated differences dy/dt(t)
-      nmse = PreTuneParameters(trees, problemData, targetVars, latentVariables, random, episodes, maxParameterOptIterations, out optTheta);
+      double nmse = PreTuneParameters(trees, problemData, targetVars, latentVariables, random, episodes, maxParameterOptIterations,
+        initialTheta, out double[] pretunedParameters);
 
       // optimize parameters using integration of f(x,y) to calculate y(t)
-      nmse = OptimizeParameters(trees, problemData, targetVars, latentVariables, episodes, maxParameterOptIterations, optTheta, numericIntegrationSteps, odeSolver, out optTheta);
-
-      if (double.IsNaN(nmse) || double.IsInfinity(nmse)) nmse = 100 * trees.Length * episodes.Sum(ep => ep.Size);
+      nmse = OptimizeParameters(trees, problemData, targetVars, latentVariables, episodes, maxParameterOptIterations, pretunedParameters, numericIntegrationSteps, odeSolver,
+        out double[] optTheta);
+
+      if (double.IsNaN(nmse) ||
+        double.IsInfinity(nmse) ||
+        nmse > 100 * trees.Length * episodes.Sum(ep => ep.Size))
+        return 100 * trees.Length * episodes.Sum(ep => ep.Size);
+
+      // update tree nodes with optimized values
+      var paramIdx = 0;
+      for (var treeIdx = 0; treeIdx < constantNodes.Length; treeIdx++) {
+        for (int i = 0; i < constantNodes[treeIdx].Length; i++)
+          constantNodes[treeIdx][i].Value = optTheta[paramIdx++];
+      }
+      return nmse;
     }
 
@@ -257,7 +267,10 @@
       IEnumerable<IntRange> episodes,
       int maxParameterOptIterations,
+      double[][] initialTheta,
       out double[] optTheta) {
       var thetas = new List<double>();
       double nmse = 0.0;
+      var maxTreeNmse = 100 * episodes.Sum(ep => ep.Size);
+
       // NOTE: the order of values in parameter matches prefix order of constant nodes in trees
       for (int treeIdx = 0; treeIdx < trees.Length; treeIdx++) {
@@ -272,33 +285,36 @@
         var paramCount = myState.nodeValueLookup.ParameterCount;
 
-        // init params randomly
-        // theta contains parameter values for trees and then the initial values for latent variables (a separate vector for each episode)
-        // inital values for latent variables are also optimized
-        var theta = new double[paramCount + latentVariables.Length * episodes.Count()];
-        for (int i = 0; i < theta.Length; i++)
-          theta[i] = random.NextDouble() * 2.0e-1 - 1.0e-1;
-
         optTheta = new double[0];
-        if (theta.Length > 0) {
-          alglib.minlmstate state;
-          alglib.minlmreport report;
-          alglib.minlmcreatevj(targetValuesDiff.Length, theta, out state);
-          alglib.minlmsetcond(state, 0.0, 0.0, 0.0, maxParameterOptIterations);
-          // alglib.minlmsetgradientcheck(state, 1.0e-3);
-          alglib.minlmoptimize(state, EvaluateObjectiveVector, EvaluateObjectiveVectorAndJacobian, null, myState);
-
-          alglib.minlmresults(state, out optTheta, out report);
-
-
-          if (report.terminationtype < 0) { throw new InvalidOperationException("there was a problem in the optimizer"); }
-
+        if (initialTheta[treeIdx].Length > 0) {
+          try {
+            alglib.minlmstate state;
+            alglib.minlmreport report;
+            var p = new double[initialTheta[treeIdx].Length];
+            var lowerBounds = Enumerable.Repeat(-100.0, p.Length).ToArray();
+            var upperBounds = Enumerable.Repeat(100.0, p.Length).ToArray();
+            Array.Copy(initialTheta[treeIdx], p, p.Length);
+            alglib.minlmcreatevj(targetValuesDiff.Length, p, out state);
+            alglib.minlmsetcond(state, 0.0, 0.0, 0.0, maxParameterOptIterations);
+            alglib.minlmsetbc(state, lowerBounds, upperBounds);
+            // alglib.minlmsetgradientcheck(state, 1.0e-3);
+            alglib.minlmoptimize(state, EvaluateObjectiveVector, EvaluateObjectiveVectorAndJacobian, null, myState);
+
+            alglib.minlmresults(state, out optTheta, out report);
+            if (report.terminationtype < 0) { optTheta = initialTheta[treeIdx]; }
+          } catch (alglib.alglibexception) {
+            optTheta = initialTheta[treeIdx];
+          }
+        }
+        var tree_nmse = EvaluateMSE(optTheta, myState);
+        if (double.IsNaN(tree_nmse) || double.IsInfinity(tree_nmse) || tree_nmse > maxTreeNmse) {
+          nmse += maxTreeNmse;
+          thetas.AddRange(initialTheta[treeIdx]);
+        } else {
+          nmse += tree_nmse;
           thetas.AddRange(optTheta);
         }
-        nmse += EvaluateMSE(optTheta, myState);
       } // foreach tree
       optTheta = thetas.ToArray();
 
-      var maxNmse = 100 * trees.Length * episodes.Sum(ep => ep.Size);
-      if (double.IsNaN(nmse) || double.IsInfinity(nmse) || nmse > maxNmse) nmse = maxNmse;
       return nmse;
     }
@@ -321,9 +337,11 @@
 
       if (initialTheta.Length > 0) {
-
+        var lowerBounds = Enumerable.Repeat(-100.0, initialTheta.Length).ToArray();
+        var upperBounds = Enumerable.Repeat(100.0, initialTheta.Length).ToArray();
         try {
           alglib.minlmstate state;
           alglib.minlmreport report;
           alglib.minlmcreatevj(rowsForDataExtraction.Length * trees.Length, initialTheta, out state);
+          alglib.minlmsetbc(state, lowerBounds, upperBounds);
           alglib.minlmsetcond(state, 0.0, 0.0, 0.0, maxParameterOptIterations);
           // alglib.minlmsetgradientcheck(state, 1.0e-3);
@@ -490,4 +508,5 @@
 
       if (OptimizeParametersForEpisodes) {
+        throw new NotSupportedException();
         var eIdx = 0;
         var trainingPredictions = new List<Tuple<double, Vector>[][]>();
@@ -527,5 +546,5 @@
         results["Models"].Value = models;
       } else {
-        var optTheta = ((DoubleArray)bestIndividualAndQuality.Item1["OptTheta"]).ToArray(); // see evaluate
+        var optTheta = Problem.ExtractParametersFromTrees(trees);
         var optimizationData = new OptimizationData(trees, targetVars, problemData, null, TrainingEpisodes.ToArray(), NumericIntegrationSteps, latentVariables, OdeSolver);
         var trainingPrediction = Integrate(optimizationData, optTheta).ToArray();
@@ -629,5 +648,6 @@
         for (int idx = 0; idx < trees.Length; idx++) {
           var tree = trees[idx];
-          optimizedTrees.Add(new SymbolicExpressionTree(FixParameters(tree.Root, optTheta.ToArray(), ref nextParIdx)));
+          // optimizedTrees.Add(new SymbolicExpressionTree(FixParameters(tree.Root, optTheta.ToArray(), ref nextParIdx)));
+          optimizedTrees.Add(tree);
         }
         var ds = problemData.Dataset;
@@ -660,14 +680,9 @@
           var tree = trees[idx];
 
-          // when we reference HeuristicLab.Problems.DataAnalysis.Symbolic we can translate symbols
-          var shownTree = new SymbolicExpressionTree(TranslateTreeNode(tree.Root, optTheta.ToArray(),
-            ref nextParIdx));
-
-
           var origTreeVar = new HeuristicLab.Core.Variable(varName + "(original)");
           origTreeVar.Value = (ISymbolicExpressionTree)tree.Clone();
           models.Add(origTreeVar);
           var simplifiedTreeVar = new HeuristicLab.Core.Variable(varName + "(simplified)");
-          simplifiedTreeVar.Value = TreeSimplifier.Simplify(shownTree);
+          simplifiedTreeVar.Value = TreeSimplifier.Simplify(tree);
           models.Add(simplifiedTreeVar);
 
@@ -678,4 +693,12 @@
       }
     }
+
+
+    public static double[] ExtractParametersFromTrees(ISymbolicExpressionTree[] trees) {
+      return trees
+        .SelectMany(t => t.IterateNodesPrefix().OfType<ConstantTreeNode>().Select(n => n.Value))
+        .ToArray();
+    }
+
 
 
@@ -1104,55 +1127,46 @@
     }
 
+    // TODO: use an existing interpreter implementation instead
     private static double InterpretRec(ISymbolicExpressionTreeNode node, NodeValueLookup nodeValues) {
-      switch (node.Symbol.Name) {
-        case "+": {
-            var f = InterpretRec(node.GetSubtree(0), nodeValues);
-            var g = InterpretRec(node.GetSubtree(1), nodeValues);
-            return f + g;
-          }
-        case "*": {
-            var f = InterpretRec(node.GetSubtree(0), nodeValues);
-            var g = InterpretRec(node.GetSubtree(1), nodeValues);
-            return f * g;
-          }
-
-        case "-": {
-            if (node.SubtreeCount == 1) {
-              var f = InterpretRec(node.GetSubtree(0), nodeValues);
-              return -f;
-            } else {
-              var f = InterpretRec(node.GetSubtree(0), nodeValues);
-              var g = InterpretRec(node.GetSubtree(1), nodeValues);
-
-              return f - g;
-            }
-          }
-        case "%": {
-            var f = InterpretRec(node.GetSubtree(0), nodeValues);
-            var g = InterpretRec(node.GetSubtree(1), nodeValues);
-
-            // protected division
-            if (g.IsAlmost(0.0)) {
-              return 0;
-            } else {
-              return f / g;
-            }
-          }
-        case "sin": {
-            var f = InterpretRec(node.GetSubtree(0), nodeValues);
-            return Math.Sin(f);
-          }
-        case "cos": {
-            var f = InterpretRec(node.GetSubtree(0), nodeValues);
-            return Math.Cos(f);
-          }
-        case "sqr": {
-            var f = InterpretRec(node.GetSubtree(0), nodeValues);
-            return f * f;
-          }
-        default: {
-            return nodeValues.NodeValue(node);
-          }
-      }
+      if (node.Symbol is Constant || node.Symbol is Variable) {
+        return nodeValues.NodeValue(node);
+      } else if (node.Symbol is Addition) {
+        var f = InterpretRec(node.GetSubtree(0), nodeValues);
+        var g = InterpretRec(node.GetSubtree(1), nodeValues);
+        return f + g;
+      } else if (node.Symbol is Multiplication) {
+        var f = InterpretRec(node.GetSubtree(0), nodeValues);
+        var g = InterpretRec(node.GetSubtree(1), nodeValues);
+        return f * g;
+      } else if (node.Symbol is Subtraction) {
+        if (node.SubtreeCount == 1) {
+          var f = InterpretRec(node.GetSubtree(0), nodeValues);
+          return -f;
+        } else {
+          var f = InterpretRec(node.GetSubtree(0), nodeValues);
+          var g = InterpretRec(node.GetSubtree(1), nodeValues);
+
+          return f - g;
+        }
+      } else if (node.Symbol is Division) {
+        var f = InterpretRec(node.GetSubtree(0), nodeValues);
+        var g = InterpretRec(node.GetSubtree(1), nodeValues);
+
+        // protected division
+        if (g.IsAlmost(0.0)) {
+          return 0;
+        } else {
+          return f / g;
+        }
+      } else if (node.Symbol is Sine) {
+        var f = InterpretRec(node.GetSubtree(0), nodeValues);
+        return Math.Sin(f);
+      } else if (node.Symbol is Cosine) {
+        var f = InterpretRec(node.GetSubtree(0), nodeValues);
+        return Math.Cos(f);
+      } else if (node.Symbol is Square) {
+        var f = InterpretRec(node.GetSubtree(0), nodeValues);
+        return f * f;
+      } else throw new NotSupportedException("unsupported symbol");
     }
 
@@ -1165,73 +1179,63 @@
       double f, g;
       Vector df, dg;
-      switch (node.Symbol.Name) {
-        case "+": {
-            InterpretRec(node.GetSubtree(0), nodeValues, out f, out df);
-            InterpretRec(node.GetSubtree(1), nodeValues, out g, out dg);
-            z = f + g;
-            dz = df + dg; // Vector.AddTo(gl, gr);
-            break;
-          }
-        case "*": {
-            InterpretRec(node.GetSubtree(0), nodeValues, out f, out df);
-            InterpretRec(node.GetSubtree(1), nodeValues, out g, out dg);
-            z = f * g;
-            dz = df * g + f * dg;  // Vector.AddTo(gl.Scale(fr), gr.Scale(fl)); // f'*g + f*g'
-            break;
-          }
-
-        case "-": {
-            if (node.SubtreeCount == 1) {
-              InterpretRec(node.GetSubtree(0), nodeValues, out f, out df);
-              z = -f;
-              dz = -df;
-            } else {
-              InterpretRec(node.GetSubtree(0), nodeValues, out f, out df);
-              InterpretRec(node.GetSubtree(1), nodeValues, out g, out dg);
-
-              z = f - g;
-              dz = df - dg;
-            }
-            break;
-          }
-        case "%": {
-            InterpretRec(node.GetSubtree(0), nodeValues, out f, out df);
-            InterpretRec(node.GetSubtree(1), nodeValues, out g, out dg);
-
-            // protected division
-            if (g.IsAlmost(0.0)) {
-              z = 0;
-              dz = Vector.Zero;
-            } else {
-              z = f / g;
-              dz = -f / (g * g) * dg + df / g; // Vector.AddTo(dg.Scale(f * -1.0 / (g * g)), df.Scale(1.0 / g)); // (f/g)' = f * (1/g)' + 1/g * f' = f * -1/g² * g' + 1/g * f'
-            }
-            break;
-          }
-        case "sin": {
-            InterpretRec(node.GetSubtree(0), nodeValues, out f, out df);
-            z = Math.Sin(f);
-            dz = Math.Cos(f) * df;
-            break;
-          }
-        case "cos": {
-            InterpretRec(node.GetSubtree(0), nodeValues, out f, out df);
-            z = Math.Cos(f);
-            dz = -Math.Sin(f) * df;
-            break;
-          }
-        case "sqr": {
-            InterpretRec(node.GetSubtree(0), nodeValues, out f, out df);
-            z = f * f;
-            dz = 2.0 * f * df;
-            break;
-          }
-        default: {
-            z = nodeValues.NodeValue(node);
-            dz = Vector.CreateNew(nodeValues.NodeGradient(node));
-            break;
-          }
-      }
-    }
+      if (node.Symbol is Constant || node.Symbol is Variable) {
+        z = nodeValues.NodeValue(node);
+        dz = Vector.CreateNew(nodeValues.NodeGradient(node));
+      } else if (node.Symbol is Addition) {
+        InterpretRec(node.GetSubtree(0), nodeValues, out f, out df);
+        InterpretRec(node.GetSubtree(1), nodeValues, out g, out dg);
+        z = f + g;
+        dz = df + dg; // Vector.AddTo(gl, gr);
+
+      } else if (node.Symbol is Multiplication) {
+        InterpretRec(node.GetSubtree(0), nodeValues, out f, out df);
+        InterpretRec(node.GetSubtree(1), nodeValues, out g, out dg);
+        z = f * g;
+        dz = df * g + f * dg;  // Vector.AddTo(gl.Scale(fr), gr.Scale(fl)); // f'*g + f*g'
+
+      } else if (node.Symbol is Subtraction) {
+        if (node.SubtreeCount == 1) {
+          InterpretRec(node.GetSubtree(0), nodeValues, out f, out df);
+          z = -f;
+          dz = -df;
+        } else {
+          InterpretRec(node.GetSubtree(0), nodeValues, out f, out df);
+          InterpretRec(node.GetSubtree(1), nodeValues, out g, out dg);
+
+          z = f - g;
+          dz = df - dg;
+        }
+
+      } else if (node.Symbol is Division) {
+        InterpretRec(node.GetSubtree(0), nodeValues, out f, out df);
+        InterpretRec(node.GetSubtree(1), nodeValues, out g, out dg);
+
+        // protected division
+        if (g.IsAlmost(0.0)) {
+          z = 0;
+          dz = Vector.Zero;
+        } else {
+          z = f / g;
+          dz = -f / (g * g) * dg + df / g; // Vector.AddTo(dg.Scale(f * -1.0 / (g * g)), df.Scale(1.0 / g)); // (f/g)' = f * (1/g)' + 1/g * f' = f * -1/g² * g' + 1/g * f'
+        }
+
+      } else if (node.Symbol is Sine) {
+        InterpretRec(node.GetSubtree(0), nodeValues, out f, out df);
+        z = Math.Sin(f);
+        dz = Math.Cos(f) * df;
+
+      } else if (node.Symbol is Cosine) {
+        InterpretRec(node.GetSubtree(0), nodeValues, out f, out df);
+        z = Math.Cos(f);
+        dz = -Math.Sin(f) * df;
+      } else if (node.Symbol is Square) {
+        InterpretRec(node.GetSubtree(0), nodeValues, out f, out df);
+        z = f * f;
+        dz = 2.0 * f * df;
+      } else {
+        throw new NotSupportedException("unsupported symbol");
+      }
+    }
+
     #endregion
 
@@ -1312,24 +1316,25 @@
 
     private void InitAllParameters() {
-      UpdateTargetVariables(); // implicitly updates the grammar and the encoding      
+      UpdateTargetVariables(); // implicitly updates the grammar and the encoding
     }
 
     private ReadOnlyCheckedItemList<StringValue> CreateFunctionSet() {
       var l = new CheckedItemList<StringValue>();
-      l.Add(new StringValue("+").AsReadOnly());
-      l.Add(new StringValue("*").AsReadOnly());
-      l.Add(new StringValue("%").AsReadOnly());
-      l.Add(new StringValue("-").AsReadOnly());
-      l.Add(new StringValue("sin").AsReadOnly());
-      l.Add(new StringValue("cos").AsReadOnly());
-      l.Add(new StringValue("sqr").AsReadOnly());
+      l.Add(new StringValue("Addition").AsReadOnly());
+      l.Add(new StringValue("Multiplication").AsReadOnly());
+      l.Add(new StringValue("Division").AsReadOnly());
+      l.Add(new StringValue("Subtraction").AsReadOnly());
+      l.Add(new StringValue("Sine").AsReadOnly());
+      l.Add(new StringValue("Cosine").AsReadOnly());
+      l.Add(new StringValue("Square").AsReadOnly());
       return l.AsReadOnly();
     }
 
     private static bool IsConstantNode(ISymbolicExpressionTreeNode n) {
-      return n.Symbol.Name[0] == 'θ';
+      // return n.Symbol.Name[0] == 'θ';
+      return n is ConstantTreeNode;
     }
     private static double GetConstantValue(ISymbolicExpressionTreeNode n) {
-      return 0.0; // TODO: needs to be updated when we write back values to the tree
+      return ((ConstantTreeNode)n).Value;
     }
     private static bool IsLatentVariableNode(ISymbolicExpressionTreeNode n) {
@@ -1340,7 +1345,6 @@
     }
     private static string GetVariableName(ISymbolicExpressionTreeNode n) {
-      return n.Symbol.Name;
-    }
-
+      return ((VariableTreeNode)n).VariableName;
+    }
 
     private void UpdateTargetVariables() {
@@ -1375,107 +1379,93 @@
 
     private ISymbolicExpressionGrammar CreateGrammar() {
-      // whenever ProblemData is changed we create a new grammar with the necessary symbols
-      var g = new SimpleSymbolicExpressionGrammar();
-      var unaryFunc = new string[] { "sin", "cos", "sqr" };
-      var binaryFunc = new string[] { "+", "-", "*", "%" };
-      foreach (var func in unaryFunc) {
-        if (FunctionSet.CheckedItems.Any(ci => ci.Value.Value == func)) g.AddSymbol(func, 1, 1);
-      }
-      foreach (var func in binaryFunc) {
-        if (FunctionSet.CheckedItems.Any(ci => ci.Value.Value == func)) g.AddSymbol(func, 2, 2);
-      }
-
-      foreach (var variableName in ProblemData.AllowedInputVariables.Union(TargetVariables.CheckedItems.Select(i => i.Value.Value)))
-        g.AddTerminalSymbol(variableName);
-
-      // generate symbols for numeric parameters for which the value is optimized using AutoDiff
-      // we generate multiple symbols to balance the probability for selecting a numeric parameter in the generation of random trees
-      var numericConstantsFactor = 2.0;
-      for (int i = 0; i < numericConstantsFactor * (ProblemData.AllowedInputVariables.Count() + TargetVariables.CheckedItems.Count()); i++) {
-        g.AddTerminalSymbol("θ" + i); // numeric parameter for which the value is optimized using AutoDiff
-      }
-
-      // generate symbols for latent variables
-      for (int i = 1; i <= NumberOfLatentVariables; i++) {
-        g.AddTerminalSymbol("λ" + i); // numeric parameter for which the value is optimized using AutoDiff
-      }
-
-      return g;
-    }
-
-
-
-
-
-    private ISymbolicExpressionTreeNode FixParameters(ISymbolicExpressionTreeNode n, double[] parameterValues, ref int nextParIdx) {
-      ISymbolicExpressionTreeNode translatedNode = null;
-      if (n.Symbol is StartSymbol) {
-        translatedNode = new StartSymbol().CreateTreeNode();
-      } else if (n.Symbol is ProgramRootSymbol) {
-        translatedNode = new ProgramRootSymbol().CreateTreeNode();
-      } else if (n.Symbol.Name == "+") {
-        translatedNode = new SimpleSymbol("+", 2).CreateTreeNode();
-      } else if (n.Symbol.Name == "-") {
-        translatedNode = new SimpleSymbol("-", 2).CreateTreeNode();
-      } else if (n.Symbol.Name == "*") {
-        translatedNode = new SimpleSymbol("*", 2).CreateTreeNode();
-      } else if (n.Symbol.Name == "%") {
-        translatedNode = new SimpleSymbol("%", 2).CreateTreeNode();
-      } else if (n.Symbol.Name == "sin") {
-        translatedNode = new SimpleSymbol("sin", 1).CreateTreeNode();
-      } else if (n.Symbol.Name == "cos") {
-        translatedNode = new SimpleSymbol("cos", 1).CreateTreeNode();
-      } else if (n.Symbol.Name == "sqr") {
-        translatedNode = new SimpleSymbol("sqr", 1).CreateTreeNode();
-      } else if (IsConstantNode(n)) {
-        translatedNode = new SimpleSymbol("c_" + nextParIdx, 0).CreateTreeNode();
-        nextParIdx++;
-      } else {
-        translatedNode = new SimpleSymbol(n.Symbol.Name, n.SubtreeCount).CreateTreeNode();
-      }
-      foreach (var child in n.Subtrees) {
-        translatedNode.AddSubtree(FixParameters(child, parameterValues, ref nextParIdx));
-      }
-      return translatedNode;
-    }
-
-
-    private ISymbolicExpressionTreeNode TranslateTreeNode(ISymbolicExpressionTreeNode n, double[] parameterValues, ref int nextParIdx) {
-      ISymbolicExpressionTreeNode translatedNode = null;
-      if (n.Symbol is StartSymbol) {
-        translatedNode = new StartSymbol().CreateTreeNode();
-      } else if (n.Symbol is ProgramRootSymbol) {
-        translatedNode = new ProgramRootSymbol().CreateTreeNode();
-      } else if (n.Symbol.Name == "+") {
-        translatedNode = new Addition().CreateTreeNode();
-      } else if (n.Symbol.Name == "-") {
-        translatedNode = new Subtraction().CreateTreeNode();
-      } else if (n.Symbol.Name == "*") {
-        translatedNode = new Multiplication().CreateTreeNode();
-      } else if (n.Symbol.Name == "%") {
-        translatedNode = new Division().CreateTreeNode();
-      } else if (n.Symbol.Name == "sin") {
-        translatedNode = new Sine().CreateTreeNode();
-      } else if (n.Symbol.Name == "cos") {
-        translatedNode = new Cosine().CreateTreeNode();
-      } else if (n.Symbol.Name == "sqr") {
-        translatedNode = new Square().CreateTreeNode();
-      } else if (IsConstantNode(n)) {
-        var constNode = (ConstantTreeNode)new Constant().CreateTreeNode();
-        constNode.Value = parameterValues[nextParIdx];
-        nextParIdx++;
-        translatedNode = constNode;
-      } else {
-        // assume a variable name
-        var varName = n.Symbol.Name;
-        var varNode = (VariableTreeNode)new Variable().CreateTreeNode();
-        varNode.Weight = 1.0;
-        varNode.VariableName = varName;
-        translatedNode = varNode;
-      }
-      foreach (var child in n.Subtrees) {
-        translatedNode.AddSubtree(TranslateTreeNode(child, parameterValues, ref nextParIdx));
-      }
-      return translatedNode;
+      var grammar = new TypeCoherentExpressionGrammar();
+      grammar.StartGrammarManipulation();
+
+      var problemData = ProblemData;
+      var ds = problemData.Dataset;
+      grammar.MaximumFunctionArguments = 0;
+      grammar.MaximumFunctionDefinitions = 0;
+      var allowedVariables = problemData.AllowedInputVariables.Concat(TargetVariables.CheckedItems.Select(chk => chk.Value.Value));
+      foreach (var varSymbol in grammar.Symbols.OfType<HeuristicLab.Problems.DataAnalysis.Symbolic.VariableBase>()) {
+        if (!varSymbol.Fixed) {
+          varSymbol.AllVariableNames = problemData.InputVariables.Select(x => x.Value).Where(x => ds.VariableHasType<double>(x));
+          varSymbol.VariableNames = allowedVariables.Where(x => ds.VariableHasType<double>(x));
+        }
+      }
+      foreach (var factorSymbol in grammar.Symbols.OfType<BinaryFactorVariable>()) {
+        if (!factorSymbol.Fixed) {
+          factorSymbol.AllVariableNames = problemData.InputVariables.Select(x => x.Value).Where(x => ds.VariableHasType<string>(x));
+          factorSymbol.VariableNames = problemData.AllowedInputVariables.Where(x => ds.VariableHasType<string>(x));
+          factorSymbol.VariableValues = factorSymbol.VariableNames
+            .ToDictionary(varName => varName, varName => ds.GetStringValues(varName).Distinct().ToList());
+        }
+      }
+      foreach (var factorSymbol in grammar.Symbols.OfType<FactorVariable>()) {
+        if (!factorSymbol.Fixed) {
+          factorSymbol.AllVariableNames = problemData.InputVariables.Select(x => x.Value).Where(x => ds.VariableHasType<string>(x));
+          factorSymbol.VariableNames = problemData.AllowedInputVariables.Where(x => ds.VariableHasType<string>(x));
+          factorSymbol.VariableValues = factorSymbol.VariableNames
+            .ToDictionary(varName => varName,
+            varName => ds.GetStringValues(varName).Distinct()
+            .Select((n, i) => Tuple.Create(n, i))
+            .ToDictionary(tup => tup.Item1, tup => tup.Item2));
+        }
+      }
+
+      grammar.ConfigureAsDefaultRegressionGrammar();
+      grammar.GetSymbol("Logarithm").Enabled = false; // not supported yet
+      grammar.GetSymbol("Exponential").Enabled = false; // not supported yet
+
+      // configure initialization of constants
+      var constSy = (Constant)grammar.GetSymbol("Constant");
+      // max and min are only relevant for initialization
+      constSy.MaxValue = +1.0e-1; // small initial values for constant opt
+      constSy.MinValue = -1.0e-1;
+      constSy.MultiplicativeManipulatorSigma = 1.0; // allow large jumps for manipulation
+      constSy.ManipulatorMu = 0.0;
+      constSy.ManipulatorSigma = 1.0; // allow large jumps
+
+      // configure initialization of variables
+      var varSy = (Variable)grammar.GetSymbol("Variable");
+      // fix variable weights to 1.0
+      varSy.WeightMu = 1.0;
+      varSy.WeightSigma = 0.0;
+      varSy.WeightManipulatorMu = 0.0;
+      varSy.WeightManipulatorSigma = 0.0;
+      varSy.MultiplicativeWeightManipulatorSigma = 0.0;
+
+      foreach (var f in FunctionSet) {
+        grammar.GetSymbol(f.Value).Enabled = FunctionSet.ItemChecked(f);
+      }
+
+      grammar.FinishedGrammarManipulation();
+      return grammar;
+      // // whenever ProblemData is changed we create a new grammar with the necessary symbols
+      // var g = new SimpleSymbolicExpressionGrammar();
+      // var unaryFunc = new string[] { "sin", "cos", "sqr" };
+      // var binaryFunc = new string[] { "+", "-", "*", "%" };
+      // foreach (var func in unaryFunc) {
+      //   if (FunctionSet.CheckedItems.Any(ci => ci.Value.Value == func)) g.AddSymbol(func, 1, 1);
+      // }
+      // foreach (var func in binaryFunc) {
+      //   if (FunctionSet.CheckedItems.Any(ci => ci.Value.Value == func)) g.AddSymbol(func, 2, 2);
+      // }
+      // 
+      // foreach (var variableName in ProblemData.AllowedInputVariables.Union(TargetVariables.CheckedItems.Select(i => i.Value.Value)))
+      //   g.AddTerminalSymbol(variableName);
+      // 
+      // // generate symbols for numeric parameters for which the value is optimized using AutoDiff
+      // // we generate multiple symbols to balance the probability for selecting a numeric parameter in the generation of random trees
+      // var numericConstantsFactor = 2.0;
+      // for (int i = 0; i < numericConstantsFactor * (ProblemData.AllowedInputVariables.Count() + TargetVariables.CheckedItems.Count()); i++) {
+      //   g.AddTerminalSymbol("θ" + i); // numeric parameter for which the value is optimized using AutoDiff
+      // }
+      // 
+      // // generate symbols for latent variables
+      // for (int i = 1; i <= NumberOfLatentVariables; i++) {
+      //   g.AddTerminalSymbol("λ" + i); // numeric parameter for which the value is optimized using AutoDiff
+      // }
+      // 
+      // return g;
     }
     #endregion
@@ -1552,5 +1542,5 @@
             }
             nodes.Add(node);
-            SetVariableValue(varName, 0.0);
+            SetVariableValue(varName, 0.0);  // this value is updated in the prediction loop
           }
         }
@@ -1574,5 +1564,5 @@
           nodes.ForEach(n => node2val[n] = Tuple.Create(val, dVal));
         } else {
-          var fakeNode = new SimpleSymbol(variableName, 0).CreateTreeNode();
+          var fakeNode = new VariableTreeNode(new Variable());
           var newNodeList = new List<ISymbolicExpressionTreeNode>();
           newNodeList.Add(fakeNode);
Index: branches/2925_AutoDiffForDynamicalModels/HeuristicLab.Problems.DynamicalSystemsModelling/3.3/Solution.cs
===================================================================
--- branches/2925_AutoDiffForDynamicalModels/HeuristicLab.Problems.DynamicalSystemsModelling/3.3/Solution.cs	(revision 16601)
+++ branches/2925_AutoDiffForDynamicalModels/HeuristicLab.Problems.DynamicalSystemsModelling/3.3/Solution.cs	(revision 16602)
@@ -8,4 +8,5 @@
 using HeuristicLab.Persistence.Default.CompositeSerializers.Storable;
 using HeuristicLab.Problems.DataAnalysis;
+using HeuristicLab.Problems.DataAnalysis.Symbolic;
 using HeuristicLab.Random;
 
@@ -88,9 +89,12 @@
       var forecastEpisode = new IntRange(episode.Start, episode.End + forecastHorizon);
 
-      double[] optL0;
       var random = new FastRandom(12345);
-      Problem.OptimizeForEpisodes(trees, problemData, targetVars, latentVariables, random, new[] { forecastEpisode }, 100, numericIntegrationSteps, odeSolver, out optL0, out snmse);
+      snmse = Problem.OptimizeForEpisodes(trees, problemData, targetVars, latentVariables, random, new[] { forecastEpisode }, 100, numericIntegrationSteps, odeSolver);
       var optimizationData = new Problem.OptimizationData(trees, targetVars, problemData, null, new[] { forecastEpisode }, numericIntegrationSteps, latentVariables, odeSolver);
-      var predictions = Problem.Integrate(optimizationData, optL0).ToArray();
+
+
+      var theta = Problem.ExtractParametersFromTrees(trees);
+
+      var predictions = Problem.Integrate(optimizationData, theta).ToArray();
       return predictions.Select(p => p.Select(pi => pi.Item1).ToArray()).ToArray();
     }
