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.ImpulseFunction(dim_input, component, sampling_time)[source]¶
Bases:
pymor.analyticalproblems.functions.Function
Interface for
Parameter
dependent analytical functions.Every
Function
is a map of the formf(μ): Ω ⊆ R^d -> R^(shape_range)
The returned values are
NumPy arrays
of arbitrary (but fixed) shape. Note that NumPy distinguishes between one-dimensional arrays of length 1 (with shape(1,)
) and zero-dimensional scalar arrays (with shape()
). In pyMOR, we usually expect scalar-valued functions to haveshape_range == ()
.While the function might raise an error if it is evaluated for an argument not in the domain Ω, the exact behavior is left undefined.
Functions are vectorized in the sense, that if
x.ndim == k
, thenf(x, μ)[i0, i1, ..., i(k-2)] == f(x[i0, i1, ..., i(k-2)], μ).
In particular,
f(x, μ).shape == x.shape[:-1] + shape_range
.Methods
Evaluate the function for given argument
x
andparameter values
mu
.- evaluate(x, mu=None)[source]¶
Evaluate the function for given argument
x
andparameter values
mu
.
- 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.
Return A, B, C, D, and E operators as matrices.
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
. SeeLTIModel
.- 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
. SeeLTIModel
.- 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
. SeeLTIModel
.- 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
. SeeLTIModel
.- 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.'pr_c_lrcf'
: low-rank Cholesky factor of the “controllability” positive real Gramian,'pr_o_lrcf'
: low-rank Cholesky factor of the “observability” positive 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. Additionally, for'pr_c_lrcf'
and'pr_o_lrcf'
, it is assumed thatD + D^T
is invertible.- 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, tol=1e-10)[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{L}_\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.- tol
Tolerance in norm computation.
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
- output
Impulse response as a 3D
NumPy array
whereoutput.shape[0]
is the number of time steps,output.shape[1]
is the number of outputs, andoutput.shape[2]
is the number of inputs.- solution
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, tol=1e-10)[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.- tol
Tolerance in norm computation.
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
- output
Step response as a 3D
NumPy array
whereoutput.shape[0]
is the number of time steps,output.shape[1]
is the number of outputs, andoutput.shape[2]
is the number of inputs.- solution
The tuple of solution
VectorArrays
for every input. Returned only whenreturn_solution
isTrue
.
- to_abcde_files(files_basename, mu=None)[source]¶
Save operators as matrices to .[ABCDE] files in Matrix Market format.
Parameters
- files_basename
The basename of files containing the operators.
- mu
The
parameter values
for which to write the operators to files.
- to_abcde_matrices(format=None, mu=None)[source]¶
Return A, B, C, D, and E operators as matrices.
Parameters
- format
Format of the resulting matrices:
NumPy array
if ‘dense’, otherwise the appropriateSciPy spmatrix
. IfNone
, a choice between dense and sparse format is automatically made.- mu
The
parameter values
for which to convert the operators.
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
).
- 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, mu=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
.- mu
The
parameter values
for which to write the operators to files.
- to_mat_file(file_name, mu=None)[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).
- mu
The
parameter values
for which to write the operators to files.
- to_matrices(format=None, mu=None)[source]¶
Return operators as matrices.
Parameters
- format
Format of the resulting matrices:
NumPy array
if ‘dense’, otherwise the appropriateSciPy spmatrix
. IfNone
, a choice between dense and sparse format is automatically made.- mu
The
parameter values
for which to convert the operators.
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, Q=None, solver_options=None, error_estimator=None, visualizer=None, name=None)[source]¶
Bases:
LTIModel
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)) Q(\mu) x(t, \mu) + (G(\mu) - P(\mu)) u(t), \\ y(t, \mu) & = (G(\mu) + P(\mu))^T Q(\mu) x(t, \mu) + (S(\mu) - N(\mu)) u(t),\end{split}\]where \(H(\mu) = Q(\mu)^T E(\mu)\),
\[\begin{split}\Gamma(\mu) = \begin{bmatrix} J(\mu) & G(\mu) \\ -G(\mu)^T & N(\mu) \end{bmatrix}, \text{ and } \mathcal{W}(\mu) = \begin{bmatrix} R(\mu) & P(\mu) \\ P(\mu)^T & S(\mu) \end{bmatrix}\end{split}\]satisfy \(H(\mu) = H(\mu)^T \succ 0\), \(\Gamma(\mu)^T = -\Gamma(\mu)\), and \(\mathcal{W}(\mu) = \mathcal{W}(\mu)^T \succcurlyeq 0\).
A dynamical system of this form, together with a given quadratic (energy) function \(\mathcal{H}(x, \mu) = \tfrac{1}{2} x^T H(\mu) x\), typically called Hamiltonian, is called a port-Hamiltonian system.
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).- Q
The
Operator
Q orNone
(then Q 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.Convert a passive
LTIModel
to aPHLTIModel
.Convert the
PHLTIModel
into its Berlin form.Return operators as matrices.
- classmethod from_matrices(J, R, G, P=None, S=None, N=None, E=None, Q=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).- Q
The
NumPy array
orSciPy spmatrix
Q orNone
(then Q 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.
- classmethod from_passive_LTIModel(model)[source]¶
Convert a passive
LTIModel
to aPHLTIModel
.Note
The method uses dense computations and converts
model
to dense matrices.Parameters
- model
The passive
LTIModel
to convert.- generalized
If
True
, the resultingPHLTIModel
will have \(Q=I\).
- to_berlin_form()[source]¶
Convert the
PHLTIModel
into its Berlin form.Returns a
PHLTIModel
with \(Q=I\), by left multiplication with \(Q^T\).Returns
- model
PHLTIModel
with \(Q=I\).
- to_matrices(format=None, mu=None)[source]¶
Return operators as matrices.
Parameters
- format
Format of the resulting matrices:
NumPy array
if ‘dense’, otherwise the appropriateSciPy spmatrix
. IfNone
, a choice between dense and sparse format is automatically made.- mu
The
parameter values
for which to convert the operators.
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 orNone
(if P is aZeroOperator
).- S
The
NumPy array
orSciPy spmatrix
S orNone
(if S is aZeroOperator
).- N
The
NumPy array
orSciPy spmatrix
N orNone
(if N is aZeroOperator
).- E
The
NumPy array
orSciPy spmatrix
E orNone
(if E is anIdentityOperator
).- Q
The
NumPy array
orSciPy spmatrix
Q orNone
(if Q is anIdentityOperator
).
- 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, tol=1e-10)[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.- tol
Tolerance in norm computation.
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, mu=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
.- mu
The
parameter values
for which to write the operators to files.
- 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(mu=None)[source]¶
Return operators as matrices.
Parameters
- mu
The
parameter values
for which to convert the operators.
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.
- class pymor.models.iosys.StepFunction(dim_input, component, sampling_time)[source]¶
Bases:
pymor.analyticalproblems.functions.Function
Interface for
Parameter
dependent analytical functions.Every
Function
is a map of the formf(μ): Ω ⊆ R^d -> R^(shape_range)
The returned values are
NumPy arrays
of arbitrary (but fixed) shape. Note that NumPy distinguishes between one-dimensional arrays of length 1 (with shape(1,)
) and zero-dimensional scalar arrays (with shape()
). In pyMOR, we usually expect scalar-valued functions to haveshape_range == ()
.While the function might raise an error if it is evaluated for an argument not in the domain Ω, the exact behavior is left undefined.
Functions are vectorized in the sense, that if
x.ndim == k
, thenf(x, μ)[i0, i1, ..., i(k-2)] == f(x[i0, i1, ..., i(k-2)], μ).
In particular,
f(x, μ).shape == x.shape[:-1] + shape_range
.Methods
Evaluate the function for given argument
x
andparameter values
mu
.- evaluate(x, mu=None)[source]¶
Evaluate the function for given argument
x
andparameter values
mu
.