From 45881823522b587d59f0e0a944e0a509cf62fa93 Mon Sep 17 00:00:00 2001 From: Mark Olesen <Mark.Olesen@esi-group.com> Date: Fri, 11 Aug 2017 14:53:24 +0200 Subject: [PATCH] ENH: add absolute weighting for surfaceFieldValue (issue #567) - can be useful either for flow-rate weighting where backflow is to be ignored in the average, or for flow-rate weighting on surfaces with inconsistent orientation. Reworked to code to make better use of Enum (the NamedEnum replacement). Enum doesn't require contiguous enumeration values, which lets us use bitmasking of similar operations to reduce duplicate code. --- .../surfaceFieldValue/surfaceFieldValue.C | 129 +++++++++-------- .../surfaceFieldValue/surfaceFieldValue.H | 130 ++++++++++++----- .../surfaceFieldValue/surfaceFieldValueI.H | 17 ++- .../surfaceFieldValueTemplates.C | 136 ++++++++++-------- .../fieldValues/volFieldValue/volFieldValue.C | 65 ++++----- .../fieldValues/volFieldValue/volFieldValue.H | 59 +++++--- .../volFieldValue/volFieldValueTemplates.C | 85 +++++------ 7 files changed, 361 insertions(+), 260 deletions(-) diff --git a/src/functionObjects/field/fieldValues/surfaceFieldValue/surfaceFieldValue.C b/src/functionObjects/field/fieldValues/surfaceFieldValue/surfaceFieldValue.C index f0eb481dd77..db92766baee 100644 --- a/src/functionObjects/field/fieldValues/surfaceFieldValue/surfaceFieldValue.C +++ b/src/functionObjects/field/fieldValues/surfaceFieldValue/surfaceFieldValue.C @@ -69,23 +69,32 @@ const Foam::Enum > Foam::functionObjects::fieldValues::surfaceFieldValue::operationTypeNames_ { + // Normal operations { operationType::opNone, "none" }, + { operationType::opMin, "min" }, + { operationType::opMax, "max" }, { operationType::opSum, "sum" }, - { operationType::opWeightedSum, "weightedSum" }, { operationType::opSumMag, "sumMag" }, { operationType::opSumDirection, "sumDirection" }, { operationType::opSumDirectionBalance, "sumDirectionBalance" }, { operationType::opAverage, "average" }, - { operationType::opWeightedAverage, "weightedAverage" }, { operationType::opAreaAverage, "areaAverage" }, - { operationType::opWeightedAreaAverage, "weightedAreaAverage" }, { operationType::opAreaIntegrate, "areaIntegrate" }, - { operationType::opWeightedAreaIntegrate, "weightedAreaIntegrate" }, - { operationType::opMin, "min" }, - { operationType::opMax, "max" }, { operationType::opCoV, "CoV" }, { operationType::opAreaNormalAverage, "areaNormalAverage" }, { operationType::opAreaNormalIntegrate, "areaNormalIntegrate" }, + + // Using weighting + { operationType::opWeightedSum, "weightedSum" }, + { operationType::opWeightedAverage, "weightedAverage" }, + { operationType::opWeightedAreaAverage, "weightedAreaAverage" }, + { operationType::opWeightedAreaIntegrate, "weightedAreaIntegrate" }, + + // Using absolute weighting + { operationType::opAbsWeightedSum, "absWeightedSum" }, + { operationType::opAbsWeightedAverage, "absWeightedAverage" }, + { operationType::opAbsWeightedAreaAverage, "absWeightedAreaAverage" }, + { operationType::opAbsWeightedAreaIntegrate, "absWeightedAreaIntegrate" }, }; const Foam::Enum @@ -244,7 +253,7 @@ void Foam::functionObjects::fieldValues::surfaceFieldValue::combineMeshGeometry { if (facePatchId_[i] != -1) { - label patchi = facePatchId_[i]; + const label patchi = facePatchId_[i]; globalFacesIs[i] += mesh_.boundaryMesh()[patchi].start(); } } @@ -279,9 +288,8 @@ void Foam::functionObjects::fieldValues::surfaceFieldValue::combineMeshGeometry // My own data first { const faceList& fcs = allFaces[Pstream::myProcNo()]; - forAll(fcs, i) + for (const face& f : fcs) { - const face& f = fcs[i]; face& newF = faces[nFaces++]; newF.setSize(f.size()); forAll(f, fp) @@ -291,9 +299,9 @@ void Foam::functionObjects::fieldValues::surfaceFieldValue::combineMeshGeometry } const pointField& pts = allPoints[Pstream::myProcNo()]; - forAll(pts, i) + for (const point& pt : pts) { - points[nPoints++] = pts[i]; + points[nPoints++] = pt; } } @@ -303,9 +311,8 @@ void Foam::functionObjects::fieldValues::surfaceFieldValue::combineMeshGeometry if (proci != Pstream::myProcNo()) { const faceList& fcs = allFaces[proci]; - forAll(fcs, i) + for (const face& f : fcs) { - const face& f = fcs[i]; face& newF = faces[nFaces++]; newF.setSize(f.size()); forAll(f, fp) @@ -315,9 +322,9 @@ void Foam::functionObjects::fieldValues::surfaceFieldValue::combineMeshGeometry } const pointField& pts = allPoints[proci]; - forAll(pts, i) + for (const point& pt : pts) { - points[nPoints++] = pts[i]; + points[nPoints++] = pt; } } } @@ -343,9 +350,9 @@ void Foam::functionObjects::fieldValues::surfaceFieldValue::combineMeshGeometry } points.transfer(newPoints); - forAll(faces, i) + for (face& f : faces) { - inplaceRenumber(oldToNew, faces[i]); + inplaceRenumber(oldToNew, f); } } } @@ -447,17 +454,19 @@ Foam::functionObjects::fieldValues::surfaceFieldValue::totalArea() const // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * // -bool Foam::functionObjects::fieldValues::surfaceFieldValue::needsSf() const +bool Foam::functionObjects::fieldValues::surfaceFieldValue::usesSf() const { - // Many operations use the Sf field + // Only a few operations do not require the Sf field switch (operation_) { case opNone: + case opMin: + case opMax: case opSum: case opSumMag: case opAverage: - case opMin: - case opMax: + case opSumDirection: + case opSumDirectionBalance: { return false; } @@ -469,26 +478,6 @@ bool Foam::functionObjects::fieldValues::surfaceFieldValue::needsSf() const } -bool Foam::functionObjects::fieldValues::surfaceFieldValue::needsWeight() const -{ - // Only a few operations use weight field - switch (operation_) - { - case opWeightedSum: - case opWeightedAverage: - case opWeightedAreaAverage: - case opWeightedAreaIntegrate: - { - return true; - } - default: - { - return false; - } - } -} - - void Foam::functionObjects::fieldValues::surfaceFieldValue::initialise ( const dictionary& dict @@ -582,19 +571,16 @@ void Foam::functionObjects::fieldValues::surfaceFieldValue::initialise weightFieldName_ = "none"; - if (needsWeight()) + if (usesWeight() && dict.readIfPresent("weightField", weightFieldName_)) { - if (dict.readIfPresent("weightField", weightFieldName_)) + if (regionType_ == stSampledSurface) { - if (regionType_ == stSampledSurface) - { - FatalIOErrorInFunction(dict) - << "Cannot use weightField for sampledSurface" - << exit(FatalIOError); - } - - Info<< " weight field = " << weightFieldName_ << nl; + FatalIOErrorInFunction(dict) + << "Cannot use weightField for sampledSurface" + << exit(FatalIOError); } + + Info<< " weight field = " << weightFieldName_ << nl; } // Backwards compatibility for v1612+ and older @@ -660,10 +646,10 @@ void Foam::functionObjects::fieldValues::surfaceFieldValue::writeFileHeader os << tab << "Area"; } - forAll(fields_, i) + for (const word& fieldName : fields_) { os << tab << operationTypeNames_[operation_] - << "(" << fields_[i] << ")"; + << "(" << fieldName << ")"; } os << endl; @@ -684,12 +670,12 @@ Foam::functionObjects::fieldValues::surfaceFieldValue::processValues { case opSumDirection: { - vector n(dict_.lookup("direction")); + const vector n(dict_.lookup("direction")); return gSum(pos(values*(Sf & n))*mag(values)); } case opSumDirectionBalance: { - vector n(dict_.lookup("direction")); + const vector n(dict_.lookup("direction")); const scalarField nv(values*(Sf & n)); return gSum(pos(nv)*mag(values) - neg(nv)*mag(values)); @@ -754,10 +740,17 @@ Foam::tmp<Foam::scalarField> Foam::functionObjects::fieldValues::surfaceFieldValue::weightingFactor ( const Field<scalar>& weightField -) +) const { - // pass through - return weightField; + if (usesMag()) + { + return mag(weightField); + } + else + { + // pass through + return weightField; + } } @@ -767,16 +760,21 @@ Foam::functionObjects::fieldValues::surfaceFieldValue::weightingFactor ( const Field<scalar>& weightField, const vectorField& Sf -) +) const { // scalar * Area if (returnReduce(weightField.empty(), andOp<bool>())) { + // No weight field - revert to unweighted form return mag(Sf); } + else if (usesMag()) + { + return mag(weightField * mag(Sf)); + } else { - return weightField * mag(Sf); + return (weightField * mag(Sf)); } } @@ -787,16 +785,21 @@ Foam::functionObjects::fieldValues::surfaceFieldValue::weightingFactor ( const Field<vector>& weightField, const vectorField& Sf -) +) const { // vector (dot) Area if (returnReduce(weightField.empty(), andOp<bool>())) { + // No weight field - revert to unweighted form return mag(Sf); } + else if (usesMag()) + { + return mag(weightField & Sf); + } else { - return weightField & Sf; + return (weightField & Sf); } } @@ -915,7 +918,7 @@ bool Foam::functionObjects::fieldValues::surfaceFieldValue::write() // Many operations use the Sf field vectorField Sf; - if (needsSf()) + if (usesSf()) { if (regionType_ == stSurface) { diff --git a/src/functionObjects/field/fieldValues/surfaceFieldValue/surfaceFieldValue.H b/src/functionObjects/field/fieldValues/surfaceFieldValue/surfaceFieldValue.H index 04697bdb342..b34a0529168 100644 --- a/src/functionObjects/field/fieldValues/surfaceFieldValue/surfaceFieldValue.H +++ b/src/functionObjects/field/fieldValues/surfaceFieldValue/surfaceFieldValue.H @@ -103,22 +103,26 @@ Usage The \c operation is one of: \plaintable none | no operation + min | minimum + max | maximum sum | sum - weightedSum | weighted sum sumMag | sum of component magnitudes - sumDirection | sum values which are positive in given direction + sumDirection | sum values that are positive in given direction sumDirectionBalance | sum of balance of values in given direction average | ensemble average - weightedAverage | weighted average areaAverage | area-weighted average - weightedAreaAverage | weighted area average areaIntegrate | area integral - weightedAreaIntegrate | weighted area integral - min | minimum - max | maximum CoV | coefficient of variation: standard deviation/mean areaNormalAverage| area-weighted average in face normal direction areaNormalIntegrate | area-weighted integral in face normal directon + weightedSum | weighted sum + weightedAverage | weighted average + weightedAreaAverage | weighted area average + weightedAreaIntegrate | weighted area integral + absWeightedSum | sum using absolute weighting + absWeightedAverage | average using absolute weighting + absWeightedAreaAverage | area average using absolute weighting + absWeightedAreaIntegrate | area integral using absolute weighting \endplaintable Note @@ -190,36 +194,73 @@ public: //- Region type enumeration enum regionTypes { - stFaceZone, - stPatch, - stSurface, - stSampledSurface + stFaceZone, //!< Calculate on a faceZone + stPatch, //!< Calculate on a patch + stSurface, //!< Calculate with fields on a surfMesh + stSampledSurface //!< Sample onto surface and calculate }; //- Region type names static const Enum<regionTypes> regionTypeNames_; + //- Bitmask values for operation variants + enum operationVariant + { + typeBase = 0, //!< Base operation + typeWeighted = 0x100, //!< Operation using weighting + typeAbsolute = 0x200, //!< Operation using mag (eg, for weighting) + }; //- Operation type enumeration enum operationType { - opNone, //!< None - opSum, //!< Sum - opWeightedSum, //!< Weighted sum - opSumMag, //!< Magnitude of sum + // Normal operations + + opNone = 0, //!< No operation + opMin, //!< Minimum value + opMax, //!< Maximum value + opSum, //!< Sum of values + opSumMag, //!< Sum of component magnitudes opSumDirection, //!< Sum in a given direction - opSumDirectionBalance, //!< Sum in a given direction for multiple - opAverage, //!< Average - opWeightedAverage, //!< Weighted average + + //! Sum of balance of values in given direction + opSumDirectionBalance, + + opAverage, //!< Ensemble average opAreaAverage, //!< Area average - opWeightedAreaAverage, //!< Weighted area average opAreaIntegrate, //!< Area integral - opWeightedAreaIntegrate, //!< Weighted area integral - opMin, //!< Minimum - opMax, //!< Maximum opCoV, //!< Coefficient of variation opAreaNormalAverage, //!< Area average in normal direction - opAreaNormalIntegrate //!< Area integral in normal direction + opAreaNormalIntegrate, //!< Area integral in normal direction + + // Weighted variants + + //! Weighted sum + opWeightedSum = (opSum | typeWeighted), + + //! Weighted average + opWeightedAverage = (opAverage | typeWeighted), + + //! Weighted area average + opWeightedAreaAverage = (opAreaAverage | typeWeighted), + + //! Weighted area integral + opWeightedAreaIntegrate = (opAreaIntegrate | typeWeighted), + + // Variants using absolute weighting + + //! Sum using abs weighting + opAbsWeightedSum = (opWeightedSum | typeAbsolute), + + //! Average using abs weighting + opAbsWeightedAverage = (opWeightedAverage | typeAbsolute), + + //! Area average using abs weighting + opAbsWeightedAreaAverage = (opWeightedAreaAverage | typeAbsolute), + + //! Area integral using abs weighting + opAbsWeightedAreaIntegrate = + (opWeightedAreaIntegrate | typeAbsolute), }; //- Operation type names @@ -229,8 +270,8 @@ public: //- Post-operation type enumeration enum postOperationType { - postOpNone, - postOpSqrt + postOpNone, //!< No additional operation after calculation + postOpSqrt //!< Perform sqrt after normal operation }; //- Operation type names @@ -325,11 +366,19 @@ protected: //- Return the local true/false list representing the face flip map inline const boolList& faceFlip() const; - //- True if the specified operation needs a surface Sf - bool needsSf() const; + //- True if the operation needs a surface Sf + bool usesSf() const; + + //- True if the operation variant uses mag + inline bool usesMag() const; - //- True if the specified operation needs a weight-field - bool needsWeight() const; + //- True if the operation variant uses a weight-field + inline bool usesWeight() const; + + //- True if operation variant uses a weight-field that is available. + // Checks for availability on any processor. + template<class WeightType> + inline bool canWeight(const Field<WeightType>& weightField) const; //- Initialise, e.g. face addressing void initialise(const dictionary& dict); @@ -365,6 +414,7 @@ protected: const Field<WeightType>& weightField ) const; + //- Filter a surface field according to faceIds template<class Type> tmp<Field<Type>> filterField @@ -379,20 +429,24 @@ protected: const GeometricField<Type, fvPatchField, volMesh>& field ) const; - //- Weighting factor + + //- Weighting factor. + // Possibly applies mag() depending on the operation type. template<class WeightType> - static tmp<scalarField> weightingFactor + tmp<scalarField> weightingFactor ( const Field<WeightType>& weightField - ); + ) const; - //- Weighting factor, weight field with the area + //- Weighting factor, weight field with the area. + // Possibly applies mag() depending on the operation type. + // Reverts to mag(Sf) if the weight field is not available. template<class WeightType> - static tmp<scalarField> weightingFactor + tmp<scalarField> weightingFactor ( const Field<WeightType>& weightField, const vectorField& Sf - ); + ) const; //- Templated helper function to output field values @@ -489,7 +543,7 @@ template<> tmp<scalarField> surfaceFieldValue::weightingFactor ( const Field<scalar>& weightField -); +) const; //- Specialisation for scalar - scalar * Area @@ -498,7 +552,7 @@ tmp<scalarField> surfaceFieldValue::weightingFactor ( const Field<scalar>& weightField, const vectorField& Sf -); +) const; //- Specialisation for vector - vector (dot) Area @@ -507,7 +561,7 @@ tmp<scalarField> surfaceFieldValue::weightingFactor ( const Field<vector>& weightField, const vectorField& Sf -); +) const; // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // diff --git a/src/functionObjects/field/fieldValues/surfaceFieldValue/surfaceFieldValueI.H b/src/functionObjects/field/fieldValues/surfaceFieldValue/surfaceFieldValueI.H index e6dacee4dfc..8623229b6e4 100644 --- a/src/functionObjects/field/fieldValues/surfaceFieldValue/surfaceFieldValueI.H +++ b/src/functionObjects/field/fieldValues/surfaceFieldValue/surfaceFieldValueI.H @@ -23,7 +23,6 @@ License \*---------------------------------------------------------------------------*/ -#include "surfaceFieldValue.H" #include "Time.H" // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * // @@ -51,6 +50,22 @@ Foam::functionObjects::fieldValues::surfaceFieldValue::faceFlip() const // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // +inline bool +Foam::functionObjects::fieldValues::surfaceFieldValue::usesMag() const +{ + // Operation specifically tagged to use mag + return (operation_ & typeAbsolute); +} + + +inline bool +Foam::functionObjects::fieldValues::surfaceFieldValue::usesWeight() const +{ + // Operation specifically tagged to require a weight field + return (operation_ & typeWeighted); +} + + inline const Foam::functionObjects::fieldValues::surfaceFieldValue::regionTypes& Foam::functionObjects::fieldValues::surfaceFieldValue::regionType() const { diff --git a/src/functionObjects/field/fieldValues/surfaceFieldValue/surfaceFieldValueTemplates.C b/src/functionObjects/field/fieldValues/surfaceFieldValue/surfaceFieldValueTemplates.C index 806a0684546..58520474880 100644 --- a/src/functionObjects/field/fieldValues/surfaceFieldValue/surfaceFieldValueTemplates.C +++ b/src/functionObjects/field/fieldValues/surfaceFieldValue/surfaceFieldValueTemplates.C @@ -33,6 +33,20 @@ License // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * // +template<class WeightType> +inline bool Foam::functionObjects::fieldValues::surfaceFieldValue::canWeight +( + const Field<WeightType>& weightField +) const +{ + return + ( + usesWeight() + && returnReduce(!weightField.empty(), orOp<bool>()) // On some processor + ); +} + + template<class Type> bool Foam::functionObjects::fieldValues::surfaceFieldValue::validField ( @@ -142,28 +156,36 @@ processSameTypeValues { break; } - case opSum: + case opMin: + { + result = gMin(values); + break; + } + case opMax: + { + result = gMax(values); + break; + } + case opSumMag: { - result = gSum(values); + result = gSum(cmptMag(values)); break; } + case opSum: case opWeightedSum: + case opAbsWeightedSum: { - if (returnReduce(weightField.empty(), andOp<bool>())) - { - result = gSum(values); - } - else + if (canWeight(weightField)) { tmp<scalarField> weight(weightingFactor(weightField)); result = gSum(weight*values); } - break; - } - case opSumMag: - { - result = gSum(cmptMag(values)); + else + { + // Unweighted form + result = gSum(values); + } break; } case opSumDirection: @@ -178,62 +200,56 @@ processSameTypeValues break; } case opAverage: - { - const label n = returnReduce(values.size(), sumOp<label>()); - result = gSum(values)/(scalar(n) + ROOTVSMALL); - break; - } case opWeightedAverage: + case opAbsWeightedAverage: { - if (returnReduce(weightField.empty(), andOp<bool>())) - { - const label n = returnReduce(values.size(), sumOp<label>()); - result = gSum(values)/(scalar(n) + ROOTVSMALL); - } - else + if (canWeight(weightField)) { const scalarField factor(weightingFactor(weightField)); result = gSum(factor*values)/(gSum(factor) + ROOTVSMALL); } + else + { + // Unweighted form + const label n = returnReduce(values.size(), sumOp<label>()); + result = gSum(values)/(scalar(n) + ROOTVSMALL); + } break; } case opAreaAverage: - { - const scalarField factor(mag(Sf)); - - result = gSum(factor*values)/gSum(factor); - break; - } case opWeightedAreaAverage: + case opAbsWeightedAreaAverage: { - const scalarField factor(weightingFactor(weightField, Sf)); + if (canWeight(weightField)) + { + const scalarField factor(weightingFactor(weightField, Sf)); - result = gSum(factor*values)/gSum(factor + ROOTVSMALL); + result = gSum(factor*values)/gSum(factor + ROOTVSMALL); + } + else + { + // Unweighted form + const scalarField factor(mag(Sf)); + result = gSum(factor*values)/gSum(factor); + } break; } case opAreaIntegrate: - { - const scalarField factor(mag(Sf)); - - result = gSum(factor*values); - break; - } case opWeightedAreaIntegrate: + case opAbsWeightedAreaIntegrate: { - const scalarField factor(weightingFactor(weightField, Sf)); - - result = gSum(factor*values); - break; - } - case opMin: - { - result = gMin(values); - break; - } - case opMax: - { - result = gMax(values); + if (canWeight(weightField)) + { + tmp<scalarField> factor(weightingFactor(weightField, Sf)); + result = gSum(factor*values); + } + else + { + // Unweighted form + tmp<scalarField> factor(mag(Sf)); + result = gSum(factor*values); + } break; } case opCoV: @@ -245,8 +261,8 @@ processSameTypeValues for (direction d=0; d < pTraits<Type>::nComponents; ++d) { - tmp<scalarField> vals = values.component(d); - scalar mean = component(meanValue, d); + tmp<scalarField> vals(values.component(d)); + const scalar mean = component(meanValue, d); scalar& res = setComponent(result, d); res = @@ -286,8 +302,10 @@ Foam::tmp<Foam::scalarField> Foam::functionObjects::fieldValues::surfaceFieldValue::weightingFactor ( const Field<WeightType>& weightField -) +) const { + // The scalar form is specialized. + // For other types always need mag() to generate a scalar field. return mag(weightField); } @@ -302,10 +320,8 @@ Foam::label Foam::functionObjects::fieldValues::surfaceFieldValue::writeAll { label nProcessed = 0; - forAll(fields_, i) + for (const word& fieldName : fields_) { - const word& fieldName = fields_[i]; - if ( writeValues<scalar>(fieldName, Sf, weightField, surfToWrite) @@ -437,8 +453,8 @@ Foam::functionObjects::fieldValues::surfaceFieldValue::filterField forAll(values, i) { - label facei = faceId_[i]; - label patchi = facePatchId_[i]; + const label facei = faceId_[i]; + const label patchi = facePatchId_[i]; if (patchi >= 0) { values[i] = field.boundaryField()[patchi][facei]; @@ -472,8 +488,8 @@ Foam::functionObjects::fieldValues::surfaceFieldValue::filterField forAll(values, i) { - label facei = faceId_[i]; - label patchi = facePatchId_[i]; + const label facei = faceId_[i]; + const label patchi = facePatchId_[i]; if (patchi >= 0) { values[i] = field.boundaryField()[patchi][facei]; diff --git a/src/functionObjects/field/fieldValues/volFieldValue/volFieldValue.C b/src/functionObjects/field/fieldValues/volFieldValue/volFieldValue.C index 912b6459bbf..eafef810961 100644 --- a/src/functionObjects/field/fieldValues/volFieldValue/volFieldValue.C +++ b/src/functionObjects/field/fieldValues/volFieldValue/volFieldValue.C @@ -49,32 +49,35 @@ const Foam::Enum > Foam::functionObjects::fieldValues::volFieldValue::operationTypeNames_ { + // Normal operations { operationType::opNone, "none" }, + { operationType::opMin, "min" }, + { operationType::opMax, "max" }, { operationType::opSum, "sum" }, - { operationType::opWeightedSum, "weightedSum" }, { operationType::opSumMag, "sumMag" }, { operationType::opAverage, "average" }, - { operationType::opWeightedAverage, "weightedAverage" }, { operationType::opVolAverage, "volAverage" }, - { operationType::opWeightedVolAverage, "weightedVolAverage" }, { operationType::opVolIntegrate, "volIntegrate" }, - { operationType::opWeightedVolIntegrate, "weightedVolIntegrate" }, - { operationType::opMin, "min" }, - { operationType::opMax, "max" }, { operationType::opCoV, "CoV" }, + + // Using weighting + { operationType::opWeightedSum, "weightedSum" }, + { operationType::opWeightedAverage, "weightedAverage" }, + { operationType::opWeightedVolAverage, "weightedVolAverage" }, + { operationType::opWeightedVolIntegrate, "weightedVolIntegrate" }, }; // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * // -bool Foam::functionObjects::fieldValues::volFieldValue::needsVol() const +bool Foam::functionObjects::fieldValues::volFieldValue::usesVol() const { - // Only some operations need the cell volume + // Only a few operations require the cell volume switch (operation_) { case opVolAverage: - case opWeightedVolAverage: case opVolIntegrate: + case opWeightedVolAverage: case opWeightedVolIntegrate: case opCoV: return true; @@ -85,20 +88,23 @@ bool Foam::functionObjects::fieldValues::volFieldValue::needsVol() const } -bool Foam::functionObjects::fieldValues::volFieldValue::needsWeight() const +bool Foam::functionObjects::fieldValues::volFieldValue::usesWeight() const { - // Only a few operations use weight field - switch (operation_) - { - case opWeightedSum: - case opWeightedAverage: - case opWeightedVolAverage: - case opWeightedVolIntegrate: - return true; + // Operation specifically tagged to require a weight field + return (operation_ & typeWeighted); +} - default: - return false; - } + +bool Foam::functionObjects::fieldValues::volFieldValue::canWeight +( + const scalarField& weightField +) const +{ + return + ( + usesWeight() + && returnReduce(!weightField.empty(), orOp<bool>()) // On some processor + ); } @@ -108,12 +114,9 @@ void Foam::functionObjects::fieldValues::volFieldValue::initialise ) { weightFieldName_ = "none"; - if (needsWeight()) + if (usesWeight() && dict.readIfPresent("weightField", weightFieldName_)) { - if (dict.readIfPresent("weightField", weightFieldName_)) - { - Info<< " weight field = " << weightFieldName_; - } + Info<< " weight field = " << weightFieldName_; } Info<< nl << endl; @@ -133,10 +136,10 @@ void Foam::functionObjects::fieldValues::volFieldValue::writeFileHeader writeCommented(os, "Time"); - forAll(fields_, fieldi) + for (const word& fieldName : fields_) { os << tab << operationTypeNames_[operation_] - << "(" << fields_[fieldi] << ")"; + << "(" << fieldName << ")"; } os << endl; @@ -151,10 +154,8 @@ Foam::label Foam::functionObjects::fieldValues::volFieldValue::writeAll { label nProcessed = 0; - forAll(fields_, i) + for (const word& fieldName : fields_) { - const word& fieldName = fields_[i]; - if ( writeValues<scalar>(fieldName, V, weightField) @@ -245,7 +246,7 @@ bool Foam::functionObjects::fieldValues::volFieldValue::write() // Only some operations need the cell volume scalarField V; - if (needsVol()) + if (usesVol()) { V = filterField(fieldValue::mesh_.V()); } diff --git a/src/functionObjects/field/fieldValues/volFieldValue/volFieldValue.H b/src/functionObjects/field/fieldValues/volFieldValue/volFieldValue.H index f2730c526f8..ae299430ae2 100644 --- a/src/functionObjects/field/fieldValues/volFieldValue/volFieldValue.H +++ b/src/functionObjects/field/fieldValues/volFieldValue/volFieldValue.H @@ -81,18 +81,18 @@ Usage The \c operation is one of: \plaintable none | No operation + min | Minimum + max | Maximum sum | Sum - weightedSum | Weighted sum sumMag | Sum of component magnitudes average | Ensemble average - weightedAverage | Weighted average volAverage | Volume weighted average - weightedVolAverage | Weighted volume average volIntegrate | Volume integral - weightedVolIntegrate | Weighted volume integral - min | Minimum - max | Maximum CoV | Coefficient of variation: standard deviation/mean + weightedSum | Weighted sum + weightedAverage | Weighted average + weightedVolAverage | Weighted volume average + weightedVolIntegrate | Weighted volume integral \endplaintable See also @@ -134,22 +134,41 @@ public: // Public data types + //- Bitmask values for operation variants + enum operationVariant + { + typeBase = 0, //!< Base operation + typeWeighted = 0x100, //!< Operation using weighting + }; + //- Operation type enumeration enum operationType { - opNone, //!< None + // Normal operations + + opNone = 0, //!< No operation + opMin, //!< Minimum + opMax, //!< Maximum opSum, //!< Sum - opWeightedSum, //!< Weighted sum opSumMag, //!< Magnitude of sum opAverage, //!< Average - opWeightedAverage, //!< Weighted average opVolAverage, //!< Volume average - opWeightedVolAverage, //!< Weighted volume average opVolIntegrate, //!< Volume integral - opWeightedVolIntegrate, //!< Weighted volume integral - opMin, //!< Minimum - opMax, //!< Maximum - opCoV //!< Coefficient of variation + opCoV, //!< Coefficient of variation + + // Weighted variants + + //! Weighted sum + opWeightedSum = (opSum | typeWeighted), + + //! Weighted average + opWeightedAverage = (opAverage | typeWeighted), + + //! Weighted volume average + opWeightedVolAverage = (opVolAverage | typeWeighted), + + //! Weighted volume integral + opWeightedVolIntegrate = (opVolIntegrate | typeWeighted), }; //- Operation type names @@ -169,11 +188,15 @@ protected: // Protected Member Functions - //- True if the specified operation needs the cell-volume - bool needsVol() const; + //- True if the operation needs the cell-volume + bool usesVol() const; + + //- True if the operation variant uses a weight-field + bool usesWeight() const; - //- True if the specified operation needs a weight-field - bool needsWeight() const; + //- True if operation variant uses a weight-field that is available. + // Checks for availability on any processor. + inline bool canWeight(const scalarField& weightField) const; //- Initialise, e.g. cell addressing void initialise(const dictionary& dict); diff --git a/src/functionObjects/field/fieldValues/volFieldValue/volFieldValueTemplates.C b/src/functionObjects/field/fieldValues/volFieldValue/volFieldValueTemplates.C index a4f1fb82a6c..09b9931808f 100644 --- a/src/functionObjects/field/fieldValues/volFieldValue/volFieldValueTemplates.C +++ b/src/functionObjects/field/fieldValues/volFieldValue/volFieldValueTemplates.C @@ -82,21 +82,18 @@ Type Foam::functionObjects::fieldValues::volFieldValue::processValues Type result = Zero; switch (operation_) { - case opSum: + case opNone: { - result = gSum(values); break; } - case opWeightedSum: + case opMin: { - if (returnReduce(weightField.empty(), andOp<bool>())) - { - result = gSum(values); - } - else - { - result = gSum(weightField*values); - } + result = gMin(values); + break; + } + case opMax: + { + result = gMax(values); break; } case opSumMag: @@ -104,71 +101,65 @@ Type Foam::functionObjects::fieldValues::volFieldValue::processValues result = gSum(cmptMag(values)); break; } - case opAverage: + case opSum: + case opWeightedSum: { - const label n = returnReduce(values.size(), sumOp<label>()); - result = gSum(values)/(scalar(n) + ROOTVSMALL); + if (canWeight(weightField)) + { + result = gSum(weightField*values); + } + else + { + // Unweighted form + result = gSum(values); + } break; } + case opAverage: case opWeightedAverage: { - if (returnReduce(weightField.empty(), andOp<bool>())) + if (canWeight(weightField)) { - const label n = returnReduce(values.size(), sumOp<label>()); - result = gSum(values)/(scalar(n) + ROOTVSMALL); + result = + gSum(weightField*values)/(gSum(weightField) + ROOTVSMALL); } else { - result = - gSum(weightField*values)/(gSum(weightField) + ROOTVSMALL); + // Unweighted form + const label n = returnReduce(values.size(), sumOp<label>()); + result = gSum(values)/(scalar(n) + ROOTVSMALL); } break; } case opVolAverage: - { - result = gSum(values*V)/gSum(V); - break; - } case opWeightedVolAverage: { - if (returnReduce(weightField.empty(), andOp<bool>())) + if (canWeight(weightField)) { - result = gSum(V*values)/(gSum(V) + ROOTVSMALL); + result = gSum(weightField*V*values) + /(gSum(weightField*V) + ROOTVSMALL); } else { - result = gSum(weightField*V*values) - /(gSum(weightField*V) + ROOTVSMALL); + // Unweighted form + result = gSum(V*values)/(gSum(V) + ROOTVSMALL); } break; } case opVolIntegrate: - { - result = gSum(V*values); - break; - } case opWeightedVolIntegrate: { - if (returnReduce(weightField.empty(), andOp<bool>())) + if (canWeight(weightField)) { - result = gSum(V*values); + result = gSum(weightField*V*values); } else { - result = gSum(weightField*V*values); + // Unweighted form + result = gSum(V*values); } break; } - case opMin: - { - result = gMin(values); - break; - } - case opMax: - { - result = gMax(values); - break; - } case opCoV: { const scalar sumV = gSum(V); @@ -177,8 +168,8 @@ Type Foam::functionObjects::fieldValues::volFieldValue::processValues for (direction d=0; d < pTraits<Type>::nComponents; ++d) { - scalarField vals(values.component(d)); - scalar mean = component(meanValue, d); + tmp<scalarField> vals(values.component(d)); + const scalar mean = component(meanValue, d); scalar& res = setComponent(result, d); res = sqrt(gSum(V*sqr(vals - mean))/sumV)/(mean + ROOTVSMALL); @@ -186,8 +177,6 @@ Type Foam::functionObjects::fieldValues::volFieldValue::processValues break; } - case opNone: - {} } return result; -- GitLab