4.4. pde.solvers package¶
Solvers define how a PDE is solved, i.e., how the initial state is advanced in time.
class controlling a simulation |
|
class for solving partial differential equations explicitly |
|
class for solving partial differential equations explicitly using MPI |
|
class for solving partial differential equations implicitly |
|
class for solving partial differential equations using scipy |
|
returns all solvers that are currently registered |
- class Controller(solver, t_range, tracker='auto')[source]¶
Bases:
object
class controlling a simulation
The controller calls a solver to advance the simulation into the future and it takes care of trackers that analyze and modify the state periodically.
- Parameters
solver (
SolverBase
) – Solver instance that is used to advance the simulation in timet_range (float or tuple) – Sets the time range for which the simulation is run. If only a single value t_end is given, the time range is assumed to be [0, t_end].
tracker (Optional[Union[Sequence[Union[TrackerBase, str]], TrackerBase, str]]) – Defines a tracker that process the state of the simulation at specified times. A tracker is either an instance of
TrackerBase
or a string, which identifies a tracker. All possible identifiers can be obtained by callingget_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’) whentqdm
is installed. More general trackers are defined intrackers
, 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.
- get_current_time()¶
process_time() -> float
Process time for profiling: sum of the kernel and user-space CPU time.
- run(initial_state, dt=None)[source]¶
run the simulation
Diagnostic information about the solver procedure are available in the diagnostics property of the instance after this function has been called.
- Parameters
state – The initial state of the simulation. This state will be copied and thus not modified by the simulation. Instead, the final state will be returned and trackers can be used to record intermediate states.
dt (float) – Time step of the chosen stepping scheme. If None, a default value based on the stepper will be chosen.
initial_state (TState) –
- Returns
The state at the final time point. If multiprocessing is used, only the main node will return the state. All other nodes return None.
- Return type
Optional[TState]
- class ExplicitSolver(pde, scheme='euler', *, backend='auto', adaptive=False, tolerance=0.0001)[source]¶
Bases:
SolverBase
class for solving partial differential equations explicitly
- Parameters
pde (
PDEBase
) – The instance describing the pde that needs to be solvedscheme (str) – Defines the explicit scheme to use. Supported values are ‘euler’ and ‘runge-kutta’ (or ‘rk’ for short).
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.
adaptive (bool) – When enabled, the time step is adjusted during the simulation using the error tolerance set with tolerance.
tolerance (float) – The error tolerance used in adaptive time stepping. This is used in adaptive time stepping to choose a time step which is small enough so the truncation error of a single step is below tolerance.
- make_stepper(state, dt=None)[source]¶
return a stepper function using an explicit scheme
- Parameters
- Returns
Function that can be called to advance the state from time t_start to time t_end. The function call signature is (state: numpy.ndarray, t_start: float, t_end: float)
- Return type
- name = 'explicit'¶
- class ImplicitSolver(pde, maxiter=100, maxerror=0.0001, backend='auto')[source]¶
Bases:
SolverBase
class for solving partial differential equations implicitly
- Parameters
pde (
PDEBase
) – The instance describing the pde that needs to be solvedmaxiter (int) – The maximal number of iterations per step
maxerror (float) – The maximal error that is permitted in each step
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.
- make_stepper(state, dt=None)[source]¶
return a stepper function using an implicit scheme
- Parameters
state (
FieldBase
) – An example for the state from which the grid and other information can be extracteddt (float) – Time step of the explicit stepping. If None, this solver specifies 1e-3 as a default value.
- Returns
Function that can be called to advance the state from time t_start to time t_end. The function call signature is (state: numpy.ndarray, t_start: float, t_end: float)
- Return type
- name = 'implicit'¶
- class ScipySolver(pde, backend='auto', **kwargs)[source]¶
Bases:
SolverBase
class for solving partial differential equations using scipy
This class is a thin wrapper around
scipy.integrate.solve_ivp()
. In particular, it supports all the methods implemented by this function.- Parameters
pde (
PDEBase
) – The instance describing the pde that needs to be solvedbackend (str) – Determines how the function is created. Accepted values are ‘numpy` and ‘numba’. Alternatively, ‘auto’ lets the code decide for the most optimal backend.
**kwargs – All extra arguments are forwarded to
scipy.integrate.solve_ivp()
.
- make_stepper(state, dt=None)[source]¶
return a stepper function
- Parameters
state (
FieldBase
) – An example for the state from which the grid and other information can be extracted.dt (float) – Initial time step for the simulation. If None, the solver will choose a suitable initial value.
- Returns
Function that can be called to advance the state from time t_start to time t_end.
- Return type
- name = 'scipy'¶