From b1045257208b98acb47b38ad5b141cedcff8b697 Mon Sep 17 00:00:00 2001 From: mattijs Date: Thu, 15 Sep 2011 16:04:25 +0100 Subject: [PATCH] ENH: execFlowFunctionObjects: added -region option. Added -noFlow option. --- .../postProcessing/foamCalc/foamCalcApp.C | 3 +- .../execFlowFunctionObjects.C | 418 ++++++++++-------- src/postProcessing/postCalc/postCalc.C | 8 +- 3 files changed, 253 insertions(+), 176 deletions(-) diff --git a/applications/utilities/postProcessing/foamCalc/foamCalcApp.C b/applications/utilities/postProcessing/foamCalc/foamCalcApp.C index e9fb697570..39d0a581d5 100644 --- a/applications/utilities/postProcessing/foamCalc/foamCalcApp.C +++ b/applications/utilities/postProcessing/foamCalc/foamCalcApp.C @@ -43,6 +43,7 @@ Description int main(int argc, char *argv[]) { Foam::timeSelector::addOptions(); +# include "addRegionOption.H" Foam::argList::addBoolOption ( "noWrite", @@ -74,7 +75,7 @@ int main(int argc, char *argv[]) # include "setRootCase.H" # include "createTime.H" Foam::instantList timeDirs = Foam::timeSelector::select0(runTime, args); -# include "createMesh.H" +# include "createNamedMesh.H" utility().tryPreCalc(args, runTime, mesh); diff --git a/applications/utilities/postProcessing/miscellaneous/execFlowFunctionObjects/execFlowFunctionObjects.C b/applications/utilities/postProcessing/miscellaneous/execFlowFunctionObjects/execFlowFunctionObjects.C index fe84a1daa0..f0487ac2a2 100644 --- a/applications/utilities/postProcessing/miscellaneous/execFlowFunctionObjects/execFlowFunctionObjects.C +++ b/applications/utilities/postProcessing/miscellaneous/execFlowFunctionObjects/execFlowFunctionObjects.C @@ -37,6 +37,11 @@ Description #include "calc.H" +#include "volFields.H" +#include "surfaceFields.H" +#include "pointFields.H" +#include "ReadFields.H" + #include "incompressible/singlePhaseTransportModel/singlePhaseTransportModel.H" #include "incompressible/RAS/RASModel/RASModel.H" @@ -84,203 +89,268 @@ namespace Foam void Foam::calc(const argList& args, const Time& runTime, const fvMesh& mesh) { - Info<< " Reading phi" << endl; - surfaceScalarField phi - ( - IOobject - ( - "phi", - runTime.timeName(), - mesh, - IOobject::MUST_READ - ), - mesh - ); - - Info<< " Reading U" << endl; - volVectorField U - ( - IOobject - ( - "U", - runTime.timeName(), - mesh, - IOobject::MUST_READ - ), - mesh - ); - - Info<< " Reading p" << endl; - volScalarField p - ( - IOobject - ( - "p", - runTime.timeName(), - mesh, - IOobject::MUST_READ - ), - mesh - ); - - if (phi.dimensions() == dimensionSet(0, 3, -1, 0, 0)) + if (args.optionFound("noFlow")) { - IOobject RASPropertiesHeader - ( - "RASProperties", - runTime.constant(), - mesh, - IOobject::MUST_READ_IF_MODIFIED, - IOobject::NO_WRITE, - false - ); + Info<< " Operating in no-flow mode; no models will be loaded." + << " All vol, surface and point fields will be loaded." << endl; - IOobject LESPropertiesHeader - ( - "LESProperties", - runTime.constant(), - mesh, - IOobject::MUST_READ_IF_MODIFIED, - IOobject::NO_WRITE, - false - ); + // Read objects in time directory + IOobjectList objects(mesh, runTime.timeName()); - if (RASPropertiesHeader.headerOk()) - { - IOdictionary RASProperties(RASPropertiesHeader); + // Read vol fields. - singlePhaseTransportModel laminarTransport(U, phi); + PtrList vsFlds; + ReadFields(mesh, objects, vsFlds); - autoPtr RASModel - ( - incompressible::RASModel::New - ( - U, - phi, - laminarTransport - ) - ); - execFlowFunctionObjects(args, runTime); - } - else if (LESPropertiesHeader.headerOk()) - { - IOdictionary LESProperties(LESPropertiesHeader); + PtrList vvFlds; + ReadFields(mesh, objects, vvFlds); - singlePhaseTransportModel laminarTransport(U, phi); + PtrList vstFlds; + ReadFields(mesh, objects, vstFlds); - autoPtr sgsModel - ( - incompressible::LESModel::New(U, phi, laminarTransport) - ); + PtrList vsymtFlds; + ReadFields(mesh, objects, vsymtFlds); - execFlowFunctionObjects(args, runTime); - } - else - { - IOdictionary transportProperties - ( - IOobject - ( - "transportProperties", - runTime.constant(), - mesh, - IOobject::MUST_READ_IF_MODIFIED, - IOobject::NO_WRITE - ) - ); + PtrList vtFlds; + ReadFields(mesh, objects, vtFlds); - dimensionedScalar nu(transportProperties.lookup("nu")); + // Read surface fields. - execFlowFunctionObjects(args, runTime); - } - } - else if (phi.dimensions() == dimensionSet(1, 0, -1, 0, 0)) - { - autoPtr thermo(basicPsiThermo::New(mesh)); + PtrList ssFlds; + ReadFields(mesh, objects, ssFlds); - volScalarField rho - ( - IOobject - ( - "rho", - runTime.timeName(), - mesh - ), - thermo->rho() - ); + PtrList svFlds; + ReadFields(mesh, objects, svFlds); - IOobject RASPropertiesHeader - ( - "RASProperties", - runTime.constant(), - mesh, - IOobject::MUST_READ_IF_MODIFIED, - IOobject::NO_WRITE, - false - ); + PtrList sstFlds; + ReadFields(mesh, objects, sstFlds); - IOobject LESPropertiesHeader - ( - "LESProperties", - runTime.constant(), - mesh, - IOobject::MUST_READ_IF_MODIFIED, - IOobject::NO_WRITE, - false - ); + PtrList ssymtFlds; + ReadFields(mesh, objects, ssymtFlds); - if (RASPropertiesHeader.headerOk()) - { - IOdictionary RASProperties(RASPropertiesHeader); + PtrList stFlds; + ReadFields(mesh, objects, stFlds); - autoPtr RASModel - ( - compressible::RASModel::New - ( - rho, - U, - phi, - thermo() - ) - ); + // Read point fields. + const pointMesh& pMesh = pointMesh::New(mesh); - execFlowFunctionObjects(args, runTime); - } - else if (LESPropertiesHeader.headerOk()) - { - IOdictionary LESProperties(LESPropertiesHeader); + PtrList psFlds; + ReadFields(pMesh, objects, psFlds); - autoPtr sgsModel - ( - compressible::LESModel::New(rho, U, phi, thermo()) - ); + PtrList pvFlds; + ReadFields(pMesh, objects, pvFlds); - execFlowFunctionObjects(args, runTime); - } - else - { - IOdictionary transportProperties - ( - IOobject - ( - "transportProperties", - runTime.constant(), - mesh, - IOobject::MUST_READ_IF_MODIFIED, - IOobject::NO_WRITE - ) - ); + PtrList pstFlds; + ReadFields(pMesh, objects, pstFlds); - dimensionedScalar mu(transportProperties.lookup("mu")); + PtrList psymtFlds; + ReadFields(pMesh, objects, psymtFlds); - execFlowFunctionObjects(args, runTime); - } + PtrList ptFlds; + ReadFields(pMesh, objects, ptFlds); + + execFlowFunctionObjects(args, runTime); } else { - FatalErrorIn(args.executable()) - << "Incorrect dimensions of phi: " << phi.dimensions() - << nl << exit(FatalError); + Info<< " Reading phi" << endl; + surfaceScalarField phi + ( + IOobject + ( + "phi", + runTime.timeName(), + mesh, + IOobject::MUST_READ + ), + mesh + ); + + Info<< " Reading U" << endl; + volVectorField U + ( + IOobject + ( + "U", + runTime.timeName(), + mesh, + IOobject::MUST_READ + ), + mesh + ); + + Info<< " Reading p" << endl; + volScalarField p + ( + IOobject + ( + "p", + runTime.timeName(), + mesh, + IOobject::MUST_READ + ), + mesh + ); + + if (phi.dimensions() == dimensionSet(0, 3, -1, 0, 0)) + { + IOobject RASPropertiesHeader + ( + "RASProperties", + runTime.constant(), + mesh, + IOobject::MUST_READ_IF_MODIFIED, + IOobject::NO_WRITE, + false + ); + + IOobject LESPropertiesHeader + ( + "LESProperties", + runTime.constant(), + mesh, + IOobject::MUST_READ_IF_MODIFIED, + IOobject::NO_WRITE, + false + ); + + if (RASPropertiesHeader.headerOk()) + { + IOdictionary RASProperties(RASPropertiesHeader); + + singlePhaseTransportModel laminarTransport(U, phi); + + autoPtr RASModel + ( + incompressible::RASModel::New + ( + U, + phi, + laminarTransport + ) + ); + execFlowFunctionObjects(args, runTime); + } + else if (LESPropertiesHeader.headerOk()) + { + IOdictionary LESProperties(LESPropertiesHeader); + + singlePhaseTransportModel laminarTransport(U, phi); + + autoPtr sgsModel + ( + incompressible::LESModel::New(U, phi, laminarTransport) + ); + + execFlowFunctionObjects(args, runTime); + } + else + { + IOdictionary transportProperties + ( + IOobject + ( + "transportProperties", + runTime.constant(), + mesh, + IOobject::MUST_READ_IF_MODIFIED, + IOobject::NO_WRITE + ) + ); + + dimensionedScalar nu(transportProperties.lookup("nu")); + + execFlowFunctionObjects(args, runTime); + } + } + else if (phi.dimensions() == dimensionSet(1, 0, -1, 0, 0)) + { + autoPtr thermo(basicPsiThermo::New(mesh)); + + volScalarField rho + ( + IOobject + ( + "rho", + runTime.timeName(), + mesh + ), + thermo->rho() + ); + + IOobject RASPropertiesHeader + ( + "RASProperties", + runTime.constant(), + mesh, + IOobject::MUST_READ_IF_MODIFIED, + IOobject::NO_WRITE, + false + ); + + IOobject LESPropertiesHeader + ( + "LESProperties", + runTime.constant(), + mesh, + IOobject::MUST_READ_IF_MODIFIED, + IOobject::NO_WRITE, + false + ); + + if (RASPropertiesHeader.headerOk()) + { + IOdictionary RASProperties(RASPropertiesHeader); + + autoPtr RASModel + ( + compressible::RASModel::New + ( + rho, + U, + phi, + thermo() + ) + ); + + execFlowFunctionObjects(args, runTime); + } + else if (LESPropertiesHeader.headerOk()) + { + IOdictionary LESProperties(LESPropertiesHeader); + + autoPtr sgsModel + ( + compressible::LESModel::New(rho, U, phi, thermo()) + ); + + execFlowFunctionObjects(args, runTime); + } + else + { + IOdictionary transportProperties + ( + IOobject + ( + "transportProperties", + runTime.constant(), + mesh, + IOobject::MUST_READ_IF_MODIFIED, + IOobject::NO_WRITE + ) + ); + + dimensionedScalar mu(transportProperties.lookup("mu")); + + execFlowFunctionObjects(args, runTime); + } + } + else + { + FatalErrorIn(args.executable()) + << "Incorrect dimensions of phi: " << phi.dimensions() + << nl << exit(FatalError); + } } } diff --git a/src/postProcessing/postCalc/postCalc.C b/src/postProcessing/postCalc/postCalc.C index f477818a4c..ff2dd174e0 100644 --- a/src/postProcessing/postCalc/postCalc.C +++ b/src/postProcessing/postCalc/postCalc.C @@ -57,11 +57,17 @@ namespace Foam int main(int argc, char *argv[]) { Foam::timeSelector::addOptions(); +# include "addRegionOption.H" Foam::argList::addBoolOption ( "noWrite", "suppress writing results" ); + Foam::argList::addBoolOption + ( + "noFlow", + "suppress creating flow models (execFlowFunctionObjects only)" + ); Foam::argList::addOption ( "dict", @@ -72,7 +78,7 @@ int main(int argc, char *argv[]) #include "setRootCase.H" #include "createTime.H" Foam::instantList timeDirs = Foam::timeSelector::select0(runTime, args); - #include "createMesh.H" + #include "createNamedMesh.H" forAll(timeDirs, timeI) {