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 "INSFVWallFunctionBC.h"
      11             : #include "NavierStokesMethods.h"
      12             : #include "NS.h"
      13             : 
      14             : registerMooseObject("NavierStokesApp", INSFVWallFunctionBC);
      15             : 
      16             : InputParameters
      17         122 : INSFVWallFunctionBC::validParams()
      18             : {
      19         122 :   InputParameters params = INSFVNaturalFreeSlipBC::validParams();
      20         122 :   params.addClassDescription("Implements a wall shear BC for the momentum equation based on "
      21             :                              "algebraic standard velocity wall functions.");
      22         244 :   params.addRequiredParam<MooseFunctorName>("u", "The velocity in the x direction.");
      23         244 :   params.addParam<MooseFunctorName>("v", "The velocity in the y direction.");
      24         244 :   params.addParam<MooseFunctorName>("w", "The velocity in the z direction.");
      25         122 :   params.addRequiredParam<MooseFunctorName>(NS::density, "fluid density");
      26         244 :   params.addRequiredParam<MooseFunctorName>("mu", "Dynamic viscosity");
      27         122 :   return params;
      28           0 : }
      29             : 
      30          72 : INSFVWallFunctionBC::INSFVWallFunctionBC(const InputParameters & params)
      31             :   : INSFVNaturalFreeSlipBC(params),
      32          72 :     _dim(_subproblem.mesh().dimension()),
      33         144 :     _u(getFunctor<ADReal>("u")),
      34         288 :     _v(isParamValid("v") ? &getFunctor<ADReal>("v") : nullptr),
      35         144 :     _w(isParamValid("w") ? &getFunctor<ADReal>("w") : nullptr),
      36          72 :     _rho(getFunctor<ADReal>(NS::density)),
      37         216 :     _mu(getFunctor<ADReal>("mu"))
      38             : {
      39          72 :   if (_rc_uo.segregated())
      40           0 :     mooseError("Wall sheer stress enforcement based wall functions are not supported with "
      41             :                "segregated solution approaches!");
      42          72 : }
      43             : 
      44             : ADReal
      45           0 : INSFVWallFunctionBC::computeSegregatedContribution()
      46             : {
      47           0 :   mooseError("Sheer-stress-based wall function enforcement not supported for segregated solvers.");
      48             :   return 0.0;
      49             : }
      50             : 
      51             : ADReal
      52      166400 : INSFVWallFunctionBC::computeStrongResidual()
      53             : {
      54             :   // Get the velocity vector
      55      166400 :   const FaceInfo & fi = *_face_info;
      56             :   const Elem & elem = fi.elem();
      57      166400 :   const Moose::ElemArg elem_arg{&elem, false};
      58      166400 :   const auto state = determineState();
      59      166400 :   ADRealVectorValue velocity(_u(elem_arg, state));
      60      166400 :   if (_v)
      61      166400 :     velocity(1) = (*_v)(elem_arg, state);
      62      166400 :   if (_w)
      63           0 :     velocity(2) = (*_w)(elem_arg, state);
      64             : 
      65             :   // Compute the velocity magnitude (parallel_speed) and
      66             :   // direction of the tangential velocity component (parallel_dir)
      67      166400 :   ADReal dist = std::abs((fi.elemCentroid() - fi.faceCentroid()) * _normal);
      68      166400 :   ADReal perpendicular_speed = velocity * _normal;
      69      166400 :   ADRealVectorValue parallel_velocity = velocity - perpendicular_speed * _normal;
      70      166400 :   ADReal parallel_speed = parallel_velocity.norm();
      71      332800 :   _a = 1 / parallel_speed;
      72             : 
      73      166400 :   if (parallel_speed.value() < 1e-7)
      74           0 :     return 0;
      75             : 
      76      166400 :   if (!std::isfinite(parallel_speed.value()))
      77           0 :     return parallel_speed;
      78             : 
      79             :   // Compute the friction velocity and the wall shear stress
      80      166400 :   const auto rho = _rho(makeElemArg(&elem), state);
      81             :   ADReal u_star =
      82      166400 :       NS::findUStar<ADReal>(_mu(makeElemArg(&elem), state), rho, parallel_speed, dist.value());
      83      166400 :   ADReal tau = u_star * u_star * rho;
      84      166400 :   _a *= tau;
      85             : 
      86             :   // Compute the shear stress component for this momentum equation
      87      166400 :   if (_index == 0)
      88             :     return _a * parallel_velocity(0);
      89       83200 :   else if (_index == 1)
      90             :     return _a * parallel_velocity(1);
      91             :   else
      92             :     return _a * parallel_velocity(2);
      93             : }
      94             : 
      95             : void
      96      166400 : INSFVWallFunctionBC::gatherRCData(const FaceInfo & fi)
      97             : {
      98      166400 :   _face_info = &fi;
      99      166400 :   _face_type = fi.faceType(std::make_pair(_var.number(), _var.sys().number()));
     100      166400 :   _normal = fi.normal();
     101             : 
     102             :   // Fill-in the coefficient _a (but without multiplication by A)
     103      166400 :   const auto strong_resid = computeStrongResidual();
     104             : 
     105      332800 :   _rc_uo.addToA((_face_type == FaceInfo::VarFaceNeighbors::ELEM) ? &fi.elem() : fi.neighborPtr(),
     106      166400 :                 _index,
     107      166400 :                 _a * (fi.faceArea() * fi.faceCoord()));
     108             : 
     109      332800 :   addResidualAndJacobian(strong_resid * (fi.faceArea() * fi.faceCoord()));
     110      166400 : }