fixed f13

Author: jacobdenobel

include/sampling.hpp

@@ -17,17 +17,17 @@ namespace sampling
     struct Sampler
     {
         Sampler(const size_t d) : d(d) {}
-        [[nodiscard]] virtual Vector operator()() = 0;
         size_t d;
-
+        
         virtual void reset(const parameters::Modules &, const size_t)
         {
         }
-
+        
         virtual Float expected_length()
         {
             return std::sqrt(static_cast<Float>(d));
         }
+        [[nodiscard]] virtual Vector operator()() = 0;
     };
 
     /**
@@ -271,6 +271,9 @@ namespace sampling
 
         [[nodiscard]] inline Vector box_muller(const Vector& u)
         {
+            static Vector u_extra;
+            static int n_extra_used = 0;
+
             size_t n = u.size();
             size_t m = n / 2;
 
@@ -282,7 +285,14 @@ namespace sampling
             }
 
             if (n % 2 != 0)
-                z(n - 1) = box_muller(u(0), u(n - 1)).first;
+            {
+                if (u_extra.size() <= n_extra_used)
+                {
+                    u_extra = (*sampler)();
+                    n_extra_used = 0;
+                }
+                z(n - 1) = box_muller(u(n - 1), u_extra(n_extra_used++)).first;
+            }
 
             return z;
         }   

modcma/modularcmaes.py

@@ -481,10 +481,13 @@ def evaluate_bbob(
     # This speeds up the import, this import is quite slow, so import it lazy here
     # pylint: disable=import-outside-toplevel
     import ioh
+    from modcma.c_maes import Settings, parameters, options, constants, utils
+    from modcma.c_maes import ModularCMAES as cModularCMAES
 
     evals, fopts = np.array([]), np.array([])
     if seed:
         np.random.seed(seed)
+        utils.set_seed(seed)
     fitness_func = ioh.get_problem(
         fid, dimension=dim, instance=instance
     )
@@ -516,13 +519,11 @@ def evaluate_bbob(
         hs = []
 
         if cpp:
-            from modcma.c_maes import Settings, parameters, options
-            from modcma.c_maes import ModularCMAES as cModularCMAES
 
             modules = parameters.Modules()
             modules.matrix_adaptation = options.COVARIANCE
             modules.ssa = options.StepSizeAdaptation.CSA
-            modules.restart_strategy = options.RestartStrategy.STOP
+            modules.restart_strategy = options.RestartStrategy.NONE
 
             settings = Settings(
                 fitness_func.meta_data.n_variables, 
@@ -534,10 +535,9 @@ def evaluate_bbob(
                 sigma0=2.0,
                 target=fitness_func.optimum.y + 1e-8,
                 budget=fitness_func.meta_data.n_variables * 10_000,
-                cs=0.4149090010980616,
-                cmu=0.0512430870135861,
             )
             optimizer = cModularCMAES(settings)
+
             while not optimizer.break_conditions():
                 optimizer.step(fitness_func)
                 ps_norm.append(np.linalg.norm(optimizer.p.adaptation.ps))
@@ -549,6 +549,7 @@ def evaluate_bbob(
                 f_values.append(optimizer.p.pop.f.mean())
                 dm.append(optimizer.p.adaptation.dm.copy())
                 hs.append(optimizer.p.adaptation.hs)
+            # print(optimizer)
             title = "modcmacpp"
         else:
             optimizer = ModularCMAES(fitness_func, dim, x0 = np.zeros(dim), target=target, **kwargs)
@@ -602,7 +603,7 @@ def evaluate_bbob(
 
             for ax in ax0, ax1, ax2, ax3, ax4:
                 ax.grid()
-                ax.set_xlim(0, 300)
+                ax.set_xlim(0, 500)
                 ax.set_yscale("log", base=10)
             plt.show()
 

src/matrix_adaptation.cpp

@@ -117,15 +117,9 @@ namespace matrix_adaptation
 			}
 			return false;
 		}
-		
-		for (int i = 0; i < d.size(); ++i)
-    		d[i] = std::max(d[i], 1e-10);
-
-		d.noalias() = d.cwiseSqrt();
-		// inv_root_C.noalias() = eigen_solver.operatorInverseSqrt();
 
-		Matrix D_inv_sqrt = d.cwiseInverse().asDiagonal();
-		inv_root_C.noalias() = B * D_inv_sqrt * B.transpose();
+		d.noalias() = d.cwiseSqrt().eval();
+		inv_root_C.noalias() = B * d.cwiseInverse().asDiagonal() * B.transpose();
 		A.noalias() = B * d.asDiagonal();
 		return true;
 	}