diff --git a/src/functionObjects/forces/forceCoeffs/forceCoeffs.C b/src/functionObjects/forces/forceCoeffs/forceCoeffs.C index 376df8038e..d9807ec976 100644 --- a/src/functionObjects/forces/forceCoeffs/forceCoeffs.C +++ b/src/functionObjects/forces/forceCoeffs/forceCoeffs.C @@ -6,7 +6,7 @@ \\/ M anipulation | ------------------------------------------------------------------------------- Copyright (C) 2011-2016 OpenFOAM Foundation - Copyright (C) 2015-2020 OpenCFD Ltd. + Copyright (C) 2015-2022 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -50,47 +50,121 @@ namespace functionObjects // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * // -void Foam::functionObjects::forceCoeffs::createFiles() +void Foam::functionObjects::forceCoeffs::initialise() { - // Note: Only possible to create bin files after bins have been initialised + if (initialised_) + { + return; + } - if (writeToFile() && !coeffFilePtr_) + initialised_ = true; +} + + +void Foam::functionObjects::forceCoeffs::reset() +{ + Cf_.reset(); + Cm_.reset(); + + forceCoeff_ == dimensionedVector(dimless, Zero); + momentCoeff_ == dimensionedVector(dimless, Zero); +} + + +Foam::HashTable +Foam::functionObjects::forceCoeffs::selectCoeffs() const +{ + HashTable coeffs(16); + + coeffs.insert("Cd", coeffDesc("Drag force", "Cd", 0, 0)); + coeffs.insert("Cs", coeffDesc("Side force", "Cs", 1, 2)); + coeffs.insert("Cl", coeffDesc("Lift force", "Cl", 2, 1)); + + // Add front/rear options + const auto frontRearCoeffs(coeffs); + forAllConstIters(frontRearCoeffs, iter) + { + const auto& fr = iter.val(); + coeffs.insert(fr.frontName(), fr.front()); + coeffs.insert(fr.rearName(), fr.rear()); + } + + // Add moments + coeffs.insert("CmRoll", coeffDesc("Roll moment", "CmRoll", 0, -1)); + coeffs.insert("CmPitch", coeffDesc("Pitch moment", "CmPitch", 1, -1)); + coeffs.insert("CmYaw", coeffDesc("Yaw moment", "CmYaw", 2, -1)); + + return coeffs; +} + + +void Foam::functionObjects::forceCoeffs::calcForceCoeffs() +{ + // Calculate scaling factors + const scalar pDyn = 0.5*rhoRef_*sqr(magUInf_); + const dimensionedScalar forceScaling + ( + dimless/dimForce, + scalar(1)/(Aref_*pDyn + SMALL) + ); + + const auto& coordSys = coordSysPtr_(); + + // Calculate force coefficients + forceCoeff_ = forceScaling*force_; + + Cf_.reset + ( + forceScaling.value()*coordSys.localVector(sumPatchForcesV_), + forceScaling.value()*coordSys.localVector(sumPatchForcesP_), + forceScaling.value()*coordSys.localVector(sumInternalForces_) + ); +} + + +void Foam::functionObjects::forceCoeffs::calcMomentCoeffs() +{ + // Calculate scaling factors + const scalar pDyn = 0.5*rhoRef_*sqr(magUInf_); + const dimensionedScalar momentScaling + ( + dimless/(dimForce*dimLength), + scalar(1)/(Aref_*pDyn*lRef_ + SMALL) + ); + + const auto& coordSys = coordSysPtr_(); + + // Calculate moment coefficients + momentCoeff_ = momentScaling*moment_; + + Cm_.reset + ( + momentScaling.value()*coordSys.localVector(sumPatchMomentsP_), + momentScaling.value()*coordSys.localVector(sumPatchMomentsV_), + momentScaling.value()*coordSys.localVector(sumInternalMoments_) + ); +} + + +void Foam::functionObjects::forceCoeffs::createIntegratedDataFile() +{ + if (!coeffFilePtr_.valid()) { coeffFilePtr_ = createFile("coefficient"); - writeIntegratedHeader("Coefficients", coeffFilePtr_()); - - if (nBin_ > 1) - { - CdBinFilePtr_ = createFile("CdBin"); - writeBinHeader("Drag coefficient bins", CdBinFilePtr_()); - - CsBinFilePtr_ = createFile("CsBin"); - writeBinHeader("Side coefficient bins", CsBinFilePtr_()); - - ClBinFilePtr_ = createFile("ClBin"); - writeBinHeader("Lift coefficient bins", ClBinFilePtr_()); - - CmRollBinFilePtr_ = createFile("CmRollBin"); - writeBinHeader("Roll moment coefficient bins", CmRollBinFilePtr_()); - - CmPitchBinFilePtr_ = createFile("CmPitchBin"); - writeBinHeader("Moment coefficient bins", CmPitchBinFilePtr_()); - - CmYawBinFilePtr_ = createFile("CmYawBin"); - writeBinHeader("Yaw moment coefficient bins", CmYawBinFilePtr_()); - } + writeIntegratedDataFileHeader("Coefficients", coeffFilePtr_()); } } -void Foam::functionObjects::forceCoeffs::writeIntegratedHeader +void Foam::functionObjects::forceCoeffs::writeIntegratedDataFileHeader ( const word& header, - Ostream& os + OFstream& os ) const { const auto& coordSys = coordSysPtr_(); - writeHeader(os, "Force coefficients"); + + writeHeader(os, "Force and moment coefficients"); writeHeaderValue(os, "dragDir", coordSys.e1()); writeHeaderValue(os, "sideDir", coordSys.e2()); writeHeaderValue(os, "liftDir", coordSys.e3()); @@ -103,128 +177,38 @@ void Foam::functionObjects::forceCoeffs::writeIntegratedHeader writeHeaderValue(os, "CofR", coordSys.origin()); writeHeader(os, ""); writeCommented(os, "Time"); - writeTabbed(os, "Cd"); - writeTabbed(os, "Cs"); - writeTabbed(os, "Cl"); - writeTabbed(os, "CmRoll"); - writeTabbed(os, "CmPitch"); - writeTabbed(os, "CmYaw"); - writeTabbed(os, "Cd(f)"); - writeTabbed(os, "Cd(r)"); - writeTabbed(os, "Cs(f)"); - writeTabbed(os, "Cs(r)"); - writeTabbed(os, "Cl(f)"); - writeTabbed(os, "Cl(r)"); - os << endl; -} - -void Foam::functionObjects::forceCoeffs::writeBinHeader -( - const word& header, - Ostream& os -) const -{ - writeHeader(os, header); - writeHeaderValue(os, "bins", nBin_); - writeHeaderValue(os, "start", binMin_); - writeHeaderValue(os, "delta", binDx_); - writeHeaderValue(os, "direction", binDir_); - - vectorField binPoints(nBin_); - writeCommented(os, "x co-ords :"); - forAll(binPoints, pointi) + for (const auto& iter : coeffs_.sorted()) { - binPoints[pointi] = (binMin_ + (pointi + 1)*binDx_)*binDir_; - os << tab << binPoints[pointi].x(); - } - os << nl; + const auto& coeff = iter.val(); - writeCommented(os, "y co-ords :"); - forAll(binPoints, pointi) - { - os << tab << binPoints[pointi].y(); - } - os << nl; + if (!coeff.active_) continue; - writeCommented(os, "z co-ords :"); - forAll(binPoints, pointi) - { - os << tab << binPoints[pointi].z(); - } - os << nl; - - writeHeader(os, ""); - writeCommented(os, "Time"); - - for (label j = 0; j < nBin_; ++j) - { - const word jn(Foam::name(j) + ':'); - writeTabbed(os, jn + "total"); - writeTabbed(os, jn + "pressure"); - writeTabbed(os, jn + "viscous"); - - if (porosity_) - { - writeTabbed(os, jn + "porous"); - } + writeTabbed(os, coeff.name_); } os << endl; } -void Foam::functionObjects::forceCoeffs::writeIntegratedData -( - const word& title, - const List>& coeff -) const +void Foam::functionObjects::forceCoeffs::writeIntegratedDataFile() { - if (!log) - { - return; - } + OFstream& os = coeffFilePtr_(); - const scalar pressure = sum(coeff[0]); - const scalar viscous = sum(coeff[1]); - const scalar porous = sum(coeff[2]); - const scalar total = pressure + viscous + porous; - - Info<< " " << title << " : " << total << token::TAB - << '(' - << "pressure: " << pressure << token::TAB - << "viscous: " << viscous; - - if (porosity_) - { - Info<< token::TAB << "porous: " << porous; - } - - Info<< ')' << endl; -} - - -void Foam::functionObjects::forceCoeffs::writeBinData -( - const List> coeffs, - Ostream& os -) const -{ writeCurrentTime(os); - for (label bini = 0; bini < nBin_; ++bini) + for (const auto& iter : coeffs_.sorted()) { - scalar total = coeffs[0][bini] + coeffs[1][bini] + coeffs[2][bini]; + const auto& coeff = iter.val(); - os << tab << total << tab << coeffs[0][bini] << tab << coeffs[1][bini]; + if (!coeff.active_) continue; - if (porosity_) - { - os << tab << coeffs[2][bini]; - } + const vector coeffValue = coeff.value(Cf_, Cm_); + + os << tab << (coeffValue.x() + coeffValue.y() + coeffValue.z()); } - os << endl; + coeffFilePtr_() << endl; } @@ -239,16 +223,39 @@ Foam::functionObjects::forceCoeffs::forceCoeffs ) : forces(name, runTime, dict, false), + Cf_(), + Cm_(), + forceCoeff_ + ( + IOobject + ( + "forceCoeff", // scopedName() is not available at ctor level + mesh_.time().timeName(), + mesh_, + IOobject::NO_READ, + IOobject::NO_WRITE + ), + mesh_, + dimensionedVector(dimless, Zero) + ), + momentCoeff_ + ( + IOobject + ( + "momentCoeff", + mesh_.time().timeName(), + mesh_, + IOobject::NO_READ, + IOobject::NO_WRITE + ), + mesh_, + dimensionedVector(dimless, Zero) + ), + coeffFilePtr_(), magUInf_(Zero), lRef_(Zero), Aref_(Zero), - coeffFilePtr_(), - CdBinFilePtr_(), - CsBinFilePtr_(), - ClBinFilePtr_(), - CmRollBinFilePtr_(), - CmPitchBinFilePtr_(), - CmYawBinFilePtr_() + initialised_(false) { if (readFields) { @@ -263,230 +270,139 @@ Foam::functionObjects::forceCoeffs::forceCoeffs bool Foam::functionObjects::forceCoeffs::read(const dictionary& dict) { - forces::read(dict); + if (!forces::read(dict)) + { + return false; + } - // Free stream velocity magnitude + initialised_ = false; + + // Reference scales dict.readEntry("magUInf", magUInf_); + dict.readEntry("lRef", lRef_); + dict.readEntry("Aref", Aref_); - // If case is compressible we must read rhoInf (store in rhoRef_) to - // calculate the reference dynamic pressure - // Note: for incompressible, rhoRef_ is already initialised + // If case is compressible we must read rhoInf + // (stored in rhoRef_) to calculate the reference dynamic pressure + // Note: for incompressible, rhoRef_ is already initialised to 1 if (rhoName_ != "rhoInf") { dict.readEntry("rhoInf", rhoRef_); } - // Reference length and area scales - dict.readEntry("lRef", lRef_); - dict.readEntry("Aref", Aref_); + Info<< " magUInf: " << magUInf_ << nl + << " lRef : " << lRef_ << nl + << " Aref : " << Aref_ << nl + << " rhoInf : " << rhoRef_ << endl; - if (writeFields_) + if (min(magUInf_, rhoRef_) < SMALL || min(lRef_, Aref_) < SMALL) { - volVectorField* forceCoeffPtr - ( - new volVectorField - ( - IOobject - ( - scopedName("forceCoeff"), - mesh_.time().timeName(), - mesh_, - IOobject::NO_READ, - IOobject::NO_WRITE - ), - mesh_, - dimensionedVector(dimless, Zero) - ) - ); + FatalIOErrorInFunction(dict) + << "Non-zero values are required for reference scales" + << exit(FatalIOError); - mesh_.objectRegistry::store(forceCoeffPtr); - - volVectorField* momentCoeffPtr - ( - new volVectorField - ( - IOobject - ( - scopedName("momentCoeff"), - mesh_.time().timeName(), - mesh_, - IOobject::NO_READ, - IOobject::NO_WRITE - ), - mesh_, - dimensionedVector(dimless, Zero) - ) - ); - - mesh_.objectRegistry::store(momentCoeffPtr); + return false; } + if (!dict.found("coefficients")) + { + Info<< " Selecting all coefficients" << nl; + + coeffs_ = selectCoeffs(); + + for (auto& iter : coeffs_.sorted()) + { + auto& coeff = iter.val(); + coeff.active_ = true; + Info<< " - " << coeff << nl; + } + } + else + { + const wordHashSet coeffs(dict.get("coefficients")); + + coeffs_ = selectCoeffs(); + + Info<< " Selecting coefficients:" << nl; + + for (const word& key : coeffs) + { + auto coeffIter = coeffs_.find(key); + + if (!coeffIter.good()) + { + FatalIOErrorInFunction(dict) + << "Unknown coefficient type " << key + << " . Valid entries are : " << coeffs_.sortedToc() + << exit(FatalIOError); + } + auto& coeff = coeffIter.val(); + coeff.active_ = true; + Info<< " - " << coeff << nl; + } + } + + forceCoeff_.rename(scopedName("forceCoeff")); + momentCoeff_.rename(scopedName("momentCoeff")); + + Info<< endl; + return true; } + bool Foam::functionObjects::forceCoeffs::execute() { - forces::calcForcesMoment(); + forces::calcForcesMoments(); - createFiles(); + initialise(); - // Storage for pressure, viscous and porous contributions to coeffs - List> dragCoeffs(3); - List> sideCoeffs(3); - List> liftCoeffs(3); - List> rollMomentCoeffs(3); - List> pitchMomentCoeffs(3); - List> yawMomentCoeffs(3); + reset(); - forAll(liftCoeffs, i) + Log << type() << " " << name() << " write:" << nl + << " " << "Coefficient" + << tab << "Total" + << tab << "Pressure" + << tab << "Viscous" + << tab << "Internal" + << nl; + + calcForceCoeffs(); + + calcMomentCoeffs(); + + auto logValues = [](const word& name, const vector& coeff, Ostream& os) { - dragCoeffs[i].setSize(nBin_); - sideCoeffs[i].setSize(nBin_); - liftCoeffs[i].setSize(nBin_); - rollMomentCoeffs[i].setSize(nBin_); - pitchMomentCoeffs[i].setSize(nBin_); - yawMomentCoeffs[i].setSize(nBin_); - } + os << " " << name << ":" + << tab << coeff.x() + coeff.y() + coeff.z() + << tab << coeff.x() + << tab << coeff.y() + << tab << coeff.z() + << nl; + }; - // Calculate coefficients - scalar CdTot = 0; - scalar CsTot = 0; - scalar ClTot = 0; - scalar CmRollTot = 0; - scalar CmPitchTot = 0; - scalar CmYawTot = 0; - - const scalar pDyn = 0.5*rhoRef_*sqr(magUInf_); - - // Avoid divide by zero in 2D cases - const scalar momentScaling = 1.0/(Aref_*pDyn*lRef_ + SMALL); - const scalar forceScaling = 1.0/(Aref_*pDyn + SMALL); - - const auto& coordSys = coordSysPtr_(); - - forAll(liftCoeffs, i) + // Always setting all results + for (const auto& iter : coeffs_.sorted()) { - const Field localForce(coordSys.localVector(force_[i])); - const Field localMoment(coordSys.localVector(moment_[i])); + const word& coeffName = iter.key(); + const auto& coeff = iter.val(); - dragCoeffs[i] = forceScaling*(localForce.component(0)); - sideCoeffs[i] = forceScaling*(localForce.component(1)); - liftCoeffs[i] = forceScaling*(localForce.component(2)); - rollMomentCoeffs[i] = momentScaling*(localMoment.component(0)); - pitchMomentCoeffs[i] = momentScaling*(localMoment.component(1)); - yawMomentCoeffs[i] = momentScaling*(localMoment.component(2)); + // Vectors for x:pressure, y:viscous, z:internal + const vector coeffValue = coeff.value(Cf_, Cm_); - CdTot += sum(dragCoeffs[i]); - CsTot += sum(sideCoeffs[i]); - ClTot += sum(liftCoeffs[i]); - CmRollTot += sum(rollMomentCoeffs[i]); - CmPitchTot += sum(pitchMomentCoeffs[i]); - CmYawTot += sum(yawMomentCoeffs[i]); - } - - // Single contributions to the front and rear - const scalar CdfTot = 0.5*CdTot + CmRollTot; - const scalar CdrTot = 0.5*CdTot - CmRollTot; - const scalar CsfTot = 0.5*CsTot + CmYawTot; - const scalar CsrTot = 0.5*CsTot - CmYawTot; - const scalar ClfTot = 0.5*ClTot + CmPitchTot; - const scalar ClrTot = 0.5*ClTot - CmPitchTot; - - Log << type() << " " << name() << " execute:" << nl - << " Coefficients" << nl; - - writeIntegratedData("Cd", dragCoeffs); - writeIntegratedData("Cs", sideCoeffs); - writeIntegratedData("Cl", liftCoeffs); - writeIntegratedData("CmRoll", rollMomentCoeffs); - writeIntegratedData("CmPitch", pitchMomentCoeffs); - writeIntegratedData("CmYaw", yawMomentCoeffs); - - Log << " Cd(f) : " << CdfTot << nl - << " Cd(r) : " << CdrTot << nl; - - Log << " Cs(f) : " << CsfTot << nl - << " Cs(r) : " << CsrTot << nl; - - Log << " Cl(f) : " << ClfTot << nl - << " Cl(r) : " << ClrTot << nl; - - if (writeToFile()) - { - writeCurrentTime(coeffFilePtr_()); - coeffFilePtr_() - << tab << CdTot << tab << CsTot << tab << ClTot - << tab << CmRollTot << tab << CmPitchTot << tab << CmYawTot - << tab << CdfTot << tab << CdrTot - << tab << CsfTot << tab << CsrTot - << tab << ClfTot << tab << ClrTot << endl; - - if (nBin_ > 1) + if (log && coeff.active_) { - if (binCumulative_) - { - forAll(liftCoeffs, i) - { - for (label bini = 1; bini < nBin_; ++bini) - { - dragCoeffs[i][bini] += dragCoeffs[i][bini-1]; - sideCoeffs[i][bini] += sideCoeffs[i][bini-1]; - liftCoeffs[i][bini] += liftCoeffs[i][bini-1]; - rollMomentCoeffs[i][bini] += - rollMomentCoeffs[i][bini-1]; - pitchMomentCoeffs[i][bini] += - pitchMomentCoeffs[i][bini-1]; - yawMomentCoeffs[i][bini] += yawMomentCoeffs[i][bini-1]; - } - } - } - - writeBinData(dragCoeffs, CdBinFilePtr_()); - writeBinData(sideCoeffs, CsBinFilePtr_()); - writeBinData(liftCoeffs, ClBinFilePtr_()); - writeBinData(rollMomentCoeffs, CmRollBinFilePtr_()); - writeBinData(pitchMomentCoeffs, CmPitchBinFilePtr_()); - writeBinData(yawMomentCoeffs, CmYawBinFilePtr_()); + logValues(coeffName, coeffValue, Info); } + + setResult(coeffName, coeffValue.x() + coeffValue.y() + coeffValue.z()); + setResult(coeffName & "pressure", coeffValue.x()); + setResult(coeffName & "viscous", coeffValue.y()); + setResult(coeffName & "internal", coeffValue.z()); } - // Write state/results information - { - setResult("Cd", CdTot); - setResult("Cs", CsTot); - setResult("Cl", ClTot); - setResult("CmRoll", CmRollTot); - setResult("CmPitch", CmPitchTot); - setResult("CmYaw", CmYawTot); - setResult("Cd(f)", CdfTot); - setResult("Cd(r)", CdrTot); - setResult("Cs(f)", CsfTot); - setResult("Cs(r)", CsrTot); - setResult("Cl(f)", ClfTot); - setResult("Cl(r)", ClrTot); - } - - if (writeFields_) - { - const volVectorField& force = - lookupObject(scopedName("force")); - - const volVectorField& moment = - lookupObject(scopedName("moment")); - - volVectorField& forceCoeff = - lookupObjectRef(scopedName("forceCoeff")); - - volVectorField& momentCoeff = - lookupObjectRef(scopedName("momentCoeff")); - - dimensionedScalar f0("f0", dimForce, Aref_*pDyn); - dimensionedScalar m0("m0", dimForce*dimLength, Aref_*lRef_*pDyn); - - forceCoeff == force/f0; - momentCoeff == moment/m0; - } + Log << endl; return true; } @@ -494,18 +410,23 @@ bool Foam::functionObjects::forceCoeffs::execute() bool Foam::functionObjects::forceCoeffs::write() { + if (writeToFile()) + { + Log << " writing force and moment coefficient files." << endl; + + createIntegratedDataFile(); + + writeIntegratedDataFile(); + } + if (writeFields_) { - const volVectorField& forceCoeff = - lookupObject(scopedName("forceCoeff")); - - const volVectorField& momentCoeff = - lookupObject(scopedName("momentCoeff")); - - forceCoeff.write(); - momentCoeff.write(); + forceCoeff_.write(); + momentCoeff_.write(); } + Log << endl; + return true; } diff --git a/src/functionObjects/forces/forceCoeffs/forceCoeffs.H b/src/functionObjects/forces/forceCoeffs/forceCoeffs.H index 0b218f51ae..85fbe23d2e 100644 --- a/src/functionObjects/forces/forceCoeffs/forceCoeffs.H +++ b/src/functionObjects/forces/forceCoeffs/forceCoeffs.H @@ -6,7 +6,7 @@ \\/ M anipulation | ------------------------------------------------------------------------------- Copyright (C) 2011-2016 OpenFOAM Foundation - Copyright (C) 2019-2020 OpenCFD Ltd. + Copyright (C) 2019-2022 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -31,142 +31,207 @@ Group grpForcesFunctionObjects Description - Extends the \c forces functionObject by providing coefficients for: - - drag, side and lift forces (Cd, Cs, and Cl) - - roll, pitch and yaw moments (CmRoll, CmPitch, and CmYaw) - - front and rear axle force contributions (C(f) and C(r)) wherein - + Computes force and moment coefficients over a given + list of patches, and optionally over given porous zones. + The following coefficients can be selected and output: \verbatim - Cd(f/r) = 0.5*Cd \pm CmRoll - Cs(f/r) = 0.5*Cs \pm CmYaw - Cl(f/r) = 0.5*Cl \pm CmPitch + Cd | Drag coefficient + Cs | Side-force coefficient + Cl | Lift coefficient + CmRoll | Roll-moment coefficient + CmPitch | Pitch-moment coefficient + CmYaw | Yaw-moment coefficient \endverbatim - The data can optionally be output into bins, defined in a given direction. + The force coefficients can also be optionally output + in terms of their front and rear axle constituents: + \verbatim + Cd{f,r} = 0.5*Cd {+,-} CmRoll + Cs{f,r} = 0.5*Cs {+,-} CmYaw + Cl{f,r} = 0.5*Cl {+,-} CmPitch + \endverbatim + where \c f and \c r represent front and rear axles, respectively. - The binned data provides the total and consitituent components per bin: - - total coefficient - - pressure coefficient contribution - - viscous coefficient contribution - - porous coefficient contribution + Force and moment coefficients are output + in their total and constituent components: + - total force and moment coefficients + - pressure contributions + - viscous contributions + - porous resistance contributions (optional) - Data is written into multiple files in the - postProcessing/\ directory: - - coefficient.dat : integrated coefficients over all geometries - - CdBin.dat : drag coefficient bins - - CsBin.dat : side coefficient bins - - ClBin.dat : lift coefficient bins - - CmRollBin.dat : roll moment coefficient bins - - CmPitchBin.dat : pitch moment coefficient bins - - CmYawBin.dat : yaw moment coefficient bins + Force and moment coefficients can be computed and output in: + - the global Cartesian coordinate system (default) + - a user-defined Cartesian coordinate system + + Operands: + \table + Operand | Type | Location + input | - | - + output file | dat | postProcessing/\/\/\s + output field | volVectorField | \/\s + \endtable + + where \c \s: + \verbatim + coefficient.dat | Integrated coefficients over all patches + \endverbatim + + where \c \s: + \verbatim + :forceCoeff | Force coefficient field + :momentCoeff | Moment coefficient field + \endverbatim Usage - Example of function object specification: + Minimal example by using \c system/controlDict.functions: \verbatim forceCoeffs1 { - type forceCoeffs; - libs (forces); + // Mandatory entries + type forceCoeffs; + libs (forces); + patches (); // (wall1 "(wall2|wall3)"); + magUInf ; + lRef ; + Aref ; + + // Optional entries + coefficients (); + directForceDensity ; + porosity ; + writeFields ; + useNamePrefix ; + + // Conditional mandatory entries + + // Cartesian coordinate system specification when evaluating + // force and moment coefficients, either of the below + + // Option-1, i.e. the centre of rotation + // by inherently using e3=(0 0 1) and e1=(1 0 0) + CofR (0 0 0); // Centre of rotation + dragDir (1 0 0); + liftDir (0 0 1); + + // Option-2, i.e. local coordinate system specification + origin (0 0 0); + e1 (1 0 0); + e3 (0 0 1); // (e1, e2) or (e2, e3) or (e3, e1) + + // Option-3, i.e. general coordinate system specification + coordinateSystem + { + type cartesian; + origin (0 0 0); + rotation + { + type axes; + e3 (0 0 1); + e1 (1 0 0); // (e1, e2) or (e2, e3) or (e3, e1) + } + } + + // Conditional optional entries + + // if directForceDensity == true + fD ; + + // if directForceDensity == false + p ; + U ; + rho ; + rhoInf ; // redundant for incompressible-flow cases + pRef ; + + // Inherited entries ... - log yes; - writeFields yes; - patches (walls); - - // input keywords for directions of force/moment coefficients - // refer below for options, and relations - - magUInf 100; - lRef 3.5; - Aref 2.2; - porosity no; - - binData - { - nBin 20; - direction (1 0 0); - cumulative yes; - } } \endverbatim - Where the entries comprise: + where the entries mean: \table - Property | Description | Required | Default - type | Type name: forceCoeffs | yes | - log | Write force data to standard output | no | no - writeFields | Write force,moment coefficient fields | no | no - patches | Patches included in the forces calculation | yes | - magUInf | Free stream velocity magnitude | yes | - rhoInf | Free stream density | for compressible cases | - lRef | Reference length scale for moment calculations | yes | - Aref | Reference area | yes | - porosity | Include porosity contributions | no | false + Property | Description | Type | Reqd | Deflt + type | Type name: forceCoeffs | word | yes | - + libs | Library name: forces | word | yes | - + patches | Names of operand patches | wordRes | yes | - + coefficients | Names of operand coefficients | wordHashSet | no | - + magUInf | Reference velocity magnitude | scalar | yes | - + lRef | Reference length scale for moment calculations | scalar | yes | - + Aref | Reference area | scalar | yes | - + directForceDensity | Flag to directly supply force density | bool | no | false + porosity | Flag to include porosity contributions | bool | no | false + writeFields | Flag to write force and moment fields | bool | no | false + useNamePrefix | Flag to include prefix for field names | bool | no | false + CofR | Centre of rotation | vector | cndtnl | - + origin | Origin of coordinate system | vector | cndtnl | - + e3 | e3 coordinate axis | vector | cndtnl | - + e1 | e1 coordinate axis | vector | cndtnl | - + coordinateSystem | Coordinate system specifier | dictionary | cndtnl | - + fD | Name of force density field | word | cndtnl-no | fD + p | Name of pressure field | word | cndtnl-no | p + U | Name of velocity field | word | cndtnl-no | U + rho | Name of density field | word | cndtnl-no | rho + rhoInf | Value of reference density | scalar | cndtnl-yes | - + pRef | Value of reference pressure | scalar | cndtnl-no | 0 \endtable - Bin data is optional, but if the dictionary is present, the entries must - be defined according to following: - \table - nBin | Number of data bins | yes | - direction | Direction along which bins are defined | yes | - cumulative | Bin data accumulated with incresing distance | yes | - \endtable - - Input of force/moment coefficient directions: - - require an origin, and two orthogonal directions; the remaining orthogonal - direction is determined accordingly. - - can be added by the three options below. - + Options for the \c coefficients entry: \verbatim - CofR (0 0 0); // Centre of rotation - dragDir (1 0 0); - liftDir (0 0 1); + Cd | Drag coefficient + Cs | Side-force coefficient + Cl | Lift coefficient + CmRoll | Roll-moment coefficient + CmPitch | Pitch-moment coefficient + CmYaw | Yaw-moment coefficient + Cdf | Front-axle drag coefficient + Csf | Front-axle side-force coefficient + Clf | Front-axle lift coefficient + Cdr | Rear-axle drag coefficient + Csr | Rear-axle side-force coefficient + Clr | Rear-axle lift coefficient \endverbatim - \verbatim - origin (0 0 0); - e1 (1 0 0); - e3 (0 0 1); // combinations: (e1, e2) or (e2, e3) or (e3, e1) - \endverbatim + The inherited entries are elaborated in: + - \link functionObject.H \endlink + - \link writeFile.H \endlink + - \link coordinateSystem.H \endlink + - \link forces.H \endlink - \verbatim - coordinateSystem - { - origin (0 0 0); - rotation - { - type axes; - e1 (1 0 0); - e3 (0 0 1); // combinations: (e1, e2) or (e2, e3) or (e3, e1) - } - } - \endverbatim - - The default direction relations are shown below: +Note + - \c rhoInf is always redundant for incompressible computations. + That is, \c rhoInf is always equal to 1 in incompressible + computations no matter which input value is assigned to \c rhoInf. + The value of \c rhoInf is only used for compressible computations. + - \c writeControl and \c writeInterval entries of function + object do control when to output force and moment files and fields. + - Input for force/moment coefficient directions + require an origin and two orthogonal directions where + the remaining orthogonal direction is automatically computed. + - The default direction relations are shown below: \table - Property | Description | Alias | Direction - dragDir | Drag direction | e1 | (1 0 0) - sideDir | Side force direction | e2 | (0 1 0) - liftDir | Lift direction | e3 | (0 0 1) - rollAxis | Roll axis | e1 | (1 0 0) - pitchAxis | Pitch axis | e2 | (0 1 0) - yawAxis | Yaw axis | e3 | (0 0 1) + Property | Description | Alias | Direction + dragDir | Drag direction | e1 | (1 0 0) + sideDir | Side force direction | e2 | (0 1 0) + liftDir | Lift direction | e3 | (0 0 1) + rollAxis | Roll axis | e1 | (1 0 0) + pitchAxis | Pitch axis | e2 | (0 1 0) + yawAxis | Yaw axis | e3 | (0 0 1) \endtable -See also - Foam::functionObject - Foam::functionObjects::timeControl - Foam::functionObjects::forces - SourceFiles forceCoeffs.C \*---------------------------------------------------------------------------*/ -#ifndef functionObjects_forceCoeffs_H -#define functionObjects_forceCoeffs_H +#ifndef Foam_functionObjects_forceCoeffs_H +#define Foam_functionObjects_forceCoeffs_H #include "forces.H" +#include "HashSet.H" +#include "Enum.H" // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // @@ -183,78 +248,261 @@ class forceCoeffs : public forces { - // Private data +public: - // Free-stream conditions + // Container for force and moment components + class forceComponents + { + // Private Data - //- Free-stream velocity magnitude - scalar magUInf_; + //- Pressure force/moment component + vector pressure_; + + //- Viscous force/moment component + vector viscous_; + + //- Internal force/moment component + vector internal_; + + + public: + + // Constructors + + //- Default construct + forceComponents() + : + pressure_(Zero), + viscous_(Zero), + internal_(Zero) + {} + + + // Member Functions + + //- Return the sum of the components + const vector total() const noexcept + { + return pressure_ + viscous_ + internal_; + } + + //- Reset the components to zero + void reset() noexcept + { + pressure_ = Zero; + viscous_ = Zero; + internal_ = Zero; + } + + //- Reset the components to given values + void reset + ( + const vector& pressure, + const vector& viscous, + const vector& internal + ) noexcept + { + pressure_ = pressure; + viscous_ = viscous; + internal_ = internal; + } + + + // Member Operators + + //- Return components in a given direction + vector operator[](const label cmpt) const + { + return vector(pressure_[cmpt], viscous_[cmpt], internal_[cmpt]); + } + }; + + + // Coefficients description + struct coeffDesc + { + enum splitType + { + stFront, + stRear, + stNone + }; + + string desc_; + word name_; + label c0_; + label c1_; + bool active_; + splitType splitType_; + + // Constructors + + //- Construct from components + coeffDesc + ( + const string& description, + const word& name, + const label c0, + const label c1 = -1 + ) + : + desc_(description), + name_(name), + c0_(c0), + c1_(c1), + active_(false), + splitType_(stNone) + {} + + //- Return name with front-name appendix + const word frontName() const noexcept + { + return name_ + "(f)"; + } + + //- Return name with rear-name appendix + const word rearName() const noexcept + { + return name_ + "(r)"; + } + + //- Return force/moment components based on the specified split type + vector value(const forceComponents& f, const forceComponents& m) const + { + if (c1_ == -1) + { + return m[c0_]; + } + else + { + switch (splitType_) + { + case stFront: + { + return 0.5*f[c0_] + m[c1_]; + } + case stRear: + { + return 0.5*f[c0_] - m[c1_]; + } + case stNone: + { + return f[c0_]; + } + } + } + + FatalErrorInFunction + << "Cannot determine value" + << abort(FatalError); + + return vector::zero; + } + + //- Return front-axle coefficient description + coeffDesc front() const + { + coeffDesc coeff(*this); + coeff.name_ = coeff.frontName(); + coeff.desc_ = coeff.desc_ + " front"; + coeff.splitType_ = stFront; + return coeff; + } + + //- Return rear-axle coefficient description + coeffDesc rear() const + { + coeffDesc coeff(*this); + coeff.name_ = coeff.rearName(); + coeff.desc_ = coeff.desc_ + " rear"; + coeff.splitType_ = stRear; + return coeff; + } + }; + + +private: + + // Private Data + + //- Force components + forceComponents Cf_; + + //- Moment components + forceComponents Cm_; + + //- Force coefficient field + volVectorField forceCoeff_; + + //- Moment coefficient field + volVectorField momentCoeff_; + + //- Operand force and moment coefficients + HashTable coeffs_; + + + // File streams + + //- File stream for integrated operand coefficients + autoPtr coeffFilePtr_; // Reference scales - //- Reference length [m] + //- Reference velocity magnitude [m/s] + scalar magUInf_; + + //- Reference length scale [m] scalar lRef_; //- Reference area [m^2] scalar Aref_; - // File streams - - //- Integrated coefficients - autoPtr coeffFilePtr_; - - //- Drag coefficient - autoPtr CdBinFilePtr_; - - //- Side coefficient - autoPtr CsBinFilePtr_; - - //- Lift coefficient - autoPtr ClBinFilePtr_; - - //- Roll moment coefficient - autoPtr CmRollBinFilePtr_; - - //- Pitch moment coefficient - autoPtr CmPitchBinFilePtr_; - - //- Yaw moment coefficient - autoPtr CmYawBinFilePtr_; - - - // Private Member Functions - - //- No copy construct - forceCoeffs(const forceCoeffs&) = delete; - - //- No copy assignment - void operator=(const forceCoeffs&) = delete; + //- Flag of initialisation (internal) + bool initialised_; protected: // Protected Member Functions - //- Create the output files - void createFiles(); + //- Initialise containers and fields + void initialise(); - //- Write header for integrated data - void writeIntegratedHeader(const word& header, Ostream& os) const; + //- Reset containers and fields + void reset(); - //- Write header for binned data - void writeBinHeader(const word& header, Ostream& os) const; - //- Write integrated data - void writeIntegratedData + // Evaluation + + //- Return the operand coefficients + HashTable selectCoeffs() const; + + //- Calculate force coefficients + void calcForceCoeffs(); + + //- Calculate moment coefficients + void calcMomentCoeffs(); + + //- Return integrated {total, pressure, viscous, porous} components + List integrateData(const List>& coeff) const; + + + // I-O + + //- Create the integrated-coefficient file + void createIntegratedDataFile(); + + //- Write header to the integrated-coefficient file + void writeIntegratedDataFileHeader ( - const word& title, - const List>& coeff + const word& header, + OFstream& os ) const; - //- Write binned data - void writeBinData(const List> coeffs, Ostream& os) const; + //- Write integrated coefficients to the integrated-coefficient file + void writeIntegratedDataFile(); public: @@ -270,10 +518,16 @@ public: ( const word& name, const Time& runTime, - const dictionary&, + const dictionary& dict, const bool readFields = true ); + //- No copy construct + forceCoeffs(const forceCoeffs&) = delete; + + //- No copy assignment + void operator=(const forceCoeffs&) = delete; + //- Destructor virtual ~forceCoeffs() = default; @@ -281,19 +535,27 @@ public: // Member Functions - //- Read the forces data - virtual bool read(const dictionary&); + //- Read the dictionary + virtual bool read(const dictionary& dict); - //- Execute + //- Execute the function object virtual bool execute(); - //- Write the forces + //- Write to data files/fields and to streams virtual bool write(); }; // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // +Ostream& operator<<(Ostream& os, const forceCoeffs::coeffDesc& coeff) +{ + os << coeff.desc_.c_str() << ": " << coeff.name_; + + return os; +} + + } // End namespace functionObjects } // End namespace Foam diff --git a/src/functionObjects/forces/forces/forces.C b/src/functionObjects/forces/forces/forces.C index 2f3648012f..5916aa2a61 100644 --- a/src/functionObjects/forces/forces/forces.C +++ b/src/functionObjects/forces/forces/forces.C @@ -6,7 +6,7 @@ \\/ M anipulation | ------------------------------------------------------------------------------- Copyright (C) 2011-2016 OpenFOAM Foundation - Copyright (C) 2015-2020 OpenCFD Ltd. + Copyright (C) 2015-2022 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -31,8 +31,8 @@ License #include "porosityModel.H" #include "turbulentTransportModel.H" #include "turbulentFluidThermoModel.H" -#include "addToRunTimeSelectionTable.H" #include "cartesianCS.H" +#include "addToRunTimeSelectionTable.H" // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // @@ -48,107 +48,6 @@ namespace functionObjects // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * // -void Foam::functionObjects::forces::createFiles() -{ - // Note: Only possible to create bin files after bins have been initialised - - if (writeToFile() && !forceFilePtr_) - { - forceFilePtr_ = createFile("force"); - writeIntegratedHeader("Force", forceFilePtr_()); - momentFilePtr_ = createFile("moment"); - writeIntegratedHeader("Moment", momentFilePtr_()); - - if (nBin_ > 1) - { - forceBinFilePtr_ = createFile("forceBin"); - writeBinHeader("Force", forceBinFilePtr_()); - momentBinFilePtr_ = createFile("momentBin"); - writeBinHeader("Moment", momentBinFilePtr_()); - } - } -} - - -void Foam::functionObjects::forces::writeIntegratedHeader -( - const word& header, - Ostream& os -) const -{ - writeHeader(os, header); - writeHeaderValue(os, "CofR", coordSysPtr_->origin()); - writeHeader(os, ""); - writeCommented(os, "Time"); - writeTabbed(os, "(total_x total_y total_z)"); - writeTabbed(os, "(pressure_x pressure_y pressure_z)"); - writeTabbed(os, "(viscous_x viscous_y viscous_z)"); - - if (porosity_) - { - writeTabbed(os, "(porous_x porous_y porous_z)"); - } - - os << endl; -} - - -void Foam::functionObjects::forces::writeBinHeader -( - const word& header, - Ostream& os -) const -{ - writeHeader(os, header + " bins"); - writeHeaderValue(os, "bins", nBin_); - writeHeaderValue(os, "start", binMin_); - writeHeaderValue(os, "end", binMax_); - writeHeaderValue(os, "delta", binDx_); - writeHeaderValue(os, "direction", binDir_); - - vectorField binPoints(nBin_); - writeCommented(os, "x co-ords :"); - forAll(binPoints, pointi) - { - binPoints[pointi] = (binMin_ + (pointi + 1)*binDx_)*binDir_; - os << tab << binPoints[pointi].x(); - } - os << nl; - - writeCommented(os, "y co-ords :"); - forAll(binPoints, pointi) - { - os << tab << binPoints[pointi].y(); - } - os << nl; - - writeCommented(os, "z co-ords :"); - forAll(binPoints, pointi) - { - os << tab << binPoints[pointi].z(); - } - os << nl; - - writeHeader(os, ""); - writeCommented(os, "Time"); - - for (label j = 0; j < nBin_; j++) - { - const word jn(Foam::name(j) + ':'); - os << tab << jn << "(total_x total_y total_z)" - << tab << jn << "(pressure_x pressure_y pressure_z)" - << tab << jn << "(viscous_x viscous_y viscous_z)"; - - if (porosity_) - { - os << tab << jn << "(porous_x porous_y porous_z)"; - } - } - - os << endl; -} - - void Foam::functionObjects::forces::setCoordinateSystem ( const dictionary& dict, @@ -189,7 +88,6 @@ void Foam::functionObjects::forces::setCoordinateSystem coordSysPtr_.reset(new coordSystem::cartesian(dict)); } } - } @@ -215,126 +113,38 @@ void Foam::functionObjects::forces::initialise() ( !foundObject(UName_) || !foundObject(pName_) - ) { FatalErrorInFunction - << "Could not find U: " << UName_ << " or p:" << pName_ - << " in database" + << "Could not find U: " << UName_ + << " or p:" << pName_ << " in database" << exit(FatalError); } if (rhoName_ != "rhoInf" && !foundObject(rhoName_)) { FatalErrorInFunction - << "Could not find rho:" << rhoName_ + << "Could not find rho:" << rhoName_ << " in database" << exit(FatalError); } } - initialiseBins(); - initialised_ = true; } -void Foam::functionObjects::forces::initialiseBins() +void Foam::functionObjects::forces::reset() { - if (nBin_ > 1) - { - const polyBoundaryMesh& pbm = mesh_.boundaryMesh(); + sumPatchForcesP_ = Zero; + sumPatchForcesV_ = Zero; + sumPatchMomentsP_ = Zero; + sumPatchMomentsV_ = Zero; - // Determine extents of patches - scalar geomMin = GREAT; - scalar geomMax = -GREAT; - for (const label patchi : patchSet_) - { - const polyPatch& pp = pbm[patchi]; - scalarField d(pp.faceCentres() & binDir_); - geomMin = min(min(d), geomMin); - geomMax = max(max(d), geomMax); - } + sumInternalForces_ = Zero; + sumInternalMoments_ = Zero; - // Include porosity - if (porosity_) - { - const HashTable models = - obr_.lookupClass(); - - const scalarField dd(mesh_.C() & binDir_); - - forAllConstIters(models, iter) - { - const porosityModel& pm = *iter(); - const labelList& cellZoneIDs = pm.cellZoneIDs(); - - for (const label zonei : cellZoneIDs) - { - const cellZone& cZone = mesh_.cellZones()[zonei]; - const scalarField d(dd, cZone); - geomMin = min(min(d), geomMin); - geomMax = max(max(d), geomMax); - } - } - } - - reduce(geomMin, minOp()); - reduce(geomMax, maxOp()); - - // Slightly boost max so that region of interest is fully within bounds - geomMax = 1.0001*(geomMax - geomMin) + geomMin; - - // Use geometry limits if not specified by the user - if (binMin_ == GREAT) - { - binMin_ = geomMin; - } - if (binMax_ == GREAT) - { - binMax_ = geomMax; - } - - binDx_ = (binMax_ - binMin_)/scalar(nBin_); - - // Create the bin mid-points used for writing - binPoints_.setSize(nBin_); - forAll(binPoints_, i) - { - binPoints_[i] = (i + 0.5)*binDir_*binDx_; - } - } - - // Allocate storage for forces and moments - forAll(force_, i) - { - force_[i].setSize(nBin_, vector::zero); - moment_[i].setSize(nBin_, vector::zero); - } -} - - -void Foam::functionObjects::forces::resetFields() -{ - force_[0] = Zero; - force_[1] = Zero; - force_[2] = Zero; - - moment_[0] = Zero; - moment_[1] = Zero; - moment_[2] = Zero; - - if (writeFields_) - { - volVectorField& force = - lookupObjectRef(scopedName("force")); - - force == dimensionedVector(force.dimensions(), Zero); - - volVectorField& moment = - lookupObjectRef(scopedName("moment")); - - moment == dimensionedVector(moment.dimensions(), Zero); - } + force_ == dimensionedVector(force_.dimensions(), Zero); + moment_ == dimensionedVector(moment_.dimensions(), Zero); } @@ -344,42 +154,45 @@ Foam::functionObjects::forces::devRhoReff() const typedef compressible::turbulenceModel cmpTurbModel; typedef incompressible::turbulenceModel icoTurbModel; - const auto& U = lookupObject(UName_); - if (foundObject(cmpTurbModel::propertiesName)) { - const cmpTurbModel& turb = + const auto& turb = lookupObject(cmpTurbModel::propertiesName); - return turb.devRhoReff(U); + return turb.devRhoReff(); } else if (foundObject(icoTurbModel::propertiesName)) { - const incompressible::turbulenceModel& turb = + const auto& turb = lookupObject(icoTurbModel::propertiesName); - return rho()*turb.devReff(U); + return rho()*turb.devReff(); } else if (foundObject(fluidThermo::dictName)) { - const fluidThermo& thermo = - lookupObject(fluidThermo::dictName); + const auto& thermo = lookupObject(fluidThermo::dictName); + + const auto& U = lookupObject(UName_); return -thermo.mu()*dev(twoSymm(fvc::grad(U))); } else if (foundObject("transportProperties")) { - const transportModel& laminarT = + const auto& laminarT = lookupObject("transportProperties"); + const auto& U = lookupObject(UName_); + return -rho()*laminarT.nu()*dev(twoSymm(fvc::grad(U))); } else if (foundObject("transportProperties")) { - const dictionary& transportProperties = + const auto& transportProperties = lookupObject("transportProperties"); - dimensionedScalar nu("nu", dimViscosity, transportProperties); + const dimensionedScalar nu("nu", dimViscosity, transportProperties); + + const auto& U = lookupObject(UName_); return -rho()*nu*dev(twoSymm(fvc::grad(U))); } @@ -398,24 +211,23 @@ Foam::tmp Foam::functionObjects::forces::mu() const { if (foundObject(basicThermo::dictName)) { - const fluidThermo& thermo = - lookupObject(basicThermo::dictName); + const auto& thermo = lookupObject(basicThermo::dictName); return thermo.mu(); } else if (foundObject("transportProperties")) { - const transportModel& laminarT = + const auto& laminarT = lookupObject("transportProperties"); return rho()*laminarT.nu(); } else if (foundObject("transportProperties")) { - const dictionary& transportProperties = + const auto& transportProperties = lookupObject("transportProperties"); - dimensionedScalar nu("nu", dimViscosity, transportProperties); + const dimensionedScalar nu("nu", dimViscosity, transportProperties); return rho()*nu; } @@ -443,11 +255,11 @@ Foam::tmp Foam::functionObjects::forces::rho() const mesh_ ), mesh_, - dimensionedScalar("rho", dimDensity, rhoRef_) + dimensionedScalar(dimDensity, rhoRef_) ); } - return(lookupObject(rhoName_)); + return (lookupObject(rhoName_)); } @@ -455,7 +267,7 @@ Foam::scalar Foam::functionObjects::forces::rho(const volScalarField& p) const { if (p.dimensions() == dimPressure) { - return 1.0; + return 1; } if (rhoName_ != "rhoInf") @@ -469,225 +281,159 @@ Foam::scalar Foam::functionObjects::forces::rho(const volScalarField& p) const } -void Foam::functionObjects::forces::applyBins -( - const vectorField& Md, - const vectorField& fN, - const vectorField& fT, - const vectorField& fP, - const vectorField& d -) -{ - if (nBin_ == 1) - { - force_[0][0] += sum(fN); - force_[1][0] += sum(fT); - force_[2][0] += sum(fP); - moment_[0][0] += sum(Md^fN); - moment_[1][0] += sum(Md^fT); - moment_[2][0] += sum(Md^fP); - } - else - { - scalarField dd((d & binDir_) - binMin_); - - forAll(dd, i) - { - label bini = min(max(floor(dd[i]/binDx_), 0), force_[0].size() - 1); - - force_[0][bini] += fN[i]; - force_[1][bini] += fT[i]; - force_[2][bini] += fP[i]; - moment_[0][bini] += Md[i]^fN[i]; - moment_[1][bini] += Md[i]^fT[i]; - moment_[2][bini] += Md[i]^fP[i]; - } - } -} - - -void Foam::functionObjects::forces::addToFields +void Foam::functionObjects::forces::addToPatchFields ( const label patchi, const vectorField& Md, - const vectorField& fN, - const vectorField& fT, - const vectorField& fP + const vectorField& fP, + const vectorField& fV ) { - if (!writeFields_) - { - return; - } + sumPatchForcesP_ += sum(fP); + sumPatchForcesV_ += sum(fV); + force_.boundaryFieldRef()[patchi] += fP + fV; - auto& force = lookupObjectRef(scopedName("force")); - vectorField& pf = force.boundaryFieldRef()[patchi]; - pf += fN + fT + fP; + const vectorField mP(Md^fP); + const vectorField mV(Md^fV); - auto& moment = lookupObjectRef(scopedName("moment")); - vectorField& pm = moment.boundaryFieldRef()[patchi]; - pm = Md^pf; + sumPatchMomentsP_ += sum(mP); + sumPatchMomentsV_ += sum(mV); + moment_.boundaryFieldRef()[patchi] += mP + mV; } -void Foam::functionObjects::forces::addToFields +void Foam::functionObjects::forces::addToInternalField ( const labelList& cellIDs, const vectorField& Md, - const vectorField& fN, - const vectorField& fT, - const vectorField& fP + const vectorField& f ) { - if (!writeFields_) - { - return; - } - - auto& force = lookupObjectRef(scopedName("force")); - auto& moment = lookupObjectRef(scopedName("moment")); - forAll(cellIDs, i) { - label celli = cellIDs[i]; - force[celli] += fN[i] + fT[i] + fP[i]; - moment[celli] = Md[i]^force[celli]; + const label celli = cellIDs[i]; + + sumInternalForces_ += f[i]; + force_[celli] += f[i]; + + const vector m(Md[i]^f[i]); + sumInternalMoments_ += m; + moment_[celli] = m; } } -void Foam::functionObjects::forces::writeIntegratedForceMoment +void Foam::functionObjects::forces::createIntegratedDataFiles() +{ + if (!forceFilePtr_.valid()) + { + forceFilePtr_ = createFile("force"); + writeIntegratedDataFileHeader("Force", forceFilePtr_()); + } + + if (!momentFilePtr_.valid()) + { + momentFilePtr_ = createFile("moment"); + writeIntegratedDataFileHeader("Moment", momentFilePtr_()); + } +} + + +void Foam::functionObjects::forces::writeIntegratedDataFileHeader ( - const string& descriptor, - const vectorField& fm0, - const vectorField& fm1, - const vectorField& fm2, - autoPtr& osPtr + const word& header, + OFstream& os ) const { - vector pressure = sum(fm0); - vector viscous = sum(fm1); - vector porous = sum(fm2); - vector total = pressure + viscous + porous; - - Log << " Sum of " << descriptor.c_str() << nl - << " Total : " << total << nl - << " Pressure : " << pressure << nl - << " Viscous : " << viscous << nl; + const auto& coordSys = coordSysPtr_(); + const auto vecDesc = [](const word& root)->string + { + return root + "_x " + root + "_y " + root + "_z"; + }; + writeHeader(os, header); + writeHeaderValue(os, "CofR", coordSys.origin()); + writeHeader(os, ""); + writeCommented(os, "Time"); + writeTabbed(os, vecDesc("total")); + writeTabbed(os, vecDesc("pressure")); + writeTabbed(os, vecDesc("viscous")); if (porosity_) { - Log << " Porous : " << porous << nl; + writeTabbed(os, vecDesc("porous")); } - if (writeToFile()) - { - Ostream& os = osPtr(); - - writeCurrentTime(os); - - os << tab << total - << tab << pressure - << tab << viscous; - - if (porosity_) - { - os << tab << porous; - } - - os << endl; - } + os << endl; } -void Foam::functionObjects::forces::writeForces() +void Foam::functionObjects::forces::writeIntegratedDataFiles() { - Log << type() << " " << name() << " write:" << nl; - const auto& coordSys = coordSysPtr_(); - writeIntegratedForceMoment + writeIntegratedDataFile ( - "forces", - coordSys.localVector(force_[0]), - coordSys.localVector(force_[1]), - coordSys.localVector(force_[2]), - forceFilePtr_ + coordSys.localVector(sumPatchForcesP_), + coordSys.localVector(sumPatchForcesV_), + coordSys.localVector(sumInternalForces_), + forceFilePtr_() ); - writeIntegratedForceMoment + writeIntegratedDataFile ( - "moments", - coordSys.localVector(moment_[0]), - coordSys.localVector(moment_[1]), - coordSys.localVector(moment_[2]), - momentFilePtr_ + coordSys.localVector(sumPatchMomentsP_), + coordSys.localVector(sumPatchMomentsV_), + coordSys.localVector(sumInternalMoments_), + momentFilePtr_() ); - - Log << endl; } -void Foam::functionObjects::forces::writeBinnedForceMoment +void Foam::functionObjects::forces::writeIntegratedDataFile ( - const List>& fm, - autoPtr& osPtr + const vector& pres, + const vector& vis, + const vector& internal, + OFstream& os ) const { - if ((nBin_ == 1) || !writeToFile()) + writeCurrentTime(os); + + writeValue(os, pres + vis + internal); + writeValue(os, pres); + writeValue(os, vis); + + if (porosity_) + { + writeValue(os, internal); + } + + os << endl; +} + + +void Foam::functionObjects::forces::logIntegratedData +( + const string& descriptor, + const vector& pres, + const vector& vis, + const vector& internal +) const +{ + if (!log) { return; } - List> f(fm); + Log << " Sum of " << descriptor.c_str() << nl + << " Total : " << (pres + vis + internal) << nl + << " Pressure : " << pres << nl + << " Viscous : " << vis << nl; - if (binCumulative_) + if (porosity_) { - for (label i = 1; i < f[0].size(); i++) - { - f[0][i] += f[0][i-1]; - f[1][i] += f[1][i-1]; - f[2][i] += f[2][i-1]; - } + Log << " Porous : " << internal << nl; } - - Ostream& os = osPtr(); - - writeCurrentTime(os); - - forAll(f[0], i) - { - vector total = f[0][i] + f[1][i] + f[2][i]; - - os << tab << total - << tab << f[0][i] - << tab << f[1][i]; - - if (porosity_) - { - os << tab << f[2][i]; - } - } - - os << nl; -} - - -void Foam::functionObjects::forces::writeBins() -{ - const auto& coordSys = coordSysPtr_(); - - List> lf(3); - List> lm(3); - lf[0] = coordSys.localVector(force_[0]); - lf[1] = coordSys.localVector(force_[1]); - lf[2] = coordSys.localVector(force_[2]); - lm[0] = coordSys.localVector(moment_[0]); - lm[1] = coordSys.localVector(moment_[1]); - lm[2] = coordSys.localVector(moment_[2]); - - writeBinnedForceMoment(lf, forceBinFilePtr_); - writeBinnedForceMoment(lm, momentBinFilePtr_); } @@ -703,29 +449,50 @@ Foam::functionObjects::forces::forces : fvMeshFunctionObject(name, runTime, dict), writeFile(mesh_, name), - force_(3), - moment_(3), + force_ + ( + IOobject + ( + "force", // scopedName() is not available at ctor level + time_.timeName(), + mesh_, + IOobject::NO_READ, + IOobject::NO_WRITE + ), + mesh_, + dimensionedVector(dimForce, Zero) + ), + moment_ + ( + IOobject + ( + "moment", + time_.timeName(), + mesh_, + IOobject::NO_READ, + IOobject::NO_WRITE + ), + mesh_, + dimensionedVector(dimForce*dimLength, Zero) + ), + sumPatchForcesP_(Zero), + sumPatchForcesV_(Zero), + sumPatchMomentsP_(Zero), + sumPatchMomentsV_(Zero), + sumInternalForces_(Zero), + sumInternalMoments_(Zero), forceFilePtr_(), momentFilePtr_(), - forceBinFilePtr_(), - momentBinFilePtr_(), + coordSysPtr_(nullptr), patchSet_(), + rhoRef_(VGREAT), + pRef_(0), pName_("p"), UName_("U"), rhoName_("rho"), - directForceDensity_(false), fDName_("fD"), - rhoRef_(VGREAT), - pRef_(0), - coordSysPtr_(nullptr), + directForceDensity_(false), porosity_(false), - nBin_(1), - binDir_(Zero), - binDx_(0), - binMin_(GREAT), - binMax_(GREAT), - binPoints_(), - binCumulative_(true), writeFields_(false), initialised_(false) { @@ -748,29 +515,50 @@ Foam::functionObjects::forces::forces : fvMeshFunctionObject(name, obr, dict), writeFile(mesh_, name), - force_(3), - moment_(3), + force_ + ( + IOobject + ( + "force", + time_.timeName(), + mesh_, + IOobject::NO_READ, + IOobject::NO_WRITE + ), + mesh_, + dimensionedVector(dimForce, Zero) + ), + moment_ + ( + IOobject + ( + "moment", + time_.timeName(), + mesh_, + IOobject::NO_READ, + IOobject::NO_WRITE + ), + mesh_, + dimensionedVector(dimForce*dimLength, Zero) + ), + sumPatchForcesP_(Zero), + sumPatchForcesV_(Zero), + sumPatchMomentsP_(Zero), + sumPatchMomentsV_(Zero), + sumInternalForces_(Zero), + sumInternalMoments_(Zero), forceFilePtr_(), momentFilePtr_(), - forceBinFilePtr_(), - momentBinFilePtr_(), + coordSysPtr_(nullptr), patchSet_(), + rhoRef_(VGREAT), + pRef_(0), pName_("p"), UName_("U"), rhoName_("rho"), - directForceDensity_(false), fDName_("fD"), - rhoRef_(VGREAT), - pRef_(0), - coordSysPtr_(nullptr), + directForceDensity_(false), porosity_(false), - nBin_(1), - binDir_(Zero), - binDx_(0), - binMin_(GREAT), - binMax_(GREAT), - binPoints_(), - binCumulative_(true), writeFields_(false), initialised_(false) { @@ -780,17 +568,6 @@ Foam::functionObjects::forces::forces setCoordinateSystem(dict); Log << endl; } - -/* - // Turn off writing to file - writeToFile_ = false; - - forAll(force_, i) - { - force_[i].setSize(nBin_, vector::zero); - moment_[i].setSize(nBin_, vector::zero); - } -*/ } @@ -798,14 +575,14 @@ Foam::functionObjects::forces::forces bool Foam::functionObjects::forces::read(const dictionary& dict) { - fvMeshFunctionObject::read(dict); - writeFile::read(dict); + if (!fvMeshFunctionObject::read(dict) || !writeFile::read(dict)) + { + return false; + } initialised_ = false; - Info<< type() << " " << name() << ":" << nl; - - directForceDensity_ = dict.getOrDefault("directForceDensity", false); + Info<< type() << " " << name() << ":" << endl; patchSet_ = mesh_.boundaryMesh().patchSet @@ -813,9 +590,10 @@ bool Foam::functionObjects::forces::read(const dictionary& dict) dict.get("patches") ); + dict.readIfPresent("directForceDensity", directForceDensity_); if (directForceDensity_) { - // Optional entry for fDName + // Optional name entry for fD if (dict.readIfPresent("fD", fDName_)) { Info<< " fD: " << fDName_ << endl; @@ -840,7 +618,7 @@ bool Foam::functionObjects::forces::read(const dictionary& dict) // Reference density needed for incompressible calculations if (rhoName_ == "rhoInf") { - rhoRef_ = dict.get("rhoInf"); + rhoRef_ = dict.getCheck("rhoInf", scalarMinMax::ge(SMALL)); Info<< " Freestream density (rhoInf) set to " << rhoRef_ << endl; } @@ -861,115 +639,43 @@ bool Foam::functionObjects::forces::read(const dictionary& dict) Info<< " Not including porosity effects" << endl; } - if (dict.found("binData")) - { - Info<< " Activated data bins" << endl; - const dictionary& binDict(dict.subDict("binData")); - nBin_ = binDict.get