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-2017 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"
33 
34 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
35 
36 namespace Foam
37 {
38 namespace functionObjects
39 {
40  defineTypeNameAndDebug(wallShearStress, 0);
41  addToRunTimeSelectionTable(functionObject, wallShearStress, dictionary);
42 }
43 }
44 
45 
46 // * * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * //
47 
49 {
50  // Add headers to output data
51  writeHeader(file(), "Wall shear stress");
52  writeCommented(file(), "Time");
53  writeTabbed(file(), "patch");
54  writeTabbed(file(), "min");
55  writeTabbed(file(), "max");
56  file() << endl;
57 }
58 
59 
61 (
62  const volSymmTensorField& Reff,
63  volVectorField& shearStress
64 )
65 {
66  shearStress.dimensions().reset(Reff.dimensions());
67 
68  forAllConstIter(labelHashSet, patchSet_, iter)
69  {
70  label patchi = iter.key();
71 
72  vectorField& ssp = shearStress.boundaryFieldRef()[patchi];
73  const vectorField& Sfp = mesh_.Sf().boundaryField()[patchi];
74  const scalarField& magSfp = mesh_.magSf().boundaryField()[patchi];
75  const symmTensorField& Reffp = Reff.boundaryField()[patchi];
76 
77  ssp = (-Sfp/magSfp) & Reffp;
78  }
79 }
80 
81 
82 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
83 
84 Foam::functionObjects::wallShearStress::wallShearStress
85 (
86  const word& name,
87  const Time& runTime,
88  const dictionary& dict
89 )
90 :
91  fvMeshFunctionObject(name, runTime, dict),
92  logFiles(obr_, name),
93  writeLocalObjects(obr_, log),
94  patchSet_()
95 {
96  volVectorField* wallShearStressPtr
97  (
98  new volVectorField
99  (
100  IOobject
101  (
102  type(),
103  mesh_.time().timeName(),
104  mesh_,
107  ),
108  mesh_,
110  (
111  "0",
113  Zero
114  )
115  )
116  );
117 
118  mesh_.objectRegistry::store(wallShearStressPtr);
119 
120  read(dict);
121  resetName(typeName);
122  resetLocalObjectName(typeName);
123 }
124 
125 
126 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
127 
129 {}
130 
131 
132 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
133 
135 {
138 
139  const polyBoundaryMesh& pbm = mesh_.boundaryMesh();
140 
141  patchSet_ =
142  mesh_.boundaryMesh().patchSet
143  (
144  wordReList(dict.lookupOrDefault("patches", wordReList()))
145  );
146 
147  Info<< type() << " " << name() << ":" << nl;
148 
149  if (patchSet_.empty())
150  {
151  forAll(pbm, patchi)
152  {
153  if (isA<wallPolyPatch>(pbm[patchi]))
154  {
155  patchSet_.insert(patchi);
156  }
157  }
158 
159  Info<< " processing all wall patches" << nl << endl;
160  }
161  else
162  {
163  Info<< " processing wall patches: " << nl;
164  labelHashSet filteredPatchSet;
165  forAllConstIter(labelHashSet, patchSet_, iter)
166  {
167  label patchi = iter.key();
168  if (isA<wallPolyPatch>(pbm[patchi]))
169  {
170  filteredPatchSet.insert(patchi);
171  Info<< " " << pbm[patchi].name() << endl;
172  }
173  else
174  {
176  << "Requested wall shear stress on non-wall boundary "
177  << "type patch: " << pbm[patchi].name() << endl;
178  }
179  }
180 
181  Info<< endl;
182 
183  patchSet_ = filteredPatchSet;
184  }
185 
186  return true;
187 }
188 
189 
191 {
192  typedef compressible::turbulenceModel cmpModel;
193  typedef incompressible::turbulenceModel icoModel;
194 
197 
199  if (mesh_.foundObject<cmpModel>(turbulenceModel::propertiesName))
200  {
201  const cmpModel& model =
202  mesh_.lookupObject<cmpModel>(turbulenceModel::propertiesName);
203 
204  Reff = model.devRhoReff();
205  }
206  else if (mesh_.foundObject<icoModel>(turbulenceModel::propertiesName))
207  {
208  const icoModel& model =
209  mesh_.lookupObject<icoModel>(turbulenceModel::propertiesName);
210 
211  Reff = model.devReff();
212  }
213  else
214  {
216  << "Unable to find turbulence model in the "
217  << "database" << exit(FatalError);
218  }
219 
220  calcShearStress(Reff(), wallShearStress);
221 
222  return true;
223 }
224 
225 
227 {
228  Log << type() << " " << name() << " write:" << nl;
229 
231 
232  logFiles::write();
233 
236 
237  const fvPatchList& patches = mesh_.boundary();
238 
239  forAllConstIter(labelHashSet, patchSet_, iter)
240  {
241  label patchi = iter.key();
242  const fvPatch& pp = patches[patchi];
243 
244  const vectorField& ssp = wallShearStress.boundaryField()[patchi];
245 
246  vector minSsp = gMin(ssp);
247  vector maxSsp = gMax(ssp);
248 
249  if (Pstream::master())
250  {
251  file() << mesh_.time().value()
252  << token::TAB << pp.name()
253  << token::TAB << minSsp
254  << token::TAB << maxSsp
255  << endl;
256  }
257 
258  Log << " min/max(" << pp.name() << ") = "
259  << minSsp << ", " << maxSsp << endl;
260  }
261 
262  Log << endl;
263 
264  return true;
265 }
266 
267 
268 // ************************************************************************* //
Foam::surfaceFields.
virtual bool write()
Write function.
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
virtual bool write()
Write function.
Definition: logFiles.C:180
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
const word & name() const
Return name.
Definition: IOobject.H:291
dimensionedScalar log(const dimensionedScalar &ds)
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
virtual bool read(const dictionary &)
Read the wallShearStress data.
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
Type gMin(const FieldField< Field, Type > &f)
const Boundary & boundaryField() const
Return const-reference to the boundary field.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:61
static bool master(const label communicator=0)
Am I the master process.
Definition: UPstream.H:412
patches[0]
Templated wrapper class to provide compressible turbulence models thermal diffusivity based thermal t...
bool insert(const Key &key)
Insert a new entry.
Definition: HashSet.H:116
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.
const ObjectType & lookupObject(const word &fieldName) const
Lookup object from the objectRegistry.
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:68
const word & name() const
Return name.
Definition: fvPatch.H:149
Macros for easy insertion into run-time selection tables.
virtual bool execute()
Calculate the wall shear-stress.
Templated abstract base class for single-phase incompressible turbulence models.
const dimensionSet & dimensions() const
Return dimensions.
virtual bool read(const dictionary &)
Read optional controls.
bool read(const char *, int32_t &)
Definition: int32IO.C:85
virtual void writeFileHeader(const label i)
File header information.
static const word propertiesName
Default name of the turbulence properties dictionary.
void calcShearStress(const volSymmTensorField &Reff, volVectorField &shearStress)
Calculate the shear-stress.
A class for handling words, derived from string.
Definition: word.H:59
Calculates and write the shear-stress at wall patches as the volVectorField field &#39;wallShearStress&#39;...
static const zero Zero
Definition: zero.H:91
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:29
Foam::polyBoundaryMesh.
static const char nl
Definition: Ostream.H:262
Type gMax(const FieldField< Field, Type > &f)
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
fileName::Type type(const fileName &, const bool followLink=true)
Return the file type: DIRECTORY or FILE.
Definition: POSIX.C:485
void reset(const dimensionSet &)
Definition: dimensionSet.C:108
virtual bool write()
Write the wall shear-stress.
label patchi
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
#define WarningInFunction
Report a warning using Foam::Warning.
Boundary & boundaryFieldRef()
Return a reference to the boundary field.
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:50
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: List.H:63
ObjectType & lookupObjectRef(const word &fieldName)
Lookup non-const object reference from the objectRegistry.
FunctionObject base class for managing a list of objects on behalf of the inheriting function object...
virtual bool read(const dictionary &)
Read the list of objects to be written.
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Definition: dimensionSets.H:51
#define Log
Report write to Foam::Info if the local log switch is true.
messageStream Info
Specialization of Foam::functionObject for an Foam::fvMesh, providing a reference to the Foam::fvMesh...
defineTypeNameAndDebug(fvMeshFunctionObject, 0)
List< wordRe > wordReList
A List of wordRe (word or regular expression)
Definition: wordReList.H:50
A class for managing temporary objects.
Definition: PtrList.H:53
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
addToRunTimeSelectionTable(functionObject, add, dictionary)
Namespace for OpenFOAM.
functionObject base class for creating, maintaining and writing log files e.g. integrated of averaged...
Definition: logFiles.H:57