fix(poly): polytrope converges to solution

first commit where the polytrope module converges to a solution. I have not yet validated if it is a correct solution
This commit is contained in:
2025-02-20 16:05:02 -05:00
parent 776174c093
commit 9925f56e34
3 changed files with 22 additions and 7 deletions

View File

@@ -9,6 +9,8 @@
#include "polyMFEMUtils.h"
#include "polyCoeff.h"
#include "warning_control.h"
// TODO: Come back to this and think of a better way to get the mesh file
const std::string SPHERICAL_MESH = std::string(getenv("MESON_SOURCE_ROOT")) + "/src/resources/mesh/sphere.msh";
@@ -111,10 +113,11 @@ void PolySolver::solve(){
// --- Extract the Solution ---
mfem::Vector u_sol_view(U.GetData(), lambdaDofOffset);
#pragma GCC diagnostic push // MFEM is using deprecated functions, I cant do anything about it so I will ignore the warning
#pragma GCC diagnostic ignored "-Wdeprecated-declarations"
u->SetData(u_sol_view);
#pragma GCC diagnostic pop
DEPRECATION_WARNING_OFF // DISABLE DEPRECATION WARNING
u->SetData(u_sol_view);
DEPRECATION_WARNING_ON // REENABLE DEPRECATION WARNING
double lambda = U[lambdaDofOffset];
std::cout << "λ = " << lambda << std::endl;

View File

@@ -6,6 +6,8 @@
#include "polyMFEMUtils.h"
#include "warning_control.h"
namespace polyMFEMUtils {
NonlinearPowerIntegrator::NonlinearPowerIntegrator(
mfem::Coefficient &coeff,
@@ -242,7 +244,11 @@ namespace polyMFEMUtils {
y[i] = F[i];
}
mfem::GridFunction u_gf(C.FESpace());
u_gf.SetData(u);
DEPRECATION_WARNING_OFF
u_gf.SetData(u);
DEPRECATION_WARNING_ON
y[lambdaDofOffset] = C.operator()(u_gf);
// add -lambda * C to the residual
@@ -275,8 +281,13 @@ namespace polyMFEMUtils {
mfem::SparseMatrix *J_aug = new mfem::SparseMatrix(height, width);
// Copy the original Jacobian into the augmented Jacobian
for (int i = 0; i < Jnfl_sparse->Height(); i++) {
for (int j = 0; j < Jnfl_sparse->Width(); j++) {
J_aug->Set(i, j, Jnfl_sparse->Elem(i, j));
const int *J_cols = Jnfl_sparse->GetRowColumns(i);
const double *J_vals = Jnfl_sparse->GetRowEntries(i);
for (int jj = 0; jj < Jnfl_sparse->RowSize(i); jj++) {
int j = J_cols[jj];
double val = J_vals[jj];
J_aug->Set(i, j, val);
}
}

View File

@@ -5,6 +5,7 @@
#include <string>
/**
* @file polyMFEMUtils.h
* @brief A collection of utilities for working with MFEM and solving the lane-emden equation.