Add elastic inversion capability using batch materials

Makefile

@@ -3,7 +3,7 @@
 ###############################################################################
 #
 # Optional Environment variables
-# MOOSE_DIR        - Root directory of the MOOSE project 
+# MOOSE_DIR        - Root directory of the MOOSE project
 #
 ###############################################################################
 # Use the MOOSE submodule if it exists and MOOSE_DIR is not set
@@ -40,6 +40,7 @@ TENSOR_MECHANICS    := yes
 WATER_STEAM_EOS     := no
 XFEM                := yes
 POROUS_FLOW         := no
+OPTIMIZATION        := yes
 
 include $(MOOSE_DIR)/modules/modules.mk
 ###############################################################################

doc/config.yml

@@ -19,6 +19,8 @@ Content:
         root_dir: ${MOOSE_DIR}/modules/stochastic_tools/doc/content
     ray_tracing:
         root_dir: ${MOOSE_DIR}/modules/ray_tracing/doc/content
+    optimization:
+        root_dir: ${MOOSE_DIR}/modules/optimization/doc/content
     modules:
         root_dir: ${MOOSE_DIR}/modules/doc/content
         content:
@@ -83,6 +85,7 @@ Extensions:
             misc: !include ${MOOSE_DIR}/modules/misc/doc/sqa_misc.yml
             xfem: !include ${MOOSE_DIR}/modules/xfem/doc/sqa_xfem.yml
             ray_tracing: !include ${MOOSE_DIR}/modules/ray_tracing/doc/sqa_ray_tracing.yml
+            optimization: !include ${MOOSE_DIR}/modules/optimization/doc/sqa_optimization.yml
             blackbear: !include ${ROOT_DIR}/doc/sqa_blackbear.yml
         requirement-groups:
             dgkernels: DGKernel Objects

doc/content/source/userobjects/BatchStressGrad.md

@@ -0,0 +1,22 @@
+# BatchStressGrad
+
+## Description
+
+This `UserObject` computes double contraction of the elastic tensor derivative and the forward mechanical strain, i.e.,
+\begin{equation}
+\underbrace{\frac{\partial \boldsymbol{C}}{\partial p}}_{\text{elastic tensor derivative}} \underbrace{(\boldsymbol{L} u)}_{\text{forward strain}}
+\end{equation}
+as a batch material.
+Here, the $\boldsymbol{C}$ is the elastic tensor, $p$ is the interested parameter, $\boldsymbol{L}$ is the differential operator for elasticity problem, and $\boldsymbol{L} u$ is strain from the forward problem. This object requires the elastic tensor derivative material property (i.e., $\frac{\partial \boldsymbol{C}}{\partial p}$) as its input.
+
+This object is used together with [AdjointStrainStressGradInnerProduct](/AdjointStrainStressGradInnerProduct.md) in a elastic inversion problem.
+
+## Example Input File Syntax
+
+!listing test/tests/bimaterial_elastic_inversion_batch_mat/grad.i block=UserObjects/stress_grad_lambda
+
+!syntax parameters /UserObjects/BatchStressGrad
+
+!syntax inputs /UserObjects/BatchStressGrad
+
+!syntax children /UserObjects/BatchStressGrad

doc/content/source/vectorpostprocessors/AdjointStrainStressGradInnerProduct.md

@@ -0,0 +1,23 @@
+# AdjointStrainStressGradInnerProduct
+
+!syntax description /VectorPostprocessors/AdjointStrainStressGradInnerProduct
+
+## Description
+
+This `VectorPostprocessor` computes the gradient of objective function with respect to interested parameter for elastic problem. Specifically,
+\begin{equation}
+  \frac{\text{d}\boldsymbol{F}}{\text{d}p} = - \int \underbrace{(\boldsymbol{L}\lambda)^T}_{\text{adjoint strain}}  \underbrace{\frac{\partial \boldsymbol{C}}{\partial p} (\boldsymbol{L} u)}_{\text{stress gradient}} \text{d}\Omega,
+\end{equation}
+where $\boldsymbol{F}$ is the objective function, $\boldsymbol{L}$ is the differential operator for the elasticity problem, $\lambda$ is the adjoint displacement, $\boldsymbol{C}$ is the elastic tensor, $p$ is the interested parameter, and $u$ is displacement from the forward problem.
+
+Specifically, this object takes the `adjoint strain` and the `stress gradient` as inputs, computes the double contraction of the two inputs, and integrate over the entire domain. Specifically, the `adjoint strain` is the strain from the adjoint problem, the `stress gradient` should be a batch material (e.g., [BatchStressGrad](/BatchStressGrad.md)).
+
+## Example Input File Syntax
+
+!listing test/tests/bimaterial_elastic_inversion_batch_mat/grad.i block=VectorPostprocessors/grad_lambda
+
+!syntax parameters /VectorPostprocessors/AdjointStrainStressGradInnerProduct
+
+!syntax inputs /VectorPostprocessors/AdjointStrainStressGradInnerProduct
+
+!syntax children /VectorPostprocessors/AdjointStrainStressGradInnerProduct

include/userobjects/BatchStressGrad.h

@@ -0,0 +1,33 @@
+//* This file is part of the MOOSE framework
+//* https://www.mooseframework.org
+//*
+//* All rights reserved, see COPYRIGHT for full restrictions
+//* https://github.com/idaholab/moose/blob/master/COPYRIGHT
+//*
+//* Licensed under LGPL 2.1, please see LICENSE for details
+//* https://www.gnu.org/licenses/lgpl-2.1.html
+
+#pragma once
+
+#include "BatchMaterial.h"
+
+typedef BatchMaterial<
+    // tuple representation
+    BatchMaterialUtils::TupleStd,
+    // output data type
+    RankTwoTensor,
+    // gathered input data types:
+    BatchMaterialUtils::GatherMatProp<RankFourTensor>,
+    BatchMaterialUtils::GatherMatProp<RankTwoTensor>>
+
+    BatchStressGradParent;
+
+class BatchStressGrad : public BatchStressGradParent
+{
+public:
+  static InputParameters validParams();
+
+  BatchStressGrad(const InputParameters & params);
+
+  void batchCompute() override;
+};

include/vectorpostprocessors/AdjointStrainStressGradInnerProduct.h

@@ -0,0 +1,31 @@
+//* This file is part of the MOOSE framework
+//* https://www.mooseframework.org
+//*
+//* All rights reserved, see COPYRIGHT for full restrictions
+//* https://github.com/idaholab/moose/blob/master/COPYRIGHT
+//*
+//* Licensed under LGPL 2.1, please see LICENSE for details
+//* https://www.gnu.org/licenses/lgpl-2.1.html
+
+#pragma once
+
+#include "ElementOptimizationFunctionInnerProduct.h"
+#include "BatchStressGrad.h"
+
+class AdjointStrainStressGradInnerProduct : public ElementOptimizationFunctionInnerProduct
+{
+public:
+  static InputParameters validParams();
+
+  AdjointStrainStressGradInnerProduct(const InputParameters & parameters);
+
+protected:
+  virtual Real computeQpInnerProduct() override;
+  /// Base name of the material system
+  const std::string _base_name;
+  /// Holds adjoint strain at current quadrature points
+  const MaterialProperty<RankTwoTensor> & _adjoint_strain;
+  /// UO that holds gradient of stress wrt material parameter at current quadrature points
+  const BatchStressGrad & _stress_grad_uo;
+  const BatchStressGrad::OutputVector & _stress_grad_uo_output;
+};

src/userobjects/BatchStressGrad.C

@@ -0,0 +1,46 @@
+//* This file is part of the MOOSE framework
+//* https://www.mooseframework.org
+//*
+//* All rights reserved, see COPYRIGHT for full restrictions
+//* https://github.com/idaholab/moose/blob/master/COPYRIGHT
+//*
+//* Licensed under LGPL 2.1, please see LICENSE for details
+//* https://www.gnu.org/licenses/lgpl-2.1.html
+
+#include "BatchStressGrad.h"
+#include "libmesh/int_range.h"
+
+registerMooseObject("BlackBearApp", BatchStressGrad);
+
+InputParameters
+BatchStressGrad::validParams()
+{
+  auto params = BatchStressGradParent::validParams();
+  params.addRequiredParam<MaterialPropertyName>(
+      "elastic_tensor_derivative", "Name of the elastic tensor derivative material property.");
+  return params;
+}
+
+BatchStressGrad::BatchStressGrad(const InputParameters & params)
+  : BatchStressGradParent(params,
+                          // here we pass the derivative of elastic tensor wrt to the parameter
+                          "elastic_tensor_derivative",
+                          // here we pass in the forward strain
+                          "forward_mechanical_strain")
+{
+}
+
+void
+BatchStressGrad::batchCompute()
+{
+  for (const auto i : index_range(_input_data))
+  {
+    const auto & input = _input_data[i];
+    auto & output = _output_data[i];
+
+    const auto & elasticity_dev = std::get<0>(input);
+    const auto & strain = std::get<1>(input);
+
+    output = elasticity_dev * strain;
+  }
+}

src/vectorpostprocessors/AdjointStrainStressGradInnerProduct.C

@@ -0,0 +1,48 @@
+//* This file is part of the MOOSE framework
+//* https://www.mooseframework.org
+//*
+//* All rights reserved, see COPYRIGHT for full restrictions
+//* https://github.com/idaholab/moose/blob/master/COPYRIGHT
+//*
+//* Licensed under LGPL 2.1, please see LICENSE for details
+//* https://www.gnu.org/licenses/lgpl-2.1.html
+
+#include "AdjointStrainStressGradInnerProduct.h"
+
+registerMooseObject("BlackBearApp", AdjointStrainStressGradInnerProduct);
+
+InputParameters
+AdjointStrainStressGradInnerProduct::validParams()
+{
+  InputParameters params = ElementOptimizationFunctionInnerProduct::validParams();
+  params.addRequiredParam<UserObjectName>(
+      "stress_grad_name",
+      "Name of the stress gradient user object with respect to the material parameter");
+
+  params.addRequiredParam<MaterialPropertyName>(
+      "adjoint_strain_name", "Name of the strain property in the adjoint problem");
+
+  params.addClassDescription("Compute the gradient for elastic material inversion");
+  return params;
+}
+AdjointStrainStressGradInnerProduct::AdjointStrainStressGradInnerProduct(
+    const InputParameters & parameters)
+  : ElementOptimizationFunctionInnerProduct(parameters),
+    _base_name(isParamValid("base_name") ? getParam<std::string>("base_name") + "_" : ""),
+    _adjoint_strain(getMaterialPropertyByName<RankTwoTensor>(
+        getParam<MaterialPropertyName>("adjoint_strain_name"))),
+    _stress_grad_uo(getUserObject<BatchStressGrad>("stress_grad_name")),
+    _stress_grad_uo_output(_stress_grad_uo.getOutputData())
+{
+}
+
+Real
+AdjointStrainStressGradInnerProduct::computeQpInnerProduct()
+{
+  if (!_stress_grad_uo.outputReady())
+    mooseError("Stress gradient batch material property is not ready to output.");
+
+  const auto index = _stress_grad_uo.getIndex(_current_elem->id());
+
+  return -_adjoint_strain[_qp].doubleContraction(_stress_grad_uo_output[index + _qp]);
+}

test/tests/bimaterial_elastic_inversion_batch_mat/gold/main_out_OptimizationReporter_0001.csv

@@ -0,0 +1,10 @@
+grad_lambda,grad_mu,lambda,measurement_time,measurement_values,measurement_xcoord,measurement_ycoord,measurement_zcoord,misfit_values,mu,simulation_values
+5.1698615078396e-11,5.2394323056988e-11,4.9999752463318,0,4.46462,-1,-1,0,9.2479240265675e-08,1.0000001055133,4.4646200924792
+5.097187284528e-12,3.5417284419815e-11,4.0000666870441,0,6.447155,-1,0,0,5.2641227377137e-08,2.0000064790597,6.4471550526412
+-2.7212788369455e-12,8.6258766757591e-12,2.9997143457064,0,8.434803,-1,1,0,-1.0044113984975e-07,2.9999940961679,8.4348028995589
+0,0,0,0,4.176264,0,-1,0,-1.14358395642e-07,0,4.1762638856416
+0,0,0,0,6.172984,0,0,0,-1.7726937251439e-07,0,6.1729838227306
+0,0,0,0,8.200859,0,1,0,1.9391855410333e-07,0,8.2008591939186
+0,0,0,0,4.049477,1,-1,0,2.2853767767117e-08,0,4.0494770228538
+0,0,0,0,6.052499,1,0,0,1.1289934853664e-07,0,6.0524991128993
+0,0,0,0,8.079385,1,1,0,-8.6613697547477e-08,0,8.0793849133863