From 07c43e9c5d4ce20d7e9a9722c426b0e5b23a02e8 Mon Sep 17 00:00:00 2001 From: Will Bainbridge Date: Thu, 6 Jan 2022 15:13:17 +0000 Subject: [PATCH 1/3] setWriter, surfaceWriter: Added lists of fields to the variadic interface --- src/sampling/sampledSet/writers/setWriter.H | 111 +++++++++++++++----- 1 file changed, 85 insertions(+), 26 deletions(-) diff --git a/src/sampling/sampledSet/writers/setWriter.H b/src/sampling/sampledSet/writers/setWriter.H index 965d5b53f7..8f4aa2d5cd 100644 --- a/src/sampling/sampledSet/writers/setWriter.H +++ b/src/sampling/sampledSet/writers/setWriter.H @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | Website: https://openfoam.org - \\ / A nd | Copyright (C) 2011-2021 OpenFOAM Foundation + \\ / A nd | Copyright (C) 2011-2022 OpenFOAM Foundation \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -69,24 +69,26 @@ public: //- Helper for variadic write template - static inline void setTypeValueSet + static inline void appendTypeValueSet ( UPtrList>& valueSets, - const label valueSeti, - const Field& valueSet - ) - {} - - //- Helper for variadic write - template - static inline void setTypeValueSet - ( - UPtrList>& valueSets, - const label valueSeti, const Field& valueSet ) { - valueSets.set(valueSeti, &valueSet); + valueSets.resize(valueSets.size() + 1); + valueSets.set(valueSets.size() - 1, nullptr); + } + + //- Helper for variadic write + template + static inline void appendTypeValueSet + ( + UPtrList>& valueSets, + const Field& valueSet + ) + { + valueSets.resize(valueSets.size() + 1); + valueSets.set(valueSets.size() - 1, &valueSet); } //- Helper for variadic write @@ -114,14 +116,11 @@ public: Args& ... args ) { - const label valueSeti = - valueSetNames.size() - 1 - sizeof...(Args)/2; - - valueSetNames[valueSeti] = valueSetName; - #define SetTypeValueSet(Type, nullArg) \ - setTypeValueSet(Type##ValueSets, valueSeti, valueSet); - FOR_ALL_FIELD_TYPES(SetTypeValueSet); - #undef SetTypeValueSet + valueSetNames.append(valueSetName); + #define AppendTypeValueSet(Type, nullArg) \ + appendTypeValueSet(Type##ValueSets, valueSet); + FOR_ALL_FIELD_TYPES(AppendTypeValueSet); + #undef AppendTypeValueSet unpackTypeValueSets ( @@ -134,6 +133,68 @@ public: ); } + //- Helper for variadic write + template + static inline void unpackTypeValueSets + ( + wordList& valueSetNames + #define TypeValueSetsNonConstArg(Type, nullArg) \ + , UPtrList>& Type##ValueSets + FOR_ALL_FIELD_TYPES(TypeValueSetsNonConstArg), + #undef TypeValueSetsNonConstArg + const wordList& valueSetNamesPart, + const UPtrList>& valueSetsPart, + Args& ... args + ) + { + forAll(valueSetNamesPart, i) + { + valueSetNames.append(valueSetNamesPart[i]); + #define AppendTypeValueSet(Type, nullArg) \ + appendTypeValueSet(Type##ValueSets, valueSetsPart[i]); + FOR_ALL_FIELD_TYPES(AppendTypeValueSet); + #undef AppendTypeValueSet + } + + unpackTypeValueSets + ( + valueSetNames + #define TypeValueSetsParameter(Type, nullArg) \ + , Type##ValueSets + FOR_ALL_FIELD_TYPES(TypeValueSetsParameter), + #undef TypeValueSetsParameter + args ... + ); + } + + //- Helper for variadic write + template + static inline void unpackTypeValueSets + ( + wordList& valueSetNames + #define TypeValueSetsNonConstArg(Type, nullArg) \ + , UPtrList>& Type##ValueSets + FOR_ALL_FIELD_TYPES(TypeValueSetsNonConstArg), + #undef TypeValueSetsNonConstArg + const wordList& valueSetNamesPart, + const PtrList>& valueSetsPart, + Args& ... args + ) + { + unpackTypeValueSets + ( + valueSetNames + #define TypeValueSetsParameter(Type, nullArg) \ + , Type##ValueSets + FOR_ALL_FIELD_TYPES(TypeValueSetsParameter), + #undef TypeValueSetsParameter + valueSetNamesPart, + reinterpret_cast>&> + (valueSetsPart), + args ... + ); + } + protected: @@ -345,11 +406,9 @@ public: const Args& ... args ) const { - const label nFields = sizeof...(Args)/2; - - wordList valueSetNames(nFields); + wordList valueSetNames; #define DeclareTypeValueSets(Type, nullArg) \ - UPtrList> Type##ValueSets(nFields); + UPtrList> Type##ValueSets; FOR_ALL_FIELD_TYPES(DeclareTypeValueSets); #undef DeclareTypeValueSets From 36c565b9bf5d0f318e0e908c2097fd18f14d1e10 Mon Sep 17 00:00:00 2001 From: Will Bainbridge Date: Fri, 7 Jan 2022 08:42:33 +0000 Subject: [PATCH 2/3] multiphaseEulerFoam: new functionObject "moments" This function calculates integral (integer moments) or mean properties (mean, variance, standard deviation) of a size distribution computed with multiphaseEulerFoam. It has to be run with multiphaseEulerFoam, either at run-time or with -postProcess. It will not work with the postProcess utility. The following function object specification for example returns the first moment of the volume-based number density function which is equivalent to the phase fraction of the particulate phase: moments { type moments; libs ("libmultiphaseEulerFoamFunctionObjects.so"); executeControl timeStep; writeControl writeTime; populationBalance bubbles; momentType integerMoment; coordinateType volume; order 1; } The same can be achieved using a packaged function: #includeFunc moments ( populationBalance=bubbles, momentType=integerMoment, coordinateType=volume, order=1, funcName=moments ) Or on the command line: multiphaseEulerFoam -postProcess -func " moments ( populationBalance=bubbles, momentType=integerMoment, coordinateType=volume, order=1, funcName=moments )" Patch contributed by Institute of Fluid Dynamics, Helmholtz-Zentrum Dresden - Rossendorf (HZDR) --- .../functionObjects/Make/files | 1 + .../functionObjects/moments/moments.C | 972 ++++++++++++++++++ .../functionObjects/moments/moments.H | 258 +++++ .../postProcessing/multiphase/moments | 34 + .../{forces => multiphase}/phaseForces | 0 .../{fields => multiphase}/phaseMap | 0 .../populationBalance/drift/Allclean | 2 +- .../drift/system/controlDict | 37 +- .../drift/validation/createGraphs | 4 +- 9 files changed, 1291 insertions(+), 17 deletions(-) create mode 100644 applications/solvers/multiphase/multiphaseEulerFoam/functionObjects/moments/moments.C create mode 100644 applications/solvers/multiphase/multiphaseEulerFoam/functionObjects/moments/moments.H create mode 100644 etc/caseDicts/postProcessing/multiphase/moments rename etc/caseDicts/postProcessing/{forces => multiphase}/phaseForces (100%) rename etc/caseDicts/postProcessing/{fields => multiphase}/phaseMap (100%) diff --git a/applications/solvers/multiphase/multiphaseEulerFoam/functionObjects/Make/files b/applications/solvers/multiphase/multiphaseEulerFoam/functionObjects/Make/files index f7eb685454..4b72b7ff33 100644 --- a/applications/solvers/multiphase/multiphaseEulerFoam/functionObjects/Make/files +++ b/applications/solvers/multiphase/multiphaseEulerFoam/functionObjects/Make/files @@ -1,3 +1,4 @@ +moments/moments.C sizeDistribution/sizeDistribution.C phaseForces/phaseForces.C phaseMap/phaseMap.C diff --git a/applications/solvers/multiphase/multiphaseEulerFoam/functionObjects/moments/moments.C b/applications/solvers/multiphase/multiphaseEulerFoam/functionObjects/moments/moments.C new file mode 100644 index 0000000000..e8e5bbb4d1 --- /dev/null +++ b/applications/solvers/multiphase/multiphaseEulerFoam/functionObjects/moments/moments.C @@ -0,0 +1,972 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | Website: https://openfoam.org + \\ / A nd | Copyright (C) 2022 OpenFOAM Foundation + \\/ M anipulation | +------------------------------------------------------------------------------- +License + This file is part of OpenFOAM. + + OpenFOAM is free software: you can redistribute it and/or modify it + under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + OpenFOAM is distributed in the hope that it will be useful, but WITHOUT + ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + for more details. + + You should have received a copy of the GNU General Public License + along with OpenFOAM. If not, see . + +\*---------------------------------------------------------------------------*/ + +#include "moments.H" +#include "addToRunTimeSelectionTable.H" + +// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // + +namespace Foam +{ +namespace functionObjects +{ + defineTypeNameAndDebug(moments, 0); + addToRunTimeSelectionTable(functionObject, moments, dictionary); +} +} + + +namespace Foam +{ + template<> + const char* NamedEnum + < + Foam::functionObjects::moments::momentType, + 4 + >::names[] = {"integerMoment", "mean", "variance", "stdDev"}; +} + + +const Foam::NamedEnum +< + Foam::functionObjects::moments::momentType, + 4 +> +Foam::functionObjects::moments::momentTypeNames_; + + +namespace Foam +{ + template<> + const char* NamedEnum + < + Foam::functionObjects::moments::coordinateType, + 3 + >::names[] = {"volume", "area", "diameter"}; +} + + +const Foam::NamedEnum +< + Foam::functionObjects::moments::coordinateType, + 3 +> +Foam::functionObjects::moments::coordinateTypeNames_; + + +namespace Foam +{ + template<> + const char* NamedEnum + < + Foam::functionObjects::moments::weightType, + 3 + >::names[] = + { + "numberConcentration", + "volumeConcentration", + "areaConcentration" + }; +} + + +const Foam::NamedEnum +< + Foam::functionObjects::moments::weightType, + 3 +> +Foam::functionObjects::moments::weightTypeNames_; + + +namespace Foam +{ + template<> + const char* NamedEnum + < + Foam::functionObjects::moments::meanType, + 3 + >::names[] = {"arithmetic", "geometric", "notApplicable"}; +} + + +const Foam::NamedEnum +< + Foam::functionObjects::moments::meanType, + 3 +> +Foam::functionObjects::moments::meanTypeNames_; + + +// * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * // + +Foam::word Foam::functionObjects::moments::coordinateTypeSymbolicName() +{ + word coordinateTypeSymbolicName(word::null); + + switch (coordinateType_) + { + case coordinateType::volume: + { + coordinateTypeSymbolicName = "v"; + + break; + } + case coordinateType::area: + { + coordinateTypeSymbolicName = "a"; + + break; + } + case coordinateType::diameter: + { + coordinateTypeSymbolicName = "d"; + + break; + } + } + + return coordinateTypeSymbolicName; +} + + +Foam::word Foam::functionObjects::moments::weightTypeSymbolicName() +{ + word weightTypeSymbolicName(word::null); + + switch (weightType_) + { + case weightType::numberConcentration: + { + weightTypeSymbolicName = "N"; + + break; + } + case weightType::volumeConcentration: + { + weightTypeSymbolicName = "V"; + + break; + } + case weightType::areaConcentration: + { + weightTypeSymbolicName = "A"; + + break; + } + } + + return weightTypeSymbolicName; +} + + +Foam::word Foam::functionObjects::moments::defaultFldName() +{ + word meanName + ( + meanType_ == meanType::geometric + ? word(meanTypeNames_[meanType_]).capitalise() + : word("") + ); + + return + word + ( + IOobject::groupName + ( + "weighted" + + meanName + + word(momentTypeNames_[momentType_]).capitalise() + + "(" + + weightTypeSymbolicName() + + "," + + coordinateTypeSymbolicName() + + ")", + popBal_.name() + ) + ); +} + + +Foam::word Foam::functionObjects::moments::integerMomentFldName() +{ + return + word + ( + IOobject::groupName + ( + word(momentTypeNames_[momentType_]) + + Foam::name(order_) + + "(" + + weightTypeSymbolicName() + + "," + + coordinateTypeSymbolicName() + + ")", + popBal_.name() + ) + ); +} + + +void Foam::functionObjects::moments::setDimensions +( + volScalarField& fld, + momentType momType +) +{ + switch (momType) + { + case momentType::integerMoment: + { + switch (coordinateType_) + { + case coordinateType::volume: + { + fld.dimensions().reset + ( + pow(dimVolume, order_)/dimVolume + ); + + break; + } + case coordinateType::area: + { + fld.dimensions().reset + ( + pow(dimArea, order_)/dimVolume + ); + + break; + } + case coordinateType::diameter: + { + fld.dimensions().reset + ( + pow(dimLength, order_)/dimVolume + ); + + break; + } + } + + switch (weightType_) + { + case weightType::volumeConcentration: + { + fld.dimensions().reset(fld.dimensions()*dimVolume); + + break; + } + case weightType::areaConcentration: + { + fld.dimensions().reset(fld.dimensions()*dimArea); + + break; + } + default: + { + break; + } + } + + break; + } + case momentType::mean: + { + switch (coordinateType_) + { + case coordinateType::volume: + { + fld.dimensions().reset(dimVolume); + + break; + } + case coordinateType::area: + { + fld.dimensions().reset(dimArea); + + break; + } + case coordinateType::diameter: + { + fld.dimensions().reset(dimLength); + + break; + } + } + + break; + } + case momentType::variance: + { + switch (coordinateType_) + { + case coordinateType::volume: + { + fld.dimensions().reset(sqr(dimVolume)); + + break; + } + case coordinateType::area: + { + fld.dimensions().reset(sqr(dimArea)); + + break; + } + case coordinateType::diameter: + { + fld.dimensions().reset(sqr(dimLength)); + + break; + } + } + + if (meanType_ == meanType::geometric) + { + fld.dimensions().reset(dimless); + } + + break; + } + case momentType::stdDev: + { + switch (coordinateType_) + { + case coordinateType::volume: + { + fld.dimensions().reset(dimVolume); + + break; + } + case coordinateType::area: + { + fld.dimensions().reset(dimArea); + + break; + } + case coordinateType::diameter: + { + fld.dimensions().reset(dimLength); + + break; + } + } + + if (meanType_ == meanType::geometric) + { + fld.dimensions().reset(dimless); + } + + break; + } + } +} + + +Foam::tmp +Foam::functionObjects::moments::totalConcentration() +{ + tmp tTotalConcentration + ( + volScalarField::New + ( + "totalConcentration", + mesh_, + dimensionedScalar(inv(dimVolume), Zero) + ) + ); + + volScalarField& totalConcentration = tTotalConcentration.ref(); + + switch (weightType_) + { + case weightType::volumeConcentration: + { + totalConcentration.dimensions().reset + ( + totalConcentration.dimensions()*dimVolume + ); + + break; + } + case weightType::areaConcentration: + { + totalConcentration.dimensions().reset + ( + totalConcentration.dimensions()*dimArea + ); + + break; + } + default: + { + break; + } + } + + forAll(popBal_.sizeGroups(), i) + { + const Foam::diameterModels::sizeGroup& fi = popBal_.sizeGroups()[i]; + + switch (weightType_) + { + case weightType::numberConcentration: + { + totalConcentration += fi*fi.phase()/fi.x(); + + break; + } + case weightType::volumeConcentration: + { + totalConcentration += fi*fi.phase(); + + break; + } + case weightType::areaConcentration: + { + totalConcentration += fi.a()*fi*fi.phase()/fi.x(); + + break; + } + } + } + + return tTotalConcentration; +} + + +Foam::tmp Foam::functionObjects::moments::mean() +{ + tmp tMean + ( + volScalarField::New + ( + "mean", + mesh_, + dimensionedScalar(dimless, Zero) + ) + ); + + volScalarField& mean = tMean.ref(); + + setDimensions(mean, momentType::mean); + + volScalarField totalConcentration(this->totalConcentration()); + + forAll(popBal_.sizeGroups(), i) + { + const Foam::diameterModels::sizeGroup& fi = popBal_.sizeGroups()[i]; + + volScalarField concentration(fi*fi.phase()/fi.x()); + + switch (weightType_) + { + case weightType::volumeConcentration: + { + concentration *= fi.x(); + + break; + } + case weightType::areaConcentration: + { + concentration *= fi.a(); + + break; + } + default: + { + break; + } + } + + switch (meanType_) + { + case meanType::geometric: + { + mean.dimensions().reset(dimless); + + switch (coordinateType_) + { + case coordinateType::volume: + { + dimensionedScalar unitVolume(dimVolume, 1); + + mean += + Foam::log(fi.x()/unitVolume) + *concentration/totalConcentration; + + break; + } + case coordinateType::area: + { + dimensionedScalar unitArea(dimArea, 1); + + mean += + Foam::log(fi.a()/unitArea) + *concentration/totalConcentration; + + break; + } + case coordinateType::diameter: + { + dimensionedScalar unitLength(dimLength, 1); + + mean += + Foam::log(fi.d()/unitLength) + *concentration/totalConcentration; + + break; + } + } + + break; + } + default: + { + switch (coordinateType_) + { + case coordinateType::volume: + { + mean += fi.x()*concentration/totalConcentration; + + break; + } + case coordinateType::area: + { + mean += fi.a()*concentration/totalConcentration; + + break; + } + case coordinateType::diameter: + { + mean += fi.d()*concentration/totalConcentration; + + break; + } + } + + break; + } + } + } + + if (meanType_ == meanType::geometric) + { + mean = exp(mean); + + setDimensions(mean, momentType::mean); + } + + return tMean; +} + + +Foam::tmp Foam::functionObjects::moments::variance() +{ + tmp tVariance + ( + volScalarField::New + ( + "variance", + mesh_, + dimensionedScalar(dimless, Zero) + ) + ); + + volScalarField& variance = tVariance.ref(); + + setDimensions(variance, momentType::variance); + + volScalarField totalConcentration(this->totalConcentration()); + volScalarField mean(this->mean()); + + forAll(popBal_.sizeGroups(), i) + { + const Foam::diameterModels::sizeGroup& fi = popBal_.sizeGroups()[i]; + + volScalarField concentration(fi*fi.phase()/fi.x()); + + switch (weightType_) + { + case weightType::volumeConcentration: + { + concentration *= fi.x(); + + break; + } + case weightType::areaConcentration: + { + concentration *= fi.a(); + + break; + } + default: + { + break; + } + } + + switch (meanType_) + { + case meanType::geometric: + { + switch (coordinateType_) + { + case coordinateType::volume: + { + variance += + sqr(Foam::log(fi.x()/mean)) + *concentration/totalConcentration; + + break; + } + case coordinateType::area: + { + variance += + sqr(Foam::log(fi.a()/mean)) + *concentration/totalConcentration; + + break; + } + case coordinateType::diameter: + { + variance += + sqr(Foam::log(fi.d()/mean)) + *concentration/totalConcentration; + + break; + } + } + + break; + } + default: + { + switch (coordinateType_) + { + case coordinateType::volume: + { + variance += + sqr(fi.x() - mean)*concentration/totalConcentration; + + + break; + } + case coordinateType::area: + { + variance += + sqr(fi.a() - mean)*concentration/totalConcentration; + + break; + } + case coordinateType::diameter: + { + variance += + sqr(fi.d() - mean)*concentration/totalConcentration; + + break; + } + } + + break; + } + } + } + + return tVariance; +} + + +Foam::tmp +Foam::functionObjects::moments::stdDev() +{ + switch (meanType_) + { + case meanType::geometric: + { + return exp(sqrt(this->variance())); + } + default: + { + return sqrt(this->variance()); + } + } +} + + +// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // + +Foam::functionObjects::moments::moments +( + const word& name, + const Time& runTime, + const dictionary& dict +) +: + fvMeshFunctionObject(name, runTime, dict), + popBal_ + ( + obr_.lookupObject + ( + dict.lookup("populationBalance") + ) + ), + momentType_(momentTypeNames_.read(dict.lookup("momentType"))), + coordinateType_(coordinateTypeNames_.read(dict.lookup("coordinateType"))), + weightType_ + ( + dict.found("weight") + ? weightTypeNames_.read(dict.lookup("weightType")) + : weightType::numberConcentration + ), + meanType_(meanType::notApplicable), + order_(-1), + fldPtr_(nullptr) +{ + read(dict); +} + + +// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * // + +Foam::functionObjects::moments::~moments() +{} + + +// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // + +bool Foam::functionObjects::moments::read(const dictionary& dict) +{ + fvMeshFunctionObject::read(dict); + + switch (momentType_) + { + case momentType::integerMoment: + { + order_ = dict.lookup("order"); + + break; + } + default: + { + meanType_ = + dict.found("meanType") + ? meanTypeNames_.read(dict.lookup("meanType")) + : meanType::arithmetic; + + break; + } + } + + switch (momentType_) + { + case momentType::integerMoment: + { + fldPtr_.set + ( + new volScalarField + ( + IOobject + ( + this->integerMomentFldName(), + mesh_.time().timeName(), + mesh_, + IOobject::NO_READ, + IOobject::NO_WRITE + ), + mesh_, + dimensionedScalar(dimless, Zero) + ) + ); + + volScalarField& integerMoment = fldPtr_(); + + setDimensions(integerMoment, momentType::integerMoment); + + break; + } + case momentType::mean: + { + fldPtr_.set + ( + new volScalarField + ( + IOobject + ( + this->defaultFldName(), + mesh_.time().timeName(), + mesh_, + IOobject::NO_READ, + IOobject::NO_WRITE + ), + this->mean() + ) + ); + + break; + } + case momentType::variance: + { + fldPtr_.set + ( + new volScalarField + ( + IOobject + ( + this->defaultFldName(), + mesh_.time().timeName(), + mesh_, + IOobject::NO_READ, + IOobject::NO_WRITE + ), + this->variance() + ) + ); + + break; + } + case momentType::stdDev: + { + fldPtr_.set + ( + new volScalarField + ( + IOobject + ( + this->defaultFldName(), + mesh_.time().timeName(), + mesh_, + IOobject::NO_READ, + IOobject::NO_WRITE + ), + this->stdDev() + ) + ); + + break; + } + } + + return true; +} + + +bool Foam::functionObjects::moments::execute() +{ + switch (momentType_) + { + case momentType::integerMoment: + { + volScalarField& integerMoment = fldPtr_(); + + integerMoment = Zero; + + forAll(popBal_.sizeGroups(), i) + { + const Foam::diameterModels::sizeGroup& fi = + popBal_.sizeGroups()[i]; + + volScalarField concentration(fi*fi.phase()/fi.x()); + + switch (weightType_) + { + case weightType::volumeConcentration: + { + concentration *= fi.x(); + + break; + } + case weightType::areaConcentration: + { + concentration *= fi.a(); + + break; + } + default: + { + break; + } + } + + switch (coordinateType_) + { + case coordinateType::volume: + { + integerMoment += + pow(fi.x(), order_)*concentration; + + break; + } + case coordinateType::area: + { + integerMoment += + pow(fi.a(), order_)*concentration; + + break; + } + case coordinateType::diameter: + { + integerMoment += + pow(fi.d(), order_)*concentration; + + break; + } + } + } + + break; + } + case momentType::mean: + { + fldPtr_() = this->mean(); + + break; + } + case momentType::variance: + { + fldPtr_() = this->variance(); + + break; + } + case momentType::stdDev: + { + fldPtr_() = sqrt(this->variance()); + + break; + } + } + + return true; +} + + +bool Foam::functionObjects::moments::write() +{ + writeObject(fldPtr_->name()); + + return true; +} + + +// ************************************************************************* // diff --git a/applications/solvers/multiphase/multiphaseEulerFoam/functionObjects/moments/moments.H b/applications/solvers/multiphase/multiphaseEulerFoam/functionObjects/moments/moments.H new file mode 100644 index 0000000000..c337a0fc85 --- /dev/null +++ b/applications/solvers/multiphase/multiphaseEulerFoam/functionObjects/moments/moments.H @@ -0,0 +1,258 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | Website: https://openfoam.org + \\ / A nd | Copyright (C) 2022 OpenFOAM Foundation + \\/ M anipulation | +------------------------------------------------------------------------------- +License + This file is part of OpenFOAM. + + OpenFOAM is free software: you can redistribute it and/or modify it + under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + OpenFOAM is distributed in the hope that it will be useful, but WITHOUT + ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + for more details. + + You should have received a copy of the GNU General Public License + along with OpenFOAM. If not, see . + +Class + Foam::functionObjects::moments + +Description + Calculates and writes out integral (integer moments) or mean properties + (mean, variance, standard deviation) of a size distribution computed with + multiphaseEulerFoam. Requires solver post-processing. + + The following function object specification for example returns the first + moment of the volume-based number density function which is equivalent to + the phase fraction of the particulate phase: + + \verbatim + moments + { + type moments; + libs ("libmultiphaseEulerFoamFunctionObjects.so"); + executeControl timeStep; + writeControl writeTime; + populationBalance bubbles; + momentType integerMoment; + coordinateType volume; + order 1; + } + \endverbatim + +Usage + \table + Property | Description | Required | Default + populationBalance | population balance name | yes | + momentType | desired moment of the distribution\\ + | yes | + coordinateType | particle property | yes | + weightType | number/volume/area concentration\\ + | no\\ + | numberConcentration + order | order of integer moment | for integer moments | + meanType | arithmetic or geometric | for non-integer moments\\ + | arithmetic + \endtable + +See also + Foam::diameterModels::populationBalanceModel + Foam::functionObjects::fvMeshFunctionObject + Foam::functionObject + +SourceFiles + moments.C + +\*---------------------------------------------------------------------------*/ + +#ifndef functionObjects_moments_H +#define functionObjects_moments_H + +#include "fvMeshFunctionObject.H" +#include "populationBalanceModel.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ +namespace functionObjects +{ + +/*---------------------------------------------------------------------------*\ + Class moments Declaration +\*---------------------------------------------------------------------------*/ + +class moments +: + public fvMeshFunctionObject +{ +public: + + //- Enumeration for the moment types + enum class momentType + { + integerMoment, + mean, + variance, + stdDev + }; + + //- Names of the moment types + static const NamedEnum momentTypeNames_; + + //- Enumeration for the coordinate types + enum class coordinateType + { + volume, + area, + diameter + }; + + //- Names of the coordinate types + static const NamedEnum coordinateTypeNames_; + + //- Enumeration for the weight types + enum class weightType + { + numberConcentration, + volumeConcentration, + areaConcentration + }; + + //- Names of the weight types + static const NamedEnum weightTypeNames_; + + //- Enumeration for the mean types + enum class meanType + { + arithmetic, + geometric, + notApplicable + }; + + //- Names of the mean types + static const NamedEnum meanTypeNames_; + + +private: + + // Private Data + + //- Reference to population balance + const Foam::diameterModels::populationBalanceModel& popBal_; + + //- Moment type + momentType momentType_; + + //- Coordinate type + coordinateType coordinateType_; + + //- Weight type + weightType weightType_; + + //- Mean type + meanType meanType_; + + //- Integer moment order + int order_; + + //- Result field + autoPtr fldPtr_; + + + // Private Member Functions + + //- Coordinate type symbolic name for shorter field names + word coordinateTypeSymbolicName(); + + //- Weight type symbolic name for shorter field names + word weightTypeSymbolicName(); + + //- Default field name + word defaultFldName(); + + //- Integer moment field name + word integerMomentFldName(); + + //- Set dimensions + void setDimensions(volScalarField& fld, momentType momType); + + //- Total concentration + tmp totalConcentration(); + + //- Mean value + tmp mean(); + + //- Variance + tmp variance(); + + //- Standard deviation + tmp stdDev(); + + +public: + + //- Runtime type information + TypeName("moments"); + + + // Constructors + + //- Construct from Time and dictionary + moments + ( + const word& name, + const Time& runTime, + const dictionary& + ); + + //- Disallow default bitwise copy construction + moments(const moments&) = delete; + + + //- Destructor + virtual ~moments(); + + + // Member Functions + + //- Read the data + virtual bool read(const dictionary&); + + //- Return the list of fields required + virtual wordList fields() const + { + return wordList::null(); + } + + //- Calculate the force fields + virtual bool execute(); + + //- Write the force fields + virtual bool write(); + + + // Member Operators + + //- Disallow default bitwise assignment + void operator=(const moments&) = delete; +}; + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace functionObjects +} // End namespace Foam + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#endif + +// ************************************************************************* // diff --git a/etc/caseDicts/postProcessing/multiphase/moments b/etc/caseDicts/postProcessing/multiphase/moments new file mode 100644 index 0000000000..912b856118 --- /dev/null +++ b/etc/caseDicts/postProcessing/multiphase/moments @@ -0,0 +1,34 @@ +/*--------------------------------*- C++ -*----------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | Website: https://openfoam.org + \\ / A nd | Version: dev + \\/ M anipulation | +------------------------------------------------------------------------------- +Description + Calculates and writes out integral (integer moments) or mean properties + (mean, variance, standard deviation) of a size distribution computed with + multiphaseEulerFoam. Requires solver post-processing. + +\*---------------------------------------------------------------------------*/ + +type moments; +libs ("libmultiphaseEulerFoamFunctionObjects.so"); + +populationBalance ; +momentType ; // integerMoment, mean, variance, + // stdDev +coordinateType ; // volume, area, diameter +weightType numberConcentration; // volumeConcentration, + // areaConcentration + // defaults to numberConcentration +order 0; // relevant for integer moments only +mean arithmetic; // geometric + // relevant for non-integer moments, + // defaults to arithmetic + +executeControl timeStep; +writeControl writeTime; + + +// ************************************************************************* // diff --git a/etc/caseDicts/postProcessing/forces/phaseForces b/etc/caseDicts/postProcessing/multiphase/phaseForces similarity index 100% rename from etc/caseDicts/postProcessing/forces/phaseForces rename to etc/caseDicts/postProcessing/multiphase/phaseForces diff --git a/etc/caseDicts/postProcessing/fields/phaseMap b/etc/caseDicts/postProcessing/multiphase/phaseMap similarity index 100% rename from etc/caseDicts/postProcessing/fields/phaseMap rename to etc/caseDicts/postProcessing/multiphase/phaseMap diff --git a/test/multiphase/multiphaseEulerFoam/populationBalance/drift/Allclean b/test/multiphase/multiphaseEulerFoam/populationBalance/drift/Allclean index fa9dc48dcc..70690e08bb 100755 --- a/test/multiphase/multiphaseEulerFoam/populationBalance/drift/Allclean +++ b/test/multiphase/multiphaseEulerFoam/populationBalance/drift/Allclean @@ -6,6 +6,6 @@ cd "${0%/*}" || exit 1 # Source clean functions . $WM_PROJECT_DIR/bin/tools/CleanFunctions -cleanCase && rm -f moments.eps numberDensity.eps +cleanCase && rm -rf *.eps 0/d.air 0/uniform/ 0/integerMoment* #------------------------------------------------------------------------------ diff --git a/test/multiphase/multiphaseEulerFoam/populationBalance/drift/system/controlDict b/test/multiphase/multiphaseEulerFoam/populationBalance/drift/system/controlDict index 780adb1fe2..ae82fc9175 100644 --- a/test/multiphase/multiphaseEulerFoam/populationBalance/drift/system/controlDict +++ b/test/multiphase/multiphaseEulerFoam/populationBalance/drift/system/controlDict @@ -28,7 +28,7 @@ deltaT 0.01; writeControl runTime; -writeInterval 6; +writeInterval 0.2; purgeWrite 0; @@ -68,21 +68,30 @@ functions densityFunction yes; } - moments - { - type sizeDistribution; - functionObjectLibs ("libmultiphaseEulerFoamFunctionObjects.so"); + #includeFunc moments + ( + populationBalance=bubbles, + momentType=integerMoment, + coordinateType=volume, + order=1 + ) - writeControl runTime; - writeInterval 0.1; + #includeFunc moments + ( + populationBalance=bubbles, + momentType=integerMoment, + coordinateType=volume, + order=0 + ) - setFormat raw; - - populationBalance bubbles; - functionType moments; - coordinateType volume; - maxOrder 1; - } + #includeFunc probes + ( + funcName=probes, + points=((0.5 0.5 0.5)), + integerMoment0(N,v).bubbles, + integerMoment1(N,v).bubbles, + writeControl=writeTime + ) } // ************************************************************************* // diff --git a/test/multiphase/multiphaseEulerFoam/populationBalance/drift/validation/createGraphs b/test/multiphase/multiphaseEulerFoam/populationBalance/drift/validation/createGraphs index 1520907a1b..924d77fbac 100755 --- a/test/multiphase/multiphaseEulerFoam/populationBalance/drift/validation/createGraphs +++ b/test/multiphase/multiphaseEulerFoam/populationBalance/drift/validation/createGraphs @@ -36,9 +36,9 @@ gnuplot< Date: Fri, 7 Jan 2022 08:47:29 +0000 Subject: [PATCH 3/3] multiphaseEulerFoam: revised sizeDistribution functionObject Following the addition of the new moments functionObject, all related functionality was removed from sizeDistribution. In its revised version, sizeDistribution allows for different kinds of weighted region averaging in case of field-dependent representative particle properties. A packaged function has also been added to allow for command line solver post-processing. For example, the following function object specification returns the volume-based number density function: numberDensity { type sizeDistribution; libs ("libmultiphaseEulerFoamFunctionObjects.so"); writeControl writeTime; populationBalance bubbles; functionType numberDensity; coordinateType volume; setFormat raw; } The same can be achieved using a packaged function: #includeFunc sizeDistribution ( populationBalance=bubbles, functionType=numberDensity, coordinateType=volume, funcName=numberDensity ) Or on the command line: multiphaseEulerFoam -postProcess -func " sizeDistribution ( populationBalance=bubbles, functionType=numberDensity, coordinateType=volume, funcName=numberDensity )" Patch contributed by Institute of Fluid Dynamics, Helmholtz-Zentrum Dresden - Rossendorf (HZDR) --- .../sizeDistribution/sizeDistribution.C | 693 +++++++++++------- .../sizeDistribution/sizeDistribution.H | 193 +++-- .../multiphase/sizeDistribution | 41 ++ .../binaryBreakup/system/controlDict | 22 +- .../binaryBreakup/validation/createGraphs | 6 +- .../breakup/system/controlDict | 22 +- .../breakup/validation/createGraphs | 6 +- .../coalescence/system/controlDict | 22 +- .../coalescence/validation/createGraphs | 6 +- .../drift/system/controlDict | 22 +- .../drift/validation/createGraphs | 2 +- .../isothermalGrowth/validation/d.bubbles | 299 -------- .../populationBalance/negativeDrift/Allclean | 2 +- .../negativeDrift/system/controlDict | 59 +- .../negativeDrift/validation/createGraphs | 6 +- .../system/controlDict | 48 +- .../validation/createGraphs | 12 +- .../RAS/bubblePipe/system/controlDict | 56 +- .../RAS/bubblePipe/validation/createGraphs | 4 +- .../system/controlDict.orig | 6 +- .../system/controlDict.orig | 6 +- .../titaniaSynthesis/system/controlDict | 32 +- .../titaniaSynthesis/validation/createGraphs | 8 +- .../system/controlDict | 32 +- .../validation/createGraphs | 8 +- 25 files changed, 731 insertions(+), 882 deletions(-) create mode 100644 etc/caseDicts/postProcessing/multiphase/sizeDistribution delete mode 100644 test/multiphase/multiphaseEulerFoam/populationBalance/isothermalGrowth/validation/d.bubbles diff --git a/applications/solvers/multiphase/multiphaseEulerFoam/functionObjects/sizeDistribution/sizeDistribution.C b/applications/solvers/multiphase/multiphaseEulerFoam/functionObjects/sizeDistribution/sizeDistribution.C index 02b761dc47..c3a6bb68a8 100644 --- a/applications/solvers/multiphase/multiphaseEulerFoam/functionObjects/sizeDistribution/sizeDistribution.C +++ b/applications/solvers/multiphase/multiphaseEulerFoam/functionObjects/sizeDistribution/sizeDistribution.C @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | Website: https://openfoam.org - \\ / A nd | Copyright (C) 2017-2021 OpenFOAM Foundation + \\ / A nd | Copyright (C) 2017-2022 OpenFOAM Foundation \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -24,10 +24,7 @@ License \*---------------------------------------------------------------------------*/ #include "sizeDistribution.H" -#include "Time.H" -#include "fvMesh.H" #include "addToRunTimeSelectionTable.H" -#include "mathematicalConstants.H" // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // @@ -45,19 +42,21 @@ const char* Foam::NamedEnum < Foam::functionObjects::sizeDistribution::functionType, - 4 + 6 >::names[] = { - "moments", - "standardDeviation", - "number", - "volume" + "numberConcentration", + "numberDensity", + "volumeConcentration", + "volumeDensity", + "areaConcentration", + "areaDensity" }; const Foam::NamedEnum < Foam::functionObjects::sizeDistribution::functionType, - 4 + 6 > Foam::functionObjects::sizeDistribution::functionTypeNames_; template<> @@ -80,10 +79,122 @@ const Foam::NamedEnum 4 > Foam::functionObjects::sizeDistribution::coordinateTypeNames_; + +namespace Foam +{ + template<> + const char* NamedEnum + < + Foam::functionObjects::sizeDistribution::weightType, + 4 + >::names[] = + { + "numberConcentration", + "volumeConcentration", + "areaConcentration", + "cellVolume" + }; +} + + +const Foam::NamedEnum +< + Foam::functionObjects::sizeDistribution::weightType, + 4 +> +Foam::functionObjects::sizeDistribution::weightTypeNames_; + using Foam::constant::mathematical::pi; -// * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * // +// * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * // + +Foam::word Foam::functionObjects::sizeDistribution::functionTypeSymbolicName() +{ + word functionTypeSymbolicName(word::null); + + switch (functionType_) + { + case functionType::numberConcentration: + { + functionTypeSymbolicName = "N"; + + break; + } + case functionType::numberDensity: + { + functionTypeSymbolicName = "n"; + + break; + } + case functionType::volumeConcentration: + { + functionTypeSymbolicName = "V"; + + break; + } + case functionType::volumeDensity: + { + functionTypeSymbolicName = "v"; + + break; + } + case functionType::areaConcentration: + { + functionTypeSymbolicName = "A"; + + break; + } + case functionType::areaDensity: + { + functionTypeSymbolicName = "a"; + + break; + } + } + + return functionTypeSymbolicName; +} + + +Foam::word Foam::functionObjects::sizeDistribution::coordinateTypeSymbolicName +( + const coordinateType& cType +) +{ + word coordinateTypeSymbolicName(word::null); + + switch (cType) + { + case coordinateType::volume: + { + coordinateTypeSymbolicName = "v"; + + break; + } + case coordinateType::area: + { + coordinateTypeSymbolicName = "a"; + + break; + } + case coordinateType::diameter: + { + coordinateTypeSymbolicName = "d"; + + break; + } + case coordinateType::projectedAreaDiameter: + { + coordinateTypeSymbolicName = "dPa"; + + break; + } + } + + return coordinateTypeSymbolicName; +} + Foam::tmp Foam::functionObjects::sizeDistribution::filterField @@ -102,252 +213,120 @@ Foam::functionObjects::sizeDistribution::filterField } -void Foam::functionObjects::sizeDistribution::correctVolAverages() +Foam::scalar Foam::functionObjects::sizeDistribution::averageCoordinateValue +( + const Foam::diameterModels::sizeGroup& fi, + const coordinateType& cType +) { - forAll(N_, i) + scalar averageCoordinateValue(Zero); + + switch (cType) { - const Foam::diameterModels::sizeGroup& fi = popBal_.sizeGroups()[i]; - - scalarField Ni(filterField(fi*fi.phase()/fi.x())); - scalarField V(filterField(mesh_.V())); - scalarField ai(filterField(fi.a())); - scalarField di(Ni.size()); - - switch (coordinateType_) + case coordinateType::volume: { - case ctProjectedAreaDiameter: - { - di = sqrt(ai/pi); - - break; - } - default: - { - di = fi.d(); - - break; - } - } - - N_[i] = gSum(V*Ni)/this->V(); - V_[i] = fi.x().value(); - a_[i] = gSum(V*ai)/this->V(); - d_[i] = gSum(V*di)/this->V(); - } -} - - -void Foam::functionObjects::sizeDistribution::writeMoments() -{ - logFiles::write(); - - Log << " writing moments of size distribution." << endl; - - if (Pstream::master()) - { - writeTime(file()); - } - - const scalarField& bin = this->bin(); - - for (label k = 0; k <= maxOrder_; k++) - { - scalar result = 0; - - forAll(N_, i) - { - result += pow(bin[i], k)*N_[i]; - } - - if (Pstream::master()) - { - file() << tab << result; - } - } - - if (Pstream::master()) - { - file() << endl; - } -} - - -void Foam::functionObjects::sizeDistribution::writeStdDev() -{ - logFiles::write(); - - Log << " writing standard deviation of size distribution." - << endl; - - if (Pstream::master()) - { - writeTime(file()); - } - - scalar stdDev = 0; - scalar mean = 0; - scalar var = 0; - - const scalarField& bin = this->bin(); - - if (sum(N_) != 0) - { - if (geometric_) - { - mean = exp(sum(Foam::log(bin)*N_/sum(N_))); - - var = - sum(sqr(Foam::log(bin) - Foam::log(mean)) - *N_/sum(N_)); - - stdDev = exp(sqrt(var)); - } - else - { - mean = sum(bin*N_/sum(N_)); - - var = sum(sqr(bin - mean)*N_/sum(N_)); - - stdDev = sqrt(var); - } - } - - if (Pstream::master()) - { - file() << tab << stdDev << tab << mean << tab << var << endl; - } -} - - -void Foam::functionObjects::sizeDistribution::writeDistribution() -{ - scalarField result(N_); - - const scalarField& bin = this->bin(); - - switch (functionType_) - { - case ftNumber: - { - Log << " writing number distribution. " - << endl; + averageCoordinateValue = fi.x().value(); break; } - case ftVolume: + case coordinateType::area: { - Log << " writing volume distribution. " - << endl; - - result *= V_; + averageCoordinateValue = + weightedAverage(fi.a(), fi); break; } - default: + case coordinateType::diameter: { + averageCoordinateValue = + weightedAverage(fi.d(), fi); + break; } + case coordinateType::projectedAreaDiameter: + { + averageCoordinateValue = + weightedAverage(sqrt(fi.a()/pi), fi); + break; + } } - if (normalise_) + return averageCoordinateValue; +} + + +Foam::scalar Foam::functionObjects::sizeDistribution::weightedAverage +( + const Foam::scalarField& fld, + const Foam::diameterModels::sizeGroup& fi +) +{ + scalar weightedAverage(Zero); + + switch (weightType_) { - if(sum(result) != 0) + case weightType::numberConcentration: { - result /= sum(result); - } + scalarField Ni(filterField(fi*fi.phase()/fi.x().value())); - } - - if (densityFunction_) - { - List bndrs(N_.size() + 1); - - bndrs.first() = bin.first(); - bndrs.last() = bin.last(); - - for (label i = 1; i < N_.size(); i++) - { - bndrs[i] = (bin[i]+ bin[i-1])/2.0; - } - - forAll(result, i) - { - if (geometric_) + if (gSum(Ni) == 0) { - result[i] /= - (Foam::log(bndrs[i+1]) - Foam::log(bndrs[i])); + weightedAverage = + gSum(filterField(mesh_.V()*fld))/this->V(); } else { - result[i] /= (bndrs[i+1] - bndrs[i]); + weightedAverage = + gSum(Ni*filterField(fld))/gSum(Ni); } + + break; } - } - - if (Pstream::master()) - { - formatterPtr_->write - ( - file_.baseTimeDir(), - name(), - coordSet(true, "volume", V_), - "area", - a_, - "diameter", - d_, - word(functionTypeNames_[functionType_]) - + (densityFunction_ ? "Density" : "Concentration"), - result - ); - } -} - - -void Foam::functionObjects::sizeDistribution::writeFileHeader -( - const label i -) -{ - volRegion::writeFileHeader(*this, file()); - - writeHeaderValue - ( - file(), - "Coordinate", - word(coordinateTypeNames_[coordinateType_]) - ); - - word str("Time"); - - switch (functionType_) - { - case ftMoments: + case weightType::volumeConcentration: { - for (label k = 0; k <= maxOrder_; k++) + scalarField Vi(filterField(fi*fi.phase())); + + if (gSum(Vi) == 0) { - str += (" k=" + std::to_string(k)); + weightedAverage = + gSum(filterField(mesh_.V()*fld))/this->V(); + } + else + { + weightedAverage = + gSum(Vi*filterField(fld))/gSum(Vi); } break; } - - case ftStdDev: + case weightType::areaConcentration: { - str += " standardDeviation mean variance"; + scalarField Ai(filterField(fi.a().ref()*fi.phase())); + + if (gSum(Ai) == 0) + { + weightedAverage = + gSum(filterField(mesh_.V()*fld))/this->V(); + } + else + { + weightedAverage = + gSum(Ai*filterField(fld))/gSum(Ai); + } break; } - - default: + case weightType::cellVolume: { + weightedAverage = + gSum(filterField(mesh_.V()*fld))/this->V(); + break; } } - writeCommented(file(), str); - - file() << endl; + return weightedAverage; } @@ -362,9 +341,8 @@ Foam::functionObjects::sizeDistribution::sizeDistribution : fvMeshFunctionObject(name, runTime, dict), volRegion(fvMeshFunctionObject::mesh_, dict), - logFiles(obr_, name), - mesh_(fvMeshFunctionObject::mesh_), file_(obr_, name), + mesh_(fvMeshFunctionObject::mesh_), popBal_ ( obr_.lookupObject @@ -374,10 +352,26 @@ Foam::functionObjects::sizeDistribution::sizeDistribution ), functionType_(functionTypeNames_.read(dict.lookup("functionType"))), coordinateType_(coordinateTypeNames_.read(dict.lookup("coordinateType"))), - N_(popBal_.sizeGroups().size(), 0), - V_(popBal_.sizeGroups().size(), 0), - a_(popBal_.sizeGroups().size(), 0), - d_(popBal_.sizeGroups().size(), 0) + allCoordinates_ + ( + dict.lookupOrDefault("allCoordinates", false) + ), + normalise_(dict.lookupOrDefault("normalise", false)), + logTransform_ + ( + dict.lookupOrDefaultBackwardsCompatible + ( + {"logTransform", "geometric"}, + false + ) + ), + weightType_ + ( + dict.found("weightType") + ? weightTypeNames_.read(dict.lookup("weightType")) + : weightType::numberConcentration + ), + formatterPtr_(nullptr) { read(dict); } @@ -393,17 +387,12 @@ Foam::functionObjects::sizeDistribution::~sizeDistribution() bool Foam::functionObjects::sizeDistribution::read(const dictionary& dict) { + Log << type() << " " << name() << ":" << nl; + fvMeshFunctionObject::read(dict); - normalise_ = dict.lookupOrDefault("normalise", false); - densityFunction_ = dict.lookupOrDefault("densityFunction", false); - geometric_ = dict.lookupOrDefault("geometric", false); - maxOrder_ = dict.lookupOrDefault("maxOrder", 3); - formatterPtr_ = setWriter::New(dict.lookup("setFormat"), dict); - resetName(name()); - return false; } @@ -414,43 +403,255 @@ bool Foam::functionObjects::sizeDistribution::execute() } -bool Foam::functionObjects::sizeDistribution::end() -{ - return true; -} - - bool Foam::functionObjects::sizeDistribution::write() { Log << type() << " " << name() << " write:" << nl; - correctVolAverages(); + const UPtrList& sizeGroups = + popBal_.sizeGroups(); + + scalarField coordinateValues(sizeGroups.size()); + scalarField boundaryValues(sizeGroups.size() + 1); + scalarField resultValues(sizeGroups.size()); + + forAll(sizeGroups, i) + { + const diameterModels::sizeGroup& fi = sizeGroups[i]; + + coordinateValues[i] = averageCoordinateValue(fi, coordinateType_); + } + + if + ( + functionType_ == functionType::numberDensity + || functionType_ == functionType::volumeDensity + || functionType_ == functionType::areaDensity + ) + { + if (logTransform_) + { + boundaryValues.first() = Foam::log(coordinateValues.first()); + boundaryValues.last() = Foam::log(coordinateValues.last()); + + for (label i = 1; i < boundaryValues.size() - 1; i++) + { + boundaryValues[i] = + 0.5 + *( + Foam::log(coordinateValues[i]) + + Foam::log(coordinateValues[i-1]) + ); + } + } + else + { + boundaryValues.first() = coordinateValues.first(); + boundaryValues.last() = coordinateValues.last(); + + for (label i = 1; i < boundaryValues.size() - 1; i++) + { + boundaryValues[i] = + (coordinateValues[i] + coordinateValues[i-1])/2; + } + } + } switch (functionType_) { - case ftMoments: + case functionType::numberConcentration: { - writeMoments(); + forAll(sizeGroups, i) + { + const diameterModels::sizeGroup& fi = sizeGroups[i]; + + resultValues[i] = + gSum(filterField(mesh_.V()*fi*fi.phase()/fi.x()))*this->V(); + } + + if (normalise_ && sum(resultValues) != 0) + { + resultValues /= sum(resultValues); + } break; } - - case ftStdDev: + case functionType::numberDensity: { - writeStdDev(); + forAll(sizeGroups, i) + { + const diameterModels::sizeGroup& fi = sizeGroups[i]; + + resultValues[i] = + gSum(filterField(mesh_.V()*fi*fi.phase()/fi.x()))*this->V(); + } + + if (normalise_ && sum(resultValues) != 0) + { + resultValues /= sum(resultValues); + } + + forAll(resultValues, i) + { + resultValues[i] /= (boundaryValues[i+1] - boundaryValues[i]); + } break; } - - default: + case functionType::volumeConcentration: { - writeDistribution(); + forAll(sizeGroups, i) + { + const diameterModels::sizeGroup& fi = sizeGroups[i]; + + resultValues[i] = + gSum(filterField(mesh_.V()*fi*fi.phase()))*this->V(); + } + + if (normalise_ && sum(resultValues) != 0) + { + resultValues /= sum(resultValues); + } + + break; + } + case functionType::volumeDensity: + { + forAll(sizeGroups, i) + { + const diameterModels::sizeGroup& fi = sizeGroups[i]; + + resultValues[i] = + gSum(filterField(mesh_.V()*fi*fi.phase()))*this->V(); + } + + if (normalise_ && sum(resultValues) != 0) + { + resultValues /= sum(resultValues); + } + + forAll(resultValues, i) + { + resultValues[i] /= (boundaryValues[i+1] - boundaryValues[i]); + } + + break; + } + case functionType::areaConcentration: + { + forAll(sizeGroups, i) + { + const diameterModels::sizeGroup& fi = sizeGroups[i]; + + resultValues[i] = + gSum + ( + filterField(mesh_.V()*fi.a().ref()*fi*fi.phase()/fi.x()) + ) + *this->V(); + } + + if (normalise_ && sum(resultValues) != 0) + { + resultValues /= sum(resultValues); + } + + break; + } + case functionType::areaDensity: + { + forAll(sizeGroups, i) + { + const diameterModels::sizeGroup& fi = sizeGroups[i]; + + resultValues[i] = + gSum + ( + filterField(mesh_.V()*fi.a().ref()*fi*fi.phase()/fi.x()) + ) + *this->V(); + } + + if (normalise_ && sum(resultValues) != 0) + { + resultValues /= sum(resultValues); + } + + forAll(resultValues, i) + { + resultValues[i] /= (boundaryValues[i+1] - boundaryValues[i]); + } break; } } - Log << endl; + + if (allCoordinates_) + { + wordList otherCoordinateSymbolicNames(coordinateTypeNames_.size()); + PtrList otherCoordinateValues(coordinateTypeNames_.size()); + typedef NamedEnum namedEnumCoordinateType; + + forAllConstIter(namedEnumCoordinateType, coordinateTypeNames_, iter) + { + const coordinateType cType = coordinateTypeNames_[iter.key()]; + + otherCoordinateSymbolicNames[cType] = + coordinateTypeSymbolicName(cType); + + otherCoordinateValues.set + ( + cType, + new scalarField(popBal_.sizeGroups().size()) + ); + + forAll(sizeGroups, i) + { + const diameterModels::sizeGroup& fi = sizeGroups[i]; + + otherCoordinateValues[cType][i] = + averageCoordinateValue(fi, cType); + } + } + + if (Pstream::master()) + { + formatterPtr_->write + ( + file_.baseTimeDir(), + name(), + coordSet + ( + true, + coordinateTypeSymbolicName(coordinateType_), + coordinateValues + ), + functionTypeSymbolicName(), + resultValues, + otherCoordinateSymbolicNames, + otherCoordinateValues + ); + } + } + else + { + if (Pstream::master()) + { + formatterPtr_->write + ( + file_.baseTimeDir(), + name(), + coordSet + ( + true, + coordinateTypeSymbolicName(coordinateType_), + coordinateValues + ), + functionTypeSymbolicName(), + resultValues + ); + } + } return true; } diff --git a/applications/solvers/multiphase/multiphaseEulerFoam/functionObjects/sizeDistribution/sizeDistribution.H b/applications/solvers/multiphase/multiphaseEulerFoam/functionObjects/sizeDistribution/sizeDistribution.H index ea684746cc..b08ec2773b 100644 --- a/applications/solvers/multiphase/multiphaseEulerFoam/functionObjects/sizeDistribution/sizeDistribution.H +++ b/applications/solvers/multiphase/multiphaseEulerFoam/functionObjects/sizeDistribution/sizeDistribution.H @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | Website: https://openfoam.org - \\ / A nd | Copyright (C) 2017-2021 OpenFOAM Foundation + \\ / A nd | Copyright (C) 2017-2022 OpenFOAM Foundation \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -25,42 +25,42 @@ Class Foam::functionObjects::sizeDistribution Description - This function object calculates and outputs volume-averaged information - about the size distribution of the dispersed phase, such as the number - density function or its moments. It is designed to be used exclusively with - the population balance modeling functionality of the multiphaseEulerFoam - solver. It can be applied to a specific cellZone or the entire domain. The - function type determines whether the density function and its moments are - based on the number of dispersed phase elements in a size group or their - total volume. + Writes out the size distribution computed with multiphaseEulerFoam for the + entire domain or a volume region. Requires solver post-processing. + + The following function object specification for example returns the volume- + based number density function: Example of function object specification: \verbatim numberDensity { - type sizeDistribution; - libs ("libmultiphaseEulerFoamFunctionObjects.so"); - ... + type sizeDistribution; + libs ("libmultiphaseEulerFoamFunctionObjects.so"); + writeControl writeTime; populationBalance bubbles; - regionType cellZone; - name zone0; - functionType number; + functionType numberDensity; coordinateType volume; - densityFunction yes; + setFormat raw; } \endverbatim Usage \table - Property | Description | Required | Default value - type | type name: sizeDistribution | yes | - populationBalance | corresponding populationBalance | yes | - functionType | number/volume/moments/stdDev | yes | - coordinateType | used for density/moment calculation | yes | - normalise | normalise concentrations | no | no - densityFunction | compute densityFunction | no | no - logBased | use log of coordinate for density | no | no - maxOrder | maxim order of moment output | no | 3 + Property | Description | Required | Default + populationBalance | population balance name | yes | + functionType | function type | yes | + coordinateType | particle property | yes | + allCoordinates | write all coordinate values | no | false + normalise | divide by total concentration | no | false + logTransform | class width based on log of coordinate\\ + | no | false + weightType | weighting in case of field-dependent particle\\ + properties | no\\ + | numberConcentration + regionType | cellZone or all | no | all + name | name of cellZone if required | no | + setFormat | output format | yes | \endtable SourceFiles @@ -73,7 +73,6 @@ SourceFiles #include "fvMeshFunctionObject.H" #include "volRegion.H" -#include "logFiles.H" #include "populationBalanceModel.H" #include "writeFile.H" #include "setWriter.H" @@ -92,8 +91,7 @@ namespace functionObjects class sizeDistribution : public fvMeshFunctionObject, - public volRegion, - public logFiles + public volRegion { public: @@ -102,111 +100,101 @@ public: //- Function type enumeration enum functionType { - ftMoments, - ftStdDev, - ftNumber, - ftVolume + numberConcentration, + numberDensity, + volumeConcentration, + volumeDensity, + areaConcentration, + areaDensity }; - //- Ordinate type names - static const NamedEnum functionTypeNames_; + //- Function type names + static const NamedEnum functionTypeNames_; //- Coordinate type enumeration enum coordinateType { - ctVolume, - ctArea, - ctDiameter, - ctProjectedAreaDiameter + volume, + area, + diameter, + projectedAreaDiameter }; //- Coordinate type names static const NamedEnum coordinateTypeNames_; + //- Enumeration for the weight types + enum class weightType + { + numberConcentration, + volumeConcentration, + areaConcentration, + cellVolume + }; -protected: + //- Names of the weight types + static const NamedEnum weightTypeNames_; - // Protected Data - //- Reference to fvMesh - const fvMesh& mesh_; +private: - //- File containing data for all functionTypes except moments + // Private Data + + //- Write file writeFile file_; - //- Output formatter - autoPtr formatterPtr_; + //- Reference to mesh + const fvMesh& mesh_; - //- Reference to populationBalanceModel + //- Reference to population balance const Foam::diameterModels::populationBalanceModel& popBal_; - //- Function to evaluate + //- Function type functionType functionType_; - //- Abscissa type + //- Coordinate type coordinateType coordinateType_; - //- List of volume-averaged number concentrations - scalarField N_; + //- Add values for all coordinate types to output + Switch allCoordinates_; - //- ??? - scalarField V_; - - //- List of volume-averaged surface areas - scalarField a_; - - //- List of volume-averaged diameters - scalarField d_; - - //- Normalise number/volume concentrations + //- Normalise result through division by sum Switch normalise_; - //- Determines whether density function is calculated - Switch densityFunction_; + //- Log transform + Switch logTransform_; - //- Geometric standard deviation/density function - Switch geometric_; + //- Weight types, relevant if particle properties are field dependent + weightType weightType_; - //- Highest moment order - label maxOrder_; + //- Set formatter + autoPtr formatterPtr_; - // Protected Member Functions + // Private Member Functions - //- Filter field according to cellIds + //- Function type symbolic name for shorter file header + word functionTypeSymbolicName(); + + //- Coordinate type symbolic name for shorter file header + word coordinateTypeSymbolicName(const coordinateType& cType); + + //- Filter a field according to cellIds tmp filterField(const scalarField& field) const; - //- Bin component used according to chosen coordinate type - inline const scalarField& bin() const - { - switch (coordinateType_) - { - case ctVolume: - return V_; - case ctArea: - return a_; - case ctDiameter: - return d_; - case ctProjectedAreaDiameter: - return d_; - } - return scalarField::null(); - } + //- Field averaged coordinate value + scalar averageCoordinateValue + ( + const diameterModels::sizeGroup& fi, + const coordinateType& cType + ); - //- Correct volume averages - void correctVolAverages(); - - //- Write moments - void writeMoments(); - - //- Write standard deviation - void writeStdDev(); - - //- Write distribution - void writeDistribution(); - - //- Output file header information for functionType moments - virtual void writeFileHeader(const label i); + //- Weighted average + scalar weightedAverage + ( + const scalarField& fld, + const diameterModels::sizeGroup& fi + ); public: @@ -235,21 +223,18 @@ public: // Member Functions - //- Read the sizeDistribution data - virtual bool read(const dictionary&); - //- Return the list of fields required virtual wordList fields() const { return wordList::null(); } + //- Read the sizeDistribution data + virtual bool read(const dictionary&); + //- Execute, currently does nothing virtual bool execute(); - //- Execute at the final time-loop, currently does nothing - virtual bool end(); - //- Calculate and write the size distribution virtual bool write(); diff --git a/etc/caseDicts/postProcessing/multiphase/sizeDistribution b/etc/caseDicts/postProcessing/multiphase/sizeDistribution new file mode 100644 index 0000000000..71e9034a92 --- /dev/null +++ b/etc/caseDicts/postProcessing/multiphase/sizeDistribution @@ -0,0 +1,41 @@ +/*--------------------------------*- C++ -*----------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | Website: https://openfoam.org + \\ / A nd | Version: dev + \\/ M anipulation | +------------------------------------------------------------------------------- +Description + Writes out the size distribution computed with multiphaseEulerFoam for the + entire domain or a volume region. Requires solver post-processing. + +\*---------------------------------------------------------------------------*/ + +type sizeDistribution; +libs ("libmultiphaseEulerFoamFunctionObjects.so"); + +populationBalance ; +functionType ; // numberConcentration, numberDensity + // volumeConcentration, volumeDensity + // areaConcentration, areaDensity +coordinateType ; // volume, area, diameter, + // projectedAreaDiameter +allCoordinates false; // defaults to false +normalise false; // defaults to false +logTransform false; // defaults to false, only relevant + // for density functions +weightType numberConcentration; // volumeConcentration, + // areaConcentration, cellVolume + // relevant for field-dependent + // particle properties, defaults to + // numberConcentration +regionType all; // cellZone + // defaults to all +name cellZoneName; // relevant for regionType all + + +setFormat raw; +writeControl writeTime; + + +// ************************************************************************* // diff --git a/test/multiphase/multiphaseEulerFoam/populationBalance/binaryBreakup/system/controlDict b/test/multiphase/multiphaseEulerFoam/populationBalance/binaryBreakup/system/controlDict index f52e70d8d2..e49e4664ff 100644 --- a/test/multiphase/multiphaseEulerFoam/populationBalance/binaryBreakup/system/controlDict +++ b/test/multiphase/multiphaseEulerFoam/populationBalance/binaryBreakup/system/controlDict @@ -52,21 +52,13 @@ maxDeltaT 1; functions { - numberDensity - { - type sizeDistribution; - functionObjectLibs ("libmultiphaseEulerFoamFunctionObjects.so"); - - writeControl outputTime; - writeInterval 1; - - setFormat raw; - - populationBalance bubbles; - functionType number; - coordinateType volume; - densityFunction yes; - } + #includeFunc sizeDistribution + ( + populationBalance=bubbles, + functionType=numberDensity, + coordinateType=volume, + funcName=numberDensity + ) } // ************************************************************************* // diff --git a/test/multiphase/multiphaseEulerFoam/populationBalance/binaryBreakup/validation/createGraphs b/test/multiphase/multiphaseEulerFoam/populationBalance/binaryBreakup/validation/createGraphs index 8b7909d56c..a3fd7df831 100755 --- a/test/multiphase/multiphaseEulerFoam/populationBalance/binaryBreakup/validation/createGraphs +++ b/test/multiphase/multiphaseEulerFoam/populationBalance/binaryBreakup/validation/createGraphs @@ -31,9 +31,9 @@ gnuplot< -50 -( -0.00300040563 -0.00300110488 -0.00300180018 -0.00300249358 -0.00300318666 -0.00300387818 -0.00300457055 -0.00300526257 -0.00300595539 -0.00300664879 -0.00300734302 -0.00300803811 -0.00300873419 -0.00300943126 -0.00301012937 -0.00301082851 -0.00301152868 -0.00301222987 -0.00301293204 -0.00301363515 -0.00301433916 -0.00301504402 -0.00301574966 -0.00301645602 -0.00301716301 -0.00301787056 -0.00301857858 -0.00301928697 -0.00301999563 -0.00302070448 -0.00302141342 -0.00302212238 -0.00302283129 -0.00302354009 -0.00302424876 -0.00302495726 -0.0030256656 -0.00302637377 -0.00302708181 -0.00302778974 -0.0030284976 -0.00302920541 -0.00302991328 -0.00303062095 -0.00303132989 -0.00303203401 -0.00303275902 -0.00303342196 -0.00303447928 -0.00303514582 -) -; - -boundaryField -{ - inlet - { - type calculated; - value uniform 0.00300000033; - } - outlet - { - type calculated; - value uniform 0.00303514582; - } - walls - { - type calculated; - value nonuniform List -200 -( -0.00300040563 -0.00300110488 -0.00300180018 -0.00300249358 -0.00300318666 -0.00300387818 -0.00300457055 -0.00300526257 -0.00300595539 -0.00300664879 -0.00300734302 -0.00300803811 -0.00300873419 -0.00300943126 -0.00301012937 -0.00301082851 -0.00301152868 -0.00301222987 -0.00301293204 -0.00301363515 -0.00301433916 -0.00301504402 -0.00301574966 -0.00301645602 -0.00301716301 -0.00301787056 -0.00301857858 -0.00301928697 -0.00301999563 -0.00302070448 -0.00302141342 -0.00302212238 -0.00302283129 -0.00302354009 -0.00302424876 -0.00302495726 -0.0030256656 -0.00302637377 -0.00302708181 -0.00302778974 -0.0030284976 -0.00302920541 -0.00302991328 -0.00303062095 -0.00303132989 -0.00303203401 -0.00303275902 -0.00303342196 -0.00303447928 -0.00303514582 -0.00300040563 -0.00300110488 -0.00300180018 -0.00300249358 -0.00300318666 -0.00300387818 -0.00300457055 -0.00300526257 -0.00300595539 -0.00300664879 -0.00300734302 -0.00300803811 -0.00300873419 -0.00300943126 -0.00301012937 -0.00301082851 -0.00301152868 -0.00301222987 -0.00301293204 -0.00301363515 -0.00301433916 -0.00301504402 -0.00301574966 -0.00301645602 -0.00301716301 -0.00301787056 -0.00301857858 -0.00301928697 -0.00301999563 -0.00302070448 -0.00302141342 -0.00302212238 -0.00302283129 -0.00302354009 -0.00302424876 -0.00302495726 -0.0030256656 -0.00302637377 -0.00302708181 -0.00302778974 -0.0030284976 -0.00302920541 -0.00302991328 -0.00303062095 -0.00303132989 -0.00303203401 -0.00303275902 -0.00303342196 -0.00303447928 -0.00303514582 -0.00300040563 -0.00300110488 -0.00300180018 -0.00300249358 -0.00300318666 -0.00300387818 -0.00300457055 -0.00300526257 -0.00300595539 -0.00300664879 -0.00300734302 -0.00300803811 -0.00300873419 -0.00300943126 -0.00301012937 -0.00301082851 -0.00301152868 -0.00301222987 -0.00301293204 -0.00301363515 -0.00301433916 -0.00301504402 -0.00301574966 -0.00301645602 -0.00301716301 -0.00301787056 -0.00301857858 -0.00301928697 -0.00301999563 -0.00302070448 -0.00302141342 -0.00302212238 -0.00302283129 -0.00302354009 -0.00302424876 -0.00302495726 -0.0030256656 -0.00302637377 -0.00302708181 -0.00302778974 -0.0030284976 -0.00302920541 -0.00302991328 -0.00303062095 -0.00303132989 -0.00303203401 -0.00303275902 -0.00303342196 -0.00303447928 -0.00303514582 -0.00300040563 -0.00300110488 -0.00300180018 -0.00300249358 -0.00300318666 -0.00300387818 -0.00300457055 -0.00300526257 -0.00300595539 -0.00300664879 -0.00300734302 -0.00300803811 -0.00300873419 -0.00300943126 -0.00301012937 -0.00301082851 -0.00301152868 -0.00301222987 -0.00301293204 -0.00301363515 -0.00301433916 -0.00301504402 -0.00301574966 -0.00301645602 -0.00301716301 -0.00301787056 -0.00301857858 -0.00301928697 -0.00301999563 -0.00302070448 -0.00302141342 -0.00302212238 -0.00302283129 -0.00302354009 -0.00302424876 -0.00302495726 -0.0030256656 -0.00302637377 -0.00302708181 -0.00302778974 -0.0030284976 -0.00302920541 -0.00302991328 -0.00303062095 -0.00303132989 -0.00303203401 -0.00303275902 -0.00303342196 -0.00303447928 -0.00303514582 -) -; - } -} - - -// ************************************************************************* // diff --git a/test/multiphase/multiphaseEulerFoam/populationBalance/negativeDrift/Allclean b/test/multiphase/multiphaseEulerFoam/populationBalance/negativeDrift/Allclean index fa9dc48dcc..70690e08bb 100755 --- a/test/multiphase/multiphaseEulerFoam/populationBalance/negativeDrift/Allclean +++ b/test/multiphase/multiphaseEulerFoam/populationBalance/negativeDrift/Allclean @@ -6,6 +6,6 @@ cd "${0%/*}" || exit 1 # Source clean functions . $WM_PROJECT_DIR/bin/tools/CleanFunctions -cleanCase && rm -f moments.eps numberDensity.eps +cleanCase && rm -rf *.eps 0/d.air 0/uniform/ 0/integerMoment* #------------------------------------------------------------------------------ diff --git a/test/multiphase/multiphaseEulerFoam/populationBalance/negativeDrift/system/controlDict b/test/multiphase/multiphaseEulerFoam/populationBalance/negativeDrift/system/controlDict index b16b25fcfd..829675d7d3 100644 --- a/test/multiphase/multiphaseEulerFoam/populationBalance/negativeDrift/system/controlDict +++ b/test/multiphase/multiphaseEulerFoam/populationBalance/negativeDrift/system/controlDict @@ -28,7 +28,7 @@ deltaT 0.01; writeControl runTime; -writeInterval 4; +writeInterval 0.2; purgeWrite 0; @@ -52,37 +52,38 @@ maxDeltaT 1; functions { - numberDensity - { - type sizeDistribution; - functionObjectLibs ("libmultiphaseEulerFoamFunctionObjects.so"); + #includeFunc sizeDistribution + ( + populationBalance=bubbles, + functionType=numberDensity, + coordinateType=volume, + funcName=numberDensity + ) - writeControl outputTime; - writeInterval 1; + #includeFunc moments + ( + populationBalance=bubbles, + momentType=integerMoment, + coordinateType=volume, + order=1 + ) - setFormat raw; + #includeFunc moments + ( + populationBalance=bubbles, + momentType=integerMoment, + coordinateType=volume, + order=0 + ) - populationBalance bubbles; - functionType number; - coordinateType volume; - densityFunction yes; - } - - moments - { - type sizeDistribution; - functionObjectLibs ("libmultiphaseEulerFoamFunctionObjects.so"); - - writeControl runTime; - writeInterval 0.1; - - setFormat raw; - - populationBalance bubbles; - functionType moments; - coordinateType volume; - maxOrder 1; - } + #includeFunc probes + ( + funcName=probes, + points=((0.5 0.5 0.5)), + integerMoment0(N,v).bubbles, + integerMoment1(N,v).bubbles, + writeControl=writeTime + ) } // ************************************************************************* // diff --git a/test/multiphase/multiphaseEulerFoam/populationBalance/negativeDrift/validation/createGraphs b/test/multiphase/multiphaseEulerFoam/populationBalance/negativeDrift/validation/createGraphs index a635a60ae5..be0d849e47 100755 --- a/test/multiphase/multiphaseEulerFoam/populationBalance/negativeDrift/validation/createGraphs +++ b/test/multiphase/multiphaseEulerFoam/populationBalance/negativeDrift/validation/createGraphs @@ -23,7 +23,7 @@ gnuplot<