pressureControl.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) 2017-2019 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 "pressureControl.H"
27 #include "findRefCell.H"
28 
29 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
30 
32 (
33  const volScalarField& p,
34  const volScalarField& rho,
35  const dictionary& dict,
36  const bool pRefRequired
37 )
38 :
39  refCell_(-1),
40  refValue_(0),
41  pMax_("pMax", dimPressure, great),
42  pMin_("pMin", dimPressure, 0),
43  limitMaxP_(false),
44  limitMinP_(false)
45 {
46  bool pLimits = false;
47  scalar pMax = -great;
48  scalar pMin = great;
49 
50  // Set the reference cell and value for closed domain simulations
51  if (pRefRequired && setRefCell(p, dict, refCell_, refValue_))
52  {
53  pLimits = true;
54 
55  pMax = refValue_;
56  pMin = refValue_;
57  }
58 
59  if (dict.found("pMax") && dict.found("pMin"))
60  {
61  pMax_.value() = readScalar(dict.lookup("pMax"));
62  limitMaxP_ = true;
63  pMin_.value() = readScalar(dict.lookup("pMin"));
64  limitMinP_ = true;
65  }
66  else
67  {
68  const volScalarField::Boundary& pbf = p.boundaryField();
69  const volScalarField::Boundary& rhobf = rho.boundaryField();
70 
71  scalar rhoRefMax = -great;
72  scalar rhoRefMin = great;
73  bool rhoLimits = false;
74 
75  forAll(pbf, patchi)
76  {
77  if (pbf[patchi].fixesValue())
78  {
79  pLimits = true;
80  rhoLimits = true;
81 
82  pMax = max(pMax, max(pbf[patchi]));
83  pMin = min(pMin, min(pbf[patchi]));
84 
85  rhoRefMax = max(rhoRefMax, max(rhobf[patchi]));
86  rhoRefMin = min(rhoRefMin, min(rhobf[patchi]));
87  }
88  }
89 
90  reduce(rhoLimits, andOp<bool>());
91  if (rhoLimits)
92  {
93  reduce(pMax, maxOp<scalar>());
94  reduce(pMin, minOp<scalar>());
95 
96  reduce(rhoRefMax, maxOp<scalar>());
97  reduce(rhoRefMin, minOp<scalar>());
98  }
99 
100  if (dict.found("pMax"))
101  {
102  pMax_.value() = readScalar(dict.lookup("pMax"));
103  limitMaxP_ = true;
104  }
105  else if (dict.found("pMaxFactor"))
106  {
107  if (!pLimits)
108  {
110  << "'pMaxFactor' specified rather than 'pMax'" << nl
111  << " but the corresponding reference pressure cannot"
112  " be evaluated from the boundary conditions." << nl
113  << " Please specify 'pMax' rather than 'pMaxFactor'"
114  << exit(FatalIOError);
115  }
116 
117  const scalar pMaxFactor(readScalar(dict.lookup("pMaxFactor")));
118  pMax_.value() = pMaxFactor*pMax;
119  limitMaxP_ = true;
120  }
121  else if (dict.found("rhoMax"))
122  {
123  // For backward-compatibility infer the pMax from rhoMax
124 
125  IOWarningInFunction(dict)
126  << "'rhoMax' specified rather than 'pMax' or 'pMaxFactor'"
127  << nl
128  << " This is supported for backward-compatibility but "
129  "'pMax' or 'pMaxFactor' are more reliable." << endl;
130 
131  if (!pLimits)
132  {
134  << "'rhoMax' specified rather than 'pMax'" << nl
135  << " but the corresponding reference pressure cannot"
136  " be evaluated from the boundary conditions." << nl
137  << " Please specify 'pMax' rather than 'rhoMax'"
138  << exit(FatalIOError);
139  }
140 
141  if (!rhoLimits)
142  {
144  << "'rhoMax' specified rather than 'pMaxFactor'" << nl
145  << " but the corresponding reference density cannot"
146  " be evaluated from the boundary conditions." << nl
147  << " Please specify 'pMaxFactor' rather than 'rhoMax'"
148  << exit(FatalIOError);
149  }
150 
151  dimensionedScalar rhoMax("rhoMax", dimDensity, dict);
152 
153  pMax_.value() = max(rhoMax.value()/rhoRefMax, 1)*pMax;
154  limitMaxP_ = true;
155  }
156 
157  if (dict.found("pMin"))
158  {
159  pMin_.value() = readScalar(dict.lookup("pMin"));
160  limitMinP_ = true;
161  }
162  else if (dict.found("pMinFactor"))
163  {
164  if (!pLimits)
165  {
167  << "'pMinFactor' specified rather than 'pMin'" << nl
168  << " but the corresponding reference pressure cannot"
169  " be evaluated from the boundary conditions." << nl
170  << " Please specify 'pMin' rather than 'pMinFactor'"
171  << exit(FatalIOError);
172  }
173 
174  const scalar pMinFactor(readScalar(dict.lookup("pMinFactor")));
175  pMin_.value() = pMinFactor*pMin;
176  limitMinP_ = true;
177  }
178  else if (dict.found("rhoMin"))
179  {
180  // For backward-compatibility infer the pMin from rhoMin
181 
182  IOWarningInFunction(dict)
183  << "'rhoMin' specified rather than 'pMin' or 'pMinFactor'" << nl
184  << " This is supported for backward-compatibility but"
185  "'pMin' or 'pMinFactor' are more reliable." << endl;
186 
187  if (!pLimits)
188  {
190  << "'rhoMin' specified rather than 'pMin'" << nl
191  << " but the corresponding reference pressure cannot"
192  " be evaluated from the boundary conditions." << nl
193  << " Please specify 'pMin' rather than 'rhoMin'"
194  << exit(FatalIOError);
195  }
196 
197  if (!rhoLimits)
198  {
200  << "'rhoMin' specified rather than 'pMinFactor'" << nl
201  << " but the corresponding reference density cannot"
202  " be evaluated from the boundary conditions." << nl
203  << " Please specify 'pMinFactor' rather than 'rhoMin'"
204  << exit(FatalIOError);
205  }
206 
207  dimensionedScalar rhoMin("rhoMin", dimDensity, dict);
208 
209  pMin_.value() = min(rhoMin.value()/rhoRefMin, 1)*pMin;
210  limitMinP_ = true;
211  }
212  }
213 
214  if (limitMaxP_ || limitMinP_)
215  {
216  Info<< "pressureControl" << nl;
217 
218  if (limitMaxP_)
219  {
220  Info<< " pMax " << pMax_.value() << nl;
221  }
222 
223  if (limitMinP_)
224  {
225  Info<< " pMin " << pMin_.value() << nl;
226  }
227 
228  Info << endl;
229  }
230 }
231 
232 
233 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
234 
236 {
237  if (limitMaxP_ || limitMinP_)
238  {
239  if (limitMaxP_)
240  {
241  const scalar pMax = max(p).value();
242 
243  if (pMax > pMax_.value())
244  {
245  Info<< "pressureControl: p max " << pMax << endl;
246  p = min(p, pMax_);
247  }
248  }
249 
250  if (limitMinP_)
251  {
252  const scalar pMin = min(p).value();
253 
254  if (pMin < pMin_.value())
255  {
256  Info<< "pressureControl: p min " << pMin << endl;
257  p = max(p, pMin_);
258  }
259  }
260 
262 
263  return true;
264  }
265  else
266  {
267  return false;
268  }
269 }
270 
271 
272 // ************************************************************************* //
#define readScalar
Definition: doubleScalar.C:38
bool found(const word &, bool recursive=false, bool patternMatch=true) const
Search dictionary for given keyword.
Definition: dictionary.C:438
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:158
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
const Boundary & boundaryField() const
Return const-reference to the boundary field.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:256
Find the reference cell nearest (in index) to the given cell but which is not on a cyclic...
const Type & value() const
Return const reference to value.
bool limit(volScalarField &p) const
Limit the pressure if necessary and return true if so.
const dimensionSet dimPressure
static const char nl
Definition: Ostream.H:265
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
dimensionedScalar pMin("pMin", dimPressure, fluid)
const dimensionSet dimDensity
label patchi
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
Definition: error.H:331
void correctBoundaryConditions()
Correct boundary field.
messageStream Info
#define IOWarningInFunction(ios)
Report an IO warning using Foam::Warning.
bool setRefCell(const volScalarField &field, const volScalarField &fieldRef, const dictionary &dict, label &refCelli, scalar &refValue, const bool forceReference=false)
If the field fieldRef needs referencing find the reference cell nearest.
Definition: findRefCell.C:31
pressureControl(const volScalarField &p, const volScalarField &rho, const dictionary &dict, const bool pRefRequired=true)
Construct from the simple/pimple sub-dictionary.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:583
IOerror FatalIOError