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 "ADSimpleTurbine1PhaseUserObject.h"
      11             : #include "SinglePhaseFluidProperties.h"
      12             : #include "VolumeJunction1Phase.h"
      13             : #include "THMIndicesVACE.h"
      14             : #include "ADNumericalFlux3EqnBase.h"
      15             : #include "Numerics.h"
      16             : 
      17             : registerMooseObject("ThermalHydraulicsApp", ADSimpleTurbine1PhaseUserObject);
      18             : 
      19             : InputParameters
      20         146 : ADSimpleTurbine1PhaseUserObject::validParams()
      21             : {
      22         146 :   InputParameters params = ADJunctionParallelChannels1PhaseUserObject::validParams();
      23         292 :   params.addRequiredParam<bool>("on", "Flag determining if turbine is operating or not");
      24         292 :   params.addRequiredParam<Real>("W_dot", "Power, [W]");
      25             : 
      26         292 :   params.addClassDescription("Computes and caches flux and residual vectors for a 1-phase turbine");
      27             : 
      28         292 :   params.declareControllable("W_dot on");
      29             : 
      30         146 :   return params;
      31           0 : }
      32             : 
      33          79 : ADSimpleTurbine1PhaseUserObject::ADSimpleTurbine1PhaseUserObject(const InputParameters & params)
      34             :   : ADJunctionParallelChannels1PhaseUserObject(params),
      35          79 :     _on(getParam<bool>("on")),
      36         237 :     _W_dot(getParam<Real>("W_dot"))
      37             : {
      38          79 : }
      39             : 
      40             : void
      41        2898 : ADSimpleTurbine1PhaseUserObject::computeFluxesAndResiduals(const unsigned int & c)
      42             : {
      43        2898 :   ADJunctionParallelChannels1PhaseUserObject::computeFluxesAndResiduals(c);
      44             : 
      45        2898 :   if ((c == 0) && _on)
      46             :   {
      47             :     const auto & rhouV = _cached_junction_var_values[VolumeJunction1Phase::RHOUV_INDEX];
      48             :     const auto & rhovV = _cached_junction_var_values[VolumeJunction1Phase::RHOVV_INDEX];
      49             :     const auto & rhowV = _cached_junction_var_values[VolumeJunction1Phase::RHOWV_INDEX];
      50             : 
      51        1214 :     const Point di = _dir[0];
      52             :     const ADRealVectorValue rhouV_vec(rhouV, rhovV, rhowV);
      53             : 
      54             :     // energy source
      55        1214 :     const ADReal S_E = _W_dot;
      56             : 
      57             :     // momentum source
      58        1214 :     const ADReal v_in = THM::v_from_rhoA_A(_rhoA[0], _A[0]);
      59             : 
      60        1214 :     const ADReal rhouA2 = _rhouA[0] * _rhouA[0];
      61        3642 :     const ADReal e_in = _rhoEA[0] / _rhoA[0] - 0.5 * rhouA2 / (_rhoA[0] * _rhoA[0]);
      62             : 
      63        1214 :     const ADReal cp = _fp.cp_from_v_e(v_in, e_in);
      64        1214 :     const ADReal cv = _fp.cv_from_v_e(v_in, e_in);
      65             :     const ADReal gamma = cp / cv;
      66        1214 :     const ADReal p_in = _fp.p_from_v_e(v_in, e_in);
      67        1214 :     const ADReal T_in = _fp.T_from_v_e(v_in, e_in);
      68        1214 :     const ADReal h_in = _fp.h_from_p_T(p_in, T_in);
      69             :     const ADReal delta_p =
      70        3642 :         p_in * (1 - std::pow((1 - _W_dot / _rhouA[0] / h_in), (gamma / (gamma - 1))));
      71             : 
      72        2428 :     const ADRealVectorValue S_M = delta_p * _A[0] * di;
      73             : 
      74        1214 :     _residual[VolumeJunction1Phase::RHOUV_INDEX] += S_M(0);
      75        1214 :     _residual[VolumeJunction1Phase::RHOVV_INDEX] += S_M(1);
      76        1214 :     _residual[VolumeJunction1Phase::RHOWV_INDEX] += S_M(2);
      77        1214 :     _residual[VolumeJunction1Phase::RHOEV_INDEX] += S_E;
      78             :   }
      79        2898 : }