14 #ifndef included_IncrementalSVD_h
15 #define included_IncrementalSVD_h
18 #include "linalg/Options.h"
19 #include "utils/Database.h"
43 const std::string& basis_file_name);
65 bool add_without_increase =
false);
131 double* u,
bool add_without_increase =
false);
bool d_save_state
If true the state of the SVD will be written to disk when the object is deleted. If there are multipl...
double checkOrthogonality(const Matrix *m)
Computes and returns the orthogonality of m.
std::vector< int > d_proc_dims
The dimension of the system on each processor.
Database * d_state_database
Pointer to the database that will hold saved state data if the state is to be saved.
bool d_update_right_SV
Whether to update the right singular vectors.
double d_linearity_tol
Tolerance to determine whether or not a sample is linearly dependent.
virtual void addNewSample(const Vector *j, const Matrix *A, const Matrix *W, Matrix *sigma)=0
Add a new, unique sample to the SVD.
static const int COMMUNICATE_U
MPI message tag.
virtual void computeBasis()=0
Computes the current basis vectors.
virtual const Vector * getSingularValues()
Returns the singular values for the current time interval.
int d_max_basis_dimension
The maximum basis dimension.
long int d_total_dim
The total dimension of the system.
std::string d_state_file_name
The name of file to which state is saved or restored from.
virtual const Matrix * getTemporalBasis()
Returns the temporal basis vectors for the current time interval as a Matrix.
void constructQ(double *&Q, const Vector *l, double k)
Construct the matrix Q whose SVD is needed.
virtual bool buildIncrementalSVD(double *u, bool add_without_increase=false)
Adds the new sampled the state vector, u, to the system.
int numSamples()
The number of samples stored.
int d_size
The total number of processors.
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....
virtual bool takeSample(double *u_in, bool add_without_increase=false)
Sample new state, u_in, at the given time.
virtual ~IncrementalSVD()
Destructor.
virtual void buildInitialSVD(double *u)=0
Constructs the first SVD.
bool d_skip_linearly_dependent
Whether to skip linearly dependent samples.
virtual const Matrix * getSnapshotMatrix()
Returns the snapshot matrix for the current time interval.
virtual const Matrix * getSpatialBasis()
Returns the basis vectors for the current time interval as a Matrix.
virtual void addLinearlyDependentSample(const Matrix *A, const Matrix *W, const Matrix *sigma)=0
Add a linearly dependent sample to the SVD.
int d_rank
The rank of the processor owning this object.
int d_num_samples
Number of samples stored for the current time interval.