Go to the documentation of this file.
43 namespace functionObjects
74 writeBinHeader(
"Roll moment coefficient bins", CmRollBinFilePtr_());
93 writeHeaderValue(os,
"dragDir", coordSys_.e1());
94 writeHeaderValue(os,
"sideDir", coordSys_.e2());
95 writeHeaderValue(os,
"liftDir", coordSys_.e3());
96 writeHeaderValue(os,
"rollAxis", coordSys_.e1());
97 writeHeaderValue(os,
"pitchAxis", coordSys_.e2());
98 writeHeaderValue(os,
"yawAxis", coordSys_.e3());
99 writeHeaderValue(os,
"magUInf", magUInf_);
100 writeHeaderValue(os,
"lRef", lRef_);
101 writeHeaderValue(os,
"Aref", Aref_);
102 writeHeaderValue(os,
"CofR", coordSys_.origin());
104 writeCommented(os,
"Time");
105 writeTabbed(os,
"Cd");
106 writeTabbed(os,
"Cs");
107 writeTabbed(os,
"Cl");
108 writeTabbed(os,
"CmRoll");
109 writeTabbed(os,
"CmPitch");
110 writeTabbed(os,
"CmYaw");
111 writeTabbed(os,
"Cd(f)");
112 writeTabbed(os,
"Cd(r)");
113 writeTabbed(os,
"Cs(f)");
114 writeTabbed(os,
"Cs(r)");
115 writeTabbed(os,
"Cl(f)");
116 writeTabbed(os,
"Cl(r)");
128 writeHeaderValue(os,
"bins", nBin_);
129 writeHeaderValue(os,
"start", binMin_);
130 writeHeaderValue(os,
"delta", binDx_);
131 writeHeaderValue(os,
"direction", binDir_);
134 writeCommented(os,
"x co-ords :");
137 binPoints[pointi] = (binMin_ + (pointi + 1)*binDx_)*binDir_;
138 os <<
tab << binPoints[pointi].x();
142 writeCommented(os,
"y co-ords :");
145 os <<
tab << binPoints[pointi].y();
149 writeCommented(os,
"z co-ords :");
152 os <<
tab << binPoints[pointi].z();
157 writeCommented(os,
"Time");
159 for (label j = 0; j < nBin_; ++j)
162 writeTabbed(os,
jn +
"total");
163 writeTabbed(os,
jn +
"pressure");
164 writeTabbed(os,
jn +
"viscous");
168 writeTabbed(os,
jn +
"porous");
188 const scalar viscous =
sum(coeff[1]);
189 const scalar porous =
sum(coeff[2]);
190 const scalar total =
pressure + viscous + porous;
195 <<
"viscous: " << viscous;
212 writeCurrentTime(os);
214 for (label bini = 0; bini < nBin_; ++bini)
216 scalar total = coeffs[0][bini] + coeffs[1][bini] + coeffs[2][bini];
218 os <<
tab << total <<
tab << coeffs[0][bini] <<
tab << coeffs[1][bini];
222 os <<
tab << coeffs[2][bini];
232 Foam::functionObjects::forceCoeffs::forceCoeffs
249 CmPitchBinFilePtr_(),
255 setCoordinateSystem(
dict,
"liftDir",
"dragDir");
268 dict.readEntry(
"magUInf", magUInf_);
273 if (rhoName_ !=
"rhoInf")
275 dict.readEntry(
"rhoInf", rhoRef_);
279 dict.readEntry(
"lRef", lRef_);
280 dict.readEntry(
"Aref", Aref_);
290 fieldName(
"forceCoeff"),
291 mesh_.time().timeName(),
301 mesh_.objectRegistry::store(forceCoeffPtr);
309 fieldName(
"momentCoeff"),
310 mesh_.time().timeName(),
320 mesh_.objectRegistry::store(momentCoeffPtr);
346 rollMomentCoeffs[i].
setSize(nBin_);
347 pitchMomentCoeffs[i].
setSize(nBin_);
348 yawMomentCoeffs[i].
setSize(nBin_);
355 scalar CmRollTot = 0;
356 scalar CmPitchTot = 0;
359 const scalar
pDyn = 0.5*rhoRef_*
sqr(magUInf_);
362 const scalar momentScaling = 1.0/(Aref_*
pDyn*lRef_ + SMALL);
363 const scalar forceScaling = 1.0/(Aref_*
pDyn + SMALL);
367 const Field<vector> localForce(coordSys_.localVector(force_[i]));
368 const Field<vector> localMoment(coordSys_.localVector(moment_[i]));
370 dragCoeffs[i] = forceScaling*(localForce.component(0));
371 sideCoeffs[i] = forceScaling*(localForce.component(1));
372 liftCoeffs[i] = forceScaling*(localForce.component(2));
373 rollMomentCoeffs[i] = momentScaling*(localMoment.component(0));
374 pitchMomentCoeffs[i] = momentScaling*(localMoment.component(1));
375 yawMomentCoeffs[i] = momentScaling*(localMoment.component(2));
378 CsTot +=
sum(sideCoeffs[i]);
379 ClTot +=
sum(liftCoeffs[i]);
380 CmRollTot +=
sum(rollMomentCoeffs[i]);
381 CmPitchTot +=
sum(pitchMomentCoeffs[i]);
382 CmYawTot +=
sum(yawMomentCoeffs[i]);
386 const scalar CdfTot = 0.5*CdTot + CmRollTot;
387 const scalar CdrTot = 0.5*CdTot - CmRollTot;
388 const scalar CsfTot = 0.5*CsTot + CmYawTot;
389 const scalar CsrTot = 0.5*CsTot - CmYawTot;
390 const scalar ClfTot = 0.5*ClTot + CmPitchTot;
391 const scalar ClrTot = 0.5*ClTot - CmPitchTot;
394 <<
" Coefficients" <<
nl;
397 writeIntegratedData(
"Cs", sideCoeffs);
398 writeIntegratedData(
"Cl", liftCoeffs);
399 writeIntegratedData(
"CmRoll", rollMomentCoeffs);
400 writeIntegratedData(
"CmPitch", pitchMomentCoeffs);
401 writeIntegratedData(
"CmYaw", yawMomentCoeffs);
403 Log <<
" Cd(f) : " << CdfTot <<
nl
404 <<
" Cd(r) : " << CdrTot <<
nl;
406 Log <<
" Cs(f) : " << CsfTot <<
nl
407 <<
" Cs(r) : " << CsrTot <<
nl;
409 Log <<
" Cl(f) : " << ClfTot <<
nl
410 <<
" Cl(r) : " << ClrTot <<
nl;
414 writeCurrentTime(coeffFilePtr_());
416 <<
tab << CdTot <<
tab << CsTot <<
tab << ClTot
417 <<
tab << CmRollTot <<
tab << CmPitchTot <<
tab << CmYawTot
418 <<
tab << CdfTot <<
tab << CdrTot
419 <<
tab << CsfTot <<
tab << CsrTot
428 for (label bini = 1; bini < nBin_; ++bini)
431 sideCoeffs[i][bini] += sideCoeffs[i][bini-1];
432 liftCoeffs[i][bini] += liftCoeffs[i][bini-1];
433 rollMomentCoeffs[i][bini] +=
434 rollMomentCoeffs[i][bini-1];
435 pitchMomentCoeffs[i][bini] +=
436 pitchMomentCoeffs[i][bini-1];
437 yawMomentCoeffs[i][bini] += yawMomentCoeffs[i][bini-1];
443 writeBinData(sideCoeffs, CsBinFilePtr_());
444 writeBinData(liftCoeffs, ClBinFilePtr_());
445 writeBinData(rollMomentCoeffs, CmRollBinFilePtr_());
446 writeBinData(pitchMomentCoeffs, CmPitchBinFilePtr_());
447 writeBinData(yawMomentCoeffs, CmYawBinFilePtr_());
453 setResult(
"Cd", CdTot);
454 setResult(
"Cs", CsTot);
455 setResult(
"Cl", ClTot);
456 setResult(
"CmRoll", CmRollTot);
457 setResult(
"CmPitch", CmPitchTot);
458 setResult(
"CmYaw", CmYawTot);
459 setResult(
"Cd(f)", CdfTot);
460 setResult(
"Cd(r)", CdrTot);
461 setResult(
"Cs(f)", CsfTot);
462 setResult(
"Cs(r)", CsrTot);
463 setResult(
"Cl(f)", ClfTot);
464 setResult(
"Cl(r)", ClrTot);
470 lookupObject<volVectorField>(fieldName(
"force"));
473 lookupObject<volVectorField>(fieldName(
"moment"));
476 lookupObjectRef<volVectorField>(fieldName(
"forceCoeff"));
479 lookupObjectRef<volVectorField>(fieldName(
"momentCoeff"));
484 forceCoeff == force/f0;
485 momentCoeff == moment/m0;
497 lookupObject<volVectorField>(fieldName(
"forceCoeff"));
500 lookupObject<volVectorField>(fieldName(
"momentCoeff"));
Defines the attributes of an object for which implicit objectRegistry management is supported,...
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Dimensionless.
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
A class for handling words, derived from Foam::string.
volScalarField pDyn(IOobject("pDyn", runTime.timeName(), mesh, IOobject::NO_READ, IOobject::AUTO_WRITE), mesh, dimensionedScalar(dimPressure, Zero))
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
Calculates the forces and moments by integrating the pressure and skin-friction forces over a given l...
static constexpr const zero Zero
Global zero (0)
dimensionedScalar jn(const int n, const dimensionedScalar &ds)
bool read(const char *buf, int32_t &val)
Same as readInt32.
Ostream & endl(Ostream &os)
Add newline and flush stream.
static void writeHeader(Ostream &os, const word &fieldName)
const dimensionSet dimForce
void writeBinData(const List< Field< scalar >> coeffs, Ostream &os) const
Write binned data.
Provides several methods to convert an input pressure field into derived forms, including:
virtual bool read(const dictionary &)
Read the forces data.
void writeBinHeader(const word &header, Ostream &os) const
Write header for binned data.
#define forAll(list, i)
Loop across all elements in list.
virtual void calcForcesMoment()
Calculate the forces and moments.
dimensioned< vector > dimensionedVector
Dimensioned vector obtained from generic dimensioned type.
void writeIntegratedHeader(const word &header, Ostream &os) const
Write header for integrated data.
label nBin_
Number of bins.
messageStream Info
Information stream (uses stdout - output is on the master only)
word name(const complex &c)
Return string representation of complex.
void writeIntegratedData(const word &title, const List< Field< scalar >> &coeff) const
Write integrated data.
autoPtr< multiphaseSystem::dragCoeffFields > dragCoeffs(fluid.dragCoeffs())
Istream and Ostream manipulators taking arguments.
virtual bool read(const dictionary &)
Read the forces data.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Macros for easy insertion into run-time selection tables.
virtual bool write()
Write the forces.
virtual bool execute()
Execute.
void createFiles()
Create the output files.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
virtual autoPtr< OFstream > createFile(const word &name, scalar timeValue) const
Return autoPtr to a new file for a given time.
addToRunTimeSelectionTable(functionObject, ObukhovLength, dictionary)
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
defineTypeNameAndDebug(ObukhovLength, 0)
Reads fields from the time directories and adds them to the mesh database for further post-processing...
Computes the natural logarithm of an input volScalarField.
virtual bool writeToFile() const
Flag to allow writing to file.
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
void setSize(const label newSize)
Alias for resize(const label)