LTIGalerkinProjection¶
Module for Galerkin projection of LTI systems.
-
class
modred.ltigalerkinproj.
LTIGalerkinProjectionArrays
(basis_vecs, adjoint_basis_vecs=None, is_basis_orthonormal=False, inner_product_weights=None, put_array=<function save_array_text>)[source]¶ Computes A, B, and C arrays of reduced-order model (ROM) of a linear time-invariant (LTI) system, from data stored in arrays.
- Args:
basis_vecs
: Array whose columns are basis vectors.- Kwargs:
adjoint_basis_vec_handles
: Array whose columns are adjoint basis vectors. If not given, then assumed to be the same asbasis_vecs
.is_basis_orthonormal
: Boolean for biorthonormality of the basis vecs.True
if the basis and adjoint basis vectors are biorthonormal. Default isFalse
.inner_product_weights
: 1D or 2D array of inner product weights. Corresponds to in inner product .put_array
: Function to put an array out of modred, e.g., write it to file.
This class projects discrete- or continuous-time dynamics onto a set of basis vectors, producing reduced-order models. For example, the basis vectors could be POD, BPOD, or DMD modes. It uses
vectorspace.VectorSpaceArrays
for low level functions.Usage:
LTI_proj = LTIGalerkinProjectionArray( basis_vecs, adjoint_basis_vecs=adjoint_basis_vecs, is_basis_orthonormal=True) A, B, C = LTI_proj.compute_model( A_on_basis_vecs, B_on_standard_basis_array, C_on_basis_vecs)
-
compute_model
(A_on_basis_vecs, B_on_standard_basis_array, C_on_basis_vecs)[source]¶ Computes the continous- or discrete-time reduced A, B, and C arrays.
- Args:
A_on_basis_vecs
: Array whose columns contain , i.e., the results of acting on individual basis vectors. can represent either discrete- or continuous-time dynamics. For continuous time systems, see alsocompute_derivs_arrays()
.B_on_standard_basis_array
: Array whose jth column contains , i.e., the action of the array on the jth standard basis vector (a vector with a 1 at the jth index and zeros elsewhere).C_on_basis_vecs
: Array whose columns contain , i.e., the results of acting on individual basis vectors.- Returns:
A_reduced
: Reduced-order A array.B_reduced
: Reduced-order B array.C_reduced
: Reduced-order C array.
-
reduce_A
(A_on_basis_vecs)[source]¶ Computes the continous- or discrete-time reduced A array.
- Args:
A_on_basis_vecs
: Array whose columns contain , i.e., the results of acting on individual basis vectors. can represent either discrete- or continuous-time dynamics. For continuous time systems, see alsocompute_derivs_arrays()
.- Returns:
A_reduced
: Reduced-order A array.
-
reduce_B
(B_on_standard_basis_array)[source]¶ Computes the continous- or discrete-time reduced B array.
- Args:
B_on_standard_basis_array
: Array whose jth column contains , i.e., the action of the array on the jth standard basis vector (a vector with a 1 at the jth index and zeros elsewhere).- Returns:
B_reduced
: Reduced-order B array.
Note that if you found the modes via sampling initial condition responses to B (and C), then your snapshots may be missing a factor of 1 / dt (where dt is the first simulation time step). This can be remedied by multiplying
B_reduced
by dt.
-
class
modred.ltigalerkinproj.
LTIGalerkinProjectionHandles
(inner_product, basis_vec_handles, adjoint_basis_vec_handles=None, is_basis_orthonormal=False, put_array=<function save_array_text>, verbosity=1, max_vecs_per_node=10000)[source]¶ Computes A, B, and C arrays of reduced-order model (ROM) of a linear time-invariant (LTI) system, using vector handles.
- Args:
inner_product
: Function that computes inner product of two vector objects.basis_vec_handles
: List of handles for basis vector objects.- Kwargs:
adjoint_basis_vec_handles
: List of handles for adjoint basis vector objects. If not given, then assumed to be the same asbasis_vec_handles
.is_basis_orthonormal
: Boolean for biorthonormality of the basis vecs.True
if the basis and adjoint basis vectors are biorthonormal. Default isFalse
.put_array
: Function to put an array out of modred, e.g., write it to file.max_vecs_per_node
: Maximum number of vectors that can be stored in memory, per node.verbosity
: 1 prints progress and warnings, 0 prints almost nothing.
This class projects discrete- or continuous-time dynamics onto a set of basis vectors, producing reduced-order models. For example, the basis vectors could be POD, BPOD, or DMD modes.
Usage:
LTI_proj = LTIGalerkinProjectionHandles( inner_product, basis_vec_handles, adjoint_basis_vec_handles=adjoint_basis_vec_handles, is_basis_orthonormal=True) A, B, C = LTI_proj.compute_model( A_on_basis_vec_handles, B_on_standard_basis_handles, C_on_basis_vecs)
-
compute_model
(A_on_basis_vec_handles, B_on_standard_basis_handles, C_on_basis_vecs)[source]¶ Computes the continous- or discrete-time reduced A, B, and C arrays.
- Args:
A_on_basis_vec_handles
: List of handles for the results of acting on individual basis vectors. can represent either discrete- or continuous-time dynamics. For continuous time systems, see alsocompute_derivs_handles()
.B_on_standard_basis_handles
: List of handles for the results of acting on standard basis vectors (vectors with a 1 at the jth index and zeros elsewhere).C_on_basis_vecs
: List of 1D arrays, each of which is the result of acting on an individual basis vector.- Returns:
A_reduced
: Reduced-order A array.B_reduced
: Reduced-order B array.C_reduced
: Reduced-order C array.
-
reduce_A
(A_on_basis_vec_handles)[source]¶ Computes the continous- or discrete-time reduced A array.
- Args:
A_on_basis_vec_handles
: List of handles for the results of acting on individual basis vectors. can represent either discrete- or continuous-time dynamics. For continuous time systems, see alsocompute_derivs_handles()
.- Returns:
A_reduced
: Reduced-order A array.
-
reduce_B
(B_on_standard_basis_handles)[source]¶ Computes the continous- or discrete-time reduced B array.
- Args:
B_on_standard_basis_handles
: List of handles for the results of acting on standard basis vectors (vectors with a 1 at the jth index and zeros elsewhere).- Returns:
B_reduced
: Reduced-order B array.
Note that if you found the modes via sampling initial condition responses to B (and C), then your snapshots may be missing a factor of 1/dt (where dt is the first simulation time step). This can be remedied by multiplying
B_reduced
by dt.
-
modred.ltigalerkinproj.
compute_derivs_arrays
(vecs, adv_vecs, dt)[source]¶ Computes 1st-order time derivatives using data stored in arrays.
- Args:
vecs
: Array whose columns are data vectors.adv_vecs
: Array whose columns are data vectors advanced in time.dt
: Time step, corresponding to time advanced betweenvec_handles
andadv_vec_handles
- Returns:
deriv_vecs
: Array whose columns are time derivatives of data vectors.
Computes d(
vec
)/dt = (vec
(t=dt) -vec
(t=0) ) / dt.
-
modred.ltigalerkinproj.
compute_derivs_handles
(vec_handles, adv_vec_handles, deriv_vec_handles, dt)[source]¶ Computes 1st-order time derivatives of vector objects, using handles.
- Args:
vec_handles
: List of handles for vector objects.adv_vec_handles
: List of handles for vector objects advanced in time.deriv_vec_handles
: List of handles for time derivatives of vector objects.dt
: Time step, corresponding to time advanced betweenvec_handles
andadv_vec_handles
Computes d(
vec
)/dt = (vec
(t=dt) -vec
(t=0) ) / dt.