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 "RDG3EqnMaterial.h"
      11             : #include "SinglePhaseFluidProperties.h"
      12             : 
      13             : registerMooseObject("ThermalHydraulicsApp", RDG3EqnMaterial);
      14             : 
      15             : InputParameters
      16           0 : RDG3EqnMaterial::validParams()
      17             : {
      18           0 :   InputParameters params = Material::validParams();
      19           0 :   params += SlopeReconstruction1DInterface<false>::validParams();
      20             : 
      21           0 :   params.addClassDescription(
      22             :       "Reconstructed solution values for the 1-D, 1-phase, variable-area Euler equations");
      23             : 
      24           0 :   params.addRequiredCoupledVar("A_elem", "Cross-sectional area, elemental");
      25           0 :   params.addRequiredCoupledVar("A_linear", "Cross-sectional area, linear");
      26           0 :   params.addRequiredCoupledVar("rhoA", "Conserved variable: rho*A");
      27           0 :   params.addRequiredCoupledVar("rhouA", "Conserved variable: rho*u*A");
      28           0 :   params.addRequiredCoupledVar("rhoEA", "Conserved variable: rho*E*A");
      29             : 
      30           0 :   params.addRequiredParam<MaterialPropertyName>("direction",
      31             :                                                 "Flow channel direction material property name");
      32             : 
      33           0 :   params.addRequiredParam<UserObjectName>("fluid_properties",
      34             :                                           "Name of fluid properties user object");
      35             : 
      36           0 :   return params;
      37           0 : }
      38             : 
      39           0 : RDG3EqnMaterial::RDG3EqnMaterial(const InputParameters & parameters)
      40             :   : Material(parameters),
      41             :     SlopeReconstruction1DInterface<false>(this),
      42             : 
      43           0 :     _A_avg(coupledValue("A_elem")),
      44           0 :     _A_linear(coupledValue("A_linear")),
      45           0 :     _rhoA_avg(coupledValue("rhoA")),
      46           0 :     _rhouA_avg(coupledValue("rhouA")),
      47           0 :     _rhoEA_avg(coupledValue("rhoEA")),
      48             : 
      49           0 :     _A_var(getVar("A_elem", 0)),
      50           0 :     _rhoA_var(getVar("rhoA", 0)),
      51           0 :     _rhouA_var(getVar("rhouA", 0)),
      52           0 :     _rhoEA_var(getVar("rhoEA", 0)),
      53             : 
      54           0 :     _dir(getMaterialProperty<RealVectorValue>("direction")),
      55             : 
      56           0 :     _rhoA(declareProperty<Real>("rhoA")),
      57           0 :     _rhouA(declareProperty<Real>("rhouA")),
      58           0 :     _rhoEA(declareProperty<Real>("rhoEA")),
      59             : 
      60           0 :     _fp(getUserObject<SinglePhaseFluidProperties>("fluid_properties"))
      61             : {
      62           0 : }
      63             : 
      64             : void
      65           0 : RDG3EqnMaterial::computeQpProperties()
      66             : {
      67             :   // Get the limited slopes of the primitive variables: {p, u, T}.
      68           0 :   const auto slopes = getElementSlopes(_current_elem);
      69           0 :   const Real p_slope = slopes[PRESSURE];
      70           0 :   const Real vel_slope = slopes[VELOCITY];
      71           0 :   const Real T_slope = slopes[TEMPERATURE];
      72             : 
      73             :   // compute primitive variables from the cell-average solution
      74           0 :   const Real rho_avg = _rhoA_avg[_qp] / _A_avg[_qp];
      75           0 :   const Real vel_avg = _rhouA_avg[_qp] / _rhoA_avg[_qp];
      76           0 :   const Real v_avg = 1.0 / rho_avg;
      77           0 :   const Real e_avg = _rhoEA_avg[_qp] / _rhoA_avg[_qp] - 0.5 * vel_avg * vel_avg;
      78           0 :   const Real p_avg = _fp.p_from_v_e(v_avg, e_avg);
      79           0 :   const Real T_avg = _fp.T_from_v_e(v_avg, e_avg);
      80             : 
      81             :   // apply slopes to primitive variables
      82           0 :   const Real delta_x = (_q_point[_qp] - _current_elem->vertex_average()) * _dir[_qp];
      83           0 :   const Real p = p_avg + p_slope * delta_x;
      84           0 :   const Real vel = vel_avg + vel_slope * delta_x;
      85           0 :   const Real T = T_avg + T_slope * delta_x;
      86             : 
      87             :   // compute reconstructed conserved variables
      88           0 :   const Real rho = _fp.rho_from_p_T(p, T);
      89           0 :   const Real e = _fp.e_from_p_rho(p, rho);
      90           0 :   const Real E = e + 0.5 * vel * vel;
      91             : 
      92           0 :   _rhoA[_qp] = rho * _A_linear[_qp];
      93           0 :   _rhouA[_qp] = _rhoA[_qp] * vel;
      94           0 :   _rhoEA[_qp] = _rhoA[_qp] * E;
      95           0 : }
      96             : 
      97             : std::vector<Real>
      98           0 : RDG3EqnMaterial::computeElementPrimitiveVariables(const Elem * elem) const
      99             : {
     100             :   // get the cell-average conserved variables
     101             :   Real A, rhoA, rhouA, rhoEA;
     102           0 :   if (_is_implicit)
     103             :   {
     104           0 :     A = _A_var->getElementalValue(elem);
     105           0 :     rhoA = _rhoA_var->getElementalValue(elem);
     106           0 :     rhouA = _rhouA_var->getElementalValue(elem);
     107           0 :     rhoEA = _rhoEA_var->getElementalValue(elem);
     108             :   }
     109             :   else
     110             :   {
     111           0 :     A = _A_var->getElementalValueOld(elem);
     112           0 :     rhoA = _rhoA_var->getElementalValueOld(elem);
     113           0 :     rhouA = _rhouA_var->getElementalValueOld(elem);
     114           0 :     rhoEA = _rhoEA_var->getElementalValueOld(elem);
     115             :   }
     116             : 
     117             :   // compute primitive variables
     118             : 
     119           0 :   const Real rho = rhoA / A;
     120           0 :   const Real vel = rhouA / rhoA;
     121           0 :   const Real v = 1.0 / rho;
     122           0 :   const Real e = rhoEA / rhoA - 0.5 * vel * vel;
     123             : 
     124           0 :   std::vector<Real> W(_n_slopes);
     125           0 :   W[PRESSURE] = _fp.p_from_v_e(v, e);
     126           0 :   W[VELOCITY] = vel;
     127           0 :   W[TEMPERATURE] = _fp.T_from_v_e(v, e);
     128             : 
     129           0 :   return W;
     130             : }