Changeset 16905 for trunk/HeuristicLab.Problems.DataAnalysis.Symbolic

Timestamp:

05/07/19 14:20:33 (5 years ago)

Author:

gkronber

Message:

#2915: Corrected calculation of cuberoot function in all interpreters and updated all formatters.
I tested a run with all interpreters and got the same results with all of them.

Location:

trunk/HeuristicLab.Problems.DataAnalysis.Symbolic/3.4

Files:

• Converters/DerivativeCalculator.cs (modified) (1 diff)
• Converters/TreeToAutoDiffTermConverter.cs (modified) (2 diffs)
• Formatters/SymbolicDataAnalysisExpressionCSharpFormatter.cs (modified) (3 diffs)
• Formatters/SymbolicDataAnalysisExpressionExcelFormatter.cs (modified) (1 diff)
• Formatters/SymbolicDataAnalysisExpressionLatexFormatter.cs (modified) (2 diffs)
• Formatters/SymbolicDataAnalysisExpressionMATLABFormatter.cs (modified) (1 diff)
• Formatters/SymbolicDataAnalysisExpressionMathematicaFormatter.cs (modified) (1 diff)
• Formatters/SymbolicDataAnalysisExpressionSmalltalkFormatter.cs (modified) (1 diff)
• Interpreter/BatchOperations.cs (modified) (1 diff)
• Interpreter/SymbolicDataAnalysisExpressionCompiledTreeInterpreter.cs (modified) (2 diffs)
• Interpreter/SymbolicDataAnalysisExpressionTreeILEmittingInterpreter.cs (modified) (1 diff)
• Interpreter/SymbolicDataAnalysisExpressionTreeInterpreter.cs (modified) (1 diff)
• Interpreter/SymbolicDataAnalysisExpressionTreeLinearInterpreter.cs (modified) (1 diff)

trunk/HeuristicLab.Problems.DataAnalysis.Symbolic/3.4/Converters/DerivativeCalculator.cs

@@ -143 +143 @@
         var f = (ISymbolicExpressionTreeNode)branch.Clone();
         var u = (ISymbolicExpressionTreeNode)branch.GetSubtree(0).Clone();
-        return Product(Div(CreateConstant(1.0), Product(CreateConstant(3.0), Square(f))), Derive(u, variableName));
+        return Product(Div(CreateConstant(1.0), Product(CreateConstant(3.0), Square(f))), Derive(u, variableName));  // 1/3 1/cbrt(f(x))^2 d/dx f(x)
       }
       if (branch.Symbol is Cube) {

trunk/HeuristicLab.Problems.DataAnalysis.Symbolic/3.4/Converters/TreeToAutoDiffTermConverter.cs

@@ -92 +92 @@
       );
 
+    private static readonly Func<Term, UnaryFunc> cbrt = UnaryFunc.Factory(
+      eval: x => x < 0 ? -Math.Pow(-x, 1.0 / 3) : Math.Pow(x, 1.0 / 3),
+      diff: x => { var cbrt_x = x < 0 ? -Math.Pow(-x, 1.0 / 3) : Math.Pow(x, 1.0 / 3); return 1.0 / (3 * cbrt_x * cbrt_x); }
+      );
+
+
+
     #endregion
 
@@ -250 +257 @@
       }
       if (node.Symbol is CubeRoot) {
-        return AutoDiff.TermBuilder.Power(
-          ConvertToAutoDiff(node.GetSubtree(0)), 1.0/3.0);
+        return cbrt(ConvertToAutoDiff(node.GetSubtree(0)));
       }
       if (node.Symbol is Sine) {

trunk/HeuristicLab.Problems.DataAnalysis.Symbolic/3.4/Formatters/SymbolicDataAnalysisExpressionCSharpFormatter.cs

@@ -140 +140 @@
           FormatPower(node, strBuilder, "3");
         } else if (node.Symbol is CubeRoot) {
-          FormatPower(node, strBuilder, "1.0/3");
+          strBuilder.Append("Cbrt(");
+          FormatRecursively(node.GetSubtree(0), strBuilder);
+          strBuilder.Append(")");
         } else if (node.Symbol is Power) {
           FormatFunction(node, "Math.Pow", strBuilder);
@@ -252 +254 @@
       strBuilder.AppendLine("public static class Model {");
       GenerateAverageSource(strBuilder);
+      GenerateCbrtSource(strBuilder);
       GenerateIfThenElseSource(strBuilder);
       GenerateFactorSource(strBuilder);
@@ -293 +296 @@
       strBuilder.AppendLine("private static double Average(params double[] values) {");
       strBuilder.AppendLine("  return values.Average();");
+      strBuilder.AppendLine("}");
+    }
+    private void GenerateCbrtSource(StringBuilder strBuilder) {
+      strBuilder.AppendLine("private static double Cbrt(double x) {");
+      strBuilder.AppendLine("  return x < 0 ? -Math.Pow(-x, 1.0 / 3.0) : Math.Pow(x, 1.0 / 3.0);");
       strBuilder.AppendLine("}");
     }

trunk/HeuristicLab.Problems.DataAnalysis.Symbolic/3.4/Formatters/SymbolicDataAnalysisExpressionExcelFormatter.cs

@@ -144 +144 @@
         stringBuilder.Append($"POWER({FormatRecursively(node.GetSubtree(0))}, 3)");
       } else if (symbol is CubeRoot) {
-        stringBuilder.Append($"POWER({FormatRecursively(node.GetSubtree(0))}, 1/3)");
+        var arg_expr = FormatRecursively(node.GetSubtree(0));
+        stringBuilder.Append($"IF({arg_expr} < 0, -POWER(-{arg_expr}, 1/3), POWER({arg_expr}, 1/3))");
       } else if (symbol is Division) {
         if (node.SubtreeCount == 1) {

trunk/HeuristicLab.Problems.DataAnalysis.Symbolic/3.4/Formatters/SymbolicDataAnalysisExpressionLatexFormatter.cs

@@ -138 +138 @@
         strBuilder.Append(@"\left(");
       } else if (node.Symbol is CubeRoot) {
-        strBuilder.Append(@"\left(");
+        strBuilder.Append(@"\operatorname{cbrt}\left(");
       } else if (node.Symbol is Sine) {
         strBuilder.Append(@"\sin \left( ");
@@ -420 +420 @@
         strBuilder.Append(@"\right)^3");
       } else if (node.Symbol is CubeRoot) {
-        strBuilder.Append(@"\right)^\frac{1}{3}");
+        strBuilder.Append(@"\right)");
       } else if (node.Symbol is Sine) {
         strBuilder.Append(@" \right) ");

trunk/HeuristicLab.Problems.DataAnalysis.Symbolic/3.4/Formatters/SymbolicDataAnalysisExpressionMATLABFormatter.cs

@@ -188 +188 @@
         stringBuilder.Append(").^3");
       } else if (symbol is CubeRoot) {
-        stringBuilder.Append("(");
-        stringBuilder.Append(FormatRecursively(node.GetSubtree(0)));
-        stringBuilder.Append(").^(1/3)");
+        stringBuilder.Append("NTHROOT(");
+        stringBuilder.Append(FormatRecursively(node.GetSubtree(0)));
+        stringBuilder.Append(", 3)");
       } else if (symbol is GreaterThan) {
         stringBuilder.Append("((");

trunk/HeuristicLab.Problems.DataAnalysis.Symbolic/3.4/Formatters/SymbolicDataAnalysisExpressionMathematicaFormatter.cs

@@ -115 +115 @@
           FormatPower(node, strBuilder, "3");
         } else if (node.Symbol is CubeRoot) {
-          FormatPower(node, strBuilder, "1/3");
+          strBuilder.Append("CubeRoot[");
+          FormatRecursively(node.GetSubtree(0), strBuilder);
+          strBuilder.Append("]");
         } else if (node.Symbol is Power) {
           FormatFunction(node, "Power", strBuilder);

trunk/HeuristicLab.Problems.DataAnalysis.Symbolic/3.4/Formatters/SymbolicDataAnalysisExpressionSmalltalkFormatter.cs

@@ -92 +92 @@
         stringBuilder.Append(FormatPower(node.GetSubtree(0), "3"));
       } else if (symbol is CubeRoot) {
-        stringBuilder.Append(FormatPower(node.GetSubtree(0), "(1/3)"));
+        stringBuilder.Append(FormatRecursively(node.GetSubtree(0)));
+        stringBuilder.Append(" cbrt");
       } else if (symbol is Division) {
         if (node.SubtreeCount == 1) {

trunk/HeuristicLab.Problems.DataAnalysis.Symbolic/3.4/Interpreter/BatchOperations.cs

@@ -101 +101 @@
     public static void CubeRoot(double[] a, double[] b) {
       for (int i = 0; i < BATCHSIZE; ++i)
-        a[i] = Math.Pow(b[i], 1d / 3d);
+        a[i] = b[i] < 0 ? -Math.Pow(-b[i], 1d / 3d) : Math.Pow(b[i], 1d / 3d);
     }
 

trunk/HeuristicLab.Problems.DataAnalysis.Symbolic/3.4/Interpreter/SymbolicDataAnalysisExpressionCompiledTreeInterpreter.cs

@@ -247 +247 @@
         case OpCodes.CubeRoot: {
             var arg = MakeExpr(node.GetSubtree(0), variableIndices, row, columns);
-            return Expression.Power(arg, Expression.Constant(1.0 / 3.0));
+            return Expression.Condition(Expression.LessThan(arg, Expression.Constant(0.0)),
+              Expression.Negate(Expression.Power(Expression.Negate(arg), Expression.Constant(1.0 / 3.0))),
+              Expression.Power(arg, Expression.Constant(1.0 / 3.0)));
           }
         case OpCodes.Root: {
@@ -514 +516 @@
             var x1 = MakeExpr(node.GetSubtree(0), variableIndices, row, columns);
             var x2 = MakeExpr(node.GetSubtree(1), variableIndices, row, columns);
-            return Expression.Divide(x1, 
-              Expression.Call(Sqrt, 
+            return Expression.Divide(x1,
+              Expression.Call(Sqrt,
               Expression.Add(
                 Expression.Constant(1.0),

trunk/HeuristicLab.Problems.DataAnalysis.Symbolic/3.4/Interpreter/SymbolicDataAnalysisExpressionTreeILEmittingInterpreter.cs

@@ -336 +336 @@
         case OpCodes.CubeRoot: {
             CompileInstructions(il, state, ds);
+            var c1 = il.DefineLabel();
+            var end = il.DefineLabel();
+
+            il.Emit(System.Reflection.Emit.OpCodes.Dup); // x
+            il.Emit(System.Reflection.Emit.OpCodes.Ldc_R8, 0.0);
+            il.Emit(System.Reflection.Emit.OpCodes.Clt); // x < 0?
+            il.Emit(System.Reflection.Emit.OpCodes.Brfalse, c1);
+            il.Emit(System.Reflection.Emit.OpCodes.Neg); // x = -x
             il.Emit(System.Reflection.Emit.OpCodes.Ldc_R8, 1.0 / 3.0);
             il.Emit(System.Reflection.Emit.OpCodes.Call, power);
+            il.Emit(System.Reflection.Emit.OpCodes.Neg); // -Math.Pow(-x, 1/3)
+            il.Emit(System.Reflection.Emit.OpCodes.Br, end);
+            il.MarkLabel(c1);
+            il.Emit(System.Reflection.Emit.OpCodes.Ldc_R8, 1.0 / 3.0);
+            il.Emit(System.Reflection.Emit.OpCodes.Call, power);
+            il.MarkLabel(end);
             return;
           }

trunk/HeuristicLab.Problems.DataAnalysis.Symbolic/3.4/Interpreter/SymbolicDataAnalysisExpressionTreeInterpreter.cs

@@ -232 +232 @@
           }
         case OpCodes.CubeRoot: {
-            return Math.Pow(Evaluate(dataset, ref row, state), 1.0 / 3.0);
+            var arg = Evaluate(dataset, ref row, state);
+            return arg < 0 ? -Math.Pow(-arg, 1.0 / 3.0) : Math.Pow(arg, 1.0 / 3.0);
           }
         case OpCodes.Root: {

trunk/HeuristicLab.Problems.DataAnalysis.Symbolic/3.4/Interpreter/SymbolicDataAnalysisExpressionTreeLinearInterpreter.cs

@@ -245 +245 @@
           instr.value = Math.Sqrt(code[instr.childIndex].value);
         } else if (instr.opCode == OpCodes.CubeRoot) {
-          instr.value = Math.Pow(code[instr.childIndex].value, 1.0 / 3.0);
+          var arg = code[instr.childIndex].value;
+          instr.value = arg < 0 ? -Math.Pow(-arg, 1.0 / 3.0) : Math.Pow(arg, 1.0 / 3.0);
         } else if (instr.opCode == OpCodes.Root) {
           double x = code[instr.childIndex].value;