From 431a47b9c734ef4858b0ee8679c89c9d281a3cb9 Mon Sep 17 00:00:00 2001 From: Emily Boudreaux Date: Mon, 21 Apr 2025 08:54:59 -0400 Subject: [PATCH] refactor(poly): updated header guards to pragma once --- src/poly/coeff/public/polyCoeff.h | 8 ++------ src/poly/solver/private/polySolver.cpp | 18 ++++++++--------- src/poly/solver/public/polySolver.h | 19 ++++++++--------- src/poly/utils/private/integrators.cpp | 18 ++++------------- src/poly/utils/public/integrators.h | 28 +++++++++++--------------- src/poly/utils/public/operator.h | 6 +----- 6 files changed, 36 insertions(+), 61 deletions(-) diff --git a/src/poly/coeff/public/polyCoeff.h b/src/poly/coeff/public/polyCoeff.h index 8cea996..6a18746 100644 --- a/src/poly/coeff/public/polyCoeff.h +++ b/src/poly/coeff/public/polyCoeff.h @@ -19,9 +19,7 @@ // // *********************************************************************** */ -#ifndef POLYCOEFF_H -#define POLYCOEFF_H - +#pragma once #include "mfem.hpp" #include @@ -55,6 +53,4 @@ namespace polycoeff constexpr static int power[51] = {0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33, 34, 35, 36, 37, 38, 39, 40, 41, 42, 43, 44, 45, 46, 47, 48, 49, 50}; constexpr static double coeff[51] = {2.449260049071003e+00, 4.340987403815444e-01, 2.592462347475860e+00, -2.283794512491552e+01, 1.135535208066129e+02, -3.460093673587010e+02, 6.948139716983651e+02, -9.564368645158952e+02, 9.184546798891624e+02, -6.130308987803503e+02, 2.747228735318193e+02, -7.416795903196430e+01, 7.318638538580859e+00, 1.749441306797260e+00, -4.240148582456829e-01, -5.818809544982156e-02, 1.514877217199105e-02, 3.228634707578998e-03, -2.862524323980516e-04, -1.622486968261819e-04, -1.253644717104076e-05, 4.334945141292894e-06, 1.296452565229763e-06, 8.634802209209870e-08, -3.337511676486084e-08, -1.094796628367775e-08, -1.228178760540410e-09, 1.416744125622751e-10, 8.513777351265677e-11, 1.624582561364811e-11, 8.377207519041114e-13, -4.363812865112836e-13, -1.535862757816461e-13, -2.485085045037669e-14, -6.566281276491033e-16, 8.405047965853478e-16, 2.673804441025638e-16, 4.176337890285142e-17, 8.150073570140493e-19, -1.531673805257016e-18, -4.746996933716653e-19, -6.976825828390195e-20, 8.807513368604331e-22, 3.307739310180278e-21, 8.593260940093030e-22, 7.385969061093440e-23, -2.211130577977291e-23, -8.557291048388455e-24, -4.169359901215994e-25, 4.609379358875657e-25, -3.590870335035984e-26}; }; -} // namespace polyCoeff - -#endif // POLYCOEFF_H \ No newline at end of file +} // namespace polyCoeff \ No newline at end of file diff --git a/src/poly/solver/private/polySolver.cpp b/src/poly/solver/private/polySolver.cpp index 1ea037e..e3e71ad 100644 --- a/src/poly/solver/private/polySolver.cpp +++ b/src/poly/solver/private/polySolver.cpp @@ -18,23 +18,24 @@ // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA // // *********************************************************************** */ -#include "mfem.hpp" +#include "polySolver.h" #include #include #include #include -#include "polySolver.h" +#include "mfem.hpp" + #include "4DSTARTypes.h" #include "config.h" #include "integrators.h" +#include "mfem.hpp" #include "operator.h" #include "polyCoeff.h" #include "probe.h" #include "resourceManager.h" #include "resourceManagerTypes.h" - #include "quill/LogMacros.h" @@ -159,8 +160,7 @@ void PolySolver::assembleBlockSystem() { // --- Assemble the NonlinearForm (f) --- auto fform = std::make_unique(m_feTheta.get()); - mfem::ConstantCoefficient oneCoeff(1.0); - fform->AddDomainIntegrator(new polyMFEMUtils::NonlinearPowerIntegrator(oneCoeff, m_polytropicIndex)); + fform->AddDomainIntegrator(new polyMFEMUtils::NonlinearPowerIntegrator(m_polytropicIndex)); // TODO: Add essential boundary conditions to the nonlinear form // -- Build the BlockOperator -- @@ -183,8 +183,8 @@ void PolySolver::solve() const { // It's safer to get the offsets directly from the operator after finalization const mfem::Array& block_offsets = m_polytropOperator->GetBlockOffsets(); // Assuming a getter exists or accessing member if public/friend mfem::BlockVector state_vector(block_offsets); - state_vector.GetBlock(0) = *m_theta; - state_vector.GetBlock(1) = *m_phi; + state_vector.GetBlock(0) = static_cast(*m_theta); + state_vector.GetBlock(1) = static_cast(*m_phi); mfem::Vector zero_rhs(block_offsets.Last()); zero_rhs = 0.0; @@ -289,7 +289,7 @@ void PolySolver::saveAndViewSolution(const mfem::BlockVector& state_vector) cons } // --- Extract the Solution --- - if (bool write11DSolution = m_config.get("Poly:Output:1D:Save", true)) { + if (m_config.get("Poly:Output:1D:Save", true)) { auto solutionPath = m_config.get("Poly:Output:1D:Path", "polytropeSolution_1D.csv"); auto derivSolPath = "d" + solutionPath; @@ -340,7 +340,7 @@ solverBundle PolySolver::setupNewtonSolver() const { LoadSolverUserParams(newtonRelTol, newtonAbsTol, newtonMaxIter, newtonPrintLevel, gmresRelTol, gmresAbsTol, gmresMaxIter, gmresPrintLevel); - solverBundle solver; + solverBundle solver; // Use this solver bundle to ensure lifetime safety solver.solver.SetRelTol(gmresRelTol); solver.solver.SetAbsTol(gmresAbsTol); solver.solver.SetMaxIter(gmresMaxIter); diff --git a/src/poly/solver/public/polySolver.h b/src/poly/solver/public/polySolver.h index b0a1a32..60af2a4 100644 --- a/src/poly/solver/public/polySolver.h +++ b/src/poly/solver/public/polySolver.h @@ -1,5 +1,4 @@ -#ifndef POLYSOLVER_H -#define POLYSOLVER_H +#pragma once #include "mfem.hpp" #include @@ -21,8 +20,8 @@ namespace laneEmden { // Struct to persist lifetime of the linear and nonlinear solvers struct solverBundle { - mfem::GMRESSolver solver; - mfem::NewtonSolver newton; + mfem::GMRESSolver solver; // Must be first so it lives longer than the newton solver + mfem::NewtonSolver newton; // Must be second so that when it is destroyed the solver is still alive preventing a double delete }; class PolySolver { @@ -32,10 +31,10 @@ public: // Public methods void solve() const; - double getN() { return m_polytropicIndex; } - double getOrder() { return m_feOrder; } - mfem::Mesh* getMesh() { return m_mesh.get(); } - mfem::GridFunction& getSolution() { return *m_theta; } + double getN() const { return m_polytropicIndex; } + double getOrder() const { return m_feOrder; } + mfem::Mesh* getMesh() const { return m_mesh.get(); } + mfem::GridFunction& getSolution() const { return *m_theta; } private: // Private Attributes Config& m_config = Config::getInstance(); @@ -70,6 +69,4 @@ private: // Private methods void LoadSolverUserParams(double &newtonRelTol, double &newtonAbsTol, int &newtonMaxIter, int &newtonPrintLevel, double &gmresRelTol, double &gmresAbsTol, int &gmresMaxIter, int &gmresPrintLevel) const; -}; - -#endif // POLYSOLVER_H \ No newline at end of file +}; \ No newline at end of file diff --git a/src/poly/utils/private/integrators.cpp b/src/poly/utils/private/integrators.cpp index 7590604..a81dbf0 100644 --- a/src/poly/utils/private/integrators.cpp +++ b/src/poly/utils/private/integrators.cpp @@ -22,15 +22,10 @@ #include #include "integrators.h" -#include "debug.h" namespace polyMFEMUtils { - NonlinearPowerIntegrator::NonlinearPowerIntegrator( - mfem::Coefficient &coeff, - double n) : - - m_coeff(coeff), + NonlinearPowerIntegrator::NonlinearPowerIntegrator(const double n) : m_polytropicIndex(n) {} void NonlinearPowerIntegrator::AssembleElementVector( @@ -60,9 +55,7 @@ namespace polyMFEMUtils { double u_safe = std::max(u_val, 0.0); double u_nl = std::pow(u_safe, m_polytropicIndex); - // double coeff_val = m_coeff.Eval(Trans, ip); - double coeff_val = 1.0; - double x2_u_nl = coeff_val * u_nl; + double x2_u_nl = u_nl; for (int i = 0; i < dof; i++){ elvect(i) += shape(i) * x2_u_nl * weight; @@ -95,13 +88,10 @@ namespace polyMFEMUtils { for (int j = 0; j < dof; j++) { u_val += elfun(j) * shape(j); } - // double coeff_val = m_coeff.Eval(Trans, ip); - double coeff_val = 1.0; - - + // Calculate the Jacobian double u_safe = std::max(u_val, 0.0); - double d_u_nl = coeff_val * m_polytropicIndex * std::pow(u_safe, m_polytropicIndex - 1); + double d_u_nl = m_polytropicIndex * std::pow(u_safe, m_polytropicIndex - 1); double x2_d_u_nl = d_u_nl; for (int i = 0; i < dof; i++) { diff --git a/src/poly/utils/public/integrators.h b/src/poly/utils/public/integrators.h index b6be026..df6893d 100644 --- a/src/poly/utils/public/integrators.h +++ b/src/poly/utils/public/integrators.h @@ -18,8 +18,7 @@ // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA // // *********************************************************************** */ -#ifndef POLYMFEMUTILS_H -#define POLYMFEMUTILS_H +#pragma once #include "mfem.hpp" #include @@ -28,7 +27,7 @@ /** - * @file polyMFEMUtils.h + * @file integrators.h * @brief A collection of utilities for working with MFEM and solving the lane-emden equation. */ @@ -42,24 +41,18 @@ namespace polyMFEMUtils { * @brief A class for nonlinear power integrator. */ class NonlinearPowerIntegrator: public mfem::NonlinearFormIntegrator { - private: - Config& m_config = Config::getInstance(); - Probe::LogManager& m_logManager = Probe::LogManager::getInstance(); - quill::Logger* m_logger = m_logManager.getLogger("log"); - mfem::Coefficient &m_coeff; - double m_polytropicIndex; public: /** * @brief Constructor for NonlinearPowerIntegrator. - * + * * @param coeff The function coefficient. * @param n The polytropic index. */ - NonlinearPowerIntegrator(mfem::Coefficient &coeff, double n); + NonlinearPowerIntegrator(double n); /** * @brief Assembles the element vector. - * + * * @param el The finite element. * @param Trans The element transformation. * @param elfun The element function. @@ -68,14 +61,17 @@ namespace polyMFEMUtils { virtual void AssembleElementVector(const mfem::FiniteElement &el, mfem::ElementTransformation &Trans, const mfem::Vector &elfun, mfem::Vector &elvect) override; /** * @brief Assembles the element gradient. - * + * * @param el The finite element. * @param Trans The element transformation. * @param elfun The element function. * @param elmat The element matrix to be assembled. */ virtual void AssembleElementGrad (const mfem::FiniteElement &el, mfem::ElementTransformation &Trans, const mfem::Vector &elfun, mfem::DenseMatrix &elmat) override; + private: + Config& m_config = Config::getInstance(); + Probe::LogManager& m_logManager = Probe::LogManager::getInstance(); + quill::Logger* m_logger = m_logManager.getLogger("log"); + double m_polytropicIndex; }; -} // namespace polyMFEMUtils - -#endif // POLYMFEMUTILS_H \ No newline at end of file +} // namespace polyMFEMUtils \ No newline at end of file diff --git a/src/poly/utils/public/operator.h b/src/poly/utils/public/operator.h index 42f5797..00c9ee8 100644 --- a/src/poly/utils/public/operator.h +++ b/src/poly/utils/public/operator.h @@ -1,5 +1,4 @@ -#ifndef POLY_UTILS_OPERATOR_H -#define POLY_UTILS_OPERATOR_H +#pragma once #include "mfem.hpp" #include "4DSTARTypes.h" @@ -81,6 +80,3 @@ private: void updateInverseNonlinearJacobian(const mfem::Operator &grad) const; void updateInverseSchurCompliment() const; }; - - -#endif // POLY_UTILS_OPERATOR_H \ No newline at end of file