diff --git a/docs/src/config/model.md b/docs/src/config/model.md index 4db0750ab..1c3c3078d 100644 --- a/docs/src/config/model.md +++ b/docs/src/config/model.md @@ -64,7 +64,7 @@ scale based on the bounding box of the computational domain. with -`"Tol" [1e-2]` : Relative error convergence tolerance for adaptive mesh refinement (AMR). +`"Tol" [1.0e-2]` : Relative error convergence tolerance for adaptive mesh refinement (AMR). `"MaxIts" [0]` : Maximum number of iterations of AMR to perform. diff --git a/docs/src/config/solver.md b/docs/src/config/solver.md index 24c006c52..e079818bf 100644 --- a/docs/src/config/solver.md +++ b/docs/src/config/solver.md @@ -434,7 +434,7 @@ eigenmode simulation type. `"DivFreeMaxIts" [1000]` : Maximum number of iterations for divergence-free cleaning use in the eigenmode simulation type. -`"EstimatorTol" [1e-6]` : Relative tolerance for flux projection used in the +`"EstimatorTol" [1.0e-6]` : Relative tolerance for flux projection used in the error estimate calculation. `"EstimatorMaxIts" [10000]` : Maximum number of iterations for flux projection use in the @@ -465,5 +465,7 @@ vectors in Krylov subspace methods or other parts of the code. - `"STRUMPACKCompressionTol" [1.0e-3]` - `"STRUMPACKLossyPrecision" [16]` - `"STRUMPACKButterflyLevels" [1]` - - `"SuperLU3D" [false]` - - `"AMSVector" [false]` + - `"SuperLU3DCommunicator" [false]` + - `"AMSVectorInterpolation" [false]` + - `"AMSSingularOperator" [false]` + - `"AMGAggressiveCoarsening" [false]` diff --git a/palace/linalg/amg.cpp b/palace/linalg/amg.cpp index ad10fca1b..cbc8e4da7 100644 --- a/palace/linalg/amg.cpp +++ b/palace/linalg/amg.cpp @@ -6,7 +6,7 @@ namespace palace { -BoomerAmgSolver::BoomerAmgSolver(int cycle_it, int smooth_it, int print) +BoomerAmgSolver::BoomerAmgSolver(int cycle_it, int smooth_it, bool agg_coarsen, int print) : mfem::HypreBoomerAMG() { HYPRE_BoomerAMGSetPrintLevel(*this, (print > 1) ? print - 1 : 0); @@ -14,7 +14,7 @@ BoomerAmgSolver::BoomerAmgSolver(int cycle_it, int smooth_it, int print) HYPRE_BoomerAMGSetTol(*this, 0.0); // Set additional BoomerAMG options. - int agg_levels = 1; // Number of aggressive coarsening levels + int agg_levels = agg_coarsen ? 1 : 0; // Number of aggressive coarsening levels double theta = 0.5; // AMG strength parameter = 0.25 is 2D optimal (0.5-0.8 for 3D) int relax_type = 8; // 8 = l1-symm. GS, 13 = l1-GS, 18 = l1-Jacobi, 16 = Chebyshev if (mfem::Device::Allows(mfem::Backend::DEVICE_MASK)) diff --git a/palace/linalg/amg.hpp b/palace/linalg/amg.hpp index e02069b70..4fb09e6e5 100644 --- a/palace/linalg/amg.hpp +++ b/palace/linalg/amg.hpp @@ -16,10 +16,12 @@ namespace palace class BoomerAmgSolver : public mfem::HypreBoomerAMG { public: - BoomerAmgSolver(int cycle_it = 1, int smooth_it = 1, int print = 0); + BoomerAmgSolver(int cycle_it = 1, int smooth_it = 1, bool agg_coarsen = true, + int print = 0); BoomerAmgSolver(const IoData &iodata, bool coarse_solver, int print) : BoomerAmgSolver(coarse_solver ? 1 : iodata.solver.linear.mg_cycle_it, - iodata.solver.linear.mg_smooth_it, print) + iodata.solver.linear.mg_smooth_it, + iodata.solver.linear.amg_agg_coarsen, print) { } }; diff --git a/palace/linalg/ams.cpp b/palace/linalg/ams.cpp index b30462c59..a50464285 100644 --- a/palace/linalg/ams.cpp +++ b/palace/linalg/ams.cpp @@ -15,7 +15,8 @@ namespace palace HypreAmsSolver::HypreAmsSolver(FiniteElementSpace &nd_fespace, FiniteElementSpace &h1_fespace, int cycle_it, int smooth_it, - bool vector_interp, bool op_pos, bool op_singular, int print) + bool vector_interp, bool singular_op, bool agg_coarsen, + int print) : mfem::HypreSolver(), // From the Hypre docs for AMS: cycles 1, 5, 8, 11, 13 are fastest, 7 yields fewest its // (MFEM default is 13). 14 is similar to 11/13 but is cheaper in that is uses additive @@ -24,12 +25,12 @@ HypreAmsSolver::HypreAmsSolver(FiniteElementSpace &nd_fespace, // When used as the coarse solver of geometric multigrid, always do only a single // V-cycle. ams_it(cycle_it), ams_smooth_it(smooth_it), - // For positive (SPD) operators, we will use aggressive coarsening but not for frequency - // domain problems when the preconditioner matrix is not SPD. - ams_pos(op_pos), // If we know the operator is singular (no mass matrix, for magnetostatic problems), // internally the AMS solver will avoid G-space corrections. - ams_singular(op_singular), print((print > 1) ? print - 1 : 0) + ams_singular(singular_op), + // For positive (SPD) operators, we will use aggressive coarsening but not for frequency + // domain problems when the preconditioner matrix is not SPD. + agg_coarsen(agg_coarsen), print((print > 1) ? print - 1 : 0) { // From MFEM: The AMS preconditioner may sometimes require inverting singular matrices // with BoomerAMG, which are handled correctly in Hypre's Solve method, but can produce @@ -151,8 +152,8 @@ void HypreAmsSolver::InitializeSolver() } // Set additional AMS options. - int coarsen_type = 10; // 10 = HMIS, 8 = PMIS, 6 = Falgout, 0 = CLJP - int amg_agg_levels = ams_pos ? 1 : 0; // Number of aggressive coarsening levels + int coarsen_type = 10; // 10 = HMIS, 8 = PMIS, 6 = Falgout, 0 = CLJP + int amg_agg_levels = agg_coarsen ? 1 : 0; // Number of aggressive coarsening levels double theta = 0.5; // AMG strength parameter = 0.25 is 2D optimal (0.5-0.8 for 3D) int amg_relax_type = 8; // 3 = GS, 6 = symm. GS, 8 = l1-symm. GS, 13 = l1-GS, // 18 = l1-Jacobi, 16 = Chebyshev diff --git a/palace/linalg/ams.hpp b/palace/linalg/ams.hpp index 2fb96d5d0..a401fd5dc 100644 --- a/palace/linalg/ams.hpp +++ b/palace/linalg/ams.hpp @@ -25,7 +25,7 @@ class HypreAmsSolver : public mfem::HypreSolver // Parameters used for preconditioner construction. const int cycle_type, space_dim, ams_it, ams_smooth_it; - const bool ams_pos, ams_singular; + const bool ams_singular, agg_coarsen; // Control print level for debugging. const int print; @@ -49,16 +49,14 @@ class HypreAmsSolver : public mfem::HypreSolver // Constructor requires the ND space, but will construct the H1 and (H1)ᵈ spaces // internally as needed. HypreAmsSolver(FiniteElementSpace &nd_fespace, FiniteElementSpace &h1_fespace, - int cycle_it, int smooth_it, bool vector_interp, bool op_pos, - bool op_singular, int print); + int cycle_it, int smooth_it, bool vector_interp, bool singular_op, + bool agg_coarsen, int print); HypreAmsSolver(const IoData &iodata, bool coarse_solver, FiniteElementSpace &nd_fespace, FiniteElementSpace &h1_fespace, int print) : HypreAmsSolver( nd_fespace, h1_fespace, coarse_solver ? 1 : iodata.solver.linear.mg_cycle_it, - iodata.solver.linear.mg_smooth_it, iodata.solver.linear.ams_vector, - (iodata.problem.type == config::ProblemData::Type::TRANSIENT || - iodata.problem.type == config::ProblemData::Type::MAGNETOSTATIC), - (iodata.problem.type == config::ProblemData::Type::MAGNETOSTATIC), print) + iodata.solver.linear.mg_smooth_it, iodata.solver.linear.ams_vector_interp, + iodata.solver.linear.ams_singular_op, iodata.solver.linear.amg_agg_coarsen, print) { } ~HypreAmsSolver() override; diff --git a/palace/linalg/divfree.cpp b/palace/linalg/divfree.cpp index 5f9aa369d..8f7031da3 100644 --- a/palace/linalg/divfree.cpp +++ b/palace/linalg/divfree.cpp @@ -83,7 +83,7 @@ DivFreeSolver::DivFreeSolver( // The system matrix for the projection is real and SPD. auto amg = std::make_unique>( - std::make_unique(1, 1, 0)); + std::make_unique(1, 1, true, 0)); std::unique_ptr> pc; if (h1_fespaces.GetNumLevels() > 1) { diff --git a/palace/linalg/errorestimator.cpp b/palace/linalg/errorestimator.cpp index 4241c9705..bbf99fdc4 100644 --- a/palace/linalg/errorestimator.cpp +++ b/palace/linalg/errorestimator.cpp @@ -76,7 +76,7 @@ auto ConfigureLinearSolver(const FiniteElementSpaceHierarchy &fespaces, double t } else { - auto amg = std::make_unique(1, 1, 0); + auto amg = std::make_unique(1, 1, true, 0); amg->SetStrengthThresh(0.8); // More coarsening to save memory if (fespaces.GetNumLevels() > 1) { diff --git a/palace/utils/configfile.cpp b/palace/utils/configfile.cpp index dfc64ce64..6b3c3155d 100644 --- a/palace/utils/configfile.cpp +++ b/palace/utils/configfile.cpp @@ -1705,8 +1705,10 @@ void LinearSolverData::SetUp(json &solver) strumpack_lossy_precision = linear->value("STRUMPACKLossyPrecision", strumpack_lossy_precision); strumpack_butterfly_l = linear->value("STRUMPACKButterflyLevels", strumpack_butterfly_l); - superlu_3d = linear->value("SuperLU3D", superlu_3d); - ams_vector = linear->value("AMSVector", ams_vector); + superlu_3d = linear->value("SuperLU3DCommunicator", superlu_3d); + ams_vector_interp = linear->value("AMSVectorInterpolation", ams_vector_interp); + ams_singular_op = linear->value("AMSSingularOperator", ams_singular_op); + amg_agg_coarsen = linear->value("AMGAggressiveCoarsening", amg_agg_coarsen); // Other linear solver options. divfree_tol = linear->value("DivFreeTol", divfree_tol); @@ -1743,8 +1745,10 @@ void LinearSolverData::SetUp(json &solver) linear->erase("STRUMPACKCompressionTol"); linear->erase("STRUMPACKLossyPrecision"); linear->erase("STRUMPACKButterflyLevels"); - linear->erase("SuperLU3D"); - linear->erase("AMSVector"); + linear->erase("SuperLU3DCommunicator"); + linear->erase("AMSVectorInterpolation"); + linear->erase("AMSSingularOperator"); + linear->erase("AMGAggressiveCoarsening"); linear->erase("DivFreeTol"); linear->erase("DivFreeMaxIts"); @@ -1785,8 +1789,11 @@ void LinearSolverData::SetUp(json &solver) std::cout << "STRUMPACKCompressionTol: " << strumpack_lr_tol << '\n'; std::cout << "STRUMPACKLossyPrecision: " << strumpack_lossy_precision << '\n'; std::cout << "STRUMPACKButterflyLevels: " << strumpack_butterfly_l << '\n'; - std::cout << "SuperLU3D: " << superlu_3d << '\n'; - std::cout << "AMSVector: " << ams_vector << '\n'; + std::cout << "SuperLU3DCommunicator: " << superlu_3d << '\n'; + std::cout << "AMSVectorInterpolation: " << ams_vector_interp << '\n'; + std::cout << "AMSSingularOperator: " << ams_singular_op << '\n'; + std::cout << "AMGAggressiveCoarsening: " << amg_agg_coarsen << '\n'; + std::cout << "DivFreeTol: " << divfree_tol << '\n'; std::cout << "DivFreeMaxIts: " << divfree_max_it << '\n'; std::cout << "EstimatorTol: " << estimator_tol << '\n'; diff --git a/palace/utils/configfile.hpp b/palace/utils/configfile.hpp index bbd783646..550165978 100644 --- a/palace/utils/configfile.hpp +++ b/palace/utils/configfile.hpp @@ -137,7 +137,7 @@ struct RefinementData { public: // Non-dimensional tolerance used to specify convergence of adaptive mesh refinement. - double tol = 1e-2; + double tol = 1.0e-2; // Maximum number of iterations to perform during adaptive mesh refinement. int max_it = 0; @@ -864,7 +864,15 @@ struct LinearSolverData bool superlu_3d = false; // Option to use vector or scalar Pi-space corrections for the AMS preconditioner. - bool ams_vector = false; + bool ams_vector_interp = false; + + // Option to tell the AMS solver that the operator is singular, like for magnetostatic + // problems. + int ams_singular_op = -1; + + // Option to use aggressive coarsening for Hypre AMG solves (with BoomerAMG or AMS). + // Typically use this when the operator is positive definite. + int amg_agg_coarsen = -1; // Relative tolerance for solving linear systems in divergence-free projector. double divfree_tol = 1.0e-12; diff --git a/palace/utils/iodata.cpp b/palace/utils/iodata.cpp index c4da2cd62..ef48debd1 100644 --- a/palace/utils/iodata.cpp +++ b/palace/utils/iodata.cpp @@ -415,6 +415,18 @@ void IoData::CheckConfiguration() { solver.linear.mg_smooth_order = std::max(2 * solver.order, 4); } + if (solver.linear.ams_singular_op < 0) + { + solver.linear.ams_singular_op = + (problem.type == config::ProblemData::Type::MAGNETOSTATIC); + } + if (solver.linear.amg_agg_coarsen < 0) + { + solver.linear.amg_agg_coarsen = + (problem.type == config::ProblemData::Type::ELECTROSTATIC || + problem.type == config::ProblemData::Type::MAGNETOSTATIC || + problem.type == config::ProblemData::Type::TRANSIENT); + } // Configure settings for quadrature rules and partial assembly. BilinearForm::pa_order_threshold = solver.pa_order_threshold; diff --git a/scripts/schema/config/solver.json b/scripts/schema/config/solver.json index 1db36129b..379d420cd 100644 --- a/scripts/schema/config/solver.json +++ b/scripts/schema/config/solver.json @@ -121,8 +121,10 @@ "STRUMPACKCompressionTol": { "type": "number", "minimum": 0.0 }, "STRUMPACKLossyPrecision": { "type": "integer", "minimum": 0 }, "STRUMPACKButterflyLevels": { "type": "integer", "minimum": 0 }, - "SuperLU3D": { "type": "boolean" }, - "AMSVector": { "type": "boolean" }, + "SuperLU3DCommunicator": { "type": "boolean" }, + "AMSVectorInterpolation": { "type": "boolean" }, + "AMSSingularOperator": { "type": "boolean" }, + "AMGAggressiveCoarsening": { "type": "boolean" }, "DivFreeTol": { "type": "number", "minimum": 0.0 }, "DivFreeMaxIts": { "type": "integer", "minimum": 0 }, "EstimatorTol": { "type": "number", "minimum": 0.0 },