libROM  v1.0
Data-driven physical simulation library
BasisGenerator.h
1 
11 // Description: The wrapper class for an SVD algorithm and
12 // sampler. This class controls all aspects of basis
13 // vector generation.
14 
15 #ifndef included_BasisGenerator_h
16 #define included_BasisGenerator_h
17 
18 #include "BasisWriter.h"
19 #include "BasisReader.h"
20 #include "Options.h"
21 #include "svd/SVD.h"
22 
23 #include "mpi.h"
24 
25 #include <cmath>
26 
27 /* Use C++11 built-in shared pointers if available; else fallback to Boost. */
28 #if __cplusplus >= 201103L
29 #include <memory>
30 #else
31 #include <boost/shared_ptr.hpp>
32 #endif
33 
34 #include <string.h>
35 
36 namespace CAROM {
37 
38 class BasisWriter;
39 class BasisReader;
40 class Matrix;
41 
48 {
49 public:
64  Options options,
65  bool incremental,
66  const std::string& basis_file_name = "",
67  Database::formats file_format = Database::formats::HDF5);
68 
73 
83  bool
85  double time);
86 
92  bool
94  {
95  return d_update_right_SV;
96  }
97 
111  bool
112  takeSample(
113  double* u_in,
114  bool add_without_increase = false);
115 
122  void
123  endSamples(const std::string& kind = "basis")
124  {
125  if (d_basis_writer) {
126  d_basis_writer->writeBasis(kind);
127  }
128  }
129 
133  void
135  {
136  if (d_basis_writer) {
137  d_basis_writer->writeBasis("snapshot");
138  }
139  }
140 
153  void
154  loadSampleRange(const std::string& base_file_name,
155  const std::string& kind = "basis",
156  int col_min = 0,
157  int col_max = 1e9,
158  Database::formats db_format = Database::formats::HDF5);
159 
170  void
171  loadSamples(const std::string& base_file_name,
172  const std::string& kind = "basis",
173  int cut_off = 1e9,
174  Database::formats db_format = Database::formats::HDF5);
175 
187  double
189  double* u_in,
190  double* rhs_in,
191  double time);
192 
199  const Matrix*
201  {
202  return d_svd->getSpatialBasis();
203  }
204 
211  const Matrix*
213  {
214  return d_svd->getTemporalBasis();
215  }
216 
223  const Vector*
225  {
226  return d_svd->getSingularValues();
227  }
228 
234  const Matrix*
236  {
237  return d_svd->getSnapshotMatrix();
238  }
239 
245  int getNumSamples() const
246  {
247  return d_svd->getNumSamples();
248  }
249 
259  void finalSummary(
260  const double energyFractionThreshold,
261  int & cutoff,
262  const std::string & cutoffOutputPath = "",
263  const int first_sv = 0);
264 
265 protected:
270 
275 
280 
284 #if __cplusplus >= 201103L
285  std::shared_ptr<SVD> d_svd;
286 #else
287  boost::shared_ptr<SVD> d_svd;
288 #endif
289 
290 private:
294  BasisGenerator();
299  const BasisGenerator& other);
300 
305  operator = (
306  const BasisGenerator& rhs);
307 
313  void
314  resetDt(
315  double new_dt);
316 
322  int
323  getDim()
324  {
325  return d_dim;
326  }
327 
331  bool d_incremental;
332 
336  bool d_update_right_SV;
337 
343  double d_tol;
344 
348  double d_max_time_between_samples;
349 
353  double d_min_sampling_time_step_scale;
354 
358  double d_sampling_time_step_scale;
359 
363  double d_max_sampling_time_step_scale;
364 
368  double d_dt;
369 
373  double d_next_sample_time;
374 
378  int d_num_procs;
379 
385  const int d_dim;
386 };
387 
388 }
389 
390 #endif
void loadSampleRange(const std::string &base_file_name, const std::string &kind="basis", int col_min=0, int col_max=1e9, Database::formats db_format=Database::formats::HDF5)
Load previously saved sample (basis or state) within a column range.
void writeSnapshot()
Write current snapshot matrix.
const Matrix * getSnapshotMatrix()
Returns the snapshot matrix for the current time interval.
int getNumSamples() const
Returns the number of samples taken.
bool d_write_snapshots
Whether to write snapshots instead of bases.
bool updateRightSV()
Check whether right basis vectors will be updated.
BasisReader * d_basis_reader
Reader of basis vectors.
bool isNextSample(double time)
Returns true if it is time for the next svd sample.
void endSamples(const std::string &kind="basis")
Signal that the final sample has been taken.
const Matrix * getTemporalBasis()
Returns the temporal basis vectors for the current time interval as a Matrix.
BasisWriter * d_basis_writer
Writer of basis vectors.
const Matrix * getSpatialBasis()
Returns the basis vectors for the current time interval as a Matrix.
void loadSamples(const std::string &base_file_name, const std::string &kind="basis", int cut_off=1e9, Database::formats db_format=Database::formats::HDF5)
Load previously saved sample (basis or state).
boost::shared_ptr< SVD > d_svd
Pointer to the abstract SVD algorithm object.
~BasisGenerator()
Destructor.
double computeNextSampleTime(double *u_in, double *rhs_in, double time)
Computes next time an svd sample is needed.
void finalSummary(const double energyFractionThreshold, int &cutoff, const std::string &cutoffOutputPath="", const int first_sv=0)
Prints the summary of recommended numbers of basis vectors.
const Vector * getSingularValues()
Returns the singular values for the current time interval as a Vector.
bool takeSample(double *u_in, bool add_without_increase=false)
Sample the new state, u_in, at the given time.
void writeBasis(const std::string &kind="basis")
Write basis or state vectors generated by d_basis_generator.
Definition: BasisWriter.cpp:74
formats
Implemented database file formats. Add to this enum each time a new database format is implemented.
Definition: Database.h:311