wallShearStress.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2013-2021 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 
62 (
63  const volSymmTensorField& tau
64 )
65 {
66  tmp<volVectorField> twallShearStress
67  (
69  (
70  type(),
71  mesh_,
73  )
74  );
75 
76  volVectorField::Boundary& wallShearStressBf =
77  twallShearStress.ref().boundaryFieldRef();
78 
79  forAllConstIter(labelHashSet, patchSet_, iter)
80  {
81  label patchi = iter.key();
82 
83  const vectorField& Sfp = mesh_.Sf().boundaryField()[patchi];
84  const scalarField& magSfp = mesh_.magSf().boundaryField()[patchi];
85 
86  wallShearStressBf[patchi] = (-Sfp/magSfp) & tau.boundaryField()[patchi];
87  }
88 
89  return twallShearStress;
90 }
91 
92 
93 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
94 
96 (
97  const word& name,
98  const Time& runTime,
99  const dictionary& dict
100 )
101 :
102  fvMeshFunctionObject(name, runTime, dict),
103  logFiles(obr_, name),
104  writeLocalObjects(obr_, log),
105  patchSet_()
106 {
107  read(dict);
108 }
109 
110 
111 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
112 
114 {}
115 
116 
117 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
118 
120 {
123 
124  const polyBoundaryMesh& pbm = mesh_.boundaryMesh();
125 
126  patchSet_ =
127  mesh_.boundaryMesh().patchSet
128  (
129  wordReList(dict.lookupOrDefault("patches", wordReList()))
130  );
131 
132  Info<< type() << " " << name() << ":" << nl;
133 
134  if (patchSet_.empty())
135  {
136  forAll(pbm, patchi)
137  {
138  if (isA<wallPolyPatch>(pbm[patchi]))
139  {
140  patchSet_.insert(patchi);
141  }
142  }
143 
144  Info<< " processing all wall patches" << nl << endl;
145  }
146  else
147  {
148  Info<< " processing wall patches: " << nl;
149  labelHashSet filteredPatchSet;
150  forAllConstIter(labelHashSet, patchSet_, iter)
151  {
152  label patchi = iter.key();
153  if (isA<wallPolyPatch>(pbm[patchi]))
154  {
155  filteredPatchSet.insert(patchi);
156  Info<< " " << pbm[patchi].name() << endl;
157  }
158  else
159  {
161  << "Requested wall shear stress on non-wall boundary "
162  << "type patch: " << pbm[patchi].name() << endl;
163  }
164  }
165 
166  Info<< endl;
167 
168  patchSet_ = filteredPatchSet;
169  }
170 
171  resetName(typeName);
172  resetLocalObjectName(typeName);
173 
174  return true;
175 }
176 
177 
179 {
180  typedef compressible::momentumTransportModel cmpModel;
182 
184  if (mesh_.foundObject<cmpModel>(momentumTransportModel::typeName))
185  {
186  const cmpModel& model =
187  mesh_.lookupObject<cmpModel>(momentumTransportModel::typeName);
188 
189  tau = model.devTau();
190  }
191  else if (mesh_.foundObject<icoModel>(momentumTransportModel::typeName))
192  {
193  const icoModel& model =
194  mesh_.lookupObject<icoModel>(momentumTransportModel::typeName);
195 
196  tau = model.devSigma();
197  }
198  else
199  {
201  << "Unable to find turbulence model in the "
202  << "database" << exit(FatalError);
203  }
204 
205  word name(type());
206 
207  return store(name, calcShearStress(tau));
208 }
209 
210 
212 {
213  Log << type() << " " << name() << " write:" << nl;
214 
216 
217  logFiles::write();
218 
221 
222  const fvPatchList& patches = mesh_.boundary();
223 
224  forAllConstIter(labelHashSet, patchSet_, iter)
225  {
226  label patchi = iter.key();
227  const fvPatch& pp = patches[patchi];
228 
229  const vectorField& ssp = wallShearStress.boundaryField()[patchi];
230 
231  vector minSsp = gMin(ssp);
232  vector maxSsp = gMax(ssp);
233 
234  if (Pstream::master())
235  {
236  file() << mesh_.time().value()
237  << tab << pp.name()
238  << tab << minSsp
239  << tab << maxSsp
240  << endl;
241  }
242 
243  Log << " min/max(" << pp.name() << ") = "
244  << minSsp << ", " << maxSsp << endl;
245  }
246 
247  Log << endl;
248 
249  return true;
250 }
251 
252 
253 // ************************************************************************* //
Foam::surfaceFields.
virtual bool write()
Write function.
Calculates the natural logarithm of the specified scalar field.
Definition: log.H:103
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
tmp< volVectorField > calcShearStress(const volSymmTensorField &tau)
Calculate the shear-stress.
virtual bool write()
Write function.
Definition: logFiles.C:167
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:303
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
virtual bool read(const dictionary &)
Read the wallShearStress data.
static const char tab
Definition: Ostream.H:259
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:323
Type gMin(const FieldField< Field, Type > &f)
Templated abstract base class for single-phase compressible turbulence models.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:181
addToRunTimeSelectionTable(functionObject, Qdot, dictionary)
static tmp< GeometricField< vector, fvPatchField, volMesh > > New(const word &name, const Internal &, const PtrList< fvPatchField< vector >> &)
Return a temporary field constructed from name,.
dimensioned< vector > dimensionedVector
Dimensioned vector obtained from generic dimensioned type.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:62
static bool master(const label communicator=0)
Am I the master process.
Definition: UPstream.H:423
wallShearStress(const word &name, const Time &runTime, const dictionary &)
Construct from Time and dictionary.
patches[0]
bool insert(const Key &key)
Insert a new entry.
Definition: HashSet.H:111
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
Macros for easy insertion into run-time selection tables.
virtual bool execute()
Calculate the wall shear-stress.
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.
void writeHeader(std::ostream &, const bool isBinary, const std::string &title)
Write header.
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:97
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:29
Foam::polyBoundaryMesh.
static const char nl
Definition: Ostream.H:260
Type gMax(const FieldField< Field, Type > &f)
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
defineTypeNameAndDebug(Qdot, 0)
virtual bool write()
Write the wall shear-stress.
Templated abstract base class for single-phase incompressible turbulence models.
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.
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: List.H:70
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.
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
#define Log
Report write to Foam::Info if the local log switch is true.
messageStream Info
virtual const word & name() const
Return name.
Definition: fvPatch.H:144
Specialisation of Foam::functionObject for an Foam::fvMesh, providing a reference to the Foam::fvMesh...
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
Namespace for OpenFOAM.
functionObject base class for creating, maintaining and writing log files e.g. integrated of averaged...
Definition: logFiles.H:57