Make required properties abstract

popt/cost_functions/quadratic.py

@@ -4,15 +4,15 @@
 from popt.cost_functions.epf import epf
 
 
-def quadratic(state, *args, **kwargs):
+def quadratic(x, *args, **kwargs):
 	r"""Quadratic objective function
 
 	$$ f(x) = ||x - b||^2_A $$
 	"""
 
 	r = kwargs.get('r', -1)
 
-	x = state[0]['vector']
+	x = x[0]['vector']
 	dim, ne = x.shape
 	A = 0.5*np.diag(np.ones(dim))
 	b = 1.0*np.ones(dim)

popt/loop/ensemble_base.py

@@ -155,20 +155,21 @@ def function(self, x, *args, **kwargs):
         self._invert_scale_state()  # ensure that state is in [lb,ub]
         self._set_multilevel_state(self.state, x)  # set multilevel state if applicable
         run_success = self.calc_prediction(save_prediction=self.save_prediction)  # calculate flow data
-        self._set_multilevel_state(self.state, x) # For some reason this has to be done again after calc_prediction
-        self._scale_state()  # scale back to [0, 1]
+        self._set_multilevel_state(self.state, x)  # toggle back after calc_prediction
 
         # Evaluate the objective function
         if run_success:
             func_values = self.obj_func(
                 self.pred_data, 
                 input_dict=self.sim.input_dict,
-                true_order=self.sim.true_order, 
+                true_order=self.sim.true_order,
+                state=self.state, # pass state for possible use in objective function
                 **kwargs
             )
         else:
             func_values = np.inf  # the simulations have crashed
 
+        self._scale_state()  # scale back to [0, 1]
         if len(x.shape) == 1: self.state_func_values = func_values
         else: self.ens_func_values = func_values
 

popt/loop/optimize.py

@@ -4,6 +4,7 @@
 import logging
 import time
 import pickle
+from abc import ABC, abstractmethod
 
 # Internal imports
 import popt.misc_tools.optim_tools as ot
@@ -20,7 +21,7 @@
 logger.addHandler(console_handler)
 
 
-class Optimize:
+class Optimize(ABC):
     """
     Class for ensemble optimization algorithms. These are classified by calculating the sensitivity or gradient using
     ensemble instead of classical derivatives. The loop is else as a classic optimization loop: a state (or control
@@ -102,19 +103,36 @@ def __init__(self, **options):
         self.epf_iteration = 0
 
         # Initialize variables (set in subclasses)
-        # TODO: these variables should be abstract properties that subclasses are forced to define 
         self.options = None
-        self.mean_state = None
         self.obj_func_values = None
-        self.fun = None  # objective function
-        self.obj_func_tol = None  # objective tolerance limit
 
         # Initialize number of function and jacobi evaluations
         self.nfev = 0
         self.njev = 0
 
         self.msg = 'Convergence was met :)'
 
+    # Abstract function that subclasses are forced to define
+    @abstractmethod
+    def fun(self, x, *args, **kwargs):  # objective function
+        pass
+
+    # Abstract properties that subclasses are forced to define
+    @property
+    @abstractmethod
+    def xk(self):  # current state
+        pass
+
+    @property
+    @abstractmethod
+    def ftol(self):  # function tolerance
+        pass
+
+    @ftol.setter
+    @abstractmethod
+    def ftol(self, value):  # setter for function tolerance
+        pass
+
     def run_loop(self):
         """
         This is the main optimization loop.
@@ -138,7 +156,7 @@ def run_loop(self):
         epf_not_converged = True
         previous_state = None
         if self.epf:
-            previous_state = self.mean_state
+            previous_state = self.xk
             logger.info(f'       -----> EPF-EnOpt: {self.epf_iteration}, {self.epf["r"]} (outer iteration, penalty factor)')  # print epf info
 
         while epf_not_converged:  # outer loop using epf
@@ -178,14 +196,14 @@ def run_loop(self):
                 if self.epf_iteration > self.epf['max_epf_iter']:  # max epf_iterations set to 10
                     logger.info(f'       -----> EPF-EnOpt: maximum epf iterations reached')  # print epf info
                     break
-                p = np.abs(previous_state-self.mean_state) / (np.abs(previous_state) + 1.0e-9)
+                p = np.abs(previous_state-self.xk) / (np.abs(previous_state) + 1.0e-9)
                 conv_crit = self.epf['conv_crit']
                 if np.any(p > conv_crit):
                     epf_not_converged = True
-                    previous_state = self.mean_state
+                    previous_state = self.xk
                     self.epf['r'] *= self.epf['r_factor']  # increase penalty factor
-                    self.obj_func_tol *= self.epf['tol_factor']  # decrease tolerance
-                    self.obj_func_values = self.fun(self.mean_state, epf = self.epf)
+                    self.ftol *= self.epf['tol_factor']  # decrease tolerance
+                    self.obj_func_values = self.fun(self.xk, epf = self.epf)
                     self.iteration = 0
                     self.epf_iteration += 1
                     optimize_result = ot.get_optimize_result(self)
@@ -196,9 +214,10 @@ def run_loop(self):
                     logger.info(f'       -----> EPF-EnOpt: {self.epf_iteration}, {r} (outer iteration, penalty factor)')  # print epf info
                 else:
                     logger.info(f'       -----> EPF-EnOpt: converged, no variables changed more than {conv_crit*100} %')  # print epf info
-                    final_obj_no_penalty = str(round(float(np.mean(self.fun(self.mean_state))),4))
+                    final_obj_no_penalty = str(round(float(np.mean(self.fun(self.xk))),4))
                     logger.info(f'       -----> EPF-EnOpt: objective value without penalty = {final_obj_no_penalty}') # print epf info
 
+
     def save(self):
         """
         We use pickle to dump all the information we have in 'self'. Can be used, e.g., if some error has occurred.
@@ -218,6 +237,7 @@ def load(self):
         # Save in 'self'
         self.__dict__.update(tmp_load)
 
+    @abstractmethod
     def calc_update(self):
         """
         This is an empty dummy function. Actual functionality must be defined by the subclasses.

popt/misc_tools/optim_tools.py

@@ -334,7 +334,7 @@ def get_optimize_result(obj):
     """
 
     # Initialize dictionary of variables to save
-    save_dict = OptimizeResult({'success': True, 'x': obj.mean_state, 'fun': np.mean(obj.obj_func_values),
+    save_dict = OptimizeResult({'success': True, 'x': obj.xk, 'fun': np.mean(obj.fk),
                                 'nit':  obj.iteration, 'nfev': obj.nfev, 'njev': obj.njev})
     if hasattr(obj, 'epf') and obj.epf:
         save_dict['epf_iteration'] = obj.epf_iteration
@@ -349,7 +349,7 @@ def get_optimize_result(obj):
 
         # Loop over variables to store in save list
         for save_typ in savedata:
-            if 'mean_state' in save_typ:
+            if 'xk' in save_typ:
                 continue  # mean_state is alwaysed saved as 'x'
             if save_typ in locals():
                 save_dict[save_typ] = eval('{}'.format(save_typ))

popt/update_schemes/enopt.py

@@ -88,7 +88,7 @@ def __set__variable(var_name=None, defalut=None):
 
         # Set input as class variables
         self.options = options  # options
-        self.fun = fun  # objective function
+        self._fun = fun  # objective function
         self.cov = args[0]  # initial covariance
         self.jac = jac  # gradient function
         self.hess = hess  # hessian function
@@ -114,7 +114,7 @@ def __set__variable(var_name=None, defalut=None):
         # Calculate objective function of startpoint
         if not self.restart:
             self.start_time = time.perf_counter()
-            self.obj_func_values = self.fun(self.mean_state, epf=self.epf)
+            self.obj_func_values = self._fun(self.mean_state, epf=self.epf)
             self.nfev += 1
             self.optimize_result = ot.get_optimize_result(self)
             ot.save_optimize_results(self.optimize_result)
@@ -142,6 +142,25 @@ def __set__variable(var_name=None, defalut=None):
         # The EnOpt class self-ignites, and it is possible to send the EnOpt class as a callale method to scipy.minimize
         self.run_loop()  # run_loop resides in the Optimization class (super)
 
+    def fun(self, x, *args, **kwargs):
+        return self._fun(x, *args, **kwargs)
+
+    @property
+    def xk(self):
+        return self.mean_state
+
+    @property
+    def fk(self):
+        return self.obj_func_values
+
+    @property
+    def ftol(self):
+        return self.obj_func_tol
+
+    @ftol.setter
+    def ftol(self, value):
+        self.obj_func_tol = value
+
     def calc_update(self):
         """
         Update using steepest descent method with ensemble gradients
@@ -152,7 +171,7 @@ def calc_update(self):
         success = False
         resampling_iter = 0
 
-        while improvement is False:  # resampling loop
+        while not improvement:  # resampling loop
 
             # Shrink covariance each time we try resampling
             shrink = self.cov_factor ** resampling_iter
@@ -179,14 +198,14 @@ def calc_update(self):
             # Initialize for this step
             alpha_iter = 0
 
-            while improvement is False:  # backtracking loop
+            while not improvement:  # backtracking loop
 
                 new_state, new_step = self.optimizer.apply_update(self.mean_state, gradient,
                                                                   hessian=hessian, iter=self.iteration)
                 new_state = ot.clip_state(new_state, self.bounds)
 
                 # Calculate new objective function
-                new_func_values = self.fun(new_state, epf=self.epf)
+                new_func_values = self._fun(new_state, epf=self.epf)
                 self.nfev += 1
 
                 if np.mean(self.obj_func_values) - np.mean(new_func_values) > self.obj_func_tol:

popt/update_schemes/genopt.py

@@ -54,7 +54,7 @@ def __set__variable(var_name=None, defalut=None):
 
         # Set input as class variables
         self.options    = options    # options
-        self.fun        = fun        # objective function
+        self.function   = fun        # objective function
         self.jac        = jac        # gradient function
         self.jac_mut    = jac_mut    # mutation function
         self.corr_adapt = corr_adapt # correlation adaption function
@@ -82,7 +82,7 @@ def __set__variable(var_name=None, defalut=None):
         # Calculate objective function of startpoint
         if not self.restart:
             self.start_time = time.perf_counter()
-            self.obj_func_values = self.fun(self.mean_state)
+            self.obj_func_values = self.function(self.mean_state)
             self.nfev += 1
             self.optimize_result = ot.get_optimize_result(self)
             ot.save_optimize_results(self.optimize_result)
@@ -110,6 +110,25 @@ def __set__variable(var_name=None, defalut=None):
         # The GenOpt class self-ignites, and it is possible to send the EnOpt class as a callale method to scipy.minimize
         self.run_loop()  # run_loop resides in the Optimization class (super)
 
+    def fun(self, x, *args, **kwargs):
+        return self.function(x, *args, **kwargs)
+
+    @property
+    def xk(self):
+        return self.mean_state
+
+    @property
+    def fk(self):
+        return self.obj_func_values
+
+    @property
+    def ftol(self):
+        return self.obj_func_tol
+
+    @ftol.setter
+    def ftol(self, value):
+        self.obj_func_tol = value
+
     def calc_update(self):
         """
         Update using steepest descent method with ensemble gradients
@@ -148,7 +167,7 @@ def calc_update(self):
                 new_state = ot.clip_state(new_state, self.bounds)
 
                 # Calculate new objective function
-                new_func_values = self.fun(new_state)
+                new_func_values = self.function(new_state)
                 self.nfev += 1
 
                 if np.mean(self.obj_func_values) - np.mean(new_func_values) > self.obj_func_tol: