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 "ConvectedMesh.h"
      11             : 
      12             : registerMooseObject("FsiApp", ConvectedMesh);
      13             : 
      14             : InputParameters
      15          82 : ConvectedMesh::validParams()
      16             : {
      17          82 :   InputParameters params = INSBase::validParams();
      18          82 :   params.addClassDescription(
      19             :       "Corrects the convective derivative for situations in which the fluid mesh is dynamic.");
      20         164 :   params.addRequiredCoupledVar("disp_x", "The x displacement");
      21         164 :   params.addCoupledVar("disp_y", "The y displacement");
      22         164 :   params.addCoupledVar("disp_z", "The z displacement");
      23         164 :   params.addParam<MaterialPropertyName>("rho_name", "rho", "The name of the density");
      24         164 :   params.addParam<bool>(
      25         164 :       "supg", false, "Whether to perform SUPG stabilization of the momentum residuals");
      26          82 :   return params;
      27           0 : }
      28             : 
      29          44 : ConvectedMesh::ConvectedMesh(const InputParameters & parameters)
      30             :   : INSBase(parameters),
      31          44 :     _disp_x_dot(coupledDot("disp_x")),
      32          44 :     _d_disp_x_dot(coupledDotDu("disp_x")),
      33          44 :     _disp_x_id(coupled("disp_x")),
      34          88 :     _disp_y_dot(isCoupled("disp_y") ? coupledDot("disp_y") : _zero),
      35          88 :     _d_disp_y_dot(isCoupled("disp_y") ? coupledDotDu("disp_y") : _zero),
      36          44 :     _disp_y_id(coupled("disp_y")),
      37          44 :     _disp_z_dot(isCoupled("disp_z") ? coupledDot("disp_z") : _zero),
      38          44 :     _d_disp_z_dot(isCoupled("disp_z") ? coupledDotDu("disp_z") : _zero),
      39          44 :     _disp_z_id(coupled("disp_z")),
      40          88 :     _rho(getMaterialProperty<Real>("rho_name")),
      41         132 :     _supg(getParam<bool>("supg"))
      42             : {
      43          44 :   if (_var.number() == _u_vel_var_number)
      44          22 :     _component = 0;
      45          22 :   else if (_var.number() == _v_vel_var_number)
      46          22 :     _component = 1;
      47           0 :   else if (_var.number() == _w_vel_var_number)
      48           0 :     _component = 2;
      49             :   else
      50           0 :     paramError("variable", "The variable must match one of the velocity variables.");
      51          44 : }
      52             : 
      53             : Real
      54     7846720 : ConvectedMesh::strongResidual()
      55             : {
      56     7846720 :   return -_rho[_qp] * RealVectorValue(_disp_x_dot[_qp], _disp_y_dot[_qp], _disp_z_dot[_qp]) *
      57     7846720 :          _grad_u[_qp];
      58             : }
      59             : 
      60             : Real
      61     3707200 : ConvectedMesh::computeQpResidual()
      62             : {
      63     3707200 :   auto test = _test[_i][_qp];
      64     3707200 :   const auto U = relativeVelocity();
      65     3707200 :   if (_supg)
      66     3707200 :     test += tau() * _grad_test[_i][_qp] * U;
      67     3707200 :   return test * strongResidual();
      68             : }
      69             : 
      70             : Real
      71     4139520 : ConvectedMesh::computePGVelocityJacobian(const unsigned short component)
      72             : {
      73     4139520 :   const auto U = relativeVelocity();
      74     4139520 :   return strongResidual() * ((dTauDUComp(component) * _grad_test[_i][_qp] * U) +
      75     4139520 :                              (tau() * _grad_test[_i][_qp](component) * _phi[_j][_qp]));
      76             : }
      77             : 
      78             : Real
      79     2069760 : ConvectedMesh::computeQpJacobian()
      80             : {
      81     2069760 :   auto test = _test[_i][_qp];
      82     2069760 :   const auto U = relativeVelocity();
      83     2069760 :   if (_supg)
      84     2069760 :     test += tau() * _grad_test[_i][_qp] * U;
      85     2069760 :   auto jac = test * -_rho[_qp] *
      86     2069760 :              RealVectorValue(_disp_x_dot[_qp], _disp_y_dot[_qp], _disp_z_dot[_qp]) *
      87     2069760 :              _grad_phi[_j][_qp];
      88     2069760 :   if (_supg)
      89     2069760 :     jac += computePGVelocityJacobian(_component);
      90             : 
      91     2069760 :   return jac;
      92             : }
      93             : 
      94             : Real
      95     8279040 : ConvectedMesh::computeQpOffDiagJacobian(unsigned int jvar)
      96             : {
      97             :   mooseAssert(jvar != _var.number(), "Making sure I understand how old hand-coded Jacobians work.");
      98             : 
      99     8279040 :   auto test = _test[_i][_qp];
     100     8279040 :   const auto U = relativeVelocity();
     101     8279040 :   if (_supg)
     102     8279040 :     test += tau() * _grad_test[_i][_qp] * U;
     103             : 
     104     8279040 :   if (jvar == _disp_x_id)
     105     2069760 :     return test * -_rho[_qp] * _phi[_j][_qp] * _d_disp_x_dot[_qp] * _grad_u[_qp](0);
     106     6209280 :   else if (jvar == _disp_y_id)
     107     2069760 :     return test * -_rho[_qp] * _phi[_j][_qp] * _d_disp_y_dot[_qp] * _grad_u[_qp](1);
     108     4139520 :   else if (jvar == _disp_z_id)
     109           0 :     return test * -_rho[_qp] * _phi[_j][_qp] * _d_disp_z_dot[_qp] * _grad_u[_qp](2);
     110     4139520 :   else if (jvar == _u_vel_var_number)
     111             :   {
     112     1034880 :     if (_supg)
     113     1034880 :       return computePGVelocityJacobian(0);
     114             :     else
     115             :       return 0;
     116             :   }
     117     3104640 :   else if (jvar == _v_vel_var_number)
     118             :   {
     119     1034880 :     if (_supg)
     120     1034880 :       return computePGVelocityJacobian(1);
     121             :     else
     122             :       return 0;
     123             :   }
     124     2069760 :   else if (jvar == _w_vel_var_number)
     125             :   {
     126           0 :     if (_supg)
     127           0 :       return computePGVelocityJacobian(2);
     128             :     else
     129             :       return 0;
     130             :   }
     131             :   else
     132             :     return 0.0;
     133             : }