4.3.2. pde.pdes.base module

Base classes

class PDEBase(*, noise: Union[int, float, complex, ndarray, Sequence[Union[int, float, complex, ndarray]], Sequence[Sequence[Any]]] = 0, rng: Optional[Generator] = None)[source]

Bases: object

base class for solving partial differential equations

Custom PDEs can be implemented by specifying their evolution rate. In the simple case of deterministic PDEs, the methods PDEBase.evolution_rate() and PDEBase._make_pde_rhs_numba() need to be overwritten for the numpy and numba backend, respectively.

Parameters
  • noise (float or ndarray) – Magnitude of the additive Gaussian white noise that is supported for all PDEs by default. If set to zero, a deterministic partial differential equation will be solved. Different noise magnitudes can be supplied for each field in coupled PDEs.

  • rng (Generator) – Random number generator (default: default_rng()). Note that this random number generator is only used for numpy function, while compiled numba code is unaffected.

Note

If more complicated noise structures are required, the methods PDEBase.noise_realization() and PDEBase._make_noise_realization_numba() need to be overwritten for the numpy and numba backend, respectively.

cache_rhs: bool = False

Flag indicating whether the right hand side of the equation should be cached. If True, the same implementation is used in subsequent calls to solve. Note that this might lead to wrong results if the parameters of the PDE were changed after the first call. This option is thus disabled by default and should be used with care.

Type

bool

check_implementation: bool = True

Flag determining whether (some) numba-compiled functions should be checked against their numpy counter-parts. This can help with implementing a correct compiled version for a PDE class.

Type

bool

check_rhs_consistency(state: FieldBase, t: float = 0, *, tol: float = 1e-07, rhs_numba: Optional[Callable] = None)[source]

check the numba compiled right hand side versus the numpy variant

Parameters
  • state (FieldBase) – The state for which the evolution rates should be compared

  • t (float) – The associated time point

  • tol (float) – Acceptance tolerance. The check passes if the evolution rates differ by less then this value

  • rhs_numba (callable) – The implementation of the numba variant that is to be checked. If omitted, an implementation is obtained by calling PDEBase._make_pde_rhs_numba_cached().

complex_valued: bool = False

Flag indicating whether the right hand side is a complex-valued PDE, which requires all involved variables to be of complex type

Type

bool

abstract evolution_rate(state: FieldBase, t: float = 0) FieldBase[source]
explicit_time_dependence: Optional[bool] = None

Flag indicating whether the right hand side of the PDE has an explicit time dependence.

Type

bool

property is_sde: bool

flag indicating whether this is a stochastic differential equation

The BasePDF class supports additive Gaussian white noise, whose magnitude is controlled by the noise property. In this case, is_sde is True if self.noise != 0.

make_modify_after_step(state: FieldBase) Callable[[ndarray], float][source]

returns a function that can be called to modify a state

This function is applied to the state after each integration step when an explicit stepper is used. The default behavior is to not change the state.

Parameters

state (FieldBase) – An example for the state from which the grid and other information can be extracted

Returns

Function that can be applied to a state to modify it and which returns a measure for the corrections applied to the state

make_pde_rhs(state: FieldBase, backend: str = 'auto') Callable[[ndarray, float], ndarray][source]

return a function for evaluating the right hand side of the PDE

Parameters
  • state (FieldBase) – An example for the state from which the grid and other information can be extracted.

  • backend (str) – Determines how the function is created. Accepted values are ‘numpy` and ‘numba’. Alternatively, ‘auto’ lets the code decide for the most optimal backend.

Returns

Function determining the right hand side of the PDE

Return type

callable

make_sde_rhs(state: FieldBase, backend: str = 'auto') Callable[[ndarray, float], Tuple[ndarray, ndarray]][source]

return a function for evaluating the right hand side of the SDE

Parameters
  • state (FieldBase) – An example for the state from which the grid and other information can be extracted

  • backend (str) – Determines how the function is created. Accepted values are ‘python` and ‘numba’. Alternatively, ‘auto’ lets the code decide for the most optimal backend.

Returns

Function determining the deterministic part of the right hand side of the PDE together with a noise realization.

noise_realization(state: FieldBase, t: float = 0, label: str = 'Noise realization') FieldBase[source]

returns a realization for the noise

Parameters
  • state (ScalarField) – The scalar field describing the concentration distribution

  • t (float) – The current time point

  • label (str) – The label for the returned field

Returns

Scalar field describing the evolution rate of the PDE

Return type

ScalarField

solve(state: FieldBase, t_range: TRangeType, dt: float = None, tracker: Optional[Union[Sequence[Union[TrackerBase, str]], TrackerBase, str]] = 'auto', method: Union[str, SolverBase] = 'auto', ret_info: bool = False, **kwargs) Union[FieldBase, Tuple[FieldBase, Dict[str, Any]]][source]

convenience method for solving the partial differential equation

The method constructs a suitable solver (SolverBase) and controller (Controller) to advance the state over the temporal range specified by t_range. To obtain full flexibility, it is advisable to construct these classes explicitly.

Parameters
  • state (FieldBase) – The initial state (which also defines the spatial grid)

  • t_range (float or tuple) – Sets the time range for which the PDE is solved. This should typically be a tuple of two numbers, (t_start, t_end), specifying the initial and final time of the simulation. If only a single value is given, it is interpreted as t_end and the time range is assumed to be (0, t_end).

  • dt (float) – Time step of the chosen stepping scheme. If None, a default value based on the stepper will be chosen. In particular, if method == ‘auto’, ScipySolver with an automatic, adaptive time step provided by scipy is used. This is a flexible choice, but can also result in unstable or slow simulations. If an adaptive stepper is used (supported by ScipySolver and ExplicitSolver), the value given here sets the initial time step.

  • tracker – Defines a tracker that process the state of the simulation at specified time intervals. A tracker is either an instance of TrackerBase or a string, which identifies a tracker. All possible identifiers can be obtained by calling get_named_trackers(). Multiple trackers can be specified as a list. The default value auto checks the state for consistency (tracker ‘consistency’) and displays a progress bar (tracker ‘progress’) when tqdm is installed. More general trackers are defined in trackers, where all options are explained in detail. In particular, the interval at which the tracker is evaluated can be chosen when creating a tracker object explicitly.

  • method (SolverBase or str) – Specifies the method for solving the differential equation. This can either be an instance of SolverBase or a descriptive name like ‘explicit’ or ‘scipy’. The valid names are given by pde.solvers.base.SolverBase.registered_solvers(). The default value ‘auto’ selects ScipySolver if dt is not specified and ExplicitSolver otherwise. Details of the solvers and additional features (like adaptive time steps) are explained in their documentation.

  • ret_info (bool) – Flag determining whether diagnostic information about the solver process should be returned.

  • **kwargs – Additional keyword arguments are forwarded to the solver class chosen with the method argument. In particular, ExplicitSolver supports several schemes and an adaptive stepper can be enabled using adaptive=True. Conversely, ScipySolver accepts the additional arguments of scipy.integrate.solve_ivp().

Returns

The state at the final time point. If ret_info == True, a tuple with the final state and a dictionary with additional information is returned.

Return type

FieldBase

expr_prod(factor: float, expression: str) str[source]

helper function for building an expression with an (optional) pre-factor

Parameters
  • factor (float) – The value of the prefactor

  • expression (str) – The remaining expression

Returns

The expression with the factor appended if necessary

Return type

str