Integration modules
The modules used to numerically integrate the model forward in time. Users may add their own recipes here.
Integrate module
Module with the function to integrate the ordinary differential equations
of the model and its linearized version.
Description of the module functions
Two main functions:
- qgs.integrators.integrate.integrate_runge_kutta(f, t0, t, dt, ic=None, forward=True, write_steps=1, b=None, c=None, a=None)[source]
Integrate the ordinary differential equations (ODEs)
\[\dot{\boldsymbol{x}} = \boldsymbol{f}(t, \boldsymbol{x})\]with a specified Runge-Kutta method. The function \(\boldsymbol{f}\) should be a Numba jitted function. This function must have a signature
f(t, x)
wherex
is the state value andt
is the time.- Parameters
f (callable) – The Numba-jitted function \(\boldsymbol{f}\). Should have the signature``f(t, x)`` where
x
is the state value andt
is the time.t0 (float) – Initial time of the time integration. Corresponds to the initial condition. Important if the ODEs are non-autonomous.
t (float) – Final time of the time integration. Corresponds to the final condition. Important if the ODEs are non-autonomous.
dt (float) – Timestep of the integration.
ic (None or ndarray(float), optional) –
Initial (or final) conditions of the system. Can be a 1D or a 2D array:
1D: Provide a single initial condition. Should be of shape (n_dim,) where n_dim = \(\mathrm{dim}(\boldsymbol{x})\).
2D: Provide an ensemble of initial condition. Should be of shape (n_traj, n_dim) where n_dim = \(\mathrm{dim}(\boldsymbol{x})\), and where n_traj is the number of initial conditions.
If None, use a zero initial condition. Default to None. If the forward argument is False, it specifies final conditions.
forward (bool, optional) – Whether to integrate the ODEs forward or backward in time. In case of backward integration, the initial condition ic becomes a final condition. Default to forward integration.
write_steps (int, optional) – Save the state of the integration in memory every write_steps steps. The other intermediary steps are lost. It determines the size of the returned objects. Default is 1. Set to 0 to return only the final state.
b (None or ndarray, optional) – Vector of coefficients \(b_i\) of the Runge-Kutta method . If None, use the classic RK4 method coefficients. Default to None.
c (None or ndarray, optional) – Matrix of coefficients \(c_{i,j}\) of the Runge-Kutta method . If None, use the classic RK4 method coefficients. Default to None.
a (None or ndarray, optional) – Vector of coefficients \(a_i\) of the Runge-Kutta method . If None, use the classic RK4 method coefficients. Default to None.
- Returns
time, traj – The result of the integration:
time: Time at which the state of the system was saved. Array of shape (n_step,) where n_step is the number of saved states of the integration.
traj: Saved dynamical system states. 3D array of shape (n_traj, n_dim, n_steps). If n_traj = 1, a 2D array of shape (n_dim, n_steps) is returned instead.
- Return type
Examples
>>> from numba import njit >>> import numpy as np >>> from qgs.integrators.integrate import integrate_runge_kutta >>> a = 0.25 >>> F = 16. >>> G = 3. >>> b = 6. >>> # Lorenz 84 example >>> @njit ... def fL84(t, x): ... xx = -x[1] ** 2 - x[2] ** 2 - a * x[0] + a * F ... yy = x[0] * x[1] - b * x[0] * x[2] - x[1] + G ... zz = b * x[0] * x[1] + x[0] * x[2] - x[2] ... return np.array([xx, yy, zz]) >>> # no ic >>> # write_steps is 1 by default >>> tt, traj = integrate_runge_kutta(fL84, t0=0., t=10., dt=0.1) # 101 steps >>> print(traj.shape) (3, 101) >>> # 1 ic >>> ic = 0.1 * np.random.randn(3) >>> tt, traj = integrate_runge_kutta(fL84, t0=0., t=10., dt=0.1, ic=ic) # 101 steps >>> print(ic.shape) (3,) >>> print(traj.shape) (3, 101) >>> # 4 ic >>> ic = 0.1 * np.random.randn(4, 3) >>> tt, traj = integrate_runge_kutta(fL84, t0=0., t=10., dt=0.1, ic=ic) # 101 steps >>> print(ic.shape) (4, 3) >>> print(traj.shape) (4, 3, 101)
- qgs.integrators.integrate.integrate_runge_kutta_tgls(f, fjac, t0, t, dt, ic=None, tg_ic=None, forward=True, adjoint=False, inverse=False, boundary=None, write_steps=1, b=None, c=None, a=None)[source]
Integrate simultaneously the ordinary differential equations (ODEs)
\[\dot{\boldsymbol{x}} = \boldsymbol{f}(t, \boldsymbol{x})\]and its tangent linear model, i.e. the linearized ODEs
\[\dot{\boldsymbol{\delta x}} = \boldsymbol{\mathrm{J}}(t, \boldsymbol{x}) \cdot \boldsymbol{\delta x}\]where \(\boldsymbol{\mathrm{J}} = \frac{\partial \boldsymbol{f}}{\partial \boldsymbol{x}}\) is the Jacobian matrix of \(\boldsymbol{f}\), with a specified Runge-Kutta method. To solve this equation, one has to integrate simultaneously both ODEs.
The function \(\boldsymbol{f}\) and \(\boldsymbol{J}\) should be Numba jitted functions. These functions must have a signature
f(t, x)
andJ(t, x)
wherex
is the state value andt
is the time.- Parameters
f (callable) – The Numba-jitted function \(\boldsymbol{f}\). Should have the signature``f(t, x)`` where
x
is the state value andt
is the time.fjac (callable) – The Numba-jitted Jacobian \(\boldsymbol{J}\). Should have the signature``J(t, x)`` where
x
is the state value andt
is the time.t0 (float) – Initial time of the time integration. Corresponds to the initial conditions.
t (float) – Final time of the time integration. Corresponds to the final conditions.
dt (float) – Timestep of the integration.
ic (None or ndarray(float), optional) –
Initial (or final) conditions of the ODEs \(\dot{\boldsymbol{x}} = \boldsymbol{f}(t, \boldsymbol{x})\). Can be a 1D or a 2D array:
1D: Provide a single initial condition. Should be of shape (n_dim,) where n_dim = \(\mathrm{dim}(\boldsymbol{x})\).
2D: Provide an ensemble of initial condition. Should be of shape (n_traj, n_dim) where n_dim = \(\mathrm{dim}(\boldsymbol{x})\), and where n_traj is the number of initial conditions.
If None, use a zero initial condition. Default to None. If the forward argument is False, it specifies final conditions.
tg_ic (None or ndarray(float), optional) –
Initial (or final) conditions of the linear ODEs \(\dot{\boldsymbol{\delta x}} = \boldsymbol{\mathrm{J}}(t, \boldsymbol{x}) \cdot \boldsymbol{\delta x}\).
Can be a 1D, a 2D or a 3D array:
1D: Provide a single initial condition. This initial condition of the linear ODEs will be the same used for each initial condition ic of the ODEs \(\dot{\boldsymbol{x}} = \boldsymbol{f}(t, \boldsymbol{x})\) Should be of shape (n_dim,) where n_dim = \(\mathrm{dim}(\boldsymbol{x})\).
2D: Two sub-cases:
If tg_ic.shape[0]`=`ic.shape[0], assumes that each initial condition ic[i] of \(\dot{\boldsymbol{x}} = \boldsymbol{f}(t, \boldsymbol{x})\), correspond to a different initial condition tg_ic[i].
Else, assumes and integrate an ensemble of n_tg_traj initial condition of the linear ODEs for each initial condition of \(\dot{\boldsymbol{x}} = \boldsymbol{f}(t, \boldsymbol{x})\).
3D: An array of shape (n_traj, n_tg_traj, n_dim) which provide an ensemble of n_tg_ic initial conditions specific to each of the n_traj initial conditions of \(\dot{\boldsymbol{x}} = \boldsymbol{f}(t, \boldsymbol{x})\).
If None, use the identity matrix as initial condition, returning the fundamental matrix of solutions of the linear ODEs. Default to None. If the forward argument is False, it specifies final conditions.
forward (bool, optional) – Whether to integrate the ODEs forward or backward in time. In case of backward integration, the initial condition ic becomes a final condition. Default to forward integration.
adjoint (bool, optional) – Wheter to integrate the tangent \(\dot{\boldsymbol{\delta x}} = \boldsymbol{\mathrm{J}}(t, \boldsymbol{x}) \cdot \boldsymbol{\delta x}\) or the adjoint linear model \(\dot{\boldsymbol{\delta x}} = \boldsymbol{\mathrm{J}}^T(t, \boldsymbol{x}) \cdot \boldsymbol{\delta x}\). Integrate the tangent model by default.
inverse (bool, optional) – Wheter or not to invert the Jacobian matrix \(\boldsymbol{\mathrm{J}}(t, \boldsymbol{x}) \rightarrow \boldsymbol{\mathrm{J}}^{-1}(t, \boldsymbol{x})\). False by default.
boundary (None or callable, optional) – Allow to add a inhomogeneous term to linear ODEs: \(\dot{\boldsymbol{\delta x}} = \boldsymbol{\mathrm{J}}(t, \boldsymbol{x}) \cdot \boldsymbol{\delta x} + \Psi(t, \boldsymbol{x})\). The boundary \(\Psi\) should have the same signature as \(\boldsymbol{\mathrm{J}}\), i.e.
func(t, x)
. If None, don’t add anything (homogeneous case). None by default.write_steps (int, optional) – Save the state of the integration in memory every write_steps steps. The other intermediary steps are lost. It determines the size of the returned objects. Default is 1. Set to 0 to return only the final state.
b (None or ndarray, optional) – Vector of coefficients \(b_i\) of the Runge-Kutta method . If None, use the classic RK4 method coefficients. Default to None.
c (None or ndarray, optional) – Matrix of coefficients \(c_{i,j}\) of the Runge-Kutta method . If None, use the classic RK4 method coefficients. Default to None.
a (None or ndarray, optional) – Vector of coefficients \(a_i\) of the Runge-Kutta method . If None, use the classic RK4 method coefficients. Default to None.
- Returns
time, traj, tg_traj – The result of the integration:
time: Time at which the state of the system was saved. Array of shape (n_step,) where n_step is the number of saved states of the integration.
traj: Saved states of the ODEs. 3D array of shape (n_traj, n_dim, n_steps). If n_traj = 1, a 2D array of shape (n_dim, n_steps) is returned instead.
tg_traj: Saved states of the linear ODEs. Depending on the input initial conditions of both ODEs, it is at maximum a 4D array of shape (n_traj, n_tg_traj `n_dim, n_steps). If one of the dimension is 1, it is squeezed.
- Return type
Examples
>>> from numba import njit >>> import numpy as np >>> from qgs.integrators.integrate import integrate_runge_kutta_tgls >>> a = 0.25 >>> F = 16. >>> G = 3. >>> b = 6. >>> # Lorenz 84 example >>> @njit ... def fL84(t, x): ... xx = -x[1] ** 2 - x[2] ** 2 - a * x[0] + a * F ... yy = x[0] * x[1] - b * x[0] * x[2] - x[1] + G ... zz = b * x[0] * x[1] + x[0] * x[2] - x[2] ... return np.array([xx, yy, zz]) >>> @njit ... def DfL84(t, x): ... return np.array([[ -a , -2. * x[1], -2. * x[2]], ... [x[1] - b * x[2], -1. + x[0], -b * x[0]], ... [b * x[1] + x[2], b * x[0], -1. + x[0]]]) >>> # 4 ic, no tg_ic (fundamental matrix computation of an ensemble of ic) >>> ic = 0.1 * np.random.randn(4, 3) >>> tt, traj, dtraj = integrate_runge_kutta_tgls(fL84, DfL84, t0=0., t=10., dt=0.1, ... ic=ic, write_steps=1) # 101 steps >>> print(ic.shape) (4, 3) >>> print(traj.shape) (4, 3, 101) >>> print(dtraj.shape) (4, 3, 3, 101) >>> # 1 ic, 1 tg_ic (one ic and its tg_ic) >>> ic = 0.1 * np.random.randn(3) >>> tg_ic = 0.001 * np.random.randn(3) >>> tt, traj, dtraj = integrate_runge_kutta_tgls(fL84, DfL84, t0=0., t=10., dt=0.1, ... ic=ic, tg_ic=tg_ic) # 101 steps >>> print(ic.shape) (3,) >>> print(tg_ic.shape) (3,) >>> print(traj.shape) (3, 101) >>> print(dtraj.shape) (3, 101) >>> # 4 ic, 1 same tg_ic (an ensemble of ic with the same tg_ic) >>> ic = 0.1 * np.random.randn(4, 3) >>> tt, traj, dtraj = integrate_runge_kutta_tgls(fL84, DfL84, t0=0., t=10., dt=0.1, ... ic=ic, tg_ic=tg_ic) # 101 steps >>> print(ic.shape) (4, 3) >>> print(tg_ic.shape) (3,) >>> print(traj.shape) (4, 3, 101) >>> print(dtraj.shape) (4, 3, 101) >>> # 1 ic, 4 tg_ic (an ic with an ensemble of tg_ic in its tangent space) >>> ic = 0.1 * np.random.randn(3) >>> tg_ic = 0.001 * np.random.randn(4, 3) >>> tt, traj, dtraj = integrate_runge_kutta_tgls(fL84, DfL84, t0=0., t=10., dt=0.1, ... ic=ic, tg_ic=tg_ic) # 101 steps >>> print(ic.shape) (3,) >>> print(tg_ic.shape) (4, 3) >>> print(traj.shape) (3, 101) >>> print(dtraj.shape) (4, 3, 101) >>> # 2 ic, same 4 tg_ic (an ensemble of 2 ic, both with the same ensemble >>> # of tg_ic in their tangent space) >>> ic = 0.1 * np.random.randn(2, 3) >>> tg_ic = 0.001 * np.random.randn(4, 3) >>> tt, traj, dtraj = integrate_runge_kutta_tgls(fL84, DfL84, t0=0., t=10., dt=0.1, ... ic=ic, tg_ic=tg_ic) # 101 steps >>> print(ic.shape) (2, 3) >>> print(tg_ic.shape) (4, 3) >>> print(traj.shape) (2, 3, 101) >>> print(dtraj.shape) (2, 4, 3, 101) >>> # 2 ic, 4 different tg_ic (an ensemble of 2 ic, with different ensemble >>> # of tg_ic in their tangent space) >>> ic = 0.1 * np.random.randn(2, 3) >>> tg_ic = 0.001 * np.random.randn(2, 4, 3) >>> tt, traj, dtraj = integrate_runge_kutta_tgls(fL84, DfL84, t0=0., t=10., dt=0.1, ... ic=ic, tg_ic=tg_ic) # 101 steps >>> print(ic.shape) (2, 3) >>> print(tg_ic.shape) (2, 4, 3) >>> print(traj.shape) (2, 3, 101) >>> print(dtraj.shape) (2, 4, 3, 101)
Integrator module
Module with the classes of integrators to multi-thread the integration of an ordinary differential equations
of the model and its linearized version.
Module classes
- class qgs.integrators.integrator.RungeKuttaIntegrator(num_threads=None, b=None, c=None, a=None, number_of_dimensions=None)[source]
Bases:
object
Class to integrate the ordinary differential equations (ODEs)
\[\dot{\boldsymbol{x}} = \boldsymbol{f}(t, \boldsymbol{x})\]with a set of
TrajectoryProcess
and a specified Runge-Kutta method.- Parameters
num_threads (None or int, optional) – Number of
TrajectoryProcess
workers (threads) to use. If None, use the number of machine’s cores available. Default to None.b (None or ndarray, optional) – Vector of coefficients \(b_i\) of the Runge-Kutta method . If None, use the classic RK4 method coefficients. Default to None.
c (None or ndarray, optional) – Matrix of coefficients \(c_{i,j}\) of the Runge-Kutta method . If None, use the classic RK4 method coefficients. Default to None.
a (None or ndarray, optional) – Vector of coefficients \(a_i\) of the Runge-Kutta method . If None, use the classic RK4 method coefficients. Default to None.
number_of_dimensions (None or int, optional) – Allow to hardcode the dynamical system dimension. If None, evaluate the dimension from the callable
func
. Default to None.
- num_threads
Number of
TrajectoryProcess
workers (threads) to use.- Type
- b
Vector of coefficients \(b_i\) of the Runge-Kutta method .
- Type
- c
Matrix of coefficients \(c_{i,j}\) of the Runge-Kutta method .
- Type
- a
Vector of coefficients \(a_i\) of the Runge-Kutta method .
- Type
- n_traj
The number of trajectories (initial conditions) computed at the last integration performed by the integrator.
- Type
- func
Last function \(\boldsymbol{f}\) used by the integrator to integrate.
- Type
callable
- get_ic()[source]
Returns the initial conditions stored in the integrator.
- Returns
The initial conditions.
- Return type
- get_trajectories()[source]
Returns the result of the previous integrator integration.
- Returns
time, traj – The result of the integration:
- Return type
- initialize(convergence_time, dt, pert_size=0.01, reconvergence_time=None, forward=True, number_of_trajectories=1, ic=None, reconverge=False)[source]
Initialize the integration on an attractor by running it for a transient time, For an ensemble of initial conditions, can do the same transient time for each, or the convergence_time for the first one, and a smaller reconvergence_time for the subsequent ones. This results into initial conditions on the attractor, stored in
ic
.- Parameters
convergence_time (float) – Transient time needed to converge to the attractor.
dt (float) – Timestep of the transient integration.
pert_size (float, optional) – If the reconvergence is activated, size of the perturbation to add to the previous ic to find the next one. Default to 0.01 .
reconvergence_time (None or float, optional) – Transient time for the subsequent trajectories after the first long transient_time.
forward (bool, optional) – Whether to integrate the ODEs forward or backward in time. In case of backward integration, the initial condition ic becomes a final condition. Default to forward integration.
number_of_trajectories (int) – Number of initial conditions to find. Default to 1. Inactive if ic is provided.
ic (None or ndarray(float), optional) –
Initial (or final) conditions of the system. Can be a 1D or a 2D array:
1D: Provide a single initial condition. Should be of shape (n_dim,) where n_dim = \(\mathrm{dim}(\boldsymbol{x})\).
2D: Provide an ensemble of initial condition. Should be of shape (n_traj, n_dim) where n_dim = \(\mathrm{dim}(\boldsymbol{x})\), and where n_traj is the number of initial conditions.
If None, use number_trajectories random initial conditions. Default to None. If the forward argument is False, it specifies final conditions.
reconverge (bool) – Use or not the smaller transient time reconvergence with a perturbation after the first initial conditions have been computed. If activated, only use the
num_threads
first initial conditions of the ic arguments. Default to False.
- integrate(t0, t, dt, ic=None, forward=True, write_steps=1)[source]
Integrate the ordinary differential equations (ODEs)
\[\dot{\boldsymbol{x}} = \boldsymbol{f}(t, \boldsymbol{x})\]with a specified Runge-Kutta method and workers. The function \(\boldsymbol{f}\) is the Numba jitted function stored in
func
. The result of the integration can be obtained afterward by callingget_trajectories()
.- Parameters
t0 (float) – Initial time of the time integration. Corresponds to the initial condition’s ic time. Important if the ODEs are non-autonomous.
t (float) – Final time of the time integration. Corresponds to the final condition. Important if the ODEs are non-autonomous.
dt (float) – Timestep of the integration.
ic (None or ndarray(float), optional) –
Initial (or final) conditions of the system. Can be a 1D or a 2D array:
1D: Provide a single initial condition. Should be of shape (n_dim,) where n_dim = \(\mathrm{dim}(\boldsymbol{x})\).
2D: Provide an ensemble of initial condition. Should be of shape (n_traj, n_dim) where n_dim = \(\mathrm{dim}(\boldsymbol{x})\), and where n_traj is the number of initial conditions.
If None, use the initial conditions stored in
ic
. If thenic
is None, use a zero initial condition. Default to None. If the forward argument is False, it specifies final conditions.forward (bool, optional) – Whether to integrate the ODEs forward or backward in time. In case of backward integration, the initial condition ic becomes a final condition. Default to forward integration.
write_steps (int, optional) – Save the state of the integration in memory every write_steps steps. The other intermediary steps are lost. It determines the size of the returned objects. Default is 1. Set to 0 to return only the final state.
- set_bca(b=None, c=None, a=None, ic_init=True)[source]
Set the coefficients of the Runge-Kutta method and restart the integrator. s
- Parameters
b (None or ndarray, optional) – Vector of coefficients \(b_i\) of the Runge-Kutta method . If None, does not reinitialize these coefficients.
c (None or ndarray, optional) – Matrix of coefficients \(c_{i,j}\) of the Runge-Kutta method . If None, does not reinitialize these coefficients.
a (None or ndarray, optional) – Vector of coefficients \(a_i\) of the Runge-Kutta method . If None, does not reinitialize these coefficients.
ic_init (bool, optional) – Re-initialize or not the initial conditions of the integrator. Default to True.
- set_func(f, ic_init=True)[source]
Set the Numba-jitted function \(\boldsymbol{f}\) to integrate.
- Parameters
Warning
This function restarts the integrator!
- set_ic(ic)[source]
Direct setter for the integrator’s initial conditions
- Parameters
Initial condition of the system. Can be a 1D or a 2D array:
1D: Provide a single initial condition. Should be of shape (n_dim,) where n_dim = \(\mathrm{dim}(\boldsymbol{x})\).
2D: Provide an ensemble of initial condition. Should be of shape (n_traj, n_dim) where n_dim = \(\mathrm{dim}(\boldsymbol{x})\), and where n_traj is the number of initial conditions.
- class qgs.integrators.integrator.RungeKuttaTglsIntegrator(num_threads=None, b=None, c=None, a=None, number_of_dimensions=None)[source]
Bases:
object
Class to integrate simultaneously the ordinary differential equations (ODEs)
\[\dot{\boldsymbol{x}} = \boldsymbol{f}(t, \boldsymbol{x})\]and its tangent linear model, i.e. the linearized ODEs
\[\dot{\boldsymbol{\delta x}} = \boldsymbol{\mathrm{J}}(t, \boldsymbol{x}) \cdot \boldsymbol{\delta x}\]where \(\boldsymbol{\mathrm{J}} = \frac{\partial \boldsymbol{f}}{\partial \boldsymbol{x}}\) is the Jacobian matrix of \(\boldsymbol{f}\), with a specified Runge-Kutta method. To solve this equation, one has to integrate simultaneously both ODEs. This class does so with a set of
TglsTrajectoryProcess
workers.The function \(\boldsymbol{f}\) and \(\boldsymbol{J}\) should be Numba jitted functions. These functions must have a signature
f(t, x)
andJ(t, x)
wherex
is the state value andt
is the time.- Parameters
num_threads (None or int, optional) – Number of
TrajectoryProcess
workers (threads) to use. If None, use the number of machine’s cores available. Default to None.b (None or ndarray, optional) – Vector of coefficients \(b_i\) of the Runge-Kutta method . If None, use the classic RK4 method coefficients. Default to None.
c (None or ndarray, optional) – Matrix of coefficients \(c_{i,j}\) of the Runge-Kutta method . If None, use the classic RK4 method coefficients. Default to None.
a (None or ndarray, optional) – Vector of coefficients \(a_i\) of the Runge-Kutta method . If None, use the classic RK4 method coefficients. Default to None.
number_of_dimensions (None or int, optional) – Allow to hardcode the dynamical system dimension. If None, evaluate the dimension from the callable
func
. Default to None.
- num_threads
Number of
TrajectoryProcess
workers (threads) to use.- Type
- b
Vector of coefficients \(b_i\) of the Runge-Kutta method .
- Type
- c
Matrix of coefficients \(c_{i,j}\) of the Runge-Kutta method .
- Type
- a
Vector of coefficients \(a_i\) of the Runge-Kutta method .
- Type
- n_traj
The number of trajectories (initial conditions) of the non-linear ODEs computed at the last integration performed by the integrator.
- Type
- n_tg_traj
The number of trajectories (initial conditions) the linear ODEs computed at the last integration performed by the integrator.
- Type
- func
Last function \(\boldsymbol{f}\) used by the integrator to integrate.
- Type
callable
- func_jac
Last Jacobian matrix function \(\boldsymbol{J}\) used by the integrator to integrate.
- Type
callable
- get_ic()[source]
Returns the initial conditions of the non-linear ODEs stored in the integrator.
- Returns
The initial conditions.
- Return type
- get_tg_ic()[source]
Returns the initial conditions of the linear ODEs stored in the integrator.
- Returns
The initial conditions.
- Return type
- get_trajectories()[source]
Returns the result of the previous integrator integration.
- Returns
time, traj, tg_traj – The result of the integration:
time: time at which the state of the system was saved. Array of shape (
n_records
,).traj: Saved dynamical system states. 3D array of shape (
n_traj
,n_dim
,n_records
). Ifn_traj
= 1, a 2D array of shape (n_dim
,n_records
) is returned instead.tg_traj: Saved states of the linear ODEs. Depending on the input initial conditions of both ODEs, it is at maximum a 4D array of shape (
n_traj
,n_tg_traj
,n_dim
,n_records
). If one of the dimension is 1, it is squeezed.
- Return type
- initialize(convergence_time, dt, pert_size=0.01, reconvergence_time=None, forward=True, number_of_trajectories=1, ic=None, reconverge=None)[source]
Initialize the integration on an attractor by running it for a transient time, For an ensemble of initial conditions, can do the same transient time for each, or the convergence_time for the first one, and a smaller reconvergence_time for the subsequent ones. This results into initial conditions on the attractor, stored in
ic
.- Parameters
convergence_time (float) – Transient time needed to converge to the attractor.
dt (float) – Timestep of the transient integration.
pert_size (float, optional) – If the reconvergence is activated, size of the perturbation to add to the previous ic to find the next one. Default to 0.01 .
reconvergence_time (None or float, optional) – Transient time for the subsequent trajectories after the first long transient_time.
forward (bool, optional) – If true, integrate the ODEs forward in time, else, integrate backward in time. In case of backward integration, the initial condition ic becomes a final condition. Default to forward integration.
number_of_trajectories (int) – Number of initial conditions to find. Default to 1. Inactive if ic is provided.
ic (None or ndarray(float), optional) –
Initial condition of the system. Can be a 1D or a 2D array:
1D: Provide a single initial condition. Should be of shape (n_dim,) where n_dim = \(\mathrm{dim}(\boldsymbol{x})\).
2D: Provide an ensemble of initial condition. Should be of shape (n_traj, n_dim) where n_dim = \(\mathrm{dim}(\boldsymbol{x})\), and where n_traj is the number of initial conditions.
If None, use number_trajectories random initial conditions. Default to None.
reconverge (bool) – Use or not the smaller transient time reconvergence with a perturbation after the first initial conditions have been computed. If activated, only use the
num_threads
first initial conditions of the ic arguments. Default to False.
- integrate(t0, t, dt, ic=None, tg_ic=None, forward=True, adjoint=False, inverse=False, boundary=None, write_steps=1)[source]
Integrate simultaneously the non-linear and linearized ordinary differential equations (ODEs)
\[\dot{\boldsymbol{x}} = \boldsymbol{f}(t, \boldsymbol{x})\]and
\[\dot{\boldsymbol{\delta x}} = \boldsymbol{\mathrm{J}}(t, \boldsymbol{x}) \cdot \boldsymbol{\delta x}\]with a specified Runge-Kutta method and workers. The function \(\boldsymbol{f}\) is the Numba jitted function stored in
func
. The function \(\boldsymbol{J}\) is the Numba jitted function stored infunc_jac
. The result of the integration can be obtained afterward by callingget_trajectories()
.- Parameters
t0 (float) – Initial time of the time integration. Corresponds to the initial condition.
t (float) – Final time of the time integration. Corresponds to the final condition.
dt (float) – Timestep of the integration.
ic (None or ndarray(float), optional) –
Initial (or final) conditions of the system. Can be a 1D or a 2D array:
1D: Provide a single initial condition. Should be of shape (n_dim,) where n_dim = \(\mathrm{dim}(\boldsymbol{x})\).
2D: Provide an ensemble of initial condition. Should be of shape (n_traj, n_dim) where n_dim = \(\mathrm{dim}(\boldsymbol{x})\), and where n_traj is the number of initial conditions.
If None, use the initial conditions stored in
ic
. If thenic
is None, use a zero initial condition. Default to None. If the forward argument is False, it specifies final conditions.tg_ic (None or ndarray(float), optional) –
Initial (or final) conditions of the linear ODEs \(\dot{\boldsymbol{\delta x}} = \boldsymbol{\mathrm{J}}(t, \boldsymbol{x}) \cdot \boldsymbol{\delta x}\).
Can be a 1D, a 2D or a 3D array:
1D: Provide a single initial condition. This initial condition of the linear ODEs will be the same used for each initial condition ic of the ODEs \(\dot{\boldsymbol{x}} = \boldsymbol{f}(t, \boldsymbol{x})\) Should be of shape (n_dim,) where n_dim = \(\mathrm{dim}(\boldsymbol{x})\).
2D: Two sub-cases:
If tg_ic.shape[0]`=`ic.shape[0], assumes that each initial condition ic[i] of \(\dot{\boldsymbol{x}} = \boldsymbol{f}(t, \boldsymbol{x})\), correspond to a different initial condition tg_ic[i].
Else, assumes and integrate an ensemble of n_tg_traj initial condition of the linear ODEs for each initial condition of \(\dot{\boldsymbol{x}} = \boldsymbol{f}(t, \boldsymbol{x})\).
3D: An array of shape (n_traj, n_tg_traj, n_dim) which provide an ensemble of n_tg_ic initial conditions specific to each of the n_traj initial conditions of \(\dot{\boldsymbol{x}} = \boldsymbol{f}(t, \boldsymbol{x})\).
If None, use the identity matrix as initial condition, returning the fundamental matrix of solutions of the linear ODEs. Default to None. If the forward argument is False, it specifies final conditions.
forward (bool, optional) – If true, integrate the ODEs forward in time, else, integrate backward in time. In case of backward integration, the initial condition ic becomes a final condition. Default to forward integration.
adjoint (bool, optional) – If true, integrate the tangent \(\dot{\boldsymbol{\delta x}} = \boldsymbol{\mathrm{J}}(t, \boldsymbol{x}) \cdot \boldsymbol{\delta x}\) , else, integrate the adjoint linear model \(\dot{\boldsymbol{\delta x}} = \boldsymbol{\mathrm{J}}^T(t, \boldsymbol{x}) \cdot \boldsymbol{\delta x}\). Integrate the tangent model by default.
inverse (bool, optional) – Wheter or not to invert the Jacobian matrix \(\boldsymbol{\mathrm{J}}(t, \boldsymbol{x}) \rightarrow \boldsymbol{\mathrm{J}}^{-1}(t, \boldsymbol{x})\). False by default.
boundary (None or callable, optional) – Allow to add a inhomogeneous term to linear ODEs: \(\dot{\boldsymbol{\delta x}} = \boldsymbol{\mathrm{J}}(t, \boldsymbol{x}) \cdot \boldsymbol{\delta x} + \Psi(t, \boldsymbol{x})\). The boundary \(\Psi\) should have the same signature as \(\boldsymbol{\mathrm{J}}\), i.e.
func(t, x)
. If None, don’t add anything (homogeneous case). None by default.write_steps (int, optional) – Save the state of the integration in memory every write_steps steps. The other intermediary steps are lost. It determines the size of the returned objects. Default is 1. Set to 0 to return only the final state.
- set_bca(b=None, c=None, a=None, ic_init=True)[source]
Set the coefficients of the Runge-Kutta method and restart the integrator. s
- Parameters
b (None or ndarray, optional) – Vector of coefficients \(b_i\) of the Runge-Kutta method . If None, does not reinitialize these coefficients.
c (None or ndarray, optional) – Matrix of coefficients \(c_{i,j}\) of the Runge-Kutta method . If None, does not reinitialize these coefficients.
a (None or ndarray, optional) – Vector of coefficients \(a_i\) of the Runge-Kutta method . If None, does not reinitialize these coefficients.
ic_init (bool, optional) – Re-initialize or not the initial conditions of the integrator. Default to True.
- set_func(f, fjac, ic_init=True)[source]
Set the Numba-jitted function \(\boldsymbol{f}\) and Jacobian matrix function \(\boldsymbol{\mathrm{J}}\) to integrate.
- Parameters
f (callable) – The Numba-jitted function \(\boldsymbol{f}\). Should have the signature
f(t, x)
wherex
is the state value andt
is the time.fjac (callable) – The Numba-jitted Jacobian matrix function \(\boldsymbol{J}\). Should have the signature
J(t, x)
wherex
is the state value andt
is the time.ic_init (bool, optional) – Re-initialize or not the initial conditions of the integrator. Default to True.
Warning
This function restarts the integrator!
- set_ic(ic)[source]
Direct setter for the integrator’s non-linear ODEs initial conditions
- Parameters
Initial condition of the system. Can be a 1D or a 2D array:
1D: Provide a single initial condition. Should be of shape (n_dim,) where n_dim = \(\mathrm{dim}(\boldsymbol{x})\).
2D: Provide an ensemble of initial condition. Should be of shape (n_traj, n_dim) where n_dim = \(\mathrm{dim}(\boldsymbol{x})\), and where n_traj is the number of initial conditions.
- set_tg_ic(tg_ic)[source]
Direct setter for the integrator’s linear ODEs initial conditions
- Parameters
Initial condition of the linear ODEs \(\dot{\boldsymbol{\delta x}} = \boldsymbol{\mathrm{J}}(t, \boldsymbol{x}) \cdot \boldsymbol{\delta x}\).
Can be a 1D, a 2D or a 3D array:
1D: Provide a single initial condition. This initial condition of the linear ODEs will be the same used for each initial condition ic of the ODEs \(\dot{\boldsymbol{x}} = \boldsymbol{f}(t, \boldsymbol{x})\) Should be of shape (n_dim,) where n_dim = \(\mathrm{dim}(\boldsymbol{x})\).
2D: Two sub-cases:
If tg_ic.shape[0]`=`ic.shape[0], assumes that each initial condition ic[i] of \(\dot{\boldsymbol{x}} = \boldsymbol{f}(t, \boldsymbol{x})\), correspond to a different initial condition tg_ic[i].
Else, assumes and integrate an ensemble of n_tg_traj initial condition of the linear ODEs for each initial condition of \(\dot{\boldsymbol{x}} = \boldsymbol{f}(t, \boldsymbol{x})\).
3D: An array of shape (n_traj, n_tg_traj, n_dim) which provide an ensemble of n_tg_ic initial conditions specific to each of the n_traj initial conditions of \(\dot{\boldsymbol{x}} = \boldsymbol{f}(t, \boldsymbol{x})\).
- class qgs.integrators.integrator.TglsTrajectoryProcess(processID, func, func_jac, b, c, a, ics_queue, traj_queue)[source]
Bases:
Process
RungeKuttaTglsIntegrator
’s workers class. Allows to multi-thread time integration.- Parameters
processID (int) – Number identifying the worker.
func (callable) – Numba-jitted function to integrate assigned to the worker.
func_jac (callable) – Numba-jitted Jacobian matrix function to integrate assigned to the worker.
b (ndarray, optional) – Vector of coefficients \(b_i\) of the Runge-Kutta method .
c (ndarray, optional) – Matrix of coefficients \(c_{i,j}\) of the Runge-Kutta method .
a (ndarray, optional) – Vector of coefficients \(a_i\) of the Runge-Kutta method .
ics_queue (multiprocessing.JoinableQueue) – Queue to which the worker ask for initial conditions input.
traj_queue (multiprocessing.Queue) – Queue to which the worker returns the integration results.
- b
Vector of coefficients \(b_i\) of the Runge-Kutta method .
- Type
- c
Matrix of coefficients \(c_{i,j}\) of the Runge-Kutta method .
- Type
- a
Vector of coefficients \(a_i\) of the Runge-Kutta method .
- Type
- class qgs.integrators.integrator.TrajectoryProcess(processID, func, b, c, a, ics_queue, traj_queue)[source]
Bases:
Process
RungeKuttaIntegrator
’s workers class. Allows to multi-thread time integration.- Parameters
processID (int) – Number identifying the worker.
func (callable) – Numba-jitted function to integrate assigned to the worker.
b (ndarray, optional) – Vector of coefficients \(b_i\) of the Runge-Kutta method .
c (ndarray, optional) – Matrix of coefficients \(c_{i,j}\) of the Runge-Kutta method .
a (ndarray, optional) – Vector of coefficients \(a_i\) of the Runge-Kutta method .
ics_queue (multiprocessing.JoinableQueue) – Queue to which the worker ask for initial conditions and parameters input.
traj_queue (multiprocessing.Queue) – Queue to which the worker returns the integration results.
- b
Vector of coefficients \(b_i\) of the Runge-Kutta method .
- Type
- c
Matrix of coefficients \(c_{i,j}\) of the Runge-Kutta method .
- Type
- a
Vector of coefficients \(a_i\) of the Runge-Kutta method .
- Type