Merge 2D magnetostatics branch into serac #63

tuckerbabcock · 2022-11-14T19:01:48Z

This merges all the changes I've made to the magnetostatic solver to make it support evaluation on 2D meshes. Includes all the EM outputs needed for a motor problem, with all of them differentiated.

Items of potentially broad interest:

RelaxedNewtonSolver: a newton solver that supports various LineSearch methods. Currently only have a backtracking linesearch implemented that fits a quadratic or cubic and calculates the step as the step that minimizes the fit polynomial.
Further developed support for reverse and forward mode output and residual sensitivities
Maybe more?

… the 3d version, and can actually converge a nonlinear residual

…n with 3D magnets

… 3D. Integrator approach follows a more straightforward AD approach than the previous integrator, using hand written AD of the MagneticEnergy integrator

…sh and test both

…d 3D. Once it is differentiated, I will rename it to ForceIntegrator and remove the old one. Set it up so that magnetic energy could be computed on a per attribute basis, so you could calculate just the energy in the airgap for example. Now only add the current density integrator to the winding attributes so that it hopefully runs more efficiently

…lize setting attributes

…fully easier to differentiate than ForceIntegrator2

…is differentiated wrt the mesh and the state

…omag integs

…n but didn't end up working)

…make sure differentiation of the 2D cases works correctly

…nce. Added LineSearch abstract class, with a BacktrackingLineSearch implementation. Added tests for RelaxedNewton using BacktrackingLineSearch on a simple 2D Rosenbrock minimization problem

…n the field named state, and differentiated it wrt the input state and the mesh

…ption

…tonal output, and included tests for it

…emoved try-catch from AbstractSolvers calcOutput and output sensitivities methods. Previously they would catch and print the what of any exception thrown, now they will not intercept. This is to allow a NotImplementedException to be thrown further up the chain

…ions to be more efficient

…lly on output sensitivity verification

…er residual

mach/mach_output.py

mach/mach_state.py

src/common/mach_residual.hpp

tuckerbabcock · 2022-11-14T19:35:08Z

src/physics/fluidflow/flow_control_residual.cpp

   double closed_double = 1.0;
   setValueFromInputs(inputs, "closed-loop", closed_double);
-   if (fabs(closed_double) < numeric_limits<double>::epsilon())
-   {
-      residual.closed_loop = false;
-   }
-   else
-   {
-      residual.closed_loop = true;
-   }
+   residual.closed_loop =
+       fabs(closed_double) >= numeric_limits<double>::epsilon();


@jehicken I think this should really be an option that is set for the flow control residual, then it can have a typed true/false value and be more expressive.

Agreed, I'm not sure why it is a data member. Would this be easy for you to fix? I only ask because I won't have time until summer to implement the change. Or, we could ignore it for now.

github-actions

⚠️ clang-tidy found issue(s) with the introduced code (1/8)

github-actions · 2022-11-14T20:33:30Z

src/physics/mfem_common_integ.hpp

@@ -1,10 +1,12 @@
 #ifndef MACH_MFEM_COMMON_INTEG
 #define MACH_MFEM_COMMON_INTEG

+#include <linalg/vector.hpp>


⚠️ clang-diagnostic-error ⚠️
linalg/vector.hpp file not found

github-actions · 2022-11-14T20:33:30Z

test/unit/test_mach_load.cpp

@@ -83,7 +83,7 @@ TEST_CASE("MachInputs Scalar Input Test",
   MachLinearForm lf(fes, fields);
   lf.addDomainIntegrator(new TestMachLoadIntegrator);

-   MachLoad ml(lf);
+   MachLoad ml(std::move(lf));


⚠️ clang-diagnostic-error ⚠️
call to implicitly-deleted copy constructor of mach::MachLinearForm

@tuckerbabcock Should we be concerned about this?

github-actions · 2022-11-14T20:33:30Z

src/physics/electromagnetics/electromag_integ.hpp

+   /// \param[in] pf_trans - the transformation between reference and physical
+   /// space
+   /// \param[out] peak_flux_bar - dJdB for the element
+   void AssembleRHSElementVect(const mfem::FiniteElement &pf_el,


⚠️ readability-inconsistent-declaration-parameter-name ⚠️
function mach::ACLossFunctionalIntegratorPeakFluxSens::AssembleRHSElementVect has a definition with different parameter names

@tuckerbabcock

github-actions · 2022-11-14T20:33:31Z

src/physics/electromagnetics/electromag_integ.hpp

+   {
+      output_scalar_sens.at("frequency")
+          .AddDomainIntegrator(
+              new SteinmetzLossIntegratorFreqSens(primal_integ), *attr_marker);


⚠️ clang-diagnostic-error ⚠️
too many arguments to function call, expected single argument nlfi, have 2 arguments

github-actions · 2022-11-14T20:33:31Z

src/physics/electromagnetics/electromag_integ.hpp

+      output_scalar_sens.at("max_flux_magnitude")
+          .AddDomainIntegrator(
+              new SteinmetzLossIntegratorMaxFluxSens(primal_integ),
+              *attr_marker);


⚠️ clang-diagnostic-error ⚠️
too many arguments to function call, expected single argument nlfi, have 2 arguments

github-actions · 2022-11-14T20:33:34Z

src/physics/fluidflow/euler_integ.hpp

@@ -803,15 +805,15 @@
   BoundaryEntropy(adept::Stack &diff_stack,
                   const mfem::FiniteElementCollection *fe_coll,
                   BCScaleFun scale,
-                   const mfem::Vector &xc,
+                   mfem::Vector xc,


⚠️ performance-unnecessary-value-param ⚠️
the parameter xc is copied for each invocation but only used as a const reference; consider making it a const reference

Suggested change

mfem::Vector xc,

const mfem::Vector xc,

github-actions · 2022-11-14T20:33:34Z

src/physics/fluidflow/euler_integ.hpp

@@ -803,15 +805,15 @@
   BoundaryEntropy(adept::Stack &diff_stack,
                   const mfem::FiniteElementCollection *fe_coll,
                   BCScaleFun scale,
-                   const mfem::Vector &xc,
+                   mfem::Vector xc,


⚠️ performance-unnecessary-value-param ⚠️
the parameter xc is copied for each invocation but only used as a const reference; consider making it a const reference

Suggested change

mfem::Vector xc,

mfem::Vector& xc,

github-actions · 2022-11-14T20:33:34Z

src/physics/fluidflow/euler_integ.hpp

                   double len = 1.0)
    : InviscidBoundaryIntegrator<BoundaryEntropy<dim, entvar>>(diff_stack,
                                                               fe_coll,
                                                               dim + 2,
                                                               1.0),
      len_scale(len),
-      x_actuator(xc),
-      control_scale(scale)
+      x_actuator(std::move(xc)),


⚠️ performance-move-const-arg ⚠️
passing result of std::move() as a const reference argument; no move will actually happen

Suggested change

x_actuator(std::move(xc)),

x_actuator(xc)),

github-actions · 2022-11-14T20:33:34Z

src/physics/fluidflow/euler_integ.hpp

                   double len = 1.0)
    : InviscidBoundaryIntegrator<BoundaryEntropy<dim, entvar>>(diff_stack,
                                                               fe_coll,
                                                               dim + 2,
                                                               1.0),
      len_scale(len),
-      x_actuator(xc),
-      control_scale(scale)
+      x_actuator(std::move(xc)),


⚠️ performance-move-const-arg ⚠️
passing result of std::move() as a const reference argument; no move will actually happen

Suggested change

x_actuator(std::move(xc)),

x_actuator(std::move(xc),

github-actions · 2022-11-14T20:33:35Z

test/unit/test_linesearch.cpp

+
+   mach::RelaxedNewton newton(MPI_COMM_SELF, newton_opts);
+   // mfem::NewtonSolver newton(MPI_COMM_SELF);
+   newton.SetPrintLevel(mfem::IterativeSolver::PrintLevel().All());


⚠️ clang-diagnostic-error ⚠️
no member named PrintLevel in mfem::IterativeSolver

github-actions

⚠️ clang-tidy found issue(s) with the introduced code (2/8)

github-actions · 2022-11-14T20:33:51Z

src/physics/common_outputs.cpp

@@ -1,6 +1,7 @@
 #include <string>
 #include <unordered_map>

+#include "mach_load.hpp"
 #include "mfem.hpp"


⚠️ cppcoreguidelines-init-variables ⚠️
variable state is not initialized

Suggested change

#include "mfem.hpp"

#include "math.h"

#include "mfem.hpp"

github-actions · 2022-11-14T20:33:51Z

src/physics/common_outputs.cpp

+                             const std::string &wrt)
+{
+   const MachInputs &inputs = *output.inputs;
+   double state = calcOutput(output.state_integ, inputs);


⚠️ cppcoreguidelines-init-variables ⚠️
variable state is not initialized

Suggested change

double state = calcOutput(output.state_integ, inputs);

double state = NAN = calcOutput(output.state_integ, inputs);

github-actions · 2022-11-14T20:33:51Z

src/physics/common_outputs.cpp

+{
+   const MachInputs &inputs = *output.inputs;
+   double state = calcOutput(output.state_integ, inputs);
+   double volume = calcOutput(output.volume, inputs);


⚠️ cppcoreguidelines-init-variables ⚠️
variable volume is not initialized

Suggested change

double volume = calcOutput(output.volume, inputs);

double volume = NAN = calcOutput(output.volume, inputs);

github-actions · 2022-11-14T20:33:52Z

src/physics/common_outputs.cpp

+                           mfem::Vector &wrt_bar)
+{
+   const MachInputs &inputs = *output.inputs;
+   double state = calcOutput(output.state_integ, inputs);


⚠️ cppcoreguidelines-init-variables ⚠️
variable state is not initialized

Suggested change

double state = calcOutput(output.state_integ, inputs);

double state = NAN = calcOutput(output.state_integ, inputs);

github-actions · 2022-11-14T20:33:52Z

src/physics/common_outputs.cpp

+{
+   const MachInputs &inputs = *output.inputs;
+   double state = calcOutput(output.state_integ, inputs);
+   double volume = calcOutput(output.volume, inputs);


⚠️ cppcoreguidelines-init-variables ⚠️
variable volume is not initialized

Suggested change

double volume = calcOutput(output.volume, inputs);

double volume = NAN = calcOutput(output.volume, inputs);

github-actions · 2022-11-14T20:33:54Z

src/physics/electromagnetics/electromag_outputs.cpp

+          output.wire_length / (strand_area * output.strands_in_hand) * rho_dot;
+
+      double loss = pow(output.rms_current, 2) * R * sqrt(2);
+      double loss_dot = pow(output.rms_current, 2) * sqrt(2) * R_dot;


⚠️ cppcoreguidelines-init-variables ⚠️
variable loss_dot is not initialized

Suggested change

double loss_dot = pow(output.rms_current, 2) * sqrt(2) * R_dot;

double loss_dot = NAN = pow(output.rms_current, 2) * sqrt(2) * R_dot;

github-actions · 2022-11-14T20:33:55Z

src/physics/electromagnetics/electromag_outputs.cpp

+      /// double loss = pow(output.rms_current, 2) * R * sqrt(2);
+      // double rms_current_bar = loss_bar * 2 * output.rms_current * R *
+      // sqrt(2);
+      double R_bar = loss_bar * pow(output.rms_current, 2) * sqrt(2);


⚠️ cppcoreguidelines-init-variables ⚠️
variable R_bar is not initialized

Suggested change

double R_bar = loss_bar * pow(output.rms_current, 2) * sqrt(2);

double R_bar = NAN = loss_bar * pow(output.rms_current, 2) * sqrt(2);

github-actions · 2022-11-14T20:33:55Z

src/physics/electromagnetics/electromag_outputs.cpp

+      // volume does not depend on any of the inputs except mesh coords
+
+      /// double loss = pow(output.rms_current, 2) * R * sqrt(2);
+      double rms_current_bar = loss_bar * 2 * output.rms_current * R * sqrt(2);


⚠️ cppcoreguidelines-init-variables ⚠️
variable rms_current_bar is not initialized

Suggested change

double rms_current_bar = loss_bar * 2 * output.rms_current * R * sqrt(2);

double rms_current_bar = NAN = loss_bar * 2 * output.rms_current * R * sqrt(2);

github-actions · 2022-11-14T20:33:55Z

src/physics/electromagnetics/electromag_outputs.cpp

+      /// double loss = pow(output.rms_current, 2) * R * sqrt(2);
+      // double rms_current_bar = loss_bar * 2 * output.rms_current * R *
+      // sqrt(2);
+      double R_bar = loss_bar * pow(output.rms_current, 2) * sqrt(2);


⚠️ cppcoreguidelines-init-variables ⚠️
variable R_bar is not initialized

Suggested change

double R_bar = loss_bar * pow(output.rms_current, 2) * sqrt(2);

double R_bar = NAN = loss_bar * pow(output.rms_current, 2) * sqrt(2);

github-actions · 2022-11-14T20:33:55Z

src/physics/electromagnetics/electromag_outputs.cpp

+      /// double loss = pow(output.rms_current, 2) * R * sqrt(2);
+      // double rms_current_bar = loss_bar * 2 * output.rms_current * R *
+      // sqrt(2);
+      double R_bar = loss_bar * pow(output.rms_current, 2) * sqrt(2);


⚠️ cppcoreguidelines-init-variables ⚠️
variable R_bar is not initialized

Suggested change

double R_bar = loss_bar * pow(output.rms_current, 2) * sqrt(2);

double R_bar = NAN = loss_bar * pow(output.rms_current, 2) * sqrt(2);

github-actions

⚠️ clang-tidy found issue(s) with the introduced code (3/8)

github-actions · 2022-11-14T20:34:12Z

src/physics/electromagnetics/electromag_outputs.cpp

+      double R =
+          output.wire_length * rho / (strand_area * output.strands_in_hand);
+
+      double loss = pow(output.rms_current, 2) * R * sqrt(2);


⚠️ cppcoreguidelines-init-variables ⚠️
variable loss is not initialized

Suggested change

double loss = pow(output.rms_current, 2) * R * sqrt(2);

double loss = NAN = pow(output.rms_current, 2) * R * sqrt(2);

github-actions · 2022-11-14T20:34:13Z

src/physics/electromagnetics/electromag_outputs.cpp

+      /// double loss = pow(output.rms_current, 2) * R * sqrt(2);
+      // double rms_current_bar =
+      //     loss_bar * 2 * output.rms_current * R * sqrt(2);
+      double R_bar = loss_bar * pow(output.rms_current, 2) * sqrt(2);


⚠️ cppcoreguidelines-init-variables ⚠️
variable R_bar is not initialized

Suggested change

double R_bar = loss_bar * pow(output.rms_current, 2) * sqrt(2);

double R_bar = NAN = loss_bar * pow(output.rms_current, 2) * sqrt(2);

github-actions · 2022-11-14T20:34:13Z

src/common/inexact_newton.cpp

@@ -130,14 +130,14 @@ void InexactNewton::Mult(const Vector &b, Vector &x) const

      if (norm <= norm_goal)
      {
-         converged = 1;
+         converged = true;


⚠️ readability-implicit-bool-conversion ⚠️
implicit conversion bool -> int

Suggested change

converged = true;

converged = 1;

github-actions · 2022-11-14T20:34:13Z

src/common/inexact_newton.cpp

         final_iter = it;
         break;
      }

      if (it >= max_iter)
      {
-         converged = 0;
+         converged = false;


⚠️ readability-implicit-bool-conversion ⚠️
implicit conversion bool -> int

Suggested change

converged = false;

converged = 0;

github-actions · 2022-11-14T20:34:13Z

src/common/inexact_newton.cpp

@@ -152,7 +152,7 @@

      if (c_scale == 0.0)
      {
-         converged = 0;
+         converged = true;


⚠️ readability-implicit-bool-conversion ⚠️
implicit conversion bool -> int

Suggested change

converged = true;

converged = 1;

github-actions · 2022-11-14T20:34:16Z

test/unit/test_l2_transfer_operator.cpp

+   op.vectorJacobianProduct(out_bar, "mesh_coords", mesh_coords_bar);
+
+   auto dout_dmesh_v_local = mesh_pert * mesh_coords_bar;
+   double dout_dmesh_v;


⚠️ cppcoreguidelines-init-variables ⚠️
variable dout_dmesh_v is not initialized

Suggested change

double dout_dmesh_v;

double dout_dmesh_v = NAN;

github-actions · 2022-11-14T20:34:16Z

test/regression/test_magnetostatic_residual.cpp

+         mfem::H1_FECollection fec(p, dim);
+         mfem::ParFiniteElementSpace fes(&mesh, &fec);
+
+         std::map<std::string, mach::FiniteElementState> fields;


⚠️ cppcoreguidelines-init-variables ⚠️
variable fields is not initialized

Suggested change

std::map<std::string, mach::FiniteElementState> fields;

std::map<std::string, mach::FiniteElementState> fields = 0;

github-actions · 2022-11-14T20:34:17Z

test/regression/test_magnetostatic_residual.cpp

+#include <random>
+
+#include "catch.hpp"
+#include "mfem.hpp"


⚠️ cppcoreguidelines-init-variables ⚠️
variable pert is not initialized

Suggested change

#include "mfem.hpp"

#include "math.h"

#include "mfem.hpp"

github-actions · 2022-11-14T20:34:17Z

test/regression/test_magnetostatic_residual.cpp

+            {"current_density:box2", -current_density}
+         };
+
+         double pert = uniform_rand(gen);


⚠️ cppcoreguidelines-init-variables ⚠️
variable pert is not initialized

Suggested change

double pert = uniform_rand(gen);

double pert = NAN = uniform_rand(gen);

github-actions · 2022-11-14T20:34:17Z

test/regression/test_magnetostatic_residual.cpp

+         jacobianVectorProduct(res, pert_vec, "current_density:box1", res_dot);
+         double drdp_fwd = res_bar * res_dot;
+
+         double drdp_rev = vectorJacobianProduct(res, res_bar, "current_density:box1") * pert;


⚠️ cppcoreguidelines-init-variables ⚠️
variable drdp_rev is not initialized

Suggested change

double drdp_rev = vectorJacobianProduct(res, res_bar, "current_density:box1") * pert;

double drdp_rev = NAN = vectorJacobianProduct(res, res_bar, "current_density:box1") * pert;

github-actions

⚠️ clang-tidy found issue(s) with the introduced code (4/8)

github-actions · 2022-11-14T20:34:33Z

test/regression/test_magnetostatic_residual.cpp

+         mfem::H1_FECollection fec(p, dim);
+         mfem::ParFiniteElementSpace fes(&mesh, &fec);
+
+         std::map<std::string, mach::FiniteElementState> fields;


⚠️ cppcoreguidelines-init-variables ⚠️
variable fields is not initialized

Suggested change

std::map<std::string, mach::FiniteElementState> fields;

std::map<std::string, mach::FiniteElementState> fields = 0;

github-actions · 2022-11-14T20:34:33Z