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