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 "PINSFVMomentumFrictionCorrection.h"
11 : #include "INSFVRhieChowInterpolator.h"
12 : #include "NS.h"
13 : #include "SystemBase.h"
14 :
15 : registerMooseObject("NavierStokesApp", PINSFVMomentumFrictionCorrection);
16 :
17 : InputParameters
18 2502 : PINSFVMomentumFrictionCorrection::validParams()
19 : {
20 2502 : auto params = INSFVFluxKernel::validParams();
21 2502 : params.addClassDescription(
22 : "Computes a correction term to avoid oscillations from average pressure interpolation in "
23 : "regions of high changes in friction coefficients.");
24 5004 : params.addParam<MooseFunctorName>("Darcy_name", "Name of the Darcy coefficients property.");
25 5004 : params.addParam<MooseFunctorName>("Forchheimer_name",
26 : "Name of the Forchheimer coefficients property.");
27 2502 : params.addRequiredParam<MooseFunctorName>(NS::density, "The density.");
28 2502 : params.addParam<MooseFunctorName>(
29 : NS::speed,
30 : "The norm of the interstitial velocity. This is required for Forchheimer calculations");
31 2502 : params.addParam<MooseFunctorName>(NS::mu, NS::mu, "The dynamic viscosity");
32 5004 : params.addRangeCheckedParam<Real>("consistent_scaling",
33 : 1,
34 : "consistent_scaling >= 0",
35 : "Smoothing scaling parameter to control "
36 : "collocated mesh oscillations");
37 2502 : return params;
38 0 : }
39 :
40 1494 : PINSFVMomentumFrictionCorrection::PINSFVMomentumFrictionCorrection(const InputParameters & params)
41 : : INSFVFluxKernel(params),
42 4482 : _D(isParamValid("Darcy_name") ? &getFunctor<ADRealVectorValue>("Darcy_name") : nullptr),
43 5976 : _F(isParamValid("Forchheimer_name") ? &getFunctor<ADRealVectorValue>("Forchheimer_name")
44 : : nullptr),
45 2988 : _use_Darcy_friction_model(isParamValid("Darcy_name")),
46 2988 : _use_Forchheimer_friction_model(isParamValid("Forchheimer_name")),
47 1494 : _rho(getFunctor<ADReal>(NS::density)),
48 1494 : _mu(getFunctor<ADReal>(NS::mu)),
49 1494 : _speed(isParamValid(NS::speed) ? &getFunctor<ADReal>(NS::speed) : nullptr),
50 4482 : _consistent_scaling(getParam<Real>("consistent_scaling"))
51 : {
52 1494 : if (!_use_Darcy_friction_model && !_use_Forchheimer_friction_model)
53 0 : mooseError("At least one friction model needs to be specified.");
54 :
55 1494 : if (_use_Forchheimer_friction_model && !_speed)
56 0 : mooseError("If using a Forchheimer friction model, then the '",
57 : NS::speed,
58 : "' parameter must be provided");
59 1494 : }
60 :
61 : void
62 11115254 : PINSFVMomentumFrictionCorrection::gatherRCData(const FaceInfo & fi)
63 : {
64 11115254 : if (skipForBoundary(fi))
65 153920 : return;
66 :
67 10961334 : _face_info = &fi;
68 10961334 : _normal = fi.normal();
69 10961334 : _face_type = fi.faceType(std::make_pair(_var.number(), _var.sys().number()));
70 :
71 : using namespace Moose::FV;
72 :
73 10961334 : const auto elem_face = elemArg();
74 10961334 : const auto neighbor_face = neighborArg();
75 10961334 : const auto state = determineState();
76 :
77 10961334 : Point _face_centroid = _face_info->faceCentroid();
78 10961334 : Point _elem_centroid = _face_info->elemCentroid();
79 :
80 10961334 : ADReal friction_term_elem = 0;
81 10961334 : ADReal friction_term_neighbor = 0;
82 :
83 : ADReal diff_face;
84 :
85 : // If we are on an internal face for the variable, we interpolate the friction terms
86 : // to the face using the cell values
87 10961334 : if (_var.isInternalFace(*_face_info))
88 : {
89 10258290 : if (_use_Darcy_friction_model)
90 : {
91 20516580 : friction_term_elem += (*_D)(elem_face, state)(_index)*_mu(elem_face, state);
92 20516580 : friction_term_neighbor += (*_D)(neighbor_face, state)(_index)*_mu(neighbor_face, state);
93 : }
94 10258290 : if (_use_Forchheimer_friction_model)
95 : {
96 : friction_term_elem +=
97 20516580 : (*_F)(elem_face, state)(_index)*_rho(elem_face, state) / 2 * (*_speed)(elem_face, state);
98 10258290 : friction_term_neighbor += (*_F)(neighbor_face, state)(_index)*_rho(neighbor_face, state) / 2 *
99 20516580 : (*_speed)(neighbor_face, state);
100 : }
101 :
102 10258290 : Point _neighbor_centroid = _face_info->neighborCentroid();
103 :
104 10258290 : Real geometric_factor = _consistent_scaling * (_neighbor_centroid - _face_centroid).norm() *
105 10258290 : (_elem_centroid - _face_centroid).norm();
106 :
107 : // Compute the diffusion driven by the velocity gradient
108 : // Interpolate viscosity divided by porosity on the face
109 10258290 : interpolate(Moose::FV::InterpMethod::Average,
110 : diff_face,
111 10258290 : friction_term_elem * geometric_factor,
112 10258290 : friction_term_neighbor * geometric_factor,
113 : *_face_info,
114 : true);
115 : }
116 : // If not, we use the extrapolated values of the functors on the face
117 : else
118 : {
119 : const auto face =
120 703044 : makeFace(*_face_info, Moose::FV::limiterType(Moose::FV::InterpMethod::Average), true);
121 703044 : if (_use_Darcy_friction_model)
122 1406088 : friction_term_elem += (*_D)(face, state)(_index)*_mu(face, state);
123 703044 : if (_use_Forchheimer_friction_model)
124 : friction_term_elem +=
125 1406088 : (*_F)(face, state)(_index)*_rho(face, state) / 2 * (*_speed)(face, state);
126 :
127 : Real geometric_factor =
128 703044 : _consistent_scaling * std::pow((_elem_centroid - _face_centroid).norm(), 2);
129 :
130 703044 : diff_face = friction_term_elem * geometric_factor;
131 : }
132 :
133 : // Compute face superficial velocity gradient
134 21922668 : auto dudn = _var.gradient(makeCDFace(*_face_info), state) * _face_info->normal();
135 :
136 10961334 : if (_face_type == FaceInfo::VarFaceNeighbors::ELEM ||
137 : _face_type == FaceInfo::VarFaceNeighbors::BOTH)
138 : {
139 10961334 : const auto dof_number = _face_info->elem().dof_number(_sys.number(), _var.number(), 0);
140 : // A gradient is a linear combination of degrees of freedom so it's safe to straight-up index
141 : // into the derivatives vector at the dof we care about
142 10961334 : ADReal ae = dudn.derivatives()[dof_number];
143 10961334 : ae *= -diff_face;
144 21922668 : _rc_uo.addToA(&fi.elem(), _index, ae * (fi.faceArea() * fi.faceCoord()));
145 : }
146 10961334 : if (_face_type == FaceInfo::VarFaceNeighbors::NEIGHBOR ||
147 : _face_type == FaceInfo::VarFaceNeighbors::BOTH)
148 : {
149 10258290 : const auto dof_number = _face_info->neighbor().dof_number(_sys.number(), _var.number(), 0);
150 10258290 : ADReal an = dudn.derivatives()[dof_number];
151 10258290 : an *= diff_face;
152 30774870 : _rc_uo.addToA(fi.neighborPtr(), _index, an * (fi.faceArea() * fi.faceCoord()));
153 : }
154 :
155 10961334 : const auto strong_resid = -diff_face * dudn;
156 :
157 21922668 : addResidualAndJacobian(strong_resid * (fi.faceArea() * fi.faceCoord()));
158 : }
|
