Line data Source code
1 : //* This file is part of the MOOSE framework
2 : //* https://mooseframework.inl.gov
3 : //*
4 : //* All rights reserved, see COPYRIGHT for full restrictions
5 : //* https://github.com/idaholab/moose/blob/master/COPYRIGHT
6 : //*
7 : //* Licensed under LGPL 2.1, please see LICENSE for details
8 : //* https://www.gnu.org/licenses/lgpl-2.1.html
9 :
10 : #include "PCNSFVStrongBC.h"
11 : #include "NS.h"
12 : #include "SinglePhaseFluidProperties.h"
13 : #include "Function.h"
14 : #include "MfrPostprocessor.h"
15 :
16 : registerMooseObject("NavierStokesApp", PCNSFVStrongBC);
17 :
18 : InputParameters
19 7912 : PCNSFVStrongBC::validParams()
20 : {
21 7912 : InputParameters params = FVFluxBC::validParams();
22 7912 : params.addClassDescription("Computes the residual of advective term using finite volume method.");
23 7912 : params.addRequiredParam<UserObjectName>(NS::fluid, "Fluid properties userobject");
24 15824 : MooseEnum eqn("mass momentum energy scalar");
25 15824 : params.addRequiredParam<MooseEnum>("eqn", eqn, "The equation you're solving.");
26 15824 : MooseEnum momentum_component("x=0 y=1 z=2");
27 15824 : params.addParam<MooseEnum>("momentum_component",
28 : momentum_component,
29 : "The component of the momentum equation that this kernel applies to.");
30 15824 : params.addParam<bool>("velocity_function_includes_rho",
31 15824 : false,
32 : "Whether the provided superficial velocity function actually includes "
33 : "multiplication by rho (e.g. the function is representative of momentum.");
34 7912 : params.addParam<FunctionName>(
35 : NS::superficial_velocity,
36 : "An optional name of a vector function for the velocity. If not provided then the "
37 : "superficial velocity will be treated implicitly (e.g. we will use the interior value");
38 7912 : params.addParam<FunctionName>(
39 : NS::pressure,
40 : "An optional name of a function for the pressure. If not provided then the pressure will be "
41 : "treated implicitly (e.g. we will use the interior value");
42 7912 : params.addParam<FunctionName>(
43 : NS::T_fluid,
44 : "An optional name of a function for the fluid temperature. If not provided then the fluid "
45 : "temperature will be treated implicitly (e.g. we will use the interior value");
46 15824 : params.addParam<FunctionName>(
47 : "scalar",
48 : "A function describing the value of the scalar at the boundary. If this function is not "
49 : "provided, then the interior value will be used.");
50 15824 : params.addParam<UserObjectName>("mfr_postprocessor",
51 : "A postprocessor used for outputting mass flow rates on the same "
52 : "boundary this object acts on");
53 7912 : return params;
54 7912 : }
55 :
56 4034 : PCNSFVStrongBC::PCNSFVStrongBC(const InputParameters & params)
57 : : FVFluxBC(params),
58 8068 : _fluid(getUserObject<SinglePhaseFluidProperties>(NS::fluid)),
59 4034 : _dim(_mesh.dimension()),
60 4034 : _sup_vel_x_elem(getADMaterialProperty<Real>(NS::superficial_velocity_x)),
61 4034 : _sup_vel_x_neighbor(getNeighborADMaterialProperty<Real>(NS::superficial_velocity_x)),
62 4034 : _grad_sup_vel_x_elem(
63 4034 : getADMaterialProperty<RealVectorValue>(NS::grad(NS::superficial_velocity_x))),
64 4034 : _grad_sup_vel_x_neighbor(
65 4034 : getNeighborADMaterialProperty<RealVectorValue>(NS::grad(NS::superficial_velocity_x))),
66 4034 : _sup_vel_y_elem(getADMaterialProperty<Real>(NS::superficial_velocity_y)),
67 4034 : _sup_vel_y_neighbor(getNeighborADMaterialProperty<Real>(NS::superficial_velocity_y)),
68 4034 : _grad_sup_vel_y_elem(
69 4034 : getADMaterialProperty<RealVectorValue>(NS::grad(NS::superficial_velocity_y))),
70 4034 : _grad_sup_vel_y_neighbor(
71 4034 : getNeighborADMaterialProperty<RealVectorValue>(NS::grad(NS::superficial_velocity_y))),
72 4034 : _sup_vel_z_elem(getADMaterialProperty<Real>(NS::superficial_velocity_z)),
73 4034 : _sup_vel_z_neighbor(getNeighborADMaterialProperty<Real>(NS::superficial_velocity_z)),
74 4034 : _grad_sup_vel_z_elem(
75 4034 : getADMaterialProperty<RealVectorValue>(NS::grad(NS::superficial_velocity_z))),
76 4034 : _grad_sup_vel_z_neighbor(
77 4034 : getNeighborADMaterialProperty<RealVectorValue>(NS::grad(NS::superficial_velocity_z))),
78 4034 : _T_fluid_elem(getADMaterialProperty<Real>(NS::T_fluid)),
79 4034 : _T_fluid_neighbor(getNeighborADMaterialProperty<Real>(NS::T_fluid)),
80 4034 : _grad_T_fluid_elem(getADMaterialProperty<RealVectorValue>(NS::grad(NS::T_fluid))),
81 4034 : _grad_T_fluid_neighbor(getNeighborADMaterialProperty<RealVectorValue>(NS::grad(NS::T_fluid))),
82 4034 : _pressure_elem(getADMaterialProperty<Real>(NS::pressure)),
83 4034 : _pressure_neighbor(getNeighborADMaterialProperty<Real>(NS::pressure)),
84 4034 : _grad_pressure_elem(getADMaterialProperty<RealVectorValue>(NS::grad(NS::pressure))),
85 4034 : _grad_pressure_neighbor(getNeighborADMaterialProperty<RealVectorValue>(NS::grad(NS::pressure))),
86 4034 : _eps_elem(getMaterialProperty<Real>(NS::porosity)),
87 4034 : _eps_neighbor(getNeighborMaterialProperty<Real>(NS::porosity)),
88 8068 : _eqn(getParam<MooseEnum>("eqn")),
89 8068 : _index(getParam<MooseEnum>("momentum_component")),
90 4034 : _sup_vel_provided(isParamValid(NS::superficial_velocity)),
91 4034 : _pressure_provided(isParamValid(NS::pressure)),
92 4034 : _T_fluid_provided(isParamValid(NS::T_fluid)),
93 4034 : _sup_vel_function(_sup_vel_provided ? &getFunction(NS::superficial_velocity) : nullptr),
94 4034 : _pressure_function(_pressure_provided ? &getFunction(NS::pressure) : nullptr),
95 4034 : _T_fluid_function(_T_fluid_provided ? &getFunction(NS::T_fluid) : nullptr),
96 4034 : _scalar_elem(_u),
97 4034 : _scalar_neighbor(_u_neighbor),
98 4034 : _grad_scalar_elem((_eqn == "scalar") ? &_var.adGradSln() : nullptr),
99 4034 : _grad_scalar_neighbor((_eqn == "scalar") ? &_var.adGradSlnNeighbor() : nullptr),
100 8068 : _scalar_function_provided(isParamValid("scalar")),
101 4034 : _scalar_function(_scalar_function_provided ? &getFunction("scalar") : nullptr),
102 8068 : _velocity_function_includes_rho(getParam<bool>("velocity_function_includes_rho")),
103 4034 : _mfr_pp(
104 4034 : isParamValid("mfr_postprocessor")
105 4034 : ? &const_cast<MfrPostprocessor &>(getUserObject<MfrPostprocessor>("mfr_postprocessor"))
106 4034 : : nullptr)
107 : {
108 8294 : if ((_eqn == "momentum") && !isParamValid("momentum_component"))
109 0 : paramError("eqn",
110 : "If 'momentum' is specified for 'eqn', then you must provide a parameter "
111 : "value for 'momentum_component'");
112 9262 : if ((_eqn != "momentum") && isParamValid("momentum_component"))
113 0 : paramError("momentum_component",
114 : "'momentum_component' should not be specified when the 'eqn' is not 'momentum'");
115 4034 : }
116 :
117 : ADReal
118 169228 : PCNSFVStrongBC::computeQpResidual()
119 : {
120 169228 : const auto ft = _face_type;
121 : const bool out_of_elem = (ft == FaceInfo::VarFaceNeighbors::ELEM);
122 169228 : const auto normal = out_of_elem ? _face_info->normal() : Point(-_face_info->normal());
123 :
124 169228 : const auto & sup_vel_x_interior = out_of_elem ? _sup_vel_x_elem[_qp] : _sup_vel_x_neighbor[_qp];
125 169228 : const auto & sup_vel_y_interior = out_of_elem ? _sup_vel_y_elem[_qp] : _sup_vel_y_neighbor[_qp];
126 169228 : const auto & sup_vel_z_interior = out_of_elem ? _sup_vel_z_elem[_qp] : _sup_vel_z_neighbor[_qp];
127 169228 : const auto & pressure_interior = out_of_elem ? _pressure_elem[_qp] : _pressure_neighbor[_qp];
128 169228 : const auto & T_fluid_interior = out_of_elem ? _T_fluid_elem[_qp] : _T_fluid_neighbor[_qp];
129 :
130 : const auto & grad_sup_vel_x_interior =
131 169228 : out_of_elem ? _grad_sup_vel_x_elem[_qp] : _grad_sup_vel_x_neighbor[_qp];
132 : const auto & grad_sup_vel_y_interior =
133 169228 : out_of_elem ? _grad_sup_vel_y_elem[_qp] : _grad_sup_vel_y_neighbor[_qp];
134 : const auto & grad_sup_vel_z_interior =
135 169228 : out_of_elem ? _grad_sup_vel_z_elem[_qp] : _grad_sup_vel_z_neighbor[_qp];
136 : const auto & grad_pressure_interior =
137 169228 : out_of_elem ? _grad_pressure_elem[_qp] : _grad_pressure_neighbor[_qp];
138 : const auto & grad_T_fluid_interior =
139 169228 : out_of_elem ? _grad_T_fluid_elem[_qp] : _grad_T_fluid_neighbor[_qp];
140 : const auto & interior_centroid =
141 169228 : out_of_elem ? _face_info->elemCentroid() : _face_info->neighborCentroid();
142 169228 : const auto dCf = _face_info->faceCentroid() - interior_centroid;
143 169228 : const auto eps_interior = out_of_elem ? _eps_elem[_qp] : _eps_neighbor[_qp];
144 :
145 : const auto pressure_boundary =
146 169228 : _pressure_provided ? ADReal(_pressure_function->value(_t, _face_info->faceCentroid()))
147 169228 : : pressure_interior + grad_pressure_interior * dCf;
148 : const auto T_fluid_boundary =
149 169228 : _T_fluid_provided ? ADReal(_T_fluid_function->value(_t, _face_info->faceCentroid()))
150 169228 : : T_fluid_interior + grad_T_fluid_interior * dCf;
151 169228 : const auto rho_boundary = _fluid.rho_from_p_T(pressure_boundary, T_fluid_boundary);
152 :
153 : ADReal sup_vel_x_boundary;
154 169228 : if (_sup_vel_provided)
155 : {
156 84614 : sup_vel_x_boundary = _sup_vel_function->vectorValue(_t, _face_info->faceCentroid())(0);
157 84614 : if (_velocity_function_includes_rho)
158 28512 : sup_vel_x_boundary /= rho_boundary;
159 : }
160 : else
161 169228 : sup_vel_x_boundary = sup_vel_x_interior + grad_sup_vel_x_interior * dCf;
162 :
163 169228 : ADReal sup_vel_y_boundary = 0;
164 169228 : if (_dim >= 2)
165 : {
166 107128 : if (_sup_vel_provided)
167 : {
168 53564 : sup_vel_y_boundary = _sup_vel_function->vectorValue(_t, _face_info->faceCentroid())(1);
169 53564 : if (_velocity_function_includes_rho)
170 28512 : sup_vel_y_boundary /= rho_boundary;
171 : }
172 : else
173 107128 : sup_vel_y_boundary = sup_vel_y_interior + grad_sup_vel_y_interior * dCf;
174 : }
175 :
176 169228 : ADReal sup_vel_z_boundary = 0;
177 169228 : if (_dim >= 3)
178 : {
179 0 : if (_sup_vel_provided)
180 : {
181 0 : sup_vel_z_boundary = _sup_vel_function->vectorValue(_t, _face_info->faceCentroid())(2);
182 0 : if (_velocity_function_includes_rho)
183 0 : sup_vel_z_boundary /= rho_boundary;
184 : }
185 : else
186 0 : sup_vel_z_boundary = sup_vel_z_interior + grad_sup_vel_z_interior * dCf;
187 : }
188 :
189 : const VectorValue<ADReal> sup_vel_boundary(
190 : sup_vel_x_boundary, sup_vel_y_boundary, sup_vel_z_boundary);
191 169228 : const auto eps_boundary = eps_interior;
192 169228 : const auto u_boundary = sup_vel_boundary / eps_boundary;
193 169228 : const auto e_boundary = _fluid.e_from_p_T(pressure_boundary, T_fluid_boundary);
194 :
195 169228 : if (_eqn == "mass")
196 : {
197 47282 : const ADReal mfr = rho_boundary * sup_vel_boundary * normal;
198 47282 : if (_mfr_pp)
199 0 : _mfr_pp->setMfr(_face_info, mfr.value(), false);
200 47282 : return mfr;
201 : }
202 121946 : else if (_eqn == "momentum")
203 : {
204 73864 : const auto rhou_boundary = u_boundary(_index) * rho_boundary;
205 73864 : return rhou_boundary * sup_vel_boundary * normal +
206 147728 : eps_boundary * pressure_boundary * normal(_index);
207 : }
208 48082 : else if (_eqn == "energy")
209 : {
210 : const auto ht_boundary =
211 94564 : e_boundary + 0.5 * u_boundary * u_boundary + pressure_boundary / rho_boundary;
212 : const auto rho_ht_boundary = rho_boundary * ht_boundary;
213 47282 : return rho_ht_boundary * sup_vel_boundary * normal;
214 : }
215 800 : else if (_eqn == "scalar")
216 : {
217 800 : const auto & scalar_interior = out_of_elem ? _scalar_elem[_qp] : _scalar_neighbor[_qp];
218 : const auto & grad_scalar_interior =
219 800 : out_of_elem ? (*_grad_scalar_elem)[_qp] : (*_grad_scalar_neighbor)[_qp];
220 : const auto scalar_boundary =
221 800 : _scalar_function_provided ? ADReal(_scalar_function->value(_t, _face_info->faceCentroid()))
222 1200 : : scalar_interior + grad_scalar_interior * dCf;
223 800 : return rho_boundary * scalar_boundary * sup_vel_boundary * normal;
224 : }
225 : else
226 0 : mooseError("Unrecognized equation type ", _eqn);
227 : }
|
