Fully biased walk among quenched broad traps
Walkers that hop always to the right through quenched traps with broadly distributed trapping times.
- Name:
biased_trap_walk - Field:
physics/soft_matter - Formalism: lattice dynamics
- References: Bertin 2021; Bouchaud and Georges 1990
Factory
@register_model(
field='physics/soft_matter',
tags=('lattice', 'anomalous diffusion', 'disorder', 'traps', 'subdiffusion'),
references=('Bertin 2021', 'Bouchaud and Georges 1990'),
defaults={'until': 2000, 'dt': 10},
observables=_trap_observables,
)
def biased_trap_walk(
n_walkers: int = 100,
*,
n_sites: int = 300,
alpha: float = 0.5,
tau0: float = 1.0,
random_state: int | np.random.Generator | None = 0,
) -> System:
"""Return walkers that hop always to the right through quenched traps with broadly distributed trapping times.
Every walker moves along its own line of sites and stays `tau_i` steps on site `i` before it hops to site
`i + 1`. The trapping times are drawn once, as `ceil(tau0 u^(-1 / alpha))` for uniform numbers `u`, so their
distribution decays as `tau0^alpha / tau^(1 + alpha)`. For `alpha < 1` their mean diverges, the time to cross `N`
sites grows as `N^(1 / alpha)`, and the displacement grows subdiffusively as `x ~ (t / tau0)^alpha`, with the
mean `<x> = (t / tau0)^alpha sin(pi alpha) / (pi alpha)` at long times. The walkers start on the first site.
Args:
n_walkers:
The number of walkers, each in its own landscape.
n_sites:
The number of sites of every landscape.
alpha:
The exponent of the distribution of trapping times.
tau0:
The smallest trapping time.
random_state:
The seed or generator of the trapping times.
Returns:
system:
The model, with the lattice `traps` of shape `(n_walkers, n_sites)` and the mean displacement `position`.
Examples:
>>> from skcomplex.models.physics.soft_matter import biased_trap_walk
>>> biased_trap_walk(4, n_sites=10).elements[0].check().prepare().size_
40
"""
shape = (n_walkers, n_sites)
tau = np.ceil(tau0 * (1 - np.random.default_rng(random_state).random(shape)) ** (-1 / alpha))
start = np.zeros(shape)
start[:, 0] = 1.0
traps = LatticeElement(
'traps',
element=TrapSiteElement,
space=LatticeSpace(shape, boundary='fixed', fill=-1.0),
values={'occupied': start, 'visited': start, 'tau': tau},
visual=VISITED,
)
hops = LocalRuleInteraction(
'hops',
target='traps',
state=('occupied', 'visited', 'clock'),
rule=_hop,
neighborhood=ALONG_ROW,
reads=('occupied', 'visited', 'clock', 'tau'),
)
return System(elements=[traps], interactions=[hops])
Components
class TrapSiteElement(Element):
"""Site of a line of traps, which holds a walker for its trapping time `tau` before the walker hops to the right.
The walker on a site counts the steps it has spent there with `clock` and hops to the next site after `tau`
steps, and every site keeps a record of whether a walker has visited it. A walker hops only onto a site of positive
trapping time, so on a lattice that reads `-1` outside its edges a walker on the last site stays there.
Args:
label:
The label of the site or the lattice group.
visual:
How the element is drawn, or `None` for `default_visual_`.
occupied:
Whether a walker is on the site.
visited:
Whether a walker has been on the site.
clock:
The number of steps the walker has spent on the site.
tau:
The trapping time of the site.
Examples:
>>> from skcomplex.models.physics.soft_matter import TrapSiteElement
>>> TrapSiteElement('site', tau=3.0).tau
3.0
"""
def __init__(
self,
label: str,
visual: VisualSpec | None = None,
*,
occupied: Annotated[float, StateSpec(domain=FiniteEnumeratedSet([0, 1]))] = 0.0,
visited: Annotated[float, StateSpec(domain=FiniteEnumeratedSet([0, 1]))] = 0.0,
clock: Annotated[float, StateSpec(domain=NN)] = 0.0,
tau: Annotated[float, StateSpec(domain=RealSet.open(0, oo))] = 1.0,
) -> None:
super().__init__(label, visual=visual)
self.occupied = occupied
self.visited = visited
self.clock = clock
self.tau = tau
class LocalRuleInteraction(Interaction):
"""Local rule of a cellular automaton, which sets the next value of every cell from its ordered neighbourhood.
The rule is a table of the next value of every configuration of the neighbourhood, `TableSpec`, a rule of the sum
of the neighbourhood, `TotalisticSpec`, a rule of the cell and the sum of its neighbours, `OuterTotalisticSpec`, an
object whose `table()` returns a `TableSpec`, or a function. A function receives one array per neighbour of every
read state, in neighbourhood order, with one value per cell, the uniform array `random` when the rule is stochastic,
and the parameters of the bound group it names, and returns the next values of the states, one array per state.
Args:
label:
The label of the interaction, unique among the parts of its system.
target:
The path of the group on the cells of a lattice, the nodes of a network, or a Cayley graph.
state:
The name of the state the rule sets, or a tuple of the names of the states it sets together.
rule:
The table, the totalistic rule, or the function.
neighborhood:
The ordered neighbourhood of every cell, such as `VonNeumannNeighborhood(center=True)`.
reads:
The names of the states whose neighbourhoods a function reads, or `None` for the states it sets.
schedule:
`'synchronous'`, `'random-sequential'`, `'random-site'`, `'checkerboard'`, or `'colored'`.
stochastic:
Whether a function receives the uniform random array `random`.
Raises:
ValidationError:
The state or the read states are not names, or the neighbourhood is not a neighbourhood.
Examples:
>>> from skcomplex.core.interactions import LocalRuleInteraction
>>> from skcomplex.core.spaces import VonNeumannNeighborhood
>>> from skcomplex.specs import TotalisticSpec
>>> parity = TotalisticSpec({1: 1, 3: 1})
>>> around = VonNeumannNeighborhood(n_dim=1, center=True)
>>> rule = LocalRuleInteraction('parity', target='cells', state='cell', rule=parity, neighborhood=around)
>>> rule.check().schedule
'synchronous'
"""
def __init__(
self,
label: str,
*,
target: Annotated[str, RoleSpec()],
state: str | tuple[str, ...],
rule: object,
neighborhood: BaseNeighborhood,
reads: Sequence[str] | None = None,
schedule: str = 'synchronous',
stochastic: bool = False,
) -> None:
super().__init__(label)
self.target = target
self.state = state
self.rule = rule
self.neighborhood = neighborhood
self.reads = reads
self.schedule = schedule
self.stochastic = stochastic
def _get_names(self) -> tuple[str, ...]:
"""Return the names of the states the rule sets."""
return (self.state,) if isinstance(self.state, str) else tuple(self.state)
def _check(self) -> None:
"""Raise an error for states that are not names, or a neighbourhood that is not one.
Raises:
ValidationError:
The state or the read states are not names, or the neighbourhood is not a neighbourhood.
"""
super()._check()
names = [*self._get_names(), *(self.reads or ())]
if not names or not all(isinstance(name, str) and name.isidentifier() for name in names):
error_msg = f'The local rule {self.label!r} sets {self.state!r} and reads {self.reads!r}. Give state names.'
raise ValidationError(error_msg)
if not isinstance(self.neighborhood, BaseNeighborhood):
error_msg = f'The neighbourhood of the local rule {self.label!r} is {self.neighborhood!r}. Give one.'
raise ValidationError(error_msg)
def _get_laws(self, namespace: SimpleNamespace) -> list[BaseLaw]:
"""Return the local rule law of the states of the target.
Args:
namespace:
The namespace of the symbols the laws can use, with `namespace.target`.
Returns:
laws:
The local rule law.
"""
targets = tuple(getattr(self._get_member(namespace, 'target'), name) for name in self._get_names())
reads = (
None
if self.reads is None
else tuple(getattr(self._get_member(namespace, 'target'), name) for name in self.reads)
)
return [
LocalRuleLaw(
targets[0] if isinstance(self.state, str) else targets,
self.neighborhood,
self.rule,
reads=reads,
schedule=self.schedule,
stochastic=self.stochastic,
label=self.label,
),
]