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             : #pragma once
      11             : 
      12             : #include "Kernel.h"
      13             : #include "Material.h"
      14             : #include "JvarMapInterface.h"
      15             : #include "DerivativeKernelInterface.h"
      16             : 
      17             : /**
      18             :  * This is the Cahn-Hilliard equation base class that implements the interfacial or gradient energy
      19             :  * term of the equation.
      20             :  * See M.R. Tonks et al. / Computational Materials Science 51 (2012) 20-29 for more information.
      21             :  */
      22             : template <typename T>
      23             : class CHInterfaceBase : public DerivativeMaterialInterface<JvarMapKernelInterface<Kernel>>
      24             : {
      25             : public:
      26             :   CHInterfaceBase(const InputParameters & parameters);
      27             : 
      28             :   static InputParameters validParams();
      29             : 
      30             : protected:
      31             :   virtual Real computeQpResidual();
      32             :   virtual Real computeQpJacobian();
      33             :   virtual Real computeQpOffDiagJacobian(unsigned int jvar);
      34             : 
      35             :   const MaterialProperty<Real> & _kappa;
      36             : 
      37             :   ///@{
      38             :   /// Mobility material property value and concentration derivatives
      39             :   const MaterialProperty<T> & _M;
      40             :   const MaterialProperty<T> & _dMdc;
      41             :   const MaterialProperty<T> & _d2Mdc2;
      42             :   ///@}
      43             : 
      44             :   ///@{
      45             :   /// Variables for second order derivatives
      46             :   const VariableSecond & _second_u;
      47             :   const VariableTestSecond & _second_test;
      48             :   const VariablePhiSecond & _second_phi;
      49             :   ///@}
      50             : 
      51             :   ///Number of variables
      52             :   unsigned int _nvar;
      53             : 
      54             :   ///@{
      55             :   /// Mobility derivatives w.r.t. its dependent variables
      56             :   std::vector<const MaterialProperty<T> *> _dMdarg;
      57             :   std::vector<const MaterialProperty<T> *> _d2Mdcdarg;
      58             :   std::vector<std::vector<const MaterialProperty<T> *>> _d2Mdargdarg;
      59             :   ///@}
      60             : 
      61             :   /// Coupled variables used in mobility
      62             :   std::vector<const VariableGradient *> _coupled_grad_vars;
      63             : };
      64             : 
      65             : template <typename T>
      66         222 : CHInterfaceBase<T>::CHInterfaceBase(const InputParameters & parameters)
      67             :   : DerivativeMaterialInterface<JvarMapKernelInterface<Kernel>>(parameters),
      68         222 :     _kappa(getMaterialProperty<Real>("kappa_name")),
      69         444 :     _M(getMaterialProperty<T>("mob_name")),
      70         222 :     _dMdc(getMaterialPropertyDerivative<T>("mob_name", _var.name())),
      71         222 :     _d2Mdc2(getMaterialPropertyDerivative<T>("mob_name", _var.name(), _var.name())),
      72         222 :     _second_u(second()),
      73         222 :     _second_test(secondTest()),
      74         222 :     _second_phi(secondPhi()),
      75         222 :     _nvar(_coupled_moose_vars.size()),
      76         222 :     _dMdarg(_nvar),
      77         222 :     _d2Mdcdarg(_nvar),
      78         222 :     _d2Mdargdarg(_nvar),
      79         444 :     _coupled_grad_vars(_nvar)
      80             : {
      81             :   // Iterate over all coupled variables
      82         234 :   for (unsigned int i = 0; i < _nvar; ++i)
      83             :   {
      84             :     // Set material property values
      85          12 :     _dMdarg[i] = &getMaterialPropertyDerivative<T>("mob_name", _coupled_moose_vars[i]->name());
      86          12 :     _d2Mdcdarg[i] =
      87          24 :         &getMaterialPropertyDerivative<T>("mob_name", _var.name(), _coupled_moose_vars[i]->name());
      88          12 :     _d2Mdargdarg[i].resize(_nvar);
      89          24 :     for (unsigned int j = 0; j < _nvar; ++j)
      90          12 :       _d2Mdargdarg[i][j] = &getMaterialPropertyDerivative<T>(
      91          12 :           "mob_name", _coupled_moose_vars[i]->name(), _coupled_moose_vars[j]->name());
      92             : 
      93             :     // Set coupled variable gradients
      94          12 :     _coupled_grad_vars[i] =
      95          24 :         isCoupled("args") ? &coupledGradient("args", i) : &coupledGradient("coupled_variables", i);
      96             :   }
      97         222 : }
      98             : 
      99             : template <typename T>
     100             : InputParameters
     101         420 : CHInterfaceBase<T>::validParams()
     102             : {
     103         420 :   InputParameters params = Kernel::validParams();
     104         420 :   params.addClassDescription("Gradient energy Cahn-Hilliard base Kernel");
     105         840 :   params.addRequiredParam<MaterialPropertyName>("kappa_name", "The kappa used with the kernel");
     106         840 :   params.addRequiredParam<MaterialPropertyName>("mob_name", "The mobility used with the kernel");
     107         840 :   params.addCoupledVar("args", "Vector of variable arguments of the mobility");
     108         840 :   params.deprecateCoupledVar("args", "coupled_variables", "02/27/2024");
     109         420 :   return params;
     110           0 : }
     111             : 
     112             : template <typename T>
     113             : Real
     114   102485152 : CHInterfaceBase<T>::computeQpResidual()
     115             : {
     116   102485152 :   RealGradient grad_M = _dMdc[_qp] * _grad_u[_qp];
     117   105110688 :   for (unsigned int i = 0; i < _nvar; ++i)
     118     2625536 :     grad_M += (*_dMdarg[i])[_qp] * (*_coupled_grad_vars[i])[_qp];
     119             : 
     120   102485152 :   return _kappa[_qp] * _second_u[_qp].tr() *
     121   102485152 :          ((_M[_qp] * _second_test[_i][_qp]).tr() + grad_M * _grad_test[_i][_qp]);
     122             : }
     123             : 
     124             : template <typename T>
     125             : Real
     126  1233086400 : CHInterfaceBase<T>::computeQpJacobian()
     127             : {
     128             :   // Set the gradient and gradient derivative values
     129  1233086400 :   RealGradient grad_M = _dMdc[_qp] * _grad_u[_qp];
     130             : 
     131             :   RealGradient dgrad_Mdc =
     132  1233086400 :       _d2Mdc2[_qp] * _phi[_j][_qp] * _grad_u[_qp] + _dMdc[_qp] * _grad_phi[_j][_qp];
     133             : 
     134  1265878976 :   for (unsigned int i = 0; i < _nvar; ++i)
     135             :   {
     136    32792576 :     grad_M += (*_dMdarg[i])[_qp] * (*_coupled_grad_vars[i])[_qp];
     137    32792576 :     dgrad_Mdc += (*_d2Mdcdarg[i])[_qp] * _phi[_j][_qp] * (*_coupled_grad_vars[i])[_qp];
     138             :   }
     139             : 
     140             :   // Jacobian value using product rule
     141  1233086400 :   Real value = _kappa[_qp] * _second_phi[_j][_qp].tr() *
     142  1233086400 :                    ((_M[_qp] * _second_test[_i][_qp]).tr() + grad_M * _grad_test[_i][_qp]) +
     143  1233086400 :                _kappa[_qp] * _second_u[_qp].tr() *
     144  1233086400 :                    ((_dMdc[_qp] * _second_test[_i][_qp]).tr() * _phi[_j][_qp] +
     145             :                     dgrad_Mdc * _grad_test[_i][_qp]);
     146             : 
     147  1233086400 :   return value;
     148             : }
     149             : 
     150             : template <typename T>
     151             : Real
     152    22962176 : CHInterfaceBase<T>::computeQpOffDiagJacobian(unsigned int jvar)
     153             : {
     154             :   // get the coupled variable jvar is referring to
     155             :   const unsigned int cvar = mapJvarToCvar(jvar);
     156             : 
     157             :   // Set the gradient derivative
     158    22962176 :   RealGradient dgrad_Mdarg = (*_d2Mdcdarg[cvar])[_qp] * _phi[_j][_qp] * _grad_u[_qp] +
     159    22962176 :                              (*_dMdarg[cvar])[_qp] * _grad_phi[_j][_qp];
     160             : 
     161    45924352 :   for (unsigned int i = 0; i < _nvar; ++i)
     162    22962176 :     dgrad_Mdarg += (*_d2Mdargdarg[cvar][i])[_qp] * _phi[_j][_qp] * (*_coupled_grad_vars[cvar])[_qp];
     163             : 
     164             :   // Jacobian value using product rule
     165    22962176 :   Real value = _kappa[_qp] * _second_u[_qp].tr() *
     166    22962176 :                (((*_dMdarg[cvar])[_qp] * _second_test[_i][_qp]).tr() * _phi[_j][_qp] +
     167    22962176 :                 dgrad_Mdarg * _grad_test[_i][_qp]);
     168             : 
     169    22962176 :   return value;
     170             : }