From 2e9de49f88105933ed7dba84670ed8ea79c0e967 Mon Sep 17 00:00:00 2001
From: Emily Boudreaux <thomas@boudreauxmail.com>
Date: Tue, 3 Jun 2025 08:11:01 -0400
Subject: [PATCH] refactor(poly): minor comment cleaning up

---
 src/poly/solver/private/polySolver.cpp      | 7 ++-----
 src/poly/utils/private/operator.cpp         | 9 +++++----
 utils/debugUtils/MeshQuality/meshQuality.py | 0
 3 files changed, 7 insertions(+), 9 deletions(-)
 create mode 100644 utils/debugUtils/MeshQuality/meshQuality.py

diff --git a/src/poly/solver/private/polySolver.cpp b/src/poly/solver/private/polySolver.cpp
index f04affe..74830fc 100644
--- a/src/poly/solver/private/polySolver.cpp
+++ b/src/poly/solver/private/polySolver.cpp
@@ -195,9 +195,7 @@ void PolySolver::solve() const {
     const solverBundle sb = setupNewtonSolver();
     sb.newton.Mult(zero_rhs, state_vector);
 
-    // Ax = b for x
-
-    // --- Save and view the solution ---
+    // --- Save and view an approximate 1D solution ---
     saveAndViewSolution(state_vector);
 }
 
@@ -405,8 +403,7 @@ solverBundle PolySolver::setupNewtonSolver() const {
     solver.solver.SetMaxIter(gmresMaxIter);
     solver.solver.SetPrintLevel(gmresPrintLevel);
 
-    // Preconditioner turned off because the polytrope operator seems *very* well conditioned without it
-    solver.solver.SetPreconditioner(m_polytropOperator->GetPreconditioner());
+    // solver.solver.SetPreconditioner(m_polytropOperator->GetPreconditioner());
     // --- Set up the Newton solver ---
     solver.newton.SetRelTol(newtonRelTol);
     solver.newton.SetAbsTol(newtonAbsTol);
diff --git a/src/poly/utils/private/operator.cpp b/src/poly/utils/private/operator.cpp
index f4f23c2..db092a0 100644
--- a/src/poly/utils/private/operator.cpp
+++ b/src/poly/utils/private/operator.cpp
@@ -154,14 +154,14 @@ void PolytropeOperator::Mult(const mfem::Vector &x, mfem::Vector &y) const {
     for (int i = 0; i < m_theta_ess_tdofs.first.Size(); i++) {
         if (int idx = m_theta_ess_tdofs.first[i]; idx >= 0 && idx < y_R0.Size()) {
             const double &targetValue = m_theta_ess_tdofs.second[i];
-            y_block.GetBlock(0)[idx] = x_theta(idx) - targetValue; // inhomogenous essential bc. This is commented out since seems it results in dramatic instabilities arising
+            y_block.GetBlock(0)[idx] = x_theta(idx) - targetValue; // inhomogenous essential bc.
         }
     }
 
     for (int i = 0; i < m_phi_ess_tdofs.first.Size(); i++) {
         if (int idx = m_phi_ess_tdofs.first[i]; idx >= 0 && idx < y_R1.Size()) {
             const double &targetValue = m_phi_ess_tdofs.second[i];
-            y_block.GetBlock(1)[idx] = x_phi(idx) - targetValue; // inhomogenous essential bc. This is commented out since seems it results in dramatic instabilities arising
+            y_block.GetBlock(1)[idx] = x_phi(idx) - targetValue; // inhomogenous essential bc.
         }
     }
 
@@ -240,8 +240,9 @@ GMRESInverter::GMRESInverter(const SchurCompliment &op) :
 mfem::Operator(op.Height(), op.Width()),
 m_op(op) {
     m_solver.SetOperator(m_op);
-    // PERF: It might be a good idea to turn down the total number of iterations and the tolerances
-    // PERF:   since we only need an approximation of the inverse
+    m_solver.SetMaxIter(100);
+    m_solver.SetRelTol(1e-1);
+    m_solver.SetAbsTol(1e-1);
 }
 
 void GMRESInverter::Mult(const mfem::Vector &x, mfem::Vector &y) const {
diff --git a/utils/debugUtils/MeshQuality/meshQuality.py b/utils/debugUtils/MeshQuality/meshQuality.py
new file mode 100644
index 0000000..e69de29