From 649645001cfeead0451a28c770953520bf423375 Mon Sep 17 00:00:00 2001 From: Guillaume Giudicelli Date: Wed, 27 Nov 2024 07:48:27 -0700 Subject: [PATCH 1/9] Add scalar advection kernels refs #29165 --- .../linearfvkernels/LinearFVScalarAdvection.h | 55 ++++++++++ .../linearfvkernels/LinearFVScalarAdvection.C | 103 ++++++++++++++++++ 2 files changed, 158 insertions(+) create mode 100644 modules/navier_stokes/include/linearfvkernels/LinearFVScalarAdvection.h create mode 100644 modules/navier_stokes/src/linearfvkernels/LinearFVScalarAdvection.C diff --git a/modules/navier_stokes/include/linearfvkernels/LinearFVScalarAdvection.h b/modules/navier_stokes/include/linearfvkernels/LinearFVScalarAdvection.h new file mode 100644 index 000000000000..8b3bedcbbff4 --- /dev/null +++ b/modules/navier_stokes/include/linearfvkernels/LinearFVScalarAdvection.h @@ -0,0 +1,55 @@ +//* 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 "LinearFVFluxKernel.h" +#include "RhieChowMassFlux.h" +#include "LinearFVAdvectionDiffusionBC.h" + +/** + * An advection kernel that implements the advection term for the enthalpy in the + * energy equation. + */ +class LinearFVScalarAdvection : public LinearFVFluxKernel +{ +public: + static InputParameters validParams(); + LinearFVScalarAdvection(const InputParameters & params); + + virtual Real computeElemMatrixContribution() override; + + virtual Real computeNeighborMatrixContribution() override; + + virtual Real computeElemRightHandSideContribution() override; + + virtual Real computeNeighborRightHandSideContribution() override; + + virtual Real computeBoundaryMatrixContribution(const LinearFVBoundaryCondition & bc) override; + + virtual Real computeBoundaryRHSContribution(const LinearFVBoundaryCondition & bc) override; + + virtual void setupFaceData(const FaceInfo * face_info) override; + +protected: + /// The Rhie-Chow user object that provides us with the face velocity + const RhieChowMassFlux & _mass_flux_provider; + +private: + /// Container for the current advected interpolation coefficients on the face to make sure + /// we don't compute it multiple times for different terms. + std::pair _advected_interp_coeffs; + + /// Container for the mass flux on the face which will be reused in the advection term's + /// matrix and right hand side contribution + Real _face_mass_flux; + + /// The interpolation method to use for the advected quantity + Moose::FV::InterpMethod _advected_interp_method; +}; diff --git a/modules/navier_stokes/src/linearfvkernels/LinearFVScalarAdvection.C b/modules/navier_stokes/src/linearfvkernels/LinearFVScalarAdvection.C new file mode 100644 index 000000000000..97c10b9b52b0 --- /dev/null +++ b/modules/navier_stokes/src/linearfvkernels/LinearFVScalarAdvection.C @@ -0,0 +1,103 @@ +//* 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 "LinearFVScalarAdvection.h" +#include "MooseLinearVariableFV.h" +#include "NSFVUtils.h" +#include "NS.h" + +registerMooseObject("NavierStokesApp", LinearFVScalarAdvection); + +InputParameters +LinearFVScalarAdvection::validParams() +{ + InputParameters params = LinearFVFluxKernel::validParams(); + params.addClassDescription("Represents the matrix and right hand side contributions of an " + "advection term for a passive scalar."); + params.addRequiredParam( + "rhie_chow_user_object", + "The rhie-chow user-object which is used to determine the face velocity."); + params += Moose::FV::advectedInterpolationParameter(); + return params; +} + +LinearFVScalarAdvection::LinearFVScalarAdvection(const InputParameters & params) + : LinearFVFluxKernel(params), + _mass_flux_provider(getUserObject("rhie_chow_user_object")), + _advected_interp_coeffs(std::make_pair(0, 0)), + _face_mass_flux(0.0) +{ + Moose::FV::setInterpolationMethod(*this, _advected_interp_method, "advected_interp_method"); +} + +Real +LinearFVScalarAdvection::computeElemMatrixContribution() +{ + return _advected_interp_coeffs.first * _face_mass_flux * _current_face_area; +} + +Real +LinearFVScalarAdvection::computeNeighborMatrixContribution() +{ + return _advected_interp_coeffs.second * _face_mass_flux * _current_face_area; +} + +Real +LinearFVScalarAdvection::computeElemRightHandSideContribution() +{ + return 0.0; +} + +Real +LinearFVScalarAdvection::computeNeighborRightHandSideContribution() +{ + return 0.0; +} + +Real +LinearFVScalarAdvection::computeBoundaryMatrixContribution(const LinearFVBoundaryCondition & bc) +{ + const auto * const adv_bc = static_cast(&bc); + mooseAssert(adv_bc, "This should be a valid BC!"); + + const auto boundary_value_matrix_contrib = adv_bc->computeBoundaryValueMatrixContribution(); + + // We support internal boundaries too so we have to make sure the normal points always outward + const auto factor = (_current_face_type == FaceInfo::VarFaceNeighbors::ELEM) ? 1.0 : -1.0; + + return boundary_value_matrix_contrib * factor * _face_mass_flux * _current_face_area; +} + +Real +LinearFVScalarAdvection::computeBoundaryRHSContribution(const LinearFVBoundaryCondition & bc) +{ + const auto * const adv_bc = static_cast(&bc); + mooseAssert(adv_bc, "This should be a valid BC!"); + + // We support internal boundaries too so we have to make sure the normal points always outward + const auto factor = (_current_face_type == FaceInfo::VarFaceNeighbors::ELEM ? 1.0 : -1.0); + + const auto boundary_value_rhs_contrib = adv_bc->computeBoundaryValueRHSContribution(); + return -boundary_value_rhs_contrib * factor * _face_mass_flux * _current_face_area; +} + +void +LinearFVScalarAdvection::setupFaceData(const FaceInfo * face_info) +{ + LinearFVFluxKernel::setupFaceData(face_info); + + // Caching the mass flux on the face which will be reused in the advection term's matrix and right + // hand side contributions + _face_mass_flux = _mass_flux_provider.getMassFlux(*face_info); + + // Caching the interpolation coefficients so they will be reused for the matrix and right hand + // side terms + _advected_interp_coeffs = + interpCoeffs(_advected_interp_method, *_current_face_info, true, _face_mass_flux); +} From 9c47414c33e829bf6093fab730ad05aea259db46 Mon Sep 17 00:00:00 2001 From: Guillaume Giudicelli Date: Wed, 27 Nov 2024 09:02:45 -0700 Subject: [PATCH 2/9] Move solve attributes for passive scalar to base SIMPLE class refs #29165 --- .../include/executioners/SIMPLESolveBase.h | 33 +++++- .../SIMPLESolveNonlinearAssembly.h | 27 +---- .../src/executioners/SIMPLESolveBase.C | 108 +++++++++++++++++- .../SIMPLESolveNonlinearAssembly.C | 108 ++---------------- 4 files changed, 145 insertions(+), 131 deletions(-) diff --git a/modules/navier_stokes/include/executioners/SIMPLESolveBase.h b/modules/navier_stokes/include/executioners/SIMPLESolveBase.h index 57437afbe3ef..3a1cc7916353 100644 --- a/modules/navier_stokes/include/executioners/SIMPLESolveBase.h +++ b/modules/navier_stokes/include/executioners/SIMPLESolveBase.h @@ -136,10 +136,31 @@ class SIMPLESolveBase : public SolveObject, public UserObjectInterface /// it needs to be scaled with a representative flux. const Real _energy_l_abs_tol; - // ************************ Iteration control **************************** // + // ************************ Passive Scalar Variables ************************ // - /// The maximum number of momentum-pressure iterations - const unsigned int _num_iterations; + /// The names of the passive scalar systems + const std::vector & _passive_scalar_system_names; + + /// Boolean for easy check if a passive scalar systems shall be solved or not + const bool _has_passive_scalar_systems; + + // The number(s) of the system(s) corresponding to the passive scalar equation(s) + std::vector _passive_scalar_system_numbers; + + /// The user-defined relaxation parameter(s) for the passive scalar equation(s) + const std::vector _passive_scalar_equation_relaxation; + + /// Options which hold the petsc settings for the passive scalar equation(s) + Moose::PetscSupport::PetscOptions _passive_scalar_petsc_options; + + /// Options for the linear solver of the passive scalar equation(s) + SIMPLESolverConfiguration _passive_scalar_linear_control; + + /// Absolute linear tolerance for the passive scalar equation(s). We need to store this, because + /// it needs to be scaled with a representative flux. + const Real _passive_scalar_l_abs_tol; + + // ************************ Iteration control **************************** // /// The user-defined absolute tolerance for determining the convergence in momentum const Real _momentum_absolute_tolerance; @@ -150,6 +171,12 @@ class SIMPLESolveBase : public SolveObject, public UserObjectInterface /// The user-defined absolute tolerance for determining the convergence in energy const Real _energy_absolute_tolerance; + /// The user-defined absolute tolerance for determining the convergence in passive scalars + const std::vector _passive_scalar_absolute_tolerance; + + /// The maximum number of momentum-pressure iterations + const unsigned int _num_iterations; + /// If solve should continue if maximum number of iterations is hit const bool _continue_on_max_its; diff --git a/modules/navier_stokes/include/executioners/SIMPLESolveNonlinearAssembly.h b/modules/navier_stokes/include/executioners/SIMPLESolveNonlinearAssembly.h index f02f3115f9a5..bb70e89ea3ba 100644 --- a/modules/navier_stokes/include/executioners/SIMPLESolveNonlinearAssembly.h +++ b/modules/navier_stokes/include/executioners/SIMPLESolveNonlinearAssembly.h @@ -78,9 +78,6 @@ class SIMPLESolveNonlinearAssembly : public SIMPLESolveBase /// Boolean for easy check if a solid energy system shall be solved or not const bool _has_solid_energy_system; - /// Boolean for easy check if a passive scalar systems shall be solved or not - const bool _has_passive_scalar_systems; - /// Boolean for easy check if turbulence systems shall be solved or not const bool _has_turbulence_systems; @@ -110,30 +107,11 @@ class SIMPLESolveNonlinearAssembly : public SIMPLESolveBase /// it needs to be scaled with a representative flux. const Real _solid_energy_l_abs_tol; - // ******************* Passive scalar Eq Variables *********************** // - - /// The names of the passive scalar systems - const std::vector & _passive_scalar_system_names; - - // The number(s) of the system(s) corresponding to the passive scalar equation(s) - std::vector _passive_scalar_system_numbers; + // ********************* Passive Scalar Eq. Variables *********************** // /// Pointer(s) to the system(s) corresponding to the passive scalar equation(s) std::vector _passive_scalar_systems; - /// The user-defined relaxation parameter(s) for the passive scalar equation(s) - const std::vector _passive_scalar_equation_relaxation; - - /// Options which hold the petsc settings for the passive scalar equation(s) - Moose::PetscSupport::PetscOptions _passive_scalar_petsc_options; - - /// Options for the linear solver of the passive scalar equation(s) - SIMPLESolverConfiguration _passive_scalar_linear_control; - - /// Absolute linear tolerance for the passive scalar equation(s). We need to store this, because - /// it needs to be scaled with a representative flux. - const Real _passive_scalar_l_abs_tol; - // ********************* Turbulence Eq Variables ************************* // /// The names of the turbulence scalar systems @@ -166,9 +144,6 @@ class SIMPLESolveNonlinearAssembly : public SIMPLESolveBase /// The user-defined absolute tolerance for determining the convergence in solid energy const Real _solid_energy_absolute_tolerance; - /// The user-defined absolute tolerance for determining the convergence in passive scalars - const std::vector _passive_scalar_absolute_tolerance; - /// The user-defined absolute tolerance for determining the convergence in turbulence equations const std::vector _turbulence_absolute_tolerance; diff --git a/modules/navier_stokes/src/executioners/SIMPLESolveBase.C b/modules/navier_stokes/src/executioners/SIMPLESolveBase.C index 393fd6f39555..9006ef8a7a1c 100644 --- a/modules/navier_stokes/src/executioners/SIMPLESolveBase.C +++ b/modules/navier_stokes/src/executioners/SIMPLESolveBase.C @@ -25,6 +25,8 @@ SIMPLESolveBase::validParams() params.addRequiredParam("pressure_system", "The solver system for the pressure equation."); params.addParam("energy_system", "The solver system for the energy equation."); + params.addParam>( + "passive_scalar_systems", {}, "The solver system for each scalar advection equation."); /* * Parameters to control the solution of the momentum equation @@ -183,6 +185,59 @@ SIMPLESolveBase::validParams() "0>("passive_scalar_equation_relaxation", + std::vector(), + "The relaxation which should be used for the passive scalar " + "equations. (=1 for no relaxation, " + "diagonal dominance will still be enforced)"); + + params.addParam("passive_scalar_petsc_options", + Moose::PetscSupport::getCommonPetscFlags(), + "Singleton PETSc options for the passive scalar equation(s)"); + params.addParam( + "passive_scalar_petsc_options_iname", + Moose::PetscSupport::getCommonPetscKeys(), + "Names of PETSc name/value pairs for the passive scalar equation(s)"); + params.addParam>( + "passive_scalar_petsc_options_value", + "Values of PETSc name/value pairs (must correspond with \"petsc_options_iname\" for the " + "passive scalar equation(s)"); + params.addParam>( + "passive_scalar_absolute_tolerance", + std::vector(), + "The absolute tolerance(s) on the normalized residual(s) of the passive scalar equation(s)."); + params.addRangeCheckedParam("passive_scalar_l_tol", + 1e-5, + "0.0<=passive_scalar_l_tol & passive_scalar_l_tol<1.0", + "The relative tolerance on the normalized residual in the " + "linear solver of the passive scalar equation(s)."); + params.addRangeCheckedParam("passive_scalar_l_abs_tol", + 1e-10, + "0.0( + "passive_scalar_l_max_its", + 10000, + "The maximum allowed iterations in the linear solver of the turbulence equation."); + + params.addParamNamesToGroup( + "passive_scalar_systems passive_scalar_equation_relaxation passive_scalar_petsc_options " + "passive_scalar_petsc_options_iname " + "passive_scalar_petsc_options_value passive_scalar_petsc_options_value " + "passive_scalar_absolute_tolerance " + "passive_scalar_l_tol passive_scalar_l_abs_tol passive_scalar_l_max_its", + "passive_scalar Equation"); + /* * SIMPLE iteration control */ @@ -214,10 +269,17 @@ SIMPLESolveBase::SIMPLESolveBase(Executioner & ex) _has_energy_system(isParamValid("energy_system")), _energy_equation_relaxation(getParam("energy_equation_relaxation")), _energy_l_abs_tol(getParam("energy_l_abs_tol")), - _num_iterations(getParam("num_iterations")), + _passive_scalar_system_names(getParam>("passive_scalar_systems")), + _has_passive_scalar_systems(!_passive_scalar_system_names.empty()), + _passive_scalar_equation_relaxation( + getParam>("passive_scalar_equation_relaxation")), + _passive_scalar_l_abs_tol(getParam("passive_scalar_l_abs_tol")), _momentum_absolute_tolerance(getParam("momentum_absolute_tolerance")), _pressure_absolute_tolerance(getParam("pressure_absolute_tolerance")), _energy_absolute_tolerance(getParam("energy_absolute_tolerance")), + _passive_scalar_absolute_tolerance( + getParam>("passive_scalar_absolute_tolerance")), + _num_iterations(getParam("num_iterations")), _continue_on_max_its(getParam("continue_on_max_its")), _print_fields(getParam("print_fields")) { @@ -274,6 +336,50 @@ SIMPLESolveBase::SIMPLESolveBase(Executioner & ex) "energy_absolute_tolerance", "energy_equation_relaxation"}, false); + + // We check for input errors with regards to the passive scalar equations. At the same time, we + // set up the corresponding system numbers + if (_has_passive_scalar_systems) + { + if (_passive_scalar_system_names.size() != _passive_scalar_equation_relaxation.size()) + paramError("passive_scalar_equation_relaxation", + "The number of equation relaxation parameters does not match the number of " + "passive scalar equations!"); + if (_passive_scalar_system_names.size() != _passive_scalar_absolute_tolerance.size()) + paramError("passive_scalar_absolute_tolerance", + "The number of absolute tolerances does not match the number of " + "passive scalar equations!"); + } + if (_has_passive_scalar_systems) + { + const auto & passive_scalar_petsc_options = + getParam("passive_scalar_petsc_options"); + const auto & passive_scalar_petsc_pair_options = getParam( + "passive_scalar_petsc_options_iname", "passive_scalar_petsc_options_value"); + Moose::PetscSupport::processPetscFlags(passive_scalar_petsc_options, + _passive_scalar_petsc_options); + Moose::PetscSupport::processPetscPairs(passive_scalar_petsc_pair_options, + _problem.mesh().dimension(), + _passive_scalar_petsc_options); + + _passive_scalar_linear_control.real_valued_data["rel_tol"] = + getParam("passive_scalar_l_tol"); + _passive_scalar_linear_control.real_valued_data["abs_tol"] = + getParam("passive_scalar_l_abs_tol"); + _passive_scalar_linear_control.int_valued_data["max_its"] = + getParam("passive_scalar_l_max_its"); + } + else + checkDependentParameterError("passive_scalar_systems", + {"passive_scalar_petsc_options", + "passive_scalar_petsc_options_iname", + "passive_scalar_petsc_options_value", + "passive_scalar_l_tol", + "passive_scalar_l_abs_tol", + "passive_scalar_l_max_its", + "passive_scalar_equation_relaxation", + "passive_scalar_absolute_tolerance"}, + false); } void diff --git a/modules/navier_stokes/src/executioners/SIMPLESolveNonlinearAssembly.C b/modules/navier_stokes/src/executioners/SIMPLESolveNonlinearAssembly.C index c7cb52af6b71..e049196173c0 100644 --- a/modules/navier_stokes/src/executioners/SIMPLESolveNonlinearAssembly.C +++ b/modules/navier_stokes/src/executioners/SIMPLESolveNonlinearAssembly.C @@ -24,23 +24,14 @@ SIMPLESolveNonlinearAssembly::validParams() /* * The names of the different systems in the segregated solver */ - params.addParam("energy_system", "The solver system for the energy equation."); params.addParam("solid_energy_system", "The solver system for the solid energy equation."); - params.addParam>( - "passive_scalar_systems", {}, "The solver system(s) for the passive scalar equation(s)."); params.addParam>( "turbulence_systems", {}, "The solver system(s) for the turbulence equation(s)."); /* * Relaxation parameters for the different system */ - - params.addParam>("passive_scalar_equation_relaxation", - std::vector(), - "The relaxation which should be used for the passive scalar " - "equations. (=1 for no relaxation, " - "diagonal dominance will still be enforced)"); params.addParam>( "turbulence_equation_relaxation", std::vector(), @@ -54,7 +45,7 @@ SIMPLESolveNonlinearAssembly::validParams() "The lower limit imposed on turbulent quantities. The recommended value for robustness " "is 1e-8."); - params.addParamNamesToGroup("energy_equation_relaxation passive_scalar_equation_relaxation " + params.addParamNamesToGroup("energy_equation_relaxation " "turbulence_equation_relaxation", "Relaxation"); @@ -72,18 +63,6 @@ SIMPLESolveNonlinearAssembly::validParams() "Values of PETSc name/value pairs (must correspond with \"petsc_options_iname\" for the " "solid energy equation"); - params.addParam("passive_scalar_petsc_options", - Moose::PetscSupport::getCommonPetscFlags(), - "Singleton PETSc options for the passive scalar equation(s)"); - params.addParam( - "passive_scalar_petsc_options_iname", - Moose::PetscSupport::getCommonPetscKeys(), - "Names of PETSc name/value pairs for the passive scalar equation(s)"); - params.addParam>( - "passive_scalar_petsc_options_value", - "Values of PETSc name/value pairs (must correspond with \"petsc_options_iname\" for the " - "passive scalar equation(s)"); - params.addParam("turbulence_petsc_options", Moose::PetscSupport::getCommonPetscFlags(), "Singleton PETSc options for the turbulence equation(s)"); @@ -97,8 +76,7 @@ SIMPLESolveNonlinearAssembly::validParams() params.addParamNamesToGroup( "solid_energy_petsc_options solid_energy_petsc_options_iname " - "solid_energy_petsc_options_value passive_scalar_petsc_options " - "passive_scalar_petsc_options_iname passive_scalar_petsc_options_value " + "solid_energy_petsc_options_value " "turbulence_petsc_options turbulence_petsc_options_iname turbulence_petsc_options_value", "PETSc Control"); @@ -110,17 +88,12 @@ SIMPLESolveNonlinearAssembly::validParams() 1e-5, "0.0>( - "passive_scalar_absolute_tolerance", - std::vector(), - "The absolute tolerance(s) on the normalized residual(s) of the passive scalar equation(s)."); params.addParam>( "turbulence_absolute_tolerance", std::vector(), "The absolute tolerance(s) on the normalized residual(s) of the turbulence equation(s)."); - params.addParamNamesToGroup("solid_energy_absolute_tolerance passive_scalar_absolute_tolerance " - "turbulence_absolute_tolerance", + params.addParamNamesToGroup("solid_energy_absolute_tolerance turbulence_absolute_tolerance", "Iteration Control"); /* * Linear iteration tolerances for the different equations @@ -140,20 +113,6 @@ SIMPLESolveNonlinearAssembly::validParams() 10000, "0("passive_scalar_l_tol", - 1e-5, - "0.0<=passive_scalar_l_tol & passive_scalar_l_tol<1.0", - "The relative tolerance on the normalized residual in the " - "linear solver of the passive scalar equation(s)."); - params.addRangeCheckedParam("passive_scalar_l_abs_tol", - 1e-10, - "0.0( - "passive_scalar_l_max_its", - 10000, - "The maximum allowed iterations in the linear solver of the turbulence equation."); params.addRangeCheckedParam("turbulence_l_tol", 1e-5, "0.0<=turbulence_l_tol & turbulence_l_tol<1.0", @@ -169,10 +128,10 @@ SIMPLESolveNonlinearAssembly::validParams() 10000, "The maximum allowed iterations in the linear solver of the turbulence equation(s)."); - params.addParamNamesToGroup("solid_energy_l_tol solid_energy_l_abs_tol solid_energy_l_max_its " - "passive_scalar_l_abs_tol passive_scalar_l_max_its turbulence_l_tol " - "turbulence_l_abs_tol turbulence_l_max_its", - "Linear Iteration Control"); + params.addParamNamesToGroup( + "solid_energy_l_tol solid_energy_l_abs_tol solid_energy_l_max_its turbulence_l_tol " + "turbulence_l_abs_tol turbulence_l_max_its", + "Linear Iteration Control"); return params; } @@ -182,8 +141,6 @@ SIMPLESolveNonlinearAssembly::SIMPLESolveNonlinearAssembly(Executioner & ex) _pressure_sys_number(_problem.nlSysNum(getParam("pressure_system"))), _pressure_system(_problem.getNonlinearSystemBase(_pressure_sys_number)), _has_solid_energy_system(_has_energy_system && isParamValid("solid_energy_system")), - _has_passive_scalar_systems( - !getParam>("passive_scalar_systems").empty()), _has_turbulence_systems(!getParam>("turbulence_systems").empty()), _energy_sys_number(_has_energy_system ? _problem.nlSysNum(getParam("energy_system")) @@ -198,17 +155,11 @@ SIMPLESolveNonlinearAssembly::SIMPLESolveNonlinearAssembly(Executioner & ex) ? &_problem.getNonlinearSystemBase(_solid_energy_sys_number) : nullptr), _solid_energy_l_abs_tol(getParam("solid_energy_l_abs_tol")), - _passive_scalar_system_names(getParam>("passive_scalar_systems")), - _passive_scalar_equation_relaxation( - getParam>("passive_scalar_equation_relaxation")), - _passive_scalar_l_abs_tol(getParam("passive_scalar_l_abs_tol")), _turbulence_system_names(getParam>("turbulence_systems")), _turbulence_equation_relaxation(getParam>("turbulence_equation_relaxation")), _turbulence_field_min_limit(getParam>("turbulence_field_min_limit")), _turbulence_l_abs_tol(getParam("turbulence_l_abs_tol")), _solid_energy_absolute_tolerance(getParam("solid_energy_absolute_tolerance")), - _passive_scalar_absolute_tolerance( - getParam>("passive_scalar_absolute_tolerance")), _turbulence_absolute_tolerance(getParam>("turbulence_absolute_tolerance")), _pressure_tag_name(getParam("pressure_gradient_tag")), _pressure_tag_id(_problem.addVectorTag(_pressure_tag_name)) @@ -240,20 +191,6 @@ SIMPLESolveNonlinearAssembly::SIMPLESolveNonlinearAssembly(Executioner & ex) &_problem.getNonlinearSystemBase(_turbulence_system_numbers[system_i])); } - // We check for input errors with regards to the passive scalar equations. At the same time, we - // set up the corresponding system numbers - if (_has_passive_scalar_systems) - { - if (_passive_scalar_system_names.size() != _passive_scalar_equation_relaxation.size()) - paramError("passive_scalar_equation_relaxation", - "The number of equation relaxation parameters does not match the number of " - "passive scalar equations!"); - if (_passive_scalar_system_names.size() != _passive_scalar_absolute_tolerance.size()) - paramError("passive_scalar_absolute_tolerance", - "The number of absolute tolerances does not match the number of " - "passive scalar equations!"); - } - // We check for input errors with regards to the turbulence equations. At the same time, we // set up the corresponding system numbers if (_has_turbulence_systems) @@ -314,37 +251,6 @@ SIMPLESolveNonlinearAssembly::SIMPLESolveNonlinearAssembly(Executioner & ex) false); } - if (_has_passive_scalar_systems) - { - const auto & passive_scalar_petsc_options = - getParam("passive_scalar_petsc_options"); - const auto & passive_scalar_petsc_pair_options = getParam( - "passive_scalar_petsc_options_iname", "passive_scalar_petsc_options_value"); - Moose::PetscSupport::processPetscFlags(passive_scalar_petsc_options, - _passive_scalar_petsc_options); - Moose::PetscSupport::processPetscPairs(passive_scalar_petsc_pair_options, - _problem.mesh().dimension(), - _passive_scalar_petsc_options); - - _passive_scalar_linear_control.real_valued_data["rel_tol"] = - getParam("passive_scalar_l_tol"); - _passive_scalar_linear_control.real_valued_data["abs_tol"] = - getParam("passive_scalar_l_abs_tol"); - _passive_scalar_linear_control.int_valued_data["max_its"] = - getParam("passive_scalar_l_max_its"); - } - else - checkDependentParameterError("passive_scalar_system", - {"passive_scalar_petsc_options", - "passive_scalar_petsc_options_iname", - "passive_scalar_petsc_options_value", - "passive_scalar_l_tol", - "passive_scalar_l_abs_tol", - "passive_scalar_l_max_its", - "passive_scalar_equation_relaxation", - "passive_scalar_absolute_tolerance"}, - false); - if (_has_turbulence_systems) { const auto & turbulence_petsc_options = getParam("turbulence_petsc_options"); From 4007e6a3f92841a3c8116bc324c8f62e067a262e Mon Sep 17 00:00:00 2001 From: Guillaume Giudicelli Date: Wed, 27 Nov 2024 18:55:45 -0700 Subject: [PATCH 3/9] Add support for passive scalar solves in SIMPLE refs #29165 --- .../include/executioners/SIMPLESolve.h | 3 ++ .../src/executioners/SIMPLESolve.C | 53 +++++++++++++++++++ 2 files changed, 56 insertions(+) diff --git a/modules/navier_stokes/include/executioners/SIMPLESolve.h b/modules/navier_stokes/include/executioners/SIMPLESolve.h index 1ab01f0e66d8..a41cef01a9a7 100644 --- a/modules/navier_stokes/include/executioners/SIMPLESolve.h +++ b/modules/navier_stokes/include/executioners/SIMPLESolve.h @@ -68,6 +68,9 @@ class SIMPLESolve : public SIMPLESolveBase /// Pointer to the nonlinear system corresponding to the fluid energy equation LinearSystem * _energy_system; + /// Pointer(s) to the system(s) corresponding to the passive scalar equation(s) + std::vector _passive_scalar_systems; + /// Pointer to the segregated RhieChow interpolation object RhieChowMassFlux * _rc_uo; }; diff --git a/modules/navier_stokes/src/executioners/SIMPLESolve.C b/modules/navier_stokes/src/executioners/SIMPLESolve.C index 9508144bb306..3c806aa916b5 100644 --- a/modules/navier_stokes/src/executioners/SIMPLESolve.C +++ b/modules/navier_stokes/src/executioners/SIMPLESolve.C @@ -34,6 +34,17 @@ SIMPLESolve::SIMPLESolve(Executioner & ex) _momentum_system_numbers.push_back(_problem.linearSysNum(_momentum_system_names[system_i])); _momentum_systems.push_back(&_problem.getLinearSystem(_momentum_system_numbers[system_i])); } + // and for the passive scalar equations + if (_has_passive_scalar_systems) + { + for (auto system_i : index_range(_passive_scalar_system_names)) + { + _passive_scalar_system_numbers.push_back( + _problem.linearSysNum(_passive_scalar_system_names[system_i])); + _passive_scalar_systems.push_back( + &_problem.getLinearSystem(_passive_scalar_system_numbers[system_i])); + } + } } void @@ -365,6 +376,48 @@ SIMPLESolve::solve() converged = NS::FV::converged(ns_residuals, ns_abs_tols); } + // If we have passive scalar equations, solve them here. We assume the material properties in the + // Navier-Stokes equations do not depend on passive scalar, as they are passive, therefore we + // solve outside of the velocity-pressure loop + if (_has_passive_scalar_systems && converged) + { + // The reason why we need more than one iteration is due to the matrix relaxation + // which can be used to stabilize the equations + bool passive_scalar_converged = false; + iteration_counter = 0; + while (iteration_counter < _num_iterations && !passive_scalar_converged) + { + std::vector> scalar_residuals( + _passive_scalar_system_names.size(), std::make_pair(0, 1.0)); + std::vector scalar_abs_tols; + for (const auto scalar_tol : _passive_scalar_absolute_tolerance) + scalar_abs_tols.push_back(scalar_tol); + + // We set the preconditioner/controllable parameters through petsc options. Linear + // tolerances will be overridden within the solver. + Moose::PetscSupport::petscSetOptions(_passive_scalar_petsc_options, solver_params); + for (const auto i : index_range(_passive_scalar_system_names)) + scalar_residuals[i] = solveAdvectedSystem(_passive_scalar_system_numbers[i], + *_passive_scalar_systems[i], + _passive_scalar_equation_relaxation[i], + _passive_scalar_linear_control, + _passive_scalar_l_abs_tol); + + if (_has_passive_scalar_systems) + _console << " Passive scalar equations:\n"; + for (const auto i : index_range(_passive_scalar_system_names)) + _console << COLOR_GREEN << " " << _passive_scalar_system_names[i] << ": " + << scalar_residuals[i].second << COLOR_DEFAULT + << " Linear its: " << scalar_residuals[i].first << std::endl; + + passive_scalar_converged = NS::FV::converged(scalar_residuals, scalar_abs_tols); + iteration_counter++; + } + + // Both flow and scalars must converge + converged = passive_scalar_converged && converged; + } + converged = _continue_on_max_its ? true : converged; return converged; From 975bdb456743a28fb9e56b0a83107b8705c538d5 Mon Sep 17 00:00:00 2001 From: Guillaume Giudicelli Date: Wed, 27 Nov 2024 08:15:47 -0700 Subject: [PATCH 4/9] Add test for 2 scalar advection refs #29165 Add docs --- .../doc/content/source/executioners/SIMPLE.md | 7 + .../LinearFVScalarAdvection.md | 23 ++ .../linear-segregated/2d-scalar/channel.i | 241 ++++++++++++++++++ .../2d-scalar/gold/channel_out.e | Bin 0 -> 95268 bytes .../linear-segregated/2d-scalar/tests | 15 ++ 5 files changed, 286 insertions(+) create mode 100644 modules/navier_stokes/doc/content/source/linearfvkernels/LinearFVScalarAdvection.md create mode 100644 modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/linear-segregated/2d-scalar/channel.i create mode 100644 modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/linear-segregated/2d-scalar/gold/channel_out.e create mode 100644 modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/linear-segregated/2d-scalar/tests diff --git a/modules/navier_stokes/doc/content/source/executioners/SIMPLE.md b/modules/navier_stokes/doc/content/source/executioners/SIMPLE.md index 5b9ba5239db2..14671b225257 100644 --- a/modules/navier_stokes/doc/content/source/executioners/SIMPLE.md +++ b/modules/navier_stokes/doc/content/source/executioners/SIMPLE.md @@ -146,6 +146,13 @@ As a last step, we add the SIMPLE executioner: !listing modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/linear-segregated/2d/2d-velocity-pressure.i block=Executioner +## Passive scalar advection + +The `SIMPLE` executioner can be used to solve coupled problems involving both flow and passive scalar advection. +Advected passive scalars do not affect the flow distribution, and therefore can be solved after the velocity and +pressure fields have been computed using the `SIMPLE` algorithm. +Several systems may be used, for each passive scalar. + !syntax parameters /Executioner/SIMPLE !syntax inputs /Executioner/SIMPLE diff --git a/modules/navier_stokes/doc/content/source/linearfvkernels/LinearFVScalarAdvection.md b/modules/navier_stokes/doc/content/source/linearfvkernels/LinearFVScalarAdvection.md new file mode 100644 index 000000000000..1a0ea66e5d9a --- /dev/null +++ b/modules/navier_stokes/doc/content/source/linearfvkernels/LinearFVScalarAdvection.md @@ -0,0 +1,23 @@ +# LinearFVScalarAdvection + +This kernel adds the contributions of the scalar advection term to the matrix and right hand side of the scalar equation system for the finite volume SIMPLE segregated solver [SIMPLE.md]. + +This term is described by $\nabla \cdot \left(\rho\vec{u} c_p T \right)$ present in the scalar equation conservation for an incompressible/weakly-compressible formulation. + +For FV, the integral of the advection term of scalar $C_i$ over a cell can be expressed as: + +\begin{equation} +\int\limits_{V_C} \nabla \cdot \left(\rho\vec{u} C_i \right) dV \approx \sum\limits_f (\rho \vec{u}\cdot \vec{n})_{RC} C_if |S_f| \, +\end{equation} + +where $C_if$ is the face value of the scalar concentration. An interpolation scheme (e.g. upwind) can be used to compute the face value. This kernel adds the face contribution for each face $f$ to the right hand side and matrix. + +The face mass flux $(\rho \vec{u}\cdot \vec{n})_{RC}$ is provided by the [RhieChowMassFlux.md] object which uses pressure +gradients and the discrete momentum equation to compute face velocities and mass fluxes. +For more information on the expression that is used, see [SIMPLE.md]. + +!syntax parameters /LinearFVKernels/LinearFVScalarAdvection + +!syntax inputs /LinearFVKernels/LinearFVScalarAdvection + +!syntax children /LinearFVKernels/LinearFVScalarAdvection diff --git a/modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/linear-segregated/2d-scalar/channel.i b/modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/linear-segregated/2d-scalar/channel.i new file mode 100644 index 000000000000..b5537906af8a --- /dev/null +++ b/modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/linear-segregated/2d-scalar/channel.i @@ -0,0 +1,241 @@ +mu = 2.6 +rho = 1.0 +advected_interp_method = 'upwind' +k1 = 0.1 +k2 = 0.2 + +[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 s1_system s2_system' + previous_nl_solution_required = true +[] + +[UserObjects] + [rc] + type = RhieChowMassFlux + u = vel_x + v = vel_y + pressure = pressure + rho = ${rho} + p_diffusion_kernel = p_diffusion + [] +[] + +[Variables] + [vel_x] + type = MooseLinearVariableFVReal + initial_condition = 0.5 + solver_sys = u_system + [] + [vel_y] + type = MooseLinearVariableFVReal + solver_sys = v_system + initial_condition = 0.0 + [] + [pressure] + type = MooseLinearVariableFVReal + solver_sys = pressure_system + initial_condition = 0.2 + [] + [scalar1] + type = MooseLinearVariableFVReal + solver_sys = s1_system + initial_condition = 1.1 + [] + [scalar2] + type = MooseLinearVariableFVReal + solver_sys = s2_system + initial_condition = 3 + [] +[] + +[LinearFVKernels] + [u_advection_stress] + type = LinearWCNSFVMomentumFlux + variable = vel_x + advected_interp_method = ${advected_interp_method} + mu = ${mu} + u = vel_x + v = vel_y + momentum_component = 'x' + rhie_chow_user_object = 'rc' + use_nonorthogonal_correction = false + [] + [v_advection_stress] + type = LinearWCNSFVMomentumFlux + variable = vel_y + advected_interp_method = ${advected_interp_method} + mu = ${mu} + u = vel_x + v = vel_y + momentum_component = 'y' + rhie_chow_user_object = 'rc' + use_nonorthogonal_correction = false + [] + [u_pressure] + type = LinearFVMomentumPressure + variable = vel_x + pressure = pressure + momentum_component = 'x' + [] + [v_pressure] + type = LinearFVMomentumPressure + variable = vel_y + pressure = pressure + momentum_component = 'y' + [] + + [p_diffusion] + type = LinearFVAnisotropicDiffusion + variable = pressure + diffusion_tensor = Ainv + use_nonorthogonal_correction = false + [] + [HbyA_divergence] + type = LinearFVDivergence + variable = pressure + face_flux = HbyA + force_boundary_execution = true + [] + + [s1_advection] + type = LinearFVScalarAdvection + variable = scalar1 + advected_interp_method = ${advected_interp_method} + rhie_chow_user_object = 'rc' + [] + [s1_diffusion] + type = LinearFVDiffusion + variable = scalar1 + diffusion_coeff = ${k1} + use_nonorthogonal_correction = false + [] + [s2_diffusion] + type = LinearFVScalarAdvection + variable = scalar2 + advected_interp_method = ${advected_interp_method} + rhie_chow_user_object = 'rc' + [] + [s2_conduction] + type = LinearFVDiffusion + variable = scalar2 + diffusion_coeff = ${k2} + use_nonorthogonal_correction = false + [] +[] + +[LinearFVBCs] + [inlet-u] + type = LinearFVAdvectionDiffusionFunctorDirichletBC + boundary = 'left' + variable = vel_x + functor = '1.1' + [] + [inlet-v] + type = LinearFVAdvectionDiffusionFunctorDirichletBC + boundary = 'left' + variable = vel_y + functor = '0.0' + [] + [walls-u] + type = LinearFVAdvectionDiffusionFunctorDirichletBC + boundary = 'top bottom' + variable = vel_x + functor = 0.0 + [] + [walls-v] + type = LinearFVAdvectionDiffusionFunctorDirichletBC + boundary = 'top bottom' + variable = vel_y + functor = 0.0 + [] + [outlet_p] + type = LinearFVAdvectionDiffusionFunctorDirichletBC + boundary = 'right' + variable = pressure + functor = 1.4 + [] + [outlet_u] + type = LinearFVAdvectionDiffusionOutflowBC + variable = vel_x + use_two_term_expansion = false + boundary = right + [] + [outlet_v] + type = LinearFVAdvectionDiffusionOutflowBC + variable = vel_y + use_two_term_expansion = false + boundary = right + [] + [inlet_wall_scalar1] + type = LinearFVAdvectionDiffusionFunctorDirichletBC + variable = scalar1 + functor = 1 + boundary = 'left' + [] + [outlet_scalar1] + type = LinearFVAdvectionDiffusionOutflowBC + variable = scalar1 + use_two_term_expansion = false + boundary = right + [] + [inlet_wall_scalar2] + type = LinearFVAdvectionDiffusionFunctorDirichletBC + variable = scalar2 + functor = 2 + boundary = 'left' + [] + [outlet_scalar2] + type = LinearFVAdvectionDiffusionOutflowBC + variable = scalar2 + use_two_term_expansion = false + boundary = right + [] +[] + +[Executioner] + type = SIMPLE + 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 + rhie_chow_user_object = 'rc' + momentum_systems = 'u_system v_system' + pressure_system = 'pressure_system' + passive_scalar_systems = 's1_system s2_system' + momentum_equation_relaxation = 0.8 + passive_scalar_equation_relaxation = '0.9 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 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' + print_fields = false + continue_on_max_its = true +[] + +[Outputs] + exodus = true + execute_on = timestep_end + hide = 'pressure vel_x vel_y' +[] diff --git a/modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/linear-segregated/2d-scalar/gold/channel_out.e b/modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/linear-segregated/2d-scalar/gold/channel_out.e new file mode 100644 index 0000000000000000000000000000000000000000..155370847a0da6b81db65a9b206bf5f8a2efb321 GIT binary patch literal 95268 zcmeHwU5q8ib>5v7CGDvo+!b?VePRi{o>U%YhrNUPOaz;ztgDO^W; zu5bHJ-<8>|W3ktF?I?B!_)&%T`@_E7@%p@F`D48IgN_@?2cAN7lzDE?~&{ zvv*rO#Bv|ho8qy1W+#&rVM`*D518F==N8%AE><5tvTVSVAh z{RjR&Cjb6!FK}Yd-?9gR=g0P)AnbrXFB3NDMZP(T>sZ(8p{?W(ncw*Py?tlb>ytmg zQ#O#>ino!Sl(&#R0pO3F1(P}PkZhn~BiVXy?8QC45P|jmP$u?6JKKRO!1wdgfixt($Q`_a|Hs_T^rg(j^$5D)o$>-Z~7non>GFS$=Z;0A z?;|hhkbU3ZBJ(Jo&%Te(POG1O-*@or``Q6GWLJPEj_(mE`46}s(a(M#S@PLE@#BPN zPy8P|LAsN?pk!EhLR!p3}G@L&0B`MSzppL3ro{~Mpy*FaF2im zdSzj0+8V;TkFbuDb_q0zv_^bJc>3{=5Z(#OW19k~?+kEXq%odZP`vpQqAu+K-Krgv zw+V>99mA#HeN4!mcyFNk9G`jLr=MHDbQTE&!uyXh`I9VjRKcT=U;i%901z4)^Rhy(odX>~6C2AP_EL zu@C7$IV1UV;U(51Uv_2YvagMl#d}Fu1pQw2H61eV4PQa#$a`li8~t09m3-cXTWjMR0zLfnaytIDc>d8f z9skd+2x~h2t}3tN@4d$G-`8~fPXb(PO~(pPnVZ3{h)?nT#rGBAH*hscbMe~OUjAQC z|LMtJ)v^9@zud3=xx%mR_`m+;r~mhF{`6n2)gAw>!nf{N_rIEB-QG=nI4Q6h#Bm@@ z8Hcm}p|`aiAI_?aTF1fQ5ZAud>bc#+-2APk|2O^rFbjag-v1AC+c$L=6rjyO+X`r# zp-qE!6WViW!=W8CYrmw;kM=p*(msrf_AJ_|Xn&$DMB6vGo<n#txPA`T&*ORs7wudRSUuG6>{ah<`ngo|b;nsI0| zJB#ZauE%hl$Mra_3%H)Z^$A?Rh-($sFX39l^(3x!Tu;y+})z({x~(GWI0VMrdf_OV;<8iNBS|3 zX_g}`na4C?k*>@my;)AvlWCSCU75!;%aO*+W18hiZ{{&gSfoAk*d{EeY0fmuJ)!b+ zdwoLX>9%8=Fi*E5+lP6i z4otHgd69WcvmE)7c}%mMg^PJi({4w8wQvyz%aKP(2c}t${K`D0S&qERJfne@5S{2T)J&^o7{&y(uw@3`@gO`>-#gf$a}gz z8gKH;`*6J<7yH67T*qpsoC%)ZM0NuDBak$0ZNwT_F=@_odCZBM$eju)b~({n;L!Ql+D{$!0q8^{;Q^faA)kBgzu3&Xayd32H(I?j~ zhgc__vHN>ude-2wu5@x6IL3w zyx3}k7KqL-ZFW{R&!0c{xE1);%F@c&)62`JSI)P-A?kVL-c=bt`~}2Uz0GTGw5=xwg z;`xQm-;(jG-;wd_vVD%+dylw3^6&#KhMQM6o|yqNVi(DQ|5L zt?R*#weq;N^w<-REy~X0WE9QTJLn8b$gHqx+_K#j166{np0KTT|d9KbHZ-`hSXU zl`2aNJ8=pnlQIX`3+_PLZw5nXm6BjiL#GFjw(o9%F?L-0r8_Q{NRkF$waf-_V4ERy zQ6E05LHl+BCPmd8Lj@2n6;NX-Ae}O&A)SOEVr>CYgGpobl;>{UNVsBZ!m`Ey(jd(N z3RkM}P|X3-l})2l60cfMuRnM5=G3a{2Wl5z8gOm#3sUVt61rh6vUUamD@UZ67q1IZ zoXN;m&$66h9KZz?-a>u{Vj}S2e4J^V3SDpiOYG>z=8JH=l$6O-Rwm7>Je9F;Ub*(n z)r~CSC#BIS?Y29BZ1CYqx(lVU>vohwxEH&j(4Np#OOt?Pe=7nUcx?xM<&|IwxQcC@ z&W_`^%lfprJfx!cd<-2_rM*=Avsug=MM-RDGqU5LSJukz>2ouau$AbPfz<0I08NJ? z3d0bt^%K^|my>xE27Qs`N zr2x7w4Mo3&V`k6Um0g%vR>rVL97`BSV7P3;*b2j*5H^K=Cox0d|A1jm;DgF-b$4i?1v|D$)XPn|Hk&K*U3SepU(3nt*LazGk%2B4IPiBwRP2Ed`MgyoJ( zOaf_H+}MSp|S)9S{(g7B8R=UYG){3}CSXns>bdhQ%R?dDj~OnD2n*I;8+G*I&@m=3Q?R zK-B@!C^L7zkyg%Os5&43&fNV*02VtSfJMu???o$9J0Aeb`De)}MD!KGHbz4c1+Zk$ zpa7Z*SOu_5{1YND)RJkYIH_`pj!&y>r0542en1FE(J&m*@&ZHSDl1@Raz<8}GG?ag z@*z}N11yRKQDf2&08R9A24GpkR1gj+hIz4|F`T`+g(~MdorH$}6@-JvFb8lBy`0c6 z4^YokOaw!*%4A3$T98ziiFv9FfEuq9W(A3@X;?m&!aODa0P8g&WEdOkGM=!HHprZQ1m>b17 z_n{5>9jNsr^D1WnmNbkxKsjfW^G7GP^*W}#>27&G`?9G##3?76L9uO9Qknuh04qLp zvT?%8&p&lBs%6DzgfZvMu&fsr&v<^%jZY72!H~^4aoJxi`lob}t1jt%q_NR8BfyEv z{^mMm%Cf(CPXWM5i%90i+S*@2B82GC1hw>#w7LYId=Mf@c%^xE5K z*QYL0r(lUjo3Mov#YsebN;^_ip`lptcfvga9|^0L*|rt}m}tDOLPK{KG7)C%IP`XP zO%qJ3R9OPdlWs16UD;R(j9da$ru)~-frY*Ppy##W(^eL$&~Oa(B=}WJ4x(m7k)^b1 zDnk%D1-^pz?yT@y0$`4l3t(5Ra66hA9DB+O@Mo+=5*obK2^Nv0Ykw!~&}aj{kV))->gr~pzS4UGyY%7T%566?UL0FxoRoSGeY6<{Vm*vgFx zz+S|aeWwCy+#SePZd3rNkcLJDjEMYV28ngxRe)JHVwK+|Dqz-%v&x521;8l?CRuw> zUGS&BK+Rq_V zuL((3`4CHNHN%-tZj*+!NNhddBEUvslL~2=B(ZM{V|-EN4wpoez>@Fu+Gr!O6Y`fQ z*eLs7z{0i2j_(9EHeupZOS=O&^k4y=xFeW)v#`vRwNc&e%V}236#x(6+Nc`UN@C4n zVj(?A!x|;_iU8P1Y?3t%jl?c?-I03|YyVqC8i}1UHB)ep5ogXdR*?~g))~T^MM>;h z91uzOS|s*Rp-6y@#3mKe&`4~}zhh$(tC83#%P>{NUGSYAV>kYXV_91Mkd@fAI3SWU zYLVF50g(V3iA^e`p^?~A*jmo@8i}2<3{wY1VmJPW6WCiwY`(}@i!bT%aR%;|S zsgQF@n(;Y7JemhMq<}8#BPz* z<{M*6iLHm&0$?MtN!Bzp61#+z^I0U;{vTsc*_C-nCHA{%jD3KP4)~NC`fhI@9v>i- z)W6=`TemwnmMq+IaVqLu{xK8|vnwe{m;CVoidJWU6}L5IG0V%Lig(1Y(R!pp8d9xy zDPe<2=7+&G?XyRa8tOiCT9$eAv#OeK5~#PjVTJw!ga#|W-5?ldMdPP4F`8T zzcU7*k!)fpKK{oAsjh`lTIea;a7&|5?fYdBFOh>?R7uv zIXo z%;?8_bh^q?=8WG9XlC?f-2{0Y95jl)QS=oQeKVtf+ZX+tee`3*md@6SKAku3pBcSW z^eJW<=Q(3a>}Ez!Dx~3%&gf@7V1pb-Vl|4sQS=oQeKVtP6#bk^YwqGg(GSxJy&c8; zUTJl0PkC66=^j{djb1wsD9q-raf=u>ik?(RgI4ta=Hm66m!H2D^j$w5_AkSGyPV!E zS-a`JS)<>|efhSM7lHUj;h6QuDgdW=-8XCW5_QWT3LK7UuNPRK^y>^LrG1OtSaAm6dsc9j_2UJrCwVcx|;UBcdyH3 zz3VAS3zQ^eqZn@URNp)=0msGGANy%qmUX8THI->DU=PISL5H^9B~yJ|s6`MdRoD1% z!0{o6tHaRM@aK<-VKY7?71E%`hw0exnTd2~$=b(ta7pKNQZBQ|W`>6bT)c`NH9g~lLh-pUEUawN?GV%-rAc-TX4 zM;;s7!vU4uFzDNUpl;DyW^ZZB5?q2~PZXw}o9`~>7>)u2Gj@B9I*HS~QM~BJ=8I4Q zr3TX~JtcBM#f}|>9XEfTo94sW=|_IV^2}`yCk?timvyv5ci@DvC(=$AhshwP>-p}? z4Wj@iP?A2ZdyGj7D+~uc4~*jV#`P3a0XuQ~9~@vaUH}w=qZx1FvHz5(hG?D5 z&brg#nEh$3z*39N8DdhX00B0Bu%tp7VEcl6g9vu~93=~+OKw?mqN>McjpX535&~f9 zCp8q3tJXtHi!0}FF4@Yt89t)e5&$#ePk~7z=#)hh90kCa)AA!J)2>GVwL9>5Q)qb3 znu&U{$%vshnTF?PY24%)0kCcwP8>U!j~pyZS`dkANvHxaa&}z%j)x&wi9OBSR0TuJ zHm0L6PX#3iLP;4799o1CWv1UFFxmC2WzqNWSFaf z3_zO(Z6|LjFdy5rgI?B>RK(JJ;zR3PFlHe@RC#vG9H}^E4nP{3E-t7^MRC%ot8U&X z?Py0AE{n|pc7kr-8I-{idFBA+`8Y1VYf=hy)w*@%+Q!XW8_(Dq*DsawaxRk|OF?Hi zj)H6y(^lr z{djNe_zrV`y1fWqjK(mG#MQjYeb4v$!@eDReOG=8W5?;G6*mVsaZ*u5z{*L5?R5n+ z@<{m@t6d4FoAcqo4ZHSM=nS%7CY-0rDL|~xkGhPQcfA5Ye4C^X$bD^u2EZkgnF>s% zNtFeLvb9r)tCs0(bbz%yd?^#F+)l_Q5#Q(SV4=?u3|gVeS%BC{>~)4tFKOfo;!9#k zE0*2fcKj{zrI&QxXU^CXAda&hS8KCM)BBJ#OaT^ZZ3d2dSppPo-3yc(P|DY!p$XGS zfWa1bb!UKNljZ<9{mw!tq|HYT7}$tUzS`Sxj0wV`$)y5@dp?f-fG{ktv^j>c+pbqh z2}vEr?w~|mgP!JWvlEZ~vN7AooQAs{r_ivu+bRGvORG6p;-&(?en!XD0niSF%Mw@A zoqlA~Wn^>DrVF9a*Q-`Cn#v8^#@%4#P`49Q`Ep)Tch025G+`V{g?!8cw%zt^u?0K< zl1-R5${o5m2<|QgZO+(mj}cy3ShX%2ZkvxBh@kd&hTsnSc_Cq5Wi@$(jE%91kR_i; zn**F^7Ox^4I9aB!8Np3vnA4B|*~s+)6rcTeZMfsn`xlL^n@F^dIooJn@qD<$j-ojy zz|k3B1sbXK&c|oA6I3zrxpMu=tt;zSbJNlZ6CX`O_(tLb zq4=D70!t9IIzV|)cfN_|shSnTIFaD{iQD!8E6fv{yln4!%={Om}DE#CI!hn z?v^G3w>^Z5gW5`&Z|r&4rW7ENG^3&gEg9Q=I5rC=8U43vZCro;%Ck4FU)#7oivx;f z^}J0n1t=R2S7IAA_}61~yA#8!Um88H-ds03g^KAH=MB0!;9Qju#8 zkUmRca3`A(zJ3*K7j$xDlm zTPtH_enISr#b5zYYvaOz2{bGenS85Pm&ZD_o3o9wED1XZy#~8OUk}*H-j!&`9W=5B zraeTbDr1;5wuYqQ$^eU^U&F23C{a6f-4au_Iluu&?RpcF#w(5C-8Y#}3co7b*S71F zZ~myz%3}>I5=-+&(E!UQ^uW;AZj&FzIS~wVY8*x!%>hPj4?}CO%VSYX#_kmW(hoz3 z8ccykojWQ39IsDgm6sRKUpW8R1@{*oU3$F0BPF+R>t%+*c~J5WN~#O{i+d*sZ(AJqSiZjRq1*HV9CO=f{$J5JuHWmuyA>nsxsEG z+FqpRCRtp~*+yVEX6OaTPyrZy%`UaEH4Hns^oHJ&&7`F9L081^|a|6tq zbR`YV43eoyQd+Vz^kWOnL&H+j_&ilE;3hf9+)XR3a$z_k;m=d$WR90@n0Ww}+S=!# zVX3WsE)2_T?ekz*YRTphwxHN1-jX-u0b)utaBDt{leRQCRb%Y2XF2d~RcM&3NE1LY zQ~@~Adk#}GL*Kbk(|Fbd%Qi8zw2uPJse~$3o@i^Y0&t>fyoz2$|m~c-P!Fd`bua-4bLk-Ykx$@z^ zf>^1_m`~gZy}EWTz%)L|*1l>zI#J6Ne`Li#JmGPpDw&U}@L|3+@-~}Ftt@wxk3A;N zyf3RkQpq$5wXf33k^9on9a+Ovq^v^2iB6YQq>b1_S>BVY4MoH;G2`aAm1sC=H2{6t z$e2dbPc+_FY7{zrl`f_TwlUpX0SqS@?<>(TG2ScQ%SM&hCIMKo8c>CX!tgxK6Sj)P zo@mJiKtntB6iarMMxkjOs*|T{!bDye;8dYut)_8u$V5wam98f=p1s}ap35?zN}JR0 z@a?@-r?J*!w};lBVrQockmtdw7fWlqx+#ttct6F?UPY(WzVWJdf9#S24YK@=`SP_J zH*RiNif+^VZl)OSNyEz<6WU;;U4-vk=!AQttMd|1a}@w*JOL2rEtPG(G6$%f?nllk zRD?8_hoQTLeI)RLl8_^Die1_)<^bU_E)JZ)s2GzT*zl!CdL$3Ez?5A0ea|e#c>bXT z2v<2A;HS60Bq+`dNwv2+Q#UC;k}z+S(KO`QI$Wv=u)5q4ZeM$ezdiQN2Ezeg zVf(!qdvFfm7Pif}HW@_#(cGg@G-qtD-*-FsWVhSf!@!suDm}j-8%=;()@?WRar$AK z)*ELFJK({~x8c$m$*myfty}~sg`aK{(D2r?&uyeMnwCne*zkoCs=2_10Qk_*mX1=^ zj!tOXyS|T>MT(!@ar;)ZqVlJWGNvIXeusm%w(^xz(VT{(0CQ9{r(puH@rQGDqga=7P*x3PTObZi597@?U5f3!n$D-ds9aAV4Wte%R zjIXVtp$UGha&q+OrW1z37ajYM9GoU8_cEK)5M00G!D9llDp3N^PcIE|4kR4Q%v%}b zLkGut$}@1etFY}3xWgx+WzwZei~x1!dgEYEIRX?;?Ink8%&Icq6X zOyx)xS{42t<}?%sduE)1GzpZ`qyhPqv=IPhcj3mXoMaF^W#;zK0O|TGc)0YRd6jXb ztRKOY?b@Y%#=Ocj1C_6Q2#n1!R1S@aOLjr2m}95`Vkh*lWb-Jz`*6I2tNC_nd;kuz zdoFTA`QeKC4#1j97n28s4`R+1f}=UnXsL!00A~Gm!v`(R;_2Z zJ@+E^OI&lJ=yGqkn zC572u?e2nU75SMLIZq*W4CdksTTz1$!#J)#at}cv%mEP7<576?U>LeyyMSKBD)WHs zZg)tZ30uxuGH|icApLoS>qawJRL?pxVA5oj@9PcVHGt7qOcPx+>~?d9_L`4iFkq5pBq(fZNkbJr4=~T(&q53_&o2x#|CSJ`*a|}x!wzlXS8&vV?e?#F?}ea=^_Ur`YuW^$w)B* ztbHibb^6}8-F0B>tLslcyMF!Y8I3{2aN;so6?v!z2kA-$-=YJAgYm)d6W7hm=PB^2 z8uJx5fg~2Wt~`5lmZKqLQ?n2qc(jQ@*EET%jBTdupt%-cB6wAt`k*=R^7GHaXrAG| zG)YW!qloYR3AbH3DyW^LWSd(!DKsrIIaLmDC-w5HQf1g{oH4sbtBt;DA8@vV{vhyS z%`VJrMN3!VL)U{9NR%%5S4EB1J|yl=Qp{DVES*Bj2nZ5XVH>f>8Y8*AGMMT3P#LEu z^L2~Oz_ASvj%0cym^VtzBZny!F#lPoQR#rVT!$8Ik!ID^sKO*$4+3n|XnyYyrxmlT zd}}(AhK(8}71FS3HELAz$C#MJYSd`TGE7y`dIw&OsV_R-^z5f%S?6@Nv=zpmoabc73=r>_LHJIr~XmKSAQVmFaNHL zzja>33->=H?@xVA#-Cc1@t01>__dFVc;UWJ$oo@!GXCU1#$Ws$8UOvSiFo1O-;npG z{*#P9wJhT=eOSh?JuTve`#vx4Pkmd)SHCIaFaLpzzZHvk;r>w4)#OZk$WZGk39SZnQ#AtjQ?1bd;P2OxszX! z@gp1Z_s^Y^@z-TLFMR(s_1x#wbAKU!|99_~@gIL$#;^a0{Qcyw$oP?yGJo^8Wc=!P SWc<2ppCk9)BcFTtf&ULv$sGRx literal 0 HcmV?d00001 diff --git a/modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/linear-segregated/2d-scalar/tests b/modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/linear-segregated/2d-scalar/tests new file mode 100644 index 000000000000..bd92441e1888 --- /dev/null +++ b/modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/linear-segregated/2d-scalar/tests @@ -0,0 +1,15 @@ +[Tests] + issues = '#28819' + design = 'SIMPLE.md LinearFVScalarAdvection.md' + [channel] + requirement = "The system shall be able to solve a " + [simple] + type = 'Exodiff' + input = channel.i + exodiff = 'channel_out.e' + recover = false # we don't support recovery for SIMPLE yet + max_threads = 1 # see libmesh issue #3808 + detail = "scalar advection problem." + [] + [] +[] From 40db0642e67e20218a0d2b77b1145e328bec6e0e Mon Sep 17 00:00:00 2001 From: Guillaume Giudicelli Date: Thu, 28 Nov 2024 18:46:54 -0700 Subject: [PATCH 5/9] Add a base class for scalar transport WCNSFV physics refs #29175 --- .../physics/WCNSFVScalarTransportPhysics.h | 56 +----- .../WCNSFVScalarTransportPhysicsBase.h | 85 ++++++++ .../physics/WCNSFVScalarTransportPhysics.C | 173 +--------------- .../WCNSFVScalarTransportPhysicsBase.C | 189 ++++++++++++++++++ 4 files changed, 290 insertions(+), 213 deletions(-) create mode 100644 modules/navier_stokes/include/physics/WCNSFVScalarTransportPhysicsBase.h create mode 100644 modules/navier_stokes/src/physics/WCNSFVScalarTransportPhysicsBase.C diff --git a/modules/navier_stokes/include/physics/WCNSFVScalarTransportPhysics.h b/modules/navier_stokes/include/physics/WCNSFVScalarTransportPhysics.h index 74066dc3e0fa..6971d3c2bc30 100644 --- a/modules/navier_stokes/include/physics/WCNSFVScalarTransportPhysics.h +++ b/modules/navier_stokes/include/physics/WCNSFVScalarTransportPhysics.h @@ -9,8 +9,7 @@ #pragma once -#include "NavierStokesPhysicsBase.h" -#include "WCNSFVCoupledAdvectionPhysicsHelper.h" +#include "WCNSFVScalarTransportPhysicsBase.h" #define registerWCNSFVScalarTransportBaseTasks(app_name, derived_name) \ registerMooseAction(app_name, derived_name, "add_variable"); \ @@ -20,44 +19,17 @@ /** * Creates all the objects needed to solve the Navier Stokes scalar transport equations + * using the nonlinear finite volume weakly-compressible discretization (WCNSFV) */ -class WCNSFVScalarTransportPhysics : public NavierStokesPhysicsBase, - public WCNSFVCoupledAdvectionPhysicsHelper +class WCNSFVScalarTransportPhysics : public WCNSFVScalarTransportPhysicsBase { public: static InputParameters validParams(); WCNSFVScalarTransportPhysics(const InputParameters & parameters); - /// Get the names of the advected scalar quantity variables - const std::vector & getAdvectedScalarNames() const - { - return _passive_scalar_names; - } - - /// Whether the physics is actually creating the scalar advection equations - bool hasScalarEquations() const { return _has_scalar_equation; } - -protected: - virtual void addFVKernels() override; - virtual void addFVBCs() override; - - /// Names of the passive scalar variables - std::vector _passive_scalar_names; - /// A boolean to help compatibility with the old Modules/NavierStokesFV syntax - /// or to deliberately skip adding the equations (for example for mixtures with a stationary phase) - const bool _has_scalar_equation; - - /// Passive scalar inlet boundary types - MultiMooseEnum _passive_scalar_inlet_types; - /// Functors describing the inlet boundary values. See passive_scalar_inlet_types for what the functors actually represent - std::vector> _passive_scalar_inlet_functors; - private: - void addNonlinearVariables() override; - void addInitialConditions() override; - - unsigned short getNumberAlgebraicGhostingLayersNeeded() const override; + virtual void addNonlinearVariables() override; /** * Functions adding kernels for the incompressible / weakly-compressible scalar transport @@ -65,22 +37,14 @@ class WCNSFVScalarTransportPhysics : public NavierStokesPhysicsBase, * If the material properties are not constant, some of these can be used for * weakly-compressible simulations as well. */ - void addScalarTimeKernels(); - void addScalarDiffusionKernels(); - void addScalarAdvectionKernels(); - virtual void setSlipVelocityParams(InputParameters & /* params */) const {} + virtual void addScalarTimeKernels() override; + virtual void addScalarDiffusionKernels() override; + virtual void addScalarAdvectionKernels() override; /// Equivalent of NSFVAction addScalarCoupledSourceKernels - void addScalarSourceKernels(); + virtual void addScalarSourceKernels() override; /// Functions adding boundary conditions for the incompressible simulation. /// These are used for weakly-compressible simulations as well. - void addScalarInletBC(); - void addScalarWallBC(); - - /// Functors for the passive scalar sources. Indexing is scalar variable index - std::vector _passive_scalar_sources; - /// Functors for the passive scalar (coupled) sources. Inner indexing is scalar variable index - std::vector> _passive_scalar_coupled_sources; - /// Coefficients multiplying for the passive scalar sources. Inner indexing is scalar variable index - std::vector> _passive_scalar_sources_coef; + virtual void addScalarInletBC() override; + virtual void addScalarWallBC() override; }; diff --git a/modules/navier_stokes/include/physics/WCNSFVScalarTransportPhysicsBase.h b/modules/navier_stokes/include/physics/WCNSFVScalarTransportPhysicsBase.h new file mode 100644 index 000000000000..dfbd48c5b3b1 --- /dev/null +++ b/modules/navier_stokes/include/physics/WCNSFVScalarTransportPhysicsBase.h @@ -0,0 +1,85 @@ +//* 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" + +#define registerWCNSFVScalarTransportBaseTasks(app_name, derived_name) \ + 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") + +/** + * Creates all the objects needed to solve the Navier Stokes scalar transport equations + */ +class WCNSFVScalarTransportPhysicsBase : public NavierStokesPhysicsBase, + public WCNSFVCoupledAdvectionPhysicsHelper +{ +public: + static InputParameters validParams(); + + WCNSFVScalarTransportPhysicsBase(const InputParameters & parameters); + + /// Get the names of the advected scalar quantity variables + const std::vector & getAdvectedScalarNames() const + { + return _passive_scalar_names; + } + + /// Whether the physics is actually creating the scalar advection equations + bool hasScalarEquations() const { return _has_scalar_equation; } + +protected: + virtual void addFVKernels() override; + virtual void addFVBCs() override; + virtual void setSlipVelocityParams(InputParameters & /* params */) const {} + + /// Names of the passive scalar variables + std::vector _passive_scalar_names; + /// A boolean to help compatibility with the old Modules/NavierStokesFV syntax + /// or to deliberately skip adding the equations (for example for mixtures with a stationary phase) + const bool _has_scalar_equation; + + /// Passive scalar inlet boundary types + MultiMooseEnum _passive_scalar_inlet_types; + /// Functors describing the inlet boundary values. See passive_scalar_inlet_types for what the functors actually represent + std::vector> _passive_scalar_inlet_functors; + + /// Functors for the passive scalar sources. Indexing is scalar variable index + std::vector _passive_scalar_sources; + /// Functors for the passive scalar (coupled) sources. Inner indexing is scalar variable index + std::vector> _passive_scalar_coupled_sources; + /// Coefficients multiplying for the passive scalar sources. Inner indexing is scalar variable index + std::vector> _passive_scalar_sources_coef; + +private: + virtual void addInitialConditions() override; + + /** + * Functions adding kernels for the incompressible / weakly-compressible scalar transport + * equation + * If the material properties are not constant, some of these can be used for + * weakly-compressible simulations as well. + */ + virtual void addScalarTimeKernels() = 0; + virtual void addScalarDiffusionKernels() = 0; + virtual void addScalarAdvectionKernels() = 0; + /// Equivalent of NSFVAction addScalarCoupledSourceKernels + virtual void addScalarSourceKernels() = 0; + + /// Functions adding boundary conditions for the incompressible simulation. + /// These are used for weakly-compressible simulations as well. + virtual void addScalarInletBC() = 0; + virtual void addScalarWallBC() = 0; + + virtual unsigned short getNumberAlgebraicGhostingLayersNeeded() const override; +}; diff --git a/modules/navier_stokes/src/physics/WCNSFVScalarTransportPhysics.C b/modules/navier_stokes/src/physics/WCNSFVScalarTransportPhysics.C index faa573c393f9..4f6a80932834 100644 --- a/modules/navier_stokes/src/physics/WCNSFVScalarTransportPhysics.C +++ b/modules/navier_stokes/src/physics/WCNSFVScalarTransportPhysics.C @@ -8,8 +8,7 @@ //* https://www.gnu.org/licenses/lgpl-2.1.html #include "WCNSFVScalarTransportPhysics.h" -#include "WCNSFVCoupledAdvectionPhysicsHelper.h" -#include "WCNSFVFlowPhysics.h" +#include "WCNSFVFlowPhysicsBase.h" #include "NSFVBase.h" #include "NS.h" @@ -19,103 +18,19 @@ registerWCNSFVScalarTransportBaseTasks("NavierStokesApp", WCNSFVScalarTransportP InputParameters WCNSFVScalarTransportPhysics::validParams() { - InputParameters params = NavierStokesPhysicsBase::validParams(); - params += WCNSFVCoupledAdvectionPhysicsHelper::validParams(); - params.addClassDescription( - "Define the Navier Stokes weakly-compressible scalar field transport equation(s)"); - - params += NSFVBase::commonScalarFieldAdvectionParams(); - - // TODO Remove the parameter once NavierStokesFV syntax has been removed - params.addParam( - "add_scalar_equation", - "Whether to add the scalar transport equation. This parameter is not necessary if " - "using the Physics syntax"); - - // These parameters are not shared because the NSFVPhysics use functors - params.addParam>>( - "passive_scalar_inlet_function", - std::vector>(), - "Functors for inlet boundaries in the passive scalar equations."); - - // New functor boundary conditions - params.deprecateParam( - "passive_scalar_inlet_function", "passive_scalar_inlet_functors", "01/01/2025"); - - // No need for the duplication - params.addParam>("passive_scalar_source", "Passive scalar sources"); - - // Spatial finite volume discretization scheme - params.transferParam(NSFVBase::validParams(), - "passive_scalar_advection_interpolation"); + InputParameters params = WCNSFVScalarTransportPhysicsBase::validParams(); + params.addClassDescription("Define the Navier Stokes weakly-compressible scalar field transport " + "equation(s) using the nonlinear finite volume discretization"); params.transferParam(NSFVBase::validParams(), "passive_scalar_face_interpolation"); params.transferParam(NSFVBase::validParams(), "passive_scalar_two_term_bc_expansion"); - - // Nonlinear equation solver scaling - params.addRangeCheckedParam>( - "passive_scalar_scaling", - "passive_scalar_scaling > 0.0", - "The scaling factor for the passive scalar field variables."); - - // Parameter groups - params.addParamNamesToGroup("passive_scalar_names initial_scalar_variables", "Variable"); params.addParamNamesToGroup( - "passive_scalar_advection_interpolation passive_scalar_face_interpolation " - "passive_scalar_two_term_bc_expansion passive_scalar_scaling", - "Numerical scheme"); - params.addParamNamesToGroup("passive_scalar_inlet_types passive_scalar_inlet_functors", - "Inlet boundary"); - + "passive_scalar_face_interpolation passive_scalar_two_term_bc_expansion", "Numerical scheme"); return params; } WCNSFVScalarTransportPhysics::WCNSFVScalarTransportPhysics(const InputParameters & parameters) - : NavierStokesPhysicsBase(parameters), - WCNSFVCoupledAdvectionPhysicsHelper(this), - _passive_scalar_names(getParam>("passive_scalar_names")), - _has_scalar_equation(isParamValid("add_scalar_equation") ? getParam("add_scalar_equation") - : !usingNavierStokesFVSyntax()), - _passive_scalar_inlet_types(getParam("passive_scalar_inlet_types")), - _passive_scalar_inlet_functors( - getParam>>("passive_scalar_inlet_functors")), - _passive_scalar_sources(getParam>("passive_scalar_source")), - _passive_scalar_coupled_sources( - getParam>>("passive_scalar_coupled_source")), - _passive_scalar_sources_coef( - getParam>>("passive_scalar_coupled_source_coeff")) + : WCNSFVScalarTransportPhysicsBase(parameters) { - if (_has_scalar_equation) - for (const auto & scalar_name : _passive_scalar_names) - saveSolverVariableName(scalar_name); - - // For compatibility with Modules/NavierStokesFV syntax - if (!_has_scalar_equation) - return; - - // These parameters must be passed for every passive scalar at a time - checkVectorParamsSameLengthIfSet( - "passive_scalar_names", "passive_scalar_diffusivity", true); - checkVectorParamsSameLengthIfSet>( - "passive_scalar_names", "passive_scalar_coupled_source", true); - checkVectorParamsSameLengthIfSet( - "passive_scalar_names", "passive_scalar_source", true); - checkVectorParamsSameLengthIfSet( - "passive_scalar_names", "passive_scalar_scaling", true); - checkVectorParamsSameLengthIfSet( - "passive_scalar_names", "initial_scalar_variables", true); - checkVectorParamsSameLengthIfSet>( - "passive_scalar_names", "passive_scalar_inlet_functors", true); - if (_passive_scalar_inlet_functors.size()) - checkTwoDVectorParamMultiMooseEnumSameLength( - "passive_scalar_inlet_functors", "passive_scalar_inlet_types", false); - - if (_passive_scalar_sources_coef.size()) - checkTwoDVectorParamsSameLength("passive_scalar_coupled_source", - "passive_scalar_coupled_source_coeff"); - - if (_porous_medium_treatment) - _flow_equations_physics->paramError("porous_medium_treatment", - "Porous media scalar advection is currently unimplemented"); } void @@ -152,22 +67,6 @@ WCNSFVScalarTransportPhysics::addNonlinearVariables() } } -void -WCNSFVScalarTransportPhysics::addFVKernels() -{ - // For compatibility with Modules/NavierStokesFV syntax - if (!_has_scalar_equation) - return; - - if (isTransient()) - addScalarTimeKernels(); - - addScalarAdvectionKernels(); - addScalarDiffusionKernels(); - if (_passive_scalar_sources.size() || _passive_scalar_coupled_sources.size()) - addScalarSourceKernels(); -} - void WCNSFVScalarTransportPhysics::addScalarTimeKernels() { @@ -258,19 +157,6 @@ WCNSFVScalarTransportPhysics::addScalarSourceKernels() } } -void -WCNSFVScalarTransportPhysics::addFVBCs() -{ - // For compatibility with Modules/NavierStokesFV syntax - if (!_has_scalar_equation) - return; - - addScalarInletBC(); - // There is typically no wall flux of passive scalars, similarly we rarely know - // their concentrations at the outlet at the beginning of the simulation - // TODO: we will know the outlet values in case of flow reversal. Implement scalar outlet -} - void WCNSFVScalarTransportPhysics::addScalarInletBC() { @@ -359,50 +245,3 @@ WCNSFVScalarTransportPhysics::addScalarInletBC() } } } - -void -WCNSFVScalarTransportPhysics::addInitialConditions() -{ - // For compatibility with Modules/NavierStokesFV syntax - if (!_has_scalar_equation) - return; - if (!_define_variables && parameters().isParamSetByUser("initial_scalar_variables")) - paramError("initial_scalar_variables", - "Scalar variables are defined externally of NavierStokesFV, so should their inital " - "conditions"); - // 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); - - // There are no default initial conditions for passive scalar variables, we however - // must obey the user-defined initial conditions, even if we are restarting - if (parameters().isParamSetByUser("initial_scalar_variables")) - { - for (unsigned int name_i = 0; name_i < _passive_scalar_names.size(); ++name_i) - { - params.set("variable") = _passive_scalar_names[name_i]; - params.set("function") = - getParam>("initial_scalar_variables")[name_i]; - - getProblem().addInitialCondition("FunctionIC", _passive_scalar_names[name_i] + "_ic", params); - } - } -} - -unsigned short -WCNSFVScalarTransportPhysics::getNumberAlgebraicGhostingLayersNeeded() const -{ - unsigned short necessary_layers = getParam("ghost_layers"); - necessary_layers = - std::max(necessary_layers, _flow_equations_physics->getNumberAlgebraicGhostingLayersNeeded()); - if (getParam("passive_scalar_face_interpolation") == "skewness-corrected") - necessary_layers = std::max(necessary_layers, (unsigned short)3); - - return necessary_layers; -} diff --git a/modules/navier_stokes/src/physics/WCNSFVScalarTransportPhysicsBase.C b/modules/navier_stokes/src/physics/WCNSFVScalarTransportPhysicsBase.C new file mode 100644 index 000000000000..fa5d5370ee49 --- /dev/null +++ b/modules/navier_stokes/src/physics/WCNSFVScalarTransportPhysicsBase.C @@ -0,0 +1,189 @@ +//* 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 "WCNSFVScalarTransportPhysicsBase.h" +#include "WCNSFVCoupledAdvectionPhysicsHelper.h" +#include "WCNSFVFlowPhysics.h" +#include "NSFVBase.h" +#include "NS.h" + +InputParameters +WCNSFVScalarTransportPhysicsBase::validParams() +{ + InputParameters params = NavierStokesPhysicsBase::validParams(); + params += WCNSFVCoupledAdvectionPhysicsHelper::validParams(); + params.addClassDescription( + "Define the Navier Stokes weakly-compressible scalar field transport equation(s)"); + + params += NSFVBase::commonScalarFieldAdvectionParams(); + + // TODO Remove the parameter once NavierStokesFV syntax has been removed + params.addParam( + "add_scalar_equation", + "Whether to add the scalar transport equation. This parameter is not necessary if " + "using the Physics syntax"); + + // These parameters are not shared because the NSFVPhysics use functors + params.addParam>>( + "passive_scalar_inlet_function", + std::vector>(), + "Functors for inlet boundaries in the passive scalar equations."); + + // New functor boundary conditions + params.deprecateParam( + "passive_scalar_inlet_function", "passive_scalar_inlet_functors", "01/01/2025"); + + // No need for the duplication + params.addParam>("passive_scalar_source", "Passive scalar sources"); + + // Spatial finite volume discretization scheme + params.transferParam(NSFVBase::validParams(), + "passive_scalar_advection_interpolation"); + + // Nonlinear equation solver scaling + params.addRangeCheckedParam>( + "passive_scalar_scaling", + "passive_scalar_scaling > 0.0", + "The scaling factor for the passive scalar field variables."); + + // Parameter groups + params.addParamNamesToGroup("passive_scalar_names initial_scalar_variables", "Variable"); + params.addParamNamesToGroup("passive_scalar_advection_interpolation passive_scalar_scaling", + "Numerical scheme"); + params.addParamNamesToGroup("passive_scalar_inlet_types passive_scalar_inlet_functors", + "Inlet boundary"); + + return params; +} + +WCNSFVScalarTransportPhysicsBase::WCNSFVScalarTransportPhysicsBase( + const InputParameters & parameters) + : NavierStokesPhysicsBase(parameters), + WCNSFVCoupledAdvectionPhysicsHelper(this), + _passive_scalar_names(getParam>("passive_scalar_names")), + _has_scalar_equation(isParamValid("add_scalar_equation") ? getParam("add_scalar_equation") + : !usingNavierStokesFVSyntax()), + _passive_scalar_inlet_types(getParam("passive_scalar_inlet_types")), + _passive_scalar_inlet_functors( + getParam>>("passive_scalar_inlet_functors")), + _passive_scalar_sources(getParam>("passive_scalar_source")), + _passive_scalar_coupled_sources( + getParam>>("passive_scalar_coupled_source")), + _passive_scalar_sources_coef( + getParam>>("passive_scalar_coupled_source_coeff")) +{ + if (_has_scalar_equation) + for (const auto & scalar_name : _passive_scalar_names) + saveSolverVariableName(scalar_name); + + // For compatibility with Modules/NavierStokesFV syntax + if (!_has_scalar_equation) + return; + + // These parameters must be passed for every passive scalar at a time + checkVectorParamsSameLengthIfSet( + "passive_scalar_names", "passive_scalar_diffusivity", true); + checkVectorParamsSameLengthIfSet>( + "passive_scalar_names", "passive_scalar_coupled_source", true); + checkVectorParamsSameLengthIfSet( + "passive_scalar_names", "passive_scalar_source", true); + checkVectorParamsSameLengthIfSet( + "passive_scalar_names", "passive_scalar_scaling", true); + checkVectorParamsSameLengthIfSet( + "passive_scalar_names", "initial_scalar_variables", true); + checkVectorParamsSameLengthIfSet>( + "passive_scalar_names", "passive_scalar_inlet_functors", true); + if (_passive_scalar_inlet_functors.size()) + checkTwoDVectorParamMultiMooseEnumSameLength( + "passive_scalar_inlet_functors", "passive_scalar_inlet_types", false); + + if (_passive_scalar_sources_coef.size()) + checkTwoDVectorParamsSameLength("passive_scalar_coupled_source", + "passive_scalar_coupled_source_coeff"); + + if (_porous_medium_treatment) + _flow_equations_physics->paramError("porous_medium_treatment", + "Porous media scalar advection is currently unimplemented"); +} + +void +WCNSFVScalarTransportPhysicsBase::addFVKernels() +{ + // For compatibility with Modules/NavierStokesFV syntax + if (!_has_scalar_equation) + return; + + if (isTransient()) + addScalarTimeKernels(); + + addScalarAdvectionKernels(); + addScalarDiffusionKernels(); + if (_passive_scalar_sources.size() || _passive_scalar_coupled_sources.size()) + addScalarSourceKernels(); +} + +void +WCNSFVScalarTransportPhysicsBase::addFVBCs() +{ + // For compatibility with Modules/NavierStokesFV syntax + if (!_has_scalar_equation) + return; + + addScalarInletBC(); + // There is typically no wall flux of passive scalars, similarly we rarely know + // their concentrations at the outlet at the beginning of the simulation + // TODO: we will know the outlet values in case of flow reversal. Implement scalar outlet +} + +void +WCNSFVScalarTransportPhysicsBase::addInitialConditions() +{ + // For compatibility with Modules/NavierStokesFV syntax + if (!_has_scalar_equation) + return; + if (!_define_variables && parameters().isParamSetByUser("initial_scalar_variables")) + paramError("initial_scalar_variables", + "Scalar variables are defined externally of NavierStokesFV, so should their inital " + "conditions"); + // 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); + + // There are no default initial conditions for passive scalar variables, we however + // must obey the user-defined initial conditions, even if we are restarting + if (parameters().isParamSetByUser("initial_scalar_variables")) + { + for (unsigned int name_i = 0; name_i < _passive_scalar_names.size(); ++name_i) + { + params.set("variable") = _passive_scalar_names[name_i]; + params.set("function") = + getParam>("initial_scalar_variables")[name_i]; + + getProblem().addInitialCondition("FunctionIC", _passive_scalar_names[name_i] + "_ic", params); + } + } +} + +unsigned short +WCNSFVScalarTransportPhysicsBase::getNumberAlgebraicGhostingLayersNeeded() const +{ + unsigned short necessary_layers = getParam("ghost_layers"); + necessary_layers = + std::max(necessary_layers, _flow_equations_physics->getNumberAlgebraicGhostingLayersNeeded()); + if (getParam("passive_scalar_face_interpolation") == "skewness-corrected") + necessary_layers = std::max(necessary_layers, (unsigned short)3); + + return necessary_layers; +} From bf692b6adb7be9f915a5ec3b171256eca4086118 Mon Sep 17 00:00:00 2001 From: Guillaume Giudicelli Date: Thu, 28 Nov 2024 19:54:52 -0700 Subject: [PATCH 6/9] Add a physics class for the linear FV discretization of the scalar advection equation refs #29175 --- .../WCNSFVCoupledAdvectionPhysicsHelper.h | 6 +- .../include/physics/WCNSFVFlowPhysics.h | 2 +- .../include/physics/WCNSFVFlowPhysicsBase.h | 2 + .../physics/WCNSFVScalarTransportPhysics.h | 3 +- .../WCNSFVScalarTransportPhysicsBase.h | 1 + .../include/physics/WCNSLinearFVFlowPhysics.h | 4 + .../WCNSLinearFVScalarTransportPhysics.h | 45 ++++ .../navier_stokes/src/base/NavierStokesApp.C | 2 + .../WCNSFVCoupledAdvectionPhysicsHelper.C | 15 +- .../physics/WCNSFVScalarTransportPhysics.C | 11 +- .../WCNSFVScalarTransportPhysicsBase.C | 8 +- .../src/physics/WCNSLinearFVFlowPhysics.C | 4 +- .../WCNSLinearFVScalarTransportPhysics.C | 232 ++++++++++++++++++ 13 files changed, 316 insertions(+), 19 deletions(-) create mode 100644 modules/navier_stokes/include/physics/WCNSLinearFVScalarTransportPhysics.h create mode 100644 modules/navier_stokes/src/physics/WCNSLinearFVScalarTransportPhysics.C diff --git a/modules/navier_stokes/include/physics/WCNSFVCoupledAdvectionPhysicsHelper.h b/modules/navier_stokes/include/physics/WCNSFVCoupledAdvectionPhysicsHelper.h index 07b652a483d3..18936fec8254 100644 --- a/modules/navier_stokes/include/physics/WCNSFVCoupledAdvectionPhysicsHelper.h +++ b/modules/navier_stokes/include/physics/WCNSFVCoupledAdvectionPhysicsHelper.h @@ -12,7 +12,7 @@ #include "PhysicsBase.h" class NavierStokesPhysicsBase; -class WCNSFVFlowPhysics; +class WCNSFVFlowPhysicsBase; class WCNSFVTurbulencePhysics; /** @@ -26,7 +26,7 @@ class WCNSFVCoupledAdvectionPhysicsHelper WCNSFVCoupledAdvectionPhysicsHelper(const NavierStokesPhysicsBase * derived_physics); - const WCNSFVFlowPhysics * getCoupledFlowPhysics() const; + const WCNSFVFlowPhysicsBase * getCoupledFlowPhysics() const; const WCNSFVTurbulencePhysics * getCoupledTurbulencePhysics() const; /// Return the porosity functor name. @@ -40,7 +40,7 @@ class WCNSFVCoupledAdvectionPhysicsHelper /// The Physics class using this helper const NavierStokesPhysicsBase * _advection_physics; /// Flow physics - const WCNSFVFlowPhysics * _flow_equations_physics; + const WCNSFVFlowPhysicsBase * _flow_equations_physics; /// Turbulence const WCNSFVTurbulencePhysics * _turbulence_physics; diff --git a/modules/navier_stokes/include/physics/WCNSFVFlowPhysics.h b/modules/navier_stokes/include/physics/WCNSFVFlowPhysics.h index 4230e424ef5d..5fd29fa45b08 100644 --- a/modules/navier_stokes/include/physics/WCNSFVFlowPhysics.h +++ b/modules/navier_stokes/include/physics/WCNSFVFlowPhysics.h @@ -23,7 +23,7 @@ class WCNSFVFlowPhysics final : public WCNSFVFlowPhysicsBase WCNSFVFlowPhysics(const InputParameters & parameters); /// Get the name of the linear friction coefficient. Returns an empty string if no friction. - MooseFunctorName getLinearFrictionCoefName() const; + virtual MooseFunctorName getLinearFrictionCoefName() const override; /// Return the name of the Rhie Chow user object UserObjectName rhieChowUOName() const override; /// Return the number of algebraic ghosting layers needed diff --git a/modules/navier_stokes/include/physics/WCNSFVFlowPhysicsBase.h b/modules/navier_stokes/include/physics/WCNSFVFlowPhysicsBase.h index a396f7bb459e..3792f633bdb8 100644 --- a/modules/navier_stokes/include/physics/WCNSFVFlowPhysicsBase.h +++ b/modules/navier_stokes/include/physics/WCNSFVFlowPhysicsBase.h @@ -70,6 +70,8 @@ class WCNSFVFlowPhysicsBase : public NavierStokesPhysicsBase const std::vector & getFluxInletDirections() const { return _flux_inlet_directions; } /// Get the inlet flux postprocessor if using a flux inlet const std::vector & getFluxInletPPs() const { return _flux_inlet_pps; } + /// Get the name of the linear friction coefficient. Returns an empty string if no friction. + virtual MooseFunctorName getLinearFrictionCoefName() const = 0; /// Return the name of the Rhie Chow user object virtual UserObjectName rhieChowUOName() const = 0; /// Return the number of algebraic ghosting layers needed diff --git a/modules/navier_stokes/include/physics/WCNSFVScalarTransportPhysics.h b/modules/navier_stokes/include/physics/WCNSFVScalarTransportPhysics.h index 6971d3c2bc30..acc77b7144f6 100644 --- a/modules/navier_stokes/include/physics/WCNSFVScalarTransportPhysics.h +++ b/modules/navier_stokes/include/physics/WCNSFVScalarTransportPhysics.h @@ -46,5 +46,6 @@ class WCNSFVScalarTransportPhysics : public WCNSFVScalarTransportPhysicsBase /// Functions adding boundary conditions for the incompressible simulation. /// These are used for weakly-compressible simulations as well. virtual void addScalarInletBC() override; - virtual void addScalarWallBC() override; + virtual void addScalarWallBC() override{}; + virtual void addScalarOutletBC() override; }; diff --git a/modules/navier_stokes/include/physics/WCNSFVScalarTransportPhysicsBase.h b/modules/navier_stokes/include/physics/WCNSFVScalarTransportPhysicsBase.h index dfbd48c5b3b1..a08fd4f42693 100644 --- a/modules/navier_stokes/include/physics/WCNSFVScalarTransportPhysicsBase.h +++ b/modules/navier_stokes/include/physics/WCNSFVScalarTransportPhysicsBase.h @@ -80,6 +80,7 @@ class WCNSFVScalarTransportPhysicsBase : public NavierStokesPhysicsBase, /// These are used for weakly-compressible simulations as well. virtual void addScalarInletBC() = 0; virtual void addScalarWallBC() = 0; + virtual void addScalarOutletBC() = 0; virtual unsigned short getNumberAlgebraicGhostingLayersNeeded() const override; }; diff --git a/modules/navier_stokes/include/physics/WCNSLinearFVFlowPhysics.h b/modules/navier_stokes/include/physics/WCNSLinearFVFlowPhysics.h index 5c6ac8cea234..66c67fe0589c 100644 --- a/modules/navier_stokes/include/physics/WCNSLinearFVFlowPhysics.h +++ b/modules/navier_stokes/include/physics/WCNSLinearFVFlowPhysics.h @@ -57,6 +57,10 @@ class WCNSLinearFVFlowPhysics final : public WCNSFVFlowPhysicsBase virtual void addRhieChowUserObjects() override; + virtual MooseFunctorName getLinearFrictionCoefName() const override + { + mooseError("Not implemented"); + }; UserObjectName rhieChowUOName() const override; unsigned short getNumberAlgebraicGhostingLayersNeeded() const override; diff --git a/modules/navier_stokes/include/physics/WCNSLinearFVScalarTransportPhysics.h b/modules/navier_stokes/include/physics/WCNSLinearFVScalarTransportPhysics.h new file mode 100644 index 000000000000..5f2b145542c4 --- /dev/null +++ b/modules/navier_stokes/include/physics/WCNSLinearFVScalarTransportPhysics.h @@ -0,0 +1,45 @@ +//* 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 "WCNSFVScalarTransportPhysicsBase.h" + +/** + * Creates all the objects needed to solve the Navier Stokes scalar transport equations + * using the nonlinear finite volume weakly-compressible discretization (WCNSFV) + */ +class WCNSLinearFVScalarTransportPhysics : public WCNSFVScalarTransportPhysicsBase +{ +public: + static InputParameters validParams(); + + WCNSLinearFVScalarTransportPhysics(const InputParameters & parameters); + +private: + // TODO Rename to linear or solver? + virtual void addNonlinearVariables() override; + + /** + * Functions adding kernels for the incompressible / weakly-compressible scalar transport + * equation + * If the material properties are not constant, some of these can be used for + * weakly-compressible simulations as well. + */ + virtual void addScalarTimeKernels() override; + virtual void addScalarDiffusionKernels() override; + virtual void addScalarAdvectionKernels() override; + virtual void addScalarSourceKernels() override; + + /// Functions adding boundary conditions for the incompressible simulation. + /// These are used for weakly-compressible simulations as well. + virtual void addScalarInletBC() override; + virtual void addScalarWallBC() override{}; + virtual void addScalarOutletBC() override; +}; diff --git a/modules/navier_stokes/src/base/NavierStokesApp.C b/modules/navier_stokes/src/base/NavierStokesApp.C index f7388ba05b5a..1b6c7ad20d00 100644 --- a/modules/navier_stokes/src/base/NavierStokesApp.C +++ b/modules/navier_stokes/src/base/NavierStokesApp.C @@ -54,6 +54,8 @@ associateSyntaxInner(Syntax & syntax, ActionFactory & /*action_factory*/) registerSyntax("WCNSLinearFVFlowPhysics", "Physics/NavierStokes/FlowSegregated/*"); registerSyntax("WCNSFVFluidHeatTransferPhysics", "Physics/NavierStokes/FluidHeatTransfer/*"); registerSyntax("WCNSFVScalarTransportPhysics", "Physics/NavierStokes/ScalarTransport/*"); + registerSyntax("WCNSLinearFVScalarTransportPhysics", + "Physics/NavierStokes/ScalarTransportSegregated/*"); registerSyntax("WCNSFVTurbulencePhysics", "Physics/NavierStokes/Turbulence/*"); registerSyntax("PNSFVSolidHeatTransferPhysics", "Physics/NavierStokes/SolidHeatTransfer/*"); registerSyntax("WCNSFVTwoPhaseMixturePhysics", "Physics/NavierStokes/TwoPhaseMixture/*"); diff --git a/modules/navier_stokes/src/physics/WCNSFVCoupledAdvectionPhysicsHelper.C b/modules/navier_stokes/src/physics/WCNSFVCoupledAdvectionPhysicsHelper.C index 6f7ccca4cd62..de095781f1fa 100644 --- a/modules/navier_stokes/src/physics/WCNSFVCoupledAdvectionPhysicsHelper.C +++ b/modules/navier_stokes/src/physics/WCNSFVCoupledAdvectionPhysicsHelper.C @@ -10,7 +10,7 @@ #include "WCNSFVCoupledAdvectionPhysicsHelper.h" #include "INSFVRhieChowInterpolator.h" #include "RelationshipManager.h" -#include "WCNSFVFlowPhysics.h" +#include "WCNSFVFlowPhysicsBase.h" #include "WCNSFVTurbulencePhysics.h" InputParameters @@ -49,17 +49,18 @@ WCNSFVCoupledAdvectionPhysicsHelper::getPorosityFunctorName(bool smoothed) const return _flow_equations_physics->getPorosityFunctorName(smoothed); } -const WCNSFVFlowPhysics * +const WCNSFVFlowPhysicsBase * WCNSFVCoupledAdvectionPhysicsHelper::getCoupledFlowPhysics() const { // User passed it, just use that if (_advection_physics->isParamValid("coupled_flow_physics")) - return _advection_physics->getCoupledPhysics( + return _advection_physics->getCoupledPhysics( _advection_physics->getParam("coupled_flow_physics")); // Look for any physics of the right type, and check the block restriction else { - const auto all_flow_physics = _advection_physics->getCoupledPhysics(); + const auto all_flow_physics = + _advection_physics->getCoupledPhysics(); for (const auto physics : all_flow_physics) if (_advection_physics->checkBlockRestrictionIdentical( physics->name(), physics->blocks(), /*error_if_not_identical=*/false)) @@ -67,9 +68,9 @@ WCNSFVCoupledAdvectionPhysicsHelper::getCoupledFlowPhysics() const return physics; } } - mooseError( - "No coupled flow Physics found of type 'WCNSFVFlowPhysics'. Use the 'coupled_flow_physics' " - "parameter to give the name of the desired WCNSFVFlowPhysics to couple with"); + mooseError("No coupled flow Physics found of type derived from 'WCNSFVFlowPhysicsBase'. Use the " + "'coupled_flow_physics' parameter to give the name of the desired " + "WCNSFVFlowPhysicsBase-derived Physics to couple with"); } const WCNSFVTurbulencePhysics * diff --git a/modules/navier_stokes/src/physics/WCNSFVScalarTransportPhysics.C b/modules/navier_stokes/src/physics/WCNSFVScalarTransportPhysics.C index 4f6a80932834..92742d12e10b 100644 --- a/modules/navier_stokes/src/physics/WCNSFVScalarTransportPhysics.C +++ b/modules/navier_stokes/src/physics/WCNSFVScalarTransportPhysics.C @@ -22,9 +22,7 @@ WCNSFVScalarTransportPhysics::validParams() params.addClassDescription("Define the Navier Stokes weakly-compressible scalar field transport " "equation(s) using the nonlinear finite volume discretization"); params.transferParam(NSFVBase::validParams(), "passive_scalar_face_interpolation"); - params.transferParam(NSFVBase::validParams(), "passive_scalar_two_term_bc_expansion"); - params.addParamNamesToGroup( - "passive_scalar_face_interpolation passive_scalar_two_term_bc_expansion", "Numerical scheme"); + params.addParamNamesToGroup("passive_scalar_face_interpolation", "Numerical scheme"); return params; } @@ -245,3 +243,10 @@ WCNSFVScalarTransportPhysics::addScalarInletBC() } } } + +void +WCNSFVScalarTransportPhysics::addScalarOutletBC() +{ + // Advection outlet is naturally handled by the advection flux kernel + return; +} diff --git a/modules/navier_stokes/src/physics/WCNSFVScalarTransportPhysicsBase.C b/modules/navier_stokes/src/physics/WCNSFVScalarTransportPhysicsBase.C index fa5d5370ee49..2926a8b58072 100644 --- a/modules/navier_stokes/src/physics/WCNSFVScalarTransportPhysicsBase.C +++ b/modules/navier_stokes/src/physics/WCNSFVScalarTransportPhysicsBase.C @@ -45,6 +45,7 @@ WCNSFVScalarTransportPhysicsBase::validParams() // Spatial finite volume discretization scheme params.transferParam(NSFVBase::validParams(), "passive_scalar_advection_interpolation"); + params.transferParam(NSFVBase::validParams(), "passive_scalar_two_term_bc_expansion"); // Nonlinear equation solver scaling params.addRangeCheckedParam>( @@ -54,7 +55,8 @@ WCNSFVScalarTransportPhysicsBase::validParams() // Parameter groups params.addParamNamesToGroup("passive_scalar_names initial_scalar_variables", "Variable"); - params.addParamNamesToGroup("passive_scalar_advection_interpolation passive_scalar_scaling", + params.addParamNamesToGroup("passive_scalar_advection_interpolation passive_scalar_scaling " + "passive_scalar_two_term_bc_expansion", "Numerical scheme"); params.addParamNamesToGroup("passive_scalar_inlet_types passive_scalar_inlet_functors", "Inlet boundary"); @@ -139,6 +141,8 @@ WCNSFVScalarTransportPhysicsBase::addFVBCs() // There is typically no wall flux of passive scalars, similarly we rarely know // their concentrations at the outlet at the beginning of the simulation // TODO: we will know the outlet values in case of flow reversal. Implement scalar outlet + addScalarWallBC(); + addScalarOutletBC(); } void @@ -182,7 +186,7 @@ WCNSFVScalarTransportPhysicsBase::getNumberAlgebraicGhostingLayersNeeded() const unsigned short necessary_layers = getParam("ghost_layers"); necessary_layers = std::max(necessary_layers, _flow_equations_physics->getNumberAlgebraicGhostingLayersNeeded()); - if (getParam("passive_scalar_face_interpolation") == "skewness-corrected") + if (getParam("passive_scalar_advection_interpolation") == "skewness-corrected") necessary_layers = std::max(necessary_layers, (unsigned short)3); return necessary_layers; diff --git a/modules/navier_stokes/src/physics/WCNSLinearFVFlowPhysics.C b/modules/navier_stokes/src/physics/WCNSLinearFVFlowPhysics.C index 555d583905fe..e8b18c78e74b 100644 --- a/modules/navier_stokes/src/physics/WCNSLinearFVFlowPhysics.C +++ b/modules/navier_stokes/src/physics/WCNSLinearFVFlowPhysics.C @@ -279,8 +279,8 @@ WCNSLinearFVFlowPhysics::addINSInletBC() { // Check the size of the BC parameters unsigned int num_velocity_functor_inlets = 0; - for (const auto & [bdy, momentum_outlet_type] : _momentum_inlet_types) - if (momentum_outlet_type == "fixed-velocity" || momentum_outlet_type == "fixed-pressure") + for (const auto & [bdy, momentum_inlet_type] : _momentum_inlet_types) + if (momentum_inlet_type == "fixed-velocity" || momentum_inlet_type == "fixed-pressure") num_velocity_functor_inlets++; if (num_velocity_functor_inlets != _momentum_inlet_functors.size()) diff --git a/modules/navier_stokes/src/physics/WCNSLinearFVScalarTransportPhysics.C b/modules/navier_stokes/src/physics/WCNSLinearFVScalarTransportPhysics.C new file mode 100644 index 000000000000..bb23899c66de --- /dev/null +++ b/modules/navier_stokes/src/physics/WCNSLinearFVScalarTransportPhysics.C @@ -0,0 +1,232 @@ +//* 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 "WCNSLinearFVScalarTransportPhysics.h" +#include "WCNSFVFlowPhysicsBase.h" +#include "NS.h" + +registerNavierStokesPhysicsBaseTasks("NavierStokesApp", WCNSLinearFVScalarTransportPhysics); +registerWCNSFVScalarTransportBaseTasks("NavierStokesApp", WCNSLinearFVScalarTransportPhysics); + +InputParameters +WCNSLinearFVScalarTransportPhysics::validParams() +{ + InputParameters params = WCNSFVScalarTransportPhysicsBase::validParams(); + params.addClassDescription("Define the Navier Stokes weakly-compressible scalar field transport " + "equation(s) using the linear finite volume discretization"); + params.addParam("use_nonorthogonal_correction", + true, + "If the nonorthogonal correction should be used when computing the normal " + "gradient, notably in the diffusion term."); + params.addParam( + "scalar_boundary_two_term_expansion", + false, + "Whether to use the two term boundary expansion for boundary conditions that allow it"); + return params; +} + +WCNSLinearFVScalarTransportPhysics::WCNSLinearFVScalarTransportPhysics( + const InputParameters & parameters) + : WCNSFVScalarTransportPhysicsBase(parameters) +{ +} + +void +WCNSLinearFVScalarTransportPhysics::addNonlinearVariables() +{ + // For compatibility with Modules/NavierStokesFV syntax + if (!_has_scalar_equation) + return; + + auto params = getFactory().getValidParams("MooseLinearVariableFVReal"); + assignBlocks(params, _blocks); + + for (const auto name_i : index_range(_passive_scalar_names)) + { + // Dont add if the user already defined the variable + if (variableExists(_passive_scalar_names[name_i], /*error_if_aux=*/true)) + { + checkBlockRestrictionIdentical( + _passive_scalar_names[name_i], + getProblem().getVariable(0, _passive_scalar_names[name_i]).blocks()); + continue; + } + + params.set("solver_sys") = getSolverSystem(name_i); + if (isParamValid("passive_scalar_scaling")) + params.set>("scaling") = { + getParam>("passive_scalar_scaling")[name_i]}; + + getProblem().addVariable("MooseLinearVariableFVReal", _passive_scalar_names[name_i], params); + } +} + +void +WCNSLinearFVScalarTransportPhysics::addScalarTimeKernels() +{ + paramError("transient", "Transient simulations are not supported at this time"); +} + +void +WCNSLinearFVScalarTransportPhysics::addScalarAdvectionKernels() +{ + const std::string kernel_type = "LinearFVScalarAdvection"; + InputParameters params = getFactory().getValidParams(kernel_type); + + assignBlocks(params, _blocks); + params.set("rhie_chow_user_object") = _flow_equations_physics->rhieChowUOName(); + params.set("advected_interp_method") = + getParam("passive_scalar_advection_interpolation"); + setSlipVelocityParams(params); + + for (const auto & vname : _passive_scalar_names) + { + params.set("variable") = vname; + getProblem().addLinearFVKernel(kernel_type, prefix() + "ins_" + vname + "_advection", params); + } +} + +void +WCNSLinearFVScalarTransportPhysics::addScalarDiffusionKernels() +{ + // Direct specification of diffusion term + const auto passive_scalar_diffusivities = + getParam>("passive_scalar_diffusivity"); + + if (passive_scalar_diffusivities.size()) + { + const std::string kernel_type = "LinearFVDiffusion"; + InputParameters params = getFactory().getValidParams(kernel_type); + assignBlocks(params, _blocks); + params.set("use_nonorthogonal_correction") = + getParam("use_nonorthogonal_correction"); + for (const auto name_i : index_range(_passive_scalar_names)) + { + params.set("variable") = _passive_scalar_names[name_i]; + params.set("diffusion_coeff") = passive_scalar_diffusivities[name_i]; + getProblem().addLinearFVKernel( + kernel_type, prefix() + "ins_" + _passive_scalar_names[name_i] + "_diffusion", params); + } + } +} + +void +WCNSLinearFVScalarTransportPhysics::addScalarSourceKernels() +{ + const std::string kernel_type = "LinearFVSource"; + InputParameters params = getFactory().getValidParams(kernel_type); + assignBlocks(params, _blocks); + + for (const auto scalar_i : index_range(_passive_scalar_names)) + { + params.set("variable") = _passive_scalar_names[scalar_i]; + + if (_passive_scalar_sources.size()) + { + // Added for backward compatibility with former Modules/NavierStokesFV syntax + params.set("source_density") = _passive_scalar_sources[scalar_i]; + getProblem().addFVKernel( + kernel_type, prefix() + "ins_" + _passive_scalar_names[scalar_i] + "_source", params); + } + + if (_passive_scalar_coupled_sources.size()) + for (const auto i : index_range(_passive_scalar_coupled_sources[scalar_i])) + { + params.set("source_density") = + _passive_scalar_coupled_sources[scalar_i][i]; + if (_passive_scalar_sources_coef.size()) + params.set("scaling_factor") = _passive_scalar_sources_coef[scalar_i][i]; + + getProblem().addLinearFVKernel(kernel_type, + prefix() + "ins_" + _passive_scalar_names[scalar_i] + + "_coupled_source_" + std::to_string(i), + params); + } + } +} + +void +WCNSLinearFVScalarTransportPhysics::addScalarInletBC() +{ + const auto & inlet_boundaries = _flow_equations_physics->getInletBoundaries(); + if (inlet_boundaries.empty()) + return; + + // Boundary checks + // TODO: once we have vectors of MooseEnum, we could use the same templated check for types and + // functors + if (inlet_boundaries.size() * _passive_scalar_names.size() != _passive_scalar_inlet_types.size()) + paramError( + "passive_scalar_inlet_types", + "The number of scalar inlet types (" + std::to_string(_passive_scalar_inlet_types.size()) + + ") is not equal to the number of inlet boundaries (" + + std::to_string(inlet_boundaries.size()) + ") times the number of passive scalars (" + + std::to_string(_passive_scalar_names.size()) + ")"); + if (_passive_scalar_names.size() != _passive_scalar_inlet_functors.size()) + paramError("passive_scalar_inlet_functors", + "The number of groups of inlet functors (" + + std::to_string(_passive_scalar_inlet_functors.size()) + + ") is not equal to the number of passive scalars (" + + std::to_string(_passive_scalar_names.size()) + ")"); + + for (const auto name_i : index_range(_passive_scalar_names)) + { + if (inlet_boundaries.size() != _passive_scalar_inlet_functors[name_i].size()) + paramError("passive_scalar_inlet_functors", + "The number of inlet boundary functors for scalar '" + + _passive_scalar_names[name_i] + + "' does not match the number of inlet boundaries (" + + std::to_string(_passive_scalar_inlet_functors[name_i].size()) + ")"); + + unsigned int num_inlets = inlet_boundaries.size(); + for (unsigned int bc_ind = 0; bc_ind < num_inlets; ++bc_ind) + { + if (_passive_scalar_inlet_types[name_i * num_inlets + bc_ind] == "fixed-value") + { + const std::string bc_type = "LinearFVAdvectionDiffusionFunctorDirichletBC"; + InputParameters params = getFactory().getValidParams(bc_type); + params.set("variable") = _passive_scalar_names[name_i]; + params.set("functor") = _passive_scalar_inlet_functors[name_i][bc_ind]; + params.set>("boundary") = {inlet_boundaries[bc_ind]}; + + getProblem().addLinearFVBC( + bc_type, _passive_scalar_names[name_i] + "_" + inlet_boundaries[bc_ind], params); + } + else if (_passive_scalar_inlet_types[name_i * num_inlets + bc_ind] == "flux-mass" || + _passive_scalar_inlet_types[name_i * num_inlets + bc_ind] == "flux-velocity") + { + mooseError("Flux boundary conditions not supported at this time using the linear finite " + "volume discretization"); + } + } + } +} + +void +WCNSLinearFVScalarTransportPhysics::addScalarOutletBC() +{ + 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("scalar_boundary_two_term_expansion"); + + for (const auto name_i : index_range(_passive_scalar_names)) + { + params.set("variable") = _passive_scalar_names[name_i]; + getProblem().addLinearFVBC(bc_type, _passive_scalar_names[name_i] + "_" + outlet_bdy, params); + } + } +} From 6dda0c846539659aac8754f29be62f0dedb5231f Mon Sep 17 00:00:00 2001 From: Guillaume Giudicelli Date: Thu, 28 Nov 2024 20:41:11 -0700 Subject: [PATCH 7/9] Add test for scalar advection physics add doco refs #29175 --- .../WCNSLinearFVScalarTransportPhysics.md | 39 +++++++ .../ScalarTransportSegregated/index.md | 11 ++ .../WCNSLinearFVScalarTransportPhysics.C | 6 +- .../2d-scalar/channel-physics.i | 104 ++++++++++++++++++ .../2d-scalar/gold/channel-physics_out.e | 1 + .../linear-segregated/2d-scalar/tests | 12 +- 6 files changed, 166 insertions(+), 7 deletions(-) create mode 100644 modules/navier_stokes/doc/content/source/physics/WCNSLinearFVScalarTransportPhysics.md create mode 100644 modules/navier_stokes/doc/content/syntax/Physics/NavierStokes/ScalarTransportSegregated/index.md create mode 100644 modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/linear-segregated/2d-scalar/channel-physics.i create mode 120000 modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/linear-segregated/2d-scalar/gold/channel-physics_out.e diff --git a/modules/navier_stokes/doc/content/source/physics/WCNSLinearFVScalarTransportPhysics.md b/modules/navier_stokes/doc/content/source/physics/WCNSLinearFVScalarTransportPhysics.md new file mode 100644 index 000000000000..58eca580a530 --- /dev/null +++ b/modules/navier_stokes/doc/content/source/physics/WCNSLinearFVScalarTransportPhysics.md @@ -0,0 +1,39 @@ +# Navier Stokes Linear Scalar Transport / WCNSLinearFVScalarTransportPhysics + +!syntax description /Physics/NavierStokes/ScalarTransportSegregated/WCNSLinearFVScalarTransportPhysics + +## Equation + +This [Physics](Physics/index.md) object creates the kernels and boundary conditions to solve the advection-diffusion-reaction +equation for several scalar quantities advected by the flow. + +!equation +\dfrac{\partial \phi_i}{\partial t} + \nabla \cdot (\phi_i \mathbf{v}) - \nabla \cdot (k_i \nabla \phi_i) - Q_i - \lambda_i \phi_i = 0 + +where: + +- $\phi_i$ is the i-th scalar quantity +- \mathbf{v} is the advecting velocity +- $k_i$ the i-th scalar diffusivity +- $Q_i$ is the i-th precursor source +- $\lambda_i$ is a reaction coefficient. It should be negative for a loss term + +The kernels created are: + +- [LinearFVScalarAdvection.md] for the scalar advection term +- [LinearFVDiffusion.md] for the scalar diffusion term +- [LinearFVSource.md] for the reaction terms +- [LinearFVSource.md] for the source terms + +## Coupling with other Physics + +Scalar advection equations can be solved concurrently with the flow equations by combining the [WCNSLinearFVFlowPhysics.md] with the `WCNSLinearFVScalarTransportPhysics`. +The following input performs this coupling for incompressible flow in a 2D channel. + +!listing test/tests/finite_volume/ins/channel-flow/linear-segregated/2d-scalar/channel-physics.i block=Physics + +!syntax parameters /Physics/NavierStokes/ScalarTransportSegregated/WCNSLinearFVScalarTransportPhysics + +!syntax inputs /Physics/NavierStokes/ScalarTransportSegregated/WCNSLinearFVScalarTransportPhysics + +!syntax children /Physics/NavierStokes/ScalarTransportSegregated/WCNSLinearFVScalarTransportPhysics diff --git a/modules/navier_stokes/doc/content/syntax/Physics/NavierStokes/ScalarTransportSegregated/index.md b/modules/navier_stokes/doc/content/syntax/Physics/NavierStokes/ScalarTransportSegregated/index.md new file mode 100644 index 000000000000..d431642f40ac --- /dev/null +++ b/modules/navier_stokes/doc/content/syntax/Physics/NavierStokes/ScalarTransportSegregated/index.md @@ -0,0 +1,11 @@ +# Navier Stokes Scalar Transport Linear Physics + +This syntax was created for the [WCNSLinearFVScalarTransportPhysics.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/ScalarTransportSegregated objects=True actions=False subsystems=False + +!syntax list /Physics/NavierStokes/ScalarTransportSegregated objects=False actions=False subsystems=True + +!syntax list /Physics/NavierStokes/ScalarTransportSegregated objects=False actions=True subsystems=False diff --git a/modules/navier_stokes/src/physics/WCNSLinearFVScalarTransportPhysics.C b/modules/navier_stokes/src/physics/WCNSLinearFVScalarTransportPhysics.C index bb23899c66de..6feb578c69f4 100644 --- a/modules/navier_stokes/src/physics/WCNSLinearFVScalarTransportPhysics.C +++ b/modules/navier_stokes/src/physics/WCNSLinearFVScalarTransportPhysics.C @@ -24,10 +24,6 @@ WCNSLinearFVScalarTransportPhysics::validParams() true, "If the nonorthogonal correction should be used when computing the normal " "gradient, notably in the diffusion term."); - params.addParam( - "scalar_boundary_two_term_expansion", - false, - "Whether to use the two term boundary expansion for boundary conditions that allow it"); return params; } @@ -221,7 +217,7 @@ WCNSLinearFVScalarTransportPhysics::addScalarOutletBC() InputParameters params = getFactory().getValidParams(bc_type); params.set>("boundary") = {outlet_bdy}; params.set("use_two_term_expansion") = - getParam("scalar_boundary_two_term_expansion"); + getParam("passive_scalar_two_term_bc_expansion"); for (const auto name_i : index_range(_passive_scalar_names)) { diff --git a/modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/linear-segregated/2d-scalar/channel-physics.i b/modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/linear-segregated/2d-scalar/channel-physics.i new file mode 100644 index 000000000000..11f7f2117f3e --- /dev/null +++ b/modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/linear-segregated/2d-scalar/channel-physics.i @@ -0,0 +1,104 @@ +mu = 2.6 +rho = 1.0 +advected_interp_method = 'upwind' +k1 = 0.1 +k2 = 0.2 + +[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 s1_system s2_system' +[] + + +[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} + + # use inlet for moving wall to match the reference input + # we could also use a noslip BC with a velocity wall functor + 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_advection_interpolation = ${advected_interp_method} + [] + [ScalarTransportSegregated/scalar] + passive_scalar_names = 'scalar1 scalar2' + system_names = 's1_system s2_system' + initial_scalar_variables = '1.1 3' + + passive_scalar_diffusivity = '${k1} ${k2}' + + passive_scalar_inlet_types = 'fixed-value fixed-value' + passive_scalar_inlet_function = '1; 2' + + passive_scalar_advection_interpolation = ${advected_interp_method} + passive_scalar_two_term_bc_expansion = false + use_nonorthogonal_correction = false + [] + [] +[] + +[Executioner] + type = SIMPLE + 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 + rhie_chow_user_object = 'ins_rhie_chow_interpolator' + momentum_systems = 'u_system v_system' + pressure_system = 'pressure_system' + passive_scalar_systems = 's1_system s2_system' + momentum_equation_relaxation = 0.8 + passive_scalar_equation_relaxation = '0.9 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 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' + print_fields = false + continue_on_max_its = true +[] + +[Outputs] + exodus = true + execute_on = timestep_end + hide = 'pressure vel_x vel_y' +[] diff --git a/modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/linear-segregated/2d-scalar/gold/channel-physics_out.e b/modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/linear-segregated/2d-scalar/gold/channel-physics_out.e new file mode 120000 index 000000000000..ffa5923ca413 --- /dev/null +++ b/modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/linear-segregated/2d-scalar/gold/channel-physics_out.e @@ -0,0 +1 @@ +channel_out.e \ No newline at end of file diff --git a/modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/linear-segregated/2d-scalar/tests b/modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/linear-segregated/2d-scalar/tests index bd92441e1888..aeef8237c5b4 100644 --- a/modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/linear-segregated/2d-scalar/tests +++ b/modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/linear-segregated/2d-scalar/tests @@ -2,14 +2,22 @@ issues = '#28819' design = 'SIMPLE.md LinearFVScalarAdvection.md' [channel] - requirement = "The system shall be able to solve a " + requirement = "The system shall be able to solve a scalar advection problem," [simple] type = 'Exodiff' input = channel.i exodiff = 'channel_out.e' recover = false # we don't support recovery for SIMPLE yet max_threads = 1 # see libmesh issue #3808 - detail = "scalar advection problem." + detail = "using the kernel-based syntax," + [] + [physics] + type = 'Exodiff' + input = channel-physics.i + exodiff = 'channel-physics_out.e' + recover = false # we don't support recovery for SIMPLE yet + max_threads = 1 # see libmesh issue #3808 + detail = "using the shorthand Physics syntax." [] [] [] From 9f5cafe412ef77e0d9e29078a84cf73360aa742e Mon Sep 17 00:00:00 2001 From: Guillaume Giudicelli Date: Sun, 1 Dec 2024 16:47:45 -0700 Subject: [PATCH 8/9] Address Peter's review: - use face velocities to advect scalars - use solverVariables for physics action on add_variable task Co-authored-by: Peter German <31662443+grmnptr@users.noreply.github.com> --- .../include/actioncomponents/ActionComponent.h | 2 +- framework/include/physics/DiffusionCG.h | 2 +- framework/include/physics/DiffusionFV.h | 2 +- framework/include/physics/PhysicsBase.h | 2 +- framework/src/actioncomponents/ActionComponent.C | 2 +- framework/src/physics/DiffusionCG.C | 2 +- framework/src/physics/DiffusionFV.C | 2 +- framework/src/physics/PhysicsBase.C | 2 +- .../include/physics/HeatConductionCG.h | 2 +- .../include/physics/HeatConductionFV.h | 2 +- .../heat_transfer/src/physics/HeatConductionCG.C | 2 +- .../heat_transfer/src/physics/HeatConductionFV.C | 2 +- .../linearfvkernels/LinearFVScalarAdvection.md | 8 ++++---- .../WCNSLinearFVScalarTransportPhysics.md | 6 +++++- .../linearfvkernels/LinearFVScalarAdvection.h | 7 +++---- .../physics/PNSFVSolidHeatTransferPhysics.h | 2 +- .../include/physics/WCNSFVFlowPhysics.h | 2 +- .../include/physics/WCNSFVFlowPhysicsBase.h | 2 +- .../physics/WCNSFVFluidHeatTransferPhysics.h | 2 +- .../physics/WCNSFVScalarTransportPhysics.h | 2 +- .../include/physics/WCNSFVTurbulencePhysics.h | 2 +- .../include/physics/WCNSLinearFVFlowPhysics.h | 2 +- .../physics/WCNSLinearFVScalarTransportPhysics.h | 5 ++--- .../include/userobjects/RhieChowMassFlux.h | 4 +++- .../navier_stokes/src/executioners/SIMPLESolve.C | 4 +--- .../linearfvkernels/LinearFVScalarAdvection.C | 16 ++++++++-------- .../src/physics/PNSFVSolidHeatTransferPhysics.C | 2 +- .../src/physics/WCNSFVFlowPhysics.C | 2 +- .../src/physics/WCNSFVFluidHeatTransferPhysics.C | 2 +- .../src/physics/WCNSFVScalarTransportPhysics.C | 2 +- .../src/physics/WCNSFVTurbulencePhysics.C | 2 +- .../src/physics/WCNSLinearFVFlowPhysics.C | 2 +- .../physics/WCNSLinearFVScalarTransportPhysics.C | 2 +- .../src/userobjects/RhieChowMassFlux.C | 14 +++++++++++++- .../include/physics/MultiSpeciesDiffusionCG.h | 2 +- .../src/physics/MultiSpeciesDiffusionCG.C | 2 +- 36 files changed, 67 insertions(+), 53 deletions(-) diff --git a/framework/include/actioncomponents/ActionComponent.h b/framework/include/actioncomponents/ActionComponent.h index 220192805a32..54681354752b 100644 --- a/framework/include/actioncomponents/ActionComponent.h +++ b/framework/include/actioncomponents/ActionComponent.h @@ -61,7 +61,7 @@ class ActionComponent : public Action virtual void setupComponent() {} // These routines can help define a component that also defines a Physics - virtual void addNonlinearVariables() {} + virtual void addSolverVariables() {} /// Used to add one or more Physics to be active on the component. /// We recommend using the PhysicsComponentInterface instead of overriding this directly diff --git a/framework/include/physics/DiffusionCG.h b/framework/include/physics/DiffusionCG.h index 1a0e9a61f242..fa8021fe5ee2 100644 --- a/framework/include/physics/DiffusionCG.h +++ b/framework/include/physics/DiffusionCG.h @@ -23,7 +23,7 @@ class DiffusionCG : public DiffusionPhysicsBase DiffusionCG(const InputParameters & parameters); private: - virtual void addNonlinearVariables() override; + virtual void addSolverVariables() override; virtual void addFEKernels() override; virtual void addFEBCs() override; diff --git a/framework/include/physics/DiffusionFV.h b/framework/include/physics/DiffusionFV.h index 52eb21920388..8e0dbbff8d50 100644 --- a/framework/include/physics/DiffusionFV.h +++ b/framework/include/physics/DiffusionFV.h @@ -23,7 +23,7 @@ class DiffusionFV : public DiffusionPhysicsBase DiffusionFV(const InputParameters & parameters); private: - virtual void addNonlinearVariables() override; + virtual void addSolverVariables() override; virtual void addFVKernels() override; virtual void addFVBCs() override; virtual void initializePhysicsAdditional() override; diff --git a/framework/include/physics/PhysicsBase.h b/framework/include/physics/PhysicsBase.h index 93213f0bd617..27cd76b9e2c0 100644 --- a/framework/include/physics/PhysicsBase.h +++ b/framework/include/physics/PhysicsBase.h @@ -198,7 +198,7 @@ class PhysicsBase : public Action, public InputParametersChecksUtils _advected_interp_coeffs; - /// Container for the mass flux on the face which will be reused in the advection term's + /// Container for the velocity on the face which will be reused in the advection term's /// matrix and right hand side contribution - Real _face_mass_flux; + Real _face_velocity; /// The interpolation method to use for the advected quantity Moose::FV::InterpMethod _advected_interp_method; diff --git a/modules/navier_stokes/include/physics/PNSFVSolidHeatTransferPhysics.h b/modules/navier_stokes/include/physics/PNSFVSolidHeatTransferPhysics.h index 7d0150ae88d6..3bf4a9faaf89 100644 --- a/modules/navier_stokes/include/physics/PNSFVSolidHeatTransferPhysics.h +++ b/modules/navier_stokes/include/physics/PNSFVSolidHeatTransferPhysics.h @@ -24,7 +24,7 @@ class PNSFVSolidHeatTransferPhysics final : public HeatConductionFV protected: private: - virtual void addNonlinearVariables() override; + virtual void addSolverVariables() override; virtual void addFVKernels() override; virtual void addMaterials() override; diff --git a/modules/navier_stokes/include/physics/WCNSFVFlowPhysics.h b/modules/navier_stokes/include/physics/WCNSFVFlowPhysics.h index 5fd29fa45b08..f693e95167b1 100644 --- a/modules/navier_stokes/include/physics/WCNSFVFlowPhysics.h +++ b/modules/navier_stokes/include/physics/WCNSFVFlowPhysics.h @@ -30,7 +30,7 @@ class WCNSFVFlowPhysics final : public WCNSFVFlowPhysicsBase unsigned short getNumberAlgebraicGhostingLayersNeeded() const override; private: - virtual void addNonlinearVariables() override; + virtual void addSolverVariables() override; virtual void addFVKernels() override; virtual void addUserObjects() override; virtual void addCorrectors() override; diff --git a/modules/navier_stokes/include/physics/WCNSFVFlowPhysicsBase.h b/modules/navier_stokes/include/physics/WCNSFVFlowPhysicsBase.h index 3792f633bdb8..24a677701c8a 100644 --- a/modules/navier_stokes/include/physics/WCNSFVFlowPhysicsBase.h +++ b/modules/navier_stokes/include/physics/WCNSFVFlowPhysicsBase.h @@ -80,7 +80,7 @@ class WCNSFVFlowPhysicsBase : public NavierStokesPhysicsBase protected: virtual void initializePhysicsAdditional() override; virtual void actOnAdditionalTasks() override; - virtual void addNonlinearVariables() override = 0; + virtual void addSolverVariables() override = 0; virtual void addInitialConditions() override; virtual void addFVKernels() override = 0; virtual void addFVBCs() override; diff --git a/modules/navier_stokes/include/physics/WCNSFVFluidHeatTransferPhysics.h b/modules/navier_stokes/include/physics/WCNSFVFluidHeatTransferPhysics.h index abfe625caaaf..9a452fbd9961 100644 --- a/modules/navier_stokes/include/physics/WCNSFVFluidHeatTransferPhysics.h +++ b/modules/navier_stokes/include/physics/WCNSFVFluidHeatTransferPhysics.h @@ -52,7 +52,7 @@ class WCNSFVFluidHeatTransferPhysics final : public NavierStokesPhysicsBase, protected: private: void actOnAdditionalTasks() override; - void addNonlinearVariables() override; + void addSolverVariables() override; void addInitialConditions() override; void addFVKernels() override; void addFVBCs() override; diff --git a/modules/navier_stokes/include/physics/WCNSFVScalarTransportPhysics.h b/modules/navier_stokes/include/physics/WCNSFVScalarTransportPhysics.h index acc77b7144f6..d5b178b77195 100644 --- a/modules/navier_stokes/include/physics/WCNSFVScalarTransportPhysics.h +++ b/modules/navier_stokes/include/physics/WCNSFVScalarTransportPhysics.h @@ -29,7 +29,7 @@ class WCNSFVScalarTransportPhysics : public WCNSFVScalarTransportPhysicsBase WCNSFVScalarTransportPhysics(const InputParameters & parameters); private: - virtual void addNonlinearVariables() override; + virtual void addSolverVariables() override; /** * Functions adding kernels for the incompressible / weakly-compressible scalar transport diff --git a/modules/navier_stokes/include/physics/WCNSFVTurbulencePhysics.h b/modules/navier_stokes/include/physics/WCNSFVTurbulencePhysics.h index b89bfd09321e..fb5fef49ac4b 100644 --- a/modules/navier_stokes/include/physics/WCNSFVTurbulencePhysics.h +++ b/modules/navier_stokes/include/physics/WCNSFVTurbulencePhysics.h @@ -49,7 +49,7 @@ class WCNSFVTurbulencePhysics final : public NavierStokesPhysicsBase, /// to add the relevant terms (turbulent diffusion notably) void retrieveCoupledPhysics(); - virtual void addNonlinearVariables() override; + virtual void addSolverVariables() override; virtual void addAuxiliaryVariables() override; virtual void addFVKernels() override; virtual void addFVBCs() override; diff --git a/modules/navier_stokes/include/physics/WCNSLinearFVFlowPhysics.h b/modules/navier_stokes/include/physics/WCNSLinearFVFlowPhysics.h index 66c67fe0589c..d44c79ffedd7 100644 --- a/modules/navier_stokes/include/physics/WCNSLinearFVFlowPhysics.h +++ b/modules/navier_stokes/include/physics/WCNSLinearFVFlowPhysics.h @@ -32,7 +32,7 @@ class WCNSLinearFVFlowPhysics final : public WCNSFVFlowPhysicsBase virtual void initializePhysicsAdditional() override; private: - virtual void addNonlinearVariables() override; + virtual void addSolverVariables() override; virtual void addFVKernels() override; virtual void addUserObjects() override; diff --git a/modules/navier_stokes/include/physics/WCNSLinearFVScalarTransportPhysics.h b/modules/navier_stokes/include/physics/WCNSLinearFVScalarTransportPhysics.h index 5f2b145542c4..b02447a91808 100644 --- a/modules/navier_stokes/include/physics/WCNSLinearFVScalarTransportPhysics.h +++ b/modules/navier_stokes/include/physics/WCNSLinearFVScalarTransportPhysics.h @@ -13,7 +13,7 @@ /** * Creates all the objects needed to solve the Navier Stokes scalar transport equations - * using the nonlinear finite volume weakly-compressible discretization (WCNSFV) + * using the linear finite volume weakly-compressible discretization (WCNSFV) */ class WCNSLinearFVScalarTransportPhysics : public WCNSFVScalarTransportPhysicsBase { @@ -23,8 +23,7 @@ class WCNSLinearFVScalarTransportPhysics : public WCNSFVScalarTransportPhysicsBa WCNSLinearFVScalarTransportPhysics(const InputParameters & parameters); private: - // TODO Rename to linear or solver? - virtual void addNonlinearVariables() override; + virtual void addSolverVariables() override; /** * Functions adding kernels for the incompressible / weakly-compressible scalar transport diff --git a/modules/navier_stokes/include/userobjects/RhieChowMassFlux.h b/modules/navier_stokes/include/userobjects/RhieChowMassFlux.h index 973e93db2e01..21fd8990a514 100644 --- a/modules/navier_stokes/include/userobjects/RhieChowMassFlux.h +++ b/modules/navier_stokes/include/userobjects/RhieChowMassFlux.h @@ -44,8 +44,10 @@ class RhieChowMassFlux : public GeneralUserObject, static InputParameters validParams(); RhieChowMassFlux(const InputParameters & params); - /// Get the face velocity (used in advection terms) + /// Get the face velocity times density (used in advection terms) Real getMassFlux(const FaceInfo & fi) const; + /// Get the face velocity (used in advection terms) + Real getFaceVelocity(const FaceInfo & fi) const; /// Initialize the container for face velocities void initFaceMassFlux(); diff --git a/modules/navier_stokes/src/executioners/SIMPLESolve.C b/modules/navier_stokes/src/executioners/SIMPLESolve.C index 3c806aa916b5..75e87778d559 100644 --- a/modules/navier_stokes/src/executioners/SIMPLESolve.C +++ b/modules/navier_stokes/src/executioners/SIMPLESolve.C @@ -36,7 +36,6 @@ SIMPLESolve::SIMPLESolve(Executioner & ex) } // and for the passive scalar equations if (_has_passive_scalar_systems) - { for (auto system_i : index_range(_passive_scalar_system_names)) { _passive_scalar_system_numbers.push_back( @@ -44,7 +43,6 @@ SIMPLESolve::SIMPLESolve(Executioner & ex) _passive_scalar_systems.push_back( &_problem.getLinearSystem(_passive_scalar_system_numbers[system_i])); } - } } void @@ -377,7 +375,7 @@ SIMPLESolve::solve() } // If we have passive scalar equations, solve them here. We assume the material properties in the - // Navier-Stokes equations do not depend on passive scalar, as they are passive, therefore we + // Navier-Stokes equations do not depend on passive scalars, as they are passive, therefore we // solve outside of the velocity-pressure loop if (_has_passive_scalar_systems && converged) { diff --git a/modules/navier_stokes/src/linearfvkernels/LinearFVScalarAdvection.C b/modules/navier_stokes/src/linearfvkernels/LinearFVScalarAdvection.C index 97c10b9b52b0..871e43f74d9b 100644 --- a/modules/navier_stokes/src/linearfvkernels/LinearFVScalarAdvection.C +++ b/modules/navier_stokes/src/linearfvkernels/LinearFVScalarAdvection.C @@ -31,7 +31,7 @@ LinearFVScalarAdvection::LinearFVScalarAdvection(const InputParameters & params) : LinearFVFluxKernel(params), _mass_flux_provider(getUserObject("rhie_chow_user_object")), _advected_interp_coeffs(std::make_pair(0, 0)), - _face_mass_flux(0.0) + _face_velocity(0.0) { Moose::FV::setInterpolationMethod(*this, _advected_interp_method, "advected_interp_method"); } @@ -39,13 +39,13 @@ LinearFVScalarAdvection::LinearFVScalarAdvection(const InputParameters & params) Real LinearFVScalarAdvection::computeElemMatrixContribution() { - return _advected_interp_coeffs.first * _face_mass_flux * _current_face_area; + return _advected_interp_coeffs.first * _face_velocity * _current_face_area; } Real LinearFVScalarAdvection::computeNeighborMatrixContribution() { - return _advected_interp_coeffs.second * _face_mass_flux * _current_face_area; + return _advected_interp_coeffs.second * _face_velocity * _current_face_area; } Real @@ -71,7 +71,7 @@ LinearFVScalarAdvection::computeBoundaryMatrixContribution(const LinearFVBoundar // We support internal boundaries too so we have to make sure the normal points always outward const auto factor = (_current_face_type == FaceInfo::VarFaceNeighbors::ELEM) ? 1.0 : -1.0; - return boundary_value_matrix_contrib * factor * _face_mass_flux * _current_face_area; + return boundary_value_matrix_contrib * factor * _face_velocity * _current_face_area; } Real @@ -84,7 +84,7 @@ LinearFVScalarAdvection::computeBoundaryRHSContribution(const LinearFVBoundaryCo const auto factor = (_current_face_type == FaceInfo::VarFaceNeighbors::ELEM ? 1.0 : -1.0); const auto boundary_value_rhs_contrib = adv_bc->computeBoundaryValueRHSContribution(); - return -boundary_value_rhs_contrib * factor * _face_mass_flux * _current_face_area; + return -boundary_value_rhs_contrib * factor * _face_velocity * _current_face_area; } void @@ -92,12 +92,12 @@ LinearFVScalarAdvection::setupFaceData(const FaceInfo * face_info) { LinearFVFluxKernel::setupFaceData(face_info); - // Caching the mass flux on the face which will be reused in the advection term's matrix and right + // Caching the velocity on the face which will be reused in the advection term's matrix and right // hand side contributions - _face_mass_flux = _mass_flux_provider.getMassFlux(*face_info); + _face_velocity = _mass_flux_provider.getFaceVelocity(*face_info); // Caching the interpolation coefficients so they will be reused for the matrix and right hand // side terms _advected_interp_coeffs = - interpCoeffs(_advected_interp_method, *_current_face_info, true, _face_mass_flux); + interpCoeffs(_advected_interp_method, *_current_face_info, true, _face_velocity); } diff --git a/modules/navier_stokes/src/physics/PNSFVSolidHeatTransferPhysics.C b/modules/navier_stokes/src/physics/PNSFVSolidHeatTransferPhysics.C index b7b41e20118c..e8ab4a785299 100644 --- a/modules/navier_stokes/src/physics/PNSFVSolidHeatTransferPhysics.C +++ b/modules/navier_stokes/src/physics/PNSFVSolidHeatTransferPhysics.C @@ -154,7 +154,7 @@ PNSFVSolidHeatTransferPhysics::PNSFVSolidHeatTransferPhysics(const InputParamete } void -PNSFVSolidHeatTransferPhysics::addNonlinearVariables() +PNSFVSolidHeatTransferPhysics::addSolverVariables() { // Dont add if the user already defined the variable if (variableExists(_solid_temperature_name, diff --git a/modules/navier_stokes/src/physics/WCNSFVFlowPhysics.C b/modules/navier_stokes/src/physics/WCNSFVFlowPhysics.C index d10855f61d1c..307147cb578e 100644 --- a/modules/navier_stokes/src/physics/WCNSFVFlowPhysics.C +++ b/modules/navier_stokes/src/physics/WCNSFVFlowPhysics.C @@ -145,7 +145,7 @@ WCNSFVFlowPhysics::WCNSFVFlowPhysics(const InputParameters & parameters) } void -WCNSFVFlowPhysics::addNonlinearVariables() +WCNSFVFlowPhysics::addSolverVariables() { if (!_has_flow_equations) return; diff --git a/modules/navier_stokes/src/physics/WCNSFVFluidHeatTransferPhysics.C b/modules/navier_stokes/src/physics/WCNSFVFluidHeatTransferPhysics.C index 0ccaee22b1ea..c292fb63de9f 100644 --- a/modules/navier_stokes/src/physics/WCNSFVFluidHeatTransferPhysics.C +++ b/modules/navier_stokes/src/physics/WCNSFVFluidHeatTransferPhysics.C @@ -115,7 +115,7 @@ WCNSFVFluidHeatTransferPhysics::WCNSFVFluidHeatTransferPhysics(const InputParame } void -WCNSFVFluidHeatTransferPhysics::addNonlinearVariables() +WCNSFVFluidHeatTransferPhysics::addSolverVariables() { // For compatibility with Modules/NavierStokesFV syntax if (!_has_energy_equation) diff --git a/modules/navier_stokes/src/physics/WCNSFVScalarTransportPhysics.C b/modules/navier_stokes/src/physics/WCNSFVScalarTransportPhysics.C index 92742d12e10b..7709d18a7baa 100644 --- a/modules/navier_stokes/src/physics/WCNSFVScalarTransportPhysics.C +++ b/modules/navier_stokes/src/physics/WCNSFVScalarTransportPhysics.C @@ -32,7 +32,7 @@ WCNSFVScalarTransportPhysics::WCNSFVScalarTransportPhysics(const InputParameters } void -WCNSFVScalarTransportPhysics::addNonlinearVariables() +WCNSFVScalarTransportPhysics::addSolverVariables() { // For compatibility with Modules/NavierStokesFV syntax if (!_has_scalar_equation) diff --git a/modules/navier_stokes/src/physics/WCNSFVTurbulencePhysics.C b/modules/navier_stokes/src/physics/WCNSFVTurbulencePhysics.C index 5bf5a7a7a1ba..15898c2aa58d 100644 --- a/modules/navier_stokes/src/physics/WCNSFVTurbulencePhysics.C +++ b/modules/navier_stokes/src/physics/WCNSFVTurbulencePhysics.C @@ -330,7 +330,7 @@ WCNSFVTurbulencePhysics::retrieveCoupledPhysics() } void -WCNSFVTurbulencePhysics::addNonlinearVariables() +WCNSFVTurbulencePhysics::addSolverVariables() { if (_turbulence_model == "mixing-length" || _turbulence_model == "none") return; diff --git a/modules/navier_stokes/src/physics/WCNSLinearFVFlowPhysics.C b/modules/navier_stokes/src/physics/WCNSLinearFVFlowPhysics.C index e8b18c78e74b..0a3046682feb 100644 --- a/modules/navier_stokes/src/physics/WCNSLinearFVFlowPhysics.C +++ b/modules/navier_stokes/src/physics/WCNSLinearFVFlowPhysics.C @@ -83,7 +83,7 @@ WCNSLinearFVFlowPhysics::initializePhysicsAdditional() } void -WCNSLinearFVFlowPhysics::addNonlinearVariables() +WCNSLinearFVFlowPhysics::addSolverVariables() { if (!_has_flow_equations) return; diff --git a/modules/navier_stokes/src/physics/WCNSLinearFVScalarTransportPhysics.C b/modules/navier_stokes/src/physics/WCNSLinearFVScalarTransportPhysics.C index 6feb578c69f4..7fd0e90fc64e 100644 --- a/modules/navier_stokes/src/physics/WCNSLinearFVScalarTransportPhysics.C +++ b/modules/navier_stokes/src/physics/WCNSLinearFVScalarTransportPhysics.C @@ -34,7 +34,7 @@ WCNSLinearFVScalarTransportPhysics::WCNSLinearFVScalarTransportPhysics( } void -WCNSLinearFVScalarTransportPhysics::addNonlinearVariables() +WCNSLinearFVScalarTransportPhysics::addSolverVariables() { // For compatibility with Modules/NavierStokesFV syntax if (!_has_scalar_equation) diff --git a/modules/navier_stokes/src/userobjects/RhieChowMassFlux.C b/modules/navier_stokes/src/userobjects/RhieChowMassFlux.C index 4e33d9722d40..991b6af20483 100644 --- a/modules/navier_stokes/src/userobjects/RhieChowMassFlux.C +++ b/modules/navier_stokes/src/userobjects/RhieChowMassFlux.C @@ -234,6 +234,18 @@ RhieChowMassFlux::getMassFlux(const FaceInfo & fi) const return _face_mass_flux.evaluate(&fi); } +Real +RhieChowMassFlux::getFaceVelocity(const FaceInfo & fi) const +{ + const Moose::FaceArg face_arg{&fi, + /*limiter_type=*/Moose::FV::LimiterType::CentralDifference, + /*elem_is_upwind=*/true, + /*correct_skewness=*/false, + &fi.elem()}; + const Real face_rho = _rho(face_arg, Moose::currentState()); + return _face_mass_flux.evaluate(&fi) / face_rho; +} + void RhieChowMassFlux::computeFaceMassFlux() { @@ -293,7 +305,7 @@ RhieChowMassFlux::computeFaceMassFlux() const auto rhs_contribution = _p_diffusion_kernel->computeBoundaryRHSContribution(*bc_pointer); - // On the boundary, only the element side has a constibution + // On the boundary, only the element side has a contribution p_grad_flux = (p_elem_value * matrix_contribution - rhs_contribution); } // Compute the new face flux diff --git a/modules/scalar_transport/include/physics/MultiSpeciesDiffusionCG.h b/modules/scalar_transport/include/physics/MultiSpeciesDiffusionCG.h index 8464c3ecf7fd..0d7ee2da158e 100644 --- a/modules/scalar_transport/include/physics/MultiSpeciesDiffusionCG.h +++ b/modules/scalar_transport/include/physics/MultiSpeciesDiffusionCG.h @@ -23,7 +23,7 @@ class MultiSpeciesDiffusionCG : public MultiSpeciesDiffusionPhysicsBase MultiSpeciesDiffusionCG(const InputParameters & parameters); protected: - virtual void addNonlinearVariables() override; + virtual void addSolverVariables() override; virtual void addFEKernels() override; virtual void addFEBCs() override; }; diff --git a/modules/scalar_transport/src/physics/MultiSpeciesDiffusionCG.C b/modules/scalar_transport/src/physics/MultiSpeciesDiffusionCG.C index 614aec8871dd..a660d7ddeb63 100644 --- a/modules/scalar_transport/src/physics/MultiSpeciesDiffusionCG.C +++ b/modules/scalar_transport/src/physics/MultiSpeciesDiffusionCG.C @@ -228,7 +228,7 @@ MultiSpeciesDiffusionCG::addFEBCs() } void -MultiSpeciesDiffusionCG::addNonlinearVariables() +MultiSpeciesDiffusionCG::addSolverVariables() { for (const auto & var_name : _species_names) { From 84790cf89acc993f1f48536a6e5c52922c835a6b Mon Sep 17 00:00:00 2001 From: Guillaume Giudicelli Date: Wed, 4 Dec 2024 09:55:26 -0700 Subject: [PATCH 9/9] Address Peter's review - face velocity -> volumetric face flux due to A term Co-authored-by: Peter German <31662443+grmnptr@users.noreply.github.com> Apply suggestions from code review Co-authored-by: Peter German <31662443+grmnptr@users.noreply.github.com> --- .../linearfvkernels/LinearFVScalarAdvection.md | 2 +- .../physics/WCNSLinearFVScalarTransportPhysics.md | 2 +- .../linearfvkernels/LinearFVScalarAdvection.h | 2 +- .../include/userobjects/RhieChowMassFlux.h | 4 ++-- .../src/linearfvkernels/LinearFVScalarAdvection.C | 14 +++++++------- .../src/userobjects/RhieChowMassFlux.C | 7 ++++--- 6 files changed, 16 insertions(+), 15 deletions(-) diff --git a/modules/navier_stokes/doc/content/source/linearfvkernels/LinearFVScalarAdvection.md b/modules/navier_stokes/doc/content/source/linearfvkernels/LinearFVScalarAdvection.md index e3215a7977bc..6b6ab638d4c5 100644 --- a/modules/navier_stokes/doc/content/source/linearfvkernels/LinearFVScalarAdvection.md +++ b/modules/navier_stokes/doc/content/source/linearfvkernels/LinearFVScalarAdvection.md @@ -12,7 +12,7 @@ For FV, the integral of the advection term of scalar $C_i$ over a cell can be ex where $C_{if}$ is the face value of the scalar concentration. An interpolation scheme (e.g. upwind) can be used to compute the face value. This kernel adds the face contribution for each face $f$ to the right hand side and matrix. -The face mass flux $(\vec{u}\cdot \vec{n})_{RC}$ is provided by the [RhieChowMassFlux.md] object which uses pressure +The volumetric face flux $(\vec{u}\cdot \vec{n})_{RC}$ is provided by the [RhieChowMassFlux.md] object which uses pressure gradients and the discrete momentum equation to compute face velocities and mass fluxes. For more information on the expression that is used, see [SIMPLE.md]. diff --git a/modules/navier_stokes/doc/content/source/physics/WCNSLinearFVScalarTransportPhysics.md b/modules/navier_stokes/doc/content/source/physics/WCNSLinearFVScalarTransportPhysics.md index 52325cd939ab..25851d427383 100644 --- a/modules/navier_stokes/doc/content/source/physics/WCNSLinearFVScalarTransportPhysics.md +++ b/modules/navier_stokes/doc/content/source/physics/WCNSLinearFVScalarTransportPhysics.md @@ -15,7 +15,7 @@ where: - $\phi_i$ is the i-th scalar quantity - \mathbf{v} is the advecting velocity - $k_i$ the i-th scalar diffusivity -- $Q_i$ is the i-th precursor source +- $Q_i$ is the i-th scalar source - $\lambda_i$ is a reaction coefficient. It should be negative for a loss term The kernels created are: diff --git a/modules/navier_stokes/include/linearfvkernels/LinearFVScalarAdvection.h b/modules/navier_stokes/include/linearfvkernels/LinearFVScalarAdvection.h index 0084904be8b1..3dffd7e69446 100644 --- a/modules/navier_stokes/include/linearfvkernels/LinearFVScalarAdvection.h +++ b/modules/navier_stokes/include/linearfvkernels/LinearFVScalarAdvection.h @@ -47,7 +47,7 @@ class LinearFVScalarAdvection : public LinearFVFluxKernel /// Container for the velocity on the face which will be reused in the advection term's /// matrix and right hand side contribution - Real _face_velocity; + Real _volumetric_face_flux; /// The interpolation method to use for the advected quantity Moose::FV::InterpMethod _advected_interp_method; diff --git a/modules/navier_stokes/include/userobjects/RhieChowMassFlux.h b/modules/navier_stokes/include/userobjects/RhieChowMassFlux.h index 21fd8990a514..b88d089a3a02 100644 --- a/modules/navier_stokes/include/userobjects/RhieChowMassFlux.h +++ b/modules/navier_stokes/include/userobjects/RhieChowMassFlux.h @@ -46,8 +46,8 @@ class RhieChowMassFlux : public GeneralUserObject, /// Get the face velocity times density (used in advection terms) Real getMassFlux(const FaceInfo & fi) const; - /// Get the face velocity (used in advection terms) - Real getFaceVelocity(const FaceInfo & fi) const; + /// Get the volumetric face flux (used in advection terms) + Real getVolumetricFaceFlux(const FaceInfo & fi) const; /// Initialize the container for face velocities void initFaceMassFlux(); diff --git a/modules/navier_stokes/src/linearfvkernels/LinearFVScalarAdvection.C b/modules/navier_stokes/src/linearfvkernels/LinearFVScalarAdvection.C index 871e43f74d9b..175c5057934b 100644 --- a/modules/navier_stokes/src/linearfvkernels/LinearFVScalarAdvection.C +++ b/modules/navier_stokes/src/linearfvkernels/LinearFVScalarAdvection.C @@ -31,7 +31,7 @@ LinearFVScalarAdvection::LinearFVScalarAdvection(const InputParameters & params) : LinearFVFluxKernel(params), _mass_flux_provider(getUserObject("rhie_chow_user_object")), _advected_interp_coeffs(std::make_pair(0, 0)), - _face_velocity(0.0) + _volumetric_face_flux(0.0) { Moose::FV::setInterpolationMethod(*this, _advected_interp_method, "advected_interp_method"); } @@ -39,13 +39,13 @@ LinearFVScalarAdvection::LinearFVScalarAdvection(const InputParameters & params) Real LinearFVScalarAdvection::computeElemMatrixContribution() { - return _advected_interp_coeffs.first * _face_velocity * _current_face_area; + return _advected_interp_coeffs.first * _volumetric_face_flux * _current_face_area; } Real LinearFVScalarAdvection::computeNeighborMatrixContribution() { - return _advected_interp_coeffs.second * _face_velocity * _current_face_area; + return _advected_interp_coeffs.second * _volumetric_face_flux * _current_face_area; } Real @@ -71,7 +71,7 @@ LinearFVScalarAdvection::computeBoundaryMatrixContribution(const LinearFVBoundar // We support internal boundaries too so we have to make sure the normal points always outward const auto factor = (_current_face_type == FaceInfo::VarFaceNeighbors::ELEM) ? 1.0 : -1.0; - return boundary_value_matrix_contrib * factor * _face_velocity * _current_face_area; + return boundary_value_matrix_contrib * factor * _volumetric_face_flux * _current_face_area; } Real @@ -84,7 +84,7 @@ LinearFVScalarAdvection::computeBoundaryRHSContribution(const LinearFVBoundaryCo const auto factor = (_current_face_type == FaceInfo::VarFaceNeighbors::ELEM ? 1.0 : -1.0); const auto boundary_value_rhs_contrib = adv_bc->computeBoundaryValueRHSContribution(); - return -boundary_value_rhs_contrib * factor * _face_velocity * _current_face_area; + return -boundary_value_rhs_contrib * factor * _volumetric_face_flux * _current_face_area; } void @@ -94,10 +94,10 @@ LinearFVScalarAdvection::setupFaceData(const FaceInfo * face_info) // Caching the velocity on the face which will be reused in the advection term's matrix and right // hand side contributions - _face_velocity = _mass_flux_provider.getFaceVelocity(*face_info); + _volumetric_face_flux = _mass_flux_provider.getVolumetricFaceFlux(*face_info); // Caching the interpolation coefficients so they will be reused for the matrix and right hand // side terms _advected_interp_coeffs = - interpCoeffs(_advected_interp_method, *_current_face_info, true, _face_velocity); + interpCoeffs(_advected_interp_method, *_current_face_info, true, _volumetric_face_flux); } diff --git a/modules/navier_stokes/src/userobjects/RhieChowMassFlux.C b/modules/navier_stokes/src/userobjects/RhieChowMassFlux.C index 991b6af20483..a0cc5e3a33cd 100644 --- a/modules/navier_stokes/src/userobjects/RhieChowMassFlux.C +++ b/modules/navier_stokes/src/userobjects/RhieChowMassFlux.C @@ -235,15 +235,16 @@ RhieChowMassFlux::getMassFlux(const FaceInfo & fi) const } Real -RhieChowMassFlux::getFaceVelocity(const FaceInfo & fi) const +RhieChowMassFlux::getVolumetricFaceFlux(const FaceInfo & fi) const { const Moose::FaceArg face_arg{&fi, /*limiter_type=*/Moose::FV::LimiterType::CentralDifference, /*elem_is_upwind=*/true, /*correct_skewness=*/false, - &fi.elem()}; + &fi.elem(), + /*state_limiter*/ nullptr}; const Real face_rho = _rho(face_arg, Moose::currentState()); - return _face_mass_flux.evaluate(&fi) / face_rho; + return libmesh_map_find(_face_mass_flux, fi.id()) / face_rho; } void