From 77b45f095bf74fb35fa8b03aab1150f95f33683e Mon Sep 17 00:00:00 2001 From: Thomas Bouquet Date: Fri, 31 Oct 2025 14:51:00 +0100 Subject: [PATCH 1/3] abc symplectic integrator --- src/simulated_bifurcation/optimizer/__init__.py | 2 +- .../optimizer/integrator/__init__.py | 1 + .../optimizer/integrator/abc_symplectic_integrator.py | 10 ++++++++++ .../{ => integrator}/symplectic_integrator.py | 4 ++-- .../optimizer/simulated_bifurcation_optimizer.py | 3 +-- 5 files changed, 15 insertions(+), 5 deletions(-) create mode 100644 src/simulated_bifurcation/optimizer/integrator/__init__.py create mode 100644 src/simulated_bifurcation/optimizer/integrator/abc_symplectic_integrator.py rename src/simulated_bifurcation/optimizer/{ => integrator}/symplectic_integrator.py (96%) diff --git a/src/simulated_bifurcation/optimizer/__init__.py b/src/simulated_bifurcation/optimizer/__init__.py index d3e4be92..91bd7156 100644 --- a/src/simulated_bifurcation/optimizer/__init__.py +++ b/src/simulated_bifurcation/optimizer/__init__.py @@ -15,10 +15,10 @@ """ from .environment import get_env, reset_env, set_env +from .integrator import SymplecticIntegrator from .simulated_bifurcation_engine import SimulatedBifurcationEngine from .simulated_bifurcation_optimizer import ( ConvergenceWarning, SimulatedBifurcationOptimizer, ) from .stop_window import StopWindow -from .symplectic_integrator import SymplecticIntegrator diff --git a/src/simulated_bifurcation/optimizer/integrator/__init__.py b/src/simulated_bifurcation/optimizer/integrator/__init__.py new file mode 100644 index 00000000..4c06c1d2 --- /dev/null +++ b/src/simulated_bifurcation/optimizer/integrator/__init__.py @@ -0,0 +1 @@ +from .symplectic_integrator import SymplecticIntegrator diff --git a/src/simulated_bifurcation/optimizer/integrator/abc_symplectic_integrator.py b/src/simulated_bifurcation/optimizer/integrator/abc_symplectic_integrator.py new file mode 100644 index 00000000..bace9e15 --- /dev/null +++ b/src/simulated_bifurcation/optimizer/integrator/abc_symplectic_integrator.py @@ -0,0 +1,10 @@ +from abc import ABC + +import torch + +from ...core.tensor_bearer import TensorBearer + + +class ABCSymplecticIntegrator(ABC, TensorBearer): + def __init__(self, dtype: torch.dtype, device: torch.device): + super().__init__(dtype=dtype, device=device) diff --git a/src/simulated_bifurcation/optimizer/symplectic_integrator.py b/src/simulated_bifurcation/optimizer/integrator/symplectic_integrator.py similarity index 96% rename from src/simulated_bifurcation/optimizer/symplectic_integrator.py rename to src/simulated_bifurcation/optimizer/integrator/symplectic_integrator.py index 92433d22..0766d676 100644 --- a/src/simulated_bifurcation/optimizer/symplectic_integrator.py +++ b/src/simulated_bifurcation/optimizer/integrator/symplectic_integrator.py @@ -3,10 +3,10 @@ import torch from numpy import minimum -from ..core.tensor_bearer import TensorBearer +from .abc_symplectic_integrator import ABCSymplecticIntegrator -class SymplecticIntegrator(TensorBearer): +class SymplecticIntegrator(ABCSymplecticIntegrator): """ Simulates the evolution of spins' momentum and position following the Hamiltonian quantum mechanics equations that drive the Simulated Bifurcation (SB) algorithm. diff --git a/src/simulated_bifurcation/optimizer/simulated_bifurcation_optimizer.py b/src/simulated_bifurcation/optimizer/simulated_bifurcation_optimizer.py index 2cd24ceb..1ae1e27d 100644 --- a/src/simulated_bifurcation/optimizer/simulated_bifurcation_optimizer.py +++ b/src/simulated_bifurcation/optimizer/simulated_bifurcation_optimizer.py @@ -3,14 +3,13 @@ from typing import Optional, Union import torch -from numpy import minimum from tqdm.auto import tqdm from ..core.tensor_bearer import TensorBearer from .environment import ENVIRONMENT +from .integrator import SymplecticIntegrator from .simulated_bifurcation_engine import SimulatedBifurcationEngine from .stop_window import StopWindow -from .symplectic_integrator import SymplecticIntegrator class ConvergenceWarning(Warning): From 0a96bd6fc5c163e979ea84f15d7bc8b9917c1fe1 Mon Sep 17 00:00:00 2001 From: Thomas Bouquet Date: Sun, 2 Nov 2025 20:49:45 +0100 Subject: [PATCH 2/3] rename euler integrator --- docs/modules/optimizer/symplectic_integrator.md | 2 +- src/simulated_bifurcation/optimizer/__init__.py | 2 +- src/simulated_bifurcation/optimizer/integrator/__init__.py | 2 +- ...plectic_integrator.py => euler_symplectic_integrator.py} | 2 +- .../optimizer/simulated_bifurcation_optimizer.py | 4 ++-- tests/optimizer/test_symplectic_integrator.py | 6 +++--- 6 files changed, 9 insertions(+), 9 deletions(-) rename src/simulated_bifurcation/optimizer/integrator/{symplectic_integrator.py => euler_symplectic_integrator.py} (98%) diff --git a/docs/modules/optimizer/symplectic_integrator.md b/docs/modules/optimizer/symplectic_integrator.md index c466304a..2b3ee80d 100644 --- a/docs/modules/optimizer/symplectic_integrator.md +++ b/docs/modules/optimizer/symplectic_integrator.md @@ -1,6 +1,6 @@ # Symplectic Integrator ```{eval-rst} -.. autoclass:: simulated_bifurcation.optimizer.SymplecticIntegrator +.. autoclass:: simulated_bifurcation.optimizer.EulerSymplecticIntegrator :members: ``` diff --git a/src/simulated_bifurcation/optimizer/__init__.py b/src/simulated_bifurcation/optimizer/__init__.py index 91bd7156..5c9b20a9 100644 --- a/src/simulated_bifurcation/optimizer/__init__.py +++ b/src/simulated_bifurcation/optimizer/__init__.py @@ -15,7 +15,7 @@ """ from .environment import get_env, reset_env, set_env -from .integrator import SymplecticIntegrator +from .integrator import EulerSymplecticIntegrator from .simulated_bifurcation_engine import SimulatedBifurcationEngine from .simulated_bifurcation_optimizer import ( ConvergenceWarning, diff --git a/src/simulated_bifurcation/optimizer/integrator/__init__.py b/src/simulated_bifurcation/optimizer/integrator/__init__.py index 4c06c1d2..71347df6 100644 --- a/src/simulated_bifurcation/optimizer/integrator/__init__.py +++ b/src/simulated_bifurcation/optimizer/integrator/__init__.py @@ -1 +1 @@ -from .symplectic_integrator import SymplecticIntegrator +from .euler_symplectic_integrator import EulerSymplecticIntegrator diff --git a/src/simulated_bifurcation/optimizer/integrator/symplectic_integrator.py b/src/simulated_bifurcation/optimizer/integrator/euler_symplectic_integrator.py similarity index 98% rename from src/simulated_bifurcation/optimizer/integrator/symplectic_integrator.py rename to src/simulated_bifurcation/optimizer/integrator/euler_symplectic_integrator.py index 0766d676..74c84b71 100644 --- a/src/simulated_bifurcation/optimizer/integrator/symplectic_integrator.py +++ b/src/simulated_bifurcation/optimizer/integrator/euler_symplectic_integrator.py @@ -6,7 +6,7 @@ from .abc_symplectic_integrator import ABCSymplecticIntegrator -class SymplecticIntegrator(ABCSymplecticIntegrator): +class EulerSymplecticIntegrator(ABCSymplecticIntegrator): """ Simulates the evolution of spins' momentum and position following the Hamiltonian quantum mechanics equations that drive the Simulated Bifurcation (SB) algorithm. diff --git a/src/simulated_bifurcation/optimizer/simulated_bifurcation_optimizer.py b/src/simulated_bifurcation/optimizer/simulated_bifurcation_optimizer.py index 1ae1e27d..99008a9a 100644 --- a/src/simulated_bifurcation/optimizer/simulated_bifurcation_optimizer.py +++ b/src/simulated_bifurcation/optimizer/simulated_bifurcation_optimizer.py @@ -7,7 +7,7 @@ from ..core.tensor_bearer import TensorBearer from .environment import ENVIRONMENT -from .integrator import SymplecticIntegrator +from .integrator import EulerSymplecticIntegrator from .simulated_bifurcation_engine import SimulatedBifurcationEngine from .stop_window import StopWindow @@ -130,7 +130,7 @@ def __init_window(self, matrix: torch.Tensor, early_stopping: bool) -> None: ) def __init_symplectic_integrator(self, matrix: torch.Tensor) -> None: - self.symplectic_integrator = SymplecticIntegrator( + self.symplectic_integrator = EulerSymplecticIntegrator( self.agents, self.time_step, self.pressure_slope, diff --git a/tests/optimizer/test_symplectic_integrator.py b/tests/optimizer/test_symplectic_integrator.py index 53afba90..10eca04b 100644 --- a/tests/optimizer/test_symplectic_integrator.py +++ b/tests/optimizer/test_symplectic_integrator.py @@ -3,7 +3,7 @@ import pytest import torch -from src.simulated_bifurcation.optimizer import SymplecticIntegrator +from src.simulated_bifurcation.optimizer import EulerSymplecticIntegrator from ..test_utils import DEVICES, DTYPES @@ -17,8 +17,8 @@ def init_integrator( device: torch.device, activation_function: Callable[[torch.Tensor], torch.Tensor], heat: bool, -) -> SymplecticIntegrator: - symplectic_integrator = SymplecticIntegrator( +) -> EulerSymplecticIntegrator: + symplectic_integrator = EulerSymplecticIntegrator( 2, 0.1, 0.01, From 00c026c927cbf9fc585425a568be63d9def2c279 Mon Sep 17 00:00:00 2001 From: Thomas Bouquet Date: Fri, 16 Jan 2026 14:51:23 +0100 Subject: [PATCH 3/3] integrator interface --- requirements.txt | 2 +- .../optimizer/integrator/__init__.py | 1 + .../integrator/abc_symplectic_integrator.py | 69 +++++++++++++++++- .../integrator/euler_symplectic_integrator.py | 55 +++----------- .../stormer_verlet_symplectic_integrator.py | 73 +++++++++++++++++++ .../simulated_bifurcation_optimizer.py | 4 +- 6 files changed, 156 insertions(+), 48 deletions(-) create mode 100644 src/simulated_bifurcation/optimizer/integrator/stormer_verlet_symplectic_integrator.py diff --git a/requirements.txt b/requirements.txt index bd277241..84f5a320 100644 --- a/requirements.txt +++ b/requirements.txt @@ -1,4 +1,4 @@ numpy>=2.0.0 sympy==1.13.1 -torch>=2.2.0 +torch==2.8.0 tqdm==4.67.1 diff --git a/src/simulated_bifurcation/optimizer/integrator/__init__.py b/src/simulated_bifurcation/optimizer/integrator/__init__.py index 71347df6..dad05fe8 100644 --- a/src/simulated_bifurcation/optimizer/integrator/__init__.py +++ b/src/simulated_bifurcation/optimizer/integrator/__init__.py @@ -1 +1,2 @@ from .euler_symplectic_integrator import EulerSymplecticIntegrator +from .stormer_verlet_symplectic_integrator import StormerVerletSymplecticIntegrator diff --git a/src/simulated_bifurcation/optimizer/integrator/abc_symplectic_integrator.py b/src/simulated_bifurcation/optimizer/integrator/abc_symplectic_integrator.py index bace9e15..5e686085 100644 --- a/src/simulated_bifurcation/optimizer/integrator/abc_symplectic_integrator.py +++ b/src/simulated_bifurcation/optimizer/integrator/abc_symplectic_integrator.py @@ -1,10 +1,75 @@ -from abc import ABC +from abc import ABC, abstractmethod +from typing import Callable, Tuple import torch +from numpy import minimum from ...core.tensor_bearer import TensorBearer class ABCSymplecticIntegrator(ABC, TensorBearer): - def __init__(self, dtype: torch.dtype, device: torch.device): + """ + Simulates the evolution of spins' momentum and position following the Hamiltonian quantum mechanics equations that + drive the Simulated Bifurcation (SB) algorithm. + """ + + def __init__( + self, + n_oscillators: int, + time_step: float, + pressure_slope: float, + heat_coefficient: float, + activation_function: Callable[[torch.Tensor], torch.Tensor], + heat: bool, + quadratic_tensor: torch.Tensor, + dtype: torch.dtype, + device: torch.device, + ): super().__init__(dtype=dtype, device=device) + n_spins = quadratic_tensor.shape[0] + self.position = self.init_oscillator((n_spins, n_oscillators)) + self.momentum = self.init_oscillator((n_spins, n_oscillators)) + self.time_step = time_step + self.pressure_slope = pressure_slope + self.heat_coefficient = heat_coefficient + self.activation_function = activation_function + self.heat = heat + self.quadratic_tensor = self._cast_tensor(quadratic_tensor) + self.quadratic_scale_parameter = ( + 0.5 * (n_spins - 1) ** 0.5 / (torch.sqrt(torch.sum(quadratic_tensor**2))) + ) + self.step = 0 + + def init_oscillator(self, shape: Tuple[int, int]) -> torch.Tensor: + return 2.0 * torch.rand(size=shape, device=self.device, dtype=self.dtype) - 1.0 + + def simulate_inelastic_walls(self) -> None: + self.momentum[torch.abs(self.position) > 1.0] = 0.0 + torch.clip(self.position, -1.0, 1.0, out=self.position) + + def simulate_heating(self, momentum_copy: torch.Tensor) -> None: + torch.add( + self.momentum, + momentum_copy, + alpha=self.time_step * self.heat_coefficient, + out=self.momentum, + ) + + def get_current_pressure(self) -> float: + return minimum(self.time_step * self.step * self.pressure_slope, 1.0) + + def integration_step(self) -> None: + if self.heat: + momentum_copy = self.momentum.clone() + self.integrate() + self.simulate_inelastic_walls() + if self.heat: + self.simulate_heating(momentum_copy) + self.step += 1 + + def sample_spins(self) -> torch.Tensor: + return torch.where(self.position >= 0.0, 1.0, -1.0) + + @abstractmethod + def integrate(self): + raise NotImplementedError() diff --git a/src/simulated_bifurcation/optimizer/integrator/euler_symplectic_integrator.py b/src/simulated_bifurcation/optimizer/integrator/euler_symplectic_integrator.py index 74c84b71..cc721183 100644 --- a/src/simulated_bifurcation/optimizer/integrator/euler_symplectic_integrator.py +++ b/src/simulated_bifurcation/optimizer/integrator/euler_symplectic_integrator.py @@ -1,7 +1,6 @@ -from typing import Callable, Tuple +from typing import Callable import torch -from numpy import minimum from .abc_symplectic_integrator import ABCSymplecticIntegrator @@ -24,23 +23,17 @@ def __init__( dtype: torch.dtype, device: torch.device, ): - super().__init__(dtype=dtype, device=device) - n_spins = quadratic_tensor.shape[0] - self.position = self.init_oscillator((n_spins, n_oscillators)) - self.momentum = self.init_oscillator((n_spins, n_oscillators)) - self.time_step = time_step - self.pressure_slope = pressure_slope - self.heat_coefficient = heat_coefficient - self.activation_function = activation_function - self.heat = heat - self.quadratic_tensor = self._cast_tensor(quadratic_tensor) - self.quadratic_scale_parameter = ( - 0.5 * (n_spins - 1) ** 0.5 / (torch.sqrt(torch.sum(quadratic_tensor**2))) + super().__init__( + n_oscillators=n_oscillators, + time_step=time_step, + pressure_slope=pressure_slope, + heat_coefficient=heat_coefficient, + activation_function=activation_function, + heat=heat, + quadratic_tensor=quadratic_tensor, + dtype=dtype, + device=device, ) - self.step = 0 - - def init_oscillator(self, shape: Tuple[int, int]) -> torch.Tensor: - return 2.0 * torch.rand(size=shape, device=self.device, dtype=self.dtype) - 1.0 def position_update(self) -> None: torch.add( @@ -67,31 +60,7 @@ def quadratic_momentum_update(self) -> None: alpha=self.time_step * self.quadratic_scale_parameter, ) - def simulate_inelastic_walls(self) -> None: - self.momentum[torch.abs(self.position) > 1.0] = 0.0 - torch.clip(self.position, -1.0, 1.0, out=self.position) - - def simulate_heating(self, momentum_copy: torch.Tensor) -> None: - torch.add( - self.momentum, - momentum_copy, - alpha=self.time_step * self.heat_coefficient, - out=self.momentum, - ) - - def get_current_pressure(self) -> float: - return minimum(self.time_step * self.step * self.pressure_slope, 1.0) - - def integration_step(self) -> None: - if self.heat: - momentum_copy = self.momentum.clone() + def integrate(self): self.momentum_update() self.quadratic_momentum_update() self.position_update() - self.simulate_inelastic_walls() - if self.heat: - self.simulate_heating(momentum_copy) - self.step += 1 - - def sample_spins(self) -> torch.Tensor: - return torch.where(self.position >= 0.0, 1.0, -1.0) diff --git a/src/simulated_bifurcation/optimizer/integrator/stormer_verlet_symplectic_integrator.py b/src/simulated_bifurcation/optimizer/integrator/stormer_verlet_symplectic_integrator.py new file mode 100644 index 00000000..e91c0a24 --- /dev/null +++ b/src/simulated_bifurcation/optimizer/integrator/stormer_verlet_symplectic_integrator.py @@ -0,0 +1,73 @@ +from typing import Callable, Tuple + +import torch +from numpy import minimum + +from .abc_symplectic_integrator import ABCSymplecticIntegrator + + +class StormerVerletSymplecticIntegrator(ABCSymplecticIntegrator): + """ + Order-2 symplectic integrator based on the Störmer-Verlet integration method. + """ + + def __init__( + self, + n_oscillators: int, + time_step: float, + pressure_slope: float, + heat_coefficient: float, + activation_function: Callable[[torch.Tensor], torch.Tensor], + heat: bool, + quadratic_tensor: torch.Tensor, + dtype: torch.dtype, + device: torch.device, + ): + super().__init__( + n_oscillators=n_oscillators, + time_step=time_step, + pressure_slope=pressure_slope, + heat_coefficient=heat_coefficient, + activation_function=activation_function, + heat=heat, + quadratic_tensor=quadratic_tensor, + dtype=dtype, + device=device, + ) + + def integrate(self): + intermediate_position = self.position.clone() + intermediate_momentum = self.momentum.clone() + n = 4 + for _ in range(n): + auxiliary_momentum = torch.add( + intermediate_momentum, + intermediate_position, + alpha=self.time_step * (self.get_current_pressure() - 1.0) / (2.0 * n), + ) + auxiliary_momentum = torch.addmm( + auxiliary_momentum, + self.quadratic_tensor, + self.activation_function(intermediate_position), + alpha=self.time_step * self.quadratic_scale_parameter / (2.0 * n), + ) + torch.add( + intermediate_position, + intermediate_momentum, + alpha=self.time_step / n, + out=intermediate_position, + ) + torch.add( + auxiliary_momentum, + intermediate_position, + alpha=self.time_step * (self.get_current_pressure() - 1.0) / (2.0 * n), + out=intermediate_momentum, + ) + intermediate_momentum = torch.addmm( + intermediate_momentum, + self.quadratic_tensor, + self.activation_function(intermediate_position), + alpha=self.time_step * self.quadratic_scale_parameter / (2.0 * n), + ) + self.position = intermediate_position.clone() + self.momentum = intermediate_momentum.clone() diff --git a/src/simulated_bifurcation/optimizer/simulated_bifurcation_optimizer.py b/src/simulated_bifurcation/optimizer/simulated_bifurcation_optimizer.py index 99008a9a..a4a84a8b 100644 --- a/src/simulated_bifurcation/optimizer/simulated_bifurcation_optimizer.py +++ b/src/simulated_bifurcation/optimizer/simulated_bifurcation_optimizer.py @@ -7,7 +7,7 @@ from ..core.tensor_bearer import TensorBearer from .environment import ENVIRONMENT -from .integrator import EulerSymplecticIntegrator +from .integrator import EulerSymplecticIntegrator, StormerVerletSymplecticIntegrator from .simulated_bifurcation_engine import SimulatedBifurcationEngine from .stop_window import StopWindow @@ -130,7 +130,7 @@ def __init_window(self, matrix: torch.Tensor, early_stopping: bool) -> None: ) def __init_symplectic_integrator(self, matrix: torch.Tensor) -> None: - self.symplectic_integrator = EulerSymplecticIntegrator( + self.symplectic_integrator = StormerVerletSymplecticIntegrator( self.agents, self.time_step, self.pressure_slope,