LCOV - code coverage report
Current view: top level - src/utils - ElasticityTensorTools.C (source / functions) Hit Total Coverage
Test: idaholab/moose solid_mechanics: f45d79 Lines: 66 68 97.1 %
Date: 2025-07-25 05:00:39 Functions: 6 6 100.0 %
Legend: Lines: hit not hit

          Line data    Source code
       1             : //* This file is part of the MOOSE framework
       2             : //* https://mooseframework.inl.gov
       3             : //*
       4             : //* All rights reserved, see COPYRIGHT for full restrictions
       5             : //* https://github.com/idaholab/moose/blob/master/COPYRIGHT
       6             : //*
       7             : //* Licensed under LGPL 2.1, please see LICENSE for details
       8             : //* https://www.gnu.org/licenses/lgpl-2.1.html
       9             : 
      10             : #include "MooseTypes.h"
      11             : #include "PermutationTensor.h"
      12             : #include "RankFourTensor.h"
      13             : #include "ElasticityTensorTools.h"
      14             : 
      15             : #include "libmesh/vector_value.h"
      16             : 
      17             : namespace ElasticityTensorTools
      18             : {
      19             : 
      20             : Real
      21  2302056944 : elasticJacobian(const RankFourTensor & r4t,
      22             :                 unsigned int i,
      23             :                 unsigned int k,
      24             :                 const RealGradient & grad_test,
      25             :                 const RealGradient & grad_phi)
      26             : {
      27             :   // d(stress_ij*d(test)/dx_j)/du_k = d(C_ijmn*du_m/dx_n dtest/dx_j)/du_k (which is nonzero for m ==
      28             :   // k)
      29             : 
      30  2302056944 :   const Real gt0 = grad_test(0);
      31  2302056944 :   const Real gt1 = grad_test(1);
      32  2302056944 :   const Real gt2 = grad_test(2);
      33  2302056944 :   const Real gp0 = grad_phi(0);
      34  2302056944 :   const Real gp1 = grad_phi(1);
      35  2302056944 :   const Real gp2 = grad_phi(2);
      36             : 
      37             :   // clang-format off
      38             :   // This is the algorithm that is unrolled below:
      39             :   //
      40             :   //    Real sum = 0.0;
      41             :   //    for (const auto j: make_range(Moose::dim))
      42             :   //      for (const auto l: make_range(Moose::dim))
      43             :   //        sum += r4t(i, j, k, l) * grad_phi(l) * grad_test(j);
      44             :   //    return sum;
      45             : 
      46             :   return
      47             :      (
      48  2302056944 :          r4t(i,0,k,0) * gp0
      49  2302056944 :        + r4t(i,0,k,1) * gp1
      50  2302056944 :        + r4t(i,0,k,2) * gp2
      51  2302056944 :      ) * gt0
      52  2302056944 :      +
      53             :      (
      54  2302056944 :          r4t(i,1,k,0) * gp0
      55  2302056944 :        + r4t(i,1,k,1) * gp1
      56  2302056944 :        + r4t(i,1,k,2) * gp2
      57  2302056944 :      ) * gt1
      58             :      +
      59             :      (
      60  2302056944 :          r4t(i,2,k,0) * gp0
      61  2302056944 :        + r4t(i,2,k,1) * gp1
      62  2302056944 :        + r4t(i,2,k,2) * gp2
      63  2302056944 :      ) * gt2;
      64             :   // clang-format on
      65             : }
      66             : 
      67             : Real
      68    13320192 : elasticJacobianWC(const RankFourTensor & r4t,
      69             :                   unsigned int i,
      70             :                   unsigned int k,
      71             :                   const RealGradient & grad_test,
      72             :                   Real phi)
      73             : {
      74             :   // d(stress_ij*d(test)/dx_j)/dw_k = d(C_ijmn*eps_mnp*w_p*dtest/dx_j)/dw_k (only nonzero for p ==
      75             :   // k)
      76             :   Real sum = 0.0;
      77    53280768 :   for (const auto j : make_range(Moose::dim))
      78   159842304 :     for (const auto m : make_range(Moose::dim))
      79   479526912 :       for (const auto n : make_range(Moose::dim))
      80   359645184 :         sum += r4t(i, j, m, n) * PermutationTensor::eps(m, n, k) * grad_test(j);
      81    13320192 :   return sum * phi;
      82             : }
      83             : 
      84             : Real
      85    13320192 : momentJacobian(const RankFourTensor & r4t,
      86             :                unsigned int i,
      87             :                unsigned int k,
      88             :                Real test,
      89             :                const RealGradient & grad_phi)
      90             : {
      91             :   // Jacobian entry: d(eps_ijm*stress_jm*test)/du_k = d(eps_ijm*C_jmln*du_l/dx_n*test)/du_k (only
      92             :   // nonzero for l == k)
      93             :   Real sum = 0.0;
      94    53280768 :   for (const auto j : make_range(Moose::dim))
      95   159842304 :     for (const auto m : make_range(Moose::dim))
      96   479526912 :       for (const auto n : make_range(Moose::dim))
      97   359645184 :         sum += PermutationTensor::eps(i, j, m) * r4t(j, m, k, n) * grad_phi(n);
      98    13320192 :   return test * sum;
      99             : }
     100             : 
     101             : Real
     102     9453568 : momentJacobianWC(const RankFourTensor & r4t, unsigned int i, unsigned int k, Real test, Real phi)
     103             : {
     104             :   // Jacobian entry: d(eps_ijm*stress_jm*test)/dw_k = d(eps_ijm*C_jmln*eps_lnp*w_p*test)/dw_k (only
     105             :   // nonzero for p ==k)
     106             :   Real sum = 0.0;
     107    37814272 :   for (const auto j : make_range(Moose::dim))
     108   113442816 :     for (const auto l : make_range(Moose::dim))
     109   340328448 :       for (const auto m : make_range(Moose::dim))
     110  1020985344 :         for (const auto n : make_range(Moose::dim))
     111   765739008 :           sum +=
     112   765739008 :               PermutationTensor::eps(i, j, m) * r4t(j, m, l, n) * PermutationTensor::eps(l, n, k);
     113             : 
     114     9453568 :   return test * phi * sum;
     115             : }
     116             : 
     117             : void
     118     7817472 : toVoigtNotationIndexConversion(int k, int & a, int & b)
     119             : {
     120     7817472 :   if (k < 3 && k >= 0)
     121     3908736 :     a = b = k;
     122     3908736 :   else if (k == 3)
     123             :   {
     124     1302912 :     a = 1;
     125     1302912 :     b = 2;
     126             :   }
     127     2605824 :   else if (k == 4)
     128             :   {
     129     1302912 :     a = 0;
     130     1302912 :     b = 2;
     131             :   }
     132     1302912 :   else if (k == 5)
     133             :   {
     134     1302912 :     a = 0;
     135     1302912 :     b = 1;
     136             :   }
     137             :   else
     138           0 :     mooseError("\nIndex out of bound while converting from tensor to voigt notation in "
     139             :                "toVoigtNotationIndexConversion");
     140     7817472 : }
     141             : 
     142             : // MOOSE uses the stress tensor with components ordered as (11, 22, 33, 12, 23, 13) in voigt form.
     143             : // This is in conflict with the voigt notation used in literature which has components ordered as
     144             : // (11, 22, 33, 23, 13, 12). Whenever an operation involving the voigt forms of stress and
     145             : // elasticity tensor has to be performed the voigt form of the elasticity tensor should be built
     146             : // taking into account this difference in component-ordering for stress. The
     147             : // toMooseVoigtNotationIndexConversion function facilitates that.
     148             : 
     149             : void
     150   529128864 : toMooseVoigtNotationIndexConversion(int k, int & a, int & b)
     151             : {
     152   529128864 :   if (k < 3 && k >= 0)
     153   264564432 :     a = b = k;
     154   264564432 :   else if (k == 3)
     155             :   {
     156    88188144 :     a = 0;
     157    88188144 :     b = 1;
     158             :   }
     159   176376288 :   else if (k == 4)
     160             :   {
     161    88188144 :     a = 1;
     162    88188144 :     b = 2;
     163             :   }
     164    88188144 :   else if (k == 5)
     165             :   {
     166    88188144 :     a = 0;
     167    88188144 :     b = 2;
     168             :   }
     169             :   else
     170           0 :     mooseError("\nIndex out of bound while converting from tensor to MOOSE voigt notation in "
     171             :                "toMooseVoigtNotationIndexConversion");
     172   529128864 : }
     173             : }

Generated by: LCOV version 1.14