From 3c44f08d5260bdb0a94a2e0cd3fd0eec30e056e0 Mon Sep 17 00:00:00 2001 From: Henry Weller Date: Sun, 22 May 2016 16:30:48 +0100 Subject: [PATCH] functionObjects::blendingFactor: simplified, standardized, rationalized --- .../functionObjects/field/Make/files | 1 + .../blendingFactor/blendingFactor.C | 35 ++++------ .../blendingFactor/blendingFactor.H | 31 ++------- .../blendingFactor/blendingFactorTemplates.C | 65 ++++--------------- .../functionObjects/utilities/Make/files | 1 - 5 files changed, 33 insertions(+), 100 deletions(-) rename src/postProcessing/functionObjects/{utilities => field}/blendingFactor/blendingFactor.C (79%) rename src/postProcessing/functionObjects/{utilities => field}/blendingFactor/blendingFactor.H (81%) rename src/postProcessing/functionObjects/{utilities => field}/blendingFactor/blendingFactorTemplates.C (57%) diff --git a/src/postProcessing/functionObjects/field/Make/files b/src/postProcessing/functionObjects/field/Make/files index 1ba35cee02..8b1788a4e6 100644 --- a/src/postProcessing/functionObjects/field/Make/files +++ b/src/postProcessing/functionObjects/field/Make/files @@ -41,5 +41,6 @@ Q/Q.C Lambda2/Lambda2.C CourantNo/CourantNo.C PecletNo/PecletNo.C +blendingFactor/blendingFactor.C LIB = $(FOAM_LIBBIN)/libfieldFunctionObjects diff --git a/src/postProcessing/functionObjects/utilities/blendingFactor/blendingFactor.C b/src/postProcessing/functionObjects/field/blendingFactor/blendingFactor.C similarity index 79% rename from src/postProcessing/functionObjects/utilities/blendingFactor/blendingFactor.C rename to src/postProcessing/functionObjects/field/blendingFactor/blendingFactor.C index a063ce82b0..50a055c0a6 100644 --- a/src/postProcessing/functionObjects/utilities/blendingFactor/blendingFactor.C +++ b/src/postProcessing/functionObjects/field/blendingFactor/blendingFactor.C @@ -24,7 +24,6 @@ License \*---------------------------------------------------------------------------*/ #include "blendingFactor.H" -#include "volFields.H" #include "addToRunTimeSelectionTable.H" // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // @@ -48,7 +47,7 @@ Foam::functionObjects::blendingFactor::blendingFactor const dictionary& dict ) : - fvMeshFunctionObject(name, runTime, dict) + fieldExpression(name, runTime, dict) { read(dict); } @@ -64,8 +63,11 @@ Foam::functionObjects::blendingFactor::~blendingFactor() bool Foam::functionObjects::blendingFactor::read(const dictionary& dict) { + fieldExpression::read(dict); + phiName_ = dict.lookupOrDefault("phi", "phi"); - dict.lookup("field") >> fieldName_; + + resultName_ = "blendingFactor:" + fieldName_; return true; } @@ -73,27 +75,18 @@ bool Foam::functionObjects::blendingFactor::read(const dictionary& dict) bool Foam::functionObjects::blendingFactor::execute(const bool postProcess) { - calc(); - calc(); + bool processed = false; - return true; -} + processed = processed || calc(); + processed = processed || calc(); + if (!processed) + { + WarningInFunction + << "Unprocessed field " << fieldName_ << endl; + } -bool Foam::functionObjects::blendingFactor::write(const bool postProcess) -{ - const word fieldName = "blendingFactor:" + fieldName_; - - const volScalarField& blendingFactor = - obr_.lookupObject(fieldName); - - Info<< type() << " " << name() << " output:" << nl - << " writing field " << blendingFactor.name() << nl - << endl; - - blendingFactor.write(); - - return true; + return processed; } diff --git a/src/postProcessing/functionObjects/utilities/blendingFactor/blendingFactor.H b/src/postProcessing/functionObjects/field/blendingFactor/blendingFactor.H similarity index 81% rename from src/postProcessing/functionObjects/utilities/blendingFactor/blendingFactor.H rename to src/postProcessing/functionObjects/field/blendingFactor/blendingFactor.H index f0d7dc61b8..7c7c66434b 100644 --- a/src/postProcessing/functionObjects/utilities/blendingFactor/blendingFactor.H +++ b/src/postProcessing/functionObjects/field/blendingFactor/blendingFactor.H @@ -33,6 +33,7 @@ Description value is calculated via the maximum blending factor for any cell face. SeeAlso + Foam::functionObjects::fieldExpression Foam::functionObjects::fvMeshFunctionObject SourceFiles @@ -43,17 +44,12 @@ SourceFiles #ifndef functionObjects_blendingFactor_H #define functionObjects_blendingFactor_H -#include "fvMeshFunctionObject.H" -#include "volFieldsFwd.H" +#include "fieldExpression.H" // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // namespace Foam { - -// Forward declaration of classes -class objectRegistry; - namespace functionObjects { @@ -63,35 +59,19 @@ namespace functionObjects class blendingFactor : - public fvMeshFunctionObject + public fieldExpression { // Private data //- Name of flux field, default is "phi" word phiName_; - //- Field name - word fieldName_; - // Private Member Functions - //- Disallow default bitwise copy construct - blendingFactor(const blendingFactor&); - - //- Disallow default bitwise assignment - void operator=(const blendingFactor&); - - //- Return the blending factor field from the database - template - volScalarField& factor - ( - const GeometricField& field - ); - //- Calculate the blending factor template - void calc(); + bool calc(); public: @@ -122,9 +102,6 @@ public: //- Calculate the blending-factor virtual bool execute(const bool postProcess = false); - - //- Write the blending-factor - virtual bool write(const bool postProcess = false); }; diff --git a/src/postProcessing/functionObjects/utilities/blendingFactor/blendingFactorTemplates.C b/src/postProcessing/functionObjects/field/blendingFactor/blendingFactorTemplates.C similarity index 57% rename from src/postProcessing/functionObjects/utilities/blendingFactor/blendingFactorTemplates.C rename to src/postProcessing/functionObjects/field/blendingFactor/blendingFactorTemplates.C index d9262f7ffa..565b77c0d7 100644 --- a/src/postProcessing/functionObjects/utilities/blendingFactor/blendingFactorTemplates.C +++ b/src/postProcessing/functionObjects/field/blendingFactor/blendingFactorTemplates.C @@ -26,64 +26,25 @@ License #include "gaussConvectionScheme.H" #include "blendedSchemeBase.H" #include "fvcCellReduce.H" -#include "zeroGradientFvPatchFields.H" // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // template -Foam::volScalarField& Foam::functionObjects::blendingFactor::factor -( - const GeometricField& field -) +bool Foam::functionObjects::blendingFactor::calc() { - const word fieldName = "blendingFactor:" + field.name(); + typedef GeometricField FieldType; - if (!mesh_.foundObject(fieldName)) + if (!foundField(fieldName_)) { - volScalarField* factorPtr = - new volScalarField - ( - IOobject - ( - fieldName, - mesh_.time().timeName(), - mesh_, - IOobject::NO_READ, - IOobject::NO_WRITE - ), - mesh_, - dimensionedScalar("0", dimless, 0.0), - zeroGradientFvPatchScalarField::typeName - ); - - obr_.store(factorPtr); + return false; } - return - const_cast - ( - mesh_.lookupObject(fieldName) - ); -} - - -template -void Foam::functionObjects::blendingFactor::calc() -{ - typedef GeometricField fieldType; - - if (!mesh_.foundObject(fieldName_)) - { - return; - } - - const fieldType& field = mesh_.lookupObject(fieldName_); + const FieldType& field = lookupField(fieldName_); const word divScheme("div(" + phiName_ + ',' + fieldName_ + ')'); ITstream& its = mesh_.divScheme(divScheme); - const surfaceScalarField& phi = - mesh_.lookupObject(phiName_); + const surfaceScalarField& phi = lookupField(phiName_); tmp> cs = fv::convectionScheme::New(mesh_, phi, its); @@ -101,15 +62,17 @@ void Foam::functionObjects::blendingFactor::calc() << exit(FatalError); } - // retrieve the face-based blending factor + // Retrieve the face-based blending factor const blendedSchemeBase& blendedScheme = refCast>(interpScheme); - const surfaceScalarField factorf(blendedScheme.blendingFactor(field)); + tmp factorf(blendedScheme.blendingFactor(field)); - // convert into vol field whose values represent the local face maxima - volScalarField& factor = this->factor(field); - factor = fvc::cellReduce(factorf, maxEqOp()); - factor.correctBoundaryConditions(); + // Convert into vol field whose values represent the local face maxima + return store + ( + resultName_, + fvc::cellReduce(factorf, maxEqOp()) + ); } diff --git a/src/postProcessing/functionObjects/utilities/Make/files b/src/postProcessing/functionObjects/utilities/Make/files index 7e7e618e3d..e9578cdb67 100644 --- a/src/postProcessing/functionObjects/utilities/Make/files +++ b/src/postProcessing/functionObjects/utilities/Make/files @@ -9,7 +9,6 @@ writeDictionary/writeDictionary.C writeRegisteredObject/writeRegisteredObject.C scalarTransport/scalarTransport.C -blendingFactor/blendingFactor.C dsmcFields/dsmcFields.C turbulenceFields/turbulenceFields.C