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             : // Navier-Stokes includes
      11             : #include "PorousPrimitiveVarMaterial.h"
      12             : #include "NS.h"
      13             : #include "NavierStokesMethods.h"
      14             : 
      15             : // FluidProperties includes
      16             : #include "SinglePhaseFluidProperties.h"
      17             : 
      18             : registerMooseObject("NavierStokesApp", PorousPrimitiveVarMaterial);
      19             : 
      20             : InputParameters
      21         628 : PorousPrimitiveVarMaterial::validParams()
      22             : {
      23         628 :   auto params = Material::validParams();
      24         628 :   params.addRequiredParam<UserObjectName>(NS::fluid, "fluid userobject");
      25         628 :   params.addRequiredCoupledVar(NS::pressure, "The pressure");
      26         628 :   params.addRequiredCoupledVar(NS::T_fluid, "The fluid temperature");
      27         628 :   params.addRequiredCoupledVar(NS::superficial_velocity_x, "The x-velocity times the porosity");
      28         628 :   params.addCoupledVar(NS::superficial_velocity_y, "The y-velocity times the porosity");
      29         628 :   params.addCoupledVar(NS::superficial_velocity_z, "The z-velocity times the porosity");
      30         628 :   params.addClassDescription("Provides access to variables for a primitive variable set "
      31             :                              "of pressure, temperature, and superficial velocity");
      32         628 :   params.addRequiredParam<MaterialPropertyName>(NS::porosity, "the porosity");
      33         628 :   return params;
      34           0 : }
      35             : 
      36         474 : PorousPrimitiveVarMaterial::PorousPrimitiveVarMaterial(const InputParameters & params)
      37             :   : Material(params),
      38         948 :     _fluid(UserObjectInterface::getUserObject<SinglePhaseFluidProperties>(NS::fluid)),
      39             :     // primitive variables
      40         474 :     _var_pressure(adCoupledValue(NS::pressure)),
      41         474 :     _grad_var_pressure(adCoupledGradient(NS::pressure)),
      42         474 :     _var_T_fluid(adCoupledValue(NS::T_fluid)),
      43         474 :     _grad_var_T_fluid(adCoupledGradient(NS::T_fluid)),
      44         474 :     _var_sup_vel_x(adCoupledValue(NS::superficial_velocity_x)),
      45         474 :     _grad_var_sup_vel_x(adCoupledGradient(NS::superficial_velocity_x)),
      46         474 :     _var_sup_vel_y(isCoupled(NS::superficial_velocity_y)
      47         474 :                        ? adCoupledValue(NS::superficial_velocity_y)
      48             :                        : _ad_zero),
      49         474 :     _grad_var_sup_vel_y(isCoupled(NS::superficial_velocity_y)
      50         474 :                             ? adCoupledGradient(NS::superficial_velocity_y)
      51             :                             : _ad_grad_zero),
      52         474 :     _var_sup_vel_z(isCoupled(NS::superficial_velocity_z)
      53         474 :                        ? adCoupledValue(NS::superficial_velocity_z)
      54             :                        : _ad_zero),
      55         474 :     _grad_var_sup_vel_z(isCoupled(NS::superficial_velocity_z)
      56         474 :                             ? adCoupledGradient(NS::superficial_velocity_z)
      57             :                             : _ad_grad_zero),
      58         474 :     _pressure_dot(_is_transient ? adCoupledDot(NS::pressure) : _ad_zero),
      59         474 :     _T_fluid_dot(_is_transient ? adCoupledDot(NS::T_fluid) : _ad_zero),
      60         474 :     _sup_vel_x_dot(_is_transient ? adCoupledDot(NS::superficial_velocity_x) : _ad_zero),
      61         498 :     _sup_vel_y_dot((isCoupled(NS::superficial_velocity_y) && _is_transient)
      62         498 :                        ? adCoupledDot(NS::superficial_velocity_y)
      63             :                        : _ad_zero),
      64         948 :     _sup_vel_z_dot((isCoupled(NS::superficial_velocity_z) && _is_transient)
      65         474 :                        ? adCoupledDot(NS::superficial_velocity_z)
      66             :                        : _ad_zero),
      67             :     // porosity
      68         474 :     _epsilon(getMaterialProperty<Real>(NS::porosity)),
      69             :     // properties: primitives
      70         474 :     _pressure(declareADProperty<Real>(NS::pressure)),
      71         474 :     _grad_pressure(declareADProperty<RealVectorValue>(NS::grad(NS::pressure))),
      72         474 :     _T_fluid(declareADProperty<Real>(NS::T_fluid)),
      73         474 :     _grad_T_fluid(declareADProperty<RealVectorValue>(NS::grad(NS::T_fluid))),
      74         474 :     _sup_vel_x(declareADProperty<Real>(NS::superficial_velocity_x)),
      75         474 :     _grad_sup_vel_x(declareADProperty<RealVectorValue>(NS::grad(NS::superficial_velocity_x))),
      76         474 :     _sup_vel_y(declareADProperty<Real>(NS::superficial_velocity_y)),
      77         474 :     _grad_sup_vel_y(declareADProperty<RealVectorValue>(NS::grad(NS::superficial_velocity_y))),
      78         474 :     _sup_vel_z(declareADProperty<Real>(NS::superficial_velocity_z)),
      79         474 :     _grad_sup_vel_z(declareADProperty<RealVectorValue>(NS::grad(NS::superficial_velocity_z))),
      80             :     // properties: for viz
      81         474 :     _rho(declareADProperty<Real>(NS::density)),
      82         474 :     _sup_rho_dot(declareADProperty<Real>(NS::time_deriv(NS::superficial_density))),
      83         474 :     _velocity(declareADProperty<RealVectorValue>(NS::velocity)),
      84         474 :     _vel_x(declareADProperty<Real>(NS::velocity_x)),
      85         474 :     _vel_y(declareADProperty<Real>(NS::velocity_y)),
      86         474 :     _vel_z(declareADProperty<Real>(NS::velocity_z)),
      87         474 :     _sup_mom_x(declareADProperty<Real>(NS::superficial_momentum_x)),
      88         474 :     _sup_mom_y(declareADProperty<Real>(NS::superficial_momentum_y)),
      89         474 :     _sup_mom_z(declareADProperty<Real>(NS::superficial_momentum_z)),
      90         474 :     _sup_mom_x_dot(declareADProperty<Real>(NS::time_deriv(NS::superficial_momentum_x))),
      91         474 :     _sup_mom_y_dot(declareADProperty<Real>(NS::time_deriv(NS::superficial_momentum_y))),
      92         474 :     _sup_mom_z_dot(declareADProperty<Real>(NS::time_deriv(NS::superficial_momentum_z))),
      93         474 :     _sup_rho_et_dot(declareADProperty<Real>(NS::time_deriv(NS::superficial_total_energy_density))),
      94         474 :     _mom(declareADProperty<RealVectorValue>(NS::momentum)),
      95         474 :     _mom_x(declareADProperty<Real>(NS::momentum_x)),
      96         474 :     _mom_y(declareADProperty<Real>(NS::momentum_y)),
      97         474 :     _mom_z(declareADProperty<Real>(NS::momentum_z)),
      98         474 :     _speed(declareADProperty<Real>(NS::speed)),
      99         474 :     _rho_et(declareADProperty<Real>(NS::total_energy_density)),
     100         474 :     _e(declareADProperty<Real>(NS::specific_internal_energy)),
     101         948 :     _ht(declareADProperty<Real>(NS::specific_total_enthalpy))
     102             : {
     103         474 :   if (_mesh.dimension() >= 2 && !isCoupled(NS::superficial_velocity_y))
     104           0 :     mooseError("You must couple in a superficial y-velocity when solving 2D or 3D problems.");
     105             : 
     106         474 :   if (_mesh.dimension() >= 3 && !isCoupled(NS::superficial_velocity_z))
     107           0 :     mooseError("You must couple in a superficial z-velocity when solving 3D problems.");
     108         474 : }
     109             : 
     110             : void
     111     3835280 : PorousPrimitiveVarMaterial::computeQpProperties()
     112             : {
     113             :   // Our primitive variable set
     114     3835280 :   _pressure[_qp] = _var_pressure[_qp];
     115     3835280 :   _grad_pressure[_qp] = _grad_var_pressure[_qp];
     116     3835280 :   _T_fluid[_qp] = _var_T_fluid[_qp];
     117     3835280 :   _grad_T_fluid[_qp] = _grad_var_T_fluid[_qp];
     118             :   const VectorValue<ADReal> superficial_velocity = {
     119     3835280 :       _var_sup_vel_x[_qp], _var_sup_vel_y[_qp], _var_sup_vel_z[_qp]};
     120     3835280 :   _sup_vel_x[_qp] = superficial_velocity(0);
     121     3835280 :   _sup_vel_y[_qp] = superficial_velocity(1);
     122     3835280 :   _sup_vel_z[_qp] = superficial_velocity(2);
     123     3835280 :   _grad_sup_vel_x[_qp] = _grad_var_sup_vel_x[_qp];
     124     3835280 :   _grad_sup_vel_y[_qp] = _grad_var_sup_vel_y[_qp];
     125     3835280 :   _grad_sup_vel_z[_qp] = _grad_var_sup_vel_z[_qp];
     126             : 
     127             :   ADReal drho_dp, drho_dT;
     128     3835280 :   _fluid.rho_from_p_T(_pressure[_qp], _T_fluid[_qp], _rho[_qp], drho_dp, drho_dT);
     129     7670560 :   const auto rho_dot = drho_dp * _pressure_dot[_qp] + drho_dT * _T_fluid_dot[_qp];
     130     7670560 :   _sup_rho_dot[_qp] = _epsilon[_qp] * rho_dot;
     131             : 
     132     7670560 :   _velocity[_qp] = superficial_velocity / _epsilon[_qp];
     133     3835280 :   _vel_x[_qp] = _velocity[_qp](0);
     134     3835280 :   _vel_y[_qp] = _velocity[_qp](1);
     135     3835280 :   _vel_z[_qp] = _velocity[_qp](2);
     136     3835280 :   const auto superficial_momentum = superficial_velocity * _rho[_qp];
     137     3835280 :   _sup_mom_x[_qp] = superficial_momentum(0);
     138     3835280 :   _sup_mom_y[_qp] = superficial_momentum(1);
     139     3835280 :   _sup_mom_z[_qp] = superficial_momentum(2);
     140    11505840 :   _sup_mom_x_dot[_qp] = _rho[_qp] * _sup_vel_x_dot[_qp] + rho_dot * _sup_vel_x[_qp];
     141    11505840 :   _sup_mom_y_dot[_qp] = _rho[_qp] * _sup_vel_y_dot[_qp] + rho_dot * _sup_vel_y[_qp];
     142    11505840 :   _sup_mom_z_dot[_qp] = _rho[_qp] * _sup_vel_z_dot[_qp] + rho_dot * _sup_vel_z[_qp];
     143             : 
     144     3835280 :   const auto v = 1. / _rho[_qp];
     145    11505840 :   const auto v_dot = -rho_dot / (_rho[_qp] * _rho[_qp]);
     146             :   ADReal de_dT, de_dv;
     147     3835280 :   _fluid.e_from_T_v(_T_fluid[_qp], v, _e[_qp], de_dT, de_dv);
     148     3835280 :   const auto e_dot = de_dT * _T_fluid_dot[_qp] + de_dv * v_dot;
     149     7670560 :   const auto et = _e[_qp] + _velocity[_qp] * _velocity[_qp] / 2.;
     150             :   const auto velocity_dot =
     151     3835280 :       VectorValue<ADReal>(_sup_vel_x_dot[_qp], _sup_vel_y_dot[_qp], _sup_vel_z_dot[_qp]) /
     152     3835280 :       _epsilon[_qp];
     153     3835280 :   const auto et_dot = e_dot + _velocity[_qp] * velocity_dot;
     154    11505840 :   _sup_rho_et_dot[_qp] = _epsilon[_qp] * (rho_dot * et + et_dot * _rho[_qp]);
     155             : 
     156     7670560 :   _mom_x[_qp] = _sup_mom_x[_qp] / _epsilon[_qp];
     157     7670560 :   _mom_y[_qp] = _sup_mom_y[_qp] / _epsilon[_qp];
     158     7670560 :   _mom_z[_qp] = _sup_mom_z[_qp] / _epsilon[_qp];
     159     3835280 :   _mom[_qp] = {_mom_x[_qp], _mom_y[_qp], _mom_z[_qp]};
     160             : 
     161     3835280 :   _speed[_qp] = NS::computeSpeed<ADReal>(_velocity[_qp]);
     162             : 
     163     7670560 :   _rho_et[_qp] = _rho[_qp] * et;
     164     7670560 :   _ht[_qp] = et + _pressure[_qp] / _rho[_qp];
     165     3835280 : }