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         100 : ADPenaltyVelocityContinuity::validParams()
      17             : {
      18         100 :   InputParameters params = InterfaceKernelBase::validParams();
      19         100 :   params.addClassDescription(
      20             :       "Enforces continuity of flux and continuity of solution via penalty across an interface.");
      21         200 :   params.addRequiredParam<Real>(
      22             :       "penalty",
      23             :       "The penalty that penalizes jump between primary and neighbor secondary variables.");
      24         200 :   params.addRequiredCoupledVar("fluid_velocity", "The fluid velocity");
      25         200 :   params.renameCoupledVar("neighbor_var", "displacements", "All the displacement variables");
      26         200 :   params.addRequiredCoupledVar("solid_velocities", "The solid velocity components");
      27         100 :   return params;
      28           0 : }
      29             : 
      30          56 : ADPenaltyVelocityContinuity::ADPenaltyVelocityContinuity(const InputParameters & parameters)
      31             :   : InterfaceKernelBase(parameters),
      32          56 :     _penalty(getParam<Real>("penalty")),
      33          56 :     _velocity_var(getVectorVar("fluid_velocity", 0)),
      34          56 :     _velocity(adCoupledVectorValue("fluid_velocity")),
      35          56 :     _ad_JxW(_assembly.adJxWFace()),
      36          56 :     _ad_coord(_assembly.adCoordTransformation())
      37             : {
      38         112 :   if (isCoupledConstant("fluid_velocity"))
      39           0 :     paramError("fluid_velocity", "The fluid velocity must be an actual variable");
      40             : 
      41         112 :   _solid_velocities.resize(coupledComponents("solid_velocities"));
      42         168 :   for (const auto i : index_range(_solid_velocities))
      43         112 :     _solid_velocities[i] = &adCoupledNeighborValue("solid_velocities", i);
      44             : 
      45         112 :   _displacements.resize(coupledComponents("displacements"));
      46         168 :   for (const auto i : index_range(_displacements))
      47         112 :     _displacements[i] = getVar("displacements", i);
      48          56 : }
      49             : 
      50             : ADRealVectorValue
      51       21744 : ADPenaltyVelocityContinuity::solidVelocity(const unsigned int qp) const
      52             : {
      53             :   ADRealVectorValue ret;
      54       65232 :   for (const auto i : index_range(_solid_velocities))
      55       43488 :     ret(i) = (*_solid_velocities[i])[qp];
      56       21744 :   return ret;
      57             : }
      58             : 
      59             : void
      60       26136 : ADPenaltyVelocityContinuity::computeResidual()
      61             : {
      62       26136 :   _qp_jumps.resize(_qrule->n_points());
      63       78408 :   for (auto & qp_jump : _qp_jumps)
      64             :     qp_jump = 0;
      65             : 
      66       78408 :   for (const auto qp : make_range(_qrule->n_points()))
      67      104544 :     _qp_jumps[qp] = _ad_JxW[qp] * _ad_coord[qp] * _penalty * (_velocity[qp] - solidVelocity(qp));
      68             : 
      69             :   // Fluid velocity residuals
      70             :   {
      71       26136 :     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       26136 :     _residuals.resize(phi.size());
      76      204696 :     for (auto & r : _residuals)
      77      178560 :       r = 0;
      78             : 
      79      204696 :     for (const auto i : index_range(phi))
      80      535680 :       for (const auto qp : make_range(_qrule->n_points()))
      81      357120 :         _residuals[i] += phi[i][qp] * _qp_jumps[qp];
      82             : 
      83       26136 :     addResidualsAndJacobian(_assembly, _residuals, dof_indices, _velocity_var->scalingFactor());
      84             :   }
      85             : 
      86             :   // Displacement residuals
      87       78408 :   for (const auto dim : index_range(_displacements))
      88             :   {
      89       52272 :     const auto * const disp_var = _displacements[dim];
      90       52272 :     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       52272 :       _residuals.resize(phi.size());
      97      230832 :       for (auto & r : _residuals)
      98      178560 :         r = 0;
      99             : 
     100      156816 :       for (const auto qp : make_range(_qrule->n_points()))
     101      461664 :         for (const auto i : index_range(phi))
     102      714240 :           _residuals[i] -= phi[i][qp] * _qp_jumps[qp](dim);
     103             : 
     104       52272 :       addResidualsAndJacobian(_assembly, _residuals, dof_indices, disp_var->scalingFactor());
     105             :     }
     106             :   }
     107       26136 : }