directionalPressureGradientExplicitSource.H
Go to the documentation of this file.
1/*---------------------------------------------------------------------------*\
2 ========= |
3 \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4 \\ / O peration |
5 \\ / A nd | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2015-2020 OpenCFD Ltd.
9-------------------------------------------------------------------------------
10License
11 This file is part of OpenFOAM.
12
13 OpenFOAM is free software: you can redistribute it and/or modify it
14 under the terms of the GNU General Public License as published by
15 the Free Software Foundation, either version 3 of the License, or
16 (at your option) any later version.
17
18 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
19 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21 for more details.
22
23 You should have received a copy of the GNU General Public License
24 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
25
26Class
27 Foam::fv::directionalPressureGradientExplicitSource
28
29Group
30 grpFvOptionsSources
31
32Description
33 Applies an explicit pressure gradient source in such a way to deflect the
34 flow towards an specific direction (\c flowDir). Alternatively add an extra
35 pressure drop in the \c flowDir direction using a model.
36
37Usage
38 Minimal example by using \c constant/fvOptions:
39 \verbatim
40 directionalPressureGradientExplicitSource1
41 {
42 // Mandatory entries (unmodifiable)
43 type directionalPressureGradientExplicitSource;
44
45 // Mandatory entries (unmodifiable)
46 model <modelName>;
47 fields (<fieldName>);
48
49 // Mandatory entries (runtime modifiable)
50 flowDir (1 1 0);
51 faceZone <faceZoneName>;
52
53 // Conditional mandatory entries (unmodifiable)
54
55 // when <timePath>/uniform/<name>Properties file exists
56 gradient <vectorField>; // reading from the aforementioned file
57
58 // when model=DarcyForchheimer
59 // deltaP = (D*mu + 0.5*rho*magUn)*magUn*length
60 D 5e7;
61 I 0;
62 length 1e-3;
63
64 // when model=constant
65 pressureDrop 40;
66
67 // when model=volumetricFlowRateTable
68 outOfBounds clamp;
69 fileName "volFlowRateTable";
70
71 // Optional entries (runtime modifiable)
72 relaxationFactor 0.3;
73
74 // Mandatory/Optional (inherited) entries
75 ...
76 }
77 \endverbatim
78
79 where the entries mean:
80 \table
81 Property | Description | Type | Reqd | Dflt
82 type | Type name: directionalPressureGradientExplicitSource <!--
83 --> | word | yes | -
84 model | Pressure drop model [Pa] | word | yes | -
85 fields | Name of operand field | word | yes | -
86 gradient | Initial pressure gradient field | vectorField | cndtnl | -
87 flowDir | Deflection flow direction | vector | yes | -
88 faceZone | Name of upstream faceZone | word | yes | -
89 relaxationFactor | Relaxation factor for flow deflection | scalar | no | 0.3
90 D | Darcy pressure loss coefficient | scalar | cndtnl | -
91 I | Inertia pressure lost coefficient | scalar | cndtnl | -
92 length | Porous media length | scalar | cndtnl | -
93 presssureDrop | Constant pressure drop | scalar | cndtnl | -
94 fileName | Interpolation table for volumetric flow rate <!--
95 --> | interpolationTable | cndtnl | -
96 \endtable
97
98 Options for the \c model entry:
99 \verbatim
100 volumetricFlowRateTable | Pressure-gradient file
101 constant | Constant pressure drop
102 DarcyForchheimer | Darcy-Forchheimer model
103 \endverbatim
104
105 The inherited entries are elaborated in:
106 - \link fvOption.H \endlink
107 - \link cellSetOption.H \endlink
108
109Note
110 - In order to obtain the upwind velocities this function loops over
111 the slaves cells of the faceZone specified in the dictionary, on the other
112 hand, the cellZone to which this source term is applied should be composed
113 of the master cells and they should be 'downwind' the faceZone.
114
115SourceFiles
116 directionalPressureGradientExplicitSource.C
117
118\*---------------------------------------------------------------------------*/
119
120#ifndef directionalPressureGradientExplicitSource_H
121#define directionalPressureGradientExplicitSource_H
122
123#include "autoPtr.H"
124#include "fvMesh.H"
125#include "volFields.H"
126#include "fvOption.H"
127#include "cellSetOption.H"
128#include "interpolationTable.H"
129
130// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
131
132namespace Foam
133{
134namespace fv
135{
136
137/*---------------------------------------------------------------------------*\
138 Class directionalPressureGradientExplicitSource Declaration
139\*---------------------------------------------------------------------------*/
140
141class directionalPressureGradientExplicitSource
142:
143 public fv::cellSetOption
144{
145public:
146
147 // Public Enumeration
148
149 //- Modes of pressure drop
151 {
153 pConstant,
155 };
156
157
158private:
159
160 // Private Data
161
162 static const Enum<pressureDropModel> pressureDropModelNames_;
163
164 //- Pressure drop model
165 pressureDropModel model_;
166
167 //- Pressure gradient before correction
168 vectorField gradP0_;
169
170 //- Change in pressure gradient
171 vectorField dGradP_;
172
173 //- Pressure drop due to porous media
174 vectorField gradPporous_;
175
176 //- Flow direction
177 vector flowDir_;
178
179 //- Matrix 1/A coefficients field pointer
180 autoPtr<volScalarField> invAPtr_;
181
182 //- Darcy pressure loss coefficient
183 scalar D_;
184
185 //- Inertia pressure lost coefficient
186 scalar I_;
187
188 //- Porous media length
189 scalar length_;
190
191 //- Constant pressure drop
192 scalar pressureDrop_;
193
194 //- Volumetric flow rate vs pressure drop table
195 interpolationTable<scalar> flowRate_;
196
197 //- Name of the faceZone at the heat exchange inlet
198 word faceZoneName_;
199
200 //- Id for the face zone
201 label zoneID_;
202
203 //- Local list of face IDs
204 labelList faceId_;
205
206 //- Local list of patch ID per face
207 labelList facePatchId_;
208
209 //- Relaxation factor
210 scalar relaxationFactor_;
211
212 //- Cells faces mapping
213 labelList cellFaceMap_;
214
215
216 // Private Member Functions
217
218 //- Init
219 void initialise();
220
221 //- Write the pressure gradient to file (for restarts etc)
222 void writeProps(const vectorField& gradP) const;
223
224 //- Correct driving force for a constant mass flow rate
225 void update(fvMatrix<vector>& eqn);
226
227
228public:
229
230 //- Runtime type information
231 TypeName("directionalPressureGradientExplicitSource");
232
234 // Constructors
235
236 //- Construct from explicit source name and mesh
238 (
239 const word& sourceName,
240 const word& modelType,
241 const dictionary& dict,
242 const fvMesh& mesh
243 );
244
245 //- No copy construct
247 (
249 ) = delete;
250
251 //- No copy assignment
252 void operator=
253 (
255 ) = delete;
256
257
258 //- Destructor
260
261
262 // Member Functions
263
264 // Evaluation
265
266 //- Correct the pressure gradient
267 virtual void correct(volVectorField& U);
268
269 //- Add explicit contribution to momentum equation
270 virtual void addSup
271 (
272 fvMatrix<vector>& eqn,
273 const label fieldI
274 );
275
276 //- Add explicit contribution to compressible momentum equation
277 virtual void addSup
278 (
279 const volScalarField& rho,
280 fvMatrix<vector>& eqn,
281 const label fieldI
282 );
283
284 //- Set 1/A coefficient
285 virtual void constrain
286 (
287 fvMatrix<vector>& eqn,
288 const label fieldI
289 );
290
291
292 // I-O
293
294 //- Write the source properties
295 virtual void writeData(Ostream& os) const;
296
297 //- Read source dictionary
298 virtual bool read(const dictionary& dict);
299};
300
301
302// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
303
304} // End namespace fv
305} // End namespace Foam
306
307// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
308
309#endif
310
311// ************************************************************************* //
Enum is a wrapper around a list of names/values that represent particular enumeration (or int) values...
Definition: Enum.H:61
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:62
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
Definition: autoPtr.H:66
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:126
A special matrix type and solver, designed for finite volume solutions of scalar equations....
Definition: fvMatrix.H:121
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:91
Intermediate abstract class for handling cell-set options for the derived fvOptions.
Applies an explicit pressure gradient source in such a way to deflect the flow towards an specific di...
virtual void writeData(Ostream &os) const
Write the source properties.
virtual bool read(const dictionary &dict)
Read source dictionary.
virtual ~directionalPressureGradientExplicitSource()=default
Destructor.
virtual void addSup(fvMatrix< vector > &eqn, const label fieldI)
Add explicit contribution to momentum equation.
virtual void constrain(fvMatrix< vector > &eqn, const label fieldI)
Set 1/A coefficient.
directionalPressureGradientExplicitSource(const directionalPressureGradientExplicitSource &)=delete
No copy construct.
TypeName("directionalPressureGradientExplicitSource")
Runtime type information.
const fvMesh & mesh() const noexcept
Return const access to the mesh database.
Definition: fvOptionI.H:37
An interpolation/look-up table of scalar vs <Type> values. The reference scalar values must be monoto...
A Vector of values with scalar precision, where scalar is float/double depending on the compilation f...
A class for handling words, derived from Foam::string.
Definition: word.H:68
U
Definition: pEqn.H:72
thermo correct()
volVectorField gradP(fvc::grad(p))
mesh update()
OBJstream os(runTime.globalPath()/outputName)
Namespace for OpenFOAM.
List< label > labelList
A List of labels.
Definition: List.H:66
Field< vector > vectorField
Specialisation of Field<T> for vector.
labelList fv(nPoints)
dictionary dict
#define TypeName(TypeNameString)
Declare a ClassName() with extra virtual type info.
Definition: typeInfo.H:73