fix compile errors and test fails

Author: dschwen

examples/degeus_mechanics/mech.i

@@ -11,6 +11,11 @@
 
 [TensorComputes]
   [Initialize]
+    [Finit]
+      type = RankTwoIdentity
+      buffer = F
+    []
+
     [phase]
       type = PhaseMechanicsTest
       buffer = phase
@@ -76,6 +81,16 @@
   []
 []
 
+[TensorSolver]
+  # no variables are integrated by this solver (FFTMechanics performs a steady state mechanics solve)
+  type = ForwardEulerSolver
+  root_compute = root
+  # deformation tensor is just forwarded Fnew -> F
+  forward_buffer = F
+  forward_buffer_new = Fnew
+  substeps = 10
+[]
+
 [TensorOutputs]
   [deformation_tensor]
     type = XDMFTensorOutput

include/base/SwiftApp.h

@@ -24,7 +24,7 @@ class SwiftApp : public MooseApp
 public:
   static InputParameters validParams();
 
-  SwiftApp(InputParameters parameters);
+  SwiftApp(const InputParameters & parameters);
   virtual ~SwiftApp();
 
   static void registerApps();

include/problems/TensorProblem.h

@@ -41,6 +41,7 @@ namespace Swift
 {
 struct ConstantBase
 {
+  virtual ~ConstantBase() = default;
   virtual std::string getType() = 0;
 };
 
@@ -185,6 +186,7 @@ class TensorProblem : public FEProblem, public DomainInterface
   /// substepping timestep
   Real _sub_dt;
   Real _sub_time;
+  unsigned int _substep;
 
   /// simulation time for the currently running output thread
   Real _output_time;

include/tensor_computes/FFTMechanics.h

@@ -69,6 +69,7 @@ class FFTMechanics : public TensorOperator<>
   const torch::Tensor * const _applied_macroscopic_strain;
 
   const bool _verbose;
+  const bool _accept_nonconverged;
 
   using TensorOperatorBase::_dim;
   using TensorOperatorBase::_domain;

include/tensor_computes/HyperElasticIsotropic.h

@@ -9,7 +9,6 @@
 #pragma once
 
 #include "TensorOperator.h"
-#include "neml2/tensors/Vec.h"
 
 /**
  * Hyperelastic isotropic material model

include/tensor_computes/NEML2TensorCompute.h

@@ -30,7 +30,7 @@ class NEML2TensorCompute : public TensorOperatorBase
 
 protected:
 #ifdef NEML2_ENABLED
-  neml2::Model & _model;
+  std::shared_ptr<neml2::Model> _model;
 
   std::vector<std::pair<const torch::Tensor *, neml2::LabeledAxisAccessor>> _input_mapping;
   std::vector<std::pair<neml2::LabeledAxisAccessor, torch::Tensor *>> _output_mapping;

moose

@@ -47,7 +47,7 @@ SwiftApp::validParams()
   return params;
 }
 
-SwiftApp::SwiftApp(InputParameters parameters) : MooseApp(parameters)
+SwiftApp::SwiftApp(const InputParameters & parameters) : MooseApp(parameters)
 {
   SwiftApp::registerAll(_factory, _action_factory, _syntax);
   MooseTensor::swift_global_settings._torch_device =

src/base/SwiftApp.C

@@ -44,6 +44,7 @@ FFTMechanics::validParams()
                                          "Applied macroscopic strain");
   params.addParam<TensorInputBufferName>("F", "F", "Deformation gradient tensor.");
   params.addParam<bool>("verbose", false, "Print non-linear residuals.");
+  params.addParam<bool>("accept_nonconverged", false, "Accept a non-converged solution");
   return params;
 }
 
@@ -71,7 +72,8 @@ FFTMechanics::FFTMechanics(const InputParameters & parameters)
     _applied_macroscopic_strain(isParamValid("applied_macroscopic_strain")
                                     ? &getInputBuffer("applied_macroscopic_strain")
                                     : nullptr),
-    _verbose(getParam<bool>("verbose"))
+    _verbose(getParam<bool>("verbose")),
+    _accept_nonconverged(getParam<bool>("accept_nonconverged"))
 {
   // Build projection tensor once
   const auto & q = _domain.getKGrid();
@@ -150,7 +152,7 @@ FFTMechanics::computeBuffer()
 
     // print nonlinear residual to the screen
     if (_verbose)
-      _console << "|R|=" << anorm << "\t|R/R0|=" << rnorm << '\n';
+      _console << "|R|=" << anorm << "\t|R/R0|=" << rnorm << std::endl;
 
     // check convergence
     if ((rnorm < _nl_rel_tol || anorm < _nl_abs_tol) && iiter > 0)
@@ -159,8 +161,18 @@ FFTMechanics::computeBuffer()
     iiter++;
 
     if (iiter > _nl_max_its)
-      paramError("nl_max_its",
-                 "Exceeded the maximum number of nonlinear iterations without converging.");
+    {
+      if (_accept_nonconverged)
+        paramWarning(
+            "nl_max_its", "Accepting non-converged solution |R|=", anorm, "\t|R/R0|=", rnorm);
+      else
+        paramError(
+            "nl_max_its",
+            "Exceeded the maximum number of nonlinear iterations without converging with |R|=",
+            anorm,
+            "\t|R/R0|=",
+            rnorm);
+    }
   }
 }
 

src/tensor_computes/FFTMechanics.C

@@ -33,7 +33,7 @@ HyperElasticIsotropic::HyperElasticIsotropic(const InputParameters & parameters)
     _tI4s((_tI4 + _tI4rt) / 2.0),
     _tII(MooseTensor::dyad22(_tI, _tI)),
     _tF(getInputBuffer("F")),
-    _pFstar(isParamValid("Fstar") ? &getInputBuffer("F") : nullptr),
+    _pFstar(isParamValid("Fstar") ? &getInputBuffer("Fstar") : nullptr),
     _tmu(getInputBuffer("mu")),
     _tK(getInputBuffer("K")),
     _tK4(getOutputBuffer("tangent_operator"))
@@ -60,8 +60,10 @@ HyperElasticIsotropic::computeBuffer()
   // First Piola-Kirchhoff stress: P = F_e * S
   _u = dot22(F_e, S);
 
+  _tK4 = dot24(S, _tI4) + ddot44(ddot44(_tI4rt, dot42(dot24(F_e, C4), trans2(F_e))), _tI4rt);
+
   // Consistent tangent: K_ijkl = F^e_{im} C_{jlmn} F^e_{kn} + delta_{ik} S_{jl}
-  const auto K_geo = torch::einsum("...im,...jlmn,...kn->...ijkl", F_e, C4, F_e);
-  const auto K_stress = torch::einsum("ik,jl,...jl->...ijkl", _ti, _ti, S);
-  _tK4 = K_geo + K_stress;
+  // const auto K_geo = torch::einsum("...im,...jlmn,...kn->...ijkl", {F_e, C4, F_e});
+  // const auto K_stress = torch::einsum("ik,jl,...jl->...ijkl", {_ti, _ti, S});
+  // _tK4 = K_geo + K_stress;
 }