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 "CNSFVHLLCBase.h"
      11             : #include "NS.h"
      12             : #include "HLLCUserObject.h"
      13             : #include "SinglePhaseFluidProperties.h"
      14             : 
      15             : namespace nms = NS;
      16             : using MetaPhysicL::raw_value;
      17             : 
      18             : InputParameters
      19        1531 : CNSFVHLLCBase::validParams()
      20             : {
      21        1531 :   InputParameters params = FVFluxKernel::validParams();
      22        1531 :   params.addRequiredParam<UserObjectName>(nms::fluid, "Fluid properties userobject");
      23        1531 :   return params;
      24           0 : }
      25             : 
      26         818 : CNSFVHLLCBase::CNSFVHLLCBase(const InputParameters & params)
      27             :   : FVFluxKernel(params),
      28         818 :     _fluid(dynamic_cast<FEProblemBase *>(&_subproblem)
      29        1636 :                ->getUserObject<SinglePhaseFluidProperties>(nms::fluid)),
      30         818 :     _specific_internal_energy_elem(getADMaterialProperty<Real>(nms::specific_internal_energy)),
      31         818 :     _specific_internal_energy_neighbor(
      32             :         getNeighborADMaterialProperty<Real>(nms::specific_internal_energy)),
      33         818 :     _rho_et_elem(getADMaterialProperty<Real>(nms::total_energy_density)),
      34         818 :     _rho_et_neighbor(getNeighborADMaterialProperty<Real>(nms::total_energy_density)),
      35         818 :     _vel_elem(getADMaterialProperty<RealVectorValue>(nms::velocity)),
      36         818 :     _vel_neighbor(getNeighborADMaterialProperty<RealVectorValue>(nms::velocity)),
      37         818 :     _speed_elem(getADMaterialProperty<Real>(nms::speed)),
      38         818 :     _speed_neighbor(getNeighborADMaterialProperty<Real>(nms::speed)),
      39         818 :     _rho_elem(getADMaterialProperty<Real>(nms::density)),
      40         818 :     _rho_neighbor(getNeighborADMaterialProperty<Real>(nms::density)),
      41        1636 :     _pressure_elem(getADMaterialProperty<Real>(nms::pressure)),
      42         818 :     _pressure_neighbor(getNeighborADMaterialProperty<Real>(nms::pressure))
      43             : {
      44         818 : }
      45             : 
      46             : HLLCData
      47     2320372 : CNSFVHLLCBase::hllcData() const
      48             : {
      49     2320372 :   return {_fluid,
      50     2320372 :           _rho_elem[_qp],
      51     2320372 :           _rho_neighbor[_qp],
      52     2320372 :           _vel_elem[_qp],
      53     2320372 :           _vel_neighbor[_qp],
      54     2320372 :           _specific_internal_energy_elem[_qp],
      55     2320372 :           _specific_internal_energy_neighbor[_qp]};
      56             : }
      57             : 
      58             : std::array<ADReal, 3>
      59     2404058 : CNSFVHLLCBase::waveSpeed(const HLLCData & hllc_data, const ADRealVectorValue & normal)
      60             : {
      61     2404058 :   const ADReal & rho1 = hllc_data.rho_elem;
      62     2404058 :   const ADReal u1 = hllc_data.vel_elem.norm();
      63     2404058 :   const ADReal q1 = normal * hllc_data.vel_elem;
      64     2404058 :   const ADReal v1 = 1.0 / rho1;
      65     2404058 :   const ADReal & e1 = hllc_data.e_elem;
      66     2404058 :   const ADReal E1 = e1 + 0.5 * u1 * u1;
      67     2404058 :   const ADReal p1 = hllc_data.fluid.p_from_v_e(v1, e1);
      68     2404058 :   const ADReal H1 = E1 + p1 / rho1;
      69     2404058 :   const ADReal c1 = hllc_data.fluid.c_from_v_e(v1, e1);
      70             : 
      71     2404058 :   const ADReal & rho2 = hllc_data.rho_neighbor;
      72     2404058 :   const ADReal u2 = hllc_data.vel_neighbor.norm();
      73     2404058 :   const ADReal q2 = normal * hllc_data.vel_neighbor;
      74     2404058 :   const ADReal v2 = 1.0 / rho2;
      75     2404058 :   const ADReal & e2 = hllc_data.e_neighbor;
      76     2404058 :   const ADReal E2 = e2 + 0.5 * u2 * u2;
      77     2404058 :   const ADReal p2 = hllc_data.fluid.p_from_v_e(v2, e2);
      78     2404058 :   const ADReal H2 = E2 + p2 / rho2;
      79     2404058 :   const ADReal c2 = hllc_data.fluid.c_from_v_e(v2, e2);
      80             : 
      81             :   // compute Roe-averaged variables
      82     2404058 :   const ADReal sqrt_rho1 = std::sqrt(rho1);
      83     2404058 :   const ADReal sqrt_rho2 = std::sqrt(rho2);
      84     2404058 :   const ADReal u_roe = (sqrt_rho1 * u1 + sqrt_rho2 * u2) / (sqrt_rho1 + sqrt_rho2);
      85     2404058 :   const ADReal q_roe = (sqrt_rho1 * q1 + sqrt_rho2 * q2) / (sqrt_rho1 + sqrt_rho2);
      86     2404058 :   const ADReal H_roe = (sqrt_rho1 * H1 + sqrt_rho2 * H2) / (sqrt_rho1 + sqrt_rho2);
      87     4808116 :   const ADReal h_roe = H_roe - 0.5 * u_roe * u_roe;
      88     2404058 :   const ADReal rho_roe = std::sqrt(rho1 * rho2);
      89     2404058 :   const ADReal v_roe = 1.0 / rho_roe;
      90     2404058 :   const ADReal e_roe = hllc_data.fluid.e_from_v_h(v_roe, h_roe);
      91     2404058 :   const ADReal c_roe = hllc_data.fluid.c_from_v_e(v_roe, e_roe);
      92             : 
      93             :   // compute wave speeds
      94     2404058 :   ADReal SL = std::min(q1 - c1, q_roe - c_roe);
      95     2404058 :   ADReal SR = std::max(q2 + c2, q_roe + c_roe);
      96     2404058 :   ADReal SM = (rho2 * q2 * (SR - q2) - rho1 * q1 * (SL - q1) + p1 - p2) /
      97     2404058 :               (rho2 * (SR - q2) - rho1 * (SL - q1));
      98             : 
      99     2404058 :   return {{std::move(SL), std::move(SM), std::move(SR)}};
     100     2404058 : }