From 1812773b2518942bcd6e7178d0121b103f90fba8 Mon Sep 17 00:00:00 2001 From: andy Date: Fri, 8 Jan 2010 17:14:53 +0000 Subject: [PATCH] Updates to the fieldValues function object - Updates to enable correct operation in parallel - Added weighted average operation for cell sources --- .../field/fieldValues/cellSource/cellSource.C | 52 ++++-- .../field/fieldValues/cellSource/cellSource.H | 35 +++- .../cellSource/cellSourceTemplates.C | 69 ++++---- .../field/fieldValues/faceSource/faceSource.C | 19 ++- .../field/fieldValues/faceSource/faceSource.H | 19 ++- .../faceSource/faceSourceTemplates.C | 150 +++++------------- .../field/fieldValues/fieldValue/fieldValue.H | 13 ++ .../fieldValue/fieldValueTemplates.C | 66 ++++++++ 8 files changed, 249 insertions(+), 174 deletions(-) create mode 100644 src/postProcessing/functionObjects/field/fieldValues/fieldValue/fieldValueTemplates.C diff --git a/src/postProcessing/functionObjects/field/fieldValues/cellSource/cellSource.C b/src/postProcessing/functionObjects/field/fieldValues/cellSource/cellSource.C index b94b687ba1..239aa4342d 100644 --- a/src/postProcessing/functionObjects/field/fieldValues/cellSource/cellSource.C +++ b/src/postProcessing/functionObjects/field/fieldValues/cellSource/cellSource.C @@ -46,10 +46,13 @@ namespace Foam fieldValues::cellSource::sourceTypeNames_; template<> - const char* NamedEnum:: - names[] = {"none", "sum", "volAverage", "volIntegrate"}; + const char* NamedEnum:: + names[] = + { + "none", "sum", "volAverage", "volIntegrate", "weightedAverage" + }; - const NamedEnum + const NamedEnum fieldValues::cellSource::operationTypeNames_; } @@ -93,7 +96,7 @@ void Foam::fieldValues::cellSource::setCellZoneCells() // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * // -void Foam::fieldValues::cellSource::initialise() +void Foam::fieldValues::cellSource::initialise(const dictionary& dict) { switch (source_) { @@ -104,7 +107,7 @@ void Foam::fieldValues::cellSource::initialise() } default: { - FatalErrorIn("cellSource::constructCellAddressing()") + FatalErrorIn("cellSource::initialise()") << "Unknown source type. Valid source types are:" << sourceTypeNames_ << nl << exit(FatalError); } @@ -114,6 +117,29 @@ void Foam::fieldValues::cellSource::initialise() << " total cells = " << cellId_.size() << nl << " total volume = " << sum(filterField(mesh().V())) << nl << endl; + + if (operation_ == opWeightedAverage) + { + dict.lookup("weightField") >> weightFieldName_; + if + ( + obr().foundObject(weightFieldName_) + ) + { + Info<< " weight field = " << weightFieldName_; + } + else + { + FatalErrorIn("cellSource::initialise()") + << type() << " " << name_ << ": " + << sourceTypeNames_[source_] << "(" << sourceName_ << "):" + << nl << " Weight field " << weightFieldName_ + << " must be a " << volScalarField::typeName + << nl << exit(FatalError); + } + } + + Info<< nl << endl; } @@ -172,7 +198,7 @@ void Foam::fieldValues::cellSource::read(const dictionary& dict) if (active_) { // no additional info to read - initialise(); + initialise(dict); } } @@ -183,9 +209,12 @@ void Foam::fieldValues::cellSource::write() if (active_) { - outputFilePtr_() - << obr_.time().value() << tab - << sum(filterField(mesh().V())); + if (Pstream::master()) + { + outputFilePtr_() + << obr_.time().value() << tab + << sum(filterField(mesh().V())); + } forAll(fields_, i) { @@ -196,7 +225,10 @@ void Foam::fieldValues::cellSource::write() writeValues(fields_[i]); } - outputFilePtr_()<< endl; + if (Pstream::master()) + { + outputFilePtr_()<< endl; + } if (log_) { diff --git a/src/postProcessing/functionObjects/field/fieldValues/cellSource/cellSource.H b/src/postProcessing/functionObjects/field/fieldValues/cellSource/cellSource.H index 89cc9ad010..a8d04f7d9a 100644 --- a/src/postProcessing/functionObjects/field/fieldValues/cellSource/cellSource.H +++ b/src/postProcessing/functionObjects/field/fieldValues/cellSource/cellSource.H @@ -39,7 +39,7 @@ Description valueOutput true; // Write values at run-time output times? source cellZone; // Type of cell source sourceName c0; - operation volAverage; // none, sum, volAverage, volIntegrate + operation volAverage; fields ( p @@ -47,6 +47,13 @@ Description ); } + where operation is one of: + - none + - sum + - volAverage + - volIntegrate + - weightedAverage + SourceFiles cellSource.C @@ -96,11 +103,12 @@ public: opNone, opSum, opVolAverage, - opVolIntegrate + opVolIntegrate, + opWeightedAverage }; //- Operation type names - static const NamedEnum operationTypeNames_; + static const NamedEnum operationTypeNames_; private: @@ -127,23 +135,34 @@ protected: //- Local list of cell IDs labelList cellId_; + //- Weight field name - only used for opWeightedAverage mode + word weightFieldName_; + // Protected member functions //- Initialise, e.g. cell addressing - void initialise(); + void initialise(const dictionary& dict); + + //- Return true if the field name is valid + template + bool validField(const word& fieldName) const; //- Insert field values into values list template - bool setFieldValues + tmp > setFieldValues ( - const word& fieldName, - List& values + const word& fieldName ) const; //- Apply the 'operation' to the values template - Type processValues(const List& values) const; + Type processValues + ( + const Field& values, + const scalarField& V, + const scalarField& weightField + ) const; //- Output file header information virtual void writeFileHeader(); diff --git a/src/postProcessing/functionObjects/field/fieldValues/cellSource/cellSourceTemplates.C b/src/postProcessing/functionObjects/field/fieldValues/cellSource/cellSourceTemplates.C index fe495fcca6..41c5e70d26 100644 --- a/src/postProcessing/functionObjects/field/fieldValues/cellSource/cellSourceTemplates.C +++ b/src/postProcessing/functionObjects/field/fieldValues/cellSource/cellSourceTemplates.C @@ -27,27 +27,16 @@ License #include "cellSource.H" #include "volFields.H" #include "IOList.H" -#include "ListListOps.H" // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * // template -bool Foam::fieldValues::cellSource::setFieldValues -( - const word& fieldName, - List& values -) const +bool Foam::fieldValues::cellSource::validField(const word& fieldName) const { - values.setSize(cellId_.size(), pTraits::zero); - typedef GeometricField vf; if (obr_.foundObject(fieldName)) { - const vf& field = obr_.lookupObject(fieldName); - - values = UIndirectList(field, cellId_); - return true; } @@ -55,10 +44,29 @@ bool Foam::fieldValues::cellSource::setFieldValues } +template +Foam::tmp > Foam::fieldValues::cellSource::setFieldValues +( + const word& fieldName +) const +{ + typedef GeometricField vf; + + if (obr_.foundObject(fieldName)) + { + return filterField(obr_.lookupObject(fieldName)); + } + + return tmp >(new Field(0.0)); +} + + template Type Foam::fieldValues::cellSource::processValues ( - const List& values + const Field& values, + const scalarField& V, + const scalarField& weightField ) const { Type result = pTraits::zero; @@ -71,13 +79,17 @@ Type Foam::fieldValues::cellSource::processValues } case opVolAverage: { - tmp V = filterField(mesh().V()); - result = sum(values*V())/sum(V()); + result = sum(values*V)/sum(V); break; } case opVolIntegrate: { - result = sum(values*filterField(mesh().V())); + result = sum(values*V); + break; + } + case opWeightedAverage: + { + result = sum(values*weightField)/sum(weightField); break; } default: @@ -95,25 +107,20 @@ Type Foam::fieldValues::cellSource::processValues template bool Foam::fieldValues::cellSource::writeValues(const word& fieldName) { - List > allValues(Pstream::nProcs()); + const bool ok = validField(fieldName); - bool validField = - setFieldValues(fieldName, allValues[Pstream::myProcNo()]); - - if (validField) + if (ok) { - Pstream::gatherList(allValues); + Field values = combineFields(setFieldValues(fieldName)); + + scalarField V = combineFields(filterField(mesh().V())); + + scalarField weightField = + combineFields(setFieldValues(weightFieldName_)); if (Pstream::master()) { - List values = - ListListOps::combine > - ( - allValues, - accessOp >() - ); - - Type result = processValues(values); + Type result = processValues(values, V, weightField); if (valueOutput_) { @@ -144,7 +151,7 @@ bool Foam::fieldValues::cellSource::writeValues(const word& fieldName) } } - return validField; + return ok; } diff --git a/src/postProcessing/functionObjects/field/fieldValues/faceSource/faceSource.C b/src/postProcessing/functionObjects/field/fieldValues/faceSource/faceSource.C index b609f4310a..e84fe07402 100644 --- a/src/postProcessing/functionObjects/field/fieldValues/faceSource/faceSource.C +++ b/src/postProcessing/functionObjects/field/fieldValues/faceSource/faceSource.C @@ -220,7 +220,7 @@ void Foam::fieldValues::faceSource::initialise(const dictionary& dict) } default: { - FatalErrorIn("faceSource::constructFaceAddressing()") + FatalErrorIn("faceSource::initiliase()") << type() << " " << name_ << ": " << sourceTypeNames_[source_] << "(" << sourceName_ << "):" << nl << " Unknown source type. Valid source types are:" @@ -245,7 +245,7 @@ void Foam::fieldValues::faceSource::initialise(const dictionary& dict) } else { - FatalErrorIn("faceSource::constructFaceAddressing()") + FatalErrorIn("faceSource::initialise()") << type() << " " << name_ << ": " << sourceTypeNames_[source_] << "(" << sourceName_ << "):" << nl << " Weight field " << weightFieldName_ @@ -326,9 +326,12 @@ void Foam::fieldValues::faceSource::write() if (active_) { - outputFilePtr_() - << obr_.time().value() << tab - << sum(filterField(mesh().magSf())); + if (Pstream::master()) + { + outputFilePtr_() + << obr_.time().value() << tab + << sum(filterField(mesh().magSf())); + } forAll(fields_, i) { @@ -339,7 +342,10 @@ void Foam::fieldValues::faceSource::write() writeValues(fields_[i]); } - outputFilePtr_()<< endl; + if (Pstream::master()) + { + outputFilePtr_()<< endl; + } if (log_) { @@ -350,4 +356,3 @@ void Foam::fieldValues::faceSource::write() // ************************************************************************* // - diff --git a/src/postProcessing/functionObjects/field/fieldValues/faceSource/faceSource.H b/src/postProcessing/functionObjects/field/fieldValues/faceSource/faceSource.H index e8e2634b12..8ef65d870b 100644 --- a/src/postProcessing/functionObjects/field/fieldValues/faceSource/faceSource.H +++ b/src/postProcessing/functionObjects/field/fieldValues/faceSource/faceSource.H @@ -153,17 +153,22 @@ protected: //- Initialise, e.g. face addressing void initialise(const dictionary& dict); - //- Insert field values into values list + //- Return true if the field name is valid template - bool setFieldValues - ( - const word& fieldName, - List& values - ) const; + bool validField(const word& fieldName) const; + + //- Return field values by looking up field name + template + tmp > setFieldValues(const word& fieldName) const; //- Apply the 'operation' to the values template - Type processValues(const List& values) const; + Type processValues + ( + const Field& values, + const scalarField& magSf, + const scalarField& weightField + ) const; //- Output file header information virtual void writeFileHeader(); diff --git a/src/postProcessing/functionObjects/field/fieldValues/faceSource/faceSourceTemplates.C b/src/postProcessing/functionObjects/field/fieldValues/faceSource/faceSourceTemplates.C index d7609573a7..a257384a50 100644 --- a/src/postProcessing/functionObjects/field/fieldValues/faceSource/faceSourceTemplates.C +++ b/src/postProcessing/functionObjects/field/fieldValues/faceSource/faceSourceTemplates.C @@ -33,70 +33,17 @@ License // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * // template -bool Foam::fieldValues::faceSource::setFieldValues -( - const word& fieldName, - List& values -) const +bool Foam::fieldValues::faceSource::validField(const word& fieldName) const { - values.setSize(faceId_.size(), pTraits::zero); - typedef GeometricField sf; typedef GeometricField vf; if (obr_.foundObject(fieldName)) { - const sf& field = obr_.lookupObject(fieldName); - - forAll(values, i) - { - label faceI = faceId_[i]; - label patchI = facePatchId_[i]; - if (patchI >= 0) - { - values[i] = field.boundaryField()[patchI][faceI]; - } - else - { - values[i] = field[faceI]; - } - - values[i] *= flipMap_[i]; - } - return true; } else if (obr_.foundObject(fieldName)) { - const vf& field = obr_.lookupObject(fieldName); - - forAll(values, i) - { - label faceI = faceId_[i]; - label patchI = facePatchId_[i]; - if (patchI >= 0) - { - values[i] = field.boundaryField()[patchI][faceI]; - } - else - { - FatalErrorIn - ( - "fieldValues::faceSource::setFieldValues" - "(" - "const word&, " - "List&" - ") const" - ) << type() << " " << name_ << ": " - << sourceTypeNames_[source_] << "(" << sourceName_ << "):" - << nl - << " Unable to process internal faces for volume field " - << fieldName << nl << abort(FatalError); - } - - values[i] *= flipMap_[i]; - } - return true; } @@ -104,10 +51,34 @@ bool Foam::fieldValues::faceSource::setFieldValues } +template +Foam::tmp > Foam::fieldValues::faceSource::setFieldValues +( + const word& fieldName +) const +{ + typedef GeometricField sf; + typedef GeometricField vf; + + if (obr_.foundObject(fieldName)) + { + return filterField(obr_.lookupObject(fieldName)); + } + else if (obr_.foundObject(fieldName)) + { + return filterField(obr_.lookupObject(fieldName)); + } + + return tmp >(new Field(0.0)); +} + + template Type Foam::fieldValues::faceSource::processValues ( - const List& values + const Field& values, + const scalarField& magSf, + const scalarField& weightField ) const { Type result = pTraits::zero; @@ -120,54 +91,17 @@ Type Foam::fieldValues::faceSource::processValues } case opAreaAverage: { - tmp magSf = filterField(mesh().magSf()); - result = sum(values*magSf())/sum(magSf()); + result = sum(values*magSf)/sum(magSf); break; } case opAreaIntegrate: { - result = sum(values*filterField(mesh().magSf())); + result = sum(values*magSf); break; } case opWeightedAverage: { - if (mesh().foundObject(weightFieldName_)) - { - tmp wField = - filterField - ( - mesh().lookupObject(weightFieldName_) - ); - result = sum(values*wField())/sum(wField()); - } - else if (mesh().foundObject(weightFieldName_)) - { - tmp wField = - filterField - ( - mesh().lookupObject - ( - weightFieldName_ - ) - ); - result = sum(values*wField())/sum(wField()); - } - else - { - FatalErrorIn - ( - "fieldValues::faceSource::processValues" - "(" - "List&" - ") const" - ) << type() << " " << name_ << ": " - << sourceTypeNames_[source_] << "(" << sourceName_ << "):" - << nl - << " Weight field " << weightFieldName_ - << " must be either a " << volScalarField::typeName - << " or " << surfaceScalarField::typeName << nl - << abort(FatalError); - } + result = sum(values*weightField)/sum(weightField); break; } default: @@ -185,25 +119,20 @@ Type Foam::fieldValues::faceSource::processValues template bool Foam::fieldValues::faceSource::writeValues(const word& fieldName) { - List > allValues(Pstream::nProcs()); + const bool ok = validField(fieldName); - bool validField = - setFieldValues(fieldName, allValues[Pstream::myProcNo()]); - - if (validField) + if (ok) { - Pstream::gatherList(allValues); + Field values = combineFields(setFieldValues(fieldName)); + + scalarField magSf = combineFields(filterField(mesh().magSf())); + + scalarField weightField = + combineFields(setFieldValues(weightFieldName_)); if (Pstream::master()) { - List values = - ListListOps::combine > - ( - allValues, - accessOp >() - ); - - Type result = processValues(values); + Type result = processValues(values, magSf, weightField); if (valueOutput_) { @@ -222,7 +151,6 @@ bool Foam::fieldValues::faceSource::writeValues(const word& fieldName) ).write(); } - outputFilePtr_()<< tab << result; if (log_) @@ -234,7 +162,7 @@ bool Foam::fieldValues::faceSource::writeValues(const word& fieldName) } } - return validField; + return ok; } diff --git a/src/postProcessing/functionObjects/field/fieldValues/fieldValue/fieldValue.H b/src/postProcessing/functionObjects/field/fieldValues/fieldValue/fieldValue.H index 1b96bb8e41..4aa6f41bd7 100644 --- a/src/postProcessing/functionObjects/field/fieldValues/fieldValue/fieldValue.H +++ b/src/postProcessing/functionObjects/field/fieldValues/fieldValue/fieldValue.H @@ -164,6 +164,13 @@ public: //- Execute the at the final time-loop, currently does nothing virtual void end(); + + //- Comnbine fields from all processor domains into single field + template + tmp > combineFields + ( + const tmp >& field + ) const; }; @@ -177,6 +184,12 @@ public: // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // +#ifdef NoRepository + #include "fieldValueTemplates.C" +#endif + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + #endif // ************************************************************************* // diff --git a/src/postProcessing/functionObjects/field/fieldValues/fieldValue/fieldValueTemplates.C b/src/postProcessing/functionObjects/field/fieldValues/fieldValue/fieldValueTemplates.C new file mode 100644 index 0000000000..23c700134f --- /dev/null +++ b/src/postProcessing/functionObjects/field/fieldValues/fieldValue/fieldValueTemplates.C @@ -0,0 +1,66 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2009-2009 OpenCFD Ltd. + \\/ 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 2 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, write to the Free Software Foundation, + Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA + +\*---------------------------------------------------------------------------*/ + +#include "fieldValue.H" +#include "ListListOps.H" +#include "Pstream.H" + +// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // + +template +Foam::tmp > Foam::fieldValue::combineFields +( + const tmp >& field +) const +{ + List > allValues(Pstream::nProcs()); + + allValues[Pstream::myProcNo()] = field(); + + Pstream::gatherList(allValues); + + if (Pstream::master()) + { + return tmp > + ( + new Field + ( + ListListOps::combine > + ( + allValues, + accessOp >() + ) + ) + ); + } + else + { + return field(); + } +} + + +// ************************************************************************* //