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 "WCNSFV2PSlipVelocityFunctorMaterial.h"
      11             : #include "INSFVVelocityVariable.h"
      12             : #include "MooseLinearVariableFV.h"
      13             : #include "Function.h"
      14             : #include "NS.h"
      15             : #include "FVKernel.h"
      16             : 
      17             : registerMooseObject("NavierStokesApp", WCNSFV2PSlipVelocityFunctorMaterial);
      18             : 
      19             : InputParameters
      20         854 : WCNSFV2PSlipVelocityFunctorMaterial::validParams()
      21             : {
      22         854 :   InputParameters params = FunctorMaterial::validParams();
      23         854 :   params.addClassDescription("Computes the slip velocity for two-phase mixture model.");
      24        1708 :   params.addRequiredCoupledVar("u", "The velocity in the x direction.");
      25        1708 :   params.addCoupledVar("v", "The velocity in the y direction.");
      26        1708 :   params.addCoupledVar("w", "The velocity in the z direction.");
      27         854 :   params.addRequiredParam<MooseFunctorName>(NS::density, "Continuous phase density.");
      28        1708 :   params.addRequiredParam<MooseFunctorName>("rho_d", "Dispersed phase density.");
      29         854 :   params.addRequiredParam<MooseFunctorName>(NS::mu, "Mixture Density");
      30         854 :   params.addParam<RealVectorValue>(
      31         854 :       "gravity", RealVectorValue(0, 0, 0), "Gravity acceleration vector");
      32        1708 :   params.addParam<Real>("force_value", 0.0, "Coefficient to multiply by the body force term");
      33        1708 :   params.addParam<FunctionName>("force_function", "0", "A function that describes the body force");
      34        1708 :   params.addParam<PostprocessorName>(
      35        1708 :       "force_postprocessor", 0, "A postprocessor whose value is multiplied by the body force");
      36         854 :   params.addParam<RealVectorValue>(
      37         854 :       "force_direction", RealVectorValue(1, 0, 0), "Gravitational acceleration vector");
      38        1708 :   params.addParam<MooseFunctorName>(
      39        1708 :       "linear_coef_name", 0.44, "Linear friction coefficient name as a material property");
      40        1708 :   params.addParam<MooseFunctorName>(
      41        1708 :       "particle_diameter", 1.0, "Diameter of particles in the dispersed phase.");
      42        1708 :   params.addParam<MooseFunctorName>("fd", 0.0, "Fraction dispersed phase.");
      43        1708 :   MooseEnum momentum_component("x=0 y=1 z=2");
      44        1708 :   params.addRequiredParam<MooseEnum>(
      45             :       "momentum_component",
      46             :       momentum_component,
      47             :       "The component of the momentum equation that this kernel applies to.");
      48        1708 :   params.addRequiredParam<MooseFunctorName>("slip_velocity_name", "the name of the slip velocity");
      49        2562 :   params.addParam<unsigned short>("ghost_layers",
      50        1708 :                                   3,
      51             :                                   "The number of layers of elements to ghost. With Rhie-Chow and "
      52             :                                   "the velocity gradient calculation below, we need 3");
      53        1708 :   params.addRelationshipManager(
      54             :       "ElementSideNeighborLayers",
      55             :       Moose::RelationshipManagerType::GEOMETRIC | Moose::RelationshipManagerType::ALGEBRAIC |
      56             :           Moose::RelationshipManagerType::COUPLING,
      57         792 :       [](const InputParameters & obj_params, InputParameters & rm_params) {
      58         792 :         rm_params.set<unsigned short>("layers") = obj_params.get<unsigned short>("ghost_layers");
      59         792 :       });
      60         854 :   return params;
      61         854 : }
      62             : 
      63         454 : WCNSFV2PSlipVelocityFunctorMaterial::WCNSFV2PSlipVelocityFunctorMaterial(
      64         454 :     const InputParameters & params)
      65             :   : FunctorMaterial(params),
      66         454 :     _dim(_subproblem.mesh().dimension()),
      67         454 :     _u_var(dynamic_cast<MooseVariableField<Real> *>(getFieldVar("u", 0))),
      68         908 :     _v_var(params.isParamValid("v") ? dynamic_cast<MooseVariableField<Real> *>(getFieldVar("v", 0))
      69             :                                     : nullptr),
      70         454 :     _w_var(params.isParamValid("w") ? dynamic_cast<MooseVariableField<Real> *>(getFieldVar("w", 0))
      71             :                                     : nullptr),
      72         454 :     _rho_mixture(getFunctor<ADReal>(NS::density)),
      73         908 :     _rho_d(getFunctor<ADReal>("rho_d")),
      74         454 :     _mu_mixture(getFunctor<ADReal>(NS::mu)),
      75         908 :     _gravity(getParam<RealVectorValue>("gravity")),
      76         908 :     _force_scale(getParam<Real>("force_value")),
      77         454 :     _force_function(getFunction("force_function")),
      78         454 :     _force_postprocessor(getPostprocessorValue("force_postprocessor")),
      79         908 :     _force_direction(getParam<RealVectorValue>("force_direction")),
      80         908 :     _linear_friction(getFunctor<ADReal>("linear_coef_name")),
      81         908 :     _particle_diameter(getFunctor<ADReal>("particle_diameter")),
      82        1362 :     _index(getParam<MooseEnum>("momentum_component"))
      83             : {
      84         454 :   if (!dynamic_cast<const INSFVVelocityVariable *>(_u_var) &&
      85         118 :       !dynamic_cast<const MooseLinearVariableFV<Real> *>(_u_var))
      86           0 :     paramError("u",
      87             :                "the u velocity must be an INSFVVelocityVariable or a MooseLinearVariableFVReal");
      88             : 
      89         454 :   if (_dim >= 2 && (!dynamic_cast<const INSFVVelocityVariable *>(_v_var) &&
      90         118 :                     !dynamic_cast<const MooseLinearVariableFV<Real> *>(_v_var)))
      91           0 :     paramError("v",
      92             :                "In two or more dimensions, the v velocity must be supplied and it must be an "
      93             :                "INSFVVelocityVariable or a MooseLinearVariableFVReal.");
      94             : 
      95         454 :   if (_dim >= 3 && (!dynamic_cast<const INSFVVelocityVariable *>(_w_var) &&
      96           0 :                     !dynamic_cast<const MooseLinearVariableFV<Real> *>(_w_var)))
      97           0 :     paramError("w",
      98             :                "In three-dimensions, the w velocity must be supplied and it must be an "
      99             :                "INSFVVelocityVariable or a MooseLinearVariableFVReal.");
     100             : 
     101             :   // Slip velocity advection term requires gradients
     102             :   // TODO: this could be set less often, keeping it false until the two phase mixture system is
     103             :   // solved
     104         454 :   if (auto u = dynamic_cast<MooseLinearVariableFV<Real> *>(_u_var))
     105             :     u->computeCellGradients();
     106         454 :   if (auto v = dynamic_cast<MooseLinearVariableFV<Real> *>(_v_var))
     107             :     v->computeCellGradients();
     108         454 :   if (auto w = dynamic_cast<MooseLinearVariableFV<Real> *>(_w_var))
     109             :     w->computeCellGradients();
     110             : 
     111        1816 :   addFunctorProperty<ADReal>(
     112         908 :       getParam<MooseFunctorName>("slip_velocity_name"),
     113    12953746 :       [this](const auto & r, const auto & t)
     114             :       {
     115             :         constexpr Real offset = 1e-15;
     116             : 
     117    12953746 :         const bool is_transient = _subproblem.isTransient();
     118    12953746 :         ADRealVectorValue term_advection(0, 0, 0);
     119    12953746 :         ADRealVectorValue term_transient(0, 0, 0);
     120             :         const ADRealVectorValue term_force(
     121    12953746 :             _force_scale * _force_postprocessor *
     122    12953746 :             _force_function.value(_t, _current_elem->vertex_average()) * _force_direction);
     123             : 
     124             :         // Adding transient term
     125             :         // TODO: add time derivative term to lienar FV variable
     126    12953746 :         if (is_transient && !dynamic_cast<MooseLinearVariableFV<Real> *>(_u_var))
     127             :         {
     128     5293734 :           term_transient(0) += _u_var->dot(r, t);
     129     5293734 :           if (_dim > 1)
     130     5293734 :             term_transient(1) += _v_var->dot(r, t);
     131     5293734 :           if (_dim > 2)
     132           0 :             term_transient(2) += _w_var->dot(r, t);
     133             :         }
     134             : 
     135             :         // Adding advection term
     136    12953746 :         const auto u_velocity = (*_u_var)(r, t);
     137    12953746 :         const auto u_grad = _u_var->gradient(r, t);
     138    12953746 :         term_advection(0) += u_velocity * u_grad(0);
     139    12953746 :         if (_dim > 1)
     140             :         {
     141    12953746 :           const auto v_velocity = (*_v_var)(r, t);
     142    12953746 :           const auto v_grad = _v_var->gradient(r, t);
     143    12953746 :           term_advection(0) += v_velocity * u_grad(1);
     144    12953746 :           term_advection(1) += u_velocity * v_grad(0) + v_velocity * v_grad(1);
     145    12953746 :           if (_dim > 2)
     146             :           {
     147           0 :             const auto w_velocity = (*_w_var)(r, t);
     148           0 :             const auto w_grad = _w_var->gradient(r, t);
     149           0 :             term_advection(0) += w_velocity * u_grad(2);
     150           0 :             term_advection(1) += w_velocity * v_grad(2);
     151           0 :             term_advection(2) +=
     152             :                 u_velocity * w_grad(0) + v_velocity * w_grad(1) + w_velocity * w_grad(2);
     153             :           }
     154             :         }
     155             : 
     156    12953746 :         const ADReal density_scaling = (_rho_d(r, t) - _rho_mixture(r, t)) / _rho_d(r, t);
     157             :         const ADReal flux_residual =
     158    12953746 :             density_scaling * (-term_transient - term_advection + _gravity + term_force)(_index);
     159             : 
     160    51814984 :         const ADReal relaxation_time =
     161    38861238 :             _rho_d(r, t) * Utility::pow<2>(_particle_diameter(r, t)) / (18.0 * _mu_mixture(r, t));
     162             : 
     163    25907492 :         const ADReal linear_friction_factor = _linear_friction(r, t) + offset;
     164             : 
     165    12953746 :         return relaxation_time / linear_friction_factor * flux_residual;
     166             :       });
     167         908 : }