fixedFluxPressureFvPatchScalarField.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 | Copyright (C) 2011-2013 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 
27 #include "fvPatchFieldMapper.H"
28 #include "volFields.H"
29 #include "surfaceFields.H"
31 
32 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
33 
35 (
36  const fvPatch& p,
38 )
39 :
40  fixedGradientFvPatchScalarField(p, iF),
41  curTimeIndex_(-1)
42 {}
43 
44 
46 (
48  const fvPatch& p,
50  const fvPatchFieldMapper& mapper
51 )
52 :
53  fixedGradientFvPatchScalarField(p, iF),
54  curTimeIndex_(-1)
55 {
56  patchType() = ptf.patchType();
57 
58  // Map gradient. Set unmapped values and overwrite with mapped ptf
59  gradient() = 0.0;
60  gradient().map(ptf.gradient(), mapper);
61 
62  // Evaluate the value field from the gradient if the internal field is valid
63  if (notNull(iF) && iF.size())
64  {
65  scalarField::operator=
66  (
67  //patchInternalField() + gradient()/patch().deltaCoeffs()
68  // ***HGW Hack to avoid the construction of mesh.deltaCoeffs
69  // which fails for AMI patches for some mapping operations
70  patchInternalField() + gradient()*(patch().nf() & patch().delta())
71  );
72  }
73 }
74 
75 
77 (
78  const fvPatch& p,
80  const dictionary& dict
81 )
82 :
83  fixedGradientFvPatchScalarField(p, iF),
84  curTimeIndex_(-1)
85 {
86  if (dict.found("value") && dict.found("gradient"))
87  {
89  (
90  scalarField("value", dict, p.size())
91  );
92  gradient() = scalarField("gradient", dict, p.size());
93  }
94  else
95  {
96  fvPatchField<scalar>::operator=(patchInternalField());
97  gradient() = 0.0;
98  }
99 }
100 
101 
103 (
105 )
106 :
107  fixedGradientFvPatchScalarField(wbppsf),
108  curTimeIndex_(-1)
109 {}
110 
111 
113 (
116 )
117 :
118  fixedGradientFvPatchScalarField(wbppsf, iF),
119  curTimeIndex_(-1)
120 {}
121 
122 
123 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
124 
126 (
127  const scalarField& snGradp
128 )
129 {
130  if (updated())
131  {
132  return;
133  }
134 
135  curTimeIndex_ = this->db().time().timeIndex();
136 
137  gradient() = snGradp;
138  fixedGradientFvPatchScalarField::updateCoeffs();
139 }
140 
141 
143 {
144  if (updated())
145  {
146  return;
147  }
148 
149  if (curTimeIndex_ != this->db().time().timeIndex())
150  {
151  FatalErrorIn("fixedFluxPressureFvPatchScalarField::updateCoeffs()")
152  << "updateCoeffs(const scalarField& snGradp) MUST be called before"
153  " updateCoeffs() or evaluate() to set the boundary gradient."
154  << exit(FatalError);
155  }
156 }
157 
158 
160 {
162  writeEntry("value", os);
163 }
164 
165 
166 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
167 
168 namespace Foam
169 {
171  (
174  );
175 }
176 
177 
178 // ************************************************************************* //
bool found(const word &, bool recursive=false, bool patternMatch=true) const
Search dictionary for given keyword.
Definition: dictionary.C:306
Foam::surfaceFields.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
label timeIndex
Definition: getTimeIndex.H:4
virtual void updateCoeffs()
Update the patch pressure gradient field.
void size(const label)
Override size to be inconsistent with allocated storage.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
Foam::fvPatchFieldMapper.
bool notNull(const T &t)
Return true if t is not a reference to the nullObject of type T.
Definition: nullObjectI.H:46
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
Namespace for OpenFOAM.
runTime write()
Macros for easy insertion into run-time selection tables.
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:65
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
This boundary condition sets the pressure gradient to the provided value such that the flux on the bo...
error FatalError
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:61
fixedFluxPressureFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
makePatchTypeField(fvPatchVectorField, SRFFreestreamVelocityFvPatchVectorField)
volScalarField scalarField(fieldObject, mesh)
virtual label size() const
Return size.
Definition: fvPatch.H:161