wallDampingModel.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2015-2023 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 \*---------------------------------------------------------------------------*/
25 
26 #include "wallDampingModel.H"
27 #include "surfaceInterpolate.H"
28 #include "wallFvPatch.H"
29 
30 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
31 
32 namespace Foam
33 {
36 }
37 
38 
39 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
40 
42 (
43  const dictionary& dict,
44  const phaseInterface& interface
45 )
46 :
47  wallDependentModel(interface.mesh()),
48  interface_
49  (
50  interface.modelCast<wallDampingModel, dispersedPhaseInterface>()
51  ),
52  Cd_("Cd", dimless, dict),
53  zeroWallDist_
54  (
55  dimensionedScalar::lookupOrDefault
56  (
57  "zeroWallDist",
58  dict,
59  dimLength,
60  0
61  )
62  ),
63  zeroInNearWallCells_
64  (
65  dict.lookupOrDefault<Switch>("zeroInNearWallCells", false)
66  )
67 {}
68 
69 
70 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
71 
73 {}
74 
75 
76 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
77 
80 {
81  tmp<volScalarField> tlimiter(limiter().ptr());
82 
83  if (zeroInNearWallCells_)
84  {
85  volScalarField& limiter = tlimiter.ref();
86 
87  const fvBoundaryMesh& bMesh = limiter.mesh().boundary();
88 
90  {
91  if (isA<wallFvPatch>(bMesh[patchi]))
92  {
93  const labelUList& faceCells = bMesh[patchi].faceCells();
94 
95  forAll(faceCells, facei)
96  {
97  limiter[faceCells[facei]] = 0;
98  }
99  }
100  }
101 
102  return tlimiter;
103  }
104  else
105  {
106  return tlimiter;
107  }
108 }
109 
110 
113 {
114  return fvc::interpolate(damping());
115 }
116 
117 
118 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
Generic GeometricField class.
A simple wrapper around bool so that it can be read as a word: true/false, on/off,...
Definition: Switch.H:61
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:162
Class to represent a interface between phases where one phase is considered dispersed within the othe...
Foam::fvBoundaryMesh.
Class to represent an interface between phases. Derivations can further specify the configuration of ...
A class for managing temporary objects.
Definition: tmp.H:55
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:181
Wall damping models can be used to filter interfacial models near the walls. This is particularly use...
virtual tmp< surfaceScalarField > dampingf() const
Return damped face coefficient.
virtual ~wallDampingModel()
Destructor.
virtual tmp< volScalarField > damping() const
Return damped coefficient.
wallDampingModel(const dictionary &dict, const phaseInterface &interface)
Construct from a dictionary and an interface.
A class which provides on-demand creation and caching of wall distance and wall normal fields for use...
const polyBoundaryMesh & bMesh
label patchi
void limiter(surfaceScalarField &lambda, const RdeltaTType &rDeltaT, const RhoType &rho, const volScalarField &psi, const surfaceScalarField &phiBD, const surfaceScalarField &phiCorr, const SpType &Sp, const SuType &Su, const PsiMaxType &psiMax, const PsiMinType &psiMin)
static tmp< SurfaceField< Type > > interpolate(const VolField< Type > &tvf, const surfaceScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
Namespace for OpenFOAM.
defineRunTimeSelectionTable(reactionRateFlameArea, dictionary)
const dimensionSet dimless
const dimensionSet dimLength
defineTypeNameAndDebug(combustionModel, 0)
dictionary dict