adjustPhi.C
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) 2011-2016 OpenFOAM Foundation
9 Copyright (C) 2019 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
12 This file is part of OpenFOAM.
13
14 OpenFOAM is free software: you can redistribute it and/or modify it
15 under the terms of the GNU General Public License as published by
16 the Free Software Foundation, either version 3 of the License, or
17 (at your option) any later version.
18
19 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22 for more details.
23
24 You should have received a copy of the GNU General Public License
25 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26
27\*---------------------------------------------------------------------------*/
28
29#include "adjustPhi.H"
30#include "volFields.H"
31#include "surfaceFields.H"
33
34// * * * * * * * * * * * * * * * Global Functions * * * * * * * * * * * * * //
35
37(
39 const volVectorField& U,
41)
42{
43 if (p.needReference())
44 {
45 scalar massIn = 0.0;
46 scalar fixedMassOut = 0.0;
47 scalar adjustableMassOut = 0.0;
48
50 phi.boundaryFieldRef();
51
52 forAll(bphi, patchi)
53 {
54 const fvPatchVectorField& Up = U.boundaryField()[patchi];
55 const fvsPatchScalarField& phip = bphi[patchi];
56
57 if (!phip.coupled())
58 {
59 if (Up.fixesValue() && !isA<inletOutletFvPatchVectorField>(Up))
60 {
61 forAll(phip, i)
62 {
63 if (phip[i] < 0.0)
64 {
65 massIn -= phip[i];
66 }
67 else
68 {
69 fixedMassOut += phip[i];
70 }
71 }
72 }
73 else
74 {
75 forAll(phip, i)
76 {
77 if (phip[i] < 0.0)
78 {
79 massIn -= phip[i];
80 }
81 else
82 {
83 adjustableMassOut += phip[i];
84 }
85 }
86 }
87 }
88 }
89
90 // Calculate the total flux in the domain, used for normalisation
91 scalar totalFlux = VSMALL + sum(mag(phi)).value();
92
93 reduce(massIn, sumOp<scalar>());
94 reduce(fixedMassOut, sumOp<scalar>());
95 reduce(adjustableMassOut, sumOp<scalar>());
96
97 scalar massCorr = 1.0;
98 scalar magAdjustableMassOut = mag(adjustableMassOut);
99
100 if
101 (
102 magAdjustableMassOut > VSMALL
103 && magAdjustableMassOut/totalFlux > SMALL
104 )
105 {
106 massCorr = (massIn - fixedMassOut)/adjustableMassOut;
107 }
108 else if (mag(fixedMassOut - massIn)/totalFlux > 1e-8)
109 {
111 << "Continuity error cannot be removed by adjusting the"
112 " outflow.\nPlease check the velocity boundary conditions"
113 " and/or run potentialFoam to initialise the outflow." << nl
114 << "Total flux : " << totalFlux << nl
115 << "Specified mass inflow : " << massIn << nl
116 << "Specified mass outflow : " << fixedMassOut << nl
117 << "Adjustable mass outflow : " << adjustableMassOut << nl
118 << exit(FatalError);
119 }
120
121 forAll(bphi, patchi)
122 {
123 const fvPatchVectorField& Up = U.boundaryField()[patchi];
124 fvsPatchScalarField& phip = bphi[patchi];
125
126 if (!phip.coupled())
127 {
128 if
129 (
130 !Up.fixesValue()
131 || isA<inletOutletFvPatchVectorField>(Up)
132 )
133 {
134 forAll(phip, i)
135 {
136 if (phip[i] > 0.0)
137 {
138 phip[i] *= massCorr;
139 }
140 }
141 }
142 }
143 }
144
145 return mag(massIn)/totalFlux < SMALL
146 && mag(fixedMassOut)/totalFlux < SMALL
147 && mag(adjustableMassOut)/totalFlux < SMALL;
148 }
149
150 return false;
151}
152
153
154// ************************************************************************* //
For cases which do no have a pressure boundary adjust the balance of fluxes to obey continuity....
reduce(hasMovingMesh, orOp< bool >())
surfaceScalarField & phi
virtual bool fixesValue() const
Return true if this patch field fixes a value.
Definition: fvPatchField.H:337
virtual bool coupled() const
Return true if this patch field is coupled.
U
Definition: pEqn.H:72
volScalarField & p
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
bool adjustPhi(surfaceScalarField &phi, const volVectorField &U, volScalarField &p)
Adjust the balance of fluxes to obey continuity.
Definition: adjustPhi.C:37
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
error FatalError
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:130
constexpr char nl
The newline '\n' character (0x0a)
Definition: Ostream.H:53
volScalarField & e
Definition: createFields.H:11
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333
Foam::surfaceFields.