Minor cleanup

palace/drivers/basesolver.hpp

@@ -18,7 +18,7 @@ class ParMesh;
 namespace palace
 {
 
-struct ErrorIndicators;
+class ErrorIndicators;
 class IoData;
 class PostOperator;
 class Timer;

palace/drivers/drivensolver.cpp

@@ -119,7 +119,6 @@ ErrorIndicators DrivenSolver::SweepUniform(SpaceOperator &spaceop, PostOperator
   // simply by setting diagonal entries of the system matrix for the corresponding dofs.
   // Because the Dirichlet BC is always homogenous, no special elimination is required on
   // the RHS. Assemble the linear system for the initial frequency (so we can call
-
   // KspSolver::SetOperators). Compute everything at the first frequency step.
   auto K = spaceop.GetStiffnessMatrix<ComplexOperator>(Operator::DIAG_ONE);
   auto C = spaceop.GetDampingMatrix<ComplexOperator>(Operator::DIAG_ZERO);
@@ -149,7 +148,7 @@ ErrorIndicators DrivenSolver::SweepUniform(SpaceOperator &spaceop, PostOperator
   // Initialize structures for storing and reducing the results of error estimation.
   ErrorIndicators indicators(spaceop.GlobalTrueVSize());
   const auto error_reducer = ErrorReductionOperator();
-  auto update_error_indicators =
+  auto UpdateErrorIndicators =
       [&timer, &estimator, &indicators, &error_reducer, &postop](const auto &E)
   {
     auto ind = estimator(E);
@@ -163,7 +162,7 @@ ErrorIndicators DrivenSolver::SweepUniform(SpaceOperator &spaceop, PostOperator
   // Main frequency sweep loop.
   double omega = omega0;
   auto t0 = timer.Now();
-  for (int step = step0; step < nstep; ++step)
+  for (int step = step0; step < nstep; step++)
   {
     const double freq = iodata.DimensionalizeValue(IoData::ValueType::FREQUENCY, omega);
     Mpi::Print("\nIt {:d}/{:d}: ω/2π = {:.3e} GHz (elapsed time = {:.2e} s)\n", step + 1,
@@ -187,12 +186,10 @@ ErrorIndicators DrivenSolver::SweepUniform(SpaceOperator &spaceop, PostOperator
     Mpi::Print("\n");
     ksp.Mult(RHS, E);
     timer.solve_time += timer.Lap();
-
-    update_error_indicators(E);
+    UpdateErrorIndicators(E);
 
     // Compute B = -1/(iω) ∇ x E on the true dofs, and set the internal GridFunctions in
     // PostOperator for all postprocessing operations.
-
     double E_elec = 0.0, E_mag = 0.0;
     Curl->Mult(E, B);
     B *= -1.0 / (1i * omega);
@@ -220,7 +217,6 @@ ErrorIndicators DrivenSolver::SweepUniform(SpaceOperator &spaceop, PostOperator
     // Increment frequency.
     omega += delta_omega;
   }
-
   SaveMetadata(ksp);
   return indicators;
 }
@@ -286,7 +282,7 @@ ErrorIndicators DrivenSolver::SweepAdaptive(SpaceOperator &spaceop, PostOperator
   // the online stage.
   ErrorIndicators indicators(spaceop.GlobalTrueVSize());
   ErrorReductionOperator error_reducer;
-  auto update_error_indicators =
+  auto UpdateErrorIndicators =
       [&local_timer, &estimator, &indicators, &error_reducer](const auto &E)
   {
     error_reducer(indicators, estimator(E));
@@ -296,12 +292,12 @@ ErrorIndicators DrivenSolver::SweepAdaptive(SpaceOperator &spaceop, PostOperator
 
   prom.SolveHDM(omega0, E);  // Print matrix stats at first HDM solve
   local_timer.solve_time += local_timer.Lap();
-  update_error_indicators(E);
+  UpdateErrorIndicators(E);
   prom.AddHDMSample(omega0, E);
   local_timer.construct_time += local_timer.Lap();
   prom.SolveHDM(omega0 + (nstep - step0 - 1) * delta_omega, E);
   local_timer.solve_time += local_timer.Lap();
-  update_error_indicators(E);
+  UpdateErrorIndicators(E);
   prom.AddHDMSample(omega0 + (nstep - step0 - 1) * delta_omega, E);
   local_timer.construct_time += local_timer.Lap();
 
@@ -327,10 +323,10 @@ ErrorIndicators DrivenSolver::SweepAdaptive(SpaceOperator &spaceop, PostOperator
         max_error);
     prom.SolveHDM(omega_star, E);
     local_timer.solve_time += local_timer.Lap();
-    update_error_indicators(E);
+    UpdateErrorIndicators(E);
     prom.AddHDMSample(omega_star, E);
     local_timer.construct_time += local_timer.Lap();
-    ++iter;
+    iter++;
   }
   Mpi::Print("\nAdaptive sampling{} {:d} frequency samples:\n"
              " n = {:d}, error = {:.3e}, tol = {:.3e}\n",
@@ -397,10 +393,9 @@ ErrorIndicators DrivenSolver::SweepAdaptive(SpaceOperator &spaceop, PostOperator
     timer.postpro_time += timer.Lap() - (timer.io_time - io_time_prev);
 
     // Increment frequency.
-    ++step;
+    step++;
     omega += delta_omega;
   }
-
   return indicators;
 }
 

palace/drivers/drivensolver.hpp

@@ -19,7 +19,7 @@ namespace palace
 {
 
 class CurlFluxErrorEstimator;
-struct ErrorIndicators;
+class ErrorIndicators;
 class IoData;
 class LumpedPortOperator;
 class PostOperator;

palace/drivers/eigensolver.cpp

@@ -265,7 +265,7 @@ ErrorIndicators EigenSolver::Solve(const std::vector<std::unique_ptr<mfem::ParMe
   // Initialize structures for storing and reducing the results of error estimation.
   ErrorIndicators indicators(spaceop.GlobalTrueVSize());
   ErrorReductionOperator error_reducer;
-  auto update_error_indicators =
+  auto UpdateErrorIndicators =
       [&timer, &estimator, &indicators, &error_reducer, &postop](const auto &E)
   {
     auto ind = estimator(E);
@@ -300,19 +300,18 @@ ErrorIndicators EigenSolver::Solve(const std::vector<std::unique_ptr<mfem::ParMe
     // Compute B = -1/(iω) ∇ x E on the true dofs, and set the internal GridFunctions in
     // PostOperator for all postprocessing operations.
     eigen->GetEigenvector(i, E);
-
     Curl->Mult(E, B);
     B *= -1.0 / (1i * omega);
+    timer.postpro_time += timer.Lap();
 
     if (i < iodata.solver.eigenmode.n)
     {
       // Only update the error indicator for targetted modes.
-      update_error_indicators(E);
+      UpdateErrorIndicators(E);
     }
 
     postop.SetEGridFunction(E);
     postop.SetBGridFunction(B);
-
     postop.UpdatePorts(spaceop.GetLumpedPortOp(), omega.real());
 
     // Postprocess the mode.
@@ -321,7 +320,6 @@ ErrorIndicators EigenSolver::Solve(const std::vector<std::unique_ptr<mfem::ParMe
                 timer);
     timer.postpro_time += timer.Lap() - (timer.io_time - io_time_prev);
   }
-
   return indicators;
 }
 

palace/drivers/eigensolver.hpp

@@ -19,7 +19,7 @@ class ParMesh;
 namespace palace
 {
 
-struct ErrorIndicators;
+class ErrorIndicators;
 class IoData;
 class LumpedPortOperator;
 class PostOperator;

palace/drivers/electrostaticsolver.cpp

@@ -26,7 +26,6 @@ ElectrostaticSolver::Solve(const std::vector<std::unique_ptr<mfem::ParMesh>> &me
   // handled eliminating the rows and columns of the system matrix for the corresponding
   // dofs. The eliminated matrix is stored in order to construct the RHS vector for nonzero
   // prescribed BC values.
-
   timer.Lap();
   LaplaceOperator laplaceop(iodata, mesh);
   auto K = laplaceop.GetStiffnessMatrix();
@@ -72,16 +71,15 @@ ElectrostaticSolver::Solve(const std::vector<std::unique_ptr<mfem::ParMesh>> &me
     timer.postpro_time += timer.Lap();
 
     // Next terminal.
-    ++step;
+    step++;
   }
 
+  // Construct estimator and reducer for error indicators.
   GradFluxErrorEstimator estimator(iodata, laplaceop.GetMaterialOp(), mesh,
                                    laplaceop.GetH1Space());
-
-  // Evaluate error estimator and reduce over all
   auto indicators = ErrorIndicators(laplaceop.GlobalTrueVSize());
   ErrorReductionOperator error_reducer;
-  auto update_error_indicators =
+  auto UpdateErrorIndicators =
       [&timer, &estimator, &indicators, &error_reducer](const auto &V)
   {
     error_reducer(indicators, estimator(V));
@@ -90,7 +88,7 @@ ElectrostaticSolver::Solve(const std::vector<std::unique_ptr<mfem::ParMesh>> &me
   timer.est_construction_time += timer.Lap();
 
   // Compute and reduce the error indicators for each solution.
-  std::for_each(V.begin(), V.end(), update_error_indicators);
+  std::for_each(V.begin(), V.end(), UpdateErrorIndicators);
 
   // Register the indicator field used to drive the overall adaptation.
   postop.SetIndicatorGridFunction(indicators.local_error_indicators);
@@ -100,7 +98,6 @@ ElectrostaticSolver::Solve(const std::vector<std::unique_ptr<mfem::ParMesh>> &me
   SaveMetadata(ksp);
   Postprocess(laplaceop, postop, V, timer);
   timer.postpro_time += timer.Lap() - (timer.io_time - io_time_prev);
-
   return indicators;
 }
 
@@ -129,7 +126,6 @@ void ElectrostaticSolver::Postprocess(LaplaceOperator &laplaceop, PostOperator &
     // for all postprocessing operations.
     E = 0.0;
     Grad->AddMult(V[i], E, -1.0);
-
     postop.SetEGridFunction(E);
     postop.SetVGridFunction(V[i]);
     double Ue = postop.GetEFieldEnergy();
@@ -146,13 +142,13 @@ void ElectrostaticSolver::Postprocess(LaplaceOperator &laplaceop, PostOperator &
 
     // Diagonal: C_ii = 2 U_e(V_i) / V_i².
     C(i, i) = Cm(i, i) = 2.0 * Ue;
-    ++i;
+    i++;
   }
 
   // Off-diagonals: C_ij = U_e(V_i + V_j) / (V_i V_j) - 1/2 (V_i/V_j C_ii + V_j/V_i C_jj).
-  for (i = 0; i < C.Height(); ++i)
+  for (i = 0; i < C.Height(); i++)
   {
-    for (int j = 0; j < C.Width(); ++j)
+    for (int j = 0; j < C.Width(); j++)
     {
       if (j < i)
       {

palace/drivers/electrostaticsolver.hpp

@@ -23,7 +23,7 @@ class ParMesh;
 namespace palace
 {
 
-struct ErrorIndicators;
+class ErrorIndicators;
 class IoData;
 class LaplaceOperator;
 class PostOperator;

palace/drivers/magnetostaticsolver.cpp

@@ -71,7 +71,7 @@ MagnetostaticSolver::Solve(const std::vector<std::unique_ptr<mfem::ParMesh>> &me
     timer.postpro_time += timer.Lap();
 
     // Next source.
-    ++step;
+    step++;
   }
 
   CurlFluxErrorEstimator estimator(iodata, curlcurlop.GetMaterialOp(), mesh,
@@ -80,7 +80,7 @@ MagnetostaticSolver::Solve(const std::vector<std::unique_ptr<mfem::ParMesh>> &me
   // Initialize structures for storing and reducing the results of error estimation.
   ErrorIndicators indicators(curlcurlop.GlobalTrueVSize());
   ErrorReductionOperator error_reducer;
-  auto update_error_indicators =
+  auto UpdateErrorIndicators =
       [&timer, &estimator, &indicators, &error_reducer](const auto &A)
   {
     error_reducer(indicators, estimator(A));
@@ -89,7 +89,7 @@ MagnetostaticSolver::Solve(const std::vector<std::unique_ptr<mfem::ParMesh>> &me
   timer.est_construction_time += timer.Lap();
 
   // Compute and reduce the error indicators for each solution.
-  std::for_each(A.begin(), A.end(), update_error_indicators);
+  std::for_each(A.begin(), A.end(), UpdateErrorIndicators);
 
   // Register the indicator field used to drive the overall adaptation.
   postop.SetIndicatorGridFunction(indicators.local_error_indicators);
@@ -99,7 +99,6 @@ MagnetostaticSolver::Solve(const std::vector<std::unique_ptr<mfem::ParMesh>> &me
   SaveMetadata(ksp);
   Postprocess(curlcurlop, postop, A, timer);
   timer.postpro_time += timer.Lap() - (timer.io_time - io_time_prev);
-
   return indicators;
 }
 
@@ -150,13 +149,13 @@ void MagnetostaticSolver::Postprocess(CurlCurlOperator &curlcurlop, PostOperator
 
     // Diagonal: M_ii = 2 U_m(A_i) / I_i².
     M(i, i) = Mm(i, i) = 2.0 * Um / (Iinc(i) * Iinc(i));
-    ++i;
+    i++;
   }
 
   // Off-diagonals: M_ij = U_m(A_i + A_j) / (I_i I_j) - 1/2 (I_i/I_j M_ii + I_j/I_i M_jj).
-  for (i = 0; i < M.Height(); ++i)
+  for (i = 0; i < M.Height(); i++)
   {
-    for (int j = 0; j < M.Width(); ++j)
+    for (int j = 0; j < M.Width(); j++)
     {
       if (j < i)
       {
@@ -223,10 +222,10 @@ void MagnetostaticSolver::PostprocessTerminals(const SurfaceCurrentOperator &sur
                      mat(i, j) * scale, table.w, table.p,
                      (idx2 == surf_j_op.rbegin()->first) ? "" : ",");
         // clang-format on
-        ++j;
+        j++;
       }
       output.print("\n");
-      ++i;
+      i++;
     }
   };
   const double H = iodata.DimensionalizeValue(IoData::ValueType::INDUCTANCE, 1.0);

palace/drivers/magnetostaticsolver.hpp

@@ -21,7 +21,7 @@ namespace palace
 {
 
 class CurlCurlOperator;
-struct ErrorIndicators;
+class ErrorIndicators;
 class IoData;
 class PostOperator;
 class SurfaceCurrentOperator;

palace/drivers/transientsolver.cpp

@@ -127,7 +127,6 @@ TransientSolver::Solve(const std::vector<std::unique_ptr<mfem::ParMesh>> &mesh,
     // Increment time step.
     step++;
   }
-
   SaveMetadata(timeop.GetLinearSolver());
   return ErrorIndicators(spaceop.GlobalTrueVSize());
 }

palace/drivers/transientsolver.hpp

@@ -19,7 +19,7 @@ class ParMesh;
 namespace palace
 {
 
-struct ErrorIndicators;
+class ErrorIndicators;
 class IoData;
 class LumpedPortOperator;
 class PostOperator;

palace/fem/coefficient.hpp

@@ -499,7 +499,6 @@ class CurlFluxCoefficient : public mfem::VectorCoefficient
             const mfem::IntegrationPoint &ip) override
   {
     V.SetSize(3);
-
     coef.Eval(muinv, T, ip);
     X.GetCurl(T, curl);
     muinv.Mult(curl, V);

palace/fem/multigrid.hpp

@@ -127,22 +127,21 @@ template <typename FECollection>
 inline mfem::ParFiniteElementSpaceHierarchy ConstructFiniteElementSpaceHierarchy(
     int mg_max_levels, bool mg_legacy_transfer, int pa_order_threshold,
     const std::vector<std::unique_ptr<mfem::ParMesh>> &mesh,
-    const std::vector<std::unique_ptr<FECollection>> &fecs, int dim = 1,
-    int ordering = mfem::Ordering::byNODES)
+    const std::vector<std::unique_ptr<FECollection>> &fecs, int dim = 1)
 {
   MFEM_VERIFY(!mesh.empty() && !fecs.empty(),
               "Empty mesh or FE collection for FE space construction!");
   int coarse_mesh_l =
       std::max(0, static_cast<int>(mesh.size() + fecs.size()) - 1 - mg_max_levels);
-  auto *fespace = new mfem::ParFiniteElementSpace(mesh[coarse_mesh_l].get(), fecs[0].get(),
-                                                  dim, ordering);
+  auto *fespace =
+      new mfem::ParFiniteElementSpace(mesh[coarse_mesh_l].get(), fecs[0].get(), dim);
   mfem::ParFiniteElementSpaceHierarchy fespaces(mesh[coarse_mesh_l].get(), fespace, false,
                                                 true);
 
   // h-refinement
   for (std::size_t l = coarse_mesh_l + 1; l < mesh.size(); l++)
   {
-    fespace = new mfem::ParFiniteElementSpace(mesh[l].get(), fecs[0].get(), dim, ordering);
+    fespace = new mfem::ParFiniteElementSpace(mesh[l].get(), fecs[0].get(), dim);
 
     auto *P = new ParOperator(
         std::make_unique<mfem::TransferOperator>(fespaces.GetFinestFESpace(), *fespace),
@@ -153,8 +152,7 @@ inline mfem::ParFiniteElementSpaceHierarchy ConstructFiniteElementSpaceHierarchy
   // p-refinement
   for (std::size_t l = 1; l < fecs.size(); l++)
   {
-    fespace =
-        new mfem::ParFiniteElementSpace(mesh.back().get(), fecs[l].get(), dim, ordering);
+    fespace = new mfem::ParFiniteElementSpace(mesh.back().get(), fecs[l].get(), dim);
 
     ParOperator *P;
     if (!mg_legacy_transfer && mfem::DeviceCanUseCeed())
@@ -185,11 +183,11 @@ mfem::ParFiniteElementSpaceHierarchy ConstructFiniteElementSpaceHierarchy(
     int mg_max_levels, bool mg_legacy_transfer, int pa_order_threshold,
     const Container<MeshT...> &mesh, const std::vector<std::unique_ptr<FECollection>> &fecs,
     const mfem::Array<int> &dbc_marker, std::vector<mfem::Array<int>> &dbc_tdof_lists,
-    int dim = 1, int ordering = mfem::Ordering::byNODES)
+    int dim = 1)
 {
   dbc_tdof_lists.clear();
-  auto fespaces = ConstructFiniteElementSpaceHierarchy(
-      mg_max_levels, mg_legacy_transfer, pa_order_threshold, mesh, fecs, dim, ordering);
+  auto fespaces = ConstructFiniteElementSpaceHierarchy(mg_max_levels, mg_legacy_transfer,
+                                                       pa_order_threshold, mesh, fecs, dim);
   for (int l = 0; l < fespaces.GetNumLevels(); l++)
   {
     fespaces.GetFESpaceAtLevel(l).GetEssentialTrueDofs(dbc_marker,