Go to the documentation of this file.
43 namespace functionObjects
74 writeBinHeader(
"Roll moment coefficient bins", CmRollBinFilePtr_());
92 const auto& coordSys = coordSysPtr_();
94 writeHeaderValue(
os,
"dragDir", coordSys.e1());
95 writeHeaderValue(
os,
"sideDir", coordSys.e2());
96 writeHeaderValue(
os,
"liftDir", coordSys.e3());
97 writeHeaderValue(
os,
"rollAxis", coordSys.e1());
98 writeHeaderValue(
os,
"pitchAxis", coordSys.e2());
99 writeHeaderValue(
os,
"yawAxis", coordSys.e3());
100 writeHeaderValue(
os,
"magUInf", magUInf_);
101 writeHeaderValue(
os,
"lRef", lRef_);
102 writeHeaderValue(
os,
"Aref", Aref_);
103 writeHeaderValue(
os,
"CofR", coordSys.origin());
105 writeCommented(
os,
"Time");
106 writeTabbed(
os,
"Cd");
107 writeTabbed(
os,
"Cs");
108 writeTabbed(
os,
"Cl");
109 writeTabbed(
os,
"CmRoll");
110 writeTabbed(
os,
"CmPitch");
111 writeTabbed(
os,
"CmYaw");
112 writeTabbed(
os,
"Cd(f)");
113 writeTabbed(
os,
"Cd(r)");
114 writeTabbed(
os,
"Cs(f)");
115 writeTabbed(
os,
"Cs(r)");
116 writeTabbed(
os,
"Cl(f)");
117 writeTabbed(
os,
"Cl(r)");
129 writeHeaderValue(
os,
"bins", nBin_);
130 writeHeaderValue(
os,
"start", binMin_);
131 writeHeaderValue(
os,
"delta", binDx_);
132 writeHeaderValue(
os,
"direction", binDir_);
135 writeCommented(
os,
"x co-ords :");
138 binPoints[pointi] = (binMin_ + (pointi + 1)*binDx_)*binDir_;
139 os <<
tab << binPoints[pointi].x();
143 writeCommented(
os,
"y co-ords :");
146 os <<
tab << binPoints[pointi].y();
150 writeCommented(
os,
"z co-ords :");
153 os <<
tab << binPoints[pointi].z();
158 writeCommented(
os,
"Time");
160 for (label j = 0; j < nBin_; ++j)
163 writeTabbed(
os,
jn +
"total");
164 writeTabbed(
os,
jn +
"pressure");
165 writeTabbed(
os,
jn +
"viscous");
169 writeTabbed(
os,
jn +
"porous");
188 const scalar pressure =
sum(coeff[0]);
189 const scalar viscous =
sum(coeff[1]);
190 const scalar porous =
sum(coeff[2]);
191 const scalar total = pressure + viscous + porous;
196 <<
"viscous: " << viscous;
213 writeCurrentTime(
os);
215 for (label bini = 0; bini < nBin_; ++bini)
217 scalar total = coeffs[0][bini] + coeffs[1][bini] + coeffs[2][bini];
219 os <<
tab << total <<
tab << coeffs[0][bini] <<
tab << coeffs[1][bini];
223 os <<
tab << coeffs[2][bini];
233 Foam::functionObjects::forceCoeffs::forceCoeffs
250 CmPitchBinFilePtr_(),
269 dict.readEntry(
"magUInf", magUInf_);
274 if (rhoName_ !=
"rhoInf")
276 dict.readEntry(
"rhoInf", rhoRef_);
280 dict.readEntry(
"lRef", lRef_);
281 dict.readEntry(
"Aref", Aref_);
291 scopedName(
"forceCoeff"),
292 mesh_.time().timeName(),
302 mesh_.objectRegistry::store(forceCoeffPtr);
310 scopedName(
"momentCoeff"),
311 mesh_.time().timeName(),
321 mesh_.objectRegistry::store(momentCoeffPtr);
347 rollMomentCoeffs[i].
setSize(nBin_);
348 pitchMomentCoeffs[i].
setSize(nBin_);
349 yawMomentCoeffs[i].
setSize(nBin_);
356 scalar CmRollTot = 0;
357 scalar CmPitchTot = 0;
360 const scalar
pDyn = 0.5*rhoRef_*
sqr(magUInf_);
363 const scalar momentScaling = 1.0/(Aref_*
pDyn*lRef_ + SMALL);
364 const scalar forceScaling = 1.0/(Aref_*
pDyn + SMALL);
366 const auto& coordSys = coordSysPtr_();
370 const Field<vector> localForce(coordSys.localVector(force_[i]));
371 const Field<vector> localMoment(coordSys.localVector(moment_[i]));
373 dragCoeffs[i] = forceScaling*(localForce.component(0));
374 sideCoeffs[i] = forceScaling*(localForce.component(1));
375 liftCoeffs[i] = forceScaling*(localForce.component(2));
376 rollMomentCoeffs[i] = momentScaling*(localMoment.component(0));
377 pitchMomentCoeffs[i] = momentScaling*(localMoment.component(1));
378 yawMomentCoeffs[i] = momentScaling*(localMoment.component(2));
381 CsTot +=
sum(sideCoeffs[i]);
382 ClTot +=
sum(liftCoeffs[i]);
383 CmRollTot +=
sum(rollMomentCoeffs[i]);
384 CmPitchTot +=
sum(pitchMomentCoeffs[i]);
385 CmYawTot +=
sum(yawMomentCoeffs[i]);
389 const scalar CdfTot = 0.5*CdTot + CmRollTot;
390 const scalar CdrTot = 0.5*CdTot - CmRollTot;
391 const scalar CsfTot = 0.5*CsTot + CmYawTot;
392 const scalar CsrTot = 0.5*CsTot - CmYawTot;
393 const scalar ClfTot = 0.5*ClTot + CmPitchTot;
394 const scalar ClrTot = 0.5*ClTot - CmPitchTot;
397 <<
" Coefficients" <<
nl;
400 writeIntegratedData(
"Cs", sideCoeffs);
401 writeIntegratedData(
"Cl", liftCoeffs);
402 writeIntegratedData(
"CmRoll", rollMomentCoeffs);
403 writeIntegratedData(
"CmPitch", pitchMomentCoeffs);
404 writeIntegratedData(
"CmYaw", yawMomentCoeffs);
406 Log <<
" Cd(f) : " << CdfTot <<
nl
407 <<
" Cd(r) : " << CdrTot <<
nl;
409 Log <<
" Cs(f) : " << CsfTot <<
nl
410 <<
" Cs(r) : " << CsrTot <<
nl;
412 Log <<
" Cl(f) : " << ClfTot <<
nl
413 <<
" Cl(r) : " << ClrTot <<
nl;
417 writeCurrentTime(coeffFilePtr_());
419 <<
tab << CdTot <<
tab << CsTot <<
tab << ClTot
420 <<
tab << CmRollTot <<
tab << CmPitchTot <<
tab << CmYawTot
421 <<
tab << CdfTot <<
tab << CdrTot
422 <<
tab << CsfTot <<
tab << CsrTot
431 for (label bini = 1; bini < nBin_; ++bini)
434 sideCoeffs[i][bini] += sideCoeffs[i][bini-1];
435 liftCoeffs[i][bini] += liftCoeffs[i][bini-1];
436 rollMomentCoeffs[i][bini] +=
437 rollMomentCoeffs[i][bini-1];
438 pitchMomentCoeffs[i][bini] +=
439 pitchMomentCoeffs[i][bini-1];
440 yawMomentCoeffs[i][bini] += yawMomentCoeffs[i][bini-1];
446 writeBinData(sideCoeffs, CsBinFilePtr_());
447 writeBinData(liftCoeffs, ClBinFilePtr_());
448 writeBinData(rollMomentCoeffs, CmRollBinFilePtr_());
449 writeBinData(pitchMomentCoeffs, CmPitchBinFilePtr_());
450 writeBinData(yawMomentCoeffs, CmYawBinFilePtr_());
456 setResult(
"Cd", CdTot);
457 setResult(
"Cs", CsTot);
458 setResult(
"Cl", ClTot);
459 setResult(
"CmRoll", CmRollTot);
460 setResult(
"CmPitch", CmPitchTot);
461 setResult(
"CmYaw", CmYawTot);
462 setResult(
"Cd(f)", CdfTot);
463 setResult(
"Cd(r)", CdrTot);
464 setResult(
"Cs(f)", CsfTot);
465 setResult(
"Cs(r)", CsrTot);
466 setResult(
"Cl(f)", ClfTot);
467 setResult(
"Cl(r)", ClrTot);
473 lookupObject<volVectorField>(scopedName(
"force"));
476 lookupObject<volVectorField>(scopedName(
"moment"));
479 lookupObjectRef<volVectorField>(scopedName(
"forceCoeff"));
482 lookupObjectRef<volVectorField>(scopedName(
"momentCoeff"));
487 forceCoeff == force/f0;
488 momentCoeff == moment/m0;
500 lookupObject<volVectorField>(scopedName(
"forceCoeff"));
503 lookupObject<volVectorField>(scopedName(
"momentCoeff"));
Defines the attributes of an object for which implicit objectRegistry management is supported,...
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
dimensionedScalar jn(const int n, const dimensionedScalar &ds)
bool read(const char *buf, int32_t &val)
Ostream & endl(Ostream &os)
static void writeHeader(Ostream &os, const word &fieldName)
const dimensionSet dimForce
void writeBinData(const List< Field< scalar >> coeffs, Ostream &os) const
Provides several methods to convert an input pressure field into derived forms, including:
virtual bool read(const dictionary &)
void setCoordinateSystem(const dictionary &dict, const word &e3Name=word::null, const word &e1Name=word::null)
void writeBinHeader(const word &header, Ostream &os) const
Abstract base-class for Time/database function objects.
virtual void calcForcesMoment()
dimensioned< vector > dimensionedVector
Dimensioned vector obtained from generic dimensioned type.
void writeIntegratedHeader(const word &header, Ostream &os) const
Generic templated field type.
void writeIntegratedData(const word &title, const List< Field< scalar >> &coeff) const
void setSize(const label n)
Extends the forces functionObject by providing coefficients for:
autoPtr< multiphaseSystem::dragCoeffFields > dragCoeffs(fluid.dragCoeffs())
Istream and Ostream manipulators taking arguments.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
virtual bool read(const dictionary &)
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
dimensionedScalar log(const dimensionedScalar &ds)
OBJstream os(runTime.globalPath()/outputName)
Macros for easy insertion into run-time selection tables.
Generic dimensioned Type class.
tmp< Field< cmptType > > component(const direction) const
GeometricField< vector, fvPatchField, volMesh > volVectorField
dimensionedSymmTensor sqr(const dimensionedVector &dv)
virtual autoPtr< OFstream > createFile(const word &name, scalar timeValue) const
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)
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...
word name(const expressions::valueTypeCode typeCode)
Computes the natural logarithm of an input volScalarField.
virtual bool writeToFile() const
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Generic GeometricField class.
const dimensionSet dimless