From ff299f8ce76082eaa85c298f0c48025d594e0801 Mon Sep 17 00:00:00 2001
From: Emily Boudreaux <emily@boudreauxmail.com>
Date: Thu, 20 Feb 2025 15:36:46 -0500
Subject: [PATCH] fix(poly): changed lambda from fespace to scalar

previously I had a lagrangian multipliers at every element; however, we are enforcing a global constraint so there need only be one lagrangian multiplier
---
 src/poly/solver/private/polySolver.cpp | 27 +++++---------------------
 src/poly/solver/public/polySolver.h    |  1 -
 2 files changed, 5 insertions(+), 23 deletions(-)

diff --git a/src/poly/solver/private/polySolver.cpp b/src/poly/solver/private/polySolver.cpp
index 84365f4..2be3947 100644
--- a/src/poly/solver/private/polySolver.cpp
+++ b/src/poly/solver/private/polySolver.cpp
@@ -20,7 +20,6 @@ PolySolver::PolySolver(double n, double order)
   mesh(meshIO.GetMesh()),
   feCollection(std::make_unique<mfem::H1_FECollection>(order, mesh.SpaceDimension())),
   feSpace(std::make_unique<mfem::FiniteElementSpace>(&mesh, feCollection.get())),
-  lambdaFeSpace(std::make_unique<mfem::FiniteElementSpace>(&mesh, feCollection.get(), 1)), // Scalar space for lambda
   compositeIntegrator(std::make_unique<polyMFEMUtils::CompositeNonlinearIntegrator>()),
   nonlinearForm(std::make_unique<mfem::NonlinearForm>(feSpace.get())),
   C(std::make_unique<mfem::LinearForm>(feSpace.get())),
@@ -73,15 +72,7 @@ void PolySolver::solve(){
 
     // --- Combine DOFs (u and λ) into a single vector ---
     int lambdaDofOffset  = feSpace->GetTrueVSize(); // Get the size of θ space
-    int totalTrueDofs = lambdaDofOffset + lambdaFeSpace->GetTrueVSize();
-
-    std::cout << "Total True Dofs: " << totalTrueDofs << std::endl;
-    std::cout << "Lambda Dof Offset: " << lambdaDofOffset << std::endl;
-    std::cout << "Lambda True Dofs: " << lambdaFeSpace->GetTrueVSize() << std::endl;
-    std::cout << "U True Dofs: " << feSpace->GetTrueVSize() << std::endl;
-    if (totalTrueDofs != lambdaDofOffset + 1) {
-        throw std::runtime_error("The total number of true dofs is not equal to the sum of the lambda offset and the lambda space size");
-    }
+    int totalTrueDofs = lambdaDofOffset + 1;
 
     mfem::Vector U(totalTrueDofs);
     U = 0.0;
@@ -113,25 +104,17 @@ void PolySolver::solve(){
     // --- Create the RHS of the augmented system ---
     mfem::Vector B(totalTrueDofs);
     B = 0.0;
-
-    // Set the constraint value (∫η(r) dΩ) in the last entry of B
-    // This sets the last entry to 1.0, this mighht be a problem later on...
-    mfem::ConstantCoefficient one(1.0);
-    mfem::LinearForm constraint_rhs(lambdaFeSpace.get());
-    constraint_rhs.AddDomainIntegrator(
-        new mfem::DomainLFIntegrator(*gaussianCoeff)
-    );
-
-    constraint_rhs.Assemble();
-    B[lambdaDofOffset] = constraint_rhs(0); // Get that single value for the rhs. Only one value because it's a scalar space
-
+    B[lambdaDofOffset] = 1.0;
 
     // --- Solve the augmented system ---
     newtonSolver.Mult(B, U);
 
     // --- 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
     double lambda = U[lambdaDofOffset];
 
     std::cout << "λ = " << lambda << std::endl;
diff --git a/src/poly/solver/public/polySolver.h b/src/poly/solver/public/polySolver.h
index 3ef60c3..854b709 100644
--- a/src/poly/solver/public/polySolver.h
+++ b/src/poly/solver/public/polySolver.h
@@ -20,7 +20,6 @@ private:
     std::unique_ptr<mfem::H1_FECollection> feCollection;
 
     std::unique_ptr<mfem::FiniteElementSpace> feSpace;
-    std::unique_ptr<mfem::FiniteElementSpace> lambdaFeSpace;
 
     std::unique_ptr<polyMFEMUtils::CompositeNonlinearIntegrator> compositeIntegrator;