Merge pull request #53 from Gromwud/main

Cumulitive update

Author: Maslyaev

epde/operators/common/fitness.py

@@ -65,22 +65,20 @@ def apply(self, objective: Equation, arguments: dict, force_out_of_place: bool =
         """        
         self_args, subop_args = self.parse_suboperator_args(arguments = arguments)
 
-        self.suboperators['sparsity'].apply(objective, subop_args['sparsity'])
+        if force_out_of_place:
+            self.suboperators['sparsity'].apply(objective, subop_args['sparsity'])
         self.suboperators['coeff_calc'].apply(objective, subop_args['coeff_calc'])
 
         _, target, features = objective.evaluate(normalize = False, return_val = False)
-        try:
-            if features is None:
-                discr_feats = 0
-            else:
-                discr_feats = np.dot(features, objective.weights_final[:-1][objective.weights_internal != 0])
+        if features is None:
+            discr_feats = 0
+        else:
+            discr_feats = np.dot(features, objective.weights_final[:-1][objective.weights_internal != 0])
 
-            discr = (discr_feats + np.full(target.shape, objective.weights_final[-1]) - target)
-            self.g_fun_vals = global_var.grid_cache.g_func.reshape(-1)
-            discr = np.multiply(discr, self.g_fun_vals)
-            rl_error = np.linalg.norm(discr, ord = 2)
-        except ValueError:
-            raise ValueError('An error in getting weights ')
+        discr = (discr_feats + np.full(target.shape, objective.weights_final[-1]) - target)
+        self.g_fun_vals = global_var.grid_cache.g_func.reshape(-1)
+        discr = np.multiply(discr, self.g_fun_vals)
+        rl_error = np.linalg.norm(discr, ord = 2)
 
         if not (self.params['penalty_coeff'] > 0. and self.params['penalty_coeff'] < 1.):
             raise ValueError('Incorrect penalty coefficient set, value shall be in (0, 1).')
@@ -124,20 +122,17 @@ def apply(self, objective: Equation, arguments: dict, force_out_of_place: bool =
         _, target, features = objective.evaluate(normalize=False, return_val=False)
 
         self.get_g_fun_vals()
+        data_shape = global_var.grid_cache.g_func.shape
 
-        try:
-            if features is None:
-                discr = target - objective.weights_final[-1]
-            else:
-                discr_feats = np.dot(features, objective.weights_final[:-1][objective.weights_internal != 0])
-                discr_feats = discr_feats + objective.weights_final[-1]
-                discr = discr_feats - target
+        if features is None:
+            discr = target - objective.weights_final[-1]
+        else:
+            discr_feats = np.dot(features, objective.weights_final[:-1][objective.weights_internal != 0])
+            discr_feats = discr_feats + objective.weights_final[-1]
+            discr = discr_feats - target
 
-            discr = np.multiply(discr, self.g_fun_vals) / np.std(target)
-            # discr = np.multiply(discr, self.g_fun_vals) / np.linalg.norm(target, 2)
-            rl_error = np.linalg.norm(discr, 2)
-        except ValueError:
-            raise ValueError('An error in getting weights ')
+        discr = np.multiply(discr, self.g_fun_vals) / np.std(target)
+        rl_error = np.sqrt(np.mean(discr ** 2))
 
         if not (self.params['penalty_coeff'] > 0. and self.params['penalty_coeff'] < 1.):
             raise ValueError('Incorrect penalty coefficient set, value shall be in (0, 1).')
@@ -153,17 +148,17 @@ def apply(self, objective: Equation, arguments: dict, force_out_of_place: bool =
         # assert objective.simplified, 'Trying to evaluate not simplified equation.'
 
         # Calculate r-loss
-        data_shape = global_var.grid_cache.g_func.shape
         target_vals = target.reshape(*data_shape)
         features_vals = []
 
         if target_vals.ndim == 1:
+            horizons_default = 30
             window_size = len(target_vals) // 2
             num_horizons = len(target_vals) - window_size + 1
-            if num_horizons < 30:
+            if num_horizons < horizons_default:
                 step_size = 1
             else:
-                step_size = num_horizons // 30
+                step_size = num_horizons // horizons_default
             eq_window_weights = []
             # Compute coefficients and collect statistics over horizons
             if features is None:
@@ -187,21 +182,23 @@ def apply(self, objective: Equation, arguments: dict, force_out_of_place: bool =
                     eq_window_weights.append(valuable_weights)
                 eq_cv = np.array([
                     np.abs(np.std(_)) if np.isclose(np.mean(_), 0)
+                    # else np.abs(np.std(_) / np.sqrt(np.mean(np.pow(_, 2))))
                     else np.abs(np.std(_) / np.mean(_))
                     for _ in zip(*eq_window_weights)
                 ])
-                lr = eq_cv.mean()
+                lr = eq_cv.sum()
 
         elif target_vals.ndim == 2:
             lr = 0
             for dim in range(target_vals.ndim):
+                horizons_default = 30
                 eq_window_weights = []
                 window_size = target_vals.shape[dim] // 2
                 num_horizons = target_vals.shape[dim] - window_size + 1
-                if num_horizons < 30:
+                if num_horizons < horizons_default:
                     step_size = 1
                 else:
-                    step_size = num_horizons // 30
+                    step_size = num_horizons // horizons_default
                 # Compute coefficients and collect statistics over horizons
                 if features is None:
                     for start_idx in range(0, num_horizons, step_size):
@@ -235,7 +232,7 @@ def apply(self, objective: Equation, arguments: dict, force_out_of_place: bool =
                         else np.abs(np.std(_) / np.mean(_))
                         for _ in zip(*eq_window_weights)
                     ])
-                    lr += eq_cv.mean()
+                    lr += eq_cv.sum()
 
         elif target_vals.ndim == 3:
             lr = 0
@@ -282,7 +279,7 @@ def apply(self, objective: Equation, arguments: dict, force_out_of_place: bool =
                         else np.abs(np.std(_) / np.mean(_))
                         for _ in zip(*eq_window_weights)
                     ])
-                    lr += eq_cv.mean()
+                    lr += eq_cv.sum()
 
         objective.fitness_calculated = True
         objective.fitness_value = fitness_value
@@ -334,8 +331,8 @@ def apply(self, objective : SoEq, arguments : dict, force_out_of_place: bool = F
             net = None
 
         self.set_adapter(net=net)
-
-        self.suboperators['sparsity'].apply(objective, subop_args['sparsity'])
+        if force_out_of_place:
+            self.suboperators['sparsity'].apply(objective, subop_args['sparsity'])
         self.suboperators['coeff_calc'].apply(objective, subop_args['coeff_calc'])
 
         print('solving equation:')
@@ -415,7 +412,8 @@ def apply(self, objective: SoEq, arguments: dict, force_out_of_place: bool = Fal
 
         self.set_adapter(net=net)
 
-        self.suboperators['sparsity'].apply(objective, subop_args['sparsity'])
+        if force_out_of_place:
+            self.suboperators['sparsity'].apply(objective, subop_args['sparsity'])
         self.suboperators['coeff_calc'].apply(objective, subop_args['coeff_calc'])
 
         print('solving equation:')
@@ -443,7 +441,7 @@ def apply(self, objective: SoEq, arguments: dict, force_out_of_place: bool = Fal
                 referential_data = global_var.tensor_cache.get((eq.main_var_to_explain, (1.0,)))
                 discr = solution[..., eq_idx] - referential_data.reshape(solution[..., eq_idx].shape)
                 discr = np.multiply(discr, self.g_fun_vals.reshape(discr.shape)) / np.std(discr)
-                rl_error = np.linalg.norm(discr, ord=2)
+                rl_error = np.sqrt(np.mean(discr ** 2))
 
                 print(f'fitness error is {rl_error}, while loss addition is {float(loss_add)}')
                 lp = rl_error + self.params['pinn_loss_mult'] * float(
@@ -493,7 +491,7 @@ def apply(self, objective: SoEq, arguments: dict, force_out_of_place: bool = Fal
                         else np.abs(np.std(_) / np.mean(_))
                         for _ in zip(*eq_window_weights)
                     ])
-                    lr = eq_cv.mean()
+                    lr = eq_cv.sum()
 
             elif target_vals.ndim == 2:
                 lr = 0
@@ -540,7 +538,7 @@ def apply(self, objective: SoEq, arguments: dict, force_out_of_place: bool = Fal
                             else np.abs(np.std(_) / np.mean(_))
                             for _ in zip(*eq_window_weights)
                         ])
-                        lr += eq_cv.mean()
+                        lr += eq_cv.sum()
 
             elif target_vals.ndim == 3:
                 lr = 0
@@ -590,7 +588,7 @@ def apply(self, objective: SoEq, arguments: dict, force_out_of_place: bool = Fal
                             else np.abs(np.std(_) / np.mean(_))
                             for _ in zip(*eq_window_weights)
                         ])
-                        lr += eq_cv.mean()
+                        lr += eq_cv.sum()
 
             eq.fitness_calculated = True
             eq.fitness_value = lp

epde/operators/multiobjective/moeadd_specific.py

@@ -366,12 +366,14 @@ def apply(self, objective: ParetoLevels, arguments: dict):
 
         while objective.unplaced_candidates:
             offspring = objective.unplaced_candidates.pop()
-            attempt = 0
+            attempt = 1
             mutation_attempt_limit = self.params['mutation_attempt_limit']
             offspring_attempt_limit = self.params['offspring_attempt_limit']
             temp_offspring = deepcopy(offspring)
             replaced = 0
             while True:
+                temp_offspring = self.suboperators['chromosome_mutation'].apply(objective=temp_offspring,
+                                                                                arguments=subop_args['chromosome_mutation'])
                 self.suboperators['right_part_selector'].apply(objective=temp_offspring,
                                                                arguments=subop_args['right_part_selector'])
                 temp_offspring.reset_state()
@@ -382,7 +384,7 @@ def apply(self, objective: ParetoLevels, arguments: dict):
                     self.suboperators['pareto_level_updater'].apply(objective=(temp_offspring, objective),
                                                                     arguments=subop_args['pareto_level_updater'])
                     objective.history.add(system)
-                    # print(temp_offspring.obj_fun)
+                    print(temp_offspring.obj_fun)
                     break
                 elif replaced == offspring_attempt_limit:
                     print("Could not generate unique offspring")
@@ -391,8 +393,6 @@ def apply(self, objective: ParetoLevels, arguments: dict):
                     temp_offspring = deepcopy(offspring)
                     replaced += 1
                     attempt = 0
-                temp_offspring = self.suboperators['chromosome_mutation'].apply(objective=temp_offspring,
-                                                                                arguments=subop_args['chromosome_mutation'])
                 attempt += 1
         return objective
 
@@ -447,7 +447,7 @@ def apply(self, objective : ParetoLevels, arguments : dict):
                 self.suboperators['chromosome_fitness'].apply(objective=candidate,
                                                               arguments=subop_args['chromosome_fitness'])
                 objective.history.add(system)
-                # print(candidate.obj_fun)
+                print(candidate.obj_fun)
             objective.initial_placing()
 
             # TODO: consider carefully, where normalizer init shall be held. If here, only the initial values are employed

epde/operators/multiobjective/mutations.py

@@ -57,13 +57,22 @@ class EquationMutation(CompoundOperator):
     def apply(self, objective : Equation, arguments : dict):
         self_args, subop_args = self.parse_suboperator_args(arguments = arguments)
 
-        for term_idx in range(objective.n_immutable, len(objective.structure)):
-            if np.random.uniform(0, 1) <= self.params['r_mutation']:
-                objective.structure[term_idx] = self.suboperators['mutation'].apply(objective = (term_idx, objective),
-                                                                                    arguments = subop_args['mutation'])
-        # term_idx = np.random.choice(len(objective.structure))
-        # objective.structure[term_idx] = self.suboperators['mutation'].apply(objective=(term_idx, objective),
-        #                                                                     arguments=subop_args['mutation'])
+        # for term_idx in range(objective.n_immutable, len(objective.structure)):
+        #     if np.random.uniform(0, 1) <= self.params['r_mutation']:
+        #         objective.structure[term_idx] = self.suboperators['mutation'].apply(objective = (term_idx, objective),
+        #                                                                             arguments = subop_args['mutation'])
+        nonzero_terms_mask = np.array([False if weight == 0 else True for weight in objective.weights_internal],
+                                      dtype=np.integer)
+        nonrs_terms_idx = [i for i, term in enumerate(objective.structure) if i != objective.target_idx]
+        nonzero_terms_idx = [item for item, keep in zip(nonrs_terms_idx, nonzero_terms_mask) if keep]
+        nonzero_terms_idx.append(objective.target_idx)
+        # term_idx = np.random.choice(nonzero_terms_idx)
+        if len(nonzero_terms_idx) > 0:
+            term_idx = np.random.choice(nonzero_terms_idx)
+        else:
+            term_idx = objective.target_idx
+        objective.structure[term_idx] = self.suboperators['mutation'].apply(objective=(term_idx, objective),
+                                                                            arguments=subop_args['mutation'])
         return objective
 
     def use_default_tags(self):

epde/operators/utils/parameters/default_parameters_multi_objective.json

@@ -10,8 +10,8 @@
         "number_of_neighbors" : 4
  	},
 	"ParetoLevelUpdater" : {
-    	"mutation_attempt_limit" : 5,
-		"offspring_attempt_limit" : 5
+    	"mutation_attempt_limit" : 3,
+		"offspring_attempt_limit" : 10
 	},
 	"InitialParetoLevelSorting" : {
 

projects/pic/data/ac/ac.py

@@ -132,10 +132,10 @@ def ac_discovery(foldername, noise_level):
     #                                     preprocessor_kwargs={'epochs_max' : 1e3})
     epde_search_obj.set_preprocessor(default_preprocessor_type='FD',
                                      preprocessor_kwargs={})
-    popsize = 8
+    popsize = 16
 
     epde_search_obj.set_moeadd_params(population_size=popsize,
-                                      training_epochs=10)
+                                      training_epochs=5)
 
     custom_grid_tokens = CacheStoredTokens(token_type='grid',
                                                 token_labels=['t', 'x'],
@@ -150,7 +150,7 @@ def ac_discovery(foldername, noise_level):
 
     factors_max_number = {'factors_num': [1, 2], 'probas': [0.65, 0.35]}
 
-    bounds = (1e-8, 1e-0)
+    bounds = (1e-12, 1e-0)
     epde_search_obj.fit(data=noised_data, variable_names=['u', ], max_deriv_order=(2, 3), derivs=None,
                         equation_terms_max_number=5, data_fun_pow=3,
                         additional_tokens=[],

projects/pic/data/kdv/kdv.py

@@ -430,7 +430,7 @@ def kdv_sindy_discovery(foldername, noise_level):
     popsize = 8
 
     epde_search_obj.set_moeadd_params(population_size=popsize,
-                                      training_epochs=5)
+                                      training_epochs=10)
 
     custom_trigonometric_eval_fun = {
         'cos(t)sin(x)': lambda *grids, **kwargs: (np.cos(grids[0]) * np.sin(grids[1])) ** kwargs['power']}
@@ -451,7 +451,7 @@ def kdv_sindy_discovery(foldername, noise_level):
 
     factors_max_number = {'factors_num': [1, 2], 'probas': [0.65, 0.35]}
 
-    bounds = (1e-8, 1e-0)
+    bounds = (1e-12, 1e-0)
     epde_search_obj.fit(data=noised_data, variable_names=['u', ], max_deriv_order=(2, 3), derivs=None,
                         equation_terms_max_number=5, data_fun_pow=3,
                         additional_tokens=[],

projects/pic/data/ode/ode.py

@@ -145,8 +145,8 @@ def ODE_discovery(foldername, noise_level):
     epde_search_obj.set_preprocessor(default_preprocessor_type='FD',
                                      preprocessor_kwargs={})
 
-    popsize = 12
-    epde_search_obj.set_moeadd_params(population_size=popsize, training_epochs=10)
+    popsize = 16
+    epde_search_obj.set_moeadd_params(population_size=popsize, training_epochs=15)
 
     factors_max_number = {'factors_num': [1, 2], 'probas': [0.65, 0.35]}
 

projects/pic/data/pde_compound/pde_compound.py

@@ -170,7 +170,7 @@ def run_test(self, noise_level: int = 0) -> bool:
         search_obj = EpdeSearch(
             use_solver=False,
             use_pic=True,
-            boundary=10,
+            boundary=(10, 25),
             coordinate_tensors=[grid[0], grid[1]],
             verbose_params={'show_iter_idx': True},
             device='cuda'
@@ -209,10 +209,10 @@ def run_discovery(self, noise_level: int = 0) -> EpdeSearch:
             preprocessor_kwargs={}
         )
 
-        popsize = 12
+        popsize = 8
         search_obj.set_moeadd_params(
             population_size=popsize,
-            training_epochs=2
+            training_epochs=5
         )
 
         # Prepare custom tokens
@@ -245,7 +245,7 @@ def run_discovery(self, noise_level: int = 0) -> EpdeSearch:
         )
 
         factors_max_number = {'factors_num': [1, 2], 'probas': [0.65, 0.35]}
-        bounds = (1e-4, 1e-3)
+        bounds = (1e-8, 1e-0)
 
         search_obj.fit(
             data=noised_data,

projects/pic/data/pde_divide/pde_divide.py

@@ -132,10 +132,10 @@ def run_discovery(self, noise_level: int = 0) -> EpdeSearch:
             preprocessor_kwargs={}
         )
 
-        popsize = 8
+        popsize = 16
         search_obj.set_moeadd_params(
             population_size=popsize,
-            training_epochs=50
+            training_epochs=20
         )
 
         grid_tokens, custom_trig_tokens = self.create_custom_tokens(grid)

projects/pic/data/vdp/vdp.py

@@ -140,7 +140,7 @@ def vdp_discovery(foldername, noise_level):
     epde_search_obj.set_preprocessor(default_preprocessor_type='FD',
                                      preprocessor_kwargs={})
 
-    popsize = 8
+    popsize = 16
     epde_search_obj.set_moeadd_params(population_size=popsize, training_epochs=10)
 
     factors_max_number = {'factors_num': [1, 2], 'probas': [0.65, 0.35]}