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 "PenaltyPeriodicSegmentalConstraint.h"
      11             : 
      12             : registerMooseObject("MooseApp", PenaltyPeriodicSegmentalConstraint);
      13             : 
      14             : InputParameters
      15       14453 : PenaltyPeriodicSegmentalConstraint::validParams()
      16             : {
      17       14453 :   InputParameters params = MortarScalarBase::validParams();
      18       14453 :   params.addClassDescription(
      19             :       "PenaltyPeriodicSegmentalConstraint enforces macro-micro periodic conditions between "
      20             :       "secondary and primary sides of a mortar interface using a penalty approach "
      21             :       "(no Lagrange multipliers needed). Must be used alongside PenaltyEqualValueConstraint.");
      22       14453 :   params.renameCoupledVar("scalar_variable", "epsilon", "Primary coupled scalar variable");
      23       14453 :   params.addRequiredCoupledVar("sigma", "Controlled scalar averaging variable");
      24       43359 :   params.addParam<Real>(
      25             :       "penalty_value",
      26       28906 :       1.0,
      27             :       "Penalty value used to impose a generalized force capturing the mortar constraint equation");
      28             : 
      29       14453 :   return params;
      30           0 : }
      31             : 
      32          94 : PenaltyPeriodicSegmentalConstraint::PenaltyPeriodicSegmentalConstraint(
      33          94 :     const InputParameters & parameters)
      34             :   : DerivativeMaterialInterface<MortarScalarBase>(parameters),
      35          94 :     _temp_jump_global(),
      36          94 :     _tau_s(),
      37          94 :     _kappa_aux_ptr(getScalarVar("sigma", 0)),
      38          94 :     _ka_order(_kappa_aux_ptr->order()),
      39          94 :     _kappa_aux(coupledScalarValue("sigma")),
      40         188 :     _pen_scale(getParam<Real>("penalty_value"))
      41             : {
      42          94 :   if (_kappa_aux_ptr->kind() != Moose::VarKindType::VAR_AUXILIARY)
      43           4 :     paramError("sigma",
      44             :                "Must assign auxiliary scalar variable to sigma, rather than nonlinear variable");
      45          90 : }
      46             : 
      47             : // Compute the stability parameters to use for all quadrature points
      48             : void
      49        3904 : PenaltyPeriodicSegmentalConstraint::precalculateResidual()
      50             : {
      51        3904 :   precalculateStability();
      52        3904 : }
      53             : void
      54           0 : PenaltyPeriodicSegmentalConstraint::precalculateJacobian()
      55             : {
      56           0 :   precalculateStability();
      57           0 : }
      58             : 
      59             : // Compute the temperature jump for current quadrature point
      60             : void
      61        7904 : PenaltyPeriodicSegmentalConstraint::initScalarQpResidual()
      62             : {
      63        7904 :   precalculateMaterial();
      64        7904 : }
      65             : 
      66             : Real
      67       97792 : PenaltyPeriodicSegmentalConstraint::computeQpResidual(const Moose::MortarType mortar_type)
      68             : {
      69             :   // Compute penalty parameter times x-jump times average heat flux
      70       97792 :   RealVectorValue dx(_phys_points_primary[_qp] - _phys_points_secondary[_qp]);
      71       97792 :   RealVectorValue kappa_vec(_kappa[0], 0, 0);
      72       97792 :   if (_k_order == 2)
      73       28672 :     kappa_vec(1) = _kappa[1];
      74       69120 :   else if (_k_order == 3)
      75             :   {
      76       69120 :     kappa_vec(1) = _kappa[1];
      77       69120 :     kappa_vec(2) = _kappa[2];
      78             :   }
      79       97792 :   Real r = _tau_s * (kappa_vec * dx);
      80             : 
      81       97792 :   switch (mortar_type)
      82             :   {
      83       48896 :     case Moose::MortarType::Secondary:
      84       48896 :       r *= _test_secondary[_i][_qp];
      85       48896 :       break;
      86       48896 :     case Moose::MortarType::Primary:
      87       48896 :       r *= -_test_primary[_i][_qp];
      88       48896 :       break;
      89           0 :     case Moose::MortarType::Lower:
      90           0 :       return 0;
      91           0 :     default:
      92           0 :       return 0;
      93             :   }
      94       97792 :   return r;
      95             : }
      96             : 
      97             : Real
      98       20128 : PenaltyPeriodicSegmentalConstraint::computeScalarQpResidual()
      99             : {
     100             :   // Stability/penalty term for residual of scalar variable
     101       20128 :   Real r = _tau_s * _temp_jump_global;
     102       20128 :   RealVectorValue dx(_phys_points_primary[_qp] - _phys_points_secondary[_qp]);
     103             : 
     104       20128 :   r *= -dx(_h);
     105             : 
     106       20128 :   RealVectorValue kappa_vec(_kappa[0], 0, 0);
     107       20128 :   RealVectorValue kappa_aux_vec(_kappa_aux[0], 0, 0);
     108       20128 :   if (_k_order == 2)
     109             :   {
     110        7168 :     kappa_vec(1) = _kappa[1];
     111        7168 :     kappa_aux_vec(1) = _kappa_aux[1];
     112             :   }
     113       12960 :   else if (_k_order == 3)
     114             :   {
     115       12960 :     kappa_vec(1) = _kappa[1];
     116       12960 :     kappa_vec(2) = _kappa[2];
     117       12960 :     kappa_aux_vec(1) = _kappa_aux[1];
     118       12960 :     kappa_aux_vec(2) = _kappa_aux[2];
     119             :   }
     120             : 
     121       20128 :   r += dx(_h) * _tau_s * (kappa_vec * dx);
     122       20128 :   r -= dx(_h) * (kappa_aux_vec * _normals[_qp]);
     123             : 
     124       20128 :   return r;
     125             : }
     126             : 
     127             : Real
     128       26400 : PenaltyPeriodicSegmentalConstraint::computeScalarQpJacobian()
     129             : {
     130             :   // Stability/penalty term for Jacobian of scalar variable
     131       26400 :   RealVectorValue dx(_phys_points_primary[_qp] - _phys_points_secondary[_qp]);
     132             : 
     133       26400 :   Real jac = dx(_h) * _tau_s * dx(_l);
     134             : 
     135       26400 :   return jac;
     136             : }
     137             : 
     138             : Real
     139      136704 : PenaltyPeriodicSegmentalConstraint::computeQpOffDiagJacobianScalar(
     140             :     const Moose::MortarType mortar_type, const unsigned int svar_num)
     141             : {
     142      136704 :   if (svar_num != _kappa_var)
     143           0 :     return 0;
     144             : 
     145             :   // Stability/penalty term for Jacobian
     146      136704 :   RealVectorValue dx(_phys_points_primary[_qp] - _phys_points_secondary[_qp]);
     147      136704 :   Real jac = _tau_s;
     148             : 
     149      136704 :   switch (mortar_type)
     150             :   {
     151             : 
     152       68352 :     case Moose::MortarType::Secondary: // Residual_sign -1  ddeltaU_ddisp sign 1;
     153       68352 :       jac *= _test_secondary[_i][_qp] * dx(_h);
     154       68352 :       break;
     155       68352 :     case Moose::MortarType::Primary: // Residual_sign -1  ddeltaU_ddisp sign -1;
     156       68352 :       jac *= -_test_primary[_i][_qp] * dx(_h);
     157       68352 :       break;
     158             : 
     159           0 :     default:
     160           0 :       return 0;
     161             :   }
     162      136704 :   return jac;
     163             : }
     164             : 
     165             : Real
     166      165760 : PenaltyPeriodicSegmentalConstraint::computeScalarQpOffDiagJacobian(
     167             :     const Moose::MortarType mortar_type, const unsigned int jvar_num)
     168             : {
     169             :   // Test if jvar is the ID of the primary variables and not some other random variable
     170      165760 :   switch (mortar_type)
     171             :   {
     172       66304 :     case Moose::MortarType::Secondary:
     173       66304 :       if (jvar_num != _secondary_var.number())
     174           0 :         return 0;
     175       66304 :       break;
     176       66304 :     case Moose::MortarType::Primary:
     177       66304 :       if (jvar_num != _primary_var.number())
     178           0 :         return 0;
     179       66304 :       break;
     180       33152 :     default:
     181       33152 :       return 0;
     182             :   }
     183             : 
     184             :   // Stability/penalty term for Jacobian
     185      132608 :   RealVectorValue dx(_phys_points_primary[_qp] - _phys_points_secondary[_qp]);
     186      132608 :   Real jac = _tau_s;
     187             : 
     188      132608 :   switch (mortar_type)
     189             :   {
     190       66304 :     case Moose::MortarType::Secondary:
     191       66304 :       jac *= (*_phi)[_j][_qp] * dx(_h);
     192       66304 :       break;
     193       66304 :     case Moose::MortarType::Primary:
     194       66304 :       jac *= -(*_phi)[_j][_qp] * dx(_h);
     195       66304 :       break;
     196             : 
     197           0 :     default:
     198           0 :       return 0;
     199             :   }
     200      132608 :   return jac;
     201             : }
     202             : 
     203             : void
     204        3904 : PenaltyPeriodicSegmentalConstraint::precalculateStability()
     205             : {
     206             :   // Example showing how the penalty could be loaded from some function
     207        3904 :   _tau_s = _pen_scale;
     208        3904 : }
     209             : 
     210             : // Compute temperature jump and flux average/jump
     211             : void
     212        7904 : PenaltyPeriodicSegmentalConstraint::precalculateMaterial()
     213             : {
     214        7904 :   _temp_jump_global = (_u_primary[_qp] - _u_secondary[_qp]);
     215        7904 : }