aphin.systems package

Submodules

aphin.systems.lti_systems module

Define and simulate dynamical systems e.g. LTI systems (LTISystem), port-Hamiltonian systems (PHSystem), …

class aphin.systems.lti_systems.DescrLTISystem(A, B=None, E=None)[source]

Bases: LTISystem

Linear time-invariant (LTI) system with descriptor formulation.

This class represents a linear time-invariant ODE system in a descriptor formulation, where the system is described by the equation: E * x’(t) = A * x(t) + B * u(t), with initial condition x(0) = x_init.

The descriptor matrix E allows for more general system representations, including differential-algebraic equations (DAEs) if E is singular.

solve(t, z_init, u=None, integrator_type='IMR', decomp_option='lu')[source]

Solve the descriptor linear time-invariant (LTI) system for given times, initial conditions, and inputs.

This method integrates the system’s differential-algebraic equations over the specified time steps.

Parameters:

  • t (array-like, shape (n_t,)) – Array of time steps for which the solution is computed.

  • z_init (array-like, shape (n,) or (n, n_s)) – Initial conditions of the system. If multiple simulations are run, this should have shape (n, n_s).

  • u (array-like, shape (n_t, n_u) or (n_t, n_u, n_s), optional) – Input signal for the system. If not provided, defaults to zero input. If multiple simulations are run, the shape should be (n_t, n_u, n_s).

  • integrator_type ({'IMR'}, optional) – The type of integrator to use. Only ‘IMR’ (Implicit Midpoint Rule) is supported for descriptor systems.

  • decomp_option ({'lu', 'linalg_solve'}, optional) – Option for decomposition in the IMR integrator. Options are ‘lu’ for LU decomposition or ‘linalg_solve’ for solving without decomposition.

Returns:

The solution of the ODE system at each time step.

Return type:

ndarray, shape (n_t, n, n_s)

class aphin.systems.lti_systems.LTISystem(A, B=None)[source]

Bases: object

General linear, time invariant ODE system (LTI system) with constant system matrices

get_system_matrix()[source]

Retrieve the system matrices A, B, and C.

This method returns the matrices that define the linear time-invariant (LTI) system.

Returns:

A tuple containing the system matrices: - A : array-like, shape (n, n)

The state matrix.

  • Barray-like, shape (n, n_u)

    The input matrix.

  • Carray-like, shape (n, n)

    The output matrix.

Return type:

tuple of array-like

static is_pos_def(M)[source]

Check if a matrix M is positive semi-definite.

This method checks whether all eigenvalues of the matrix M are non-negative. The matrix is considered positive semi-definite if the smallest eigenvalue is greater than or equal to a small threshold determined by the machine precision.

Parameters:

M (array-like, shape (n, n)) – The matrix to be checked for positive semi-definiteness.

Returns:

True if the matrix M is positive semi-definite, False otherwise.

Return type:

bool

is_regular(M)[source]

Check if a matrix M is regular.

A matrix is considered regular if it has full rank and a condition number below a specified threshold.

Parameters:

M (array-like, shape (n, n)) – The matrix to be checked for regularity.

Returns:

True if the matrix M is regular, False otherwise.

Return type:

bool

static is_skew_sym(M)[source]

Check if a matrix M is skew-symmetric.

This method checks whether the matrix M is equal to the negative of its transpose.

Parameters:

M (array-like, shape (n, n)) – The matrix to be checked for skew-symmetry.

Returns:

True if the matrix M is skew-symmetric, False otherwise.

Return type:

bool

static is_sym(M)[source]

Check if a matrix M is symmetric.

This method checks whether the matrix M is equal to its transpose.

Parameters:

M (array-like, shape (n, n)) – The matrix to be checked for symmetry.

Returns:

True if the matrix M is symmetric, False otherwise.

Return type:

bool

static quad(xl, M, xr=None)[source]

Calculate the quadratic form ( xl^T @ M @ xr ).

This method computes the quadratic form given by the expression: xl^T @ M @ xr, where xl and xr are vectors and M is a matrix.

Parameters:

  • xl (array-like, shape (n,)) – Vector on the left side of the matrix M.

  • M (array-like, shape (n, n)) – Matrix representing the quadratic form.

  • xr (array-like, shape (n,), optional) – Vector on the right side of the matrix M. If not provided, xr is assumed to be equal to xl.

Returns:

The result of the quadratic form.

Return type:

float

solve(t, z_init, u=None, integrator_type='IMR', decomp_option='lu')[source]

Solve the system for given times, initial conditions, and inputs.

Parameters:

  • t (array-like, shape (n_t,)) – Array with time steps for the solution.

  • z_init (array-like, shape (n,) or (n, n_sim)) – Initial conditions of the system. If the shape is (n,), it represents a single simulation; if the shape is (n, n_sim), it represents multiple simulations.

  • u (array-like, shape (n_t, n_u) or (n_sim, n_t, n_u), optional) – Input signal for the system. If not provided, it defaults to zeros. The shape should match the dimensions of the input matrix B.

  • integrator_type ({'IMR', 'lsim'}, optional) – The integration method to use. ‘IMR’ stands for implicit midpoint rule, and ‘lsim’ uses the matrix exponential. Default is ‘IMR’.

  • decomp_option ({'lu', 'linalg_solve'}, optional) – The decomposition option for the IMR integrator. ‘lu’ uses LU decomposition, and ‘linalg_solve’ solves without decomposition. Default is ‘lu’.

Returns:

z_sol – The solution of the ODE system. If multiple simulations are performed, the result is a 3D array where each slice corresponds to a simulation.

Return type:

ndarray, shape (n_sim, n_t, n_f)

solve_dt(t, z_init, U=None, integrator_type='IMR', decomp_option='lu')[source]

Solves the system of ordinary differential equations (ODEs) for given time steps, initial conditions, and inputs. Also computes the exact time derivative of the solution.

Parameters:

tndarray

Array of time steps with shape (n_t,).

z_initndarray

Initial conditions with shape (n_f, n_sim) or (n_f,).

Undarray, optional

Input signal with shape (n_t, n_u) or (n_sim, n_t, n_u). Default is None. If None, a zero input is used.

integrator_typestr, optional

Type of integrator to use. Choices are: - ‘IMR’ : Implicit Midpoint Rule - ‘LSIM’ : Matrix Exponential (Default is ‘IMR’).

decomp_optionstr, optional

Decomposition option for the Implicit Midpoint Rule. Choices are: - ‘lu’ : LU decomposition - ‘linalg_solve’ : Solve without decomposition (Default is ‘lu’).

Returns:

xndarray

Solution of the ODE system with shape (n_sim, n_t, n_f).

dx_dtndarray

Time derivative of the solution with shape (n_sim, n_t, n_f).

property stable

Check if the system is stable.

This property calculates the eigenvalues of the matrix E_inv@A and checks if the real part of the largest eigenvalue is less than a negative threshold, indicating stability.

Parameters:

eps (float, optional) – A small positive threshold to determine stability. Default is 1e-4.

Returns:

True if the system is stable, False otherwise.

Return type:

bool

aphin.systems.ph_systems module

Define and simulate dynamical systems e.g. LTI systems (LTISystem), port-Hamiltonian systems (PHSystem), …

class aphin.systems.ph_systems.CheckPHProperties[source]

Bases: object

A class for checking properties of matrices related to port-Hamiltonian (pH) systems.

check_pH_properties(J, R, Q, E=None, rtol=1e-05, atol=1e-08)[source]

Check if the matrices J, R, and Q (with optional descriptor E) satisfy port-Hamiltonian system properties.

This method verifies the following properties: - J is skew-symmetric. - R is symmetric positive definite. - transpose(E)@Q is symmetric positive definite

Parameters:

  • J (array-like, shape (n, n)) – Skew-symmetric matrix related to the port-Hamiltonian system.

  • R (array-like, shape (n, n)) – Symmetric positive definite matrix.

  • Q (array-like, shape (n, n), optional) – Matrix to be checked for positive definiteness.

  • E (array-like, shape (n, n), optional) – Descriptor matrix. If not provided, defaults to the identity matrix.

  • rtol (float, optional) – Relative tolerance for numerical comparison (default is 1e-05).

  • atol (float, optional) – Absolute tolerance for numerical comparison (default is 1e-08).

Returns:

True if all properties are satisfied, False otherwise.

Return type:

bool

check_spd(A, rtol=1e-06, atol=1e-08)[source]

Check if a matrix is symmetric positive definite.

This method verifies if the matrix A is symmetric and positive definite. It does so by attempting a Cholesky decomposition. If the decomposition fails due to the matrix not being positive definite, a small regularization term is added to the matrix to try and make it positive definite.

Parameters:

  • A (array-like, shape (n, n)) – The matrix to be checked for symmetric positive definiteness.

  • rtol (float, optional) – Relative tolerance for checking symmetry of the matrix (default is 1e-06).

  • atol (float, optional) – Absolute tolerance for checking symmetry of the matrix (default is 1e-08).

Returns:

True if the matrix is symmetric positive definite, False otherwise.

Return type:

bool

References

Adapted from [MorandinNicodemusUnger22].

class aphin.systems.ph_systems.DescrPHSystem(J_ph, R_ph, E, B=None, Q_ph=None)[source]

Bases: DescrLTISystem, CheckPHProperties

Linear port-Hamiltonian (pH) system in descriptor form with constant system matrices.

This class represents a linear port-Hamiltonian system described by the following differential equation: E * x’(t) = (J - R) * x(t) + B * u(t), with x(0) = x_init, where: - J is skew-symmetric: J = -J.T - R is symmetric positive semi-definite: R = R.T >= 0 - E is a descriptor matrix

H(x)[source]

Calculate the Hamiltonian for a time series of states.

Computes the Hamiltonian function H(x) = 0.5 * x.T @ (E.T @ Q_ph) @ x for each time step in the time series x.

Parameters:

x (ndarray, shape (n_t, n)) – Time series of states where n_t is the number of time steps and n is the dimension of the state vector.

Returns:

The Hamiltonian values computed for each time step in the time series x.

Return type:

ndarray, shape (n_t,)

solve(t, z_init, u=None, integrator_type='IMR', decomp_option='lu')[source]

Solve the descriptor port-Hamiltonian system using the specified integration method.

Computes the solution to the differential equation system E * x’(t) = (J - R) * x(t) + B * u(t) over the given time steps with the provided initial conditions and inputs.

Parameters:

  • t (ndarray, shape (n_t,)) – Array of time steps for the solution.

  • z_init (ndarray, shape (n,) or (n, n_s)) – Initial conditions for the system state.

  • u (ndarray, shape (n_t, n_u) or (n_t, n_u, n_s), optional) – Input signal to the system. If not provided, assumes zero input.

  • integrator_type (str, optional) – Integration method to use. Default is “IMR” (Implicit Midpoint Rule).

  • decomp_option (str, optional) – Decomposition method for the integrator. Default is “lu”. Other options may be available depending on the integrator.

Returns:

x – Solution of the ODE system at each time step.

Return type:

ndarray, shape (n_t, n, n_s)

class aphin.systems.ph_systems.PHSystem(J_ph, R_ph, B=None, Q_ph=None)[source]

Bases: LTISystem, CheckPHProperties

General linear, port-Hamiltonian (pH) system with constant system matrices.

This class represents a linear port-Hamiltonian system, defined by the differential equation:

x’(t) = (J - R) * x(t) + B * u(t), x(0) = x_init

where: - J is a skew-symmetric matrix (J = -J.T). - R is a symmetric positive semi-definite matrix (R = R.T >= 0).

H(x)[source]

Calculate the Hamiltonian for a time series of states.

Computes the Hamiltonian for a time series of states x, given by the quadratic form:

H(x) = 0.5 * x.T @ Q_ph @ x

Parameters:

x (array-like, shape (n_t, n)) – Time series of states where n_t is the number of time steps and n is the dimension of the state vector.

Returns:

  • ndarray, shape (n_t,)

  • The Hamiltonian values corresponding to each time step.

get_system_matrix()[source]

Retrieve the system matrices A, B, and C.

This method returns the matrices that define the linear time-invariant (LTI) system.

Returns:

A tuple containing the system matrices: - A : array-like, shape (n, n)

The state matrix.

  • Barray-like, shape (n, n_u)

    The input matrix.

  • Carray-like, shape (n, n)

    The output matrix.

Return type:

tuple of array-like

transform_pH_to_Q_identity(solver='Q', seed=1)[source]

Transform the current port-Hamiltonian (pH) system to a Q = I (identity matrix) form, if possible.

This method attempts to transform the port-Hamiltonian system described by the matrices J_ph, R_ph, Q_ph, B_ph, and C_ph to a new system where Q is replaced by the identity matrix. If the transformation is successful, it returns a new PHSystem instance with the transformed matrices. If the transformation fails, it logs a message and returns the current instance.

Parameters:

  • solver ({'scipy', 'pymor', 'Q'}, optional) – Method for solving the transformation: - ‘scipy’ or ‘pymor’: Use the Riccati equation solver. - ‘Q’: Use the provided Q matrix directly.

  • seed (int, optional) – Random seed for generating a dummy input matrix when solving the Riccati equation.

Returns:

  • PHSystem or self – If the transformation to Q = I is successful, returns a new PHSystem instance with the transformed matrices. Otherwise, returns the current instance.

  • T (ndarray, shape (n, n)) – The transformation matrix, only returned if the transformation is successful.

  • T_inv (ndarray, shape (n, n)) – The inverse of the transformation matrix, only returned if the transformation is successful.

Module contents