From 1d5c203677809cd817403ebfbecded0314a308d9 Mon Sep 17 00:00:00 2001 From: Dylan Copeland Date: Fri, 14 Jul 2023 14:53:30 -0700 Subject: [PATCH 01/11] Adding QDEIM option to some examples. --- .gitignore | 3 ++ examples/prom/mixed_nonlinear_diffusion.cpp | 28 ++++++++++++++ .../prom/nonlinear_elasticity_global_rom.cpp | 37 ++++++++++++++----- 3 files changed, 58 insertions(+), 10 deletions(-) diff --git a/.gitignore b/.gitignore index 34dc5c0a6..016626da8 100644 --- a/.gitignore +++ b/.gitignore @@ -12,6 +12,9 @@ dependencies *.log *.status +# Do not track backup files +*~ + # Do not track autogenerated config headers CAROM_config.h FCMangle.h diff --git a/examples/prom/mixed_nonlinear_diffusion.cpp b/examples/prom/mixed_nonlinear_diffusion.cpp index 2388d6cfc..816fdb728 100644 --- a/examples/prom/mixed_nonlinear_diffusion.cpp +++ b/examples/prom/mixed_nonlinear_diffusion.cpp @@ -91,6 +91,7 @@ #include "linalg/BasisReader.h" #include "linalg/NNLS.h" #include "hyperreduction/DEIM.h" +#include "hyperreduction/QDEIM.h" #include "hyperreduction/GNAT.h" #include "hyperreduction/S_OPT.h" #include "mfem/SampleMesh.hpp" @@ -550,6 +551,7 @@ int main(int argc, char *argv[]) bool offline = false; bool merge = false; bool online = false; + bool use_qdeim = false; bool use_sopt = false; bool use_eqp = false; bool writeSampleMesh = false; @@ -619,6 +621,8 @@ int main(int argc, char *argv[]) "Enable or disable the online phase."); args.AddOption(&merge, "-merge", "--merge", "-no-merge", "--no-merge", "Enable or disable the merge phase."); + args.AddOption(&use_qdeim, "-qdeim", "--qdeim", "-no-qdeim", "--no-qdeim", + "Use Q-DEIM sampling instead of DEIM for the hyperreduction."); args.AddOption(&use_sopt, "-sopt", "--sopt", "-no-sopt", "--no-sopt", "Use S-OPT sampling instead of DEIM for the hyperreduction."); args.AddOption(&num_samples_req, "-nsr", "--nsr", @@ -1032,6 +1036,19 @@ int main(int argc, char *argv[]) num_procs, nsamp_R); } + else if (use_qdeim) + { + if (myid == 0) + printf("Using QDEIM sampling\n"); + CAROM::QDEIM(FR_librom, + nldim, + sample_dofs, + num_sample_dofs_per_proc, + *Bsinv, + myid, + num_procs, + nsamp_R); + } else if (nsamp_R != nldim) { if (myid == 0) @@ -1106,6 +1123,17 @@ int main(int argc, char *argv[]) num_procs, nsamp_S); } + else if (use_qdeim) + { + CAROM::QDEIM(S_librom, + nsdim, + sample_dofs_S, + num_sample_dofs_per_proc_S, + *Ssinv, + myid, + num_procs, + nsamp_S); + } else if (nsamp_S != nsdim) { CAROM::GNAT(S_librom, diff --git a/examples/prom/nonlinear_elasticity_global_rom.cpp b/examples/prom/nonlinear_elasticity_global_rom.cpp index e004086a9..c7871f8f7 100644 --- a/examples/prom/nonlinear_elasticity_global_rom.cpp +++ b/examples/prom/nonlinear_elasticity_global_rom.cpp @@ -30,25 +30,25 @@ // and nonlinear term basis, with velocity initial condition: // // Offline phase: -// ./nonlinear_elasticity_global_rom -offline -dt 0.01 -tf 5.0 -s 14 -vs 100 -sc 0.90 -id 0 +// ./nonlinear_elasticity_global_rom -offline -dt 0.01 -tf 5.0 -s 14 -vs 100 -sc 0.9 -id 0 // -// ./nonlinear_elasticity_global_rom -offline -dt 0.01 -tf 5.0 -s 14 -vs 100 -sc 1.10 -id 1 +// ./nonlinear_elasticity_global_rom -offline -dt 0.01 -tf 5.0 -s 14 -vs 100 -sc 1.1 -id 1 // // Merge phase: // ./nonlinear_elasticity_global_rom -merge -ns 2 -dt 0.01 -tf 5.0 // // Create FOM comparison data: -// ./nonlinear_elasticity_global_rom -offline -dt 0.01 -tf 5.0 -s 14 -vs 100 -sc 1.00 -id 2 +// ./nonlinear_elasticity_global_rom -offline -dt 0.01 -tf 5.0 -s 14 -vs 100 -sc 1.0 -id 2 // // Online phase with full sampling: -// ./nonlinear_elasticity_global_rom -online -dt 0.01 -tf 5.0 -s 14 -vs 100 -hyp -rvdim 40 -rxdim 10 -hdim 71 -nsr 1170 -sc 1.00 +// ./nonlinear_elasticity_global_rom -online -dt 0.01 -tf 5.0 -s 14 -vs 100 -hyp -rvdim 40 -rxdim 10 -hdim 71 -nsr 1170 -sc 1.0 // Output message: // Elapsed time for time integration loop 1.80759 // Relative error of ROM position (x) at t_final: 5 is 0.000231698 // Relative error of ROM velocity (v) at t_final: 5 is 0.466941 // // Online phase with strong hyper-reduction: -// ./nonlinear_elasticity_global_rom -online -dt 0.01 -tf 5.0 -s 14 -vs 100 -hyp -rvdim 40 -rxdim 10 -hdim 71 -nsr 100 -sc 1.00 +// ./nonlinear_elasticity_global_rom -online -dt 0.01 -tf 5.0 -s 14 -vs 100 -hyp -rvdim 40 -rxdim 10 -hdim 71 -nsr 100 -sc 1.0 // Output message: // Elapsed time for time integration loop 1.08048 // Relative error of ROM position (x) at t_final: 5 is 0.00209877 @@ -59,24 +59,24 @@ // Sample runs and results for parametric ROM using only displacement basis // and nonlinear term basis: // Offline phase: -// ./nonlinear_elasticity_global_rom -offline -dt 0.01 -tf 5.0 -s 14 -vs 100 -sc 0.90 -xbo -def-ic -id 0 -// ./nonlinear_elasticity_global_rom -offline -dt 0.01 -tf 5.0 -s 14 -vs 100 -sc 1.10 -xbo -def-ic -id 1 +// ./nonlinear_elasticity_global_rom -offline -dt 0.01 -tf 5.0 -s 14 -vs 100 -sc 0.9 -xbo -def-ic -id 0 +// ./nonlinear_elasticity_global_rom -offline -dt 0.01 -tf 5.0 -s 14 -vs 100 -sc 1.1 -xbo -def-ic -id 1 // // Merge phase: // ./nonlinear_elasticity_global_rom -merge -ns 2 -dt 0.01 -tf 5.0 -xbo // // Create FOM comparison data: -// ./nonlinear_elasticity_global_rom -offline -dt 0.01 -tf 5.0 -s 14 -vs 100 -sc 1.00 -xbo -def-ic -id 2 +// ./nonlinear_elasticity_global_rom -offline -dt 0.01 -tf 5.0 -s 14 -vs 100 -sc 1.0 -xbo -def-ic -id 2 // // Online phase with full sampling: -// ./nonlinear_elasticity_global_rom -online -dt 0.01 -tf 5.0 -s 14 -vs 100 -hyp -rxdim 57 -hdim 183 -nsr 1170 -sc 1.00 -xbo -def-ic +// ./nonlinear_elasticity_global_rom -online -dt 0.01 -tf 5.0 -s 14 -vs 100 -hyp -rxdim 57 -hdim 183 -nsr 1170 -sc 1.0 -xbo -def-ic // Output message: // Elapsed time for time integration loop 18.9874 // Relative error of ROM position (x) at t_final: 5 is 7.08272e-05 // Relative error of ROM velocity (v) at t_final: 5 is 0.00387647 // // Online phase with strong hyper reduction: -// ./nonlinear_elasticity_global_rom -online -dt 0.01 -tf 5.0 -s 14 -vs 100 -hyp -rxdim 2 -hdim 4 -nsr 10 -sc 1.00 -xbo -def-ic +// ./nonlinear_elasticity_global_rom -online -dt 0.01 -tf 5.0 -s 14 -vs 100 -hyp -rxdim 2 -hdim 4 -nsr 10 -sc 1.0 -xbo -def-ic // Output message: // Elapsed time for time integration loop 1.01136 // Relative error of ROM position (x) at t_final: 5 is 0.0130818 @@ -90,6 +90,7 @@ #include "linalg/BasisGenerator.h" #include "linalg/BasisReader.h" #include "hyperreduction/DEIM.h" +#include "hyperreduction/QDEIM.h" #include "hyperreduction/GNAT.h" #include "hyperreduction/S_OPT.h" #include "mfem/SampleMesh.hpp" @@ -388,6 +389,7 @@ int main(int argc, char* argv[]) bool offline = false; bool merge = false; bool online = false; + bool use_qdeim = false; bool use_sopt = false; bool hyperreduce = true; bool x_base_only = false; @@ -444,6 +446,8 @@ int main(int argc, char* argv[]) "Enable or disable the online phase."); args.AddOption(&merge, "-merge", "--merge", "-no-merge", "--no-merge", "Enable or disable the merge phase."); + args.AddOption(&use_qdeim, "-qdeim", "--qdeim", "-no-qdeim", "--no-qdeim", + "Use QDEIM sampling instead of DEIM for the hyperreduction."); args.AddOption(&use_sopt, "-sopt", "--sopt", "-no-sopt", "--no-sopt", "Use S-OPT sampling instead of DEIM for the hyperreduction."); args.AddOption(&num_samples_req, "-nsr", "--nsr", @@ -851,6 +855,19 @@ int main(int argc, char* argv[]) num_procs, nsamp_H); } + else if (use_qdeim) + { + if (myid == 0) + printf("Using QDEIM sampling\n"); + CAROM::QDEIM(H_librom, + hdim, + sample_dofs, + num_sample_dofs_per_proc, + *Hsinv, + myid, + num_procs, + nsamp_H); + } else if (nsamp_H != hdim) { if (myid == 0) From a18d70f1ea0fb81488a6dc5f292e43625d7d5609 Mon Sep 17 00:00:00 2001 From: Dylan Copeland Date: Tue, 12 Sep 2023 16:03:21 -0700 Subject: [PATCH 02/11] Adding new Hyperreduction class to simplify calls to different types of hyperreduction sampling functions. --- examples/prom/mixed_nonlinear_diffusion.cpp | 147 ++++-------------- .../prom/nonlinear_elasticity_global_rom.cpp | 115 +++++--------- lib/CMakeLists.txt | 1 + lib/hyperreduction/Hyperreduction.cpp | 87 +++++++++++ lib/hyperreduction/Hyperreduction.h | 73 +++++++++ lib/hyperreduction/Utilities.h | 4 - unit_tests/test_IncrementalSVDBrand.cpp | 22 +-- 7 files changed, 240 insertions(+), 209 deletions(-) create mode 100644 lib/hyperreduction/Hyperreduction.cpp create mode 100644 lib/hyperreduction/Hyperreduction.h diff --git a/examples/prom/mixed_nonlinear_diffusion.cpp b/examples/prom/mixed_nonlinear_diffusion.cpp index 816fdb728..95277f313 100644 --- a/examples/prom/mixed_nonlinear_diffusion.cpp +++ b/examples/prom/mixed_nonlinear_diffusion.cpp @@ -31,8 +31,9 @@ // Analytic test (reproductive) // mpirun -n 1 ./mixed_nonlinear_diffusion -offline // mpirun -n 1 ./mixed_nonlinear_diffusion -merge -ns 1 -// mpirun -n 1 ./mixed_nonlinear_diffusion -online -rrdim 8 -rwdim 8 -nldim 20 -nsdim 20 -// mpirun -n 1 ./mixed_nonlinear_diffusion -online -rrdim 8 -rwdim 8 -nldim 20 -nsdim 20 -sopt +// mpirun -n 1 ./mixed_nonlinear_diffusion -online -rrdim 8 -rwdim 8 -nldim 20 -nsdim 20 -hrtype deim +// mpirun -n 1 ./mixed_nonlinear_diffusion -online -rrdim 8 -rwdim 8 -nldim 20 -nsdim 20 -hrtype qdeim +// mpirun -n 1 ./mixed_nonlinear_diffusion -online -rrdim 8 -rwdim 8 -nldim 20 -nsdim 20 -hrtype sopt // mpirun -n 1 ./mixed_nonlinear_diffusion -online -rrdim 8 -rwdim 8 -nldim 20 -nsdim 20 -ns 1 -eqp // // Relative l2 error of ROM solution at final timestep using DEIM sampling: 1.096776797994166e-08 @@ -48,8 +49,8 @@ // Initial step test (reproductive) // mpirun -n 1 ./mixed_nonlinear_diffusion -offline -p 1 // mpirun -n 1 ./mixed_nonlinear_diffusion -merge -ns 1 -p 1 -// mpirun -n 1 ./mixed_nonlinear_diffusion -online -rrdim 8 -rwdim 8 -nldim 20 -p 1 -// mpirun -n 1 ./mixed_nonlinear_diffusion -online -rrdim 8 -rwdim 8 -nldim 20 -p 1 -sopt +// mpirun -n 1 ./mixed_nonlinear_diffusion -online -rrdim 8 -rwdim 8 -nldim 20 -p 1 -hrtype deim +// mpirun -n 1 ./mixed_nonlinear_diffusion -online -rrdim 8 -rwdim 8 -nldim 20 -p 1 -hrtype sopt // mpirun -n 1 ./mixed_nonlinear_diffusion -online -rrdim 8 -rwdim 8 -nldim 20 -p 1 -ns 1 -eqp // // Relative l2 error of ROM solution at final timestep using DEIM sampling: 0.0003712362376412496 @@ -65,8 +66,8 @@ // mpirun -n 1 ./mixed_nonlinear_diffusion -p 1 -offline -id 2 -sh 0.35 // mpirun -n 1 ./mixed_nonlinear_diffusion -p 1 -merge -ns 3 // mpirun -n 1 ./mixed_nonlinear_diffusion -p 1 -offline -id 3 -sh 0.3 -// mpirun -n 1 ./mixed_nonlinear_diffusion -p 1 -online -rrdim 8 -rwdim 8 -nldim 20 -sh 0.3 -id 3 -// mpirun -n 1 ./mixed_nonlinear_diffusion -p 1 -online -rrdim 8 -rwdim 8 -nldim 20 -sh 0.3 -id 3 -sopt +// mpirun -n 1 ./mixed_nonlinear_diffusion -p 1 -online -rrdim 8 -rwdim 8 -nldim 20 -sh 0.3 -id 3 -hrtype deim +// mpirun -n 1 ./mixed_nonlinear_diffusion -p 1 -online -rrdim 8 -rwdim 8 -nldim 20 -sh 0.3 -id 3 -hrtype sopt // mpirun -n 1 ./mixed_nonlinear_diffusion -p 1 -online -rrdim 8 -rwdim 8 -nldim 20 -sh 0.3 -id 3 -ns 3 -eqp -maxnnls 30 // // Relative l2 error of ROM solution at final timestep using DEIM sampling: 0.002681387312231006 @@ -90,10 +91,7 @@ #include "linalg/BasisGenerator.h" #include "linalg/BasisReader.h" #include "linalg/NNLS.h" -#include "hyperreduction/DEIM.h" -#include "hyperreduction/QDEIM.h" -#include "hyperreduction/GNAT.h" -#include "hyperreduction/S_OPT.h" +#include "hyperreduction/Hyperreduction.h" #include "mfem/SampleMesh.hpp" #include "mfem/PointwiseSnapshot.hpp" @@ -551,8 +549,7 @@ int main(int argc, char *argv[]) bool offline = false; bool merge = false; bool online = false; - bool use_qdeim = false; - bool use_sopt = false; + const char *samplingType = "deim"; bool use_eqp = false; bool writeSampleMesh = false; int num_samples_req = -1; @@ -621,10 +618,8 @@ int main(int argc, char *argv[]) "Enable or disable the online phase."); args.AddOption(&merge, "-merge", "--merge", "-no-merge", "--no-merge", "Enable or disable the merge phase."); - args.AddOption(&use_qdeim, "-qdeim", "--qdeim", "-no-qdeim", "--no-qdeim", - "Use Q-DEIM sampling instead of DEIM for the hyperreduction."); - args.AddOption(&use_sopt, "-sopt", "--sopt", "-no-sopt", "--no-sopt", - "Use S-OPT sampling instead of DEIM for the hyperreduction."); + args.AddOption(&samplingType, "-hrtype", "--hrsamplingtype", + "Sampling type for hyperreduction."); args.AddOption(&num_samples_req, "-nsr", "--nsr", "Number of samples for the sampling algorithm to select."); args.AddOption(&use_eqp, "-eqp", "--eqp", "-no-eqp", "--no-eqp", @@ -954,8 +949,6 @@ int main(int argc, char *argv[]) CAROM::BasisReader readerFR("basisFR"); FR_librom = readerFR.getSpatialBasis(0.0); - // Compute sample points using DEIM, for hyperreduction - if (nldim == -1) { nldim = FR_librom->numColumns(); @@ -1009,6 +1002,7 @@ int main(int argc, char *argv[]) else { // Setup hyperreduction using DEIM, GNAT, or S-OPT + CAROM::Hyperreduction hr(samplingType); vector num_sample_dofs_per_proc(num_procs); if (num_samples_req != -1) @@ -1023,57 +1017,15 @@ int main(int argc, char *argv[]) // Now execute the chosen sampling algorithm to get the sampling information. Bsinv = new CAROM::Matrix(nsamp_R, nldim, false); vector sample_dofs(nsamp_R); // Indices of the sampled rows - if (use_sopt) - { - if (myid == 0) - printf("Using S_OPT sampling\n"); - CAROM::S_OPT(FR_librom, - nldim, - sample_dofs, - num_sample_dofs_per_proc, - *Bsinv, - myid, - num_procs, - nsamp_R); - } - else if (use_qdeim) - { - if (myid == 0) - printf("Using QDEIM sampling\n"); - CAROM::QDEIM(FR_librom, - nldim, - sample_dofs, - num_sample_dofs_per_proc, - *Bsinv, - myid, - num_procs, - nsamp_R); - } - else if (nsamp_R != nldim) - { - if (myid == 0) - printf("Using GNAT sampling\n"); - CAROM::GNAT(FR_librom, - nldim, - sample_dofs, - num_sample_dofs_per_proc, - *Bsinv, - myid, - num_procs, - nsamp_R); - } - else - { - if (myid == 0) - printf("Using DEIM sampling\n"); - CAROM::DEIM(FR_librom, - nldim, - sample_dofs, - num_sample_dofs_per_proc, - *Bsinv, - myid, - num_procs); - } + + hr.ComputeSamples(FR_librom, + nldim, + sample_dofs, + num_sample_dofs_per_proc, + *Bsinv, + myid, + num_procs, + nsamp_R); vector sample_dofs_S; // Indices of the sampled rows vector num_sample_dofs_per_proc_S(num_procs); @@ -1085,8 +1037,6 @@ int main(int argc, char *argv[]) readerS = new CAROM::BasisReader("basisS"); S_librom = readerS->getSpatialBasis(0.0); - // Compute sample points using DEIM - if (nsdim == -1) { nsdim = S_librom->numColumns(); @@ -1100,7 +1050,7 @@ int main(int argc, char *argv[]) if (myid == 0) printf("reduced S dim = %d\n",nsdim); - // Now execute the DEIM algorithm to get the sampling information. + // Now use a hyperreduction method to compute samples. if (num_samples_req != -1) { nsamp_S = num_samples_req; @@ -1112,49 +1062,16 @@ int main(int argc, char *argv[]) Ssinv = new CAROM::Matrix(nsamp_S, nsdim, false); sample_dofs_S.resize(nsamp_S); - if (use_sopt) - { - CAROM::S_OPT(S_librom, - nsdim, - sample_dofs_S, - num_sample_dofs_per_proc_S, - *Ssinv, - myid, - num_procs, - nsamp_S); - } - else if (use_qdeim) - { - CAROM::QDEIM(S_librom, - nsdim, - sample_dofs_S, - num_sample_dofs_per_proc_S, - *Ssinv, - myid, - num_procs, - nsamp_S); - } - else if (nsamp_S != nsdim) - { - CAROM::GNAT(S_librom, - nsdim, - sample_dofs_S, - num_sample_dofs_per_proc_S, - *Ssinv, - myid, - num_procs, - nsamp_S); - } - else - { - CAROM::DEIM(S_librom, - nsdim, - sample_dofs_S, - num_sample_dofs_per_proc_S, - *Ssinv, - myid, - num_procs); - } + + + hr.ComputeSamples(S_librom, + nsdim, + sample_dofs_S, + num_sample_dofs_per_proc_S, + *Ssinv, + myid, + num_procs, + nsamp_S); } // Construct sample mesh diff --git a/examples/prom/nonlinear_elasticity_global_rom.cpp b/examples/prom/nonlinear_elasticity_global_rom.cpp index c7871f8f7..bc8505377 100644 --- a/examples/prom/nonlinear_elasticity_global_rom.cpp +++ b/examples/prom/nonlinear_elasticity_global_rom.cpp @@ -89,10 +89,7 @@ #include "linalg/Vector.h" #include "linalg/BasisGenerator.h" #include "linalg/BasisReader.h" -#include "hyperreduction/DEIM.h" -#include "hyperreduction/QDEIM.h" -#include "hyperreduction/GNAT.h" -#include "hyperreduction/S_OPT.h" +#include "hyperreduction/Hyperreduction.h" #include "mfem/SampleMesh.hpp" #include "mfem/Utilities.hpp" @@ -167,7 +164,7 @@ class RomOperator : public TimeDependentOperator bool hyperreduce; bool x_base_only; - CAROM::Vector* pfom_librom, * pfom_v_librom; + CAROM::Vector *pfom_librom, *pfom_v_librom; Vector* pfom; Vector* pfom_x; Vector* pfom_v; @@ -263,7 +260,8 @@ CAROM::Matrix* GetFirstColumns(const int N, const CAROM::Matrix* A) // TODO: move this to the library? void BasisGeneratorFinalSummary(CAROM::BasisGenerator* bg, - const double energyFraction, int& cutoff, const std::string cutoffOutputPath) + const double energyFraction, int& cutoff, + const std::string cutoffOutputPath) { const int rom_dim = bg->getSpatialBasis()->numColumns(); const CAROM::Vector* sing_vals = bg->getSingularValues(); @@ -389,16 +387,15 @@ int main(int argc, char* argv[]) bool offline = false; bool merge = false; bool online = false; - bool use_qdeim = false; - bool use_sopt = false; bool hyperreduce = true; bool x_base_only = false; int num_samples_req = -1; + const char *samplingType = "gnat"; int nsets = 0; int id_param = 0; - // number of basis vectors to use + // Number of basis vectors to use int rxdim = -1; int rvdim = -1; int hdim = -1; @@ -446,10 +443,8 @@ int main(int argc, char* argv[]) "Enable or disable the online phase."); args.AddOption(&merge, "-merge", "--merge", "-no-merge", "--no-merge", "Enable or disable the merge phase."); - args.AddOption(&use_qdeim, "-qdeim", "--qdeim", "-no-qdeim", "--no-qdeim", - "Use QDEIM sampling instead of DEIM for the hyperreduction."); - args.AddOption(&use_sopt, "-sopt", "--sopt", "-no-sopt", "--no-sopt", - "Use S-OPT sampling instead of DEIM for the hyperreduction."); + args.AddOption(&samplingType, "-hrtype", "--hrsamplingtype", + "Sampling type for hyperreduction."); args.AddOption(&num_samples_req, "-nsr", "--nsr", "number of samples we want to select for the sampling algorithm."); args.AddOption(&rxdim, "-rxdim", "--rxdim", @@ -573,8 +568,8 @@ int main(int argc, char* argv[]) BlockVector vx(true_offset); ParGridFunction v_gf, x_gf; - v_gf.MakeTRef(&fespace, vx, - true_offset[0]); // Associate a new FiniteElementSpace and new true-dof data with the GridFunction. + // Associate a FiniteElementSpace and true-dof data with the GridFunctions. + v_gf.MakeTRef(&fespace, vx, true_offset[0]); x_gf.MakeTRef(&fespace, vx, true_offset[1]); ParGridFunction x_ref(&fespace); @@ -692,7 +687,7 @@ int main(int argc, char* argv[]) vis_v.precision(8); visualize(vis_v, pmesh, &x_gf, &v_gf, "Velocity", true); // Make sure all ranks have sent their 'v' solution before initiating - // another set of GLVis connections (one from each rank): + // another set of GLVis connections (one from each rank) MPI_Barrier(pmesh->GetComm()); vis_w.open(vishost, visport); if (vis_w) @@ -842,57 +837,17 @@ int main(int argc, char* argv[]) CAROM::Matrix* Hsinv = new CAROM::Matrix(nsamp_H, hdim, false); vector sample_dofs(nsamp_H); - if (use_sopt) - { - if (myid == 0) - printf("Using S_OPT sampling\n"); - CAROM::S_OPT(H_librom, - hdim, - sample_dofs, - num_sample_dofs_per_proc, - *Hsinv, - myid, - num_procs, - nsamp_H); - } - else if (use_qdeim) - { - if (myid == 0) - printf("Using QDEIM sampling\n"); - CAROM::QDEIM(H_librom, - hdim, - sample_dofs, - num_sample_dofs_per_proc, - *Hsinv, - myid, - num_procs, - nsamp_H); - } - else if (nsamp_H != hdim) - { - if (myid == 0) - printf("Using GNAT sampling\n"); - CAROM::GNAT(H_librom, - hdim, - sample_dofs, - num_sample_dofs_per_proc, - *Hsinv, - myid, - num_procs, - nsamp_H); - } - else - { - if (myid == 0) - printf("Using DEIM sampling\n"); - CAROM::DEIM(H_librom, - hdim, - sample_dofs, - num_sample_dofs_per_proc, - *Hsinv, - myid, - num_procs); - } + + // Setup hyperreduction using DEIM, GNAT, or S-OPT + CAROM::Hyperreduction hr(samplingType); + hr.ComputeSamples(H_librom, + hdim, + sample_dofs, + num_sample_dofs_per_proc, + *Hsinv, + myid, + num_procs, + nsamp_H); // Construct sample mesh const int nspaces = 1; @@ -921,7 +876,8 @@ int main(int argc, char* argv[]) w_x = new CAROM::Vector(rxdim, false); *w = 0.0; - // Note that some of this could be done only on the ROM solver process, but it is tricky, since RomOperator assembles Bsp in parallel. + // Note that some of this could be done only on the ROM solver process, + // but it is tricky, since RomOperator assembles Bsp in parallel. wMFEM = new Vector(&((*w)(0)), rxdim + rvdim); // Initial conditions @@ -947,8 +903,8 @@ int main(int argc, char* argv[]) // 12. Set the initial conditions for v_gf, x_gf and vx, and define the // boundary conditions on a beam-like mesh (see description above). - sp_v_gf.MakeTRef(sp_XV_space, sp_vx, - sp_offset[0]); // Associate a new FiniteElementSpace and new true-dof data with the GridFunction. + // Associate a FiniteElementSpace and true-dof data with the GridFunctions. + sp_v_gf.MakeTRef(sp_XV_space, sp_vx, sp_offset[0]); sp_x_gf.MakeTRef(sp_XV_space, sp_vx, sp_offset[1]); VectorFunctionCoefficient* velo = 0; @@ -1134,16 +1090,16 @@ int main(int argc, char* argv[]) if (x_base_only == false && basis_generator_v->isNextSample(t)) { basis_generator_v->takeSample(vx_diff.GetBlock(0), t, dt); - basis_generator_v->computeNextSampleTime(vx_diff.GetBlock(0), dvdt.GetData(), - t); + basis_generator_v->computeNextSampleTime(vx_diff.GetBlock(0), + dvdt.GetData(), t); basis_generator_H->takeSample(oper.H_sp.GetData(), t, dt); } if (basis_generator_x->isNextSample(t)) { basis_generator_x->takeSample(vx_diff.GetBlock(1), t, dt); - basis_generator_x->computeNextSampleTime(vx_diff.GetBlock(1), dxdt.GetData(), - t); + basis_generator_x->computeNextSampleTime(vx_diff.GetBlock(1), + dxdt.GetData(), t); if (x_base_only == true) { @@ -1224,7 +1180,8 @@ int main(int argc, char* argv[]) if (offline) { - // Sample final solution, to prevent extrapolation in ROM between the last sample and the end of the simulation. + // Sample final solution, to prevent extrapolation in ROM between the + // last sample and the end of the simulation. dvxdt = oper.dvxdt_sp.GetData(); vx_diff = BlockVector(vx); vx_diff -= vx0; @@ -1320,10 +1277,10 @@ int main(int argc, char* argv[]) if (myid == 0) { - cout << "Relative error of ROM position (x) at t_final: " << t_final << - " is " << tot_diff_norm_x / tot_x_fom_norm << endl; - cout << "Relative error of ROM velocity (v) at t_final: " << t_final << - " is " << tot_diff_norm_v / tot_v_fom_norm << endl; + cout << "Relative error of ROM position (x) at t_final: " << t_final + << " is " << tot_diff_norm_x / tot_x_fom_norm << endl; + cout << "Relative error of ROM velocity (v) at t_final: " << t_final + << " is " << tot_diff_norm_v / tot_v_fom_norm << endl; } } diff --git a/lib/CMakeLists.txt b/lib/CMakeLists.txt index 677d1a799..a523dd67d 100644 --- a/lib/CMakeLists.txt +++ b/lib/CMakeLists.txt @@ -50,6 +50,7 @@ set(module_list hyperreduction/S_OPT hyperreduction/STSampling hyperreduction/Utilities + hyperreduction/Hyperreduction utils/Database utils/HDFDatabase utils/CSVDatabase diff --git a/lib/hyperreduction/Hyperreduction.cpp b/lib/hyperreduction/Hyperreduction.cpp new file mode 100644 index 000000000..b110dc326 --- /dev/null +++ b/lib/hyperreduction/Hyperreduction.cpp @@ -0,0 +1,87 @@ +/****************************************************************************** + * + * Copyright (c) 2013-2023, Lawrence Livermore National Security, LLC + * and other libROM project developers. See the top-level COPYRIGHT + * file for details. + * + * SPDX-License-Identifier: (Apache-2.0 OR MIT) + * + *****************************************************************************/ + +#include "Hyperreduction.h" +#include "DEIM.h" +#include "QDEIM.h" +#include "GNAT.h" +#include "S_OPT.h" + +#include "linalg/Matrix.h" +#include "utils/Utilities.h" + +namespace CAROM { + +Hyperreduction::Hyperreduction(const char* sampling_type) +{ + auto iter = SamplingTypeMap.find(sampling_type); + CAROM_VERIFY(iter != std::end(SamplingTypeMap)); + + samplingType = iter->second; +} + +void Hyperreduction::ComputeSamples(const Matrix* f_basis, + int num_f_basis_vectors_used, + std::vector& f_sampled_row, + std::vector& f_sampled_rows_per_proc, + Matrix& f_basis_sampled_inv, + int myid, + int num_procs, + const int num_samples_req, + std::vector *init_samples, + bool qr_factorize) +{ + switch (samplingType) + { + case deim: + CAROM_VERIFY(num_samples_req == f_basis->numColumns()); + DEIM(f_basis, + num_f_basis_vectors_used, + f_sampled_row, + f_sampled_rows_per_proc, + f_basis_sampled_inv, + myid, num_procs); + return; + case gnat: + GNAT(f_basis, + num_f_basis_vectors_used, + f_sampled_row, + f_sampled_rows_per_proc, + f_basis_sampled_inv, + myid, num_procs, + num_samples_req, + init_samples); + return; + case qdeim: + QDEIM(f_basis, + num_f_basis_vectors_used, + f_sampled_row, + f_sampled_rows_per_proc, + f_basis_sampled_inv, + myid, num_procs, + num_samples_req); + return; + case sopt: + S_OPT(f_basis, + num_f_basis_vectors_used, + f_sampled_row, + f_sampled_rows_per_proc, + f_basis_sampled_inv, + myid, num_procs, + num_samples_req, + init_samples, + qr_factorize); + return; + default: + CAROM_ERROR("Sampling type not supported"); + } +} + +} diff --git a/lib/hyperreduction/Hyperreduction.h b/lib/hyperreduction/Hyperreduction.h new file mode 100644 index 000000000..4fd4fc91a --- /dev/null +++ b/lib/hyperreduction/Hyperreduction.h @@ -0,0 +1,73 @@ +/****************************************************************************** + * + * Copyright (c) 2013-2023, Lawrence Livermore National Security, LLC + * and other libROM project developers. See the top-level COPYRIGHT + * file for details. + * + * SPDX-License-Identifier: (Apache-2.0 OR MIT) + * + *****************************************************************************/ + +// Description: Interface to hyperreduction algorithms. + +#ifndef included_Hyperreduction_h +#define included_Hyperreduction_h + + +#include +#include + +#include + +namespace CAROM { + +enum SamplingType +{ + deim, // Default, DEIM + gnat, // GNAT + qdeim, // QDEIM + sopt // S-OPT +}; + +static std::unordered_map SamplingTypeMap = +{ + {"deim", deim}, + {"gnat", gnat}, + {"qdeim", qdeim}, + {"sopt", sopt} +}; + +class Matrix; + +class Hyperreduction +{ +public: + + Hyperreduction(SamplingType stype) : + samplingType(stype) + { } + + Hyperreduction(const char* sampling_type); + + void SetSamplingType(SamplingType stype) + { + samplingType = stype; + } + + void ComputeSamples(const Matrix* f_basis, + int num_f_basis_vectors_used, + std::vector& f_sampled_row, + std::vector& f_sampled_rows_per_proc, + Matrix& f_basis_sampled_inv, + int myid, + int num_procs, + const int num_samples_req = -1, + std::vector *init_samples=NULL, + bool qr_factorize = false); + +private: + SamplingType samplingType; +}; + +} +#endif diff --git a/lib/hyperreduction/Utilities.h b/lib/hyperreduction/Utilities.h index 2fa04304d..2572a47c5 100644 --- a/lib/hyperreduction/Utilities.h +++ b/lib/hyperreduction/Utilities.h @@ -17,8 +17,6 @@ namespace CAROM { -#ifndef DOXYGEN_IGNORE - /** * @brief Struct to hold the local maximum absolute value of a basis vector, * the row it is in, and the processor that owns it. We will reduce this @@ -37,8 +35,6 @@ typedef struct */ void RowInfoMax(RowInfo* a, RowInfo* b, int* len, MPI_Datatype* type); -#endif /* DOXYGEN_IGNORE */ - } #endif diff --git a/unit_tests/test_IncrementalSVDBrand.cpp b/unit_tests/test_IncrementalSVDBrand.cpp index eeaddd1d5..5cc287ca2 100644 --- a/unit_tests/test_IncrementalSVDBrand.cpp +++ b/unit_tests/test_IncrementalSVDBrand.cpp @@ -72,9 +72,9 @@ TEST(IncrementalSVDBrandTest, Test_IncrementalSVDBrand) double* basis_right_true_ans = new double[9] { -1.78651649346571794741E-01, 5.44387957786310106023E-01, -8.19588518467042281834E-01, - -9.49719639253861602768E-01, -3.13100149275943651084E-01, -9.50441422536040881122E-04, - -2.57130696341890396805E-01, 7.78209514167382598870E-01, 5.72951792961765460355E-01 - }; + -9.49719639253861602768E-01, -3.13100149275943651084E-01, -9.50441422536040881122E-04, + -2.57130696341890396805E-01, 7.78209514167382598870E-01, 5.72951792961765460355E-01 + }; double* sv_true_ans = new double[3] { 4.84486375065219387892E+00, 3.66719976398777269821E+00, 2.69114625366671811335E+00 @@ -83,18 +83,18 @@ TEST(IncrementalSVDBrandTest, Test_IncrementalSVDBrand) bool fast_update = true; bool fast_update_brand = true; CAROM::Options incremental_svd_options = CAROM::Options(d_num_rows, 3, -1, true) - .setMaxBasisDimension(num_total_rows) + .setMaxBasisDimension(num_total_rows) .setIncrementalSVD(1e-1, - 1e-1, - 1e-1, - 1e-1, - fast_update, - fast_update_brand, - false); + 1e-1, + 1e-1, + 1e-1, + fast_update, + fast_update_brand, + false); CAROM::BasisGenerator sampler( incremental_svd_options, - true, + true, "irrelevant.txt"); sampler.takeSample(&sample1[row_offset[d_rank]], 0, 1e-1); sampler.takeSample(&sample2[row_offset[d_rank]], 0, 1e-1); From d47c5582ac91f63efe4c0d90490fe5b364e0bda0 Mon Sep 17 00:00:00 2001 From: Dylan Copeland Date: Tue, 12 Sep 2023 16:16:26 -0700 Subject: [PATCH 03/11] Minor edits. --- lib/hyperreduction/GNAT.cpp | 20 +++++++++------- lib/hyperreduction/QDEIM.cpp | 45 +++++++++++++++++++----------------- lib/hyperreduction/S_OPT.cpp | 6 ++--- 3 files changed, 39 insertions(+), 32 deletions(-) diff --git a/lib/hyperreduction/GNAT.cpp b/lib/hyperreduction/GNAT.cpp index 83c47b8f5..3c247787f 100644 --- a/lib/hyperreduction/GNAT.cpp +++ b/lib/hyperreduction/GNAT.cpp @@ -97,8 +97,8 @@ void GNAT(const Matrix* f_basis, const int ns_mod_nr = num_samples % num_basis_vectors; int ns = 0; - // The small matrix inverted by the algorithm. We'll allocate the largest - // matrix we'll need and set its size at each step in the algorithm. + // The small matrix inverted by the algorithm. We allocate the largest + // size needed and set the size at each step in the algorithm. Matrix M(num_samples, std::max(num_basis_vectors-1, 1), false); // Scratch space used throughout the algorithm. @@ -125,8 +125,8 @@ void GNAT(const Matrix* f_basis, std::vector all_num_init_samples(num_procs); std::vector all_init_samples(total_num_init_samples); - MPI_Allgather(&num_init_samples, 1, MPI_INT, all_num_init_samples.data(), 1, - MPI_INT, MPI_COMM_WORLD); + MPI_Allgather(&num_init_samples, 1, MPI_INT, all_num_init_samples.data(), + 1, MPI_INT, MPI_COMM_WORLD); for (int i = 0; i < myid; ++i) { @@ -134,7 +134,7 @@ void GNAT(const Matrix* f_basis, } } - // Figure out the 1st sampled rows of the RHS. + // Figure out the first sampled rows of the RHS. RowInfo f_bv_max_local, f_bv_max_global; const int ns0 = 0 < ns_mod_nr ? (num_samples / num_basis_vectors) + 1 : @@ -215,10 +215,14 @@ void GNAT(const Matrix* f_basis, double tmp = 0.0; for (int minv_col = 0; minv_col < ns; ++minv_col) { if (ns == i) + { tmp += M.item(minv_row, minv_col)*tmp_fs.item(minv_col, i); + } else - tmp += M.item(minv_col, minv_row)*tmp_fs.item(minv_col, - i); // Transposing M^+, which is stored as its transpose. + { + // Transposing M^+, which is stored as its transpose. + tmp += M.item(minv_col, minv_row)*tmp_fs.item(minv_col, i); + } } c[minv_row] = tmp; } @@ -288,7 +292,7 @@ void GNAT(const Matrix* f_basis, CAROM_ASSERT(num_samples == ns); - // Fill f_sampled_row, and f_sampled_rows_per_proc. Unscramble tmp_fs into + // Fill f_sampled_row, and f_sampled_rows_per_proc. Unscramble tmp_fs into // f_basis_sampled_inv. int idx = 0; for (int i = 0; i < num_procs; ++i) { diff --git a/lib/hyperreduction/QDEIM.cpp b/lib/hyperreduction/QDEIM.cpp index a07e66c65..e11eea942 100644 --- a/lib/hyperreduction/QDEIM.cpp +++ b/lib/hyperreduction/QDEIM.cpp @@ -37,8 +37,8 @@ QDEIM(const Matrix* f_basis, // processor that owns each sampled row, and fills f_basis_sampled_inv with the // sampled rows of the basis of the RHS. - CAROM_VERIFY(num_f_basis_vectors_used == - f_basis->numColumns()); // The QR implementation uses the entire matrix. + // The QR implementation uses the entire matrix. + CAROM_VERIFY(num_f_basis_vectors_used == f_basis->numColumns()); CAROM_VERIFY(f_basis->numColumns() <= num_samples_req && num_samples_req <= f_basis->numDistributedRows()); CAROM_VERIFY(num_samples_req == f_basis_sampled_inv.numRows() @@ -99,7 +99,8 @@ QDEIM(const Matrix* f_basis, ns[owner]++; } - // Reorder f_sampled_row and f_sampled_row_owner to match the order of f_basis_sampled_inv + // Reorder f_sampled_row and f_sampled_row_owner to match the order + // of f_basis_sampled_inv for (int i=0; imult(V); // distributed CAROM_VERIFY(Ubt->distributed() && Ubt->numRows() == f_basis->numRows() @@ -254,9 +256,9 @@ QDEIM(const Matrix* f_basis, for (int i=0; iitem(i, j) * Ubt->item(i, - j); // column sums of Ub.^2, which are row sums of Ubt.^2 + r[i] += Ubt->item(i, j) * Ubt->item(i,j); r[i] -= sqrt((r[i] * r[i]) - (4.0 * g * Ubt->item(i, n-1) * Ubt->item(i, n-1))); } @@ -311,8 +313,9 @@ QDEIM(const Matrix* f_basis, globalSamples.insert(f_sampled_row[s]); } - // Send one row of f_basis, corresponding to sample s, to root process for f_basis_sampled_inv. - // First, scatter from root to tell the owning process the sample index. + // Send one row of f_basis, corresponding to sample s, to the root + // process for f_basis_sampled_inv. First, scatter from root to tell + // the owning process the sample index. int sample = -1; MPI_Scatter(ns.data(), 1, MPI_INT, &sample, 1, MPI_INT, 0, MPI_COMM_WORLD); @@ -321,22 +324,22 @@ QDEIM(const Matrix* f_basis, if (sample > -1) { CAROM_VERIFY(sample >= row_offset[myid]); - MPI_Send(f_basis->getData() + ((sample - row_offset[myid]) * numCol), numCol, - MPI_DOUBLE, 0, tagSendRecv, MPI_COMM_WORLD); + MPI_Send(f_basis->getData() + ((sample - row_offset[myid]) * numCol), + numCol, MPI_DOUBLE, 0, tagSendRecv, MPI_COMM_WORLD); } if (myid == 0) { MPI_Status status; - MPI_Recv(sampled_row_data.data() + (s*numCol), numCol, MPI_DOUBLE, owner, - tagSendRecv, MPI_COMM_WORLD, &status); + MPI_Recv(sampled_row_data.data() + (s*numCol), numCol, MPI_DOUBLE, + owner, tagSendRecv, MPI_COMM_WORLD, &status); } delete Ubt; } // loop s over samples - // Subtract row_offset to convert f_sampled_row from global to local indices - // Also, reorder f_sampled_row by process + // Subtract row_offset to convert f_sampled_row from global to local indices. + // Also, reorder f_sampled_row by process. if (myid == 0) { ns[0] = 0; @@ -414,8 +417,8 @@ QDEIM(const Matrix* f_basis, if (!f_basis->distributed()) { - CAROM_VERIFY(numCol == - num_samples_req); // GappyPOD+E not implemented for oversampling if not distributed + // GappyPOD+E not implemented for oversampling if not distributed + CAROM_VERIFY(numCol == num_samples_req); } // Now invert f_basis_sampled_inv, storing its transpose. diff --git a/lib/hyperreduction/S_OPT.cpp b/lib/hyperreduction/S_OPT.cpp index 7efa656b6..8166febe7 100644 --- a/lib/hyperreduction/S_OPT.cpp +++ b/lib/hyperreduction/S_OPT.cpp @@ -380,8 +380,8 @@ S_OPT(const Matrix* f_basis, { tmp += g1.item(j, k) * GG.item(j, k); } - A->item(j) = std::max(0.0, ata + (Vo->item(j, i - 1) * Vo->item(j, - i - 1)) - tmp); + A->item(j) = std::max(0.0, ata + + (Vo->item(j, i - 1) * Vo->item(j, i - 1)) - tmp); } nV.setSize(i); @@ -491,7 +491,7 @@ S_OPT(const Matrix* f_basis, delete noM; } - // Fill f_sampled_row, and f_sampled_rows_per_proc. Unscramble tmp_fs into + // Fill f_sampled_row, and f_sampled_rows_per_proc. Unscramble tmp_fs into // f_basis_sampled_inv. int idx = 0; for (int i = 0; i < num_procs; ++i) { From 604f53e9f2059fa6c6e245c5e864af92c6b1bc95 Mon Sep 17 00:00:00 2001 From: Dylan Copeland Date: Tue, 12 Sep 2023 16:20:47 -0700 Subject: [PATCH 04/11] Style --- .../prom/nonlinear_elasticity_global_rom.cpp | 32 +++++++++---------- lib/hyperreduction/GNAT.cpp | 12 +++---- lib/hyperreduction/QDEIM.cpp | 20 ++++++------ lib/hyperreduction/S_OPT.cpp | 2 +- 4 files changed, 33 insertions(+), 33 deletions(-) diff --git a/examples/prom/nonlinear_elasticity_global_rom.cpp b/examples/prom/nonlinear_elasticity_global_rom.cpp index bc8505377..217d18588 100644 --- a/examples/prom/nonlinear_elasticity_global_rom.cpp +++ b/examples/prom/nonlinear_elasticity_global_rom.cpp @@ -261,7 +261,7 @@ CAROM::Matrix* GetFirstColumns(const int N, const CAROM::Matrix* A) // TODO: move this to the library? void BasisGeneratorFinalSummary(CAROM::BasisGenerator* bg, const double energyFraction, int& cutoff, - const std::string cutoffOutputPath) + const std::string cutoffOutputPath) { const int rom_dim = bg->getSpatialBasis()->numColumns(); const CAROM::Vector* sing_vals = bg->getSingularValues(); @@ -838,16 +838,16 @@ int main(int argc, char* argv[]) CAROM::Matrix* Hsinv = new CAROM::Matrix(nsamp_H, hdim, false); vector sample_dofs(nsamp_H); - // Setup hyperreduction using DEIM, GNAT, or S-OPT - CAROM::Hyperreduction hr(samplingType); - hr.ComputeSamples(H_librom, - hdim, - sample_dofs, - num_sample_dofs_per_proc, - *Hsinv, - myid, - num_procs, - nsamp_H); + // Setup hyperreduction using DEIM, GNAT, or S-OPT + CAROM::Hyperreduction hr(samplingType); + hr.ComputeSamples(H_librom, + hdim, + sample_dofs, + num_sample_dofs_per_proc, + *Hsinv, + myid, + num_procs, + nsamp_H); // Construct sample mesh const int nspaces = 1; @@ -877,7 +877,7 @@ int main(int argc, char* argv[]) *w = 0.0; // Note that some of this could be done only on the ROM solver process, - // but it is tricky, since RomOperator assembles Bsp in parallel. + // but it is tricky, since RomOperator assembles Bsp in parallel. wMFEM = new Vector(&((*w)(0)), rxdim + rvdim); // Initial conditions @@ -903,7 +903,7 @@ int main(int argc, char* argv[]) // 12. Set the initial conditions for v_gf, x_gf and vx, and define the // boundary conditions on a beam-like mesh (see description above). - // Associate a FiniteElementSpace and true-dof data with the GridFunctions. + // Associate a FiniteElementSpace and true-dof data with the GridFunctions. sp_v_gf.MakeTRef(sp_XV_space, sp_vx, sp_offset[0]); sp_x_gf.MakeTRef(sp_XV_space, sp_vx, sp_offset[1]); @@ -1091,7 +1091,7 @@ int main(int argc, char* argv[]) { basis_generator_v->takeSample(vx_diff.GetBlock(0), t, dt); basis_generator_v->computeNextSampleTime(vx_diff.GetBlock(0), - dvdt.GetData(), t); + dvdt.GetData(), t); basis_generator_H->takeSample(oper.H_sp.GetData(), t, dt); } @@ -1099,7 +1099,7 @@ int main(int argc, char* argv[]) { basis_generator_x->takeSample(vx_diff.GetBlock(1), t, dt); basis_generator_x->computeNextSampleTime(vx_diff.GetBlock(1), - dxdt.GetData(), t); + dxdt.GetData(), t); if (x_base_only == true) { @@ -1181,7 +1181,7 @@ int main(int argc, char* argv[]) if (offline) { // Sample final solution, to prevent extrapolation in ROM between the - // last sample and the end of the simulation. + // last sample and the end of the simulation. dvxdt = oper.dvxdt_sp.GetData(); vx_diff = BlockVector(vx); vx_diff -= vx0; diff --git a/lib/hyperreduction/GNAT.cpp b/lib/hyperreduction/GNAT.cpp index 3c247787f..faa65a243 100644 --- a/lib/hyperreduction/GNAT.cpp +++ b/lib/hyperreduction/GNAT.cpp @@ -126,7 +126,7 @@ void GNAT(const Matrix* f_basis, std::vector all_init_samples(total_num_init_samples); MPI_Allgather(&num_init_samples, 1, MPI_INT, all_num_init_samples.data(), - 1, MPI_INT, MPI_COMM_WORLD); + 1, MPI_INT, MPI_COMM_WORLD); for (int i = 0; i < myid; ++i) { @@ -215,14 +215,14 @@ void GNAT(const Matrix* f_basis, double tmp = 0.0; for (int minv_col = 0; minv_col < ns; ++minv_col) { if (ns == i) - { + { tmp += M.item(minv_row, minv_col)*tmp_fs.item(minv_col, i); - } + } else - { - // Transposing M^+, which is stored as its transpose. + { + // Transposing M^+, which is stored as its transpose. tmp += M.item(minv_col, minv_row)*tmp_fs.item(minv_col, i); - } + } } c[minv_row] = tmp; } diff --git a/lib/hyperreduction/QDEIM.cpp b/lib/hyperreduction/QDEIM.cpp index e11eea942..34fbf3619 100644 --- a/lib/hyperreduction/QDEIM.cpp +++ b/lib/hyperreduction/QDEIM.cpp @@ -100,7 +100,7 @@ QDEIM(const Matrix* f_basis, } // Reorder f_sampled_row and f_sampled_row_owner to match the order - // of f_basis_sampled_inv + // of f_basis_sampled_inv for (int i=0; imult(V); // distributed CAROM_VERIFY(Ubt->distributed() && Ubt->numRows() == f_basis->numRows() @@ -256,7 +256,7 @@ QDEIM(const Matrix* f_basis, for (int i=0; iitem(i, j) * Ubt->item(i,j); @@ -314,7 +314,7 @@ QDEIM(const Matrix* f_basis, } // Send one row of f_basis, corresponding to sample s, to the root - // process for f_basis_sampled_inv. First, scatter from root to tell + // process for f_basis_sampled_inv. First, scatter from root to tell // the owning process the sample index. int sample = -1; @@ -325,7 +325,7 @@ QDEIM(const Matrix* f_basis, { CAROM_VERIFY(sample >= row_offset[myid]); MPI_Send(f_basis->getData() + ((sample - row_offset[myid]) * numCol), - numCol, MPI_DOUBLE, 0, tagSendRecv, MPI_COMM_WORLD); + numCol, MPI_DOUBLE, 0, tagSendRecv, MPI_COMM_WORLD); } if (myid == 0) @@ -417,7 +417,7 @@ QDEIM(const Matrix* f_basis, if (!f_basis->distributed()) { - // GappyPOD+E not implemented for oversampling if not distributed + // GappyPOD+E not implemented for oversampling if not distributed CAROM_VERIFY(numCol == num_samples_req); } diff --git a/lib/hyperreduction/S_OPT.cpp b/lib/hyperreduction/S_OPT.cpp index 8166febe7..f21534dc0 100644 --- a/lib/hyperreduction/S_OPT.cpp +++ b/lib/hyperreduction/S_OPT.cpp @@ -381,7 +381,7 @@ S_OPT(const Matrix* f_basis, tmp += g1.item(j, k) * GG.item(j, k); } A->item(j) = std::max(0.0, ata + - (Vo->item(j, i - 1) * Vo->item(j, i - 1)) - tmp); + (Vo->item(j, i - 1) * Vo->item(j, i - 1)) - tmp); } nV.setSize(i); From e6cfaa588c9fc971fd9875a56ea2f8afd8f91cbb Mon Sep 17 00:00:00 2001 From: Dylan Copeland Date: Wed, 13 Sep 2023 16:06:28 -0700 Subject: [PATCH 05/11] Update regression tests. --- regression_tests/tests/prom/mixed_nonlinear_diffusion.sh | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/regression_tests/tests/prom/mixed_nonlinear_diffusion.sh b/regression_tests/tests/prom/mixed_nonlinear_diffusion.sh index bf8af6063..d1a01f4a5 100755 --- a/regression_tests/tests/prom/mixed_nonlinear_diffusion.sh +++ b/regression_tests/tests/prom/mixed_nonlinear_diffusion.sh @@ -3,8 +3,8 @@ source $GITHUB_WORKSPACE/regression_tests/common.sh CMDS=( "./mixed_nonlinear_diffusion -offline -tf 0.01" "./mixed_nonlinear_diffusion -merge -ns 1 -tf 0.01" - "./mixed_nonlinear_diffusion -online -tf 0.01 -rrdim 2 -rwdim 2 -nldim 3" - "./mixed_nonlinear_diffusion -online -tf 0.01 -rrdim 2 -rwdim 2 -nldim 3 -sopt" + "./mixed_nonlinear_diffusion -online -tf 0.01 -rrdim 2 -rwdim 2 -nldim 3 -hrtype deim" + "./mixed_nonlinear_diffusion -online -tf 0.01 -rrdim 2 -rwdim 2 -nldim 3 -hrtype sopt" "./mixed_nonlinear_diffusion -online -tf 0.01 -rrdim 2 -rwdim 2 -nldim 3 -ns 1 -eqp -maxnnls 4" ) TYPE="PROM" From 9c9bf28308d6dd6425b6e989d168d0017713e7b8 Mon Sep 17 00:00:00 2001 From: Axel Larsson <65452706+axla-io@users.noreply.github.com> Date: Wed, 20 Sep 2023 10:09:25 -0400 Subject: [PATCH 06/11] Added hrtype to command line examples --- examples/prom/nonlinear_elasticity_global_rom.cpp | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/examples/prom/nonlinear_elasticity_global_rom.cpp b/examples/prom/nonlinear_elasticity_global_rom.cpp index 217d18588..188be1e94 100644 --- a/examples/prom/nonlinear_elasticity_global_rom.cpp +++ b/examples/prom/nonlinear_elasticity_global_rom.cpp @@ -41,14 +41,14 @@ // ./nonlinear_elasticity_global_rom -offline -dt 0.01 -tf 5.0 -s 14 -vs 100 -sc 1.0 -id 2 // // Online phase with full sampling: -// ./nonlinear_elasticity_global_rom -online -dt 0.01 -tf 5.0 -s 14 -vs 100 -hyp -rvdim 40 -rxdim 10 -hdim 71 -nsr 1170 -sc 1.0 +// ./nonlinear_elasticity_global_rom -online -dt 0.01 -tf 5.0 -s 14 -vs 100 -hyp -hrtype gnat -rvdim 40 -rxdim 10 -hdim 71 -nsr 1170 -sc 1.0 // Output message: // Elapsed time for time integration loop 1.80759 // Relative error of ROM position (x) at t_final: 5 is 0.000231698 // Relative error of ROM velocity (v) at t_final: 5 is 0.466941 // // Online phase with strong hyper-reduction: -// ./nonlinear_elasticity_global_rom -online -dt 0.01 -tf 5.0 -s 14 -vs 100 -hyp -rvdim 40 -rxdim 10 -hdim 71 -nsr 100 -sc 1.0 +// ./nonlinear_elasticity_global_rom -online -dt 0.01 -tf 5.0 -s 14 -vs 100 -hyp -hrtype gnat -rvdim 40 -rxdim 10 -hdim 71 -nsr 100 -sc 1.0 // Output message: // Elapsed time for time integration loop 1.08048 // Relative error of ROM position (x) at t_final: 5 is 0.00209877 @@ -69,14 +69,14 @@ // ./nonlinear_elasticity_global_rom -offline -dt 0.01 -tf 5.0 -s 14 -vs 100 -sc 1.0 -xbo -def-ic -id 2 // // Online phase with full sampling: -// ./nonlinear_elasticity_global_rom -online -dt 0.01 -tf 5.0 -s 14 -vs 100 -hyp -rxdim 57 -hdim 183 -nsr 1170 -sc 1.0 -xbo -def-ic +// ./nonlinear_elasticity_global_rom -online -dt 0.01 -tf 5.0 -s 14 -vs 100 -hyp -hrtype gnat -rxdim 57 -hdim 183 -nsr 1170 -sc 1.0 -xbo -def-ic // Output message: // Elapsed time for time integration loop 18.9874 // Relative error of ROM position (x) at t_final: 5 is 7.08272e-05 // Relative error of ROM velocity (v) at t_final: 5 is 0.00387647 // // Online phase with strong hyper reduction: -// ./nonlinear_elasticity_global_rom -online -dt 0.01 -tf 5.0 -s 14 -vs 100 -hyp -rxdim 2 -hdim 4 -nsr 10 -sc 1.0 -xbo -def-ic +// ./nonlinear_elasticity_global_rom -online -dt 0.01 -tf 5.0 -s 14 -vs 100 -hyp -hrtype gnat -rxdim 2 -hdim 4 -nsr 10 -sc 1.0 -xbo -def-ic // Output message: // Elapsed time for time integration loop 1.01136 // Relative error of ROM position (x) at t_final: 5 is 0.0130818 From 97cf09c32ef6a856aea4bf25fe7bc174947f72c7 Mon Sep 17 00:00:00 2001 From: Axel Larsson <65452706+axla-io@users.noreply.github.com> Date: Wed, 20 Sep 2023 16:10:41 -0400 Subject: [PATCH 07/11] Increased max iterations for M_solver and M_hat_solver --- examples/prom/nonlinear_elasticity_global_rom.cpp | 5 ++++- 1 file changed, 4 insertions(+), 1 deletion(-) diff --git a/examples/prom/nonlinear_elasticity_global_rom.cpp b/examples/prom/nonlinear_elasticity_global_rom.cpp index 188be1e94..e3deb484e 100644 --- a/examples/prom/nonlinear_elasticity_global_rom.cpp +++ b/examples/prom/nonlinear_elasticity_global_rom.cpp @@ -1355,7 +1355,7 @@ HyperelasticOperator::HyperelasticOperator(ParFiniteElementSpace& f, M_solver.iterative_mode = false; M_solver.SetRelTol(rel_tol); M_solver.SetAbsTol(0.0); - M_solver.SetMaxIter(30); + M_solver.SetMaxIter(1000); M_solver.SetPrintLevel(0); M_prec.SetType(HypreSmoother::Jacobi); M_solver.SetPreconditioner(M_prec); @@ -1510,6 +1510,9 @@ RomOperator::RomOperator(HyperelasticOperator* fom_, M_hat = new CAROM::Matrix(rvdim, rvdim, false); M_hat_inv = new CAROM::Matrix(rvdim, rvdim, false); + // Set the max iterations for the mass matrix solver + M_hat_solver->SetMaxIter(1000); + // Create S_hat ComputeCtAB(fom->Smat, V_v, V_v, *S_hat); From c05bda78a19a6fc87200cbd62223144a6700e2c2 Mon Sep 17 00:00:00 2001 From: Axel Larsson <65452706+axla-io@users.noreply.github.com> Date: Wed, 20 Sep 2023 16:13:56 -0400 Subject: [PATCH 08/11] Small bug fix --- examples/prom/nonlinear_elasticity_global_rom.cpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/examples/prom/nonlinear_elasticity_global_rom.cpp b/examples/prom/nonlinear_elasticity_global_rom.cpp index e3deb484e..23cc0512b 100644 --- a/examples/prom/nonlinear_elasticity_global_rom.cpp +++ b/examples/prom/nonlinear_elasticity_global_rom.cpp @@ -1511,7 +1511,7 @@ RomOperator::RomOperator(HyperelasticOperator* fom_, M_hat_inv = new CAROM::Matrix(rvdim, rvdim, false); // Set the max iterations for the mass matrix solver - M_hat_solver->SetMaxIter(1000); + M_hat_solver.SetMaxIter(1000); // Create S_hat ComputeCtAB(fom->Smat, V_v, V_v, *S_hat); From 893d7cb7b0306252f6685887047bd54079f42995 Mon Sep 17 00:00:00 2001 From: Axel Larsson <65452706+axla-io@users.noreply.github.com> Date: Wed, 20 Sep 2023 17:43:32 -0400 Subject: [PATCH 09/11] Revert changes made in error --- examples/prom/nonlinear_elasticity_global_rom.cpp | 3 +++ 1 file changed, 3 insertions(+) diff --git a/examples/prom/nonlinear_elasticity_global_rom.cpp b/examples/prom/nonlinear_elasticity_global_rom.cpp index 188be1e94..b0e66197f 100644 --- a/examples/prom/nonlinear_elasticity_global_rom.cpp +++ b/examples/prom/nonlinear_elasticity_global_rom.cpp @@ -1510,6 +1510,9 @@ RomOperator::RomOperator(HyperelasticOperator* fom_, M_hat = new CAROM::Matrix(rvdim, rvdim, false); M_hat_inv = new CAROM::Matrix(rvdim, rvdim, false); + // Set the max iterations for the mass matrix solver + M_hat_solver.SetMaxIter(1000); + // Create S_hat ComputeCtAB(fom->Smat, V_v, V_v, *S_hat); From 16a76797504d47ae8735ab72c56bd92bd6f95a53 Mon Sep 17 00:00:00 2001 From: Axel Larsson <65452706+axla-io@users.noreply.github.com> Date: Wed, 20 Sep 2023 17:44:57 -0400 Subject: [PATCH 10/11] Actually reverting --- examples/prom/nonlinear_elasticity_global_rom.cpp | 5 +---- 1 file changed, 1 insertion(+), 4 deletions(-) diff --git a/examples/prom/nonlinear_elasticity_global_rom.cpp b/examples/prom/nonlinear_elasticity_global_rom.cpp index 23cc0512b..188be1e94 100644 --- a/examples/prom/nonlinear_elasticity_global_rom.cpp +++ b/examples/prom/nonlinear_elasticity_global_rom.cpp @@ -1355,7 +1355,7 @@ HyperelasticOperator::HyperelasticOperator(ParFiniteElementSpace& f, M_solver.iterative_mode = false; M_solver.SetRelTol(rel_tol); M_solver.SetAbsTol(0.0); - M_solver.SetMaxIter(1000); + M_solver.SetMaxIter(30); M_solver.SetPrintLevel(0); M_prec.SetType(HypreSmoother::Jacobi); M_solver.SetPreconditioner(M_prec); @@ -1510,9 +1510,6 @@ RomOperator::RomOperator(HyperelasticOperator* fom_, M_hat = new CAROM::Matrix(rvdim, rvdim, false); M_hat_inv = new CAROM::Matrix(rvdim, rvdim, false); - // Set the max iterations for the mass matrix solver - M_hat_solver.SetMaxIter(1000); - // Create S_hat ComputeCtAB(fom->Smat, V_v, V_v, *S_hat); From 568c86a28f6d4825c0dc8c59c2c56bfe2490de64 Mon Sep 17 00:00:00 2001 From: Dylan Copeland Date: Wed, 18 Oct 2023 20:22:27 -0700 Subject: [PATCH 11/11] Added qdeim regression test. Added results for QDEIM sample runs. --- examples/prom/mixed_nonlinear_diffusion.cpp | 28 +++++++++++-------- .../prom/nonlinear_elasticity_global_rom.cpp | 22 ++++++++++++--- .../tests/prom/mixed_nonlinear_diffusion.sh | 1 + 3 files changed, 36 insertions(+), 15 deletions(-) diff --git a/examples/prom/mixed_nonlinear_diffusion.cpp b/examples/prom/mixed_nonlinear_diffusion.cpp index 95277f313..8b6f7916f 100644 --- a/examples/prom/mixed_nonlinear_diffusion.cpp +++ b/examples/prom/mixed_nonlinear_diffusion.cpp @@ -37,11 +37,13 @@ // mpirun -n 1 ./mixed_nonlinear_diffusion -online -rrdim 8 -rwdim 8 -nldim 20 -nsdim 20 -ns 1 -eqp // // Relative l2 error of ROM solution at final timestep using DEIM sampling: 1.096776797994166e-08 -// Elapsed time for time integration loop using DEIM sampling: 0.6351594580000001 +// Elapsed time for time integration loop using DEIM sampling: 0.5686423310000001 +// Relative l2 error of ROM solution at final timestep using QDEIM sampling: 1.227055339859034e-08 +// Elapsed time for time integration loop using QDEIM sampling: 0.575228571 // Relative l2 error of ROM solution at final timestep using S_OPT sampling: 1.01945081122054e-08 -// Elapsed time for time integration loop using S_OPT sampling: 0.6669736559999999 -// Relative l2 error of ROM solution at final timestep using EQP: 1.46205848438194e-07 -// Elapsed time for time integration loop using EQP: 0.431521853 +// Elapsed time for time integration loop using S_OPT sampling: 0.5311964850000001 +// Relative l2 error of ROM solution at final timestep using EQP: 1.052559143494963e-08 +// Elapsed time for time integration loop using EQP: 0.346165296 // // Note that the timing of the time integration loop does not include setup, // which can be greater for S_OPT and EQP than for DEIM. @@ -54,11 +56,13 @@ // mpirun -n 1 ./mixed_nonlinear_diffusion -online -rrdim 8 -rwdim 8 -nldim 20 -p 1 -ns 1 -eqp // // Relative l2 error of ROM solution at final timestep using DEIM sampling: 0.0003712362376412496 -// Elapsed time for time integration loop using DEIM sampling: 0.364855569 +// Elapsed time for time integration loop using DEIM sampling: 0.339288686 +// Relative l2 error of ROM solution at final timestep using QDEIM sampling: 0.000375066186385267 +// Elapsed time for time integration loop using QDEIM sampling: 0.32320409 // Relative l2 error of ROM solution at final timestep using S_OPT sampling: 0.0003797338657417907 -// Elapsed time for time integration loop using S_OPT sampling: 0.300462563 -// Relative l2 error of ROM solution at final timestep using EQP sampling: 0.0003710336208386964 -// Elapsed time for time integration loop using EQP sampling: 0.481740662 +// Elapsed time for time integration loop using S_OPT sampling: 0.32422266 +// Relative l2 error of ROM solution at final timestep using EQP sampling: 0.0003709977143117392 +// Elapsed time for time integration loop using EQP sampling: 0.457011311 // // Initial step parametric test (predictive) // mpirun -n 1 ./mixed_nonlinear_diffusion -p 1 -offline -id 0 -sh 0.25 @@ -71,11 +75,13 @@ // mpirun -n 1 ./mixed_nonlinear_diffusion -p 1 -online -rrdim 8 -rwdim 8 -nldim 20 -sh 0.3 -id 3 -ns 3 -eqp -maxnnls 30 // // Relative l2 error of ROM solution at final timestep using DEIM sampling: 0.002681387312231006 -// Elapsed time for time integration loop using DEIM sampling: 0.355846074 +// Elapsed time for time integration loop using DEIM sampling: 0.359206509 +// Relative l2 error of ROM solution at final timestep using QDEIM sampling: 0.00266355670204476 +// Elapsed time for time integration loop using QDEIM sampling: 0.309432016 // Relative l2 error of ROM solution at final timestep using S_OPT sampling: 0.002701713369494112 -// Elapsed time for time integration loop using S_OPT sampling: 0.348985935 +// Elapsed time for time integration loop using S_OPT sampling: 0.41139788 // Relative l2 error of ROM solution at final timestep using EQP: 0.002659978000520714 -// Elapsed time for time integration loop using EQP sampling: 0.176821221 +// Elapsed time for time integration loop using EQP sampling: 0.164961644 // // Pointwise snapshots for DMD input // mpirun -n 1 ./mixed_nonlinear_diffusion -pwsnap -pwx 101 -pwy 101 diff --git a/examples/prom/nonlinear_elasticity_global_rom.cpp b/examples/prom/nonlinear_elasticity_global_rom.cpp index 188be1e94..984f1624d 100644 --- a/examples/prom/nonlinear_elasticity_global_rom.cpp +++ b/examples/prom/nonlinear_elasticity_global_rom.cpp @@ -47,13 +47,20 @@ // Relative error of ROM position (x) at t_final: 5 is 0.000231698 // Relative error of ROM velocity (v) at t_final: 5 is 0.466941 // -// Online phase with strong hyper-reduction: +// Online phase with strong hyper-reduction, using GNAT (over-sampled DEIM): // ./nonlinear_elasticity_global_rom -online -dt 0.01 -tf 5.0 -s 14 -vs 100 -hyp -hrtype gnat -rvdim 40 -rxdim 10 -hdim 71 -nsr 100 -sc 1.0 // Output message: -// Elapsed time for time integration loop 1.08048 +// Elapsed time for time integration loop 1.0111 // Relative error of ROM position (x) at t_final: 5 is 0.00209877 // Relative error of ROM velocity (v) at t_final: 5 is 1.39472 // +// Online phase with strong hyper-reduction, using QDEIM: +// ./nonlinear_elasticity_global_rom -online -dt 0.01 -tf 5.0 -s 14 -vs 100 -hyp -hrtype qdeim -rvdim 40 -rxdim 10 -hdim 71 -nsr 100 -sc 1.0 +// Output message: +// Elapsed time for time integration loop 1.02559 +// Relative error of ROM position (x) at t_final: 5 is 0.00188458 +// Relative error of ROM velocity (v) at t_final: 5 is 0.978726 +// // ================================================================================= // // Sample runs and results for parametric ROM using only displacement basis @@ -75,13 +82,20 @@ // Relative error of ROM position (x) at t_final: 5 is 7.08272e-05 // Relative error of ROM velocity (v) at t_final: 5 is 0.00387647 // -// Online phase with strong hyper reduction: +// Online phase with strong hyper reduction, using GNAT (over-sampled DEIM): // ./nonlinear_elasticity_global_rom -online -dt 0.01 -tf 5.0 -s 14 -vs 100 -hyp -hrtype gnat -rxdim 2 -hdim 4 -nsr 10 -sc 1.0 -xbo -def-ic // Output message: -// Elapsed time for time integration loop 1.01136 +// Elapsed time for time integration loop 0.120194 // Relative error of ROM position (x) at t_final: 5 is 0.0130818 // Relative error of ROM velocity (v) at t_final: 5 is 0.979978 // +// Online phase with strong hyper reduction, using QDEIM: +// ./nonlinear_elasticity_global_rom -online -dt 0.01 -tf 5.0 -s 14 -vs 100 -hyp -hrtype qdeim -rxdim 2 -hdim 4 -nsr 10 -sc 1.0 -xbo -def-ic +// Output message: +// Elapsed time for time integration loop 0.10806 +// Relative error of ROM position (x) at t_final: 5 is 0.0108709 +// Relative error of ROM velocity (v) at t_final: 5 is 1.30704 +// // This example runs in parallel with MPI, by using the same number of MPI ranks // in all phases (offline, merge, online). diff --git a/regression_tests/tests/prom/mixed_nonlinear_diffusion.sh b/regression_tests/tests/prom/mixed_nonlinear_diffusion.sh index d1a01f4a5..035ac6f32 100755 --- a/regression_tests/tests/prom/mixed_nonlinear_diffusion.sh +++ b/regression_tests/tests/prom/mixed_nonlinear_diffusion.sh @@ -4,6 +4,7 @@ CMDS=( "./mixed_nonlinear_diffusion -offline -tf 0.01" "./mixed_nonlinear_diffusion -merge -ns 1 -tf 0.01" "./mixed_nonlinear_diffusion -online -tf 0.01 -rrdim 2 -rwdim 2 -nldim 3 -hrtype deim" + "./mixed_nonlinear_diffusion -online -tf 0.01 -rrdim 2 -rwdim 2 -nldim 3 -hrtype qdeim" "./mixed_nonlinear_diffusion -online -tf 0.01 -rrdim 2 -rwdim 2 -nldim 3 -hrtype sopt" "./mixed_nonlinear_diffusion -online -tf 0.01 -rrdim 2 -rwdim 2 -nldim 3 -ns 1 -eqp -maxnnls 4" )