libROM  v1.0
Data-driven physical simulation library
Public Member Functions | Protected Member Functions | Protected Attributes | Static Protected Attributes | List of all members
CAROM::IncrementalSVD Class Referenceabstract

#include <IncrementalSVD.h>

Inheritance diagram for CAROM::IncrementalSVD:
CAROM::SVD CAROM::IncrementalSVDBrand CAROM::IncrementalSVDFastUpdate CAROM::IncrementalSVDStandard

Public Member Functions

 IncrementalSVD (Options options, const std::string &basis_file_name)
 Constructor. More...
 
virtual ~IncrementalSVD ()
 Destructor.
 
virtual bool takeSample (double *u_in, bool add_without_increase=false)
 Sample new state, u_in, at the given time. More...
 
virtual const MatrixgetSpatialBasis ()
 Returns the basis vectors for the current time interval as a Matrix. More...
 
virtual const MatrixgetTemporalBasis ()
 Returns the temporal basis vectors for the current time interval as a Matrix. More...
 
virtual const VectorgetSingularValues ()
 Returns the singular values for the current time interval. More...
 
virtual const MatrixgetSnapshotMatrix ()
 Returns the snapshot matrix for the current time interval. More...
 
- Public Member Functions inherited from CAROM::SVD
 SVD (Options options)
 Constructor. More...
 
 ~SVD ()
 
int getDim () const
 Returns the dimension of the system on this processor. More...
 
int getNumSamples () const
 Get the number of samples taken.
 
int getMaxNumSamples () const
 Get the maximum number of samples that can be taken. SVD class will return an error if the number of samples exceeds the maximum.
 

Protected Member Functions

virtual void buildInitialSVD (double *u)=0
 Constructs the first SVD. More...
 
virtual bool buildIncrementalSVD (double *u, bool add_without_increase=false)
 Adds the new sampled the state vector, u, to the system. More...
 
virtual void computeBasis ()=0
 Computes the current basis vectors.
 
void constructQ (double *&Q, const Vector *l, double k)
 Construct the matrix Q whose SVD is needed. More...
 
bool svd (double *A, Matrix *&U, Matrix *&S, Matrix *&V)
 Given a matrix, A, returns 2 of the 3 components of its singular value decomposition. The right singular vectors are not needed and therefore not returned. More...
 
virtual void addLinearlyDependentSample (const Matrix *A, const Matrix *W, const Matrix *sigma)=0
 Add a linearly dependent sample to the SVD. More...
 
virtual void addNewSample (const Vector *j, const Matrix *A, const Matrix *W, Matrix *sigma)=0
 Add a new, unique sample to the SVD. More...
 
int numSamples ()
 The number of samples stored. More...
 
double checkOrthogonality (const Matrix *m)
 Computes and returns the orthogonality of m. More...
 
- Protected Member Functions inherited from CAROM::SVD
bool isFirstSample () const
 Returns true if the next sample will result in a new time interval. More...
 

Protected Attributes

double d_linearity_tol
 Tolerance to determine whether or not a sample is linearly dependent.
 
bool d_skip_linearly_dependent
 Whether to skip linearly dependent samples.
 
int d_max_basis_dimension
 The maximum basis dimension.
 
int d_size
 The total number of processors.
 
int d_rank
 The rank of the processor owning this object.
 
std::vector< int > d_proc_dims
 The dimension of the system on each processor.
 
long int d_total_dim
 The total dimension of the system.
 
bool d_save_state
 If true the state of the SVD will be written to disk when the object is deleted. If there are multiple time intervals then the state will not be saved as restoring such a state makes no sense.
 
bool d_update_right_SV
 Whether to update the right singular vectors.
 
Databased_state_database
 Pointer to the database that will hold saved state data if the state is to be saved.
 
std::string d_state_file_name
 The name of file to which state is saved or restored from.
 
- Protected Attributes inherited from CAROM::SVD
const int d_dim
 Dimension of the system.
 
int d_num_samples
 Number of samples stored for the current time interval.
 
int d_num_rows_of_W
 Number of rows in right singular matrix.
 
const int d_max_num_samples
 The maximum number of samples.
 
Matrixd_basis
 The globalized basis vectors for the current time interval. More...
 
Matrixd_basis_right
 The globalized right basis vectors for the current time interval. More...
 
Matrixd_U
 The matrix U which is large. More...
 
Matrixd_W
 The matrix U which is large. More...
 
Vectord_S
 The vector S which is small. More...
 
Matrixd_snapshots
 The globalized snapshot vectors for the current time interval. More...
 
bool d_debug_algorithm
 Flag to indicate if results of algorithm should be printed for debugging purposes.
 

Static Protected Attributes

static const int COMMUNICATE_U = 666
 MPI message tag.
 

Detailed Description

Abstract class IncrementalSVD defines the internal API of the incremental SVD algorithm.

Definition at line 27 of file IncrementalSVD.h.

Constructor & Destructor Documentation

◆ IncrementalSVD()

CAROM::IncrementalSVD::IncrementalSVD ( Options  options,
const std::string &  basis_file_name 
)

Constructor.

Parameters
[in]optionsThe struct containing the options for this abstract SVD class.
[in]basis_file_nameThe base part of the name of the file containing the basis vectors. Each process will append its process ID to this base name.
See also
Options

Definition at line 38 of file IncrementalSVD.cpp.

Member Function Documentation

◆ addLinearlyDependentSample()

virtual void CAROM::IncrementalSVD::addLinearlyDependentSample ( const Matrix A,
const Matrix W,
const Matrix sigma 
)
protectedpure virtual

Add a linearly dependent sample to the SVD.

Precondition
A != 0
sigma != 0
Parameters
[in]AThe left singular vectors.
[in]WThe right singular vectors.
[in]sigmaThe singular values.

◆ addNewSample()

virtual void CAROM::IncrementalSVD::addNewSample ( const Vector j,
const Matrix A,
const Matrix W,
Matrix sigma 
)
protectedpure virtual

Add a new, unique sample to the SVD.

Precondition
j != 0
A != 0
sigma != 0
Parameters
[in]jThe new column of d_U.
[in]AThe left singular vectors.
[in]WThe right singular vectors.
[in]sigmaThe singular values.

◆ buildIncrementalSVD()

bool CAROM::IncrementalSVD::buildIncrementalSVD ( double *  u,
bool  add_without_increase = false 
)
protectedvirtual

Adds the new sampled the state vector, u, to the system.

Precondition
u != 0
Parameters
[in]uThe new state.
[in]add_without_increaseIf true, addLinearlyDependent is invoked.
Returns
True if building the incremental SVD was successful.

Definition at line 281 of file IncrementalSVD.cpp.

◆ buildInitialSVD()

virtual void CAROM::IncrementalSVD::buildInitialSVD ( double *  u)
protectedpure virtual

Constructs the first SVD.

Precondition
u != 0
Parameters
[in]uThe first state.

◆ checkOrthogonality()

double CAROM::IncrementalSVD::checkOrthogonality ( const Matrix m)
protected

Computes and returns the orthogonality of m.

Precondition
m != 0
Parameters
[in]mThe matrix to check.
Returns
The orthogonality of m.

Definition at line 511 of file IncrementalSVD.cpp.

◆ constructQ()

void CAROM::IncrementalSVD::constructQ ( double *&  Q,
const Vector l,
double  k 
)
protected

Construct the matrix Q whose SVD is needed.

Precondition
l != 0
l.dim() == numSamples()
Parameters
[out]QThe matrix to be constructed. [d_S,l; 0,k]
[in]lThe last column of Q.
[in]kThe lower right element of Q.

Definition at line 393 of file IncrementalSVD.cpp.

◆ getSingularValues()

const Vector * CAROM::IncrementalSVD::getSingularValues ( )
virtual

Returns the singular values for the current time interval.

Returns
The singular values for the current time interval.

Implements CAROM::SVD.

Definition at line 266 of file IncrementalSVD.cpp.

◆ getSnapshotMatrix()

const Matrix * CAROM::IncrementalSVD::getSnapshotMatrix ( )
virtual

Returns the snapshot matrix for the current time interval.

Returns
The snapshot matrix for the current time interval.

Implements CAROM::SVD.

Definition at line 273 of file IncrementalSVD.cpp.

◆ getSpatialBasis()

const Matrix * CAROM::IncrementalSVD::getSpatialBasis ( )
virtual

Returns the basis vectors for the current time interval as a Matrix.

Returns
The basis vectors for the current time interval.

Implements CAROM::SVD.

Reimplemented in CAROM::IncrementalSVDBrand.

Definition at line 252 of file IncrementalSVD.cpp.

◆ getTemporalBasis()

const Matrix * CAROM::IncrementalSVD::getTemporalBasis ( )
virtual

Returns the temporal basis vectors for the current time interval as a Matrix.

Returns
The temporal basis vectors for the current time interval.

Implements CAROM::SVD.

Reimplemented in CAROM::IncrementalSVDBrand.

Definition at line 259 of file IncrementalSVD.cpp.

◆ numSamples()

int CAROM::IncrementalSVD::numSamples ( )
inlineprotected

The number of samples stored.

Returns
The number of samples stored.

Definition at line 220 of file IncrementalSVD.h.

◆ svd()

bool CAROM::IncrementalSVD::svd ( double *  A,
Matrix *&  U,
Matrix *&  S,
Matrix *&  V 
)
protected

Given a matrix, A, returns 2 of the 3 components of its singular value decomposition. The right singular vectors are not needed and therefore not returned.

Precondition
A != 0
Parameters
[in]AThe matrix whose SVD is needed.
[out]UThe left singular vectors of A.
[out]SThe singular values of A.
[out]VThe right singular vectors of A.
Returns
True if the SVD succeeded.

Definition at line 430 of file IncrementalSVD.cpp.

◆ takeSample()

bool CAROM::IncrementalSVD::takeSample ( double *  u_in,
bool  add_without_increase = false 
)
virtual

Sample new state, u_in, at the given time.

Precondition
u_in != 0
Parameters
[in]u_inThe state at the specified time.
[in]add_without_increaseIf true, addLinearlyDependent is invoked.
Returns
True if the sampling was successful.

Implements CAROM::SVD.

Definition at line 174 of file IncrementalSVD.cpp.


The documentation for this class was generated from the following files: