Pre-solve

src/simulated_bifurcation/core/ising.py

@@ -23,7 +23,12 @@
 import torch
 from numpy import ndarray
 
-from ..optimizer import SimulatedBifurcationEngine, SimulatedBifurcationOptimizer
+from ..optimizer import (
+    Postprocessing,
+    Preprocessing,
+    SimulatedBifurcationEngine,
+    SimulatedBifurcationOptimizer,
+)
 
 # Workaround because `Self` type is only available in Python >= 3.11
 SelfIsing = TypeVar("SelfIsing", bound="Ising")
@@ -249,6 +254,7 @@ def minimize(
         sampling_period: int = 50,
         convergence_threshold: int = 50,
         timeout: Optional[float] = None,
+        presolve: bool = False
     ) -> None:
         """
         Minimize the energy of the Ising model using the Simulated Bifurcation
@@ -400,6 +406,30 @@ def minimize(
             sampling_period,
             convergence_threshold,
         )
+        if presolve:
+            presolved_spins, reduced_J, reduced_h = Preprocessing(
+                self.J, self.h
+            ).presolve()
+            if reduced_J.shape[1] == 0:
+                self.computed_spins = presolved_spins.repeat(agents, 1).t()
+                return
+            reduced_model = Ising(reduced_J, reduced_h, self.dtype, self.device)
+            reduced_model.minimize(
+                agents=agents,
+                max_steps=max_steps,
+                ballistic=ballistic,
+                heated=heated,
+                verbose=verbose,
+                use_window=use_window,
+                sampling_period=sampling_period,
+                convergence_threshold=convergence_threshold,
+                timeout=timeout,
+                presolve=False,
+            )
+            self.computed_spins = Postprocessing.reconstruct_spins(
+                reduced_model.computed_spins, presolved_spins
+            )
+            return
         tensor = self.as_simulated_bifurcation_tensor()
         spins = optimizer.run_integrator(tensor, use_window)
         if self.linear_term:

src/simulated_bifurcation/core/quadratic_polynomial.py

@@ -322,6 +322,7 @@ def optimize(
         sampling_period: int = 50,
         convergence_threshold: int = 50,
         timeout: Optional[float] = None,
+        presolve: bool = False,
     ) -> Tuple[torch.Tensor, torch.Tensor]:
         """
         Computes a local extremum of the model by optimizing
@@ -395,6 +396,10 @@ def optimize(
         minimize : bool, optional
             if `True` the optimization direction is minimization, otherwise it
             is maximization (default is True)
+        presolve : bool, optional
+            if `True`, the model will be preprocessed to find spins that can be
+            optimized before entering the SB algorithm, reducing the
+            computation size (default is False)
 
         Returns
         -------
@@ -414,6 +419,7 @@ def optimize(
             sampling_period=sampling_period,
             convergence_threshold=convergence_threshold,
             timeout=timeout,
+            presolve=presolve,
         )
         self.sb_result = self.convert_spins(ising_equivalent, domain)
         result = self.sb_result.t().to(dtype=self.dtype)
@@ -438,6 +444,7 @@ def minimize(
         sampling_period: int = 50,
         convergence_threshold: int = 50,
         timeout: Optional[float] = None,
+        presolve: bool = False,
     ) -> Tuple[torch.Tensor, torch.Tensor]:
         """
         Computes a local minimum of the model by optimizing
@@ -508,6 +515,10 @@ def minimize(
             if `True` only the best found solution to the optimization problem
             is returned, otherwise all the solutions found by the simulated
             bifurcation algorithm.
+        presolve : bool, optional
+            if `True`, the model will be preprocessed to find spins that can be
+            optimized before entering the SB algorithm, reducing the
+            computation size (default is False)
 
         Returns
         -------
@@ -526,6 +537,7 @@ def minimize(
             sampling_period=sampling_period,
             convergence_threshold=convergence_threshold,
             timeout=timeout,
+            presolve=presolve,
         )
 
     def maximize(
@@ -542,6 +554,7 @@ def maximize(
         sampling_period: int = 50,
         convergence_threshold: int = 50,
         timeout: Optional[float] = None,
+        presolve: bool = False,
     ) -> Tuple[torch.Tensor, torch.Tensor]:
         """
         Computes a local maximum of the model by optimizing
@@ -611,6 +624,10 @@ def maximize(
             if `True` only the best found solution to the optimization problem
             is returned, otherwise all the solutions found by the simulated
             bifurcation algorithm.
+        presolve : bool, optional
+            if `True`, the model will be preprocessed to find spins that can be
+            optimized before entering the SB algorithm, reducing the
+            computation size (default is False)
 
         Returns
         -------
@@ -629,6 +646,7 @@ def maximize(
             sampling_period=sampling_period,
             convergence_threshold=convergence_threshold,
             timeout=timeout,
+            presolve=presolve,
         )
 
     @staticmethod

src/simulated_bifurcation/models/abc_model.py

@@ -34,7 +34,8 @@ def optimize(
         use_window: bool = True,
         sampling_period: int = 50,
         convergence_threshold: int = 50,
-        timeout: Optional[float] = None
+        timeout: Optional[float] = None,
+        presolve: bool = False,
     ) -> Tuple[torch.Tensor, torch.Tensor]:
         return super().optimize(
             self.domain,
@@ -49,6 +50,7 @@ def optimize(
             sampling_period=sampling_period,
             convergence_threshold=convergence_threshold,
             timeout=timeout,
+            presolve=presolve,
         )
 
     def minimize(
@@ -63,7 +65,8 @@ def minimize(
         use_window: bool = True,
         sampling_period: int = 50,
         convergence_threshold: int = 50,
-        timeout: Optional[float] = None
+        timeout: Optional[float] = None,
+        presolve: bool = False,
     ) -> Tuple[torch.Tensor, torch.Tensor]:
         return self.optimize(
             agents,
@@ -77,6 +80,7 @@ def minimize(
             sampling_period=sampling_period,
             convergence_threshold=convergence_threshold,
             timeout=timeout,
+            presolve=presolve,
         )
 
     def maximize(
@@ -91,7 +95,8 @@ def maximize(
         use_window: bool = True,
         sampling_period: int = 50,
         convergence_threshold: int = 50,
-        timeout: Optional[float] = None
+        timeout: Optional[float] = None,
+        presolve: bool = False,
     ) -> Tuple[torch.Tensor, torch.Tensor]:
         return self.optimize(
             agents,
@@ -105,4 +110,5 @@ def maximize(
             sampling_period=sampling_period,
             convergence_threshold=convergence_threshold,
             timeout=timeout,
+            presolve=presolve,
         )

src/simulated_bifurcation/optimizer/__init__.py

@@ -15,6 +15,8 @@
 """
 
 from .environment import get_env, reset_env, set_env
+from .postprocessing import Postprocessing
+from .preprocessing import Preprocessing
 from .simulated_bifurcation_engine import SimulatedBifurcationEngine
 from .simulated_bifurcation_optimizer import (
     ConvergenceWarning,

src/simulated_bifurcation/optimizer/postprocessing.py

@@ -0,0 +1,46 @@
+import torch
+
+
+class Postprocessing:
+    """
+    Utilility class to reconstrcut the spin agents from presolved spins
+    on the original Ising model and spins optimized by the Simulated
+    Bifurcation algorithm on the reduced Ising model.
+    """
+
+    @staticmethod
+    def reconstruct_spins(
+        optimized_spins: torch.Tensor, pre_solved_spins: torch.Tensor
+    ) -> torch.Tensor:
+        """
+        Reconstruct the spin vectors by merging the presolved spins and
+        the spins optimized by the Simulated Bifurcation algorithm on the
+        reduced Ising model.
+
+        Parameters
+        ----------
+        optimized_spins : torch.Tensor
+            Spins returned by the Simulated Bifurcation algorithm
+            that correspond to the presolved reduced Ising model.
+        pre_solved_spins : torch.Tensor
+            Presolved spins of the original Ising model.
+
+        Returns
+        -------
+        torch.Tensor
+            Reconstructed spins taking in account both the SB-optimized
+            spins and the presolved spins.
+        """
+        original_dimension = pre_solved_spins.shape[0]
+        agents = optimized_spins.shape[1]
+        reconstructed_spins = torch.zeros(original_dimension, agents, dtype=torch.int32)
+        reconstructed_spins[pre_solved_spins == 0] = optimized_spins.clone().to(
+            dtype=torch.int32
+        )
+        reconstructed_spins[torch.abs(pre_solved_spins) == 1] = (
+            pre_solved_spins[torch.abs(pre_solved_spins) == 1]
+            .repeat(agents, 1)
+            .t()
+            .to(dtype=torch.int32)
+        )
+        return reconstructed_spins