waveSurfacePressureFvPatchScalarField.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) 2018-2020 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
31#include "fvPatchFieldMapper.H"
32#include "volFields.H"
33#include "surfaceFields.H"
34#include "gravityMeshObject.H"
35#include "EulerDdtScheme.H"
37#include "backwardDdtScheme.H"
38
39// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
40
41const Foam::Enum
42<
44>
45Foam::waveSurfacePressureFvPatchScalarField::ddtSchemeTypeNames_
46({
47 {
48 ddtSchemeType::tsEuler,
49 fv::EulerDdtScheme<scalar>::typeName_()
50 },
51 {
52 ddtSchemeType::tsCrankNicolson,
53 fv::CrankNicolsonDdtScheme<scalar>::typeName_()
54 },
55 {
56 ddtSchemeType::tsBackward,
57 fv::backwardDdtScheme<scalar>::typeName_()
58 },
59});
60
61
62// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
63
64Foam::waveSurfacePressureFvPatchScalarField::
65waveSurfacePressureFvPatchScalarField
66(
67 const fvPatch& p,
69)
70:
71 fixedValueFvPatchScalarField(p, iF),
72 phiName_("phi"),
73 zetaName_("zeta"),
74 rhoName_("rho")
75{}
76
77
78Foam::waveSurfacePressureFvPatchScalarField::
79waveSurfacePressureFvPatchScalarField
80(
81 const fvPatch& p,
83 const dictionary& dict
84)
85:
86 fixedValueFvPatchScalarField(p, iF, dict),
87 phiName_(dict.getOrDefault<word>("phi", "phi")),
88 zetaName_(dict.getOrDefault<word>("zeta", "zeta")),
89 rhoName_(dict.getOrDefault<word>("rho", "rho"))
90{}
91
92
93Foam::waveSurfacePressureFvPatchScalarField::
94waveSurfacePressureFvPatchScalarField
95(
97 const fvPatch& p,
99 const fvPatchFieldMapper& mapper
100)
101:
102 fixedValueFvPatchScalarField(ptf, p, iF, mapper),
103 phiName_(ptf.phiName_),
104 zetaName_(ptf.zetaName_),
105 rhoName_(ptf.rhoName_)
106{}
107
108
109Foam::waveSurfacePressureFvPatchScalarField::
110waveSurfacePressureFvPatchScalarField
111(
113)
114:
115 fixedValueFvPatchScalarField(wspsf),
116 phiName_(wspsf.phiName_),
117 zetaName_(wspsf.zetaName_),
118 rhoName_(wspsf.rhoName_)
119{}
120
121
122Foam::waveSurfacePressureFvPatchScalarField::
123waveSurfacePressureFvPatchScalarField
124(
127)
128:
129 fixedValueFvPatchScalarField(wspsf, iF),
130 phiName_(wspsf.phiName_),
131 zetaName_(wspsf.zetaName_),
132 rhoName_(wspsf.rhoName_)
133{}
134
135
136// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
137
139{
140 if (updated())
141 {
142 return;
143 }
144
145 const label patchi = patch().index();
146
147 const scalar dt = db().time().deltaTValue();
148
149 // Retrieve non-const access to zeta field from the database
150 volVectorField& zeta = db().lookupObjectRef<volVectorField>(zetaName_);
151 vectorField& zetap = zeta.boundaryFieldRef()[patchi];
152
153 // Lookup d/dt scheme from database for zeta
154 const word ddtSchemeName(zeta.mesh().ddtScheme(zeta.name()));
155 ddtSchemeType ddtScheme(ddtSchemeTypeNames_[ddtSchemeName]);
156
157 // Retrieve the flux field from the database
158 const surfaceScalarField& phi =
159 db().lookupObject<surfaceScalarField>(phiName_);
160
161 // Cache the patch face-normal vectors
162 tmp<vectorField> nf(patch().nf());
163
164 // Change in zeta due to flux
165 vectorField dZetap(dt*nf()*phi.boundaryField()[patchi]/patch().magSf());
166
168 {
169 const scalarField& rhop =
170 patch().lookupPatchField<volScalarField, scalar>(rhoName_);
171
172 dZetap /= rhop;
173 }
174
175 const volVectorField& zeta0 = zeta.oldTime();
176
177 switch (ddtScheme)
178 {
179 case tsEuler:
180 case tsCrankNicolson:
181 {
182 zetap = zeta0.boundaryField()[patchi] + dZetap;
183
184 break;
185 }
186 case tsBackward:
187 {
188 scalar dt0 = db().time().deltaT0Value();
189
190 scalar c = 1.0 + dt/(dt + dt0);
191 scalar c00 = dt*dt/(dt0*(dt + dt0));
192 scalar c0 = c + c00;
193
194 zetap =
195 (
196 c0*zeta0.boundaryField()[patchi]
197 - c00*zeta0.oldTime().boundaryField()[patchi]
198 + dZetap
199 )/c;
200
201 break;
202 }
203 default:
204 {
206 << ddtSchemeName << nl
207 << " on patch " << this->patch().name()
208 << " of field " << this->internalField().name()
209 << " in file " << this->internalField().objectPath()
210 << abort(FatalError);
211 }
212 }
213
214
215 Info<< "min/max zetap = " << gMin(zetap & nf()) << ", "
216 << gMax(zetap & nf()) << endl;
217
218 // Update the surface pressure
220 meshObjects::gravity::New(db().time());
221
222 operator==(-g.value() & zetap);
223
224 fixedValueFvPatchScalarField::updateCoeffs();
225}
226
227
229{
231 os.writeEntryIfDifferent<word>("phi", "phi", phiName_);
232 os.writeEntryIfDifferent<word>("zeta", "zeta", zetaName_);
233 os.writeEntryIfDifferent<word>("rho", "rho", rhoName_);
234 writeEntry("value", os);
235}
236
237
238// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
239
240namespace Foam
241{
243 (
246 );
247}
248
249// ************************************************************************* //
Macros for easy insertion into run-time selection tables.
const uniformDimensionedVectorField & g
surfaceScalarField & phi
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
const dimensionSet & dimensions() const
Return dimensions.
const Mesh & mesh() const
Return mesh.
Enum is a wrapper around a list of names/values that represent particular enumeration (or int) values...
Definition: Enum.H:61
const GeometricField< Type, PatchField, GeoMesh > & oldTime() const
Return old time field.
Boundary & boundaryFieldRef(const bool updateAccessTime=true)
Return a reference to the boundary field.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
const word & name() const noexcept
Return the object name.
Definition: IOobjectI.H:65
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:62
Ostream & writeEntryIfDifferent(const word &key, const T &value1, const T &value2)
Write a keyword/value entry only when the two values differ.
Definition: Ostream.H:251
static autoPtr< Time > New()
Construct (dummy) Time - no functionObjects or libraries.
Definition: Time.C:717
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:126
const Type & value() const
Return const reference to value.
virtual bool write()
Write the output fields.
A FieldMapper for finite-volume patch fields.
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:71
A class for managing temporary objects.
Definition: tmp.H:65
This is a pressure boundary condition, whose value is calculated as the hydrostatic pressure based on...
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
ddtSchemeType
Enumeration defining the available ddt schemes.
A class for handling words, derived from Foam::string.
Definition: word.H:68
volScalarField & p
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
OBJstream os(runTime.globalPath()/outputName)
#define makePatchTypeField(PatchTypeField, typePatchTypeField)
Definition: fvPatchField.H:676
Namespace for OpenFOAM.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:82
const dimensionSet dimArea(sqr(dimLength))
Definition: dimensionSets.H:59
const dimensionSet dimVelocity
messageStream Info
Information stream (stdout output on master, null elsewhere)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
tmp< faMatrix< Type > > operator==(const faMatrix< Type > &, const faMatrix< Type > &)
errorManip< error > abort(error &err)
Definition: errorManip.H:144
error FatalError
Type gMin(const FieldField< Field, Type > &f)
const dimensionSet dimDensity
Type gMax(const FieldField< Field, Type > &f)
constexpr char nl
The newline '\n' character (0x0a)
Definition: Ostream.H:53
dictionary dict
Foam::surfaceFields.