POD¶
-
class
modred.pod.
PODHandles
(inner_product=None, get_array=<function load_array_text>, put_array=<function save_array_text>, max_vecs_per_node=None, verbosity=1)[source]¶ Proper Orthogonal Decomposition implemented for large datasets.
- Kwargs:
inner_product
: Function that computes inner product of two vector objects.put_array
: Function to put an array out of modred, e.g., write it to file.get_array
: Function to get an array into modred, e.g., load it from 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.
Computes POD modes from vector objects (or handles). Uses
vectorspace.VectorSpaceHandles
for low level functions.Usage:
myPOD = POD(inner_product=my_inner_product) myPOD.compute_decomp(vec_handles) myPOD.compute_modes(range(10), modes)
See also
compute_POD_arrays_snaps_method()
,compute_POD_arrays_direct_method()
, andvectors
.-
compute_decomp
(vec_handles, atol=1e-13, rtol=None)[source]¶ Computes correlation array and its eigendecomposition.
- Args:
vec_handles
: List of handles for vector objects.- Kwargs:
atol
: Level below which eigenvalues of correlation array are truncated.rtol
: Maximum relative difference between largest and smallest eigenvalues of correlation array. Smaller ones are truncated.- Returns:
eigvals
: 1D array of eigenvalues of correlation array.eigvecs
: Array whose columns are eigenvectors of correlation array.
-
compute_eigendecomp
(atol=1e-13, rtol=None)[source]¶ Computes eigendecomposition of correlation array.
- Kwargs:
atol
: Level below which eigenvalues of correlation array are truncated.rtol
: Maximum relative difference between largest and smallest eigenvalues of correlation array. Smaller ones are truncated.
Useful if you already have the correlation array and to want to avoid recomputing it.
Usage:
POD.correlation_array = pre_existing_correlation_array POD.compute_eigendecomp() POD.compute_modes(range(10), mode_handles, vec_handles=vec_handles)
-
compute_modes
(mode_indices, mode_handles, vec_handles=None)[source]¶ Computes POD modes and calls
put
on them using mode handles.- Args:
mode_indices
: List of indices describing which modes to compute, e.g.range(10)
or[3, 0, 5]
.mode_handles
: List of handles for modes to compute.- Kwargs:
vec_handles
: List of handles for vector objects. Optional if when callingcompute_decomp()
.
-
compute_proj_coeffs
()[source]¶ Computes orthogonal projection of vector objects onto POD modes.
- Returns:
proj_coeffs
: Array of projection coefficients for vector objects, expressed as a linear combination of POD modes. Columns correspond to vector objects, rows correspond to POD modes.
-
get_correlation_array
(src)[source]¶ Gets the correlation array from source (memory or file).
- Args:
src
: Source from which to retrieve correlation array.
-
get_decomp
(eigvals_src, eigvecs_src)[source]¶ Gets the decomposition arrays from sources (memory or file).
- Args:
eigvals_src
: Source from which to retrieve eigenvalues of correlation array.eigvecs_src
: Source from which to retrieve eigenvectors of correlation array.
-
get_proj_coeffs
(src)[source]¶ Gets projection coefficients from source (memory or file).
- Args:
src
: Source from which to retrieve projection coefficients.
-
modred.pod.
compute_POD_arrays_direct_method
(vecs, mode_indices=None, inner_product_weights=None, atol=1e-13, rtol=None)[source]¶ Computes POD modes using data stored in an array, using direct method.
- Args:
vecs
: Array whose columns are data vectors ().- Kwargs:
mode_indices
: List of indices describing which modes to compute. Examples arerange(10)
or[3, 0, 6, 8]
. If no mode indices are specified, then all modes will be computed.inner_product_weights
: 1D or 2D array of inner product weights. Corresponds to in inner product .atol
: Level below which eigenvalues of correlation array are truncated.rtol
: Maximum relative difference between largest and smallest eigenvalues of correlation array. Smaller ones are truncated.return_all
: Return more objects; see below. Default is false.- Returns:
res
: Results of POD computation, stored in a namedtuple with the following attributes:eigvals
: 1D array of eigenvalues of correlation array ().modes
: Array whose columns are POD modes.proj_coeffs
: Array of projection coefficients for vector objects, expressed as a linear combination of POD modes. Columns correspond to vector objects, rows correspond to POD modes.eigvecs
: Array wholse columns are eigenvectors of correlation array ().
Attributes can be accessed using calls like
res.modes
. To see all available attributes, useprint(res)
.
The algorithm is
- SVD
- Modes are
where , , , , correspond to
vecs
,inner_product_weights
,eigvals**0.5
, andeigvecs
, respectively.Since this method does not square the vectors and singular values, it is more accurate than taking the eigendecomposition of the correlation array , as in the method of snapshots (
compute_POD_arrays_snaps_method()
). However, this method is slower when has more rows than columns, i.e. there are fewer vectors than elements in each vector.
-
modred.pod.
compute_POD_arrays_snaps_method
(vecs, mode_indices=None, inner_product_weights=None, atol=1e-13, rtol=None)[source]¶ Computes POD modes using data stored in an array, using the method of snapshots.
- Args:
vecs
: Array whose columns are data vectors ().- Kwargs:
mode_indices
: List of indices describing which modes to compute. Examples arerange(10)
or[3, 0, 6, 8]
. If no mode indices are specified, then all modes will be computed.inner_product_weights
: 1D or 2D array of inner product weights. Corresponds to in inner product .atol
: Level below which eigenvalues of correlation array are truncated.rtol
: Maximum relative difference between largest and smallest eigenvalues of correlation array. Smaller ones are truncated.- Returns:
res
: Results of POD computation, stored in a namedtuple with the following attributes:eigvals
: 1D array of eigenvalues of correlation array ().modes
: Array whose columns are POD modes.proj_coeffs
: Array of projection coefficients for vector objects, expressed as a linear combination of POD modes. Columns correspond to vector objects, rows correspond to POD modes.eigvecs
: Array wholse columns are eigenvectors of correlation array ().correlation_array
: Correlation array ().
Attributes can be accessed using calls like
res.modes
. To see all available attributes, useprint(res)
.
The algorithm is
- Solve eigenvalue problem
- Coefficient array
- Modes are
where , , , and correspond to
vecs
,inner_product_weights
,correlation_array
, andbuild_coeffs
, respectively.Since this method “squares” the vector array and thus its singular values, it is slightly less accurate than taking the SVD of directly, as in
compute_POD_arrays_direct_method()
. However, this method is faster when has more rows than columns, i.e. there are more elements in each vector than there are vectors.