pymor.models.iosys
¶
Module Contents¶
- class pymor.models.iosys.BilinearModel(A, N, B, C, D, E=None, sampling_time=0, error_estimator=None, visualizer=None, name=None)[source]¶
Bases:
pymor.models.interface.Model
Class for bilinear systems.
This class describes input-output systems given by
\[\begin{split}E x'(t) & = A x(t) + \sum_{i = 1}^m{N_i x(t) u_i(t)} + B u(t), \\ y(t) & = C x(t) + D u(t),\end{split}\]if continuous-time, or
\[\begin{split}E x(k + 1) & = A x(k) + \sum_{i = 1}^m{N_i x(k) u_i(k)} + B u(k), \\ y(k) & = C x(k) + D u(t),\end{split}\]if discrete-time, where \(E\), \(A\), \(N_i\), \(B\), \(C\), and \(D\) are linear operators and \(m\) is the number of inputs.
Parameters
- A
The
Operator
A.- N
The tuple of
Operators
N_i.- B
The
Operator
B.- C
The
Operator
C.- D
The
Operator
D orNone
(then D is assumed to be zero).- E
The
Operator
E orNone
(then E is assumed to be identity).- sampling_time
0
if the system is continuous-time, otherwise a positive number that denotes the sampling time (in seconds).- error_estimator
An error estimator for the problem. This can be any object with an
estimate_error(U, mu, model)
method. Iferror_estimator
is notNone
, anestimate_error(U, mu)
method is added to the model which will callerror_estimator.estimate_error(U, mu, self)
.- visualizer
A visualizer for the problem. This can be any object with a
visualize(U, model, ...)
method. Ifvisualizer
is notNone
, avisualize(U, *args, **kwargs)
method is added to the model which forwards its arguments to the visualizer’svisualize
method.- name
Name of the system.
- class pymor.models.iosys.LTIModel(A, B, C, D=None, E=None, sampling_time=0, T=None, initial_data=None, time_stepper=None, num_values=None, presets=None, solver_options=None, ast_pole_data=None, error_estimator=None, visualizer=None, name=None)[source]¶
Bases:
pymor.models.interface.Model
Class for linear time-invariant systems.
This class describes input-state-output systems given by
\[\begin{split}E(\mu) \dot{x}(t, \mu) & = A(\mu) x(t, \mu) + B(\mu) u(t), \\ y(t, \mu) & = C(\mu) x(t, \mu) + D(\mu) u(t),\end{split}\]if continuous-time, or
\[\begin{split}E(\mu) x(k + 1, \mu) & = A(\mu) x(k, \mu) + B(\mu) u(k), \\ y(k, \mu) & = C(\mu) x(k, \mu) + D(\mu) u(k),\end{split}\]if discrete-time, where \(A\), \(B\), \(C\), \(D\), and \(E\) are linear operators.
All methods related to the transfer function (e.g., frequency response calculation and Bode plots) are attached to the
transfer_function
attribute.Parameters
- A
The
Operator
A.- B
The
Operator
B.- C
The
Operator
C.- D
The
Operator
D orNone
(then D is assumed to be zero).- E
The
Operator
E orNone
(then E is assumed to be identity).- sampling_time
0
if the system is continuous-time, otherwise a positive number that denotes the sampling time (in seconds).- T
The final time T.
- initial_data
The initial data
x_0
. Either aVectorArray
of length 1 or (for theParameter
-dependent case) a vector-likeOperator
(i.e. a linearOperator
withsource.dim == 1
) which applied toNumpyVectorArray(np.array([1]))
will yield the initial data for givenparameter values
. IfNone
, it is assumed to be zero.- time_stepper
The
time-stepper
to be used bysolve
.- num_values
The number of returned vectors of the solution trajectory. If
None
, each intermediate vector that is calculated is returned.- presets
A
dict
of preset attributes orNone
. The dict must only contain keys that correspond to attributes ofLTIModel
such aspoles
,c_lrcf
,o_lrcf
,c_dense
,o_dense
,hsv
,h2_norm
,hinf_norm
,l2_norm
andlinf_norm
. Additionally, the frequency at which the \(\mathcal{H}_\infty/\mathcal{L}_\infty\) norm is attained can be preset withfpeak
.- solver_options
The solver options to use to solve matrix equations.
- ast_pole_data
Used in
get_ast_spectrum
. Can be:dictionary of parameters for
eigs
,list of anti-stable eigenvalues (scalars),
- tuple
(lev, ew, rev)
whereew
contains the sorted anti-stable eigenvalues and
lev
andrev
areVectorArrays
representing the eigenvectors.
- tuple
None
if anti-stable eigenvalues should be computed via dense methods.
- error_estimator
An error estimator for the problem. This can be any object with an
estimate_error(U, mu, model)
method. Iferror_estimator
is notNone
, anestimate_error(U, mu)
method is added to the model which will callerror_estimator.estimate_error(U, mu, self)
.- visualizer
A visualizer for the problem. This can be any object with a
visualize(U, model, ...)
method. Ifvisualizer
is notNone
, avisualize(U, *args, **kwargs)
method is added to the model which forwards its arguments to the visualizer’svisualize
method.- name
Name of the system.
Methods
Create
LTIModel
from matrices stored in .[ABCDE] files.Create
LTIModel
from matrices stored in separate files.Create
LTIModel
from matrices stored in a .mat file.Create
LTIModel
from matrices.Compute anti-stable subset of the poles of the
LTIModel
.Compute a Gramian.
Compute the \(\mathcal{H}_2\)-norm of the
LTIModel
.Compute the Hankel-norm of the
LTIModel
.Compute the \(\mathcal{H}_\infty\)-norm of the
LTIModel
.Hankel singular values.
Compute impulse response from all inputs.
Compute the \(\mathcal{L}_2\)-norm of the
LTIModel
.Compute the \(\mathcal{L}_\infty\)-norm of the
LTIModel
.Create a transformed
LTIModel
by applying an arbitrary Moebius transformation.Compute system poles.
Compute step response from all inputs.
Save operators as matrices to .[ABCDE] files in Matrix Market format.
Converts a discrete-time
LTIModel
to a continuous-timeLTIModel
.Converts a continuous-time
LTIModel
to a discrete-timeLTIModel
.Write operators to files as matrices.
Save operators as matrices to .mat file.
Return operators as matrices.
- classmethod from_abcde_files(files_basename, sampling_time=0, T=None, time_stepper=None, num_values=None, presets=None, state_id='STATE', solver_options=None, error_estimator=None, visualizer=None, name=None)[source]¶
Create
LTIModel
from matrices stored in .[ABCDE] files.Parameters
- files_basename
The basename of files containing A, B, C, and optionally D and E.
- sampling_time
0
if the system is continuous-time, otherwise a positive number that denotes the sampling time (in seconds).- T
The final time T.
- time_stepper
The
time-stepper
to be used bysolve
.- num_values
The number of returned vectors of the solution trajectory. If
None
, each intermediate vector that is calculated is returned.- presets
A
dict
of preset attributes orNone
. See__init__
.- state_id
Id of the state space.
- solver_options
The solver options to use to solve the Lyapunov equations.
- error_estimator
An error estimator for the problem. This can be any object with an
estimate_error(U, mu, model)
method. Iferror_estimator
is notNone
, anestimate_error(U, mu)
method is added to the model which will callerror_estimator.estimate_error(U, mu, self)
.- visualizer
A visualizer for the problem. This can be any object with a
visualize(U, model, ...)
method. Ifvisualizer
is notNone
, avisualize(U, *args, **kwargs)
method is added to the model which forwards its arguments to the visualizer’svisualize
method.- name
Name of the system.
Returns
- lti
The
LTIModel
with operators A, B, C, D, and E.
- classmethod from_files(A_file, B_file, C_file, D_file=None, E_file=None, sampling_time=0, T=None, initial_data_file=None, time_stepper=None, num_values=None, presets=None, state_id='STATE', solver_options=None, error_estimator=None, visualizer=None, name=None)[source]¶
Create
LTIModel
from matrices stored in separate files.Parameters
- A_file
The name of the file (with extension) containing A.
- B_file
The name of the file (with extension) containing B.
- C_file
The name of the file (with extension) containing C.
- D_file
None
or the name of the file (with extension) containing D.- E_file
None
or the name of the file (with extension) containing E.- sampling_time
0
if the system is continuous-time, otherwise a positive number that denotes the sampling time (in seconds).- T
The final time T.
- initial_data_file
None
or the name of the file (with extension) containing the initial data.- time_stepper
The
time-stepper
to be used bysolve
.- num_values
The number of returned vectors of the solution trajectory. If
None
, each intermediate vector that is calculated is returned.- presets
A
dict
of preset attributes orNone
. See__init__
.- state_id
Id of the state space.
- solver_options
The solver options to use to solve the Lyapunov equations.
- error_estimator
An error estimator for the problem. This can be any object with an
estimate_error(U, mu, model)
method. Iferror_estimator
is notNone
, anestimate_error(U, mu)
method is added to the model which will callerror_estimator.estimate_error(U, mu, self)
.- visualizer
A visualizer for the problem. This can be any object with a
visualize(U, model, ...)
method. Ifvisualizer
is notNone
, avisualize(U, *args, **kwargs)
method is added to the model which forwards its arguments to the visualizer’svisualize
method.- name
Name of the system.
Returns
- lti
The
LTIModel
with operators A, B, C, D, and E.
- classmethod from_mat_file(file_name, sampling_time=0, T=None, time_stepper=None, num_values=None, presets=None, state_id='STATE', solver_options=None, error_estimator=None, visualizer=None, name=None)[source]¶
Create
LTIModel
from matrices stored in a .mat file.Supports the format used in the SLICOT benchmark collection.
Parameters
- file_name
The name of the .mat file (extension .mat does not need to be included) containing A, B, and optionally C, D, and E.
- sampling_time
0
if the system is continuous-time, otherwise a positive number that denotes the sampling time (in seconds).- T
The final time T.
- time_stepper
The
time-stepper
to be used bysolve
.- num_values
The number of returned vectors of the solution trajectory. If
None
, each intermediate vector that is calculated is returned.- presets
A
dict
of preset attributes orNone
. See__init__
.- state_id
Id of the state space.
- solver_options
The solver options to use to solve the Lyapunov equations.
- error_estimator
An error estimator for the problem. This can be any object with an
estimate_error(U, mu, model)
method. Iferror_estimator
is notNone
, anestimate_error(U, mu)
method is added to the model which will callerror_estimator.estimate_error(U, mu, self)
.- visualizer
A visualizer for the problem. This can be any object with a
visualize(U, model, ...)
method. Ifvisualizer
is notNone
, avisualize(U, *args, **kwargs)
method is added to the model which forwards its arguments to the visualizer’svisualize
method.- name
Name of the system.
Returns
- lti
The
LTIModel
with operators A, B, C, D, and E.
- classmethod from_matrices(A, B, C, D=None, E=None, sampling_time=0, T=None, initial_data=None, time_stepper=None, num_values=None, presets=None, state_id='STATE', solver_options=None, error_estimator=None, visualizer=None, name=None)[source]¶
Create
LTIModel
from matrices.Parameters
- A
The
NumPy array
orSciPy spmatrix
A.- B
The
NumPy array
orSciPy spmatrix
B.- C
The
NumPy array
orSciPy spmatrix
C.- D
The
NumPy array
orSciPy spmatrix
D orNone
(then D is assumed to be zero).- E
The
NumPy array
orSciPy spmatrix
E orNone
(then E is assumed to be identity).- sampling_time
0
if the system is continuous-time, otherwise a positive number that denotes the sampling time (in seconds).- T
The final time T.
- initial_data
The initial data
x_0
as aNumPy array
. IfNone
, it is assumed to be zero.- time_stepper
The
time-stepper
to be used bysolve
.- num_values
The number of returned vectors of the solution trajectory. If
None
, each intermediate vector that is calculated is returned.- presets
A
dict
of preset attributes orNone
. See__init__
.- state_id
Id of the state space.
- solver_options
The solver options to use to solve the Lyapunov equations.
- error_estimator
An error estimator for the problem. This can be any object with an
estimate_error(U, mu, model)
method. Iferror_estimator
is notNone
, anestimate_error(U, mu)
method is added to the model which will callerror_estimator.estimate_error(U, mu, self)
.- visualizer
A visualizer for the problem. This can be any object with a
visualize(U, model, ...)
method. Ifvisualizer
is notNone
, avisualize(U, *args, **kwargs)
method is added to the model which forwards its arguments to the visualizer’svisualize
method.- name
Name of the system.
Returns
- lti
The
LTIModel
with operators A, B, C, D, and E.
- get_ast_spectrum(mu=None)[source]¶
Compute anti-stable subset of the poles of the
LTIModel
.Parameters
- mu
Returns
- lev
VectorArray
of left eigenvectors.- ew
One-dimensional
NumPy array
of anti-stable eigenvalues sorted from smallest to largest.- rev
VectorArray
of right eigenvectors.
- gramian(typ, mu=None)[source]¶
Compute a Gramian.
Parameters
- typ
The type of the Gramian:
'c_lrcf'
: low-rank Cholesky factor of the controllability Gramian,'o_lrcf'
: low-rank Cholesky factor of the observability Gramian,'c_dense'
: dense controllability Gramian,'o_dense'
: dense observability Gramian,'bs_c_lrcf'
: low-rank Cholesky factor of the Bernoulli stabilized controllability Gramian,'bs_o_lrcf'
: low-rank Cholesky factor of the Bernoulli stabilized observability Gramian,'lqg_c_lrcf'
: low-rank Cholesky factor of the “controllability” LQG Gramian,'lqg_o_lrcf'
: low-rank Cholesky factor of the “observability” LQG Gramian,('br_c_lrcf', gamma)
: low-rank Cholesky factor of the “controllability” bounded real Gramian,('br_o_lrcf', gamma)
: low-rank Cholesky factor of the “observability” bounded real Gramian.
Note
For
'*_lrcf'
types, the method assumes the system is asymptotically stable. For'*_dense'
types, the method assumes that the underlying Lyapunov equation has a unique solution, i.e. no pair of system poles adds to zero in the continuous-time case and no pair of system poles multiplies to one in the discrete-time case.- mu
Returns
If typ ends with
'_lrcf'
, then the Gramian factor as aVectorArray
fromself.A.source
. If typ ends with'_dense'
, then the Gramian as aNumPy array
.
- h2_norm(mu=None)[source]¶
Compute the \(\mathcal{H}_2\)-norm of the
LTIModel
.Note
Assumes the system is asymptotically stable.
Parameters
Returns
- norm
\(\mathcal{H}_2\)-norm.
- hankel_norm(mu=None)[source]¶
Compute the Hankel-norm of the
LTIModel
.Note
Assumes the system is asymptotically stable.
Parameters
Returns
- norm
Hankel-norm.
- hinf_norm(mu=None, return_fpeak=False, ab13dd_equilibrate=False)[source]¶
Compute the \(\mathcal{H}_\infty\)-norm of the
LTIModel
.Note
Assumes the system is asymptotically stable. Under this is assumption the \(\mathcal{H}_\infty\)-norm is equal to the \(\mathcal{H}_\infty\)-norm. Accordingly, this method calls
linf_norm
.Parameters
- mu
- return_fpeak
Whether to return the frequency at which the maximum is achieved.
- ab13dd_equilibrate
Whether
slycot.ab13dd
should use equilibration.
Returns
- norm
\(\mathcal{H}_\infty\)-norm.
- fpeak
Frequency at which the maximum is achieved (if
return_fpeak
isTrue
).
- hsv(mu=None)[source]¶
Hankel singular values.
Note
Assumes the system is asymptotically stable.
Parameters
Returns
- sv
One-dimensional
NumPy array
of singular values.
- impulse_resp(mu=None, return_solution=False)[source]¶
Compute impulse response from all inputs.
Parameters
- mu
Parameter values
for which to solve.- return_solution
If
True
, the modelsolution
for the givenparameter values
mu
is returned.
Returns
- y
Impulse response as a 3D
NumPy array
wherey.shape[0]
is the number of time steps,y.shape[1]
is the number of outputs, andy.shape[2]
is the number of inputs.- Xs
The tuple of solution
VectorArrays
for every input. Returned only whenreturn_solution
isTrue
.
- l2_norm(mu=None)[source]¶
Compute the \(\mathcal{L}_2\)-norm of the
LTIModel
.The \(\mathcal{L}_2\)-norm of an
LTIModel
is defined via the integral\[\lVert H \rVert_{\mathcal{L}_2} = \left( \frac{1}{2 \pi} \int_{-\infty}^{\infty} \lVert H(\boldsymbol{\imath} \omega) \rVert_{\operatorname{F}}^2 \operatorname{d}\!\omega \right)^{\frac{1}{2}}.\]Parameters
- mu
Returns
- norm
\(\mathcal{L}_2\)-norm.
- linf_norm(mu=None, return_fpeak=False, ab13dd_equilibrate=False)[source]¶
Compute the \(\mathcal{L}_\infty\)-norm of the
LTIModel
.The \(\mathcal{L}_\infty\)-norm of an
LTIModel
is defined via\[\lVert H \rVert_{\mathcal{L}_\infty} = \sup_{\omega \in \mathbb{R}} \lVert H(\boldsymbol{\imath} \omega) \rVert_2.\]Parameters
- mu
- return_fpeak
Whether to return the frequency at which the maximum is achieved.
- ab13dd_equilibrate
Whether
slycot.ab13dd
should use equilibration.
Returns
- norm
\(\mathcal{L}_\infty\)-norm.
- fpeak
Frequency at which the maximum is achieved (if
return_fpeak
isTrue
).
- moebius_substitution(M, sampling_time=0)[source]¶
Create a transformed
LTIModel
by applying an arbitrary Moebius transformation.This method returns a transformed
LTIModel
such that the transfer function of the original and transformedLTIModel
are related by a Moebius substitution of the frequency argument:\[H(s)=\tilde{H}(M(s)),\]where
\[M(s) = \frac{as+b}{cs+d}\]is a Moebius transformation. See [CCA96] for details.
Parameters
- M
The
MoebiusTransformation
that defines the frequency mapping.- sampling_time
The sampling time of the transformed system (in seconds).
0
if the system is continuous-time, otherwise a positive number. Defaults to zero.
Returns
- sys
The transformed
LTIModel
.
- poles(mu=None)[source]¶
Compute system poles.
Note
Assumes the systems is small enough to use a dense eigenvalue solver.
Parameters
- mu
Parameter values
for which to compute the systems poles.
Returns
One-dimensional
NumPy array
of system poles.
- step_resp(mu=None, return_solution=False)[source]¶
Compute step response from all inputs.
Parameters
- mu
Parameter values
for which to solve.- return_solution
If
True
, the model solution for the givenparameter values
mu
is returned.
Returns
- y
Step response as a 3D
NumPy array
wherey.shape[0]
is the number of time steps,y.shape[1]
is the number of outputs, andy.shape[2]
is the number of inputs.- Xs
The tuple of solution
VectorArrays
for every input. Returned only whenreturn_solution
isTrue
.
- to_abcde_files(files_basename)[source]¶
Save operators as matrices to .[ABCDE] files in Matrix Market format.
Parameters
- files_basename
The basename of files containing the operators.
- to_continuous(method='Tustin', w0=0)[source]¶
Converts a discrete-time
LTIModel
to a continuous-timeLTIModel
.Parameters
- method
A string that defines the transformation method. At the moment only Tustin’s method is supported.
- w0
If
method=='Tustin'
, this parameter can be used to specify the prewarping-frequency. Defaults to zero.
Returns
- sys
Continuous-time
LTIModel
.
- to_discrete(sampling_time, method='Tustin', w0=0)[source]¶
Converts a continuous-time
LTIModel
to a discrete-timeLTIModel
.Parameters
- sampling_time
A positive number that denotes the sampling time of the resulting system (in seconds).
- method
A string that defines the transformation method. At the moment only Tustin’s method is supported.
- w0
If
method=='Tustin'
, this parameter can be used to specify the prewarping-frequency. Defaults to zero.
Returns
- sys
Discrete-time
LTIModel
.
- to_files(A_file, B_file, C_file, D_file=None, E_file=None)[source]¶
Write operators to files as matrices.
Parameters
- A_file
The name of the file (with extension) containing A.
- B_file
The name of the file (with extension) containing B.
- C_file
The name of the file (with extension) containing C.
- D_file
The name of the file (with extension) containing D or
None
if D is aZeroOperator
.- E_file
The name of the file (with extension) containing E or
None
if E is anIdentityOperator
.
- to_mat_file(file_name)[source]¶
Save operators as matrices to .mat file.
Parameters
- file_name
The name of the .mat file (extension .mat does not need to be included).
- to_matrices()[source]¶
Return operators as matrices.
Returns
- A
The
NumPy array
orSciPy spmatrix
A.- B
The
NumPy array
orSciPy spmatrix
B.- C
The
NumPy array
orSciPy spmatrix
C.- D
The
NumPy array
orSciPy spmatrix
D orNone
(if D is aZeroOperator
).- E
The
NumPy array
orSciPy spmatrix
E orNone
(if E is anIdentityOperator
).
- class pymor.models.iosys.LinearDelayModel(A, Ad, tau, B, C, D=None, E=None, sampling_time=0, error_estimator=None, visualizer=None, name=None)[source]¶
Bases:
pymor.models.interface.Model
Class for linear delay systems.
This class describes input-state-output systems given by
\[\begin{split}E x'(t) & = A x(t) + \sum_{i = 1}^q{A_i x(t - \tau_i)} + B u(t), \\ y(t) & = C x(t) + D u(t),\end{split}\]if continuous-time, or
\[\begin{split}E x(k + 1) & = A x(k) + \sum_{i = 1}^q{A_i x(k - \tau_i)} + B u(k), \\ y(k) & = C x(k) + D u(k),\end{split}\]if discrete-time, where \(E\), \(A\), \(A_i\), \(B\), \(C\), and \(D\) are linear operators.
All methods related to the transfer function (e.g., frequency response calculation and Bode plots) are attached to the
transfer_function
attribute.Parameters
- A
The
Operator
A.- Ad
The tuple of
Operators
A_i.- tau
The tuple of delay times (positive floats or ints).
- B
The
Operator
B.- C
The
Operator
C.- D
The
Operator
D orNone
(then D is assumed to be zero).- E
The
Operator
E orNone
(then E is assumed to be identity).- sampling_time
0
if the system is continuous-time, otherwise a positive number that denotes the sampling time (in seconds).- error_estimator
An error estimator for the problem. This can be any object with an
estimate_error(U, mu, model)
method. Iferror_estimator
is notNone
, anestimate_error(U, mu)
method is added to the model which will callerror_estimator.estimate_error(U, mu, self)
.- visualizer
A visualizer for the problem. This can be any object with a
visualize(U, model, ...)
method. Ifvisualizer
is notNone
, avisualize(U, *args, **kwargs)
method is added to the model which forwards its arguments to the visualizer’svisualize
method.- name
Name of the system.
- class pymor.models.iosys.LinearStochasticModel(A, As, B, C, D=None, E=None, sampling_time=0, error_estimator=None, visualizer=None, name=None)[source]¶
Bases:
pymor.models.interface.Model
Class for linear stochastic systems.
This class describes input-state-output systems given by
\[\begin{split}E \mathrm{d}x(t) & = A x(t) \mathrm{d}t + \sum_{i = 1}^q{A_i x(t) \mathrm{d}\omega_i(t)} + B u(t) \mathrm{d}t, \\ y(t) & = C x(t) + D u(t),\end{split}\]if continuous-time, or
\[\begin{split}E x(k + 1) & = A x(k) + \sum_{i = 1}^q{A_i x(k) \omega_i(k)} + B u(k), \\ y(k) & = C x(k) + D u(t),\end{split}\]if discrete-time, where \(E\), \(A\), \(A_i\), \(B\), \(C\), and \(D\) are linear operators and \(\omega_i\) are stochastic processes.
Parameters
- A
The
Operator
A.- As
The tuple of
Operators
A_i.- B
The
Operator
B.- C
The
Operator
C.- D
The
Operator
D orNone
(then D is assumed to be zero).- E
The
Operator
E orNone
(then E is assumed to be identity).- sampling_time
0
if the system is continuous-time, otherwise a positive number that denotes the sampling time (in seconds).- error_estimator
An error estimator for the problem. This can be any object with an
estimate_error(U, mu, model)
method. Iferror_estimator
is notNone
, anestimate_error(U, mu)
method is added to the model which will callerror_estimator.estimate_error(U, mu, self)
.- visualizer
A visualizer for the problem. This can be any object with a
visualize(U, model, ...)
method. Ifvisualizer
is notNone
, avisualize(U, *args, **kwargs)
method is added to the model which forwards its arguments to the visualizer’svisualize
method.- name
Name of the system.
- class pymor.models.iosys.PHLTIModel(J, R, G, P=None, S=None, N=None, E=None, solver_options=None, error_estimator=None, visualizer=None, name=None)[source]¶
Bases:
pymor.models.interface.Model
Class for (continuous) port-Hamiltonian linear time-invariant systems.
This class describes input-state-output systems given by
\[\begin{split}E(\mu) \dot{x}(t, \mu) & = (J(\mu) - R(\mu)) x(t, \mu) + (G(\mu) - P(\mu)) u(t), \\ y(t, \mu) & = (G(\mu) + P(\mu))^T x(t, \mu) + (S(\mu) - N(\mu)) u(t),\end{split}\]with \(E(\mu) \succeq 0\), \(J(\mu) = -J(\mu)^T\), \(N(\mu) = -N(\mu)^T\) and
\[\begin{split}\mathcal{R}(\mu) = \begin{bmatrix} R(\mu) & P(\mu) \\ P(\mu)^T & S(\mu) \end{bmatrix} \succeq 0.\end{split}\]All methods related to the transfer function (e.g., frequency response calculation and Bode plots) are attached to the
transfer_function
attribute.Parameters
- J
The
Operator
J.- R
The
Operator
R.- G
The
Operator
G.- P
The
Operator
P orNone
(then P is assumed to be zero).- S
The
Operator
S orNone
(then S is assumed to be zero).- N
The
Operator
N orNone
(then N is assumed to be zero).- E
The
Operator
E orNone
(then E is assumed to be identity).- solver_options
The solver options to use to solve the Lyapunov equations.
- name
Name of the system.
Methods
Create
PHLTIModel
from matrices.Compute a Gramian.
Compute the H2-norm.
Compute the Hankel-norm.
Compute the H_infinity-norm.
Hankel singular values.
Compute system poles.
Return a standard linear time-invariant system representation.
Return operators as matrices.
- classmethod from_matrices(J, R, G, P=None, S=None, N=None, E=None, state_id='STATE', solver_options=None, error_estimator=None, visualizer=None, name=None)[source]¶
Create
PHLTIModel
from matrices.Parameters
- J
The
NumPy array
orSciPy spmatrix
J.- R
The
NumPy array
orSciPy spmatrix
R.- G
The
NumPy array
orSciPy spmatrix
G.- P
The
NumPy array
orSciPy spmatrix
P orNone
(then P is assumed to be zero).- S
The
NumPy array
orSciPy spmatrix
S orNone
(then S is assumed to be zero).- N
The
NumPy array
orSciPy spmatrix
N orNone
(then N is assumed to be zero).- E
The
NumPy array
orSciPy spmatrix
E orNone
(then E is assumed to be identity).- state_id
Id of the state space.
- solver_options
The solver options to use to solve the Lyapunov equations.
- error_estimator
An error estimator for the problem. This can be any object with an
estimate_error(U, mu, model)
method. Iferror_estimator
is notNone
, anestimate_error(U, mu)
method is added to the model which will callerror_estimator.estimate_error(U, mu, self)
.- visualizer
A visualizer for the problem. This can be any object with a
visualize(U, model, ...)
method. Ifvisualizer
is notNone
, avisualize(U, *args, **kwargs)
method is added to the model which forwards its arguments to the visualizer’svisualize
method.- name
Name of the system.
Returns
- phlti
The
PHLTIModel
with operators J, R, G, P, S, N, and E.
- gramian(typ, mu=None)[source]¶
Compute a Gramian.
Parameters
- typ
The type of the Gramian:
'c_lrcf'
: low-rank Cholesky factor of the controllability Gramian,'o_lrcf'
: low-rank Cholesky factor of the observability Gramian,'c_dense'
: dense controllability Gramian,'o_dense'
: dense observability Gramian.
Note
For
'*_lrcf'
types, the method assumes the system is asymptotically stable. For'*_dense'
types, the method assumes that the underlying Lyapunov equation has a unique solution, i.e. no pair of system poles adds to zero in the continuous-time case and no pair of system poles multiplies to one in the discrete-time case.- mu
Returns
If typ is
'c_lrcf'
or'o_lrcf'
, then the Gramian factor as aVectorArray
fromself.A.source
. If typ is'c_dense'
or'o_dense'
, then the Gramian as aNumPy array
.
- h2_norm(mu=None)[source]¶
Compute the H2-norm.
Note
Assumes the system is asymptotically stable.
Parameters
Returns
- norm
H_2-norm.
- hankel_norm(mu=None)[source]¶
Compute the Hankel-norm.
Note
Assumes the system is asymptotically stable.
Parameters
Returns
- norm
Hankel-norm.
- hinf_norm(mu=None, return_fpeak=False, ab13dd_equilibrate=False)[source]¶
Compute the H_infinity-norm.
Note
Assumes the system is asymptotically stable.
Parameters
- mu
- return_fpeak
Should the frequency at which the maximum is achieved should be returned.
- ab13dd_equilibrate
Should
slycot.ab13dd
use equilibration.
Returns
- norm
H_infinity-norm.
- fpeak
Frequency at which the maximum is achieved (if
return_fpeak
isTrue
).
- hsv(mu=None)[source]¶
Hankel singular values.
Note
Assumes the system is asymptotically stable.
Parameters
Returns
- sv
One-dimensional
NumPy array
of singular values.
- poles(mu=None)[source]¶
Compute system poles.
Note
Assumes the systems is small enough to use a dense eigenvalue solver.
Parameters
Returns
One-dimensional
NumPy array
of system poles.
- to_lti()[source]¶
Return a standard linear time-invariant system representation.
The representation
\[A = J - R,\qquad B = G - P,\qquad C = (G + P)^T,\qquad D = S - N,\qquad E = E\]is returned.
Returns
- lti
LTIModel
equivalent to the port-Hamiltonian model.
- to_matrices()[source]¶
Return operators as matrices.
Returns
- J
The
NumPy array
orSciPy spmatrix
J.- R
The
NumPy array
orSciPy spmatrix
R.- G
The
NumPy array
orSciPy spmatrix
G.- P
The
NumPy array
orSciPy spmatrix
P.- S
The
NumPy array
orSciPy spmatrix
S orNone
(if Cv is aZeroOperator
).- N
The
NumPy array
orSciPy spmatrix
N orNone
(if Cv is aZeroOperator
).- E
The
NumPy array
orSciPy spmatrix
E.
- class pymor.models.iosys.SecondOrderModel(M, E, K, B, Cp, Cv=None, D=None, sampling_time=0, solver_options=None, error_estimator=None, visualizer=None, name=None)[source]¶
Bases:
pymor.models.interface.Model
Class for linear second order systems.
This class describes input-output systems given by
\[\begin{split}M(\mu) \ddot{x}(t, \mu) + E(\mu) \dot{x}(t, \mu) + K(\mu) x(t, \mu) & = B(\mu) u(t), \\ y(t, \mu) & = C_p(\mu) x(t, \mu) + C_v(\mu) \dot{x}(t, \mu) + D(\mu) u(t),\end{split}\]if continuous-time, or
\[\begin{split}M(\mu) x(k + 2, \mu) + E(\mu) x(k + 1, \mu) + K(\mu) x(k, \mu) & = B(\mu) u(k), \\ y(k, \mu) & = C_p(\mu) x(k, \mu) + C_v(\mu) x(k + 1, \mu) + D(\mu) u(k),\end{split}\]if discrete-time, where \(M\), \(E\), \(K\), \(B\), \(C_p\), \(C_v\), and \(D\) are linear operators.
All methods related to the transfer function (e.g., frequency response calculation and Bode plots) are attached to the
transfer_function
attribute.Parameters
- M
The
Operator
M.- E
The
Operator
E.- K
The
Operator
K.- B
The
Operator
B.- Cp
The
Operator
Cp.- Cv
The
Operator
Cv orNone
(then Cv is assumed to be zero).- D
The
Operator
D orNone
(then D is assumed to be zero).- sampling_time
0
if the system is continuous-time, otherwise a positive number that denotes the sampling time (in seconds).- solver_options
The solver options to use to solve the Lyapunov equations.
- error_estimator
An error estimator for the problem. This can be any object with an
estimate_error(U, mu, model)
method. Iferror_estimator
is notNone
, anestimate_error(U, mu)
method is added to the model which will callerror_estimator.estimate_error(U, mu, self)
.- visualizer
A visualizer for the problem. This can be any object with a
visualize(U, model, ...)
method. Ifvisualizer
is notNone
, avisualize(U, *args, **kwargs)
method is added to the model which forwards its arguments to the visualizer’svisualize
method.- name
Name of the system.
Methods
Create
LTIModel
from matrices stored in separate files.Create a second order system from matrices.
Compute a second-order Gramian.
Compute the \(\mathcal{H}_2\)-norm.
Compute the Hankel-norm.
Compute the \(\mathcal{H}_\infty\)-norm.
Compute system poles.
Position singular values.
Position-velocity singular values.
Write operators to files as matrices.
Return a first order representation.
Return operators as matrices.
Velocity-position singular values.
Velocity singular values.
- classmethod from_files(M_file, E_file, K_file, B_file, Cp_file, Cv_file=None, D_file=None, sampling_time=0, state_id='STATE', solver_options=None, error_estimator=None, visualizer=None, name=None)[source]¶
Create
LTIModel
from matrices stored in separate files.Parameters
- M_file
The name of the file (with extension) containing A.
- E_file
The name of the file (with extension) containing E.
- K_file
The name of the file (with extension) containing K.
- B_file
The name of the file (with extension) containing B.
- Cp_file
The name of the file (with extension) containing Cp.
- Cv_file
None
or the name of the file (with extension) containing Cv.- D_file
None
or the name of the file (with extension) containing D.- sampling_time
0
if the system is continuous-time, otherwise a positive number that denotes the sampling time (in seconds).- state_id
Id of the state space.
- solver_options
The solver options to use to solve the Lyapunov equations.
- error_estimator
An error estimator for the problem. This can be any object with an
estimate_error(U, mu, model)
method. Iferror_estimator
is notNone
, anestimate_error(U, mu)
method is added to the model which will callerror_estimator.estimate_error(U, mu, self)
.- visualizer
A visualizer for the problem. This can be any object with a
visualize(U, model, ...)
method. Ifvisualizer
is notNone
, avisualize(U, *args, **kwargs)
method is added to the model which forwards its arguments to the visualizer’svisualize
method.- name
Name of the system.
Returns
- some
The
SecondOrderModel
with operators M, E, K, B, Cp, Cv, and D.
- classmethod from_matrices(M, E, K, B, Cp, Cv=None, D=None, sampling_time=0, state_id='STATE', solver_options=None, error_estimator=None, visualizer=None, name=None)[source]¶
Create a second order system from matrices.
Parameters
- M
The
NumPy array
orSciPy spmatrix
M.- E
The
NumPy array
orSciPy spmatrix
E.- K
The
NumPy array
orSciPy spmatrix
K.- B
The
NumPy array
orSciPy spmatrix
B.- Cp
The
NumPy array
orSciPy spmatrix
Cp.- Cv
The
NumPy array
orSciPy spmatrix
Cv orNone
(then Cv is assumed to be zero).- D
The
NumPy array
orSciPy spmatrix
D orNone
(then D is assumed to be zero).- sampling_time
0
if the system is continuous-time, otherwise a positive number that denotes the sampling time (in seconds).- solver_options
The solver options to use to solve the Lyapunov equations.
- error_estimator
An error estimator for the problem. This can be any object with an
estimate_error(U, mu, model)
method. Iferror_estimator
is notNone
, anestimate_error(U, mu)
method is added to the model which will callerror_estimator.estimate_error(U, mu, self)
.- visualizer
A visualizer for the problem. This can be any object with a
visualize(U, model, ...)
method. Ifvisualizer
is notNone
, avisualize(U, *args, **kwargs)
method is added to the model which forwards its arguments to the visualizer’svisualize
method.- name
Name of the system.
Returns
- lti
The SecondOrderModel with operators M, E, K, B, Cp, Cv, and D.
- gramian(typ, mu=None)[source]¶
Compute a second-order Gramian.
Parameters
- typ
The type of the Gramian:
'pc_lrcf'
: low-rank Cholesky factor of the position controllability Gramian,'vc_lrcf'
: low-rank Cholesky factor of the velocity controllability Gramian,'po_lrcf'
: low-rank Cholesky factor of the position observability Gramian,'vo_lrcf'
: low-rank Cholesky factor of the velocity observability Gramian,'pc_dense'
: dense position controllability Gramian,'vc_dense'
: dense velocity controllability Gramian,'po_dense'
: dense position observability Gramian,'vo_dense'
: dense velocity observability Gramian.
Note
For
'*_lrcf'
types, the method assumes the system is asymptotically stable. For'*_dense'
types, the method assumes that the underlying Lyapunov equation has a unique solution, i.e. no pair of system poles adds to zero in the continuous-time case and no pair of system poles multiplies to one in the discrete-time case.- mu
Returns
If typ is
'pc_lrcf'
,'vc_lrcf'
,'po_lrcf'
or'vo_lrcf'
, then the Gramian factor as aVectorArray
fromself.M.source
. If typ is'pc_dense'
,'vc_dense'
,'po_dense'
or'vo_dense'
, then the Gramian as aNumPy array
.
- h2_norm(mu=None)[source]¶
Compute the \(\mathcal{H}_2\)-norm.
Note
Assumes the system is asymptotically stable.
Parameters
Returns
- norm
\(\mathcal{H}_2\)-norm.
- hankel_norm(mu=None)[source]¶
Compute the Hankel-norm.
Note
Assumes the system is asymptotically stable.
Parameters
Returns
- norm
Hankel-norm.
- hinf_norm(mu=None, return_fpeak=False, ab13dd_equilibrate=False)[source]¶
Compute the \(\mathcal{H}_\infty\)-norm.
Note
Assumes the system is asymptotically stable.
Parameters
- mu
- return_fpeak
Should the frequency at which the maximum is achieved should be returned.
- ab13dd_equilibrate
Should
slycot.ab13dd
use equilibration.
Returns
- norm
\(\mathcal{H}_\infty\).
- fpeak
Frequency at which the maximum is achieved (if
return_fpeak
isTrue
).
- poles(mu=None)[source]¶
Compute system poles.
Note
Assumes the systems is small enough to use a dense eigenvalue solver.
Parameters
Returns
One-dimensional
NumPy array
of system poles.
- psv(mu=None)[source]¶
Position singular values.
Note
Assumes the system is asymptotically stable.
Parameters
Returns
One-dimensional
NumPy array
of singular values.
- pvsv(mu=None)[source]¶
Position-velocity singular values.
Note
Assumes the system is asymptotically stable.
Parameters
Returns
One-dimensional
NumPy array
of singular values.
- to_files(M_file, E_file, K_file, B_file, Cp_file, Cv_file=None, D_file=None)[source]¶
Write operators to files as matrices.
Parameters
- M_file
The name of the file (with extension) containing M.
- E_file
The name of the file (with extension) containing E.
- K_file
The name of the file (with extension) containing K.
- B_file
The name of the file (with extension) containing B.
- Cp_file
The name of the file (with extension) containing Cp.
- Cv_file
The name of the file (with extension) containing Cv or
None
if D is aZeroOperator
.- D_file
The name of the file (with extension) containing D or
None
if D is aZeroOperator
.
- to_lti()[source]¶
Return a first order representation.
The first order representation
\[\begin{split}\begin{bmatrix} I & 0 \\ 0 & M \end{bmatrix} \frac{\mathrm{d}}{\mathrm{d}t}\! \begin{bmatrix} x(t) \\ \dot{x}(t) \end{bmatrix} & = \begin{bmatrix} 0 & I \\ -K & -E \end{bmatrix} \begin{bmatrix} x(t) \\ \dot{x}(t) \end{bmatrix} + \begin{bmatrix} 0 \\ B \end{bmatrix} u(t), \\ y(t) & = \begin{bmatrix} C_p & C_v \end{bmatrix} \begin{bmatrix} x(t) \\ \dot{x}(t) \end{bmatrix} + D u(t)\end{split}\]is returned.
Returns
- lti
LTIModel
equivalent to the second-order model.
- to_matrices()[source]¶
Return operators as matrices.
Returns
- M
The
NumPy array
orSciPy spmatrix
M.- E
The
NumPy array
orSciPy spmatrix
E.- K
The
NumPy array
orSciPy spmatrix
K.- B
The
NumPy array
orSciPy spmatrix
B.- Cp
The
NumPy array
orSciPy spmatrix
Cp.- Cv
The
NumPy array
orSciPy spmatrix
Cv orNone
(if Cv is aZeroOperator
).- D
The
NumPy array
orSciPy spmatrix
D orNone
(if D is aZeroOperator
).
- vpsv(mu=None)[source]¶
Velocity-position singular values.
Note
Assumes the system is asymptotically stable.
Parameters
Returns
One-dimensional
NumPy array
of singular values.
- vsv(mu=None)[source]¶
Velocity singular values.
Note
Assumes the system is asymptotically stable.
Parameters
Returns
One-dimensional
NumPy array
of singular values.