20 #ifndef included_ParametricDMD_h
21 #define included_ParametricDMD_h
23 #include "manifold_interp/MatrixInterpolator.h"
24 #include "linalg/Matrix.h"
25 #include "linalg/Vector.h"
53 std::vector<Vector*>& parameter_points,
54 std::vector<T*>& dmds,
56 std::string rbf =
"G",
57 std::string interp_method =
"LS",
58 double closest_rbf_val = 0.9,
59 bool reorthogonalize_W =
false)
61 CAROM_VERIFY(parameter_points.size() == dmds.size());
62 CAROM_VERIFY(dmds.size() > 1);
63 for (
int i = 0; i < dmds.size() - 1; i++)
65 CAROM_VERIFY(dmds[i]->d_dt == dmds[i + 1]->d_dt);
66 CAROM_VERIFY(dmds[i]->d_k == dmds[i + 1]->d_k);
68 CAROM_VERIFY(closest_rbf_val >= 0.0 && closest_rbf_val <= 1.0);
71 MPI_Initialized(&mpi_init);
73 MPI_Init(
nullptr,
nullptr);
76 MPI_Comm_rank(MPI_COMM_WORLD, &rank);
78 std::vector<CAROM::Matrix*> bases;
79 std::vector<CAROM::Matrix*> A_tildes;
80 for (
int i = 0; i < dmds.size(); i++)
82 bases.push_back(dmds[i]->d_basis);
83 A_tildes.push_back(dmds[i]->d_A_tilde);
92 rotation_matrices, bases, ref_point,
"B", rbf, interp_method, closest_rbf_val);
97 rotation_matrices, A_tildes, ref_point,
"R", rbf, interp_method,
103 std::vector<std::complex<double>> eigs = eigenpair.
eigs;
109 parametric_dmd =
new T(eigs, phi_real, phi_imaginary, dmds[0]->d_k,
110 dmds[0]->d_dt, dmds[0]->d_t_offset,
111 dmds[0]->d_state_offset);
118 for (
auto m : rotation_matrices)
142 std::vector<Vector*>& parameter_points,
143 std::vector<std::string>& dmd_paths,
145 std::string rbf =
"G",
146 std::string interp_method =
"LS",
147 double closest_rbf_val = 0.9,
148 bool reorthogonalize_W =
false)
150 std::vector<T*> dmds;
151 for (
int i = 0; i < dmd_paths.size(); i++)
153 T* dmd =
new T(dmd_paths[i]);
158 rbf, interp_method, closest_rbf_val,
160 for (
int i = 0; i < dmds.size(); i++)
Matrix * mult(const Matrix &other) const
Multiplies this Matrix with other and returns the product, reference version.
Matrix * interpolate(Vector *point, bool orthogonalize=false)
Obtain the interpolated reduced matrix of the unsampled parameter point.
void getParametricDMD(T *¶metric_dmd, std::vector< Vector * > ¶meter_points, std::vector< T * > &dmds, Vector *desired_point, std::string rbf="G", std::string interp_method="LS", double closest_rbf_val=0.9, bool reorthogonalize_W=false)
Constructor.
struct ComplexEigenPair NonSymmetricRightEigenSolve(Matrix *A)
Computes the eigenvectors/eigenvalues of an NxN real nonsymmetric matrix.
int getClosestPoint(std::vector< Vector * > &points, Vector *test_point)
Get closest point to a test point among a group of points.
std::vector< Matrix * > obtainRotationMatrices(std::vector< Vector * > parameter_points, std::vector< Matrix * > bases, int ref_point)
Obtain the rotation matrices for all the parameter points using the basis of the reference point.
Matrix * ev_real
The real component of the eigenvectors in matrix form.
std::vector< std::complex< double > > eigs
The corresponding complex eigenvalues.
Matrix * ev_imaginary
The imaginary component of the eigenvectors in matrix form.