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 "VolumeJunction1Phase.h"
      11             : #include "FlowModelSinglePhase.h"
      12             : #include "THMMesh.h"
      13             : 
      14             : registerMooseObject("ThermalHydraulicsApp", VolumeJunction1Phase);
      15             : 
      16             : const unsigned int VolumeJunction1Phase::N_EQ = 5;
      17             : 
      18             : InputParameters
      19        2198 : VolumeJunction1Phase::validParams()
      20             : {
      21        2198 :   InputParameters params = FlowJunction1Phase::validParams();
      22             : 
      23        4396 :   params.addDeprecatedParam<bool>(
      24             :       "use_scalar_variables",
      25             :       "True if the junction variables are scalar variables",
      26             :       "Please remove this parameter; it no longer has any effect. The behavior corresponding to "
      27             :       "'use_scalar_variables = false' is now the only option.");
      28             : 
      29        4396 :   params.addRequiredParam<Real>("volume", "Volume of the junction [m^3]");
      30        4396 :   params.addRequiredParam<Point>("position", "Spatial position of the center of the junction [m]");
      31             : 
      32        4396 :   params.addParam<FunctionName>("initial_p", "Initial pressure [Pa]");
      33        4396 :   params.addParam<FunctionName>("initial_T", "Initial temperature [K]");
      34        4396 :   params.addParam<FunctionName>("initial_vel_x", "Initial velocity in x-direction [m/s]");
      35        4396 :   params.addParam<FunctionName>("initial_vel_y", "Initial velocity in y-direction [m/s]");
      36        4396 :   params.addParam<FunctionName>("initial_vel_z", "Initial velocity in z-direction [m/s]");
      37             : 
      38        4396 :   params.addParam<Real>("scaling_factor_rhoV", 1.0, "Scaling factor for rho*V [-]");
      39        4396 :   params.addParam<Real>("scaling_factor_rhouV", 1.0, "Scaling factor for rho*u*V [-]");
      40        4396 :   params.addParam<Real>("scaling_factor_rhovV", 1.0, "Scaling factor for rho*v*V [-]");
      41        4396 :   params.addParam<Real>("scaling_factor_rhowV", 1.0, "Scaling factor for rho*w*V [-]");
      42        4396 :   params.addParam<Real>("scaling_factor_rhoEV", 1.0, "Scaling factor for rho*E*V [-]");
      43             : 
      44        4396 :   params.addParam<Real>("K", 0., "Form loss factor [-]");
      45        4396 :   params.addParam<Real>("A_ref", "Reference area [m^2]");
      46             : 
      47        6594 :   params.addParam<bool>("apply_velocity_scaling",
      48        4396 :                         false,
      49             :                         "Set to true to apply the scaling to the normal velocity. See "
      50             :                         "documentation for more information.");
      51             : 
      52        4396 :   params.declareControllable("K");
      53        2198 :   params.addClassDescription("Junction between 1-phase flow channels that has a non-zero volume");
      54             : 
      55        2198 :   return params;
      56           0 : }
      57             : 
      58        1099 : VolumeJunction1Phase::VolumeJunction1Phase(const InputParameters & params)
      59             :   : FlowJunction1Phase(params),
      60             : 
      61        1099 :     _volume(getParam<Real>("volume")),
      62        2198 :     _position(getParam<Point>("position")),
      63             : 
      64        2198 :     _scaling_factor_rhoV(getParam<Real>("scaling_factor_rhoV")),
      65        2198 :     _scaling_factor_rhouV(getParam<Real>("scaling_factor_rhouV")),
      66        2198 :     _scaling_factor_rhovV(getParam<Real>("scaling_factor_rhovV")),
      67        2198 :     _scaling_factor_rhowV(getParam<Real>("scaling_factor_rhowV")),
      68        2198 :     _scaling_factor_rhoEV(getParam<Real>("scaling_factor_rhoEV")),
      69             : 
      70        1099 :     _rhoV_var_name(junctionVariableName("rhoV")),
      71        1099 :     _rhouV_var_name(junctionVariableName("rhouV")),
      72        1099 :     _rhovV_var_name(junctionVariableName("rhovV")),
      73        1099 :     _rhowV_var_name(junctionVariableName("rhowV")),
      74        1099 :     _rhoEV_var_name(junctionVariableName("rhoEV")),
      75        1099 :     _pressure_var_name(junctionVariableName("p")),
      76        1099 :     _temperature_var_name(junctionVariableName("T")),
      77        1099 :     _velocity_var_name(junctionVariableName("vel")),
      78             : 
      79        2198 :     _K(getParam<Real>("K")),
      80        4127 :     _A_ref(isParamValid("A_ref") ? getParam<Real>("A_ref") : _zero)
      81             : {
      82             :   // Note: 'A_ref' can be required by child classes
      83        2249 :   if (!params.isParamRequired("A_ref") && params.isParamSetByUser("A_ref") &&
      84        1201 :       !params.isParamSetByUser("K"))
      85           2 :     logWarning("Parameter 'A_ref' is specified, but 'K' is not specified, so the junction will "
      86             :                "behave as if there were no form loss.");
      87        1099 : }
      88             : 
      89             : void
      90        1099 : VolumeJunction1Phase::setupMesh()
      91             : {
      92        1099 :   FlowJunction1Phase::setupMesh();
      93             : 
      94             :   // Add a NodeElem to the mesh
      95        1099 :   auto * node = addNode(_position);
      96        1099 :   auto * elem = addNodeElement(node->id());
      97        1099 :   _junction_subdomain_id = mesh().getNextSubdomainId();
      98        1099 :   elem->subdomain_id() = _junction_subdomain_id;
      99        1099 :   setSubdomainInfo(_junction_subdomain_id, name());
     100             : 
     101             :   // Add coupling between the flow channel end elements and the NodeElem
     102             :   const auto & elem_ids = getConnectedElementIDs();
     103        3541 :   for (unsigned int i = 0; i < elem_ids.size(); i++)
     104        2442 :     getTHMProblem().augmentSparsity(elem_ids[i], elem->id());
     105        1099 : }
     106             : 
     107             : void
     108        1060 : VolumeJunction1Phase::check() const
     109             : {
     110        1060 :   FlowJunction1Phase::check();
     111             : 
     112             :   bool ics_set =
     113        2112 :       getTHMProblem().hasInitialConditionsFromFile() ||
     114        7272 :       (isParamValid("initial_p") && isParamValid("initial_T") && isParamValid("initial_vel_x") &&
     115        5188 :        isParamValid("initial_vel_y") && isParamValid("initial_vel_z"));
     116             : 
     117        1060 :   if (!ics_set && !_app.isRestarting())
     118             :   {
     119             :     // create a list of the missing IC parameters
     120             :     const std::vector<std::string> ic_params{
     121          24 :         "initial_p", "initial_T", "initial_vel_x", "initial_vel_y", "initial_vel_z"};
     122           4 :     std::ostringstream oss;
     123          24 :     for (const auto & ic_param : ic_params)
     124          20 :       if (!isParamValid(ic_param))
     125          20 :         oss << " " << ic_param;
     126             : 
     127           4 :     logError("The following initial condition parameters are missing:", oss.str());
     128           4 :   }
     129             : 
     130             :   // https://github.com/idaholab/moose/issues/28670
     131        1060 :   if (getTHMProblem().hasInitialConditionsFromFile() && libMesh::n_threads() > 1 &&
     132           0 :       _app.n_processors() > 1)
     133           0 :     mooseDocumentedError("moose",
     134             :                          28670,
     135             :                          "Using initial conditions from a file for VolumeJunction1Phase is "
     136             :                          "currently not tested for parallel threading.");
     137        1068 : }
     138             : 
     139             : void
     140        1074 : VolumeJunction1Phase::addVariables()
     141             : {
     142        1074 :   addJunctionVariable(true, _rhoV_var_name, _scaling_factor_rhoV);
     143        1074 :   addJunctionVariable(true, _rhouV_var_name, _scaling_factor_rhouV);
     144        1074 :   addJunctionVariable(true, _rhovV_var_name, _scaling_factor_rhovV);
     145        1074 :   addJunctionVariable(true, _rhowV_var_name, _scaling_factor_rhowV);
     146        1074 :   addJunctionVariable(true, _rhoEV_var_name, _scaling_factor_rhoEV);
     147             : 
     148        1074 :   addJunctionVariable(false, _pressure_var_name);
     149        1074 :   addJunctionVariable(false, _temperature_var_name);
     150        1074 :   addJunctionVariable(false, _velocity_var_name);
     151             : 
     152        7398 :   if (isParamValid("initial_p") && isParamValid("initial_T") && isParamValid("initial_vel_x") &&
     153        6324 :       isParamValid("initial_vel_y") && isParamValid("initial_vel_z"))
     154             :   {
     155        1050 :     addVolumeJunctionIC(_rhoV_var_name, "rhoV");
     156        1050 :     addVolumeJunctionIC(_rhouV_var_name, "rhouV");
     157        1050 :     addVolumeJunctionIC(_rhovV_var_name, "rhovV");
     158        1050 :     addVolumeJunctionIC(_rhowV_var_name, "rhowV");
     159        1050 :     addVolumeJunctionIC(_rhoEV_var_name, "rhoEV");
     160             : 
     161        1050 :     addVolumeJunctionIC(_pressure_var_name, "p");
     162        1050 :     addVolumeJunctionIC(_temperature_var_name, "T");
     163        2100 :     addVolumeJunctionIC(_velocity_var_name, "vel");
     164             :   }
     165        1074 : }
     166             : 
     167             : void
     168         419 : VolumeJunction1Phase::buildVolumeJunctionUserObject()
     169             : {
     170         419 :   ExecFlagEnum execute_on(MooseUtils::getDefaultExecFlagEnum());
     171        2095 :   execute_on = {EXEC_INITIAL, EXEC_LINEAR, EXEC_NONLINEAR};
     172             : 
     173             :   {
     174         419 :     const std::string class_name = "ADVolumeJunction1PhaseUserObject";
     175         419 :     InputParameters params = _factory.getValidParams(class_name);
     176         419 :     params.set<bool>("use_scalar_variables") = false;
     177         419 :     params.set<subdomain_id_type>("junction_subdomain_id") = _junction_subdomain_id;
     178         419 :     params.set<std::vector<BoundaryName>>("boundary") = _boundary_names;
     179         838 :     params.set<std::vector<Real>>("normals") = _normals;
     180         419 :     params.set<std::vector<processor_id_type>>("processor_ids") = getConnectedProcessorIDs();
     181         419 :     params.set<std::vector<UserObjectName>>("numerical_flux_names") = _numerical_flux_names;
     182         419 :     params.set<Real>("volume") = _volume;
     183        1257 :     params.set<std::vector<VariableName>>("A") = {FlowModel::AREA};
     184        1257 :     params.set<std::vector<VariableName>>("rhoA") = {FlowModelSinglePhase::RHOA};
     185        1257 :     params.set<std::vector<VariableName>>("rhouA") = {FlowModelSinglePhase::RHOUA};
     186        1257 :     params.set<std::vector<VariableName>>("rhoEA") = {FlowModelSinglePhase::RHOEA};
     187        1257 :     params.set<std::vector<VariableName>>("rhoV") = {_rhoV_var_name};
     188        1257 :     params.set<std::vector<VariableName>>("rhouV") = {_rhouV_var_name};
     189        1257 :     params.set<std::vector<VariableName>>("rhovV") = {_rhovV_var_name};
     190        1257 :     params.set<std::vector<VariableName>>("rhowV") = {_rhowV_var_name};
     191        1257 :     params.set<std::vector<VariableName>>("rhoEV") = {_rhoEV_var_name};
     192         419 :     params.set<Real>("K") = _K;
     193         419 :     params.set<Real>("A_ref") = _A_ref;
     194         419 :     params.set<UserObjectName>("fp") = _fp_name;
     195         838 :     params.set<bool>("apply_velocity_scaling") = getParam<bool>("apply_velocity_scaling");
     196         419 :     params.set<ExecFlagEnum>("execute_on") = execute_on;
     197         419 :     getTHMProblem().addUserObject(class_name, _junction_uo_name, params);
     198         419 :     connectObject(params, _junction_uo_name, "K");
     199         419 :   }
     200         838 : }
     201             : 
     202             : void
     203        1074 : VolumeJunction1Phase::addMooseObjects()
     204             : {
     205        1074 :   buildVolumeJunctionUserObject();
     206             : 
     207             :   // Add BC to each of the connected flow channels
     208        3453 :   for (std::size_t i = 0; i < _boundary_names.size(); i++)
     209             :   {
     210             :     const std::vector<NonlinearVariableName> var_names = {
     211        9516 :         FlowModelSinglePhase::RHOA, FlowModelSinglePhase::RHOUA, FlowModelSinglePhase::RHOEA};
     212        9516 :     for (std::size_t j = 0; j < var_names.size(); j++)
     213             :     {
     214        7137 :       const std::string class_name = "ADVolumeJunction1PhaseBC";
     215        7137 :       InputParameters params = _factory.getValidParams(class_name);
     216       21411 :       params.set<std::vector<BoundaryName>>("boundary") = {_boundary_names[i]};
     217        7137 :       params.set<Real>("normal") = _normals[i];
     218        7137 :       params.set<NonlinearVariableName>("variable") = var_names[j];
     219       14274 :       params.set<UserObjectName>("volume_junction_uo") = _junction_uo_name;
     220        7137 :       params.set<unsigned int>("connection_index") = i;
     221       21411 :       params.set<std::vector<VariableName>>("A_elem") = {FlowModel::AREA};
     222       21411 :       params.set<std::vector<VariableName>>("A_linear") = {FlowModel::AREA_LINEAR};
     223       21411 :       params.set<std::vector<VariableName>>("rhoA") = {FlowModelSinglePhase::RHOA};
     224       21411 :       params.set<std::vector<VariableName>>("rhouA") = {FlowModelSinglePhase::RHOUA};
     225       21411 :       params.set<std::vector<VariableName>>("rhoEA") = {FlowModelSinglePhase::RHOEA};
     226        7137 :       params.set<bool>("implicit") = getTHMProblem().getImplicitTimeIntegrationFlag();
     227       14274 :       getTHMProblem().addBoundaryCondition(
     228        7137 :           class_name, genName(name(), i, var_names[j] + ":" + class_name), params);
     229        7137 :     }
     230        2379 :   }
     231             : 
     232             :   // Add scalar kernels for the junction
     233        1074 :   std::vector<NonlinearVariableName> var_names(N_EQ);
     234        1074 :   var_names[RHOV_INDEX] = _rhoV_var_name;
     235        1074 :   var_names[RHOUV_INDEX] = _rhouV_var_name;
     236        1074 :   var_names[RHOVV_INDEX] = _rhovV_var_name;
     237        1074 :   var_names[RHOWV_INDEX] = _rhowV_var_name;
     238        1074 :   var_names[RHOEV_INDEX] = _rhoEV_var_name;
     239        6444 :   for (std::size_t i = 0; i < N_EQ; i++)
     240             :   {
     241             :     {
     242        5370 :       const std::string class_name = "ADTimeDerivative";
     243        5370 :       InputParameters params = _factory.getValidParams(class_name);
     244        5370 :       params.set<NonlinearVariableName>("variable") = var_names[i];
     245        5370 :       const std::string obj_name = genName(name(), var_names[i], "td");
     246        5370 :       params.set<std::vector<SubdomainName>>("block") = getSubdomainNames();
     247        5370 :       getTHMProblem().addKernel(class_name, obj_name, params);
     248        5370 :     }
     249             :     {
     250        5370 :       const std::string class_name = "ADVolumeJunctionAdvectionKernel";
     251        5370 :       InputParameters params = _factory.getValidParams(class_name);
     252        5370 :       params.set<NonlinearVariableName>("variable") = var_names[i];
     253       10740 :       params.set<UserObjectName>("volume_junction_uo") = _junction_uo_name;
     254        5370 :       params.set<unsigned int>("equation_index") = i;
     255        5370 :       const std::string obj_name = genName(name(), var_names[i], "vja_sk");
     256        5370 :       params.set<std::vector<SubdomainName>>("block") = getSubdomainNames();
     257        5370 :       getTHMProblem().addKernel(class_name, obj_name, params);
     258        5370 :     }
     259             :   }
     260             : 
     261             :   const std::vector<std::pair<std::string, VariableName>> quantities = {
     262        1074 :       {"pressure", _pressure_var_name},
     263        1074 :       {"temperature", _temperature_var_name},
     264        4296 :       {"speed", _velocity_var_name}};
     265        4296 :   for (const auto & quantity_and_name : quantities)
     266             :   {
     267        3222 :     const std::string class_name = "VolumeJunction1PhaseAux";
     268        3222 :     InputParameters params = _factory.getValidParams(class_name);
     269        6444 :     params.set<AuxVariableName>("variable") = quantity_and_name.second;
     270        3222 :     params.set<MooseEnum>("quantity") = quantity_and_name.first;
     271        3222 :     params.set<Real>("volume") = _volume;
     272        9666 :     params.set<std::vector<VariableName>>("rhoV") = {_rhoV_var_name};
     273        9666 :     params.set<std::vector<VariableName>>("rhouV") = {_rhouV_var_name};
     274        9666 :     params.set<std::vector<VariableName>>("rhovV") = {_rhovV_var_name};
     275        9666 :     params.set<std::vector<VariableName>>("rhowV") = {_rhowV_var_name};
     276        9666 :     params.set<std::vector<VariableName>>("rhoEV") = {_rhoEV_var_name};
     277        3222 :     params.set<UserObjectName>("fp") = _fp_name;
     278        6444 :     const std::string obj_name = genName(name(), quantity_and_name.first + "_aux");
     279        3222 :     params.set<std::vector<SubdomainName>>("block") = getSubdomainNames();
     280        3222 :     getTHMProblem().addAuxKernel(class_name, obj_name, params);
     281        3222 :   }
     282             : 
     283             :   // An error message results if there is any block without a material, so
     284             :   // until this restriction is removed, we must add a dummy material that
     285             :   // computes no material properties.
     286             :   {
     287        1074 :     const std::string class_name = "GenericConstantMaterial";
     288        1074 :     InputParameters params = _factory.getValidParams(class_name);
     289        1074 :     params.set<std::vector<SubdomainName>>("block") = getSubdomainNames();
     290        1074 :     params.set<std::vector<std::string>>("prop_names") = {};
     291        1074 :     params.set<std::vector<Real>>("prop_values") = {};
     292        2148 :     getTHMProblem().addMaterial(class_name, genName(name(), "dummy_mat"), params);
     293        1074 :   }
     294        2148 : }
     295             : 
     296             : std::string
     297        9859 : VolumeJunction1Phase::junctionVariableName(const std::string & var_base) const
     298             : {
     299        9859 :   return var_base;
     300             : }
     301             : 
     302             : void
     303        9629 : VolumeJunction1Phase::addJunctionVariable(bool is_nonlinear,
     304             :                                           const VariableName & var,
     305             :                                           Real scaling_factor)
     306             : {
     307        9629 :   auto & problem = getTHMProblem();
     308             : 
     309             :   const libMesh::FEType fe_type(CONSTANT, MONOMIAL);
     310        9629 :   const auto & subdomains = getSubdomainNames();
     311             : 
     312        9629 :   if (is_nonlinear)
     313        5370 :     problem.addSimVariable(is_nonlinear, var, fe_type, subdomains, scaling_factor);
     314             :   else
     315        4259 :     problem.addSimVariable(is_nonlinear, var, fe_type, subdomains);
     316        9629 : }
     317             : 
     318             : void
     319         970 : VolumeJunction1Phase::addJunctionIC(const VariableName & var, Real value)
     320             : {
     321         970 :   getTHMProblem().addConstantIC(var, value, getSubdomainNames());
     322         970 : }
     323             : 
     324             : void
     325        8400 : VolumeJunction1Phase::addVolumeJunctionIC(const VariableName & var, const std::string & quantity)
     326             : {
     327        8400 :   const std::string class_name = "VolumeJunction1PhaseIC";
     328        8400 :   InputParameters params = _factory.getValidParams(class_name);
     329        8400 :   params.set<std::vector<SubdomainName>>("block") = getSubdomainNames();
     330        8400 :   params.set<VariableName>("variable") = var;
     331        8400 :   params.set<MooseEnum>("quantity") = quantity;
     332       67200 :   params.applySpecificParameters(parameters(),
     333             :                                  {"initial_p",
     334             :                                   "initial_T",
     335             :                                   "initial_vel_x",
     336             :                                   "initial_vel_y",
     337             :                                   "initial_vel_z",
     338             :                                   "volume",
     339             :                                   "position"});
     340        8400 :   params.set<UserObjectName>("fluid_properties") = _fp_name;
     341       16800 :   getTHMProblem().addSimInitialCondition(class_name, genName(name(), var, "ic"), params);
     342       33600 : }