From 107df1c51d70df624a6f94e1995c5a2a14fefc9e Mon Sep 17 00:00:00 2001 From: Guillaume Giudicelli Date: Sat, 30 Nov 2024 18:21:04 -0700 Subject: [PATCH 01/16] Add a base class for fluid heat transfer physics refs #29175 --- .../physics/WCNSFVFluidHeatTransferPhysics.h | 87 +----- .../WCNSFVFluidHeatTransferPhysicsBase.h | 116 ++++++++ .../physics/WCNSFVFluidHeatTransferPhysics.C | 238 +--------------- .../WCNSFVFluidHeatTransferPhysicsBase.C | 253 ++++++++++++++++++ 4 files changed, 384 insertions(+), 310 deletions(-) create mode 100644 modules/navier_stokes/include/physics/WCNSFVFluidHeatTransferPhysicsBase.h create mode 100644 modules/navier_stokes/src/physics/WCNSFVFluidHeatTransferPhysicsBase.C diff --git a/modules/navier_stokes/include/physics/WCNSFVFluidHeatTransferPhysics.h b/modules/navier_stokes/include/physics/WCNSFVFluidHeatTransferPhysics.h index 9a452fbd9961..3e4143eb3a32 100644 --- a/modules/navier_stokes/include/physics/WCNSFVFluidHeatTransferPhysics.h +++ b/modules/navier_stokes/include/physics/WCNSFVFluidHeatTransferPhysics.h @@ -9,102 +9,37 @@ #pragma once -#include "NavierStokesPhysicsBase.h" -#include "WCNSFVCoupledAdvectionPhysicsHelper.h" +#include "WCNSFVFluidHeatTransferPhysicsBase.h" #include "NS.h" /** * Creates all the objects needed to solve the Navier Stokes energy equation */ -class WCNSFVFluidHeatTransferPhysics final : public NavierStokesPhysicsBase, - public WCNSFVCoupledAdvectionPhysicsHelper +class WCNSFVFluidHeatTransferPhysics final : public WCNSFVFluidHeatTransferPhysicsBase { public: static InputParameters validParams(); WCNSFVFluidHeatTransferPhysics(const InputParameters & parameters); - /// Get the name of the fluid temperature variable - const NonlinearVariableName & getFluidTemperatureName() const { return _fluid_temperature_name; } - - /// Get the name of the specific heat material property - const MooseFunctorName & getSpecificHeatName() const { return _specific_heat_name; } - MooseFunctorName getSpecificEnthalpyName() const { return NS::specific_enthalpy; } - const std::vector & getThermalConductivityName() const - { - return _thermal_conductivity_name; - } - - /// Get the ambient convection parameters for parameter checking - const std::vector> & getAmbientConvectionBlocks() const - { - return _ambient_convection_blocks; - } - /// Name of the ambient convection heat transfer coefficients for each block-group - const std::vector & getAmbientConvectionHTCs() const - { - return _ambient_convection_alpha; - } - - /// Whether the physics is actually creating the heat equation - bool hasEnergyEquation() const { return _has_energy_equation; } - protected: private: - void actOnAdditionalTasks() override; - void addSolverVariables() override; - void addInitialConditions() override; - void addFVKernels() override; - void addFVBCs() override; - void addMaterials() override; - - unsigned short getNumberAlgebraicGhostingLayersNeeded() const override; + virtual void addSolverVariables() override; /** * Functions adding kernels for the incompressible / weakly compressible energy equation * If the material properties are not constant, some of these can be used for * weakly-compressible simulations as well. */ - void addINSEnergyTimeKernels(); - void addWCNSEnergyTimeKernels(); - void addINSEnergyHeatConductionKernels(); - void addINSEnergyAdvectionKernels(); - void addINSEnergyAmbientConvection(); - void addINSEnergyExternalHeatSource(); + void addINSEnergyTimeKernels() override; + void addWCNSEnergyTimeKernels() override; + void addINSEnergyHeatConductionKernels() override; + void addINSEnergyAdvectionKernels() override; + void addINSEnergyAmbientConvection() override; + void addINSEnergyExternalHeatSource() override; /// Functions adding boundary conditions for the incompressible simulation. /// These are used for weakly-compressible simulations as well. - void addINSEnergyInletBC(); - void addINSEnergyWallBC(); - - /// Process thermal conductivity (multiple functor input options are available). - /// Return true if we have vector thermal conductivity and false if scalar - bool processThermalConductivity(); - - /// A boolean to help compatibility with the old Modules/NavierStokesFV syntax - const bool _has_energy_equation; - /// Fluid temperature name - NonlinearVariableName _fluid_temperature_name; - /// Name of the specific heat material property - MooseFunctorName _specific_heat_name; - /// Vector of subdomain groups where we want to have different thermal conduction - std::vector> _thermal_conductivity_blocks; - /// Name of the thermal conductivity functor for each block-group - std::vector _thermal_conductivity_name; - - /// Vector of subdomain groups where we want to have different ambient convection - std::vector> _ambient_convection_blocks; - /// Name of the ambient convection heat transfer coefficients for each block-group - std::vector _ambient_convection_alpha; - /// Name of the solid domain temperature for each block-group - std::vector _ambient_temperature; - - /// Energy inlet boundary types - MultiMooseEnum _energy_inlet_types; - /// Functors describing the inlet boundary values. See energy_inlet_types for what the functors actually represent - std::vector _energy_inlet_functors; - /// Energy wall boundary types - MultiMooseEnum _energy_wall_types; - /// Functors describing the wall boundary values. See energy_wall_types for what the functors actually represent - std::vector _energy_wall_functors; + void addINSEnergyInletBC() override; + void addINSEnergyWallBC() override; }; diff --git a/modules/navier_stokes/include/physics/WCNSFVFluidHeatTransferPhysicsBase.h b/modules/navier_stokes/include/physics/WCNSFVFluidHeatTransferPhysicsBase.h new file mode 100644 index 000000000000..2465dc6cb392 --- /dev/null +++ b/modules/navier_stokes/include/physics/WCNSFVFluidHeatTransferPhysicsBase.h @@ -0,0 +1,116 @@ +//* This file is part of the MOOSE framework +//* https://www.mooseframework.org +//* +//* All rights reserved, see COPYRIGHT for full restrictions +//* https://github.com/idaholab/moose/blob/master/COPYRIGHT +//* +//* Licensed under LGPL 2.1, please see LICENSE for details +//* https://www.gnu.org/licenses/lgpl-2.1.html + +#pragma once + +#include "NavierStokesPhysicsBase.h" +#include "WCNSFVCoupledAdvectionPhysicsHelper.h" +#include "NS.h" + +#define registerWCNSFVFluidHeatTransferPhysicsBaseTasks(app_name, derived_name) \ + registerMooseAction(app_name, derived_name, "get_turbulence_physics"); \ + registerMooseAction(app_name, derived_name, "add_variable"); \ + registerMooseAction(app_name, derived_name, "add_ic"); \ + registerMooseAction(app_name, derived_name, "add_fv_kernel"); \ + registerMooseAction(app_name, derived_name, "add_fv_bc"); \ + registerMooseAction(app_name, derived_name, "add_material") + +/** + * Creates all the objects needed to solve the Navier Stokes energy equation + */ +class WCNSFVFluidHeatTransferPhysicsBase : public NavierStokesPhysicsBase, + public WCNSFVCoupledAdvectionPhysicsHelper +{ +public: + static InputParameters validParams(); + + WCNSFVFluidHeatTransferPhysicsBase(const InputParameters & parameters); + + /// Get the name of the fluid temperature variable + const NonlinearVariableName & getFluidTemperatureName() const { return _fluid_temperature_name; } + + /// Get the name of the specific heat material property + const MooseFunctorName & getSpecificHeatName() const { return _specific_heat_name; } + MooseFunctorName getSpecificEnthalpyName() const { return NS::specific_enthalpy; } + const std::vector & getThermalConductivityName() const + { + return _thermal_conductivity_name; + } + + /// Get the ambient convection parameters for parameter checking + const std::vector> & getAmbientConvectionBlocks() const + { + return _ambient_convection_blocks; + } + /// Name of the ambient convection heat transfer coefficients for each block-group + const std::vector & getAmbientConvectionHTCs() const + { + return _ambient_convection_alpha; + } + + /// Whether the physics is actually creating the heat equation + bool hasEnergyEquation() const { return _has_energy_equation; } + +protected: + void actOnAdditionalTasks() override; + void addInitialConditions() override; + void addFVKernels() override; + void addFVBCs() override; + void addMaterials() override; + + unsigned short getNumberAlgebraicGhostingLayersNeeded() const override; + + /** + * Functions adding kernels for the incompressible / weakly compressible energy equation + * If the material properties are not constant, some of these can be used for + * weakly-compressible simulations as well. + */ + virtual void addINSEnergyTimeKernels() = 0; + virtual void addWCNSEnergyTimeKernels() = 0; + virtual void addINSEnergyHeatConductionKernels() = 0; + virtual void addINSEnergyAdvectionKernels() = 0; + virtual void addINSEnergyAmbientConvection() = 0; + virtual void addINSEnergyExternalHeatSource() = 0; + + /// Functions adding boundary conditions for the incompressible simulation. + /// These are used for weakly-compressible simulations as well. + virtual void addINSEnergyInletBC() = 0; + virtual void addINSEnergyWallBC() = 0; + + /// Process thermal conductivity (multiple functor input options are available). + /// Return true if we have vector thermal conductivity and false if scalar + bool processThermalConductivity(); + + /// A boolean to help compatibility with the old Modules/NavierStokesFV syntax + const bool _has_energy_equation; + /// Fluid temperature name + NonlinearVariableName _fluid_temperature_name; + /// Name of the specific heat material property + MooseFunctorName _specific_heat_name; + /// Vector of subdomain groups where we want to have different thermal conduction + std::vector> _thermal_conductivity_blocks; + /// Name of the thermal conductivity functor for each block-group + std::vector _thermal_conductivity_name; + + /// Vector of subdomain groups where we want to have different ambient convection + std::vector> _ambient_convection_blocks; + /// Name of the ambient convection heat transfer coefficients for each block-group + std::vector _ambient_convection_alpha; + /// Name of the solid domain temperature for each block-group + std::vector _ambient_temperature; + + /// Energy inlet boundary types + MultiMooseEnum _energy_inlet_types; + /// Functors describing the inlet boundary values. See energy_inlet_types for what the functors actually represent + std::vector _energy_inlet_functors; + /// Energy wall boundary types + MultiMooseEnum _energy_wall_types; + /// Functors describing the wall boundary values. See energy_wall_types for what the functors actually represent + std::vector _energy_wall_functors; +}; diff --git a/modules/navier_stokes/src/physics/WCNSFVFluidHeatTransferPhysics.C b/modules/navier_stokes/src/physics/WCNSFVFluidHeatTransferPhysics.C index c292fb63de9f..69981cb7aaeb 100644 --- a/modules/navier_stokes/src/physics/WCNSFVFluidHeatTransferPhysics.C +++ b/modules/navier_stokes/src/physics/WCNSFVFluidHeatTransferPhysics.C @@ -14,104 +14,19 @@ #include "NSFVBase.h" registerNavierStokesPhysicsBaseTasks("NavierStokesApp", WCNSFVFluidHeatTransferPhysics); -registerMooseAction("NavierStokesApp", WCNSFVFluidHeatTransferPhysics, "get_turbulence_physics"); -registerMooseAction("NavierStokesApp", WCNSFVFluidHeatTransferPhysics, "add_variable"); -registerMooseAction("NavierStokesApp", WCNSFVFluidHeatTransferPhysics, "add_ic"); -registerMooseAction("NavierStokesApp", WCNSFVFluidHeatTransferPhysics, "add_fv_kernel"); -registerMooseAction("NavierStokesApp", WCNSFVFluidHeatTransferPhysics, "add_fv_bc"); -registerMooseAction("NavierStokesApp", WCNSFVFluidHeatTransferPhysics, "add_material"); +registerWCNSFVFluidHeatTransferPhysicsBaseTasks("NavierStokesApp", WCNSFVFluidHeatTransferPhysics); InputParameters WCNSFVFluidHeatTransferPhysics::validParams() { - InputParameters params = NavierStokesPhysicsBase::validParams(); - params += WCNSFVCoupledAdvectionPhysicsHelper::validParams(); - params.addClassDescription("Define the Navier Stokes weakly-compressible energy equation"); - - params += NSFVBase::commonFluidEnergyEquationParams(); - params.transferParam(PINSFVEnergyAnisotropicDiffusion::validParams(), - "effective_conductivity"); - - // TODO Remove the parameter once NavierStokesFV syntax has been removed - params.addParam("add_energy_equation", - "Whether to add the energy equation. This parameter is not necessary if " - "using the Physics syntax"); - params.addParam( - "fluid_temperature_variable", NS::T_fluid, "Name of the fluid temperature variable"); - - // New functor boundary conditions - params.deprecateParam("energy_inlet_function", "energy_inlet_functors", "01/01/2025"); - params.deprecateParam("energy_wall_function", "energy_wall_functors", "01/01/2025"); - - // Spatial finite volume discretization scheme - params.transferParam(NSFVBase::validParams(), "energy_advection_interpolation"); - params.transferParam(NSFVBase::validParams(), "energy_face_interpolation"); - params.transferParam(NSFVBase::validParams(), "energy_two_term_bc_expansion"); - - // Nonlinear equation solver scaling - params.transferParam(NSFVBase::validParams(), "energy_scaling"); - - params.addParamNamesToGroup("specific_heat thermal_conductivity thermal_conductivity_blocks " - "use_external_enthalpy_material", - "Material properties"); - params.addParamNamesToGroup("energy_advection_interpolation energy_face_interpolation " - "energy_two_term_bc_expansion energy_scaling", - "Numerical scheme"); - params.addParamNamesToGroup("energy_inlet_types energy_inlet_functors", - "Inlet boundary conditions"); - params.addParamNamesToGroup("energy_wall_types energy_wall_functors", "Wall boundary conditions"); - + InputParameters params = WCNSFVFluidHeatTransferPhysicsBase::validParams(); return params; } WCNSFVFluidHeatTransferPhysics::WCNSFVFluidHeatTransferPhysics(const InputParameters & parameters) - : NavierStokesPhysicsBase(parameters), - WCNSFVCoupledAdvectionPhysicsHelper(this), - _has_energy_equation( - isParamValid("add_energy_equation") - ? getParam("add_energy_equation") - : (usingNavierStokesFVSyntax() ? isParamSetByUser("energy_inlet_function") : true)), - _fluid_temperature_name(getParam("fluid_temperature_variable")), - _specific_heat_name(getParam("specific_heat")), - _thermal_conductivity_blocks( - parameters.isParamValid("thermal_conductivity_blocks") - ? getParam>>("thermal_conductivity_blocks") - : std::vector>()), - _thermal_conductivity_name(getParam>("thermal_conductivity")), - _ambient_convection_blocks( - getParam>>("ambient_convection_blocks")), - _ambient_convection_alpha(getParam>("ambient_convection_alpha")), - _ambient_temperature(getParam>("ambient_temperature")), - _energy_inlet_types(getParam("energy_inlet_types")), - _energy_inlet_functors(getParam>("energy_inlet_functors")), - _energy_wall_types(getParam("energy_wall_types")), - _energy_wall_functors(getParam>("energy_wall_functors")) -{ - // For compatibility with Modules/NavierStokesFV syntax - if (!_has_energy_equation) - return; + : WCNSFVFluidHeatTransferPhysicsBase(parameters) - saveSolverVariableName(_fluid_temperature_name); - - // set block restrictions if not set by user - // This should probably be done for all the coupled physics, tbd - if (!isParamSetByUser("block")) - _blocks = _flow_equations_physics->blocks(); - - // Parameter checks - checkVectorParamsSameLengthIfSet("ambient_convection_alpha", - "ambient_temperature"); - checkSecondParamSetOnlyIfFirstOneSet("external_heat_source", "external_heat_source_coeff"); - - // Check boundary parameters if provided. - // The boundaries are checked again when the boundary conditions are added as we want - // to be able to more boundary conditions to a Physics dynamically - if (isParamValid("energy_inlet_types")) - checkVectorParamAndMultiMooseEnumLength("energy_inlet_functors", - "energy_inlet_types"); - if (isParamValid("energy_wall_types")) - checkVectorParamAndMultiMooseEnumLength("energy_wall_functors", - "energy_wall_types"); +{ } void @@ -143,29 +58,6 @@ WCNSFVFluidHeatTransferPhysics::addSolverVariables() ") supplied to the WCNSFVFluidHeatTransferPhysics does not exist!"); } -void -WCNSFVFluidHeatTransferPhysics::addFVKernels() -{ - // For compatibility with Modules/NavierStokesFV syntax - if (!_has_energy_equation) - return; - - if (isTransient()) - { - if (_compressibility == "incompressible") - addINSEnergyTimeKernels(); - else - addWCNSEnergyTimeKernels(); - } - - addINSEnergyAdvectionKernels(); - addINSEnergyHeatConductionKernels(); - if (getParam>("ambient_temperature").size()) - addINSEnergyAmbientConvection(); - if (isParamValid("external_heat_source")) - addINSEnergyExternalHeatSource(); -} - void WCNSFVFluidHeatTransferPhysics::addINSEnergyTimeKernels() { @@ -348,17 +240,6 @@ WCNSFVFluidHeatTransferPhysics::addINSEnergyExternalHeatSource() getProblem().addFVKernel(kernel_type, prefix() + "external_heat_source", params); } -void -WCNSFVFluidHeatTransferPhysics::addFVBCs() -{ - // For compatibility with Modules/NavierStokesFV syntax - if (!_has_energy_equation) - return; - - addINSEnergyInletBC(); - addINSEnergyWallBC(); -} - void WCNSFVFluidHeatTransferPhysics::addINSEnergyInletBC() { @@ -437,14 +318,6 @@ WCNSFVFluidHeatTransferPhysics::addINSEnergyInletBC() } } -void -WCNSFVFluidHeatTransferPhysics::actOnAdditionalTasks() -{ - // Turbulence physics would not be initialized before this task - if (_current_task == "get_turbulence_physics") - _turbulence_physics = getCoupledTurbulencePhysics(); -} - void WCNSFVFluidHeatTransferPhysics::addINSEnergyWallBC() { @@ -521,106 +394,3 @@ WCNSFVFluidHeatTransferPhysics::addINSEnergyWallBC() } } } - -bool -WCNSFVFluidHeatTransferPhysics::processThermalConductivity() -{ - checkBlockwiseConsistency("thermal_conductivity_blocks", - {"thermal_conductivity"}); - bool have_scalar = false; - bool have_vector = false; - - for (unsigned int i = 0; i < _thermal_conductivity_name.size(); ++i) - { - // First, check if the name is just a number (only in case of isotropic conduction) - if (MooseUtils::parsesToReal(_thermal_conductivity_name[i])) - have_scalar = true; - // Now we determine what kind of functor we are dealing with - else - { - if (getProblem().hasFunctorWithType(_thermal_conductivity_name[i], - /*thread_id=*/0)) - have_scalar = true; - else - { - if (getProblem().hasFunctorWithType(_thermal_conductivity_name[i], - /*thread_id=*/0)) - have_vector = true; - else - paramError("thermal_conductivity", - "We only allow functor of type ADReal or ADRealVectorValue for thermal " - "conductivity! Functor '" + - _thermal_conductivity_name[i] + "' is not of the requested type."); - } - } - } - - if (have_vector && !_porous_medium_treatment) - paramError("thermal_conductivity", "Cannot use anisotropic diffusion with non-porous flows!"); - - if (have_vector == have_scalar) - paramError("thermal_conductivity", - "The entries on thermal conductivity shall either be scalars of vectors, mixing " - "them is not supported!"); - return have_vector; -} - -void -WCNSFVFluidHeatTransferPhysics::addInitialConditions() -{ - // For compatibility with Modules/NavierStokesFV syntax - if (!_has_energy_equation) - return; - if (!_define_variables && parameters().isParamSetByUser("initial_temperature")) - paramError( - "initial_temperature", - "T_fluid is defined externally of WCNSFVFluidHeatTransferPhysics, so should the inital " - "condition"); - // do not set initial conditions if we load from file - if (getParam("initialize_variables_from_mesh_file")) - return; - // do not set initial conditions if we are not defining variables - if (!_define_variables) - return; - - InputParameters params = getFactory().getValidParams("FunctionIC"); - assignBlocks(params, _blocks); - - if (!_app.isRestarting() || parameters().isParamSetByUser("initial_temperature")) - { - params.set("variable") = _fluid_temperature_name; - params.set("function") = getParam("initial_temperature"); - - getProblem().addInitialCondition("FunctionIC", _fluid_temperature_name + "_ic", params); - } -} - -void -WCNSFVFluidHeatTransferPhysics::addMaterials() -{ - // For compatibility with Modules/NavierStokesFV syntax - if (!_has_energy_equation) - return; - - InputParameters params = getFactory().getValidParams("INSFVEnthalpyFunctorMaterial"); - assignBlocks(params, _blocks); - - params.set(NS::density) = _density_name; - params.set(NS::cp) = _specific_heat_name; - params.set("temperature") = _fluid_temperature_name; - - getProblem().addMaterial( - "INSFVEnthalpyFunctorMaterial", prefix() + "ins_enthalpy_material", params); -} - -unsigned short -WCNSFVFluidHeatTransferPhysics::getNumberAlgebraicGhostingLayersNeeded() const -{ - unsigned short necessary_layers = getParam("ghost_layers"); - necessary_layers = - std::max(necessary_layers, _flow_equations_physics->getNumberAlgebraicGhostingLayersNeeded()); - if (getParam("energy_face_interpolation") == "skewness-corrected") - necessary_layers = std::max(necessary_layers, (unsigned short)3); - - return necessary_layers; -} diff --git a/modules/navier_stokes/src/physics/WCNSFVFluidHeatTransferPhysicsBase.C b/modules/navier_stokes/src/physics/WCNSFVFluidHeatTransferPhysicsBase.C new file mode 100644 index 000000000000..14582550dc5f --- /dev/null +++ b/modules/navier_stokes/src/physics/WCNSFVFluidHeatTransferPhysicsBase.C @@ -0,0 +1,253 @@ +//* This file is part of the MOOSE framework +//* https://www.mooseframework.org +//* +//* All rights reserved, see COPYRIGHT for full restrictions +//* https://github.com/idaholab/moose/blob/master/COPYRIGHT +//* +//* Licensed under LGPL 2.1, please see LICENSE for details +//* https://www.gnu.org/licenses/lgpl-2.1.html + +#include "WCNSFVFluidHeatTransferPhysicsBase.h" +#include "WCNSFVCoupledAdvectionPhysicsHelper.h" +#include "WCNSFVFlowPhysics.h" +#include "PINSFVEnergyAnisotropicDiffusion.h" +#include "NSFVBase.h" + +InputParameters +WCNSFVFluidHeatTransferPhysicsBase::validParams() +{ + InputParameters params = NavierStokesPhysicsBase::validParams(); + params += WCNSFVCoupledAdvectionPhysicsHelper::validParams(); + params.addClassDescription("Define the Navier Stokes weakly-compressible energy equation"); + + params += NSFVBase::commonFluidEnergyEquationParams(); + params.transferParam(PINSFVEnergyAnisotropicDiffusion::validParams(), + "effective_conductivity"); + + // TODO Remove the parameter once NavierStokesFV syntax has been removed + params.addParam("add_energy_equation", + "Whether to add the energy equation. This parameter is not necessary if " + "using the Physics syntax"); + params.addParam( + "fluid_temperature_variable", NS::T_fluid, "Name of the fluid temperature variable"); + + // New functor boundary conditions + params.deprecateParam("energy_inlet_function", "energy_inlet_functors", "01/01/2025"); + params.deprecateParam("energy_wall_function", "energy_wall_functors", "01/01/2025"); + + // Spatial finite volume discretization scheme + params.transferParam(NSFVBase::validParams(), "energy_advection_interpolation"); + params.transferParam(NSFVBase::validParams(), "energy_face_interpolation"); + params.transferParam(NSFVBase::validParams(), "energy_two_term_bc_expansion"); + + // Nonlinear equation solver scaling + params.transferParam(NSFVBase::validParams(), "energy_scaling"); + + params.addParamNamesToGroup("specific_heat thermal_conductivity thermal_conductivity_blocks " + "use_external_enthalpy_material", + "Material properties"); + params.addParamNamesToGroup("energy_advection_interpolation energy_face_interpolation " + "energy_two_term_bc_expansion energy_scaling", + "Numerical scheme"); + params.addParamNamesToGroup("energy_inlet_types energy_inlet_functors", + "Inlet boundary conditions"); + params.addParamNamesToGroup("energy_wall_types energy_wall_functors", "Wall boundary conditions"); + + return params; +} + +WCNSFVFluidHeatTransferPhysicsBase::WCNSFVFluidHeatTransferPhysicsBase( + const InputParameters & parameters) + : NavierStokesPhysicsBase(parameters), + WCNSFVCoupledAdvectionPhysicsHelper(this), + _has_energy_equation( + isParamValid("add_energy_equation") + ? getParam("add_energy_equation") + : (usingNavierStokesFVSyntax() ? isParamSetByUser("energy_inlet_function") : true)), + _fluid_temperature_name(getParam("fluid_temperature_variable")), + _specific_heat_name(getParam("specific_heat")), + _thermal_conductivity_blocks( + parameters.isParamValid("thermal_conductivity_blocks") + ? getParam>>("thermal_conductivity_blocks") + : std::vector>()), + _thermal_conductivity_name(getParam>("thermal_conductivity")), + _ambient_convection_blocks( + getParam>>("ambient_convection_blocks")), + _ambient_convection_alpha(getParam>("ambient_convection_alpha")), + _ambient_temperature(getParam>("ambient_temperature")), + _energy_inlet_types(getParam("energy_inlet_types")), + _energy_inlet_functors(getParam>("energy_inlet_functors")), + _energy_wall_types(getParam("energy_wall_types")), + _energy_wall_functors(getParam>("energy_wall_functors")) +{ + // For compatibility with Modules/NavierStokesFV syntax + if (!_has_energy_equation) + return; + + saveSolverVariableName(_fluid_temperature_name); + + // set block restrictions if not set by user + // This should probably be done for all the coupled physics, tbd + if (!isParamSetByUser("block")) + _blocks = _flow_equations_physics->blocks(); + + // Parameter checks + checkVectorParamsSameLengthIfSet("ambient_convection_alpha", + "ambient_temperature"); + checkSecondParamSetOnlyIfFirstOneSet("external_heat_source", "external_heat_source_coeff"); + + // Check boundary parameters if provided. + // The boundaries are checked again when the boundary conditions are added as we want + // to be able to more boundary conditions to a Physics dynamically + if (isParamValid("energy_inlet_types")) + checkVectorParamAndMultiMooseEnumLength("energy_inlet_functors", + "energy_inlet_types"); + if (isParamValid("energy_wall_types")) + checkVectorParamAndMultiMooseEnumLength("energy_wall_functors", + "energy_wall_types"); +} + +void +WCNSFVFluidHeatTransferPhysicsBase::addFVKernels() +{ + // For compatibility with Modules/NavierStokesFV syntax + if (!_has_energy_equation) + return; + + if (isTransient()) + { + if (_compressibility == "incompressible") + addINSEnergyTimeKernels(); + else + addWCNSEnergyTimeKernels(); + } + + addINSEnergyAdvectionKernels(); + addINSEnergyHeatConductionKernels(); + if (getParam>("ambient_temperature").size()) + addINSEnergyAmbientConvection(); + if (isParamValid("external_heat_source")) + addINSEnergyExternalHeatSource(); +} + +void +WCNSFVFluidHeatTransferPhysicsBase::addFVBCs() +{ + // For compatibility with Modules/NavierStokesFV syntax + if (!_has_energy_equation) + return; + + addINSEnergyInletBC(); + addINSEnergyWallBC(); +} + +void +WCNSFVFluidHeatTransferPhysicsBase::actOnAdditionalTasks() +{ + // Turbulence physics would not be initialized before this task + if (_current_task == "get_turbulence_physics") + _turbulence_physics = getCoupledTurbulencePhysics(); +} + +bool +WCNSFVFluidHeatTransferPhysicsBase::processThermalConductivity() +{ + checkBlockwiseConsistency("thermal_conductivity_blocks", + {"thermal_conductivity"}); + bool have_scalar = false; + bool have_vector = false; + + for (unsigned int i = 0; i < _thermal_conductivity_name.size(); ++i) + { + // First, check if the name is just a number (only in case of isotropic conduction) + if (MooseUtils::parsesToReal(_thermal_conductivity_name[i])) + have_scalar = true; + // Now we determine what kind of functor we are dealing with + else + { + if (getProblem().hasFunctorWithType(_thermal_conductivity_name[i], + /*thread_id=*/0)) + have_scalar = true; + else + { + if (getProblem().hasFunctorWithType(_thermal_conductivity_name[i], + /*thread_id=*/0)) + have_vector = true; + else + paramError("thermal_conductivity", + "We only allow functor of type ADReal or ADRealVectorValue for thermal " + "conductivity! Functor '" + + _thermal_conductivity_name[i] + "' is not of the requested type."); + } + } + } + + if (have_vector && !_porous_medium_treatment) + paramError("thermal_conductivity", "Cannot use anisotropic diffusion with non-porous flows!"); + + if (have_vector == have_scalar) + paramError("thermal_conductivity", + "The entries on thermal conductivity shall either be scalars of vectors, mixing " + "them is not supported!"); + return have_vector; +} + +void +WCNSFVFluidHeatTransferPhysicsBase::addInitialConditions() +{ + // For compatibility with Modules/NavierStokesFV syntax + if (!_has_energy_equation) + return; + if (!_define_variables && parameters().isParamSetByUser("initial_temperature")) + paramError( + "initial_temperature", + "T_fluid is defined externally of WCNSFVFluidHeatTransferPhysicsBase, so should the inital " + "condition"); + // do not set initial conditions if we load from file + if (getParam("initialize_variables_from_mesh_file")) + return; + // do not set initial conditions if we are not defining variables + if (!_define_variables) + return; + + InputParameters params = getFactory().getValidParams("FunctionIC"); + assignBlocks(params, _blocks); + + if (!_app.isRestarting() || parameters().isParamSetByUser("initial_temperature")) + { + params.set("variable") = _fluid_temperature_name; + params.set("function") = getParam("initial_temperature"); + + getProblem().addInitialCondition("FunctionIC", _fluid_temperature_name + "_ic", params); + } +} + +void +WCNSFVFluidHeatTransferPhysicsBase::addMaterials() +{ + // For compatibility with Modules/NavierStokesFV syntax + if (!_has_energy_equation) + return; + + InputParameters params = getFactory().getValidParams("INSFVEnthalpyFunctorMaterial"); + assignBlocks(params, _blocks); + + params.set(NS::density) = _density_name; + params.set(NS::cp) = _specific_heat_name; + params.set("temperature") = _fluid_temperature_name; + + getProblem().addMaterial( + "INSFVEnthalpyFunctorMaterial", prefix() + "ins_enthalpy_material", params); +} + +unsigned short +WCNSFVFluidHeatTransferPhysicsBase::getNumberAlgebraicGhostingLayersNeeded() const +{ + unsigned short necessary_layers = getParam("ghost_layers"); + necessary_layers = + std::max(necessary_layers, _flow_equations_physics->getNumberAlgebraicGhostingLayersNeeded()); + if (getParam("energy_face_interpolation") == "skewness-corrected") + necessary_layers = std::max(necessary_layers, (unsigned short)3); + + return necessary_layers; +} From 6454223bd203047ce79e4405788e6c4ba99a32dd Mon Sep 17 00:00:00 2001 From: Guillaume Giudicelli Date: Sun, 1 Dec 2024 12:58:28 -0700 Subject: [PATCH 02/16] Add a linear FV physics class for energy transfer refs #29175 --- .../WCNSLinearFVFluidHeatTransferPhysics.h | 46 ++++ .../navier_stokes/src/base/NavierStokesApp.C | 2 + .../physics/WCNSFVFluidHeatTransferPhysics.C | 4 +- .../WCNSFVFluidHeatTransferPhysicsBase.C | 3 +- .../WCNSLinearFVFluidHeatTransferPhysics.C | 260 ++++++++++++++++++ 5 files changed, 311 insertions(+), 4 deletions(-) create mode 100644 modules/navier_stokes/include/physics/WCNSLinearFVFluidHeatTransferPhysics.h create mode 100644 modules/navier_stokes/src/physics/WCNSLinearFVFluidHeatTransferPhysics.C diff --git a/modules/navier_stokes/include/physics/WCNSLinearFVFluidHeatTransferPhysics.h b/modules/navier_stokes/include/physics/WCNSLinearFVFluidHeatTransferPhysics.h new file mode 100644 index 000000000000..834e1e6c2f42 --- /dev/null +++ b/modules/navier_stokes/include/physics/WCNSLinearFVFluidHeatTransferPhysics.h @@ -0,0 +1,46 @@ +//* This file is part of the MOOSE framework +//* https://www.mooseframework.org +//* +//* All rights reserved, see COPYRIGHT for full restrictions +//* https://github.com/idaholab/moose/blob/master/COPYRIGHT +//* +//* Licensed under LGPL 2.1, please see LICENSE for details +//* https://www.gnu.org/licenses/lgpl-2.1.html + +#pragma once + +#include "WCNSFVFluidHeatTransferPhysicsBase.h" +#include "NS.h" + +/** + * Creates all the objects needed to solve the Navier Stokes energy equation using a linear finite + * volume discretization + */ +class WCNSLinearFVFluidHeatTransferPhysics final : public WCNSFVFluidHeatTransferPhysicsBase +{ +public: + static InputParameters validParams(); + + WCNSLinearFVFluidHeatTransferPhysics(const InputParameters & parameters); + +protected: +private: + virtual void addSolverVariables() override; + + /** + * Functions adding kernels for the incompressible / weakly compressible energy equation + * If the material properties are not constant, some of these can be used for + * weakly-compressible simulations as well. + */ + void addINSEnergyTimeKernels() override; + void addWCNSEnergyTimeKernels() override; + void addINSEnergyHeatConductionKernels() override; + void addINSEnergyAdvectionKernels() override; + void addINSEnergyAmbientConvection() override; + void addINSEnergyExternalHeatSource() override; + + /// Functions adding boundary conditions for the incompressible simulation. + /// These are used for weakly-compressible simulations as well. + void addINSEnergyInletBC() override; + void addINSEnergyWallBC() override; +}; diff --git a/modules/navier_stokes/src/base/NavierStokesApp.C b/modules/navier_stokes/src/base/NavierStokesApp.C index 1b6c7ad20d00..9542cc833591 100644 --- a/modules/navier_stokes/src/base/NavierStokesApp.C +++ b/modules/navier_stokes/src/base/NavierStokesApp.C @@ -53,6 +53,8 @@ associateSyntaxInner(Syntax & syntax, ActionFactory & /*action_factory*/) registerSyntax("WCNSFVFlowPhysics", "Physics/NavierStokes/Flow/*"); registerSyntax("WCNSLinearFVFlowPhysics", "Physics/NavierStokes/FlowSegregated/*"); registerSyntax("WCNSFVFluidHeatTransferPhysics", "Physics/NavierStokes/FluidHeatTransfer/*"); + registerSyntax("WCNSLinearFVFluidHeatTransferPhysics", + "Physics/NavierStokes/FluidHeatTransferSegregated/*"); registerSyntax("WCNSFVScalarTransportPhysics", "Physics/NavierStokes/ScalarTransport/*"); registerSyntax("WCNSLinearFVScalarTransportPhysics", "Physics/NavierStokes/ScalarTransportSegregated/*"); diff --git a/modules/navier_stokes/src/physics/WCNSFVFluidHeatTransferPhysics.C b/modules/navier_stokes/src/physics/WCNSFVFluidHeatTransferPhysics.C index 69981cb7aaeb..491329e11570 100644 --- a/modules/navier_stokes/src/physics/WCNSFVFluidHeatTransferPhysics.C +++ b/modules/navier_stokes/src/physics/WCNSFVFluidHeatTransferPhysics.C @@ -8,9 +8,7 @@ //* https://www.gnu.org/licenses/lgpl-2.1.html #include "WCNSFVFluidHeatTransferPhysics.h" -#include "WCNSFVCoupledAdvectionPhysicsHelper.h" #include "WCNSFVFlowPhysics.h" -#include "PINSFVEnergyAnisotropicDiffusion.h" #include "NSFVBase.h" registerNavierStokesPhysicsBaseTasks("NavierStokesApp", WCNSFVFluidHeatTransferPhysics); @@ -20,6 +18,8 @@ InputParameters WCNSFVFluidHeatTransferPhysics::validParams() { InputParameters params = WCNSFVFluidHeatTransferPhysicsBase::validParams(); + params.transferParam(NSFVBase::validParams(), "energy_face_interpolation"); + params.addParamNamesToGroup("energy_face_interpolation", "Numerical scheme"); return params; } diff --git a/modules/navier_stokes/src/physics/WCNSFVFluidHeatTransferPhysicsBase.C b/modules/navier_stokes/src/physics/WCNSFVFluidHeatTransferPhysicsBase.C index 14582550dc5f..bffce2fb646a 100644 --- a/modules/navier_stokes/src/physics/WCNSFVFluidHeatTransferPhysicsBase.C +++ b/modules/navier_stokes/src/physics/WCNSFVFluidHeatTransferPhysicsBase.C @@ -37,7 +37,6 @@ WCNSFVFluidHeatTransferPhysicsBase::validParams() // Spatial finite volume discretization scheme params.transferParam(NSFVBase::validParams(), "energy_advection_interpolation"); - params.transferParam(NSFVBase::validParams(), "energy_face_interpolation"); params.transferParam(NSFVBase::validParams(), "energy_two_term_bc_expansion"); // Nonlinear equation solver scaling @@ -46,7 +45,7 @@ WCNSFVFluidHeatTransferPhysicsBase::validParams() params.addParamNamesToGroup("specific_heat thermal_conductivity thermal_conductivity_blocks " "use_external_enthalpy_material", "Material properties"); - params.addParamNamesToGroup("energy_advection_interpolation energy_face_interpolation " + params.addParamNamesToGroup("energy_advection_interpolation " "energy_two_term_bc_expansion energy_scaling", "Numerical scheme"); params.addParamNamesToGroup("energy_inlet_types energy_inlet_functors", diff --git a/modules/navier_stokes/src/physics/WCNSLinearFVFluidHeatTransferPhysics.C b/modules/navier_stokes/src/physics/WCNSLinearFVFluidHeatTransferPhysics.C new file mode 100644 index 000000000000..7b254bfc8099 --- /dev/null +++ b/modules/navier_stokes/src/physics/WCNSLinearFVFluidHeatTransferPhysics.C @@ -0,0 +1,260 @@ +//* This file is part of the MOOSE framework +//* https://www.mooseframework.org +//* +//* All rights reserved, see COPYRIGHT for full restrictions +//* https://github.com/idaholab/moose/blob/master/COPYRIGHT +//* +//* Licensed under LGPL 2.1, please see LICENSE for details +//* https://www.gnu.org/licenses/lgpl-2.1.html + +#include "WCNSLinearFVFluidHeatTransferPhysics.h" +#include "WCNSFVFlowPhysics.h" +#include "PINSFVEnergyAnisotropicDiffusion.h" +#include "NSFVBase.h" + +registerNavierStokesPhysicsBaseTasks("NavierStokesApp", WCNSLinearFVFluidHeatTransferPhysics); +registerWCNSFVFluidHeatTransferPhysicsBaseTasks("NavierStokesApp", + WCNSLinearFVFluidHeatTransferPhysics); + +InputParameters +WCNSLinearFVFluidHeatTransferPhysics::validParams() +{ + InputParameters params = WCNSFVFluidHeatTransferPhysicsBase::validParams(); + params.addParam( + "use_nonorthogonal_correction", + true, + "Whether to use a non-orthogonal correction. This will add off-diagonal terms potentially " + "slowing down convergence, but reducing numerical dispersion"); + params.set>("system_names") = {"energy_system"}; + + // Not implemented + params.suppressParameter("effective_conductivity"); + return params; +} + +WCNSLinearFVFluidHeatTransferPhysics::WCNSLinearFVFluidHeatTransferPhysics( + const InputParameters & parameters) + : WCNSFVFluidHeatTransferPhysicsBase(parameters) +{ + if (_porous_medium_treatment) + paramError("porous_medium_treatment", "Porous media not supported at this time"); +} + +void +WCNSLinearFVFluidHeatTransferPhysics::addSolverVariables() +{ + // For compatibility with Modules/NavierStokesFV syntax + if (!_has_energy_equation) + return; + + // Dont add if the user already defined the variable + if (variableExists(_fluid_temperature_name, + /*error_if_aux=*/true)) + checkBlockRestrictionIdentical(_fluid_temperature_name, + getProblem().getVariable(0, _fluid_temperature_name).blocks()); + else if (_define_variables) + { + const auto var_type = "MooseLinearVariableFVReal"; + auto params = getFactory().getValidParams(var_type); + assignBlocks(params, _blocks); + params.set>("scaling") = {getParam("energy_scaling")}; + params.set("solver_sys") = getSolverSystem(_fluid_temperature_name); + getProblem().addVariable(var_type, _fluid_temperature_name, params); + } + else + paramError("fluid_temperature_variable", + "Variable (" + _fluid_temperature_name + + ") supplied to the WCNSLinearFVFluidHeatTransferPhysics does not exist!"); +} + +void +WCNSLinearFVFluidHeatTransferPhysics::addINSEnergyTimeKernels() +{ + paramError("transient", + "Transient simulations not supported at this time with the linear FV discretization"); +} + +void +WCNSLinearFVFluidHeatTransferPhysics::addWCNSEnergyTimeKernels() +{ + paramError("transient", + "Transient simulations not supported at this time with the linear FV discretization"); +} + +void +WCNSLinearFVFluidHeatTransferPhysics::addINSEnergyAdvectionKernels() +{ + std::string kernel_type = "LinearFVEnergyAdvection"; + std::string kernel_name = prefix() + "ins_energy_advection"; + + InputParameters params = getFactory().getValidParams(kernel_type); + params.set("variable") = _fluid_temperature_name; + assignBlocks(params, _blocks); + params.set("rhie_chow_user_object") = _flow_equations_physics->rhieChowUOName(); + params.set("advected_interp_method") = + getParam("energy_advection_interpolation"); + + getProblem().addFVKernel(kernel_type, kernel_name, params); +} + +void +WCNSLinearFVFluidHeatTransferPhysics::addINSEnergyHeatConductionKernels() +{ + const auto vector_conductivity = processThermalConductivity(); + const auto num_blocks = _thermal_conductivity_blocks.size(); + const auto num_used_blocks = num_blocks ? num_blocks : 1; + + for (const auto block_i : make_range(num_used_blocks)) + { + std::string block_name = ""; + if (num_blocks) + block_name = Moose::stringify(_thermal_conductivity_blocks[block_i]); + else + block_name = "all"; + + const std::string kernel_type = "LinearFVDiffusion"; + InputParameters params = getFactory().getValidParams(kernel_type); + params.set("variable") = _fluid_temperature_name; + std::vector block_names = + num_blocks ? _thermal_conductivity_blocks[block_i] : _blocks; + assignBlocks(params, block_names); + params.set("diffusion_coeff") = _thermal_conductivity_name[block_i]; + params.set("use_nonorthogonal_correction") = + getParam("use_nonorthogonal_correction"); + + getProblem().addFVKernel(kernel_type, prefix() + "ins_energy_diffusion_" + block_name, params); + } +} + +void +WCNSLinearFVFluidHeatTransferPhysics::addINSEnergyAmbientConvection() +{ + unsigned int num_convection_blocks = _ambient_convection_blocks.size(); + unsigned int num_used_blocks = num_convection_blocks ? num_convection_blocks : 1; + + const std::string kernel_type = "LinearFVVolumetricHeatTransfer"; + InputParameters params = getFactory().getValidParams(kernel_type); + params.set("variable") = _fluid_temperature_name; + params.set(NS::T_fluid) = _fluid_temperature_name; + params.set("is_solid") = false; + + for (unsigned int block_i = 0; block_i < num_used_blocks; ++block_i) + { + std::string block_name = ""; + if (num_convection_blocks) + { + params.set>("block") = _ambient_convection_blocks[block_i]; + block_name = Moose::stringify(_ambient_convection_blocks[block_i]); + } + else + { + assignBlocks(params, _blocks); + block_name = std::to_string(block_i); + } + + params.set("h_solid_fluid") = _ambient_convection_alpha[block_i]; + params.set(NS::T_solid) = _ambient_temperature[block_i]; + + getProblem().addFVKernel(kernel_type, prefix() + "ambient_convection_" + block_name, params); + } +} + +void +WCNSLinearFVFluidHeatTransferPhysics::addINSEnergyExternalHeatSource() +{ + const std::string kernel_type = "LinearFVSource"; + InputParameters params = getFactory().getValidParams(kernel_type); + params.set("variable") = _fluid_temperature_name; + assignBlocks(params, _blocks); + params.set("source_density") = + getParam("external_heat_source"); + params.set("scaling_factor") = getParam("external_heat_source_coeff"); + + getProblem().addFVKernel(kernel_type, prefix() + "external_heat_source", params); +} + +void +WCNSLinearFVFluidHeatTransferPhysics::addINSEnergyInletBC() +{ + const auto & inlet_boundaries = _flow_equations_physics->getInletBoundaries(); + // These are parameter errors for now. If Components add boundaries to Physics, the error + // may not be due to parameters anymore. + if (inlet_boundaries.size() != _energy_inlet_types.size()) + paramError("energy_inlet_types", + "Energy inlet types (size " + std::to_string(_energy_inlet_types.size()) + + ") should be the same size as inlet_boundaries (size " + + std::to_string(inlet_boundaries.size()) + ")"); + if (inlet_boundaries.size() != _energy_inlet_functors.size()) + paramError("energy_inlet_functors", + "Energy inlet functors (size " + std::to_string(_energy_inlet_functors.size()) + + ") should be the same size as inlet_boundaries (size " + + std::to_string(inlet_boundaries.size()) + ")"); + + unsigned int flux_bc_counter = 0; + for (const auto bc_ind : index_range(_energy_inlet_types)) + { + if (_energy_inlet_types[bc_ind] == "fixed-temperature") + { + const std::string bc_type = "LinearFVAdvectionDiffusionFunctorDirichletBC"; + InputParameters params = getFactory().getValidParams(bc_type); + params.set("variable") = _fluid_temperature_name; + params.set("functor") = _energy_inlet_functors[bc_ind]; + params.set>("boundary") = {inlet_boundaries[bc_ind]}; + + getProblem().addFVBC( + bc_type, _fluid_temperature_name + "_" + inlet_boundaries[bc_ind], params); + } + else if (_energy_inlet_types[bc_ind] == "heatflux") + { + paramError("energy_inlet_types", "Heat flux inlet boundary conditions not yet supported"); + } + else if (_energy_inlet_types[bc_ind] == "flux-mass" || + _energy_inlet_types[bc_ind] == "flux-velocity") + paramError("energy_inlet_types", "Flux inlet boundary conditions not yet supported"); + } +} + +void +WCNSLinearFVFluidHeatTransferPhysics::addINSEnergyWallBC() +{ + const auto & wall_boundaries = _flow_equations_physics->getWallBoundaries(); + if (wall_boundaries.size() != _energy_wall_types.size()) + paramError("energy_wall_types", + "Energy wall types (size " + std::to_string(_energy_wall_types.size()) + + ") should be the same size as wall_boundaries (size " + + std::to_string(wall_boundaries.size()) + ")"); + if (wall_boundaries.size() != _energy_wall_functors.size()) + paramError("energy_wall_functors", + "Energy wall functors (size " + std::to_string(_energy_wall_functors.size()) + + ") should be the same size as wall_boundaries (size " + + std::to_string(wall_boundaries.size()) + ")"); + + for (unsigned int bc_ind = 0; bc_ind < _energy_wall_types.size(); ++bc_ind) + { + if (_energy_wall_types[bc_ind] == "fixed-temperature") + { + const std::string bc_type = "LinearFVAdvectionDiffusionFunctorDirichletBC"; + InputParameters params = getFactory().getValidParams(bc_type); + params.set("variable") = _fluid_temperature_name; + params.set("functor") = _energy_wall_functors[bc_ind]; + params.set>("boundary") = {wall_boundaries[bc_ind]}; + + getProblem().addFVBC( + bc_type, _fluid_temperature_name + "_" + wall_boundaries[bc_ind], params); + } + else if (_energy_wall_types[bc_ind] == "heatflux") + { + paramError("energy_inlet_types", "Heat flux wall boundary conditions not yet supported"); + } + // We add this boundary condition here to facilitate the input of wall boundaries / functors for + // energy. If there are too many turbulence options and this gets out of hand we will have to + // move this to the turbulence Physics + else if (_energy_wall_types[bc_ind] == "wallfunction") + mooseError("Wall-function boundary condition not supported at this time."); + } +} + +void +WCNSLinearFVFluidHeatTransferPhysics::addINSEnergyOutletBC() +{ +} From 79ecaa187b9922e3d15292a5b61a35429a043c30 Mon Sep 17 00:00:00 2001 From: Guillaume Giudicelli Date: Sun, 1 Dec 2024 17:41:39 -0700 Subject: [PATCH 03/16] Add testing for the linear energy transfer physics refs #29175 --- .../WCNSLinearFVFluidHeatTransferPhysics.C | 2 +- .../finite_volume/ins/action/errors/tests | 2 +- .../2d-heated/fluid-physics.i | 150 ++++++++++++++++++ .../2d-heated/gold/fluid-physics_out.e | 1 + .../2d-heated/gold/fluid-physics_out_solid0.e | 1 + .../linear-segregated/2d-heated/tests | 9 ++ 6 files changed, 163 insertions(+), 2 deletions(-) create mode 100644 modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/linear-segregated/2d-heated/fluid-physics.i create mode 120000 modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/linear-segregated/2d-heated/gold/fluid-physics_out.e create mode 120000 modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/linear-segregated/2d-heated/gold/fluid-physics_out_solid0.e diff --git a/modules/navier_stokes/src/physics/WCNSLinearFVFluidHeatTransferPhysics.C b/modules/navier_stokes/src/physics/WCNSLinearFVFluidHeatTransferPhysics.C index 7b254bfc8099..cc5c301c3254 100644 --- a/modules/navier_stokes/src/physics/WCNSLinearFVFluidHeatTransferPhysics.C +++ b/modules/navier_stokes/src/physics/WCNSLinearFVFluidHeatTransferPhysics.C @@ -100,7 +100,6 @@ WCNSLinearFVFluidHeatTransferPhysics::addINSEnergyAdvectionKernels() void WCNSLinearFVFluidHeatTransferPhysics::addINSEnergyHeatConductionKernels() { - const auto vector_conductivity = processThermalConductivity(); const auto num_blocks = _thermal_conductivity_blocks.size(); const auto num_used_blocks = num_blocks ? num_blocks : 1; @@ -118,6 +117,7 @@ WCNSLinearFVFluidHeatTransferPhysics::addINSEnergyHeatConductionKernels() std::vector block_names = num_blocks ? _thermal_conductivity_blocks[block_i] : _blocks; assignBlocks(params, block_names); + // NOTE: vector conductivities not supported at this time params.set("diffusion_coeff") = _thermal_conductivity_name[block_i]; params.set("use_nonorthogonal_correction") = getParam("use_nonorthogonal_correction"); diff --git a/modules/navier_stokes/test/tests/finite_volume/ins/action/errors/tests b/modules/navier_stokes/test/tests/finite_volume/ins/action/errors/tests index 21a9e3ea69c4..2bf317b64cdc 100644 --- a/modules/navier_stokes/test/tests/finite_volume/ins/action/errors/tests +++ b/modules/navier_stokes/test/tests/finite_volume/ins/action/errors/tests @@ -359,7 +359,7 @@ "Modules/NavierStokesFV/define_variables=false Modules/velocity_variable='vel_x " "vel_y' Modules/pressure_variable='pressure'" expect_err = "(T_fluid is defined externally of NavierStokesFV, so should the inital " - "condition|T_fluid is defined externally of WCNSFVFluidHeatTransferPhysics, so " + "condition|T_fluid is defined externally of WCNSFVFluidHeatTransferPhysics(Base|), so " "should the inital condition)" [] [] diff --git a/modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/linear-segregated/2d-heated/fluid-physics.i b/modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/linear-segregated/2d-heated/fluid-physics.i new file mode 100644 index 000000000000..012436fc879a --- /dev/null +++ b/modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/linear-segregated/2d-heated/fluid-physics.i @@ -0,0 +1,150 @@ +mu = 2.6 +rho = 1.0 +advected_interp_method = 'upwind' +cp = 1000 +k = 1 +h_cv = 100 + +[Mesh] + [mesh] + type = CartesianMeshGenerator + dim = 2 + dx = '0.25 0.25' + dy = '0.2' + ix = '5 5' + iy = '5' + subdomain_id = '0 1' + [] +[] + +[Problem] + linear_sys_names = 'u_system v_system pressure_system energy_system' + previous_nl_solution_required = true +[] + +[Physics] + [NavierStokes] + [FlowSegregated/flow] + velocity_variable = 'vel_x vel_y' + pressure_variable = 'pressure' + + initial_velocity = '0.5 0 0' + initial_pressure = '0.2' + + density = ${rho} + dynamic_viscosity = ${mu} + + inlet_boundaries = 'left' + momentum_inlet_types = 'fixed-velocity' + momentum_inlet_functors = '1.1 0' + + wall_boundaries = 'top bottom' + momentum_wall_types = 'noslip noslip' + + outlet_boundaries = 'right' + momentum_outlet_types = 'fixed-pressure' + pressure_functors = '1.4' + + orthogonality_correction = false + pressure_two_term_bc_expansion = false + momentum_two_term_bc_expansion = false + momentum_advection_interpolation = ${advected_interp_method} + [] + [FluidHeatTransferSegregated/energy] + thermal_conductivity = ${k} + specific_heat = ${cp} + + initial_temperature = 300 + + # Left is an inlet + energy_inlet_types = 'fixed-temperature' + energy_inlet_functors = '300' + + # Top and bottom are fixed temperature 'inlets' in the reference kernel-based input + energy_wall_types = 'fixed-temperature fixed-temperature' + energy_wall_functors = '300 300' + + external_heat_source = 0 + + ambient_convection_alpha = ${h_cv} + ambient_temperature = 'T_solid' + ambient_convection_blocks = '1' + + energy_advection_interpolation = ${advected_interp_method} + energy_two_term_bc_expansion = false + use_nonorthogonal_correction = false + [] + [] +[] + +[AuxVariables] + [T_solid] + type = MooseLinearVariableFVReal + initial_condition = 300 + block = 1 + [] +[] + +[MultiApps] + inactive = 'solid' + [solid] + type = FullSolveMultiApp + input_files = solid.i + execute_on = timestep_begin + no_restore = true + [] +[] + +[Transfers] + inactive = 'from_solid to_solid' + [from_solid] + type = MultiAppGeneralFieldShapeEvaluationTransfer + from_multi_app = solid + source_variable = 'T_solid' + variable = 'T_solid' + execute_on = timestep_begin + from_blocks = 1 + [] + [to_solid] + type = MultiAppGeneralFieldShapeEvaluationTransfer + to_multi_app = solid + source_variable = 'T_fluid' + variable = 'T_fluid' + execute_on = timestep_begin + to_blocks = 1 + [] +[] + +[Executioner] + type = SIMPLE + momentum_l_abs_tol = 1e-13 + pressure_l_abs_tol = 1e-13 + energy_l_abs_tol = 1e-13 + momentum_l_tol = 0 + pressure_l_tol = 0 + energy_l_tol = 0 + rhie_chow_user_object = 'ins_rhie_chow_interpolator' + momentum_systems = 'u_system v_system' + pressure_system = 'pressure_system' + energy_system = 'energy_system' + momentum_equation_relaxation = 0.8 + energy_equation_relaxation = 0.9 + pressure_variable_relaxation = 0.3 + num_iterations = 20 + pressure_absolute_tolerance = 1e-10 + momentum_absolute_tolerance = 1e-10 + energy_absolute_tolerance = 1e-10 + momentum_petsc_options_iname = '-pc_type -pc_hypre_type' + momentum_petsc_options_value = 'hypre boomeramg' + pressure_petsc_options_iname = '-pc_type -pc_hypre_type' + pressure_petsc_options_value = 'hypre boomeramg' + energy_petsc_options_iname = '-pc_type -pc_hypre_type' + energy_petsc_options_value = 'hypre boomeramg' + print_fields = false + continue_on_max_its = true +[] + +[Outputs] + exodus = true + execute_on = timestep_end +[] diff --git a/modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/linear-segregated/2d-heated/gold/fluid-physics_out.e b/modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/linear-segregated/2d-heated/gold/fluid-physics_out.e new file mode 120000 index 000000000000..9b3925cf922b --- /dev/null +++ b/modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/linear-segregated/2d-heated/gold/fluid-physics_out.e @@ -0,0 +1 @@ +fluid_out.e \ No newline at end of file diff --git a/modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/linear-segregated/2d-heated/gold/fluid-physics_out_solid0.e b/modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/linear-segregated/2d-heated/gold/fluid-physics_out_solid0.e new file mode 120000 index 000000000000..c45f3a1eea2a --- /dev/null +++ b/modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/linear-segregated/2d-heated/gold/fluid-physics_out_solid0.e @@ -0,0 +1 @@ +fluid_out_solid0.e \ No newline at end of file diff --git a/modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/linear-segregated/2d-heated/tests b/modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/linear-segregated/2d-heated/tests index 098bc4b2fd10..f14a44afaa37 100644 --- a/modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/linear-segregated/2d-heated/tests +++ b/modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/linear-segregated/2d-heated/tests @@ -21,5 +21,14 @@ detail = "with the solid temperature problem as the main application." cli_args = "Executioner/fixed_point_max_its=10 Executioner/fixed_point_rel_tol=1e-12 MultiApps/inactive='' Transfers/inactive=''" [] + [simple-main-physics] + type = 'Exodiff' + input = fluid-physics.i + exodiff = 'fluid-physics_out.e fluid-physics_out_solid0.e' + recover = false # we don't support recovery for SIMPLE yet + max_threads = 1 # see libmesh issue #3808 + detail = "with the fluid flow problem as the main application using a Physics shorthand syntax, " + cli_args = "Executioner/fixed_point_max_its=10 Executioner/fixed_point_rel_tol=1e-12 MultiApps/inactive='' Transfers/inactive=''" + [] [] [] From eb1a96679c3b79a03cf26ea2796737045f858fe6 Mon Sep 17 00:00:00 2001 From: Guillaume Giudicelli Date: Sun, 1 Dec 2024 18:04:55 -0700 Subject: [PATCH 04/16] Add support for creating outlet BCs for fluid energy physics --- .../physics/WCNSFVFluidHeatTransferPhysics.h | 1 + .../WCNSFVFluidHeatTransferPhysicsBase.h | 1 + .../WCNSLinearFVFluidHeatTransferPhysics.h | 1 + .../WCNSFVFluidHeatTransferPhysicsBase.C | 1 + .../WCNSLinearFVFluidHeatTransferPhysics.C | 32 ++++++++++++++----- 5 files changed, 28 insertions(+), 8 deletions(-) diff --git a/modules/navier_stokes/include/physics/WCNSFVFluidHeatTransferPhysics.h b/modules/navier_stokes/include/physics/WCNSFVFluidHeatTransferPhysics.h index 3e4143eb3a32..c8f6733d9c38 100644 --- a/modules/navier_stokes/include/physics/WCNSFVFluidHeatTransferPhysics.h +++ b/modules/navier_stokes/include/physics/WCNSFVFluidHeatTransferPhysics.h @@ -42,4 +42,5 @@ class WCNSFVFluidHeatTransferPhysics final : public WCNSFVFluidHeatTransferPhysi /// These are used for weakly-compressible simulations as well. void addINSEnergyInletBC() override; void addINSEnergyWallBC() override; + void addINSEnergyOutletBC() override {} }; diff --git a/modules/navier_stokes/include/physics/WCNSFVFluidHeatTransferPhysicsBase.h b/modules/navier_stokes/include/physics/WCNSFVFluidHeatTransferPhysicsBase.h index 2465dc6cb392..35805fe7559b 100644 --- a/modules/navier_stokes/include/physics/WCNSFVFluidHeatTransferPhysicsBase.h +++ b/modules/navier_stokes/include/physics/WCNSFVFluidHeatTransferPhysicsBase.h @@ -82,6 +82,7 @@ class WCNSFVFluidHeatTransferPhysicsBase : public NavierStokesPhysicsBase, /// These are used for weakly-compressible simulations as well. virtual void addINSEnergyInletBC() = 0; virtual void addINSEnergyWallBC() = 0; + virtual void addINSEnergyOutletBC() = 0; /// Process thermal conductivity (multiple functor input options are available). /// Return true if we have vector thermal conductivity and false if scalar diff --git a/modules/navier_stokes/include/physics/WCNSLinearFVFluidHeatTransferPhysics.h b/modules/navier_stokes/include/physics/WCNSLinearFVFluidHeatTransferPhysics.h index 834e1e6c2f42..a90889d98548 100644 --- a/modules/navier_stokes/include/physics/WCNSLinearFVFluidHeatTransferPhysics.h +++ b/modules/navier_stokes/include/physics/WCNSLinearFVFluidHeatTransferPhysics.h @@ -43,4 +43,5 @@ class WCNSLinearFVFluidHeatTransferPhysics final : public WCNSFVFluidHeatTransfe /// These are used for weakly-compressible simulations as well. void addINSEnergyInletBC() override; void addINSEnergyWallBC() override; + void addINSEnergyOutletBC() override; }; diff --git a/modules/navier_stokes/src/physics/WCNSFVFluidHeatTransferPhysicsBase.C b/modules/navier_stokes/src/physics/WCNSFVFluidHeatTransferPhysicsBase.C index bffce2fb646a..a2ea67a86c22 100644 --- a/modules/navier_stokes/src/physics/WCNSFVFluidHeatTransferPhysicsBase.C +++ b/modules/navier_stokes/src/physics/WCNSFVFluidHeatTransferPhysicsBase.C @@ -138,6 +138,7 @@ WCNSFVFluidHeatTransferPhysicsBase::addFVBCs() addINSEnergyInletBC(); addINSEnergyWallBC(); + addINSEnergyOutletBC(); } void diff --git a/modules/navier_stokes/src/physics/WCNSLinearFVFluidHeatTransferPhysics.C b/modules/navier_stokes/src/physics/WCNSLinearFVFluidHeatTransferPhysics.C index cc5c301c3254..bd349a718fd3 100644 --- a/modules/navier_stokes/src/physics/WCNSLinearFVFluidHeatTransferPhysics.C +++ b/modules/navier_stokes/src/physics/WCNSLinearFVFluidHeatTransferPhysics.C @@ -24,7 +24,8 @@ WCNSLinearFVFluidHeatTransferPhysics::validParams() "use_nonorthogonal_correction", true, "Whether to use a non-orthogonal correction. This will add off-diagonal terms potentially " - "slowing down convergence, but reducing numerical dispersion"); + "slowing down convergence, but reducing numerical dispersion on non-orthogonal meshes. Can " + "be safely turned off on orthogonal meshes."); params.set>("system_names") = {"energy_system"}; // Not implemented @@ -94,7 +95,7 @@ WCNSLinearFVFluidHeatTransferPhysics::addINSEnergyAdvectionKernels() params.set("advected_interp_method") = getParam("energy_advection_interpolation"); - getProblem().addFVKernel(kernel_type, kernel_name, params); + getProblem().addLinearFVKernel(kernel_type, kernel_name, params); } void @@ -122,7 +123,8 @@ WCNSLinearFVFluidHeatTransferPhysics::addINSEnergyHeatConductionKernels() params.set("use_nonorthogonal_correction") = getParam("use_nonorthogonal_correction"); - getProblem().addFVKernel(kernel_type, prefix() + "ins_energy_diffusion_" + block_name, params); + getProblem().addLinearFVKernel( + kernel_type, prefix() + "ins_energy_diffusion_" + block_name, params); } } @@ -155,7 +157,8 @@ WCNSLinearFVFluidHeatTransferPhysics::addINSEnergyAmbientConvection() params.set("h_solid_fluid") = _ambient_convection_alpha[block_i]; params.set(NS::T_solid) = _ambient_temperature[block_i]; - getProblem().addFVKernel(kernel_type, prefix() + "ambient_convection_" + block_name, params); + getProblem().addLinearFVKernel( + kernel_type, prefix() + "ambient_convection_" + block_name, params); } } @@ -170,7 +173,7 @@ WCNSLinearFVFluidHeatTransferPhysics::addINSEnergyExternalHeatSource() getParam("external_heat_source"); params.set("scaling_factor") = getParam("external_heat_source_coeff"); - getProblem().addFVKernel(kernel_type, prefix() + "external_heat_source", params); + getProblem().addLinearFVKernel(kernel_type, prefix() + "external_heat_source", params); } void @@ -190,7 +193,6 @@ WCNSLinearFVFluidHeatTransferPhysics::addINSEnergyInletBC() ") should be the same size as inlet_boundaries (size " + std::to_string(inlet_boundaries.size()) + ")"); - unsigned int flux_bc_counter = 0; for (const auto bc_ind : index_range(_energy_inlet_types)) { if (_energy_inlet_types[bc_ind] == "fixed-temperature") @@ -201,7 +203,7 @@ WCNSLinearFVFluidHeatTransferPhysics::addINSEnergyInletBC() params.set("functor") = _energy_inlet_functors[bc_ind]; params.set>("boundary") = {inlet_boundaries[bc_ind]}; - getProblem().addFVBC( + getProblem().addLinearFVBC( bc_type, _fluid_temperature_name + "_" + inlet_boundaries[bc_ind], params); } else if (_energy_inlet_types[bc_ind] == "heatflux") @@ -239,7 +241,7 @@ WCNSLinearFVFluidHeatTransferPhysics::addINSEnergyWallBC() params.set("functor") = _energy_wall_functors[bc_ind]; params.set>("boundary") = {wall_boundaries[bc_ind]}; - getProblem().addFVBC( + getProblem().addLinearFVBC( bc_type, _fluid_temperature_name + "_" + wall_boundaries[bc_ind], params); } else if (_energy_wall_types[bc_ind] == "heatflux") @@ -257,4 +259,18 @@ WCNSLinearFVFluidHeatTransferPhysics::addINSEnergyWallBC() void WCNSLinearFVFluidHeatTransferPhysics::addINSEnergyOutletBC() { + const auto & outlet_boundaries = _flow_equations_physics->getOutletBoundaries(); + if (outlet_boundaries.empty()) + return; + + for (const auto & outlet_bdy : outlet_boundaries) + { + const std::string bc_type = "LinearFVAdvectionDiffusionOutflowBC"; + InputParameters params = getFactory().getValidParams(bc_type); + params.set>("boundary") = {outlet_bdy}; + params.set("use_two_term_expansion") = getParam("energy_two_term_bc_expansion"); + + params.set("variable") = _fluid_temperature_name; + getProblem().addLinearFVBC(bc_type, _fluid_temperature_name + "_" + outlet_bdy, params); + } } From 2661a2267218bd5d8dda41a12aed6f705e4fe851 Mon Sep 17 00:00:00 2001 From: Guillaume Giudicelli Date: Sun, 1 Dec 2024 18:34:21 -0700 Subject: [PATCH 05/16] Add documentation for new physics --- .../WCNSLinearFVFluidHeatTransferPhysics.md | 59 +++++++++++++++++++ .../NavierStokes/FluidHeatTransfer/index.md | 2 +- .../FluidHeatTransferSegregated/index.md | 11 ++++ .../WCNSFVFluidHeatTransferPhysicsBase.C | 3 +- 4 files changed, 73 insertions(+), 2 deletions(-) create mode 100644 modules/navier_stokes/doc/content/source/physics/WCNSLinearFVFluidHeatTransferPhysics.md create mode 100644 modules/navier_stokes/doc/content/syntax/Physics/NavierStokes/FluidHeatTransferSegregated/index.md diff --git a/modules/navier_stokes/doc/content/source/physics/WCNSLinearFVFluidHeatTransferPhysics.md b/modules/navier_stokes/doc/content/source/physics/WCNSLinearFVFluidHeatTransferPhysics.md new file mode 100644 index 000000000000..811ecd563a63 --- /dev/null +++ b/modules/navier_stokes/doc/content/source/physics/WCNSLinearFVFluidHeatTransferPhysics.md @@ -0,0 +1,59 @@ +# Navier Stokes Fluid Heat Transfer / WCNSLinearFVFluidHeatTransferPhysics + +!syntax description /Physics/NavierStokes/FluidHeatTransferSegregated/WCNSLinearFVFluidHeatTransferPhysics + +## Equation + +This [Physics](Physics/index.md) object creates the kernels and boundary conditions to solve the advection-diffusion +equation for the fluid temperature. +For free flow in a non-porous media: + +!equation +\dfrac{\partial \rho h}{\partial t} + \nabla \cdot (\rho h \mathbf{v}) - \nabla \cdot (k_f \nabla T_f) - Q + \alpha (T_f - T_{ambient}) = 0 + +!alert note +Porous medium treatment is not implemented for the linear finite volume discretization. + +where: + +- $h$ is the fluid enthalpy, computed from the specific heat $c_p$ +- $\rho$ is the fluid density +- $\epsilon$ is the porosity +- $T_f$ is the fluid temperature +- \mathbf{v} is the advecting velocity (clean flow) +- \mathbf{v}_D is the advecting superficial velocity (porous media flow) +- $kappa_f$ the fluid effective thermal conductivity +- $Q$ is the source term, corresponding to energy deposited directly in the fluid +- $\alpha$ is the ambient convection volumetric heat transfer coefficient +- $T_{ambient}$ is the ambient temperature + +The kernels created for flow in a non-porous medium are: + +- [LinearFVEnergyAdvection.md] for advection +- [LinearFVDiffusion.md] for diffusion +- [LinearFVSource.md] for the energy source term +- [LinearFVVolumetricHeatTransfer.md] for the volumetric ambient convection term, if present + + +## Automatically defined variables + +The `WCNSLinearFVFluidHeatTransferPhysics` automatically sets up the variables which are +necessary for solving the energy transport equation: + +- Fluid temperature: + + !listing include/base/NS.h line=std::string T_fluid + +For the default names of other variables used in this action, visit [this site](include/base/NS.h). + + +## Coupling with other Physics + +The heat advection equation can be solved concurrently with the flow equations by combining both the `WCNSLinearFVFluidHeatTransferPhysics` +and the [WCNSLinearFVFlowPhysics.md] using the [!param](/Physics/NavierStokes/FluidHeatTransferSegregated/WCNSLinearFVFluidHeatTransferPhysics/coupled_flow_physics) parameter. + +!syntax parameters /Physics/NavierStokes/FluidHeatTransferSegregated/WCNSLinearFVFluidHeatTransferPhysics + +!syntax inputs /Physics/NavierStokes/FluidHeatTransferSegregated/WCNSLinearFVFluidHeatTransferPhysics + +!syntax children /Physics/NavierStokes/FluidHeatTransferSegregated/WCNSLinearFVFluidHeatTransferPhysics diff --git a/modules/navier_stokes/doc/content/syntax/Physics/NavierStokes/FluidHeatTransfer/index.md b/modules/navier_stokes/doc/content/syntax/Physics/NavierStokes/FluidHeatTransfer/index.md index 9d9cca50dfa1..8096c2748764 100644 --- a/modules/navier_stokes/doc/content/syntax/Physics/NavierStokes/FluidHeatTransfer/index.md +++ b/modules/navier_stokes/doc/content/syntax/Physics/NavierStokes/FluidHeatTransfer/index.md @@ -1,6 +1,6 @@ # Navier Stokes Fluid Heat Transfer Physics -This syntax was created for the [WCNSLinearFVFlowPhysics.md] `Physics`. +This syntax was created for the [WCNSFVFluidHeatTransferPhysics.md] `Physics`. The additional nesting is intended to allow the definition of multiple instances of this `Physics`, with different block restrictions. diff --git a/modules/navier_stokes/doc/content/syntax/Physics/NavierStokes/FluidHeatTransferSegregated/index.md b/modules/navier_stokes/doc/content/syntax/Physics/NavierStokes/FluidHeatTransferSegregated/index.md new file mode 100644 index 000000000000..e4b76eedf980 --- /dev/null +++ b/modules/navier_stokes/doc/content/syntax/Physics/NavierStokes/FluidHeatTransferSegregated/index.md @@ -0,0 +1,11 @@ +# Navier Stokes Linear finite volume Fluid Heat Transfer Physics + +This syntax was created for the [WCNSLinearFVFluidHeatTransferPhysics.md] `Physics`. +The additional nesting is intended to allow the definition of multiple instances of this `Physics`, +with different block restrictions. + +!syntax list /Physics/NavierStokes/FluidHeatTransferSegregated objects=True actions=False subsystems=False + +!syntax list /Physics/NavierStokes/FluidHeatTransferSegregated objects=False actions=False subsystems=True + +!syntax list /Physics/NavierStokes/FluidHeatTransferSegregated objects=False actions=True subsystems=False diff --git a/modules/navier_stokes/src/physics/WCNSFVFluidHeatTransferPhysicsBase.C b/modules/navier_stokes/src/physics/WCNSFVFluidHeatTransferPhysicsBase.C index a2ea67a86c22..71255de881b6 100644 --- a/modules/navier_stokes/src/physics/WCNSFVFluidHeatTransferPhysicsBase.C +++ b/modules/navier_stokes/src/physics/WCNSFVFluidHeatTransferPhysicsBase.C @@ -37,6 +37,7 @@ WCNSFVFluidHeatTransferPhysicsBase::validParams() // Spatial finite volume discretization scheme params.transferParam(NSFVBase::validParams(), "energy_advection_interpolation"); + params.transferParam(NSFVBase::validParams(), "energy_face_interpolation"); params.transferParam(NSFVBase::validParams(), "energy_two_term_bc_expansion"); // Nonlinear equation solver scaling @@ -45,7 +46,7 @@ WCNSFVFluidHeatTransferPhysicsBase::validParams() params.addParamNamesToGroup("specific_heat thermal_conductivity thermal_conductivity_blocks " "use_external_enthalpy_material", "Material properties"); - params.addParamNamesToGroup("energy_advection_interpolation " + params.addParamNamesToGroup("energy_advection_interpolation energy_face_interpolation " "energy_two_term_bc_expansion energy_scaling", "Numerical scheme"); params.addParamNamesToGroup("energy_inlet_types energy_inlet_functors", From 9d7953e73be24e76052849297b3c533899d6e7e8 Mon Sep 17 00:00:00 2001 From: Guillaume Giudicelli Date: Mon, 2 Dec 2024 08:36:00 -0700 Subject: [PATCH 06/16] Add multi-system and linear system support to DumpObjectProblem --- .../include/problems/DumpObjectsProblem.h | 16 +++++--- framework/include/problems/FEProblemBase.h | 10 ++--- .../include/systems/DumpObjectsLinearSystem.h | 37 +++++++++++++++++++ framework/src/problems/DumpObjectsProblem.C | 29 +++++++++++++-- .../src/systems/DumpObjectsLinearSystem.C | 16 ++++++++ 5 files changed, 95 insertions(+), 13 deletions(-) create mode 100644 framework/include/systems/DumpObjectsLinearSystem.h create mode 100644 framework/src/systems/DumpObjectsLinearSystem.C diff --git a/framework/include/problems/DumpObjectsProblem.h b/framework/include/problems/DumpObjectsProblem.h index 2b4207eabf27..2ebf7de7bfc1 100644 --- a/framework/include/problems/DumpObjectsProblem.h +++ b/framework/include/problems/DumpObjectsProblem.h @@ -41,6 +41,8 @@ class DumpObjectsProblem : public FEProblemBase /// output data in solve virtual void solve(unsigned int nl_sys_num) override; + virtual void solveLinearSystem(unsigned int linear_sys_num, + const Moose::PetscSupport::PetscOptions * po) override; virtual void initialSetup() override {} virtual void advanceState() override {} @@ -77,8 +79,6 @@ class DumpObjectsProblem : public FEProblemBase /// store input syntax to build objects generated by a specific action std::map> _generated_syntax; - std::shared_ptr _nl_sys; - /// Whether to include all user-specified parameters in the dump or only parameters that differ from the default value const bool _include_all_user_specified_params; @@ -88,20 +88,25 @@ class DumpObjectsProblem : public FEProblemBase captureDump(addAuxScalarKernel, "AuxScalarKernels") captureDump(addAuxVariable, "AuxVariables") captureDump(addBoundaryCondition, "BCs") - captureDump(addFVBC, "FVBCs") captureDump(addConstraint, "Constraints") captureDump(addDamper, "Dampers") captureDump(addDGKernel, "DGKernels") captureDump(addDiracKernel, "DiracKernels") captureDump(addDistribution, "Distributions") captureDump(addFunction, "Functions") + captureDump(addFunctorMaterial, "FunctorMaterials") + captureDump(addFVBC, "FVBCs") + captureDump(addFVInitialCondition, "FVICs") + captureDump(addFVInterfaceKernel, "FVInterfaceKernels") captureDump(addFVKernel, "FVKernels") + captureDump(addHDGIntegratedBC, "HDGBCs") + captureDump(addHDGKernel, "HDGKernels") captureDump(addIndicator, "Adaptivity/Indicators") captureDump(addInitialCondition, "ICs") - captureDump(addFVInitialCondition, "FVICs") captureDump(addInterfaceKernel, "InterfaceKernels") - captureDump(addFVInterfaceKernel, "FVInterfaceKernels") captureDump(addKernel, "Kernels") + captureDump(addLinearFVBC, "LinearFVBCs") + captureDump(addLinearFVKernel, "LinearFVKernels") captureDump(addMarker, "Adaptivity/Markers") captureDump(addMaterial, "Materials") captureDump(addMultiApp, "MultiApps") @@ -111,6 +116,7 @@ class DumpObjectsProblem : public FEProblemBase captureDump(addSampler, "Samplers") captureDump(addScalarKernel, "ScalarKernels") captureDump(addTransfer, "Transfers") + captureDump(addTimeIntegrator, "TimeIntegrators") captureDump(addVariable, "Variables") captureDump(addVectorPostprocessor, "VectorPostprocessors") // clang-format off diff --git a/framework/include/problems/FEProblemBase.h b/framework/include/problems/FEProblemBase.h index 334be7ff6787..28ba2851a640 100644 --- a/framework/include/problems/FEProblemBase.h +++ b/framework/include/problems/FEProblemBase.h @@ -392,8 +392,8 @@ class FEProblemBase : public SubProblem, public Restartable * @param linear_sys_num The number of the linear system (1,..,num. of lin. systems) * @param po The petsc options for the solve, if not supplied, the defaults are used */ - void solveLinearSystem(const unsigned int linear_sys_num, - const Moose::PetscSupport::PetscOptions * po = nullptr); + virtual void solveLinearSystem(const unsigned int linear_sys_num, + const Moose::PetscSupport::PetscOptions * po = nullptr); ///@{ /** @@ -885,9 +885,9 @@ class FEProblemBase : public SubProblem, public Restartable virtual void addInterfaceMaterial(const std::string & material_name, const std::string & name, InputParameters & parameters); - void addFunctorMaterial(const std::string & functor_material_name, - const std::string & name, - InputParameters & parameters); + virtual void addFunctorMaterial(const std::string & functor_material_name, + const std::string & name, + InputParameters & parameters); /** * Add the MooseVariables and the material properties that the current materials depend on to the diff --git a/framework/include/systems/DumpObjectsLinearSystem.h b/framework/include/systems/DumpObjectsLinearSystem.h new file mode 100644 index 000000000000..0a3310d4e8b7 --- /dev/null +++ b/framework/include/systems/DumpObjectsLinearSystem.h @@ -0,0 +1,37 @@ +//* This file is part of the MOOSE framework +//* https://www.mooseframework.org +//* +//* All rights reserved, see COPYRIGHT for full restrictions +//* https://github.com/idaholab/moose/blob/master/COPYRIGHT +//* +//* Licensed under LGPL 2.1, please see LICENSE for details +//* https://www.gnu.org/licenses/lgpl-2.1.html + +#pragma once + +#include "LinearSystem.h" +#include "libmesh/nonlinear_solver.h" + +namespace libMesh +{ +template +class NumericVector; +} + +/** + * Linear system for dumping objects + * TODO: consider creating a base class for LinearSystems to be able to create a more minimal linear + * system. For now the additional code complexity is not worth it + */ +class DumpObjectsLinearSystem : public LinearSystem +{ +public: + DumpObjectsLinearSystem(FEProblemBase & problem, const std::string & name); + + virtual void solve() override {} + virtual void stopSolve(const ExecFlagType &, const std::set &) override {} + virtual bool converged() override { return true; } + +protected: + NumericVector * _dummy; +}; diff --git a/framework/src/problems/DumpObjectsProblem.C b/framework/src/problems/DumpObjectsProblem.C index 58b2f6cae955..c852443fd234 100644 --- a/framework/src/problems/DumpObjectsProblem.C +++ b/framework/src/problems/DumpObjectsProblem.C @@ -9,6 +9,7 @@ #include "DumpObjectsProblem.h" #include "DumpObjectsNonlinearSystem.h" +#include "DumpObjectsLinearSystem.h" #include "AuxiliarySystem.h" #include @@ -37,12 +38,27 @@ DumpObjectsProblem::validParams() DumpObjectsProblem::DumpObjectsProblem(const InputParameters & parameters) : FEProblemBase(parameters), - _nl_sys(std::make_shared(*this, "nl0")), _include_all_user_specified_params(getParam("include_all_user_specified_params")) { - _nl[0] = _nl_sys; - _solver_systems[0] = std::dynamic_pointer_cast(_nl[0]); + // Make dummy systems based on parameters passed + _solver_systems.resize(0); + for (const auto i : index_range(_nl_sys_names)) + { + const auto & sys_name = _nl_sys_names[i]; + const auto & new_sys = std::make_shared(*this, sys_name); + _solver_systems.push_back(new_sys); + _nl[i] = new_sys; + } + for (const auto i : index_range(_linear_sys_names)) + { + const auto & sys_name = _linear_sys_names[i]; + const auto & new_sys = std::make_shared(*this, sys_name); + _solver_systems.push_back(new_sys); + _linear_systems[i] = new_sys; + } _aux = std::make_shared(*this, "aux0"); + + // Create a dummy assembly for the systems at hand newAssemblyArray(_solver_systems); // Create extra vectors and matrices if any @@ -153,6 +169,13 @@ DumpObjectsProblem::solve(unsigned int) dumpAllGeneratedSyntax(); } +void +DumpObjectsProblem::solveLinearSystem(unsigned int linear_sys_num, + const Moose::PetscSupport::PetscOptions *) +{ + solve(linear_sys_num); +} + void DumpObjectsProblem::dumpGeneratedSyntax(const std::string path) { diff --git a/framework/src/systems/DumpObjectsLinearSystem.C b/framework/src/systems/DumpObjectsLinearSystem.C new file mode 100644 index 000000000000..7851c6b352a3 --- /dev/null +++ b/framework/src/systems/DumpObjectsLinearSystem.C @@ -0,0 +1,16 @@ +//* This file is part of the MOOSE framework +//* https://www.mooseframework.org +//* +//* All rights reserved, see COPYRIGHT for full restrictions +//* https://github.com/idaholab/moose/blob/master/COPYRIGHT +//* +//* Licensed under LGPL 2.1, please see LICENSE for details +//* https://www.gnu.org/licenses/lgpl-2.1.html + +#include "DumpObjectsLinearSystem.h" +#include "FEProblemBase.h" + +DumpObjectsLinearSystem::DumpObjectsLinearSystem(FEProblemBase & problem, const std::string & name) + : LinearSystem(problem, name), _dummy(nullptr) +{ +} From 38cc1ff0ad32a771203878afaac838c65ef0e1a3 Mon Sep 17 00:00:00 2001 From: Guillaume Giudicelli Date: Mon, 2 Dec 2024 09:21:09 -0700 Subject: [PATCH 07/16] Output objects in a different manner: use a task at the end of the setup phase rather than relying on the problem solve calls --- framework/include/actions/DumpObjectsAction.h | 25 ++++++++++++ .../include/problems/DumpObjectsProblem.h | 5 ++- framework/src/actions/DumpObjectsAction.C | 30 ++++++++++++++ framework/src/base/Moose.C | 3 +- framework/src/problems/DumpObjectsProblem.C | 39 +++++++++++++++---- .../src/executioners/SIMPLESolve.C | 4 ++ 6 files changed, 97 insertions(+), 9 deletions(-) create mode 100644 framework/include/actions/DumpObjectsAction.h create mode 100644 framework/src/actions/DumpObjectsAction.C diff --git a/framework/include/actions/DumpObjectsAction.h b/framework/include/actions/DumpObjectsAction.h new file mode 100644 index 000000000000..7c441d4d7c21 --- /dev/null +++ b/framework/include/actions/DumpObjectsAction.h @@ -0,0 +1,25 @@ +//* This file is part of the MOOSE framework +//* https://www.mooseframework.org +//* +//* All rights reserved, see COPYRIGHT for full restrictions +//* https://github.com/idaholab/moose/blob/master/COPYRIGHT +//* +//* Licensed under LGPL 2.1, please see LICENSE for details +//* https://www.gnu.org/licenses/lgpl-2.1.html + +#pragma once + +#include "Action.h" + +/** + * Uses the DumpObjectProblem to output the syntax + */ +class DumpObjectsAction : public Action +{ +public: + static InputParameters validParams(); + + DumpObjectsAction(const InputParameters & params); + + virtual void act() override; +}; diff --git a/framework/include/problems/DumpObjectsProblem.h b/framework/include/problems/DumpObjectsProblem.h index 2ebf7de7bfc1..c54e04e3f41a 100644 --- a/framework/include/problems/DumpObjectsProblem.h +++ b/framework/include/problems/DumpObjectsProblem.h @@ -39,11 +39,14 @@ class DumpObjectsProblem : public FEProblemBase /// output input blocks for all paths void dumpAllGeneratedSyntax() const; - /// output data in solve + /// output data in solve (if ever called) virtual void solve(unsigned int nl_sys_num) override; virtual void solveLinearSystem(unsigned int linear_sys_num, const Moose::PetscSupport::PetscOptions * po) override; + // output data (we expect to call this from the DumpObjectsAction) + void printObjects(); + virtual void initialSetup() override {} virtual void advanceState() override {} virtual void timestepSetup() override {} diff --git a/framework/src/actions/DumpObjectsAction.C b/framework/src/actions/DumpObjectsAction.C new file mode 100644 index 000000000000..aa3b7460c38d --- /dev/null +++ b/framework/src/actions/DumpObjectsAction.C @@ -0,0 +1,30 @@ +//* This file is part of the MOOSE framework +//* https://www.mooseframework.org +//* +//* All rights reserved, see COPYRIGHT for full restrictions +//* https://github.com/idaholab/moose/blob/master/COPYRIGHT +//* +//* Licensed under LGPL 2.1, please see LICENSE for details +//* https://www.gnu.org/licenses/lgpl-2.1.html + +#include "DumpObjectsAction.h" +#include "DumpObjectsProblem.h" + +registerMooseAction("MooseApp", DumpObjectsAction, "dump_objects"); + +InputParameters +DumpObjectsAction::validParams() +{ + InputParameters params = Action::validParams(); + return params; +} + +DumpObjectsAction::DumpObjectsAction(const InputParameters & params) : Action(params) {} + +void +DumpObjectsAction::act() +{ + if (!dynamic_cast(_problem.get())) + return; + dynamic_cast(_problem.get())->printObjects(); +} diff --git a/framework/src/base/Moose.C b/framework/src/base/Moose.C index 91250cd67442..376f764460b7 100644 --- a/framework/src/base/Moose.C +++ b/framework/src/base/Moose.C @@ -449,7 +449,8 @@ registerActions(Syntax & syntax, { Registry::registerActionsTo(action_factory, obj_labels); - // TODO: Why is this here? + // Add these actions here so they are always executed last, without setting any dependency + registerTask("dump_objects", false); registerTask("finish_input_file_output", false); } diff --git a/framework/src/problems/DumpObjectsProblem.C b/framework/src/problems/DumpObjectsProblem.C index c852443fd234..90ad64ed0e3a 100644 --- a/framework/src/problems/DumpObjectsProblem.C +++ b/framework/src/problems/DumpObjectsProblem.C @@ -8,9 +8,11 @@ //* https://www.gnu.org/licenses/lgpl-2.1.html #include "DumpObjectsProblem.h" +#include "DumpObjectsAction.h" #include "DumpObjectsNonlinearSystem.h" #include "DumpObjectsLinearSystem.h" #include "AuxiliarySystem.h" +#include "InputParameters.h" #include #include "libmesh/string_to_enum.h" @@ -33,6 +35,14 @@ DumpObjectsProblem::validParams() true, "Whether to include all parameters that have been specified by a user in the dump, even if " "they match the default value of the parameter in the Factory"); + + // Change the default because any complex solve or executioners needs the problem to perform its + // setup duties (all the calls in initialSetup()) which are skipped by the DumpObjectsProblem + params.addParam( + "solve", + false, + "Whether to attempt to solve the Problem. This will only cause additional outputs of the " + "objects and their parameters. This is unlikely to succeed with more complex executioners."); return params; } @@ -66,6 +76,16 @@ DumpObjectsProblem::DumpObjectsProblem(const InputParameters & parameters) // Create extra solution vectors if any createTagSolutions(); + + // Add an action to call printObjects at the end of the action/tasks phase + // NOTE: We previously relied on problem.solve() but some executioners (SIMPLE in NavierStokes) do + // not support this + auto action_params = _app.getActionFactory().getValidParams("DumpObjectsAction"); + action_params.applyParameters(parameters); + // action_params.finalize(""); + auto dump_objects_action = + _app.getActionFactory().create("DumpObjectsAction", "dump_objects", action_params); + _app.actionWarehouse().addActionBlock(dump_objects_action); } std::string @@ -161,6 +181,18 @@ DumpObjectsProblem::dumpVariableHelper(const std::string & system, void DumpObjectsProblem::solve(unsigned int) +{ + printObjects(); +} + +void +DumpObjectsProblem::solveLinearSystem(unsigned int, const Moose::PetscSupport::PetscOptions *) +{ + printObjects(); +} + +void +DumpObjectsProblem::printObjects() { const auto path = getParam("dump_path"); if (path != "all") @@ -169,13 +201,6 @@ DumpObjectsProblem::solve(unsigned int) dumpAllGeneratedSyntax(); } -void -DumpObjectsProblem::solveLinearSystem(unsigned int linear_sys_num, - const Moose::PetscSupport::PetscOptions *) -{ - solve(linear_sys_num); -} - void DumpObjectsProblem::dumpGeneratedSyntax(const std::string path) { diff --git a/modules/navier_stokes/src/executioners/SIMPLESolve.C b/modules/navier_stokes/src/executioners/SIMPLESolve.C index c6157b85d6eb..5c94e83bcda3 100644 --- a/modules/navier_stokes/src/executioners/SIMPLESolve.C +++ b/modules/navier_stokes/src/executioners/SIMPLESolve.C @@ -271,6 +271,10 @@ SIMPLESolve::solveAdvectedSystem(const unsigned int system_num, bool SIMPLESolve::solve() { + // Do not solve if problem is set not to + if (!_problem.shouldSolve()) + return true; + // Dummy solver parameter file which is needed for switching petsc options SolverParams solver_params; solver_params._type = Moose::SolveType::ST_LINEAR; From 62d88d7b03ff85b970e2ac45168dba77ec98dd65 Mon Sep 17 00:00:00 2001 From: Guillaume Giudicelli Date: Mon, 2 Dec 2024 09:45:55 -0700 Subject: [PATCH 08/16] Add boolean for two term BCs for momentum outflow bcs --- modules/navier_stokes/src/physics/WCNSFVFlowPhysics.C | 10 ++++------ .../navier_stokes/src/physics/WCNSFVFlowPhysicsBase.C | 8 +++++--- .../src/physics/WCNSLinearFVFlowPhysics.C | 1 + .../segregated-comparison/segregated-linear-physics.i | 1 + 4 files changed, 11 insertions(+), 9 deletions(-) diff --git a/modules/navier_stokes/src/physics/WCNSFVFlowPhysics.C b/modules/navier_stokes/src/physics/WCNSFVFlowPhysics.C index 307147cb578e..748c9dd1b507 100644 --- a/modules/navier_stokes/src/physics/WCNSFVFlowPhysics.C +++ b/modules/navier_stokes/src/physics/WCNSFVFlowPhysics.C @@ -63,7 +63,6 @@ WCNSFVFlowPhysics::validParams() params.transferParam(NSFVBase::validParams(), "mass_advection_interpolation"); params.transferParam(NSFVBase::validParams(), "pressure_allow_expansion_on_bernoulli_faces"); - params.transferParam(NSFVBase::validParams(), "momentum_two_term_bc_expansion"); // Nonlinear solver parameters params.transferParam(NSFVBase::validParams(), "mass_scaling"); @@ -75,11 +74,10 @@ WCNSFVFlowPhysics::validParams() "porosity_interface_pressure_treatment pressure_allow_expansion_on_bernoulli_faces " "porosity_smoothing_layers use_friction_correction consistent_scaling", "Flow medium discontinuity treatment"); - params.addParamNamesToGroup( - "pressure_face_interpolation momentum_face_interpolation " - "mass_advection_interpolation momentum_advection_interpolation " - "momentum_two_term_bc_expansion mass_scaling momentum_scaling characteristic_speed", - "Numerical scheme"); + params.addParamNamesToGroup("pressure_face_interpolation momentum_face_interpolation " + "mass_advection_interpolation momentum_advection_interpolation " + "mass_scaling momentum_scaling characteristic_speed", + "Numerical scheme"); // TODO Add default preconditioning and move scaling parameters to a preconditioning group diff --git a/modules/navier_stokes/src/physics/WCNSFVFlowPhysicsBase.C b/modules/navier_stokes/src/physics/WCNSFVFlowPhysicsBase.C index a92c6ff9e560..29137c348d12 100644 --- a/modules/navier_stokes/src/physics/WCNSFVFlowPhysicsBase.C +++ b/modules/navier_stokes/src/physics/WCNSFVFlowPhysicsBase.C @@ -67,6 +67,7 @@ WCNSFVFlowPhysicsBase::validParams() // Specify the numerical schemes for interpolations of velocity and pressure params.transferParam(NSFVBase::validParams(), "velocity_interpolation"); params.transferParam(NSFVBase::validParams(), "momentum_advection_interpolation"); + params.transferParam(NSFVBase::validParams(), "momentum_two_term_bc_expansion"); params.transferParam(NSFVBase::validParams(), "pressure_two_term_bc_expansion"); MooseEnum coeff_interp_method("average harmonic", "harmonic"); params.addParam("mu_interp_method", @@ -83,9 +84,10 @@ WCNSFVFlowPhysicsBase::validParams() "Outlet boundary conditions"); params.addParamNamesToGroup("wall_boundaries momentum_wall_types momentum_wall_functors", "Wall boundary conditions"); - params.addParamNamesToGroup("velocity_interpolation momentum_advection_interpolation " - "pressure_two_term_bc_expansion mu_interp_method", - "Numerical scheme"); + params.addParamNamesToGroup( + "velocity_interpolation momentum_advection_interpolation " + "momentum_two_term_bc_expansion pressure_two_term_bc_expansion mu_interp_method", + "Numerical scheme"); params.addParamNamesToGroup("thermal_expansion", "Gravity treatment"); return params; diff --git a/modules/navier_stokes/src/physics/WCNSLinearFVFlowPhysics.C b/modules/navier_stokes/src/physics/WCNSLinearFVFlowPhysics.C index 0a3046682feb..55948b7b264e 100644 --- a/modules/navier_stokes/src/physics/WCNSLinearFVFlowPhysics.C +++ b/modules/navier_stokes/src/physics/WCNSLinearFVFlowPhysics.C @@ -370,6 +370,7 @@ WCNSLinearFVFlowPhysics::addINSOutletBC() const std::string bc_type = "LinearFVAdvectionDiffusionOutflowBC"; InputParameters params = getFactory().getValidParams(bc_type); params.set>("boundary") = {outlet_bdy}; + params.set("use_two_term_expansion") = getParam("momentum_two_term_bc_expansion"); for (const auto d : make_range(dimension())) { diff --git a/modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/segregated-comparison/segregated-linear-physics.i b/modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/segregated-comparison/segregated-linear-physics.i index d87fb6f42a71..a121c7c58d5b 100644 --- a/modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/segregated-comparison/segregated-linear-physics.i +++ b/modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/segregated-comparison/segregated-linear-physics.i @@ -44,6 +44,7 @@ advected_interp_method = 'average' pressure_functors = '1.4' orthogonality_correction = false + momentum_two_term_bc_expansion = false pressure_two_term_bc_expansion = false momentum_advection_interpolation = ${advected_interp_method} [] From 9b997057491488a22a575119efbd5f463ec57985 Mon Sep 17 00:00:00 2001 From: Guillaume Giudicelli Date: Mon, 2 Dec 2024 10:40:33 -0700 Subject: [PATCH 09/16] Handle specific heat better in kernel and physics --- .../src/linearfvkernels/LinearFVEnergyAdvection.C | 2 +- .../src/physics/WCNSLinearFVFluidHeatTransferPhysics.C | 3 +++ 2 files changed, 4 insertions(+), 1 deletion(-) diff --git a/modules/navier_stokes/src/linearfvkernels/LinearFVEnergyAdvection.C b/modules/navier_stokes/src/linearfvkernels/LinearFVEnergyAdvection.C index a8fb51986640..aa4f98c455b5 100644 --- a/modules/navier_stokes/src/linearfvkernels/LinearFVEnergyAdvection.C +++ b/modules/navier_stokes/src/linearfvkernels/LinearFVEnergyAdvection.C @@ -21,7 +21,7 @@ LinearFVEnergyAdvection::validParams() params.addClassDescription("Represents the matrix and right hand side contributions of an " "advection term for the energy e.g. h=int(cp dT). A user may still " "override what quantity is advected, but the default is temperature."); - params.addParam("cp", 1, "Specific heat value"); + params.addRequiredParam("cp", "Specific heat value"); params.addRequiredParam( "rhie_chow_user_object", "The rhie-chow user-object which is used to determine the face velocity."); diff --git a/modules/navier_stokes/src/physics/WCNSLinearFVFluidHeatTransferPhysics.C b/modules/navier_stokes/src/physics/WCNSLinearFVFluidHeatTransferPhysics.C index bd349a718fd3..97f961018234 100644 --- a/modules/navier_stokes/src/physics/WCNSLinearFVFluidHeatTransferPhysics.C +++ b/modules/navier_stokes/src/physics/WCNSLinearFVFluidHeatTransferPhysics.C @@ -91,6 +91,9 @@ WCNSLinearFVFluidHeatTransferPhysics::addINSEnergyAdvectionKernels() InputParameters params = getFactory().getValidParams(kernel_type); params.set("variable") = _fluid_temperature_name; assignBlocks(params, _blocks); + if (!MooseUtils::isFloat(_specific_heat_name)) + paramError(NS::cp, "Must be a Real number. Functors not supported at this time"); + params.set("cp") = std::atof(_specific_heat_name.c_str()); params.set("rhie_chow_user_object") = _flow_equations_physics->rhieChowUOName(); params.set("advected_interp_method") = getParam("energy_advection_interpolation"); From 523316adbbdd41dec344c2363eb45e6640b3859b Mon Sep 17 00:00:00 2001 From: Guillaume Giudicelli Date: Mon, 16 Dec 2024 17:27:48 -0700 Subject: [PATCH 10/16] Address Alex-s review - more compact syntax - remove extra header - remove any form of solve from DumpObjectsProblem but dont both printing again Co-authored-by: Alex Lindsay --- framework/include/problems/DumpObjectsProblem.h | 8 +++++--- framework/include/systems/DumpObjectsLinearSystem.h | 1 - framework/src/actions/DumpObjectsAction.C | 7 ++++--- framework/src/problems/DumpObjectsProblem.C | 13 ------------- 4 files changed, 9 insertions(+), 20 deletions(-) diff --git a/framework/include/problems/DumpObjectsProblem.h b/framework/include/problems/DumpObjectsProblem.h index c54e04e3f41a..468f9c9416e5 100644 --- a/framework/include/problems/DumpObjectsProblem.h +++ b/framework/include/problems/DumpObjectsProblem.h @@ -40,9 +40,11 @@ class DumpObjectsProblem : public FEProblemBase void dumpAllGeneratedSyntax() const; /// output data in solve (if ever called) - virtual void solve(unsigned int nl_sys_num) override; - virtual void solveLinearSystem(unsigned int linear_sys_num, - const Moose::PetscSupport::PetscOptions * po) override; + virtual void solve(unsigned int /*nl_sys_num*/) override {} + virtual void solveLinearSystem(unsigned int /*linear_sys_num*/, + const Moose::PetscSupport::PetscOptions *) override + { + } // output data (we expect to call this from the DumpObjectsAction) void printObjects(); diff --git a/framework/include/systems/DumpObjectsLinearSystem.h b/framework/include/systems/DumpObjectsLinearSystem.h index 0a3310d4e8b7..b5aa56056a30 100644 --- a/framework/include/systems/DumpObjectsLinearSystem.h +++ b/framework/include/systems/DumpObjectsLinearSystem.h @@ -10,7 +10,6 @@ #pragma once #include "LinearSystem.h" -#include "libmesh/nonlinear_solver.h" namespace libMesh { diff --git a/framework/src/actions/DumpObjectsAction.C b/framework/src/actions/DumpObjectsAction.C index aa3b7460c38d..856ee60ab807 100644 --- a/framework/src/actions/DumpObjectsAction.C +++ b/framework/src/actions/DumpObjectsAction.C @@ -24,7 +24,8 @@ DumpObjectsAction::DumpObjectsAction(const InputParameters & params) : Action(pa void DumpObjectsAction::act() { - if (!dynamic_cast(_problem.get())) - return; - dynamic_cast(_problem.get())->printObjects(); + mooseAssert(dynamic_cast(_problem.get()), + "Problem should be DumpObjectProblem"); + if (auto * const dop = dynamic_cast(_problem.get()); dop) + dop->printObjects(); } diff --git a/framework/src/problems/DumpObjectsProblem.C b/framework/src/problems/DumpObjectsProblem.C index 90ad64ed0e3a..ac3fe6580a51 100644 --- a/framework/src/problems/DumpObjectsProblem.C +++ b/framework/src/problems/DumpObjectsProblem.C @@ -82,7 +82,6 @@ DumpObjectsProblem::DumpObjectsProblem(const InputParameters & parameters) // not support this auto action_params = _app.getActionFactory().getValidParams("DumpObjectsAction"); action_params.applyParameters(parameters); - // action_params.finalize(""); auto dump_objects_action = _app.getActionFactory().create("DumpObjectsAction", "dump_objects", action_params); _app.actionWarehouse().addActionBlock(dump_objects_action); @@ -179,18 +178,6 @@ DumpObjectsProblem::dumpVariableHelper(const std::string & system, // clang-format on } -void -DumpObjectsProblem::solve(unsigned int) -{ - printObjects(); -} - -void -DumpObjectsProblem::solveLinearSystem(unsigned int, const Moose::PetscSupport::PetscOptions *) -{ - printObjects(); -} - void DumpObjectsProblem::printObjects() { From 44a41aeb4b5cc5641d57168c920da543e57e2c1f Mon Sep 17 00:00:00 2001 From: Guillaume Giudicelli Date: Tue, 24 Dec 2024 02:19:51 +0100 Subject: [PATCH 11/16] Add a Physics for linear finite volume for two phase flow --- .../WCNSLinearFVTwoPhaseMixturePhysics.h | 81 +++++ .../navier_stokes/src/base/NavierStokesApp.C | 2 + .../WCNSLinearFVTwoPhaseMixturePhysics.C | 312 ++++++++++++++++++ 3 files changed, 395 insertions(+) create mode 100644 modules/navier_stokes/include/physics/WCNSLinearFVTwoPhaseMixturePhysics.h create mode 100644 modules/navier_stokes/src/physics/WCNSLinearFVTwoPhaseMixturePhysics.C diff --git a/modules/navier_stokes/include/physics/WCNSLinearFVTwoPhaseMixturePhysics.h b/modules/navier_stokes/include/physics/WCNSLinearFVTwoPhaseMixturePhysics.h new file mode 100644 index 000000000000..a48d6b192213 --- /dev/null +++ b/modules/navier_stokes/include/physics/WCNSLinearFVTwoPhaseMixturePhysics.h @@ -0,0 +1,81 @@ +//* This file is part of the MOOSE framework +//* https://www.mooseframework.org +//* +//* All rights reserved, see COPYRIGHT for full restrictions +//* https://github.com/idaholab/moose/blob/master/COPYRIGHT +//* +//* Licensed under LGPL 2.1, please see LICENSE for details +//* https://www.gnu.org/licenses/lgpl-2.1.html + +#pragma once + +#include "WCNSLinearFVScalarTransportPhysics.h" + +class WCNSLinearFVFluidHeatTransferPhysics; + +/** + * Creates all the objects needed to solve the mixture terms for the weakly-compressible + * and incompressible two-phase equations. + * Can also add a phase transport equation + */ +class WCNSLinearFVTwoPhaseMixturePhysics final : public WCNSLinearFVScalarTransportPhysics +{ +public: + static InputParameters validParams(); + + WCNSLinearFVTwoPhaseMixturePhysics(const InputParameters & parameters); + +private: + virtual void addFVKernels() override; + virtual void addMaterials() override; + + /// Adds the slip velocity parameters + virtual void setSlipVelocityParams(InputParameters & params) const override; + + /** + * Functions adding kernels for the other physics + */ + void addPhaseInterfaceTerm(); + void addPhaseChangeEnergySource(); + void addPhaseDriftFluxTerm(); + void addAdvectionSlipTerm(); + + /// Fluid heat transfer physics + const WCNSLinearFVFluidHeatTransferPhysics * _fluid_energy_physics; + + /// Convenience boolean to keep track of whether the phase transport equation is requested + const bool _add_phase_equation; + /// Convenience boolean to keep track of whether the fluid energy equation is present + bool _has_energy_equation; + + /// Name of the first phase fraction (usually, liquid) + const MooseFunctorName _phase_1_fraction_name; + /// Name of the second phase fraction (usually, dispersed or advected by the liquid) + const MooseFunctorName _phase_2_fraction_name; + + /// Name of the density of the other phase + const MooseFunctorName _phase_1_density; + /// Name of the dyanmic viscosity of the other phase + const MooseFunctorName _phase_1_viscosity; + /// Name of the specific heat of the other phase + const MooseFunctorName _phase_1_specific_heat; + /// Name of the thermal conductivity of the other phase + const MooseFunctorName _phase_1_thermal_conductivity; + + /// Name of the density of the other phase + const MooseFunctorName _phase_2_density; + /// Name of the dynamic viscosity of the other phase + const MooseFunctorName _phase_2_viscosity; + /// Name of the specific heat of the other phase + const MooseFunctorName _phase_2_specific_heat; + /// Name of the thermal conductivity of the other phase + const MooseFunctorName _phase_2_thermal_conductivity; + + /// Whether to define the mixture model internally or use fluid properties instead + const bool _use_external_mixture_properties; + + /// Whether to add the drift flux momentum terms to each component momentum equation + const bool _use_drift_flux; + /// Whether to add the advection slip term to each component of the momentum equation + const bool _use_advection_slip; +}; diff --git a/modules/navier_stokes/src/base/NavierStokesApp.C b/modules/navier_stokes/src/base/NavierStokesApp.C index 9542cc833591..bf0c679fe47b 100644 --- a/modules/navier_stokes/src/base/NavierStokesApp.C +++ b/modules/navier_stokes/src/base/NavierStokesApp.C @@ -61,6 +61,8 @@ associateSyntaxInner(Syntax & syntax, ActionFactory & /*action_factory*/) registerSyntax("WCNSFVTurbulencePhysics", "Physics/NavierStokes/Turbulence/*"); registerSyntax("PNSFVSolidHeatTransferPhysics", "Physics/NavierStokes/SolidHeatTransfer/*"); registerSyntax("WCNSFVTwoPhaseMixturePhysics", "Physics/NavierStokes/TwoPhaseMixture/*"); + registerSyntax("WCNSLinearFVTwoPhaseMixturePhysics", + "Physics/NavierStokes/TwoPhaseMixtureSegregated/*"); // Create the Action syntax registerSyntax("CNSAction", "Modules/CompressibleNavierStokes"); diff --git a/modules/navier_stokes/src/physics/WCNSLinearFVTwoPhaseMixturePhysics.C b/modules/navier_stokes/src/physics/WCNSLinearFVTwoPhaseMixturePhysics.C new file mode 100644 index 000000000000..b784aa0624c8 --- /dev/null +++ b/modules/navier_stokes/src/physics/WCNSLinearFVTwoPhaseMixturePhysics.C @@ -0,0 +1,312 @@ +//* This file is part of the MOOSE framework +//* https://www.mooseframework.org +//* +//* All rights reserved, see COPYRIGHT for full restrictions +//* https://github.com/idaholab/moose/blob/master/COPYRIGHT +//* +//* Licensed under LGPL 2.1, please see LICENSE for details +//* https://www.gnu.org/licenses/lgpl-2.1.html + +#include "WCNSLinearFVTwoPhaseMixturePhysics.h" +#include "WCNSFVTwoPhaseMixturePhysics.h" +#include "WCNSFVFluidHeatTransferPhysics.h" +#include "WCNSFVFlowPhysics.h" + +registerNavierStokesPhysicsBaseTasks("NavierStokesApp", WCNSLinearFVTwoPhaseMixturePhysics); +registerWCNSFVScalarTransportBaseTasks("NavierStokesApp", WCNSLinearFVTwoPhaseMixturePhysics); +registerMooseAction("NavierStokesApp", WCNSLinearFVTwoPhaseMixturePhysics, "add_material"); + +InputParameters +WCNSLinearFVTwoPhaseMixturePhysics::validParams() +{ + // The parameters are mostly the same being the linear and nonlinear version + InputParameters params = WCNSLinearFVScalarTransportPhysics::validParams(); + WCNSFVTwoPhaseMixturePhysics::renamePassiveScalarToMixtureParams(params); + // The flow physics is obtained from the scalar transport base class + // The fluid heat transfer physics is retrieved even if unspecified + params.addParam( + "fluid_heat_transfer_physics", + "NavierStokesFV", + "WCNSLinearFVFluidHeatTransferPhysics generating the fluid energy equation"); + params += WCNSFVTwoPhaseMixturePhysics::commonMixtureParams(); + params.addParamNamesToGroup("fluid_heat_transfer_physics", "Phase change"); + params.addClassDescription("Define the additional terms for a mixture model for the two phase " + "weakly-compressible Navier Stokes equations using the linearized " + "segregated finite volume discretization"); + + return params; +} + +WCNSLinearFVTwoPhaseMixturePhysics::WCNSLinearFVTwoPhaseMixturePhysics( + const InputParameters & parameters) + : WCNSLinearFVScalarTransportPhysics(parameters), + _add_phase_equation(_has_scalar_equation), + _phase_1_fraction_name(getParam("phase_1_fraction_name")), + _phase_2_fraction_name(_passive_scalar_names[0]), + _phase_1_density(getParam("phase_1_density_name")), + _phase_1_viscosity(getParam("phase_1_viscosity_name")), + _phase_1_specific_heat(getParam("phase_1_specific_heat_name")), + _phase_1_thermal_conductivity(getParam("phase_1_thermal_conductivity_name")), + _phase_2_density(getParam("phase_2_density_name")), + _phase_2_viscosity(getParam("phase_2_viscosity_name")), + _phase_2_specific_heat(getParam("phase_2_specific_heat_name")), + _phase_2_thermal_conductivity(getParam("phase_2_thermal_conductivity_name")), + _use_external_mixture_properties(getParam("use_external_mixture_properties")), + _use_drift_flux(getParam("add_drift_flux_momentum_terms")), + _use_advection_slip(getParam("add_advection_slip_term")) +{ + // Check that only one scalar was passed, as we are using vector parameters + if (_passive_scalar_names.size() > 1) + paramError("phase_fraction_name", "Only one phase fraction currently supported."); + if (_passive_scalar_inlet_functors.size() > 1) + paramError("phase_fraction_inlet_functors", "Only one phase fraction currently supported"); + + // Retrieve the fluid energy equation if it exists + if (isParamValid("fluid_heat_transfer_physics")) + { + _fluid_energy_physics = getCoupledPhysics( + getParam("fluid_heat_transfer_physics"), true); + // Check for a missing parameter / do not support isolated physics for now + if (!_fluid_energy_physics && + !getCoupledPhysics(true).empty()) + paramError( + "fluid_heat_transfer_physics", + "We currently do not support creating both a phase transport equation and fluid heat " + "transfer physics that are not coupled together"); + if (_fluid_energy_physics && _fluid_energy_physics->hasEnergyEquation()) + _has_energy_equation = true; + else + _has_energy_equation = false; + } + else + { + _has_energy_equation = false; + _fluid_energy_physics = nullptr; + } + + // Check that the mixture parameters are correctly in use in the other physics + if (_has_energy_equation) + { + if (_fluid_energy_physics->densityName() != "rho_mixture") + mooseError("Density name should for Physics '", + _fluid_energy_physics->name(), + "' should be 'rho_mixture'"); + if (_fluid_energy_physics->getSpecificHeatName() != "cp_mixture") + mooseError("Specific heat name should for Physics '", + _fluid_energy_physics->name(), + "' should be 'cp_mixture'"); + } + if (_flow_equations_physics) + { + if (_flow_equations_physics->densityName() != "rho_mixture") + mooseError("Density name should for Physics ,", + _flow_equations_physics->name(), + "' should be 'rho_mixture'"); + } + + if (_verbose) + { + if (_flow_equations_physics) + mooseInfoRepeated("Coupled to fluid flow physics " + _flow_equations_physics->name()); + if (_has_energy_equation) + mooseInfoRepeated("Coupled to fluid heat transfer physics " + _fluid_energy_physics->name()); + } + + // Parameter checking + // The two models are not consistent + if (isParamSetByUser("alpha_exchange") && getParam("add_phase_change_energy_term")) + paramError("alpha_exchange", + "A phase exchange coefficient cannot be specified if the phase change is handled " + "with a phase change heat loss model"); + if (_phase_1_fraction_name == _phase_2_fraction_name) + paramError("phase_1_fraction_name", + "First phase fraction name should be different from second phase fraction name"); + if (_use_drift_flux && _use_advection_slip) + paramError("add_drift_flux_momentum_terms", + "Drift flux model cannot be used at the same time as the advection slip model"); + if (!getParam("add_drift_flux_momentum_terms")) + errorDependentParameter("add_drift_flux_momentum_terms", "true", {"density_interp_method"}); + if (!getParam("use_dispersed_phase_drag_model")) + errorDependentParameter("use_dispersed_phase_drag_model", "true", {"particle_diameter"}); +} + +void +WCNSLinearFVTwoPhaseMixturePhysics::addFVKernels() +{ + WCNSLinearFVScalarTransportPhysics::addFVKernels(); + + if (_add_phase_equation && isParamSetByUser("alpha_exchange")) + addPhaseInterfaceTerm(); + + if (_fluid_energy_physics && _fluid_energy_physics->hasEnergyEquation() && + getParam("add_phase_change_energy_term")) + addPhaseChangeEnergySource(); + + if (_flow_equations_physics && _flow_equations_physics->hasFlowEquations() && _use_drift_flux) + addPhaseDriftFluxTerm(); + if (_flow_equations_physics && _flow_equations_physics->hasFlowEquations() && _use_advection_slip) + addAdvectionSlipTerm(); +} + +void +WCNSLinearFVTwoPhaseMixturePhysics::setSlipVelocityParams(InputParameters & params) const +{ + params.set("u_slip") = "vel_slip_x"; + if (dimension() >= 2) + params.set("v_slip") = "vel_slip_y"; + if (dimension() >= 3) + params.set("w_slip") = "vel_slip_z"; +} + +void +WCNSLinearFVTwoPhaseMixturePhysics::addPhaseInterfaceTerm() +{ + mooseError("Phase interface term not implemented at this time for linear finite volume"); +} + +void +WCNSLinearFVTwoPhaseMixturePhysics::addPhaseChangeEnergySource() +{ + mooseError("Phase change energy source not implemented at this time for linear finite volume"); +} + +void +WCNSLinearFVTwoPhaseMixturePhysics::addPhaseDriftFluxTerm() +{ + mooseError("Phase drift flux not implemented at this time for linear finite volume"); +} + +void +WCNSLinearFVTwoPhaseMixturePhysics::addAdvectionSlipTerm() +{ + mooseError("Phase advection slip not implemented at this time for linear finite volume"); +} + +void +WCNSLinearFVTwoPhaseMixturePhysics::addMaterials() +{ + // Add the phase fraction variable, for output purposes mostly + if (!getProblem().hasFunctor(_phase_1_fraction_name, /*thread_id=*/0)) + { + auto params = getFactory().getValidParams("ADParsedFunctorMaterial"); + assignBlocks(params, _blocks); + params.set("expression") = "1 - " + _phase_2_fraction_name; + params.set>("functor_names") = {_phase_2_fraction_name}; + params.set("property_name") = _phase_1_fraction_name; + params.set>("output_properties") = {_phase_1_fraction_name}; + params.set>("outputs") = {"all"}; + getProblem().addMaterial("ADParsedFunctorMaterial", prefix() + "phase_1_fraction", params); + + // One of the phase fraction should exist though (either as a variable or set by a + // NSLiquidFractionAux) + if (!getProblem().hasFunctor(_phase_2_fraction_name, /*thread_id=*/0)) + paramError("Phase 2 fraction should be defined as a variable or auxiliary variable"); + } + if (!getProblem().hasFunctor(_phase_2_fraction_name, /*thread_id=*/0)) + { + auto params = getFactory().getValidParams("ADParsedFunctorMaterial"); + assignBlocks(params, _blocks); + params.set("expression") = "1 - " + _phase_1_fraction_name; + params.set>("functor_names") = {_phase_1_fraction_name}; + params.set("property_name") = _phase_2_fraction_name; + params.set>("output_properties") = {_phase_2_fraction_name}; + params.set>("outputs") = {"all"}; + getProblem().addMaterial("ADParsedFunctorMaterial", prefix() + "phase_2_fraction", params); + } + + // Compute mixture properties + if (!_use_external_mixture_properties) + { + auto params = getFactory().getValidParams("NSFVMixtureFunctorMaterial"); + assignBlocks(params, _blocks); + params.set>("prop_names") = { + "rho_mixture", "mu_mixture", "cp_mixture", "k_mixture"}; + // The phase_1 and phase_2 assignments are only local to this object. + // We use the phase 2 variable to save a functor evaluation as we expect + // the phase 2 variable to be a nonlinear variable in the phase transport equation + params.set>("phase_2_names") = {_phase_1_density, + _phase_1_viscosity, + _phase_1_specific_heat, + _phase_1_thermal_conductivity}; + params.set>("phase_1_names") = {_phase_2_density, + _phase_2_viscosity, + _phase_2_specific_heat, + _phase_2_thermal_conductivity}; + params.set("phase_1_fraction") = _phase_2_fraction_name; + if (getParam("output_all_properties")) + params.set>("outputs") = {"all"}; + getProblem().addMaterial("NSFVMixtureFunctorMaterial", prefix() + "mixture_material", params); + } + + // Compute slip terms as functors, used by the drift flux kernels + if (_use_advection_slip || _use_drift_flux || _add_phase_equation) + { + const std::vector vel_components = {"u", "v", "w"}; + const std::vector components = {"x", "y", "z"}; + for (const auto dim : make_range(dimension())) + { + auto params = getFactory().getValidParams("WCNSFV2PSlipVelocityFunctorMaterial"); + assignBlocks(params, _blocks); + params.set("slip_velocity_name") = "vel_slip_" + components[dim]; + params.set("momentum_component") = components[dim]; + for (const auto j : make_range(dimension())) + params.set>(vel_components[j]) = { + _flow_equations_physics->getVelocityNames()[j]}; + params.set(NS::density) = _phase_1_density; + params.set(NS::mu) = "mu_mixture"; + params.set("rho_d") = _phase_2_density; + params.set("gravity") = _flow_equations_physics->gravityVector(); + if (isParamValid("slip_linear_friction_name")) + params.set("linear_coef_name") = + getParam("slip_linear_friction_name"); + else if (getParam("use_dispersed_phase_drag_model")) + params.set("linear_coef_name") = "Darcy_coefficient"; + else if (_flow_equations_physics) + { + if (!_flow_equations_physics->getLinearFrictionCoefName().empty()) + params.set("linear_coef_name") = + _flow_equations_physics->getLinearFrictionCoefName(); + else + params.set("linear_coef_name") = "0"; + } + else + paramError("slip_linear_friction_name", + "WCNSFV2PSlipVelocityFunctorMaterial created by this Physics required a scalar " + "field linear friction factor."); + params.set("particle_diameter") = + getParam("particle_diameter"); + if (getParam("output_all_properties")) + { + if (!isTransient()) + params.set>("outputs") = {"all"}; + else + paramInfo("output_all_properties", + "Slip velocity functor material output currently unsupported in Physics " + "in transient conditions."); + } + getProblem().addMaterial( + "WCNSFV2PSlipVelocityFunctorMaterial", prefix() + "slip_" + components[dim], params); + } + } + + // Add a default drag model for a dispersed phase + if (getParam("use_dispersed_phase_drag_model")) + { + const std::vector vel_components = {"u", "v", "w"}; + + auto params = getFactory().getValidParams("NSFVDispersePhaseDragFunctorMaterial"); + assignBlocks(params, _blocks); + params.set("drag_coef_name") = "Darcy_coefficient"; + for (const auto j : make_range(dimension())) + params.set(vel_components[j]) = { + _flow_equations_physics->getVelocityNames()[j]}; + params.set(NS::density) = "rho_mixture"; + params.set(NS::mu) = "mu_mixture"; + params.set("particle_diameter") = + getParam("particle_diameter"); + if (getParam("output_all_properties")) + params.set>("outputs") = {"all"}; + getProblem().addMaterial( + "NSFVDispersePhaseDragFunctorMaterial", prefix() + "dispersed_drag", params); + } +} From 9031694f7c1d771d6509f42f95ce98047b2931be Mon Sep 17 00:00:00 2001 From: Guillaume Giudicelli Date: Tue, 24 Dec 2024 19:07:37 +0100 Subject: [PATCH 12/16] Adapt slip material and use of slip velocity to linear variables --- .../WCNSFV2PSlipVelocityFunctorMaterial.h | 6 ++-- .../linearfvkernels/LinearFVScalarAdvection.h | 10 +++++++ .../WCNSFV2PSlipVelocityFunctorMaterial.C | 18 +++++++----- .../linearfvkernels/LinearFVScalarAdvection.C | 28 ++++++++++++++++++- 4 files changed, 51 insertions(+), 11 deletions(-) diff --git a/modules/navier_stokes/include/functormaterials/WCNSFV2PSlipVelocityFunctorMaterial.h b/modules/navier_stokes/include/functormaterials/WCNSFV2PSlipVelocityFunctorMaterial.h index 0ecb6ac816f7..01ba0dc4ee4a 100644 --- a/modules/navier_stokes/include/functormaterials/WCNSFV2PSlipVelocityFunctorMaterial.h +++ b/modules/navier_stokes/include/functormaterials/WCNSFV2PSlipVelocityFunctorMaterial.h @@ -28,11 +28,11 @@ class WCNSFV2PSlipVelocityFunctorMaterial : public FunctorMaterial const unsigned int _dim; /// x-velocity - const INSFVVelocityVariable * const _u_var; + const MooseVariableField * const _u_var; /// y-velocity - const INSFVVelocityVariable * const _v_var; + const MooseVariableField * const _v_var; /// z-velocity - const INSFVVelocityVariable * const _w_var; + const MooseVariableField * const _w_var; /// Continuous phase density const Moose::Functor & _rho_mixture; diff --git a/modules/navier_stokes/include/linearfvkernels/LinearFVScalarAdvection.h b/modules/navier_stokes/include/linearfvkernels/LinearFVScalarAdvection.h index 3dffd7e69446..ea797f102578 100644 --- a/modules/navier_stokes/include/linearfvkernels/LinearFVScalarAdvection.h +++ b/modules/navier_stokes/include/linearfvkernels/LinearFVScalarAdvection.h @@ -49,6 +49,16 @@ class LinearFVScalarAdvection : public LinearFVFluxKernel /// matrix and right hand side contribution Real _volumetric_face_flux; + /// slip velocity in direction x + const Moose::Functor * const _u_slip; + /// slip velocity in direction y + const Moose::Functor * const _v_slip; + /// slip velocity in direction z + const Moose::Functor * const _w_slip; + + /// Whether to use an additional slip velocity to compute the face flux + bool _add_slip_model; + /// The interpolation method to use for the advected quantity Moose::FV::InterpMethod _advected_interp_method; }; diff --git a/modules/navier_stokes/src/functormaterials/WCNSFV2PSlipVelocityFunctorMaterial.C b/modules/navier_stokes/src/functormaterials/WCNSFV2PSlipVelocityFunctorMaterial.C index 688ee5e1a716..00d6ecbff1dd 100644 --- a/modules/navier_stokes/src/functormaterials/WCNSFV2PSlipVelocityFunctorMaterial.C +++ b/modules/navier_stokes/src/functormaterials/WCNSFV2PSlipVelocityFunctorMaterial.C @@ -9,6 +9,7 @@ #include "WCNSFV2PSlipVelocityFunctorMaterial.h" #include "INSFVVelocityVariable.h" +#include "MooseLinearVariableFV.h" #include "Function.h" #include "NS.h" #include "FVKernel.h" @@ -63,12 +64,12 @@ WCNSFV2PSlipVelocityFunctorMaterial::WCNSFV2PSlipVelocityFunctorMaterial( const InputParameters & params) : FunctorMaterial(params), _dim(_subproblem.mesh().dimension()), - _u_var(dynamic_cast(getFieldVar("u", 0))), + _u_var(dynamic_cast *>(getFieldVar("u", 0))), _v_var(params.isParamValid("v") - ? dynamic_cast(getFieldVar("v", 0)) + ? dynamic_cast *>(getFieldVar("v", 0)) : nullptr), _w_var(params.isParamValid("w") - ? dynamic_cast(getFieldVar("w", 0)) + ? dynamic_cast *>(getFieldVar("w", 0)) : nullptr), _rho_mixture(getFunctor(NS::density)), _rho_d(getFunctor("rho_d")), @@ -82,15 +83,18 @@ WCNSFV2PSlipVelocityFunctorMaterial::WCNSFV2PSlipVelocityFunctorMaterial( _particle_diameter(getFunctor("particle_diameter")), _index(getParam("momentum_component")) { - if (!_u_var) - paramError("u", "the u velocity must be an INSFVVelocityVariable."); + if (!dynamic_cast(_u_var) && + !dynamic_cast *>(_u_var)) + paramError("u", "the u velocity must be an INSFVVelocityVariable or a ."); - if (_dim >= 2 && !_v_var) + if (_dim >= 2 && (!dynamic_cast(_v_var) && + !dynamic_cast *>(_v_var))) paramError("v", "In two or more dimensions, the v velocity must be supplied and it must be an " "INSFVVelocityVariable."); - if (_dim >= 3 && !_w_var) + if (_dim >= 3 && (!dynamic_cast(_w_var) && + !dynamic_cast *>(_w_var))) paramError("w", "In three-dimensions, the w velocity must be supplied and it must be an " "INSFVVelocityVariable."); diff --git a/modules/navier_stokes/src/linearfvkernels/LinearFVScalarAdvection.C b/modules/navier_stokes/src/linearfvkernels/LinearFVScalarAdvection.C index 175c5057934b..bd98e38d4883 100644 --- a/modules/navier_stokes/src/linearfvkernels/LinearFVScalarAdvection.C +++ b/modules/navier_stokes/src/linearfvkernels/LinearFVScalarAdvection.C @@ -24,6 +24,9 @@ LinearFVScalarAdvection::validParams() "rhie_chow_user_object", "The rhie-chow user-object which is used to determine the face velocity."); params += Moose::FV::advectedInterpolationParameter(); + params.addParam("u_slip", "The slip-velocity in the x direction."); + params.addParam("v_slip", "The slip-velocity in the y direction."); + params.addParam("w_slip", "The slip-velocity in the z direction."); return params; } @@ -31,7 +34,11 @@ LinearFVScalarAdvection::LinearFVScalarAdvection(const InputParameters & params) : LinearFVFluxKernel(params), _mass_flux_provider(getUserObject("rhie_chow_user_object")), _advected_interp_coeffs(std::make_pair(0, 0)), - _volumetric_face_flux(0.0) + _volumetric_face_flux(0.0), + _u_slip(isParamValid("u_slip") ? &getFunctor("u_slip") : nullptr), + _v_slip(isParamValid("v_slip") ? &getFunctor("v_slip") : nullptr), + _w_slip(isParamValid("w_slip") ? &getFunctor("w_slip") : nullptr), + _add_slip_model(isParamValid("u_slip") ? true : false) { Moose::FV::setInterpolationMethod(*this, _advected_interp_method, "advected_interp_method"); } @@ -96,6 +103,25 @@ LinearFVScalarAdvection::setupFaceData(const FaceInfo * face_info) // hand side contributions _volumetric_face_flux = _mass_flux_provider.getVolumetricFaceFlux(*face_info); + // Adjust volumetric face flux using the slip velocity + if (_u_slip) + { + const auto state = determineState(); + Moose::FaceArg face_arg; + // TODO Add boundary treatment to be able select two-term expansion if desired + face_arg = Moose::FaceArg{ + face_info, Moose::FV::LimiterType::CentralDifference, true, false, nullptr, nullptr}; + + RealVectorValue velocity_slip_vel_vec; + if (_u_slip) + velocity_slip_vel_vec(0) = (*_u_slip)(face_arg, state).value(); + if (_v_slip) + velocity_slip_vel_vec(1) = (*_v_slip)(face_arg, state).value(); + if (_w_slip) + velocity_slip_vel_vec(2) = (*_w_slip)(face_arg, state).value(); + _volumetric_face_flux += velocity_slip_vel_vec * face_info->normal(); + } + // Caching the interpolation coefficients so they will be reused for the matrix and right hand // side terms _advected_interp_coeffs = From 46a6512cde0f975ab8ecfd206f9884f4ce88a6be Mon Sep 17 00:00:00 2001 From: Guillaume Giudicelli Date: Tue, 24 Dec 2024 22:26:32 +0100 Subject: [PATCH 13/16] Fix gravity kernel creation in Linear FV flow physics --- modules/navier_stokes/src/physics/WCNSLinearFVFlowPhysics.C | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/modules/navier_stokes/src/physics/WCNSLinearFVFlowPhysics.C b/modules/navier_stokes/src/physics/WCNSLinearFVFlowPhysics.C index 55948b7b264e..73e99f434767 100644 --- a/modules/navier_stokes/src/physics/WCNSLinearFVFlowPhysics.C +++ b/modules/navier_stokes/src/physics/WCNSLinearFVFlowPhysics.C @@ -261,9 +261,9 @@ WCNSLinearFVFlowPhysics::addINSMomentumGravityKernels() if (gravity_vector(d) != 0) { params.set("source_density") = std::to_string(gravity_vector(d)); - params.set("variable") = _velocity_names[d]; + params.set("variable") = _velocity_names[d]; - getProblem().addFVKernel(kernel_type, kernel_name + NS::directions[d], params); + getProblem().addLinearFVKernel(kernel_type, kernel_name + NS::directions[d], params); } } } From f4b39779cd9c7ba238218ce88376e680a2fd76ce Mon Sep 17 00:00:00 2001 From: Guillaume Giudicelli Date: Tue, 24 Dec 2024 22:27:35 +0100 Subject: [PATCH 14/16] Add parameter setting routines for two phase mixtures to be shared between linear and nonlinear implementations --- .../physics/WCNSFVTwoPhaseMixturePhysics.h | 2 + .../physics/WCNSFVTwoPhaseMixturePhysics.C | 92 +++++++++++-------- 2 files changed, 56 insertions(+), 38 deletions(-) diff --git a/modules/navier_stokes/include/physics/WCNSFVTwoPhaseMixturePhysics.h b/modules/navier_stokes/include/physics/WCNSFVTwoPhaseMixturePhysics.h index b72e6a96d5c0..d14d5de6718b 100644 --- a/modules/navier_stokes/include/physics/WCNSFVTwoPhaseMixturePhysics.h +++ b/modules/navier_stokes/include/physics/WCNSFVTwoPhaseMixturePhysics.h @@ -22,6 +22,8 @@ class WCNSFVTwoPhaseMixturePhysics final : public WCNSFVScalarTransportPhysics { public: static InputParameters validParams(); + static InputParameters commonMixtureParams(); + static void renamePassiveScalarToMixtureParams(InputParameters & params); WCNSFVTwoPhaseMixturePhysics(const InputParameters & parameters); diff --git a/modules/navier_stokes/src/physics/WCNSFVTwoPhaseMixturePhysics.C b/modules/navier_stokes/src/physics/WCNSFVTwoPhaseMixturePhysics.C index 8e20c7ba6309..e8ad4ffda993 100644 --- a/modules/navier_stokes/src/physics/WCNSFVTwoPhaseMixturePhysics.C +++ b/modules/navier_stokes/src/physics/WCNSFVTwoPhaseMixturePhysics.C @@ -19,22 +19,32 @@ InputParameters WCNSFVTwoPhaseMixturePhysics::validParams() { InputParameters params = WCNSFVScalarTransportPhysics::validParams(); - params.addClassDescription("Define the additional terms for a mixture model for the two phase " - "weakly-compressible Navier Stokes equations"); - // It can be useful to define the mixture materials with a fixed phase fraction instead - // of solving the equations - params.addParam("add_scalar_equation", true, ""); - params.renameParam("add_scalar_equation", - "add_phase_transport_equation", - "Whether to add the phase transport equation."); + // First rename the parameters from passive scalar to mixture + renamePassiveScalarToMixtureParams(params); + params.renameParam("passive_scalar_face_interpolation", + "phase_face_interpolation", + "The numerical scheme to interpolate the phase fraction variable to the " + "face (separate from the advected quantity interpolation)"); + // Then add parameters specific to mixtures // The flow physics is obtained from the scalar transport base class // The fluid heat transfer physics is retrieved even if unspecified params.addParam( "fluid_heat_transfer_physics", "NavierStokesFV", "WCNSFVFluidHeatTransferPhysics generating the fluid energy equation"); + params += commonMixtureParams(); + params.addParamNamesToGroup("fluid_heat_transfer_physics", "Phase change"); + params.addClassDescription("Define the additional terms for a mixture model for the two phase " + "weakly-compressible Navier Stokes equations"); + return params; +} + +InputParameters +WCNSFVTwoPhaseMixturePhysics::commonMixtureParams() +{ + InputParameters params = emptyInputParameters(); params.addParam( "use_external_mixture_properties", @@ -45,13 +55,6 @@ WCNSFVTwoPhaseMixturePhysics::validParams() false, "Whether to output every functor material property defined to Exodus"); - params.renameParam("initial_scalar_variables", - "initial_phase_fraction", - "Initial value of the main phase fraction variable"); - params.renameParam("passive_scalar_diffusivity", - "phase_fraction_diffusivity", - "Functor names for the diffusivities used for the main phase fraction."); - // Phase change parameters params.addParam( NS::alpha_exchange, 0, "Name of the volumetric phase exchange coefficient"); @@ -89,9 +92,6 @@ WCNSFVTwoPhaseMixturePhysics::validParams() "Name of the thermal conductivity functor for phase 1"); // Properties of phase 2 (can be solid, another liquid, or gaseous) - params.renameParam("passive_scalar_names", - "phase_2_fraction_name", - "Name of the second phase fraction variable (can be a dispersed phase)"); params.addRequiredParam("phase_2_density_name", "Name of the density functor for phase 2"); params.addRequiredParam("phase_2_viscosity_name", @@ -108,6 +108,42 @@ WCNSFVTwoPhaseMixturePhysics::validParams() false, "Adds a linear friction term with the dispersed phase drag model"); + // Parameter groups + params.addParamNamesToGroup("phase_1_density_name phase_1_viscosity_name " + "phase_1_specific_heat_name phase_1_thermal_conductivity_name " + "phase_2_density_name phase_2_viscosity_name " + "phase_2_specific_heat_name phase_2_thermal_conductivity_name " + "use_external_mixture_properties", + "Mixture material properties"); + + params.addParamNamesToGroup(NS::alpha_exchange + " add_phase_change_energy_term", "Phase change"); + params.addParamNamesToGroup("add_drift_flux_momentum_terms density_interp_method", + "Drift flux model"); + params.addParamNamesToGroup("add_advection_slip_term", "Advection slip model"); + return params; +} + +void +WCNSFVTwoPhaseMixturePhysics::renamePassiveScalarToMixtureParams(InputParameters & params) +{ + // It can be useful to define the mixture materials with a fixed phase fraction instead + // of solving the equations + params.addParam("add_scalar_equation", true, ""); + params.renameParam("add_scalar_equation", + "add_phase_transport_equation", + "Whether to add the phase transport equation."); + + params.renameParam("initial_scalar_variables", + "initial_phase_fraction", + "Initial value of the main phase fraction variable"); + params.renameParam("passive_scalar_diffusivity", + "phase_fraction_diffusivity", + "Functor names for the diffusivities used for the main phase fraction."); + + params.renameParam("passive_scalar_names", + "phase_2_fraction_name", + "Name of the second phase fraction variable (can be a dispersed phase)"); + // Not applicable currently params.suppressParameter>("passive_scalar_source"); params.suppressParameter>>( @@ -127,10 +163,6 @@ WCNSFVTwoPhaseMixturePhysics::validParams() "phase_advection_interpolation", "The numerical scheme to use for interpolating the phase fraction variable, " "as an advected quantity, to the face."); - params.renameParam("passive_scalar_face_interpolation", - "phase_face_interpolation", - "The numerical scheme to interpolate the phase fraction variable to the " - "face (separate from the advected quantity interpolation)"); params.renameParam( "passive_scalar_two_term_bc_expansion", "phase_two_term_bc_expansion", @@ -142,23 +174,7 @@ WCNSFVTwoPhaseMixturePhysics::validParams() "phase_scaling", "The scaling factor for the phase transport equation"); - // Parameter groups params.renameParameterGroup("Passive scalar control", "Mixture transport control"); - params.addParamNamesToGroup("phase_1_density_name phase_1_viscosity_name " - "phase_1_specific_heat_name phase_1_thermal_conductivity_name " - "phase_2_density_name phase_2_viscosity_name " - "phase_2_specific_heat_name phase_2_thermal_conductivity_name " - "use_external_mixture_properties", - "Mixture material properties"); - - params.addParamNamesToGroup("fluid_heat_transfer_physics " + NS::alpha_exchange + - " add_phase_change_energy_term", - "Phase change"); - params.addParamNamesToGroup("add_drift_flux_momentum_terms density_interp_method", - "Drift flux model"); - params.addParamNamesToGroup("add_advection_slip_term", "Advection slip model"); - - return params; } WCNSFVTwoPhaseMixturePhysics::WCNSFVTwoPhaseMixturePhysics(const InputParameters & parameters) From b438a24c9ed4af2399dab18ec7be1e863a1cd52a Mon Sep 17 00:00:00 2001 From: Guillaume Giudicelli Date: Tue, 24 Dec 2024 22:40:19 +0100 Subject: [PATCH 15/16] Add an input for lid-driven + linear FV Currently a steady solve with only velocity-pressure, not converging. Try transient WIP --- .../segregated/lid-driven-two-phase-physics.i | 235 ++++++++++++++++++ 1 file changed, 235 insertions(+) create mode 100644 modules/navier_stokes/test/tests/finite_volume/two_phase/mixture_model/segregated/lid-driven-two-phase-physics.i diff --git a/modules/navier_stokes/test/tests/finite_volume/two_phase/mixture_model/segregated/lid-driven-two-phase-physics.i b/modules/navier_stokes/test/tests/finite_volume/two_phase/mixture_model/segregated/lid-driven-two-phase-physics.i new file mode 100644 index 000000000000..fde10a450537 --- /dev/null +++ b/modules/navier_stokes/test/tests/finite_volume/two_phase/mixture_model/segregated/lid-driven-two-phase-physics.i @@ -0,0 +1,235 @@ +mu = 1.0 +rho = 1.0e3 +mu_d = 0.3 +rho_d = 1.0 +dp = 0.01 +U_lid = 0.1 +g = -9.81 +advected_interp_method = 'upwind' + +k = 1 +k_d = 1 +cp = 1 +cp_d = 1 + +[Mesh] + [gen] + type = GeneratedMeshGenerator + dim = 2 + xmin = 0 + xmax = .1 + ymin = 0 + ymax = .1 + nx = 10 + ny = 10 + [] +[] + +[Problem] + linear_sys_names = 'u_system v_system pressure_system phi_system' +[] + +[Physics] + [NavierStokes] + [FlowSegregated] + [flow] + compressibility = 'incompressible' + + density = ${rho} #'rho_mixture' + dynamic_viscosity = ${mu} #'mu_mixture' + + # Initial conditions + initial_velocity = '1e-12 1e-12 0' + initial_pressure = 0 + + # Pressure pin + pin_pressure = true + pinned_pressure_type = 'point-value' + pinned_pressure_point = '0.01 0.099 0.0' + pinned_pressure_value = '0' + + # Gravity + gravity = '0 ${g} 0' + + # Boundary conditions are defined outside of the Physics + # Moving walls are not that common of a problem + + momentum_advection_interpolation = '${advected_interp_method}' + [] + [] + # [TwoPhaseMixtureSegregated] + # [mixture] + # system_names = 'phi_system' + # phase_1_fraction_name = 'phase_1' + # phase_2_fraction_name = 'phase_2' + + # add_phase_transport_equation = true + # phase_advection_interpolation = '${advected_interp_method}' + # phase_fraction_diffusivity = 1e-3 + + # # We could consider adding fixed-value-yet-not-an-inlet + # # boundary conditions to the TwoPhaseMixture physics + + # # Base phase material properties + # phase_1_density_name = ${rho} + # phase_1_viscosity_name = ${mu} + # phase_1_specific_heat_name = ${cp} + # phase_1_thermal_conductivity_name = ${k} + + # # Other phase material properties + # phase_2_density_name = ${rho_d} + # phase_2_viscosity_name = ${mu_d} + # phase_2_specific_heat_name = ${cp_d} + # phase_2_thermal_conductivity_name = ${k_d} + # output_all_properties = true + + # # Friction model, not actually used! + # use_dispersed_phase_drag_model = true + # particle_diameter = ${dp} + # add_advection_slip_term = false + # [] + # [] + [] +[] + +[LinearFVBCs] + [moving-lid-x] + type = LinearFVAdvectionDiffusionFunctorDirichletBC + boundary = 'top' + variable = vel_x + functor = '${U_lid}' + [] + [no-slip-wall-u] + type = LinearFVAdvectionDiffusionFunctorDirichletBC + boundary = 'left right bottom' + variable = vel_x + functor = '0' + [] + [no-slip-wall-v] + type = LinearFVAdvectionDiffusionFunctorDirichletBC + boundary = 'left right top bottom' + variable = vel_y + functor = '0' + [] + # [botttom-phase-2] + # type = LinearFVAdvectionDiffusionFunctorDirichletBC + # boundary = 'bottom' + # variable = phase_2 + # functor = '0' + # [] + # [top-phase-2] + # type = LinearFVAdvectionDiffusionFunctorDirichletBC + # boundary = 'top' + # variable = phase_2 + # functor = '0' + # [] +[] + +[AuxVariables] + [drag_coefficient] + type = MooseVariableFVReal + [] +[] + +[AuxKernels] + [populate_cd] + type = FunctorAux + variable = drag_coefficient + functor = 'Darcy_coefficient' + execute_on = 'TIMESTEP_END' + [] +[] + +[AuxVariables] + [phase_2] + [] + [Darcy_coefficient] + [] + [vel_slip_x] + [] + [vel_slip_y] + [] +[] + +[Postprocessors] + [average_void] + type = ElementAverageValue + variable = 'phase_2' + [] + [max_y_velocity] + type = ElementExtremeValue + variable = 'vel_y' + value_type = max + [] + [min_y_velocity] + type = ElementExtremeValue + variable = 'vel_y' + value_type = min + [] + [max_x_velocity] + type = ElementExtremeValue + variable = 'vel_x' + value_type = max + [] + [min_x_velocity] + type = ElementExtremeValue + variable = 'vel_x' + value_type = min + [] + [max_x_slip_velocity] + type = ElementExtremeFunctorValue + functor = 'vel_slip_x' + value_type = max + [] + [max_y_slip_velocity] + type = ElementExtremeFunctorValue + functor = 'vel_slip_y' + value_type = max + [] + [max_drag_coefficient] + type = ElementExtremeFunctorValue + functor = 'drag_coefficient' + value_type = max + [] +[] + +[Executioner] + type = SIMPLE + rhie_chow_user_object = 'ins_rhie_chow_interpolator' + + # Systems + momentum_systems = 'u_system v_system' + pressure_system = 'pressure_system' + # passive_scalar_systems = 'phi_system' + momentum_equation_relaxation = 0.8 + # passive_scalar_equation_relaxation = '0.9' + pressure_variable_relaxation = 0.3 + + # We need to converge the problem to show conservation + num_iterations = 200 + pressure_absolute_tolerance = 1e-10 + momentum_absolute_tolerance = 1e-10 + # passive_scalar_absolute_tolerance = '1e-10' + momentum_petsc_options_iname = '-pc_type -pc_hypre_type' + momentum_petsc_options_value = 'hypre boomeramg' + pressure_petsc_options_iname = '-pc_type -pc_hypre_type' + pressure_petsc_options_value = 'hypre boomeramg' + # passive_scalar_petsc_options_iname = '-pc_type -pc_hypre_type' + # passive_scalar_petsc_options_value = 'hypre boomeramg' + momentum_l_abs_tol = 1e-13 + pressure_l_abs_tol = 1e-13 + # passive_scalar_l_abs_tol = 1e-13 + momentum_l_tol = 0 + pressure_l_tol = 0 + # passive_scalar_l_tol = 0 + print_fields = false + continue_on_max_its = true +[] + +[Outputs] + exodus = false + [CSV] + type = CSV + execute_on = 'FINAL' + [] +[] From c7c56d29c3834c1fac9322824fe8920ee129747f Mon Sep 17 00:00:00 2001 From: Guillaume Giudicelli Date: Wed, 25 Dec 2024 00:11:28 +0100 Subject: [PATCH 16/16] Add support for active scalars in SIMPLE --- .../include/executioners/SIMPLESolve.h | 30 ++++++ .../src/executioners/SIMPLESolve.C | 97 ++++++++++++++++++- 2 files changed, 125 insertions(+), 2 deletions(-) diff --git a/modules/navier_stokes/include/executioners/SIMPLESolve.h b/modules/navier_stokes/include/executioners/SIMPLESolve.h index 3931ef4d898c..8d66a93b5515 100644 --- a/modules/navier_stokes/include/executioners/SIMPLESolve.h +++ b/modules/navier_stokes/include/executioners/SIMPLESolve.h @@ -71,6 +71,36 @@ class SIMPLESolve : public SIMPLESolveBase /// Pointer(s) to the system(s) corresponding to the passive scalar equation(s) std::vector _passive_scalar_systems; + /// Pointer(s) to the system(s) corresponding to the active scalar equation(s) + std::vector _active_scalar_systems; + /// Pointer to the segregated RhieChow interpolation object RhieChowMassFlux * _rc_uo; + + // ************************ Active Scalar Variables ************************ // + + /// The names of the active scalar systems + const std::vector & _active_scalar_system_names; + + /// Boolean for easy check if a active scalar systems shall be solved or not + const bool _has_active_scalar_systems; + + // The number(s) of the system(s) corresponding to the active scalar equation(s) + std::vector _active_scalar_system_numbers; + + /// The user-defined relaxation parameter(s) for the active scalar equation(s) + const std::vector _active_scalar_equation_relaxation; + + /// Options which hold the petsc settings for the active scalar equation(s) + Moose::PetscSupport::PetscOptions _active_scalar_petsc_options; + + /// Options for the linear solver of the active scalar equation(s) + SIMPLESolverConfiguration _active_scalar_linear_control; + + /// Absolute linear tolerance for the active scalar equation(s). We need to store this, because + /// it needs to be scaled with a representative flux. + const Real _active_scalar_l_abs_tol; + + /// The user-defined absolute tolerance for determining the convergence in active scalars + const std::vector _active_scalar_absolute_tolerance; }; diff --git a/modules/navier_stokes/src/executioners/SIMPLESolve.C b/modules/navier_stokes/src/executioners/SIMPLESolve.C index 5c94e83bcda3..871ffb97c3e4 100644 --- a/modules/navier_stokes/src/executioners/SIMPLESolve.C +++ b/modules/navier_stokes/src/executioners/SIMPLESolve.C @@ -18,6 +18,54 @@ InputParameters SIMPLESolve::validParams() { InputParameters params = SIMPLESolveBase::validParams(); + + /* + * Parameters to control the solution of each scalar advection system + */ + params.addParam>("active_scalar_equation_relaxation", + std::vector(), + "The relaxation which should be used for the active scalar " + "equations. (=1 for no relaxation, " + "diagonal dominance will still be enforced)"); + + params.addParam("active_scalar_petsc_options", + Moose::PetscSupport::getCommonPetscFlags(), + "Singleton PETSc options for the active scalar equation(s)"); + params.addParam( + "active_scalar_petsc_options_iname", + Moose::PetscSupport::getCommonPetscKeys(), + "Names of PETSc name/value pairs for the active scalar equation(s)"); + params.addParam>( + "active_scalar_petsc_options_value", + "Values of PETSc name/value pairs (must correspond with \"petsc_options_iname\" for the " + "active scalar equation(s)"); + params.addParam>( + "active_scalar_absolute_tolerance", + std::vector(), + "The absolute tolerance(s) on the normalized residual(s) of the active scalar equation(s)."); + params.addRangeCheckedParam("active_scalar_l_tol", + 1e-5, + "0.0<=active_scalar_l_tol & active_scalar_l_tol<1.0", + "The relative tolerance on the normalized residual in the " + "linear solver of the active scalar equation(s)."); + params.addRangeCheckedParam("active_scalar_l_abs_tol", + 1e-10, + "0.0( + "active_scalar_l_max_its", + 10000, + "The maximum allowed iterations in the linear solver of the turbulence equation."); + + params.addParamNamesToGroup( + "active_scalar_systems active_scalar_equation_relaxation active_scalar_petsc_options " + "active_scalar_petsc_options_iname " + "active_scalar_petsc_options_value active_scalar_petsc_options_value " + "active_scalar_absolute_tolerance " + "active_scalar_l_tol active_scalar_l_abs_tol active_scalar_l_max_its", + "Active Scalars Equations"); + return params; } @@ -28,7 +76,14 @@ SIMPLESolve::SIMPLESolve(Executioner & ex) _energy_sys_number(_has_energy_system ? _problem.linearSysNum(getParam("energy_system")) : libMesh::invalid_uint), - _energy_system(_has_energy_system ? &_problem.getLinearSystem(_energy_sys_number) : nullptr) + _energy_system(_has_energy_system ? &_problem.getLinearSystem(_energy_sys_number) : nullptr), + _active_scalar_system_names(getParam>("active_scalar_systems")), + _has_active_scalar_systems(!_active_scalar_system_names.empty()), + _active_scalar_equation_relaxation( + getParam>("active_scalar_equation_relaxation")), + _active_scalar_l_abs_tol(getParam("active_scalar_l_abs_tol")), + _active_scalar_absolute_tolerance( + getParam>("active_scalar_absolute_tolerance")) { // We fetch the systems and their numbers for the momentum equations. for (auto system_i : index_range(_momentum_system_names)) @@ -45,6 +100,15 @@ SIMPLESolve::SIMPLESolve(Executioner & ex) _passive_scalar_systems.push_back( &_problem.getLinearSystem(_passive_scalar_system_numbers[system_i])); } + // and for the active scalar equations + if (_has_active_scalar_systems) + for (auto system_i : index_range(_active_scalar_system_names)) + { + _active_scalar_system_numbers.push_back( + _problem.linearSysNum(_active_scalar_system_names[system_i])); + _active_scalar_systems.push_back( + &_problem.getLinearSystem(_active_scalar_system_numbers[system_i])); + } } void @@ -284,12 +348,16 @@ SIMPLESolve::solve() unsigned int iteration_counter = 0; // Assign residuals to general residual vector - const unsigned int no_systems = _momentum_systems.size() + 1 + _has_energy_system; + const unsigned int no_systems = + _momentum_systems.size() + 1 + _has_energy_system + _active_scalar_absolute_tolerance.size(); std::vector> ns_residuals(no_systems, std::make_pair(0, 1.0)); std::vector ns_abs_tols(_momentum_systems.size(), _momentum_absolute_tolerance); ns_abs_tols.push_back(_pressure_absolute_tolerance); if (_has_energy_system) ns_abs_tols.push_back(_energy_absolute_tolerance); + if (_has_active_scalar_systems) + for (const auto scalar_tol : _active_scalar_absolute_tolerance) + ns_abs_tols.push_back(scalar_tol); bool converged = false; // Loop until converged or hit the maximum allowed iteration number @@ -356,6 +424,22 @@ SIMPLESolve::solve() _energy_linear_control, _energy_l_abs_tol); } + + // If we have active scalar equations, solve them here in case they depend on temperature + if (_has_active_scalar_systems && converged) + { + // We set the preconditioner/controllable parameters through petsc options. Linear + // tolerances will be overridden within the solver. + Moose::PetscSupport::petscSetOptions(_active_scalar_petsc_options, solver_params); + for (const auto i : index_range(_active_scalar_system_names)) + ns_residuals[momentum_residual.size() + 1 + _has_energy_system + i] = + solveAdvectedSystem(_active_scalar_system_numbers[i], + *_active_scalar_systems[i], + _active_scalar_equation_relaxation[i], + _active_scalar_linear_control, + _active_scalar_l_abs_tol); + } + _problem.execute(EXEC_NONLINEAR); // Printing residuals @@ -377,6 +461,15 @@ SIMPLESolve::solve() << ns_residuals[momentum_residual.size() + 1].second << COLOR_DEFAULT << " Linear its: " << ns_residuals[momentum_residual.size() + 1].first << std::endl; + if (_has_active_scalar_systems) + _console << " Active scalar equations:\n"; + for (const auto i : index_range(_active_scalar_system_names)) + _console << COLOR_GREEN << " " << _active_scalar_system_names[i] << ": " + << ns_residuals[momentum_residual.size() + 1 + _has_energy_system + i].second + << COLOR_DEFAULT << " Linear its: " + << ns_residuals[momentum_residual.size() + 1 + _has_energy_system + i].first + << std::endl; + converged = NS::FV::converged(ns_residuals, ns_abs_tols); }