Fixes in Division_Node operation

* There were problems in both MATLAB and Python implementations added
  in the previous commit
* Added a bug folder including some scripts to reproduce the precision
  problems reported in previous issues #27

Author: Ezetowers

MLC/Common/Lisp_Tree_Expr/Operation_Nodes.py

@@ -1,4 +1,3 @@
-import math
 import importlib
 
 import MLC.Log.log as lg
@@ -108,11 +107,11 @@ def formal(self):
 
     def _process_division(self, dividend, divisor):
         if type(divisor) == np.ndarray:
-            divisor = [Division_Node.PROTECTION if abs(x) < Division_Node.PROTECTION else x for x in divisor]
-            return np.sign(divisor) * dividend / np.asarray(divisor)
+            new_divisor = [Division_Node.PROTECTION if np.abs(x) < Division_Node.PROTECTION else np.abs(x) for x in divisor]
+            return np.sign(divisor) * dividend / np.asarray(new_divisor)
         else:
             if abs(divisor) < Division_Node.PROTECTION:
-                return dividend / Division_Node.PROTECTION
+                return np.sign(divisor) * dividend / Division_Node.PROTECTION
 
         return dividend / divisor
 
@@ -149,7 +148,7 @@ def formal(self):
 
     def op_simplify(self):
         if not self._nodes[0].is_sensor():
-            arg = math.sin(float(self._nodes[0].to_string()))
+            arg = np.sin(float(self._nodes[0].to_string()))
             return Leaf_Node(process_float(arg))
         else:
             return self
@@ -168,7 +167,7 @@ def formal(self):
 
     def op_simplify(self):
         if not self._nodes[0].is_sensor():
-            arg = math.cos(float(self._nodes[0].to_string()))
+            arg = np.cos(float(self._nodes[0].to_string()))
             return Leaf_Node(process_float(arg))
         else:
             return self
@@ -189,7 +188,7 @@ def formal(self):
 
     def _process_arg(self, arg):
         if type(arg) == np.ndarray:
-            return [Logarithm_Node.PROTECTION if abs(x) < Logarithm_Node.PROTECTION else abs(x) for x in arg]
+            return [Logarithm_Node.PROTECTION if np.abs(x) < Logarithm_Node.PROTECTION else np.abs(x) for x in arg]
         else:
             if abs(arg) < Logarithm_Node.PROTECTION:
                 return Logarithm_Node.PROTECTION
@@ -199,9 +198,9 @@ def _process_arg(self, arg):
     def op_simplify(self):
         if not self._nodes[0].is_sensor():
             if float(self._nodes[0].to_string()) < Logarithm_Node.SIMPLIFY_PROTECTION:
-                arg = math.log(Logarithm_Node.SIMPLIFY_PROTECTION)
+                arg = np.log(Logarithm_Node.SIMPLIFY_PROTECTION)
             else:
-                arg = math.log(float(self._nodes[0].to_string()))
+                arg = np.log(float(self._nodes[0].to_string()))
 
             return Leaf_Node(process_float(arg))
         else:
@@ -223,9 +222,9 @@ def op_simplify(self):
         if not self._nodes[0].is_sensor():
             lg.logger_.debug("[EXP NODE] Value: " + self._nodes[0].to_string())
             try:
-                arg = math.exp(float(self._nodes[0].to_string()))
+                arg = np.exp(float(self._nodes[0].to_string()))
             except OverflowError:
-                # FIXME: See what to do with this expression, because there are problems with
+                # FIXME: See what to do with this expression, because there are problems when
                 # an infinite value is the argumento of a sinusoidal function
                 return Leaf_Node(process_float(float("inf")))
 
@@ -247,7 +246,7 @@ def formal(self):
 
     def op_simplify(self):
         if not self._nodes[0].is_sensor():
-            arg = math.tanh(float(self._nodes[0].to_string()))
+            arg = np.tanh(float(self._nodes[0].to_string()))
             return Leaf_Node(process_float(arg))
         else:
             return self

MLC/Scripts/Evaluation/toy_problem.py

@@ -8,11 +8,7 @@
 
 
 def individual_data(indiv):
-    x = np.arange(-10, 10 + 0.1, 0.1)
-    # FIXME: Numpy precision it's driving me up the wall. It doesn't put a zero
-    # in the middle element of the array (put a extremely small number instead).
-    # That generates problems with the division operation
-    x[len(x) / 2] = 0
+    x = np.linspace(-10.0, 10.0, num=201)
     y = np.tanh(x**3 - x**2 - 1)
 
     eng = MatlabEngine.engine()

MLC/Scripts/Evaluation/toy_problem_python_ev.py

@@ -7,11 +7,7 @@
 
 
 def individual_data(indiv):
-    x = np.arange(-10, 10 + 0.1, 0.1)
-    # FIXME: Numpy precision it's driving me up the wall. It doesn't put a zero
-    # in the middle element of the array (put a extremely small number instead).
-    # That generates problems with the division operation
-    x[len(x) / 2] = 0
+    x = np.linspace(-10.0, 10.0, num=201)
     y = np.tanh(x**3 - x**2 - 1)
 
     config = Config.get_instance()

bugs/precision_bug/expr_test.py

@@ -0,0 +1,47 @@
+import sys
+sys.path.append("../..")
+import binascii
+import numpy as np
+import struct
+import MLC.Log.log as lg
+
+from MLC.Common.Lisp_Tree_Expr.Lisp_Tree_Expr import Lisp_Tree_Expr
+from MLC.Log.log import set_logger
+from MLC.mlc_parameters.mlc_parameters import Config
+
+
+def initialize_config():
+    config = Config.get_instance()
+    config.read('../conf/configuration.ini')
+    return config
+
+# Set printable resolution (don't alter numpy interval resolution)
+np.set_printoptions(precision=9)
+# Show full arrays, no matter what size do they have
+np.set_printoptions(threshold=np.inf)
+# Don't show scientific notation
+np.set_printoptions(suppress=True)
+
+initialize_config()
+set_logger('console')
+
+# Full expression
+expr6 = "(root (cos (exp (- -6.3726 (* -7.1746 S0)))))"
+expr61 = "(root (exp (- -6.3726 (* -7.1746 S0))))"
+expr612 = "(root (- -6.3726 (* -7.1746 S0)))"
+
+tree = Lisp_Tree_Expr(expr612)
+x = np.linspace(-10.0, 10.0, num=201)
+mlc_y = tree.calculate_expression([x])
+
+# Calculate the Mean squared error
+y = np.tanh(x**3 - x**2 - 1)
+evaluation = float(np.sum((mlc_y - y)**2))
+
+# print mlc_y
+with open("./costs_python.txt", "w") as f:
+    for elem in mlc_y:
+        f.write("%50.50f\n" % elem)
+
+print evaluation
+print np.sum(mlc_y)

bugs/precision_bug/matlab_expr_test.m

@@ -0,0 +1,19 @@
+close all;
+clear all;
+
+%--------------------- %
+
+x = -10:0.1:10;
+% m = 'exp(((-6.3726) - ((-7.1746) .* S0)))';
+m = '((-6.3726) - ((-7.1746) .* S0))';
+y = tanh(x.^3-x.^2-1);
+m = strrep(m, 'S0', 'x'); y2=eval(m); J=sum((y2-y).^2); fprintf('%.10f\n', J);
+
+filename = './costs_matlab.txt';
+file = fopen(filename, 'w');
+
+for i=1:length(x)
+    fprintf(file, '%50.50f\n', y2(i));
+end
+
+fclose(file);

main.py

@@ -22,7 +22,7 @@ def main():
 
     eng = MatlabEngine.engine()
     config = initialize_config()
-    # MatlabEngine.load_random_values("./tests/integration_tests/matlab_randoms.txt")
+    MatlabEngine.load_random_values("./tests/integration_tests/matlab_randoms.txt")
 
     # Create the MLC2 object and store it in the workspace. With this
     # feature we will be able to call every function of the MATLAB code

matlab_code/MLC_tools/my_div.m

@@ -8,7 +8,7 @@
         b=arg1 ./ arg2;
     else
         inverse_protect_array = ~protect_array;
-        valid_elements = arg2 .* inverse_protect_array;
+        valid_elements = abs(arg2) .* inverse_protect_array;
         non_valid_elements = protect_array .* protection;
         new_arg2 = valid_elements + non_valid_elements;
         b = sign_array .* arg1 ./ new_arg2;

tests/integration_tests/test_basic/configuration.ini

@@ -51,8 +51,8 @@ cascade = 1,1
 # evaluation_method = standalone_function
 # evaluation_method = standalone_files
 evaluation_method = mfile_standalone
-# evaluation_function = toy_problem
-evaluation_function = toy_problem_python_ev
+evaluation_function = toy_problem
+# evaluation_function = toy_problem_python_ev
 indfile = ind.dat
 Jfile = J.dat
 # exchangedir = fullfile(pwd,evaluator0)