wallShearStress.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) 2013-2015 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 "wallShearStress.H"
27 #include "volFields.H"
28 #include "surfaceFields.H"
31 #include "wallPolyPatch.H"
32 
33 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
34 
35 namespace Foam
36 {
37  defineTypeNameAndDebug(wallShearStress, 0);
38 }
39 
40 
41 // * * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * //
42 
44 {
45  // Add headers to output data
46  writeHeader(file(), "Wall shear stress");
47  writeCommented(file(), "Time");
48  writeTabbed(file(), "patch");
49  writeTabbed(file(), "min");
50  writeTabbed(file(), "max");
51  file() << endl;
52 }
53 
54 
56 (
57  const fvMesh& mesh,
58  const volSymmTensorField& Reff,
59  volVectorField& shearStress
60 )
61 {
62  forAllConstIter(labelHashSet, patchSet_, iter)
63  {
64  label patchI = iter.key();
65  const polyPatch& pp = mesh.boundaryMesh()[patchI];
66 
67  vectorField& ssp = shearStress.boundaryField()[patchI];
68  const vectorField& Sfp = mesh.Sf().boundaryField()[patchI];
69  const scalarField& magSfp = mesh.magSf().boundaryField()[patchI];
70  const symmTensorField& Reffp = Reff.boundaryField()[patchI];
71 
72  ssp = (-Sfp/magSfp) & Reffp;
73 
74  vector minSsp = gMin(ssp);
75  vector maxSsp = gMax(ssp);
76 
77  if (Pstream::master())
78  {
79  file() << mesh.time().value()
80  << token::TAB << pp.name()
81  << token::TAB << minSsp
82  << token::TAB << maxSsp
83  << endl;
84  }
85 
86  if (log_) Info<< " min/max(" << pp.name() << ") = "
87  << minSsp << ", " << maxSsp << endl;
88  }
89 }
90 
91 
92 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
93 
95 (
96  const word& name,
97  const objectRegistry& obr,
98  const dictionary& dict,
99  const bool loadFromFiles
100 )
101 :
102  functionObjectFile(obr, name, typeName),
103  name_(name),
104  obr_(obr),
105  active_(true),
106  log_(true),
107  patchSet_()
108 {
109  // Check if the available mesh is an fvMesh, otherwise deactivate
110  if (!isA<fvMesh>(obr_))
111  {
112  active_ = false;
113  WarningIn
114  (
115  "wallShearStress::wallShearStress"
116  "("
117  "const word&, "
118  "const objectRegistry&, "
119  "const dictionary&, "
120  "const bool"
121  ")"
122  ) << "No fvMesh available, deactivating " << name_ << nl
123  << endl;
124  }
125 
126  if (active_)
127  {
128  const fvMesh& mesh = refCast<const fvMesh>(obr_);
129 
130  volVectorField* wallShearStressPtr
131  (
132  new volVectorField
133  (
134  IOobject
135  (
136  type(),
137  mesh.time().timeName(),
138  mesh,
141  ),
142  mesh,
144  (
145  "0",
148  )
149  )
150  );
151 
152  mesh.objectRegistry::store(wallShearStressPtr);
153  }
154 
155  read(dict);
156 }
157 
158 
159 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
160 
162 {}
163 
164 
165 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
166 
168 {
169  if (active_)
170  {
171  log_ = dict.lookupOrDefault<Switch>("log", true);
172 
173  const fvMesh& mesh = refCast<const fvMesh>(obr_);
174  const polyBoundaryMesh& pbm = mesh.boundaryMesh();
175 
176  patchSet_ =
177  mesh.boundaryMesh().patchSet
178  (
179  wordReList(dict.lookupOrDefault("patches", wordReList()))
180  );
181 
182  Info<< type() << " " << name_ << ":" << nl;
183 
184  if (patchSet_.empty())
185  {
186  forAll(pbm, patchI)
187  {
188  if (isA<wallPolyPatch>(pbm[patchI]))
189  {
190  patchSet_.insert(patchI);
191  }
192  }
193 
194  Info<< " processing all wall patches" << nl << endl;
195  }
196  else
197  {
198  Info<< " processing wall patches: " << nl;
199  labelHashSet filteredPatchSet;
200  forAllConstIter(labelHashSet, patchSet_, iter)
201  {
202  label patchI = iter.key();
203  if (isA<wallPolyPatch>(pbm[patchI]))
204  {
205  filteredPatchSet.insert(patchI);
206  Info<< " " << pbm[patchI].name() << endl;
207  }
208  else
209  {
210  WarningIn("void wallShearStress::read(const dictionary&)")
211  << "Requested wall shear stress on non-wall boundary "
212  << "type patch: " << pbm[patchI].name() << endl;
213  }
214  }
215 
216  Info<< endl;
217 
218  patchSet_ = filteredPatchSet;
219  }
220  }
221 }
222 
223 
225 {
226  typedef compressible::turbulenceModel cmpModel;
227  typedef incompressible::turbulenceModel icoModel;
228 
229  if (active_)
230  {
232 
233  const fvMesh& mesh = refCast<const fvMesh>(obr_);
234 
236  const_cast<volVectorField&>
237  (
239  );
240 
241  if (log_) Info<< type() << " " << name_ << " output:" << nl;
242 
243 
245  if (mesh.foundObject<cmpModel>(turbulenceModel::propertiesName))
246  {
247  const cmpModel& model =
249 
250  Reff = model.devRhoReff();
251  }
252  else if (mesh.foundObject<icoModel>(turbulenceModel::propertiesName))
253  {
254  const icoModel& model =
256 
257  Reff = model.devReff();
258  }
259  else
260  {
261  FatalErrorIn("void Foam::wallShearStress::execute()")
262  << "Unable to find turbulence model in the "
263  << "database" << exit(FatalError);
264  }
265 
266  calcShearStress(mesh, Reff(), wallShearStress);
267  }
268 }
269 
270 
272 {
273  if (active_)
274  {
275  execute();
276  }
277 }
278 
279 
281 {
282  // Do nothing
283 }
284 
285 
287 {
288  if (active_)
289  {
291 
293  obr_.lookupObject<volVectorField>(type());
294 
295  if (log_) Info<< type() << " " << name_ << " output:" << nl
296  << " writing field " << wallShearStress.name() << nl
297  << endl;
298 
299  wallShearStress.write();
300  }
301 }
302 
303 
304 // ************************************************************************* //
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
Foam::surfaceFields.
virtual void read(const dictionary &)
Read the wallShearStress data.
const word & name() const
Return name.
GeometricBoundaryField & boundaryField()
Return reference to GeometricBoundaryField.
wallShearStress(const wallShearStress &)
Disallow default bitwise copy construct.
const surfaceVectorField & Sf() const
Return cell face area vectors.
bool foundObject(const word &name) const
Is the named Type found?
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Definition: dimensionSets.H:51
Templated abstract base class for single-phase incompressible turbulence models.
List< wordRe > wordReList
A List of wordRe (word or regular expression)
Definition: wordReList.H:50
A class for handling words, derived from string.
Definition: word.H:59
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
static word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
Definition: Time.C:741
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
A simple wrapper around bool so that it can be read as a word: true/false, on/off, yes/no, y/n, t/f, or none.
Definition: Switch.H:60
messageStream Info
dynamicFvMesh & mesh
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
virtual ~wallShearStress()
Destructor.
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:50
This function object evaluates and outputs the shear stress at wall patches. The result is written as...
Namespace for OpenFOAM.
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:66
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:243
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:39
virtual bool write() const
Write using setting from DB.
static const char nl
Definition: Ostream.H:260
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
const surfaceScalarField & magSf() const
Return cell face area magnitudes.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:91
#define WarningIn(functionName)
Report a warning using Foam::Warning.
Foam::polyBoundaryMesh.
#define forAll(list, i)
Definition: UList.H:421
const Type & lookupObject(const word &name) const
Lookup and return the object of the given Type.
void calcShearStress(const fvMesh &mesh, const volSymmTensorField &Reff, volVectorField &shearStress)
Calculate the shear stress.
virtual void end()
Execute at the final time-loop, currently does nothing.
const word & name() const
Return name.
Definition: IOobject.H:260
Type gMin(const FieldField< Field, Type > &f)
fileName::Type type(const fileName &)
Return the file type: DIRECTORY or FILE.
Definition: POSIX.C:589
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:421
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
virtual void write()
Write function.
labelHashSet patchSet(const UList< wordRe > &patchNames, const bool warnNotFound=true, const bool usePatchGroups=true) const
Return the set of patch IDs corresponding to the given names.
virtual void execute()
Execute, currently does nothing.
virtual void write()
Calculate the wallShearStress and write.
error FatalError
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
Registry of regIOobjects.
static const Vector zero
Definition: Vector.H:80
static const word propertiesName
Default name of the turbulence properties dictionary.
Base class for output file data handling.
static bool master(const label communicator=0)
Am I the master process.
Definition: UPstream.H:398
bool read(const char *, int32_t &)
Definition: int32IO.C:87
dimensionedSymmTensor sqr(const dimensionedVector &dv)
Type gMax(const FieldField< Field, Type > &f)
Templated wrapper class to provide compressible turbulence models thermal diffusivity based thermal t...
A class for managing temporary objects.
Definition: PtrList.H:118
const Type & value() const
Return const reference to value.
defineTypeNameAndDebug(combustionModel, 0)
virtual void timeSet()
Called when time was set at the end of the Time::operator++.
virtual void writeFileHeader(const label i)
File header information.
bool insert(const Key &key)
Insert a new entry.
Definition: HashSet.H:116