Skip to content

On-off flashing ratchet

The on-off flashing ratchet, Brownian particles in a sawtooth potential switched on and off at random.

  • Name: flashing_ratchet
  • Field: physics/soft_matter
  • Formalism: stochastic differential equations
  • References: Bertin 2021; Ajdari and Prost 1992; Reimann 2002

Factory

@register_model(
    field='physics/soft_matter',
    tags=('brownian motion', 'sde', 'ratchet', 'nonequilibrium steady state', 'current'),
    references=('Bertin 2021', 'Ajdari and Prost 1992', 'Reimann 2002'),
    defaults={'until': 100.0, 'dt': 0.5, 'solver': EulerMaruyamaSolver(dt=0.002)},
)
def flashing_ratchet(
    n_particles: int = 500,
    *,
    barrier: float = 1.0,
    asymmetry: float = 0.1,
    L: float = 1.0,
    rate: float = 1.0,
    duration: float = 100.0,
    gamma: float = 1.0,
    kT: float = 0.05,
    random_state: int | np.random.Generator | None = 0,
) -> System:
    """Return the on-off flashing ratchet, Brownian particles in a sawtooth potential switched on and off at random.

    The potential switches at the rate `rate`, so on and off phases last `1 / rate` on average. While it is off the
    particles diffuse freely, and while it is on they slide to the minima of the sawtooth. A potential that rises over
    a short fraction `asymmetry` of its period sends more particles over its near maximum than over its far one, so
    the particles drift with a mean velocity of the sign of `0.5 - asymmetry`, and a symmetric potential gives no
    current. The switching times are drawn once, as a Poisson process over `duration` shared by every particle. The
    particles start at the minimum at the origin.

    Args:
        n_particles:
            The number of independent particles.

        barrier:
            The height of the sawtooth.

        asymmetry:
            The fraction of the period over which the potential rises, 0.5 for a symmetric potential.

        L:
            The period of the sawtooth.

        rate:
            The rate at which the potential switches.

        duration:
            The time span of the drawn switching times.

        gamma:
            The friction coefficient.

        kT:
            The thermal energy.

        random_state:
            The seed or generator of the switching times.

    Returns:
        system:
            The model, with the positions `particles.x`.

    Examples:
        >>> from skcomplex.models.physics.soft_matter import flashing_ratchet
        >>> len(flashing_ratchet(10, rate=2.0, duration=5.0).interactions[0].switches) > 0
        True
    """
    rng = np.random.default_rng(random_state)
    times = np.cumsum(rng.exponential(1 / rate, int(2 * rate * duration) + 10))
    particles = GroupElement(
        'particles',
        element=BrownianParticleElement,
        size=n_particles,
        values={'x': np.zeros((n_particles, 1)), 'gamma': gamma, 'temperature': kT},
    )
    sawtooth = FlashingSawtoothInteraction(
        'sawtooth',
        target='particles',
        barrier=barrier,
        L=L,
        asymmetry=asymmetry,
        switches=tuple(float(time) for time in times[times < duration]),
    )
    bath = ThermalNoiseInteraction('bath', target='particles')
    return System(elements=[particles], interactions=[sawtooth, bath], space=ContinuousSpace(n_dim=1))

Components

class BrownianParticleElement(Element):
    """Overdamped Brownian particle at the position `x`, with the friction `gamma` and the thermal energy `temperature`.

    The particle moves by `dx = F / gamma dt + sqrt(2 temperature / gamma) dW`, with the thermal noise of
    `ThermalNoiseInteraction` and the forces of other interactions, which divide by `gamma`.

    Args:
        label:
            The label of the particle or the group.

        visual:
            How the element is drawn, or `None` for `default_visual_`.

        x:
            The initial position.

        gamma:
            The friction coefficient.

        temperature:
            The thermal energy.

    Examples:
        >>> from skcomplex.models.physics.soft_matter import BrownianParticleElement
        >>> BrownianParticleElement('bead', temperature=2.0).temperature
        2.0
    """

    def __init__(
        self,
        label: str,
        visual: VisualSpec | None = POINTS,
        *,
        x: Annotated[npt.ArrayLike, StateSpec(vector=True, role='position')] = 0.0,
        gamma: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 1.0,
        temperature: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 1.0,
    ) -> None:
        super().__init__(label, visual=visual)
        self.x = x
        self.gamma = gamma
        self.temperature = temperature

class FlashingSawtoothInteraction(RateInteraction):
    """Sawtooth potential on overdamped particles on a line, turned on and off at given times.

    Its force is the rate of the position, divided by the friction `gamma` of the target.

    The potential of period `L` has a minimum at every multiple of `L`, rises linearly by `barrier` over the fraction
    `asymmetry` of the period, and falls back over the rest. It is on from time 0 to the first of the `switches`, off
    until the next one, and so on, and it keeps its state after the last one.

    Args:
        label:
            The label of the interaction.

        target:
            The path of the particle or the group.

        barrier:
            The height of the sawtooth.

        L:
            The period of the sawtooth.

        asymmetry:
            The fraction of the period over which the potential rises, 0.5 for a symmetric potential.

        switches:
            The increasing times at which the potential turns off and on in turn.

    Examples:
        >>> from skcomplex.models.physics.soft_matter import FlashingSawtoothInteraction
        >>> FlashingSawtoothInteraction('ratchet', target='beads', switches=(1.0, 2.5)).switches
        (1.0, 2.5)
    """

    def __init__(
        self,
        label: str,
        *,
        target: Annotated[str, RoleSpec()],
        barrier: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 1.0,
        L: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 1.0,
        asymmetry: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 0.1,
        switches: tuple[float, ...] = (),
    ) -> None:
        super().__init__(label, target=target, state='x', divisor='gamma')
        self.barrier = barrier
        self.L = L
        self.asymmetry = asymmetry
        self.switches = switches

    def _get_rate(self, namespace: SimpleNamespace) -> object:
        """Return the force of the sawtooth while it is on.

        Args:
            namespace:
                The namespace of the symbols the rate can use.

        Returns:
            rate:
                The force, before the division by the friction.
        """
        s = namespace
        times = (0.0, *self.switches)
        on = series(s.t, times, [float(index % 2 == 0) for index in range(len(times))], 'previous')
        rise = s.asymmetry * s.L
        return on * where(mod(s.target.x[0], s.L) < rise, -s.barrier / rise, s.barrier / (s.L - rise))

class ThermalNoiseInteraction(RateInteraction):
    """Thermal noise of overdamped particles, the random force `sqrt(2 gamma temperature) dW` divided by `gamma`.

    The friction `gamma` and the thermal energy `temperature` are parameters of the target, so every component of the
    position diffuses with the coefficient `temperature / gamma`.

    Args:
        label:
            The label of the interaction.

        target:
            The path of the particle or the group.

        state:
            The name of the state that holds the positions.

    Examples:
        >>> from skcomplex.models.physics.soft_matter import ThermalNoiseInteraction
        >>> ThermalNoiseInteraction('bath', target='beads').divisor
        'gamma'
    """

    def __init__(self, label: str, *, target: Annotated[str, RoleSpec()], state: str = 'x') -> None:
        super().__init__(label, target=target, state=state, divisor='gamma')

    def _get_noise(self, namespace: SimpleNamespace) -> object:
        """Return the amplitude `sqrt(2 gamma temperature)` of the random force.

        Args:
            namespace:
                The namespace of the symbols the amplitude can use.

        Returns:
            amplitude:
                The noise amplitude, before the division by the friction.
        """
        target = namespace.target
        return sqrt(2 * target.gamma * target.temperature)

Equations

particles.x[i, 0]: diff(x_particles_0(t), t) == where(t < 0.0, 1.0, where(t < 50.36077023991866, where(t < 27.228537619829638, where(t < 17.491739929907155, where(t < 2.271947867344106, where(t < 1.6995290054347743, where(t < 0.6799319039689096, 1.0, 0.0), where(t < 1.7193356680238296, 1.0, where(t < 1.7216049947050578, 0.0, 1.0))), where(t < 5.330772612495114, where(t < 3.901888302002491, 0.0, where(t < 4.575471254669723, 1.0, 0.0)), where(t < 8.147558591570839, 1.0, where(t < 14.205311672013412, 0.0, 1.0)))), where(t < 21.75295336114063, where(t < 19.762122344508885, where(t < 17.493027680241976, 0.0, 1.0), where(t < 19.834620029428777, 0.0, where(t < 20.904020331596094, 1.0, 0.0))), where(t < 25.563983195700718, where(t < 24.902862013879123, 1.0, where(t < 25.25687205607237, 0.0, 1.0)), where(t < 27.05620260882674, 0.0, where(t < 27.09320680761243, 1.0, 0.0))))), where(t < 37.505988124250315, where(t < 31.63284330367876, where(t < 29.033283944455277, where(t < 28.25873572648208, 1.0, 0.0), where(t < 30.77096245449593, 1.0, where(t < 31.184632540330632, 0.0, 1.0))), where(t < 35.533608540332665, where(t < 33.484935656659054, 0.0, where(t < 35.21014926678207, 1.0, 0.0)), where(t < 35.74436996121096, 1.0, where(t < 36.96766392097993, 0.0, 1.0)))), where(t < 45.189229059222264, where(t < 40.64156607213524, where(t < 38.36067321796768, 0.0, where(t < 39.40939542823748, 1.0, 0.0)), where(t < 43.569735385448304, 1.0, where(t < 43.72342414403644, 0.0, 1.0))), where(t < 48.81050898076727, where(t < 46.31346658829151, 0.0, where(t < 48.211517006621065, 1.0, 0.0)), where(t < 49.279193645087155, 1.0, where(t < 50.09630096156472, 0.0, 1.0)))))), where(t < 67.85221741644159, where(t < 61.23747007603508, where(t < 55.69094282069004, where(t < 51.212654784902504, where(t < 51.10659272892388, 0.0, 1.0), where(t < 52.859602178572835, 0.0, where(t < 55.138840346152136, 1.0, 0.0))), where(t < 59.240919195574776, where(t < 58.853673773400516, 1.0, where(t < 58.96137640484926, 0.0, 1.0)), where(t < 59.53816787860719, 0.0, where(t < 60.30172481553574, 1.0, 0.0)))), where(t < 65.27822853514715, where(t < 63.895958277490514, where(t < 61.77341288681851, 1.0, where(t < 61.89796692092076, 0.0, 1.0)), where(t < 64.39342219294886, 0.0, where(t < 64.50699210760939, 1.0, 0.0))), where(t < 66.24252117886115, where(t < 65.57472817772877, 1.0, where(t < 65.99065030284662, 0.0, 1.0)), where(t < 66.75768758786631, 0.0, where(t < 67.75539224721531, 1.0, 0.0))))), where(t < 82.55214583146152, where(t < 76.98401784230862, where(t < 71.91484960985242, where(t < 69.03474684710383, 1.0, 0.0), where(t < 72.77959298957964, 1.0, where(t < 76.47712067583673, 0.0, 1.0))), where(t < 81.13468146048186, where(t < 78.26986561981286, 0.0, where(t < 80.2886897475347, 1.0, 0.0)), where(t < 81.8787462217468, 1.0, where(t < 82.47946863539438, 0.0, 1.0)))), where(t < 91.97897448506934, where(t < 89.19698229371365, where(t < 87.00472991542622, 0.0, where(t < 88.67505015088898, 1.0, 0.0)), where(t < 90.95201468264467, 1.0, where(t < 91.96692511932557, 0.0, 1.0))), where(t < 95.26811756539416, where(t < 93.3492982191213, 0.0, where(t < 94.76297135582796, 1.0, 0.0)), where(t < 95.28145762643209, 1.0, where(t < 98.01406095409945, 0.0, 1.0))))))))*where(mod(x_particles_0(t), L_sawtooth) < L_sawtooth*asymmetry_sawtooth, -barrier_sawtooth/(L_sawtooth*asymmetry_sawtooth), -barrier_sawtooth/(L_sawtooth*asymmetry_sawtooth - L_sawtooth))/gamma_particles