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 "NestedKKSMultiSplitCHCRes.h"
      11             : 
      12             : registerMooseObject("PhaseFieldApp", NestedKKSMultiSplitCHCRes);
      13             : 
      14             : InputParameters
      15          30 : NestedKKSMultiSplitCHCRes::validParams()
      16             : {
      17          30 :   InputParameters params = JvarMapKernelInterface<Kernel>::validParams();
      18          30 :   params.addClassDescription(
      19             :       "KKS model kernel for the split Bulk Cahn-Hilliard term. This kernel operates on the "
      20             :       "physical concentration 'c' as the non-linear variable.");
      21          60 :   params.addCoupledVar("all_etas", "Phase parameters for all phases.");
      22          60 :   params.addRequiredCoupledVar("global_cs", "The interpolated concentrations c, b, etc.");
      23          60 :   params.addCoupledVar("w", "Chemical potential non-linear helper variable for the split solve.");
      24          60 :   params.addParam<std::vector<MaterialPropertyName>>(
      25             :       "c1_names",
      26             :       "Phase concentrations in the frist phase of all_etas. The order must match global_cs, for "
      27             :       "example, c1, b1, etc.");
      28          60 :   params.addParam<MaterialPropertyName>("F1_name", "Free energy of the first phase in all_etas.");
      29          30 :   return params;
      30           0 : }
      31             : 
      32          16 : NestedKKSMultiSplitCHCRes::NestedKKSMultiSplitCHCRes(const InputParameters & parameters)
      33             :   : DerivativeMaterialInterface<JvarMapKernelInterface<Kernel>>(parameters),
      34          32 :     _eta_names(coupledNames("all_etas")),
      35          16 :     _num_j(_eta_names.size()),
      36          16 :     _eta_map(getParameterJvarMap("all_etas")),
      37          16 :     _c_names(coupledNames("global_cs")),
      38          16 :     _num_c(coupledComponents("global_cs")),
      39          16 :     _c_map(getParameterJvarMap("global_cs")),
      40          16 :     _o(-1),
      41          16 :     _w_var(coupled("w")),
      42          16 :     _w(coupledValue("w")),
      43          32 :     _c1_names(getParam<std::vector<MaterialPropertyName>>("c1_names")),
      44          16 :     _F1_name(getParam<MaterialPropertyName>("F1_name")),
      45          16 :     _dF1dc1(_num_c),
      46          16 :     _d2F1dc1db1(_num_c),
      47          16 :     _dc1db(_num_c),
      48          16 :     _dc1detaj(_num_c),
      49          32 :     _d2F1dc1darg(_n_args)
      50             : 
      51             : {
      52          32 :   for (const auto i : make_range(_num_c))
      53             :   {
      54             :     // Set _o to the position of the nonlinear variable in the list of global_cs
      55          16 :     if (coupled("global_cs", i) == _var.number())
      56          16 :       _o = i;
      57             :   }
      58             : 
      59             :   // _dcideta and _dcidb are computed in KKSPhaseConcentrationDerivatives
      60          32 :   for (const auto m : make_range(_num_c))
      61             :   {
      62          16 :     _dc1detaj[m].resize(_num_j);
      63          64 :     for (const auto n : make_range(_num_j))
      64          48 :       _dc1detaj[m][n] = &getMaterialPropertyDerivative<Real>(_c1_names[m], _eta_names[n]);
      65             : 
      66          16 :     _dc1db[m].resize(_num_c);
      67          32 :     for (const auto n : make_range(_num_c))
      68          16 :       _dc1db[m][n] = &getMaterialPropertyDerivative<Real>(_c1_names[m], _c_names[n]);
      69             :   }
      70             : 
      71             :   // _dF1dc1 and _d2F1dc1db1 are computed in KKSPhaseConcentrationMaterial
      72          32 :   for (const auto m : make_range(_num_c))
      73             :   {
      74          16 :     _dF1dc1[m] = &getMaterialPropertyDerivative<Real>(_F1_name, _c1_names[m]);
      75          16 :     _d2F1dc1db1[m] = &getMaterialPropertyDerivative<Real>(_F1_name, _c1_names[_o], _c1_names[m]);
      76             :   }
      77             : 
      78             :   // _d2F1dc1darg is computed in KKSPhaseConcentrationMaterial
      79         160 :   for (const auto m : make_range(_n_args))
      80         144 :     _d2F1dc1darg[m] = &getMaterialPropertyDerivative<Real>(_F1_name, _c1_names[_o], m);
      81          16 : }
      82             : 
      83             : Real
      84     3795200 : NestedKKSMultiSplitCHCRes::computeQpResidual()
      85             : {
      86     3795200 :   return ((*_dF1dc1[_o])[_qp] - _w[_qp]) * _test[_i][_qp];
      87             : }
      88             : 
      89             : Real
      90     2252800 : NestedKKSMultiSplitCHCRes::computeQpJacobian()
      91             : {
      92             :   Real sum = 0.0;
      93             : 
      94     4505600 :   for (const auto m : make_range(_num_c))
      95     2252800 :     sum += (*_d2F1dc1db1[m])[_qp] * (*_dc1db[m][_o])[_qp];
      96             : 
      97     2252800 :   return sum * _phi[_j][_qp] * _test[_i][_qp];
      98             : }
      99             : 
     100             : Real
     101     9011200 : NestedKKSMultiSplitCHCRes::computeQpOffDiagJacobian(unsigned int jvar)
     102             : {
     103             :   Real sum = 0.0;
     104             : 
     105             :   // treat w variable explicitly
     106     9011200 :   if (jvar == _w_var)
     107     2252800 :     return -_phi[_j][_qp] * _test[_i][_qp];
     108             : 
     109             :   // if b is the coupled variable
     110     6758400 :   auto compvar = mapJvarToCvar(jvar, _c_map);
     111     6758400 :   if (compvar >= 0)
     112             :   {
     113           0 :     for (const auto m : make_range(_num_c))
     114           0 :       sum += (*_d2F1dc1db1[m])[_qp] * (*_dc1db[m][compvar])[_qp];
     115             : 
     116           0 :     return sum * _phi[_j][_qp] * _test[_i][_qp];
     117             :   }
     118             : 
     119             :   // if order parameters are the coupled variables
     120     6758400 :   auto etavar = mapJvarToCvar(jvar, _eta_map);
     121     6758400 :   if (etavar >= 0)
     122             :   {
     123    13516800 :     for (const auto m : make_range(_num_c))
     124     6758400 :       sum += (*_d2F1dc1db1[m])[_qp] * (*_dc1detaj[m][etavar])[_qp];
     125             : 
     126     6758400 :     return sum * _phi[_j][_qp] * _test[_i][_qp];
     127             :   }
     128             : 
     129             :   // for all other vars get the coupled variable jvar is referring to
     130             :   const unsigned int cvar = mapJvarToCvar(jvar);
     131             : 
     132           0 :   return (*_d2F1dc1darg[cvar])[_qp] * _phi[_j][_qp] * _test[_i][_qp];
     133             : }