forceCoeffs.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-2016 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 "forceCoeffs.H"
28 
29 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
30 
31 namespace Foam
32 {
33 namespace functionObjects
34 {
35  defineTypeNameAndDebug(forceCoeffs, 0);
36  addToRunTimeSelectionTable(functionObject, forceCoeffs, dictionary);
37 }
38 }
39 
40 
41 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
42 
44 {
45  switch (fileID(i))
46  {
47  case MAIN_FILE:
48  {
49  // force coeff data
50 
51  writeHeader(file(i), "Force coefficients");
52  writeHeaderValue(file(i), "liftDir", liftDir_);
53  writeHeaderValue(file(i), "dragDir", dragDir_);
54  writeHeaderValue(file(i), "pitchAxis", pitchAxis_);
55  writeHeaderValue(file(i), "magUInf", magUInf_);
56  writeHeaderValue(file(i), "lRef", lRef_);
57  writeHeaderValue(file(i), "Aref", Aref_);
58  writeHeaderValue(file(i), "CofR", coordSys_.origin());
59  writeCommented(file(i), "Time");
60  writeTabbed(file(i), "Cm");
61  writeTabbed(file(i), "Cd");
62  writeTabbed(file(i), "Cl");
63  writeTabbed(file(i), "Cl(f)");
64  writeTabbed(file(i), "Cl(r)");
65 
66  break;
67  }
68  case BINS_FILE:
69  {
70  // bin coeff data
71 
72  writeHeader(file(i), "Force coefficient bins");
73  writeHeaderValue(file(i), "bins", nBin_);
74  writeHeaderValue(file(i), "start", binMin_);
75  writeHeaderValue(file(i), "delta", binDx_);
76  writeHeaderValue(file(i), "direction", binDir_);
77 
78  vectorField binPoints(nBin_);
79  writeCommented(file(i), "x co-ords :");
80  forAll(binPoints, pointi)
81  {
82  binPoints[pointi] = (binMin_ + (pointi + 1)*binDx_)*binDir_;
83  file(i) << tab << binPoints[pointi].x();
84  }
85  file(i) << nl;
86 
87  writeCommented(file(i), "y co-ords :");
88  forAll(binPoints, pointi)
89  {
90  file(i) << tab << binPoints[pointi].y();
91  }
92  file(i) << nl;
93 
94  writeCommented(file(i), "z co-ords :");
95  forAll(binPoints, pointi)
96  {
97  file(i) << tab << binPoints[pointi].z();
98  }
99  file(i) << nl;
100 
101  writeCommented(file(i), "Time");
102 
103  for (label j = 0; j < nBin_; j++)
104  {
105  const word jn('(' + Foam::name(j) + ')');
106  writeTabbed(file(i), "Cm" + jn);
107  writeTabbed(file(i), "Cd" + jn);
108  writeTabbed(file(i), "Cl" + jn);
109  }
110 
111  break;
112  }
113  default:
114  {
116  << "Unhandled file index: " << i
117  << abort(FatalError);
118  }
119  }
120 
121  file(i)<< endl;
122 }
123 
124 
125 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
126 
127 Foam::functionObjects::forceCoeffs::forceCoeffs
128 (
129  const word& name,
130  const Time& runTime,
131  const dictionary& dict
132 )
133 :
134  forces(name, runTime, dict),
135  liftDir_(Zero),
136  dragDir_(Zero),
137  pitchAxis_(Zero),
138  magUInf_(0.0),
139  lRef_(0.0),
140  Aref_(0.0)
141 {
142  read(dict);
143 }
144 
145 
146 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
147 
149 {}
150 
151 
152 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
153 
155 {
156  forces::read(dict);
157 
158  // Directions for lift and drag forces, and pitch moment
159  dict.lookup("liftDir") >> liftDir_;
160  dict.lookup("dragDir") >> dragDir_;
161  dict.lookup("pitchAxis") >> pitchAxis_;
162 
163  // Free stream velocity magnitude
164  dict.lookup("magUInf") >> magUInf_;
165 
166  // Reference (free stream) density
167  dict.lookup("rhoInf") >> rhoRef_;
168 
169  // Reference length and area scales
170  dict.lookup("lRef") >> lRef_;
171  dict.lookup("Aref") >> Aref_;
172 
173  return true;
174 }
175 
176 
178 {
179  return true;
180 }
181 
182 
184 {
186 
187  if (Pstream::master())
188  {
189  logFiles::write();
190 
191  scalar pDyn = 0.5*rhoRef_*magUInf_*magUInf_;
192 
193  Field<vector> totForce(force_[0] + force_[1] + force_[2]);
194  Field<vector> totMoment(moment_[0] + moment_[1] + moment_[2]);
195 
196  List<Field<scalar>> coeffs(3);
197  coeffs[0].setSize(nBin_);
198  coeffs[1].setSize(nBin_);
199  coeffs[2].setSize(nBin_);
200 
201  // lift, drag and moment
202  coeffs[0] = (totForce & liftDir_)/(Aref_*pDyn);
203  coeffs[1] = (totForce & dragDir_)/(Aref_*pDyn);
204  coeffs[2] = (totMoment & pitchAxis_)/(Aref_*lRef_*pDyn);
205 
206  scalar Cl = sum(coeffs[0]);
207  scalar Cd = sum(coeffs[1]);
208  scalar Cm = sum(coeffs[2]);
209 
210  scalar Clf = Cl/2.0 + Cm;
211  scalar Clr = Cl/2.0 - Cm;
212 
213  writeTime(file(MAIN_FILE));
214  file(MAIN_FILE)
215  << tab << Cm << tab << Cd
216  << tab << Cl << tab << Clf << tab << Clr << endl;
217 
218  Log << type() << " " << name() << " write:" << nl
219  << " Cm = " << Cm << nl
220  << " Cd = " << Cd << nl
221  << " Cl = " << Cl << nl
222  << " Cl(f) = " << Clf << nl
223  << " Cl(r) = " << Clr << endl;
224 
225  if (nBin_ > 1)
226  {
227  if (binCumulative_)
228  {
229  for (label i = 1; i < coeffs[0].size(); i++)
230  {
231  coeffs[0][i] += coeffs[0][i-1];
232  coeffs[1][i] += coeffs[1][i-1];
233  coeffs[2][i] += coeffs[2][i-1];
234  }
235  }
236 
237  writeTime(file(BINS_FILE));
238 
239  forAll(coeffs[0], i)
240  {
241  file(BINS_FILE)
242  << tab << coeffs[2][i]
243  << tab << coeffs[1][i]
244  << tab << coeffs[0][i];
245  }
246 
247  file(BINS_FILE) << endl;
248  }
249 
250  Log << endl;
251  }
252 
253  return true;
254 }
255 
256 
257 // ************************************************************************* //
#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
virtual bool write()
Write the forces.
Definition: forceCoeffs.C:183
static const char tab
Definition: Ostream.H:261
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
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
Definition: HashTable.H:60
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:163
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
static bool master(const label communicator=0)
Am I the master process.
Definition: UPstream.H:412
virtual bool read(const dictionary &)
Read the forces data.
Definition: forces.C:612
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.
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
virtual bool read(const dictionary &)
Read the forces data.
Definition: forceCoeffs.C:154
bool read(const char *, int32_t &)
Definition: int32IO.C:85
A class for handling words, derived from string.
Definition: word.H:59
dimensionedScalar jn(const int n, const dimensionedScalar &ds)
virtual void calcForcesMoment()
Calculate the forces and moments.
Definition: forces.C:746
static const zero Zero
Definition: zero.H:91
errorManip< error > abort(error &err)
Definition: errorManip.H:131
static const char nl
Definition: Ostream.H:262
virtual void writeFileHeader(const label i)
Output file header information.
Definition: forceCoeffs.C:43
fileName::Type type(const fileName &, const bool followLink=true)
Return the file type: DIRECTORY or FILE.
Definition: POSIX.C:485
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
void setSize(const label)
Reset size of List.
Definition: List.C:281
virtual ~forceCoeffs()
Destructor.
Definition: forceCoeffs.C:148
#define Log
Report write to Foam::Info if the local log switch is true.
volScalarField pDyn(IOobject("pDyn", runTime.timeName(), mesh, IOobject::NO_READ, IOobject::AUTO_WRITE), mesh, dimensionedScalar("zero", dimPressure, 0.0))
defineTypeNameAndDebug(fvMeshFunctionObject, 0)
Calculates the forces and moments by integrating the pressure and skin-friction forces over a given l...
Definition: forces.H:198
fileID
Enumeration for ensuring the right file is accessed.
Definition: forces.H:209
virtual bool execute()
Execute, currently does nothing.
Definition: forceCoeffs.C:177
addToRunTimeSelectionTable(functionObject, add, dictionary)
Namespace for OpenFOAM.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:576