nutopy.ivp.jexp

jexp(df, tf, t0, x0, dx0, pars=None, *, options=None, time_steps=None, args=(), kwargs={})[source]

The jexp function computes a numerical approximation of the solution of the linearization of the initial value problem

\[\frac{d}{d t} x = \dot{x} = f(t, x, pars), \quad x(t_0) = x_0,\]

where \(pars \in \mathbf{R}^k\) is a vector of parameters, \(t\) is the time and \(x\) is the state variable. The time \(t_0\) is refered as the initial time, \(x_0\) as the initial condition and \(f\) as the dynamics.

Hence, dexp approximates the solution of

\[\begin{split}\dot{x} &= f(t, x, pars), \\ \dot{\delta x} &= \frac{\partial f}{\partial x}(t, x, pars) \cdot \delta x, \\ x(t_0) &= x_0, \\ \delta x(t_0) &= \delta x_0,\end{split}\]

where \(\delta x_0\) may be a vector of \(\mathbf{R}^n\) (\(n\) being the dimension of the state \(x\)) or a real matrix of \(n\) rows and \(k\) columns.

The solution at the final time \(t_f\) is denoted \((x(t_f, t_0, x_0, pars), \delta x(t_f, t_0, x_0, \delta x_0, pars))\).

Parameters
  • df (callable) –

    A vector function

    y, dy = df(t, x, dx) if pars is None, where y = f(t, x) and dy = df/dx(t, x) . dx

    or

    y, dy = df(t, x, dx, pars) if pars is not None, where y = f(t, x, pars) and dy = df/dx(t, x, pars) . dx.

    Remark. In the definition of df, the argument dx is a vector of the same shape as x.

  • tf (float) – Final time

  • t0 (float) – Initial time`

  • x0 (float vector) – Initial condition

  • dx0 (float vector or matrix) – Initial condition \(\delta x_0\) of the linear part of the equation

  • pars (float vector, optional, positional) – Parameters

  • options (nutopy.ivp.Options, optional, key only) – A dictionnary of solver options. The integrator scheme is chosen setting SolverMethod.

  • time_steps (float vector, optional, key only) –

    If time_steps is not None, then it is a vector of time steps, that is sol.tout = time_steps. If time_steps is None, then sol.tout = [t0, tf] if a fixed-step integrator scheme is used while sol.tout is given by the integrator if a variable step-size scheme is used.

    Remark. Note that in the case of a variable step-size integrator scheme, time_steps is used only for output, that is the steps are computed by the integrator but dense output is used to provide by interpolation the solution x at times in time_steps.

  • args (tuple, key only) – optional arguments for df

  • kwargs (dictionnary, key only) – keywords optional arguments for df

Returns

sol – A dictionary/struct of outputs:

  • sol.xf or sol.get('xf') : final point x(tf, t0, x0, pars);

  • sol.dxf : principal linear part dx(tf, t0, dx0, x0, pars);

  • sol.tout : the integration time-steps: tout=time_steps if provided;

  • sol.xout : the integration points: xout[i] = x(tout[i], t0, x0, pars);

  • sol.dxout : principal linear part at integration time-steps: dxout[i] = dx(tout[i], t0, x0, dx0, pars);

  • sol.success : whether or not the integrator exited successfully;

  • sol.status : termination status of the integrator;

  • sol.message : description of the cause of the termination;

  • sol.nfev : number of evaluations of the dynamics;

  • sol.nsteps : number of integration time-steps.

Notes

  • Use numpy.ndarray to define float vectors.

  • If dx0 is not a matrix but a vector, then jexp(df, tf, t0, x0, dx0, pars) = dexp(df, tf, 0, t0, 0, x0, dx0, pars, 0).

Examples

Example: dx/dx0

>>> import numpy as np
>>> from nutopy import ivp
>>> def df(t, x, dx):
...     y     = np.zeros(x.size)
...     y[0]  = -x[0] + x[1]
...     y[1]  =  x[1]
...     dy    = np.zeros(x.size)
...     dy[0] = -dx[0] + dx[1]
...     dy[1] =  dx[1]
...     return (y, dy)
>>> tf = 1.0
>>> t0 = 0.0
>>> x0 = np.array([-1.5, 0.5])
>>> dx0 = np.eye(x0.size)
>>> sol = ivp.jexp(df, tf, t0, x0, dx0)
>>> print('dx/dx0(tf, t0, x0) = ', sol.dxf)
dx/dx0(tf, t0, x0) =  [[0.36787944 1.17520119]
[0.         2.71828183]]
Raises