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 "ADPenaltyVelocityContinuity.h"
      11             : #include "Assembly.h"
      12             : 
      13             : registerMooseObject("FsiApp", ADPenaltyVelocityContinuity);
      14             : 
      15             : InputParameters
      16          60 : ADPenaltyVelocityContinuity::validParams()
      17             : {
      18          60 :   InputParameters params = InterfaceKernelBase::validParams();
      19          60 :   params.addClassDescription(
      20             :       "Enforces continuity of flux and continuity of solution via penalty across an interface.");
      21         120 :   params.addRequiredParam<Real>(
      22             :       "penalty",
      23             :       "The penalty that penalizes jump between primary and neighbor secondary variables.");
      24         120 :   params.addRequiredCoupledVar("fluid_velocity", "The fluid velocity");
      25         120 :   params.renameCoupledVar("neighbor_var", "displacements", "All the displacement variables");
      26         120 :   params.addRequiredCoupledVar("solid_velocities", "The solid velocity components");
      27          60 :   return params;
      28           0 : }
      29             : 
      30          32 : ADPenaltyVelocityContinuity::ADPenaltyVelocityContinuity(const InputParameters & parameters)
      31             :   : InterfaceKernelBase(parameters),
      32          32 :     _penalty(getParam<Real>("penalty")),
      33          32 :     _velocity_var(getVectorVar("fluid_velocity", 0)),
      34          32 :     _velocity(adCoupledVectorValue("fluid_velocity")),
      35          32 :     _ad_JxW(_assembly.adJxWFace()),
      36          32 :     _ad_coord(_assembly.adCoordTransformation())
      37             : {
      38          64 :   if (isCoupledConstant("fluid_velocity"))
      39           0 :     paramError("fluid_velocity", "The fluid velocity must be an actual variable");
      40             : 
      41          64 :   _solid_velocities.resize(coupledComponents("solid_velocities"));
      42          96 :   for (const auto i : index_range(_solid_velocities))
      43          64 :     _solid_velocities[i] = &adCoupledNeighborValue("solid_velocities", i);
      44             : 
      45          64 :   _displacements.resize(coupledComponents("displacements"));
      46          96 :   for (const auto i : index_range(_displacements))
      47          64 :     _displacements[i] = getVar("displacements", i);
      48          32 : }
      49             : 
      50             : ADRealVectorValue
      51       18612 : ADPenaltyVelocityContinuity::solidVelocity(const unsigned int qp) const
      52             : {
      53             :   ADRealVectorValue ret;
      54       55836 :   for (const auto i : index_range(_solid_velocities))
      55       37224 :     ret(i) = (*_solid_velocities[i])[qp];
      56       18612 :   return ret;
      57             : }
      58             : 
      59             : void
      60       22042 : ADPenaltyVelocityContinuity::computeResidual()
      61             : {
      62       22042 :   _qp_jumps.resize(_qrule->n_points());
      63       66126 :   for (auto & qp_jump : _qp_jumps)
      64             :     qp_jump = 0;
      65             : 
      66       66126 :   for (const auto qp : make_range(_qrule->n_points()))
      67       88168 :     _qp_jumps[qp] = _ad_JxW[qp] * _ad_coord[qp] * _penalty * (_velocity[qp] - solidVelocity(qp));
      68             : 
      69             :   // Fluid velocity residuals
      70             :   {
      71       22042 :     const auto & phi = _velocity_var->phiFace();
      72             :     const auto & dof_indices = _velocity_var->dofIndices();
      73             :     mooseAssert(phi.size() == dof_indices.size(), "These should be the same");
      74             : 
      75       22042 :     _residuals.resize(phi.size());
      76      172906 :     for (auto & r : _residuals)
      77      150864 :       r = 0;
      78             : 
      79      172906 :     for (const auto i : index_range(phi))
      80      452592 :       for (const auto qp : make_range(_qrule->n_points()))
      81      301728 :         _residuals[i] += phi[i][qp] * _qp_jumps[qp];
      82             : 
      83       22042 :     addResidualsAndJacobian(_assembly, _residuals, dof_indices, _velocity_var->scalingFactor());
      84             :   }
      85             : 
      86             :   // Displacement residuals
      87       66126 :   for (const auto dim : index_range(_displacements))
      88             :   {
      89       44084 :     const auto * const disp_var = _displacements[dim];
      90       44084 :     if (disp_var)
      91             :     {
      92             :       const auto & phi = disp_var->phiFaceNeighbor();
      93             :       const auto & dof_indices = disp_var->dofIndicesNeighbor();
      94             :       mooseAssert(phi.size() == dof_indices.size(), "These should be the same");
      95             : 
      96       44084 :       _residuals.resize(phi.size());
      97      194948 :       for (auto & r : _residuals)
      98      150864 :         r = 0;
      99             : 
     100      132252 :       for (const auto qp : make_range(_qrule->n_points()))
     101      389896 :         for (const auto i : index_range(phi))
     102      603456 :           _residuals[i] -= phi[i][qp] * _qp_jumps[qp](dim);
     103             : 
     104       44084 :       addResidualsAndJacobian(_assembly, _residuals, dof_indices, disp_var->scalingFactor());
     105             :     }
     106             :   }
     107       22042 : }