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 "BoundaryFlux3EqnBC.h"
      11             : #include "MooseVariable.h"
      12             : #include "THMIndicesVACE.h"
      13             : 
      14             : registerMooseObject("ThermalHydraulicsApp", BoundaryFlux3EqnBC);
      15             : 
      16             : InputParameters
      17           0 : BoundaryFlux3EqnBC::validParams()
      18             : {
      19           0 :   InputParameters params = OneDIntegratedBC::validParams();
      20             : 
      21           0 :   params.addClassDescription(
      22             :       "Boundary conditions 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<UserObjectName>("boundary_flux", "Name of boundary flux user object");
      31             : 
      32           0 :   return params;
      33           0 : }
      34             : 
      35           0 : BoundaryFlux3EqnBC::BoundaryFlux3EqnBC(const InputParameters & parameters)
      36             :   : OneDIntegratedBC(parameters),
      37             : 
      38           0 :     _A_elem(coupledValue("A_elem")),
      39           0 :     _A_linear(coupledValue("A_linear")),
      40             : 
      41           0 :     _rhoA(coupledValue("rhoA")),
      42           0 :     _rhouA(coupledValue("rhouA")),
      43           0 :     _rhoEA(coupledValue("rhoEA")),
      44             : 
      45           0 :     _rhoA_var(coupled("rhoA")),
      46           0 :     _rhouA_var(coupled("rhouA")),
      47           0 :     _rhoEA_var(coupled("rhoEA")),
      48             : 
      49           0 :     _jmap(getIndexMapping()),
      50           0 :     _equation_index(_jmap.at(_var.number())),
      51             : 
      52           0 :     _flux(getUserObject<BoundaryFluxBase>("boundary_flux"))
      53             : {
      54           0 : }
      55             : 
      56             : Real
      57           0 : BoundaryFlux3EqnBC::computeQpResidual()
      58             : {
      59           0 :   const std::vector<Real> U = {_rhoA[_qp], _rhouA[_qp], _rhoEA[_qp], _A_elem[_qp]};
      60           0 :   const auto & flux = _flux.getFlux(_current_side, _current_elem->id(), U, _normals[_qp]);
      61             : 
      62             :   // Note that the ratio A_linear / A_elem is necessary because A_elem is passed
      63             :   // to the flux function, but A_linear is to be used on the boundary.
      64           0 :   return flux[_equation_index] * _A_linear[_qp] / _A_elem[_qp] * _normal * _test[_i][_qp];
      65             : }
      66             : 
      67             : Real
      68           0 : BoundaryFlux3EqnBC::computeQpJacobian()
      69             : {
      70           0 :   return computeQpOffDiagJacobian(_var.number());
      71             : }
      72             : 
      73             : Real
      74           0 : BoundaryFlux3EqnBC::computeQpOffDiagJacobian(unsigned int jvar)
      75             : {
      76           0 :   const std::vector<Real> U = {_rhoA[_qp], _rhouA[_qp], _rhoEA[_qp], _A_elem[_qp]};
      77           0 :   const auto & J = _flux.getJacobian(_current_side, _current_elem->id(), U, {_normal, 0, 0});
      78             : 
      79           0 :   return J(_equation_index, _jmap.at(jvar)) * _A_linear[_qp] / _A_elem[_qp] * _normal *
      80           0 :          _phi[_j][_qp] * _test[_i][_qp];
      81             : }
      82             : 
      83             : std::map<unsigned int, unsigned int>
      84           0 : BoundaryFlux3EqnBC::getIndexMapping() const
      85             : {
      86             :   std::map<unsigned int, unsigned int> jmap;
      87           0 :   jmap.insert(std::pair<unsigned int, unsigned int>(_rhoA_var, THMVACE1D::MASS));
      88           0 :   jmap.insert(std::pair<unsigned int, unsigned int>(_rhouA_var, THMVACE1D::MOMENTUM));
      89           0 :   jmap.insert(std::pair<unsigned int, unsigned int>(_rhoEA_var, THMVACE1D::ENERGY));
      90             : 
      91           0 :   return jmap;
      92             : }