ROMHandler: projection routines for vector.

LLNL · Jan 12, 2024 · 4f30164 · 4f30164
1 parent 88981ac
commit 4f30164
Show file tree

Hide file tree

Showing 2 changed files with 45 additions and 17 deletions.
diff --git a/include/rom_handler.hpp b/include/rom_handler.hpp
@@ -148,9 +148,11 @@ class ROMHandlerBase
    virtual void ProjectInterfaceToDomainBasis(const int &c1, const int &c2, const Array2D<Operator*> &mats, Array2D<SparseMatrix *> &rom_mats) = 0;
    virtual void ProjectVariableToDomainBasis(const int &vi, const int &vj, const Array2D<Operator*> &mats, Array2D<SparseMatrix *> &rom_mats) = 0;
    virtual void ProjectGlobalToDomainBasis(const Array2D<Operator*> &mats, Array2D<SparseMatrix *> &rom_mats) = 0;
-
    virtual void ProjectOperatorOnReducedBasis(const Array2D<Operator*> &mats) = 0;
-   virtual void ProjectVectorOnReducedBasis(const BlockVector* vec, mfem::BlockVector*& rom_vec) = 0;
+
+   virtual void ProjectToRefBasis(const int &i, const Vector &vec, Vector &rom_vec) = 0;
+   virtual void ProjectToDomainBasis(const int &i, const Vector &vec, Vector &rom_vec) = 0;
+   virtual void ProjectGlobalToDomainBasis(const BlockVector* vec, mfem::BlockVector*& rom_vec) = 0;
    virtual void ProjectRHSOnReducedBasis(const BlockVector* RHS) = 0;
 
    virtual void Solve(BlockVector* U) = 0;
@@ -224,11 +226,14 @@ class MFEMROMHandler : public ROMHandlerBase
    virtual void ProjectInterfaceToDomainBasis(const int &c1, const int &c2, const Array2D<Operator*> &mats, Array2D<SparseMatrix *> &rom_mats);
    virtual void ProjectVariableToDomainBasis(const int &vi, const int &vj, const Array2D<Operator*> &mats, Array2D<SparseMatrix *> &rom_mats);
    virtual void ProjectGlobalToDomainBasis(const Array2D<Operator*> &mats, Array2D<SparseMatrix *> &rom_mats);
-
    virtual void ProjectOperatorOnReducedBasis(const Array2D<Operator*> &mats);
-   virtual void ProjectVectorOnReducedBasis(const BlockVector* vec, mfem::BlockVector*& rom_vec);
+
+   virtual void ProjectToRefBasis(const int &i, const Vector &vec, Vector &rom_vec);
+   virtual void ProjectToDomainBasis(const int &i, const Vector &vec, Vector &rom_vec);
+   virtual void ProjectGlobalToDomainBasis(const BlockVector* vec, mfem::BlockVector*& rom_vec);
    virtual void ProjectRHSOnReducedBasis(const BlockVector* RHS) override
-   { ProjectVectorOnReducedBasis(RHS, reduced_rhs); }
+   { ProjectGlobalToDomainBasis(RHS, reduced_rhs); }
+
    virtual void Solve(BlockVector* U);
    virtual void NonlinearSolve(Operator &oper, BlockVector* U, Solver *prec=NULL) override;
 

diff --git a/src/rom_handler.cpp b/src/rom_handler.cpp
@@ -405,24 +405,47 @@ void MFEMROMHandler::SaveOperator(const std::string input_prefix)
    assert(errf >= 0);
 }
 
-void MFEMROMHandler::ProjectVectorOnReducedBasis(const BlockVector* vec, BlockVector*& rom_vec)
+void MFEMROMHandler::ProjectToRefBasis(const int &i, const Vector &vec, Vector &rom_vec)
 {
-   assert(vec->NumBlocks() == numSub);
+   assert(basis_loaded);
+   assert((i >= 0) && (i < num_rom_ref_blocks));
+   DenseMatrix* basis_i;
+
+   GetReferenceBasis(i, basis_i);
+   assert(vec.Size() == basis_i->NumRows());
+   assert(rom_vec.Size() == basis_i->NumCols());
+
+   basis_i->MultTranspose(vec, rom_vec);
+}
+
+void MFEMROMHandler::ProjectToDomainBasis(const int &i, const Vector &vec, Vector &rom_vec)
+{
+   assert(basis_loaded);
+   assert((i >= 0) && (i < num_rom_blocks));
+   DenseMatrix* basis_i;
+
+   GetDomainBasis(i, basis_i);
+   assert(vec.Size() == basis_i->NumRows());
+   assert(rom_vec.Size() == basis_i->NumCols());
+
+   basis_i->MultTranspose(vec, rom_vec);
+}
+
+void MFEMROMHandler::ProjectGlobalToDomainBasis(const BlockVector* vec, BlockVector*& rom_vec)
+{
+   assert(vec->NumBlocks() == num_rom_blocks);
    // reset rom_vec if initiated a priori.
    if (rom_vec) delete rom_vec;
 
-   rom_vec = new BlockVector(rom_block_offsets);
-
-   if (!basis_loaded) LoadReducedBasis();
+   if (separate_variable)
+      rom_vec = new BlockVector(rom_varblock_offsets);
+   else
+      rom_vec = new BlockVector(rom_block_offsets);
 
-   // Each basis is applied to the same column blocks.
-   for (int i = 0; i < numSub; i++)
+   for (int i = 0; i < num_rom_blocks; i++)
    {
-      assert(vec->GetBlock(i).Size() == fom_num_vdofs[i]);
-
-      DenseMatrix* basis_i;
-      GetBasisOnSubdomain(i, basis_i);
-      basis_i->MultTranspose(vec->GetBlock(i).GetData(), rom_vec->GetBlock(i).GetData());
+      int idx = (separate_variable) ? (i % num_var) * numSub + (i / num_var) : i;
+      ProjectToDomainBasis(i, vec->GetBlock(i), rom_vec->GetBlock(idx));
    }
 }