diff --git a/applications/solvers/lagrangian/icoUncoupledKinematicParcelFoam/icoUncoupledKinematicParcelFoam.C b/applications/solvers/lagrangian/icoUncoupledKinematicParcelFoam/icoUncoupledKinematicParcelFoam.C index 6e9c2bfe7e..2bf1238452 100644 --- a/applications/solvers/lagrangian/icoUncoupledKinematicParcelFoam/icoUncoupledKinematicParcelFoam.C +++ b/applications/solvers/lagrangian/icoUncoupledKinematicParcelFoam/icoUncoupledKinematicParcelFoam.C @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | - \\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation + \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -22,7 +22,7 @@ License along with OpenFOAM. If not, see . Application - uncoupledKinematicParcelFoam + icoUncoupledKinematicParcelFoam Description Transient solver for the passive transport of a single kinematic diff --git a/applications/solvers/lagrangian/uncoupledKinematicParcelFoam/uncoupledKinematicParcelFoam.C b/applications/solvers/lagrangian/uncoupledKinematicParcelFoam/uncoupledKinematicParcelFoam.C index 661353669b..352678f853 100644 --- a/applications/solvers/lagrangian/uncoupledKinematicParcelFoam/uncoupledKinematicParcelFoam.C +++ b/applications/solvers/lagrangian/uncoupledKinematicParcelFoam/uncoupledKinematicParcelFoam.C @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | - \\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation + \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -26,7 +26,7 @@ Application Description Transient solver for the passive transport of a single kinematic - particle could. + particle cloud. Uses a pre- calculated velocity field to evolve the cloud. diff --git a/applications/test/dictionary/testDict b/applications/test/dictionary/testDict index 9aee1a35bf..6e52395a48 100644 --- a/applications/test/dictionary/testDict +++ b/applications/test/dictionary/testDict @@ -47,6 +47,10 @@ x 5; varName x; +//Indirection for keys +key inlet_9; + + boundaryField { Default_Boundary_Region @@ -63,6 +67,7 @@ boundaryField inlet_6 { $.inactive } // Test scoping inlet_7 { ${inactive}} // Test variable expansion inlet_8 { $inactive } + ${key} { $inactive } #include "testDictInc" diff --git a/applications/utilities/mesh/generation/snappyHexMesh/snappyHexMeshDict b/applications/utilities/mesh/generation/snappyHexMesh/snappyHexMeshDict index 4e5b20e4ab..f448faa40c 100644 --- a/applications/utilities/mesh/generation/snappyHexMesh/snappyHexMeshDict +++ b/applications/utilities/mesh/generation/snappyHexMesh/snappyHexMeshDict @@ -402,9 +402,10 @@ addLayersControls nLayerIter 50; // Max number of iterations after which relaxed meshQuality controls - // get used. Up to nRelaxIter it uses the settings in + // get used. Up to nRelaxedIter it uses the settings in // meshQualityControls, - // after nRelaxIter it uses the values in meshQualityControls::relaxed. + // after nRelaxedIter it uses the values in + // meshQualityControls::relaxed. nRelaxedIter 20; // Additional reporting: if there are just a few faces where there diff --git a/src/finiteVolume/cfdTools/general/porosityModel/DarcyForchheimer/DarcyForchheimer.C b/src/finiteVolume/cfdTools/general/porosityModel/DarcyForchheimer/DarcyForchheimer.C index 3806c56fec..472fc4019b 100644 --- a/src/finiteVolume/cfdTools/general/porosityModel/DarcyForchheimer/DarcyForchheimer.C +++ b/src/finiteVolume/cfdTools/general/porosityModel/DarcyForchheimer/DarcyForchheimer.C @@ -52,8 +52,8 @@ Foam::porosityModels::DarcyForchheimer::DarcyForchheimer ) : porosityModel(name, modelType, mesh, dict, cellZoneName), - D_(cellZoneIds_.size()), - F_(cellZoneIds_.size()), + D_(cellZoneIDs_.size()), + F_(cellZoneIDs_.size()), rhoName_(coeffs_.lookupOrDefault("rho", "rho")), muName_(coeffs_.lookupOrDefault("mu", "thermo:mu")), nuName_(coeffs_.lookupOrDefault("nu", "nu")) @@ -67,7 +67,7 @@ Foam::porosityModels::DarcyForchheimer::DarcyForchheimer if (coordSys_.R().uniform()) { - forAll (cellZoneIds_, zoneI) + forAll (cellZoneIDs_, zoneI) { D_[zoneI].setSize(1, tensor::zero); F_[zoneI].setSize(1, tensor::zero); @@ -89,9 +89,9 @@ Foam::porosityModels::DarcyForchheimer::DarcyForchheimer } else { - forAll(cellZoneIds_, zoneI) + forAll(cellZoneIDs_, zoneI) { - const labelList& cells = mesh_.cellZones()[cellZoneIds_[zoneI]]; + const labelList& cells = mesh_.cellZones()[cellZoneIDs_[zoneI]]; D_[zoneI].setSize(cells.size(), tensor::zero); F_[zoneI].setSize(cells.size(), tensor::zero); @@ -122,6 +122,24 @@ Foam::porosityModels::DarcyForchheimer::~DarcyForchheimer() // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // +void Foam::porosityModels::DarcyForchheimer::calcForce +( + const volVectorField& U, + const volScalarField& rho, + const volScalarField& mu, + vectorField& force +) const +{ + scalarField Udiag(U.size()); + vectorField Usource(U.size()); + const scalarField& V = mesh_.V(); + + apply(Udiag, Usource, V, rho, mu, U); + + force = Udiag*U - Usource; +} + + void Foam::porosityModels::DarcyForchheimer::correct ( fvVectorMatrix& UEqn @@ -194,10 +212,12 @@ void Foam::porosityModels::DarcyForchheimer::correct } -void Foam::porosityModels::DarcyForchheimer::writeData(Ostream& os) const +bool Foam::porosityModels::DarcyForchheimer::writeData(Ostream& os) const { os << indent << name_ << endl; dict_.write(os); + + return true; } diff --git a/src/finiteVolume/cfdTools/general/porosityModel/DarcyForchheimer/DarcyForchheimer.H b/src/finiteVolume/cfdTools/general/porosityModel/DarcyForchheimer/DarcyForchheimer.H index 5648fb0900..7b6ea1dd4e 100644 --- a/src/finiteVolume/cfdTools/general/porosityModel/DarcyForchheimer/DarcyForchheimer.H +++ b/src/finiteVolume/cfdTools/general/porosityModel/DarcyForchheimer/DarcyForchheimer.H @@ -143,6 +143,15 @@ public: // Member Functions + //- Calculate the porosity force + virtual void calcForce + ( + const volVectorField& U, + const volScalarField& rho, + const volScalarField& mu, + vectorField& force + ) const; + //- Add resistance virtual void correct(fvVectorMatrix& UEqn) const; @@ -165,7 +174,7 @@ public: // I-O //- Write - void writeData(Ostream& os) const; + bool writeData(Ostream& os) const; }; // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // diff --git a/src/finiteVolume/cfdTools/general/porosityModel/DarcyForchheimer/DarcyForchheimerTemplates.C b/src/finiteVolume/cfdTools/general/porosityModel/DarcyForchheimer/DarcyForchheimerTemplates.C index f785f67345..4ec9d4ecac 100644 --- a/src/finiteVolume/cfdTools/general/porosityModel/DarcyForchheimer/DarcyForchheimerTemplates.C +++ b/src/finiteVolume/cfdTools/general/porosityModel/DarcyForchheimer/DarcyForchheimerTemplates.C @@ -36,12 +36,12 @@ void Foam::porosityModels::DarcyForchheimer::apply const vectorField& U ) const { - forAll(cellZoneIds_, zoneI) + forAll(cellZoneIDs_, zoneI) { const tensorField& dZones = D_[zoneI]; const tensorField& fZones = F_[zoneI]; - const labelList& cells = mesh_.cellZones()[cellZoneIds_[zoneI]]; + const labelList& cells = mesh_.cellZones()[cellZoneIDs_[zoneI]]; forAll(cells, i) { @@ -68,12 +68,12 @@ void Foam::porosityModels::DarcyForchheimer::apply const vectorField& U ) const { - forAll(cellZoneIds_, zoneI) + forAll(cellZoneIDs_, zoneI) { const tensorField& dZones = D_[zoneI]; const tensorField& fZones = F_[zoneI]; - const labelList& cells = mesh_.cellZones()[cellZoneIds_[zoneI]]; + const labelList& cells = mesh_.cellZones()[cellZoneIDs_[zoneI]]; forAll(cells, i) { diff --git a/src/finiteVolume/cfdTools/general/porosityModel/fixedCoeff/fixedCoeff.C b/src/finiteVolume/cfdTools/general/porosityModel/fixedCoeff/fixedCoeff.C index 2ba2d99b55..ee19f0dc45 100644 --- a/src/finiteVolume/cfdTools/general/porosityModel/fixedCoeff/fixedCoeff.C +++ b/src/finiteVolume/cfdTools/general/porosityModel/fixedCoeff/fixedCoeff.C @@ -50,12 +50,12 @@ void Foam::porosityModels::fixedCoeff::apply const scalar rho ) const { - forAll(cellZoneIds_, zoneI) + forAll(cellZoneIDs_, zoneI) { const tensorField& alphaZones = alpha_[zoneI]; const tensorField& betaZones = beta_[zoneI]; - const labelList& cells = mesh_.cellZones()[cellZoneIds_[zoneI]]; + const labelList& cells = mesh_.cellZones()[cellZoneIDs_[zoneI]]; forAll(cells, i) { @@ -79,12 +79,12 @@ void Foam::porosityModels::fixedCoeff::apply ) const { - forAll(cellZoneIds_, zoneI) + forAll(cellZoneIDs_, zoneI) { const tensorField& alphaZones = alpha_[zoneI]; const tensorField& betaZones = beta_[zoneI]; - const labelList& cells = mesh_.cellZones()[cellZoneIds_[zoneI]]; + const labelList& cells = mesh_.cellZones()[cellZoneIDs_[zoneI]]; forAll(cells, i) { @@ -111,8 +111,8 @@ Foam::porosityModels::fixedCoeff::fixedCoeff ) : porosityModel(name, modelType, mesh, dict, cellZoneName), - alpha_(cellZoneIds_.size()), - beta_(cellZoneIds_.size()) + alpha_(cellZoneIDs_.size()), + beta_(cellZoneIDs_.size()) { dimensionedVector alpha(coeffs_.lookup("alpha")); dimensionedVector beta(coeffs_.lookup("beta")); @@ -122,7 +122,7 @@ Foam::porosityModels::fixedCoeff::fixedCoeff if (coordSys_.R().uniform()) { - forAll (cellZoneIds_, zoneI) + forAll (cellZoneIDs_, zoneI) { alpha_[zoneI].setSize(1, tensor::zero); beta_[zoneI].setSize(1, tensor::zero); @@ -140,9 +140,9 @@ Foam::porosityModels::fixedCoeff::fixedCoeff } else { - forAll(cellZoneIds_, zoneI) + forAll(cellZoneIDs_, zoneI) { - const labelList& cells = mesh_.cellZones()[cellZoneIds_[zoneI]]; + const labelList& cells = mesh_.cellZones()[cellZoneIDs_[zoneI]]; alpha_[zoneI].setSize(cells.size(), tensor::zero); beta_[zoneI].setSize(cells.size(), tensor::zero); @@ -175,6 +175,25 @@ Foam::porosityModels::fixedCoeff::~fixedCoeff() // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // +void Foam::porosityModels::fixedCoeff::calcForce +( + const volVectorField& U, + const volScalarField& rho, + const volScalarField& mu, + vectorField& force +) const +{ + scalarField Udiag(U.size()); + vectorField Usource(U.size()); + const scalarField& V = mesh_.V(); + scalar rhoRef = readScalar(coeffs_.lookup("rhoRef")); + + apply(Udiag, Usource, V, U, rhoRef); + + force = Udiag*U - Usource; +} + + void Foam::porosityModels::fixedCoeff::correct ( fvVectorMatrix& UEqn @@ -235,10 +254,12 @@ void Foam::porosityModels::fixedCoeff::correct } -void Foam::porosityModels::fixedCoeff::writeData(Ostream& os) const +bool Foam::porosityModels::fixedCoeff::writeData(Ostream& os) const { os << indent << name_ << endl; dict_.write(os); + + return true; } diff --git a/src/finiteVolume/cfdTools/general/porosityModel/fixedCoeff/fixedCoeff.H b/src/finiteVolume/cfdTools/general/porosityModel/fixedCoeff/fixedCoeff.H index ad00c5ff8c..487478030d 100644 --- a/src/finiteVolume/cfdTools/general/porosityModel/fixedCoeff/fixedCoeff.H +++ b/src/finiteVolume/cfdTools/general/porosityModel/fixedCoeff/fixedCoeff.H @@ -119,6 +119,15 @@ public: // Member Functions + //- Calculate the porosity force + virtual void calcForce + ( + const volVectorField& U, + const volScalarField& rho, + const volScalarField& mu, + vectorField& force + ) const; + //- Add resistance virtual void correct(fvVectorMatrix& UEqn) const; @@ -141,7 +150,7 @@ public: // I-O //- Write - void writeData(Ostream& os) const; + bool writeData(Ostream& os) const; }; // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // diff --git a/src/finiteVolume/cfdTools/general/porosityModel/porosityModel/porosityModel.C b/src/finiteVolume/cfdTools/general/porosityModel/porosityModel/porosityModel.C index 56395ee245..ad83fb0cda 100644 --- a/src/finiteVolume/cfdTools/general/porosityModel/porosityModel/porosityModel.C +++ b/src/finiteVolume/cfdTools/general/porosityModel/porosityModel/porosityModel.C @@ -88,13 +88,14 @@ Foam::porosityModel::porosityModel const word& cellZoneName ) : + MeshObject(mesh), name_(name), mesh_(mesh), dict_(dict), coeffs_(dict.subDict(modelType + "Coeffs")), active_(true), zoneName_(cellZoneName), - cellZoneIds_(), + cellZoneIDs_(), coordSys_(coordinateSystem::New(mesh, coeffs_)) { if (zoneName_ == word::null) @@ -103,11 +104,11 @@ Foam::porosityModel::porosityModel dict_.lookup("cellZone") >> zoneName_; } - cellZoneIds_ = mesh_.cellZones().findIndices(zoneName_); + cellZoneIDs_ = mesh_.cellZones().findIndices(zoneName_); Info<< " creating porous zone: " << zoneName_ << endl; - bool foundZone = !cellZoneIds_.empty(); + bool foundZone = !cellZoneIDs_.empty(); reduce(foundZone, orOp()); if (!foundZone && Pstream::master()) @@ -136,12 +137,30 @@ Foam::porosityModel::~porosityModel() // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // +Foam::tmp Foam::porosityModel::porosityModel::force +( + const volVectorField& U, + const volScalarField& rho, + const volScalarField& mu +) const +{ + tmp tforce(new vectorField(U.size(), vector::zero)); + + if (!cellZoneIDs_.empty()) + { + this->calcForce(U, rho, mu, tforce()); + } + + return tforce; +} + + void Foam::porosityModel::porosityModel::addResistance ( fvVectorMatrix& UEqn ) const { - if (cellZoneIds_.empty()) + if (cellZoneIDs_.empty()) { return; } @@ -157,7 +176,7 @@ void Foam::porosityModel::porosityModel::addResistance const volScalarField& mu ) const { - if (cellZoneIds_.empty()) + if (cellZoneIDs_.empty()) { return; } @@ -173,7 +192,7 @@ void Foam::porosityModel::porosityModel::addResistance bool correctAUprocBC ) const { - if (cellZoneIds_.empty()) + if (cellZoneIDs_.empty()) { return; } @@ -190,12 +209,30 @@ void Foam::porosityModel::porosityModel::addResistance } +bool Foam::porosityModel::porosityModel::movePoints() +{ + // no updates necessary; all member data independent of mesh + return true; +} + + +void Foam::porosityModel::porosityModel::updateMesh(const mapPolyMesh& mpm) +{ + // no updates necessary; all member data independent of mesh +} + + +bool Foam::porosityModel::writeData(Ostream& os) const +{ + return true; +} + bool Foam::porosityModel::read(const dictionary& dict) { active_ = readBool(dict.lookup("active")); coeffs_ = dict.subDict(type() + "Coeffs"); dict.lookup("cellZone") >> zoneName_; - cellZoneIds_ = mesh_.cellZones().findIndices(zoneName_); + cellZoneIDs_ = mesh_.cellZones().findIndices(zoneName_); return true; } diff --git a/src/finiteVolume/cfdTools/general/porosityModel/porosityModel/porosityModel.H b/src/finiteVolume/cfdTools/general/porosityModel/porosityModel/porosityModel.H index 0217ebd584..36a5a9eafb 100644 --- a/src/finiteVolume/cfdTools/general/porosityModel/porosityModel/porosityModel.H +++ b/src/finiteVolume/cfdTools/general/porosityModel/porosityModel/porosityModel.H @@ -36,6 +36,7 @@ SourceFiles #ifndef porosityModel_H #define porosityModel_H +#include "MeshObject.H" #include "fvMesh.H" #include "dictionary.H" #include "fvMatricesFwd.H" @@ -54,6 +55,8 @@ namespace Foam \*---------------------------------------------------------------------------*/ class porosityModel +: + public MeshObject { private: @@ -88,8 +91,8 @@ protected: //- Name(s) of cell-zone keyType zoneName_; - //- Cell zone Ids - labelList cellZoneIds_; + //- Cell zone IDs + labelList cellZoneIDs_; //- Local co-ordinate system coordinateSystem coordSys_; @@ -100,6 +103,15 @@ protected: //- Adjust negative resistance values to be multiplier of max value void adjustNegativeResistance(dimensionedVector& resist); + //- Calculate the porosity force + virtual void calcForce + ( + const volVectorField& U, + const volScalarField& rho, + const volScalarField& mu, + vectorField& force + ) const = 0; + virtual void correct(fvVectorMatrix& UEqn) const = 0; virtual void correct @@ -207,6 +219,17 @@ public: //- Return const access to the porosity active flag inline bool active() const; + //- Return const access to the cell zone IDs + inline const labelList& cellZoneIDs() const; + + //- Return the force over the cell zone(s) + virtual tmp force + ( + const volVectorField& U, + const volScalarField& rho, + const volScalarField& mu + ) const; + //- Add resistance virtual void addResistance(fvVectorMatrix& UEqn) const; @@ -227,10 +250,19 @@ public: ) const; + // Topology change + + //- Move points + virtual bool movePoints(); + + //- Update on meshUpdate + virtual void updateMesh(const mapPolyMesh& mpm); + + // I-O //- Write - virtual void writeData(Ostream& os) const = 0; + virtual bool writeData(Ostream& os) const; //- Read porosity dictionary virtual bool read(const dictionary& dict); diff --git a/src/finiteVolume/cfdTools/general/porosityModel/porosityModel/porosityModelI.H b/src/finiteVolume/cfdTools/general/porosityModel/porosityModel/porosityModelI.H index 1df2d11972..63f28091b3 100644 --- a/src/finiteVolume/cfdTools/general/porosityModel/porosityModel/porosityModelI.H +++ b/src/finiteVolume/cfdTools/general/porosityModel/porosityModel/porosityModelI.H @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | - \\ / A nd | Copyright (C) 2012 OpenFOAM Foundation + \\ / A nd | Copyright (C) 2012-2013 OpenFOAM Foundation \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -35,4 +35,10 @@ inline bool Foam::porosityModel::active() const } +inline const Foam::labelList& Foam::porosityModel::cellZoneIDs() const +{ + return cellZoneIDs_; +} + + // ************************************************************************* // diff --git a/src/finiteVolume/cfdTools/general/porosityModel/powerLaw/powerLaw.C b/src/finiteVolume/cfdTools/general/porosityModel/powerLaw/powerLaw.C index 2516c8c6b2..ee3157e43e 100644 --- a/src/finiteVolume/cfdTools/general/porosityModel/powerLaw/powerLaw.C +++ b/src/finiteVolume/cfdTools/general/porosityModel/powerLaw/powerLaw.C @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | - \\ / A nd | Copyright (C) 2012 OpenFOAM Foundation + \\ / A nd | Copyright (C) 2012-2013 OpenFOAM Foundation \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -66,6 +66,23 @@ Foam::porosityModels::powerLaw::~powerLaw() // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // +void Foam::porosityModels::powerLaw::calcForce +( + const volVectorField& U, + const volScalarField& rho, + const volScalarField& mu, + vectorField& force +) const +{ + scalarField Udiag(U.size()); + const scalarField& V = mesh_.V(); + + apply(Udiag, V, rho, U); + + force = Udiag*U; +} + + void Foam::porosityModels::powerLaw::correct ( fvVectorMatrix& UEqn @@ -126,10 +143,12 @@ void Foam::porosityModels::powerLaw::correct } -void Foam::porosityModels::powerLaw::writeData(Ostream& os) const +bool Foam::porosityModels::powerLaw::writeData(Ostream& os) const { os << indent << name_ << endl; dict_.write(os); + + return true; } diff --git a/src/finiteVolume/cfdTools/general/porosityModel/powerLaw/powerLaw.H b/src/finiteVolume/cfdTools/general/porosityModel/powerLaw/powerLaw.H index fb0797aaa7..299335b729 100644 --- a/src/finiteVolume/cfdTools/general/porosityModel/powerLaw/powerLaw.H +++ b/src/finiteVolume/cfdTools/general/porosityModel/powerLaw/powerLaw.H @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | - \\ / A nd | Copyright (C) 2012 OpenFOAM Foundation + \\ / A nd | Copyright (C) 2012-2013 OpenFOAM Foundation \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -127,6 +127,15 @@ public: // Member Functions + //- Calculate the porosity force + virtual void calcForce + ( + const volVectorField& U, + const volScalarField& rho, + const volScalarField& mu, + vectorField& force + ) const; + //- Add resistance virtual void correct(fvVectorMatrix& UEqn) const; @@ -149,7 +158,7 @@ public: // I-O //- Write - void writeData(Ostream& os) const; + bool writeData(Ostream& os) const; }; // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // diff --git a/src/finiteVolume/cfdTools/general/porosityModel/powerLaw/powerLawTemplates.C b/src/finiteVolume/cfdTools/general/porosityModel/powerLaw/powerLawTemplates.C index 0e1d7d2343..086085c565 100644 --- a/src/finiteVolume/cfdTools/general/porosityModel/powerLaw/powerLawTemplates.C +++ b/src/finiteVolume/cfdTools/general/porosityModel/powerLaw/powerLawTemplates.C @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | - \\ / A nd | Copyright (C) 2012 OpenFOAM Foundation + \\ / A nd | Copyright (C) 2012-2013 OpenFOAM Foundation \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -37,9 +37,9 @@ void Foam::porosityModels::powerLaw::apply const scalar C0 = C0_; const scalar C1m1b2 = (C1_ - 1.0)/2.0; - forAll(cellZoneIds_, zoneI) + forAll(cellZoneIDs_, zoneI) { - const labelList& cells = mesh_.cellZones()[cellZoneIds_[zoneI]]; + const labelList& cells = mesh_.cellZones()[cellZoneIDs_[zoneI]]; forAll(cells, i) { @@ -63,9 +63,9 @@ void Foam::porosityModels::powerLaw::apply const scalar C0 = C0_; const scalar C1m1b2 = (C1_ - 1.0)/2.0; - forAll(cellZoneIds_, zoneI) + forAll(cellZoneIDs_, zoneI) { - const labelList& cells = mesh_.cellZones()[cellZoneIds_[zoneI]]; + const labelList& cells = mesh_.cellZones()[cellZoneIDs_[zoneI]]; forAll(cells, i) { diff --git a/src/postProcessing/functionObjects/forces/forceCoeffs/forceCoeffs.C b/src/postProcessing/functionObjects/forces/forceCoeffs/forceCoeffs.C index 5a6d8ed4e7..9537f6a28a 100644 --- a/src/postProcessing/functionObjects/forces/forceCoeffs/forceCoeffs.C +++ b/src/postProcessing/functionObjects/forces/forceCoeffs/forceCoeffs.C @@ -120,8 +120,8 @@ void Foam::forceCoeffs::write() scalar pDyn = 0.5*rhoRef_*magUInf_*magUInf_; - Field totForce(force_[0] + force_[1]); - Field totMoment(moment_[0] + moment_[1]); + Field totForce(force_[0] + force_[1] + force_[2]); + Field totMoment(moment_[0] + moment_[1] + moment_[2]); List > coeffs(3); coeffs[0].setSize(nBin_); diff --git a/src/postProcessing/functionObjects/forces/forces/forces.C b/src/postProcessing/functionObjects/forces/forces/forces.C index 4f272cdb35..c159b652e8 100644 --- a/src/postProcessing/functionObjects/forces/forces/forces.C +++ b/src/postProcessing/functionObjects/forces/forces/forces.C @@ -37,6 +37,7 @@ License #include "compressible/RAS/RASModel/RASModel.H" #include "compressible/LES/LESModel/LESModel.H" +#include "porosityModel.H" // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // @@ -52,14 +53,15 @@ void Foam::forces::writeFileHeader(const label i) { file() << "# Time" << tab - << "forces(pressure, viscous) moment(pressure, viscous)"; + << "forces(pressure, viscous, porous) " + << "moment(pressure, viscous, porous)"; if (localSystem_) { file() << tab - << "local forces(pressure, viscous) " - << "local moment(pressure, viscous)"; + << "local forces(pressure, viscous, porous) " + << "local moment(pressure, viscous, porous)"; } file()<< endl; @@ -132,7 +134,7 @@ Foam::tmp Foam::forces::devRhoReff() const else { FatalErrorIn("forces::devRhoReff()") - << "No valid model for viscous stress calculation." + << "No valid model for viscous stress calculation" << exit(FatalError); return volSymmTensorField::null(); @@ -140,6 +142,46 @@ Foam::tmp Foam::forces::devRhoReff() const } +Foam::tmp Foam::forces::mu() const +{ + if (obr_.foundObject("thermophysicalProperties")) + { + const fluidThermo& thermo = + obr_.lookupObject("thermophysicalProperties"); + + return thermo.mu(); + } + else if + ( + obr_.foundObject("transportProperties") + ) + { + const singlePhaseTransportModel& laminarT = + obr_.lookupObject + ("transportProperties"); + + return rho()*laminarT.nu(); + } + else if (obr_.foundObject("transportProperties")) + { + const dictionary& transportProperties = + obr_.lookupObject("transportProperties"); + + dimensionedScalar nu(transportProperties.lookup("nu")); + + return rho()*nu; + } + else + { + FatalErrorIn("forces::mu()") + << "No valid model for dynamic viscosity calculation" + << exit(FatalError); + + return volScalarField::null(); + } +} + + Foam::tmp Foam::forces::rho() const { if (rhoName_ == "rhoInf") @@ -190,32 +232,35 @@ Foam::scalar Foam::forces::rho(const volScalarField& p) const void Foam::forces::applyBins ( - const label patchI, - const vectorField fN, - const vectorField Md, - const vectorField fT + const vectorField& Md, + const vectorField& fN, + const vectorField& fT, + const vectorField& fP, + const vectorField& d ) { if (nBin_ == 1) { force_[0][0] += sum(fN); force_[1][0] += sum(fT); - moment_[0][0] += sum(Md ^ fN); - moment_[1][0] += sum(Md ^ fT); + force_[2][0] += sum(fP); + moment_[0][0] += sum(Md^fN); + moment_[1][0] += sum(Md^fT); + moment_[2][0] += sum(Md^fP); } else { - const fvMesh& mesh = refCast(obr_); - const vectorField& Cf = mesh.C().boundaryField()[patchI]; - scalarField d((Cf & binDir_) - binMin_); + scalarField dd((d & binDir_) - binMin_); - forAll(d, i) + forAll(dd, i) { - label binI = floor(d[i]/binDx_); + label binI = floor(dd[i]/binDx_); force_[0][binI] += fN[i]; force_[1][binI] += fT[i]; + force_[2][binI] += fP[i]; moment_[0][binI] += Md[i]^fN[i]; moment_[1][binI] += Md[i]^fT[i]; + moment_[2][binI] += Md[i]^fP[i]; } } } @@ -239,7 +284,7 @@ void Foam::forces::writeBins() const Info<< " Writing bins to " << forcesDir << endl; } - wordList fieldNames(IStringStream("(pressure viscous)")()); + wordList fieldNames(IStringStream("(pressure viscous porous)")()); List > f(force_); List > m(moment_); @@ -250,8 +295,11 @@ void Foam::forces::writeBins() const { f[0][i] += f[0][i-1]; f[1][i] += f[1][i-1]; + f[2][i] += f[2][i-1]; + m[0][i] += m[0][i-1]; m[1][i] += m[1][i-1]; + m[2][i] += m[2][i-1]; } } @@ -264,12 +312,14 @@ void Foam::forces::writeBins() const if (localSystem_) { - List > localForce(2); - List > localMoment(2); + List > localForce(3); + List > localMoment(3); localForce[0] = coordSys_.localVector(force_[0]); localForce[1] = coordSys_.localVector(force_[1]); + localForce[2] = coordSys_.localVector(force_[2]); localMoment[0] = coordSys_.localVector(moment_[0]); localMoment[1] = coordSys_.localVector(moment_[1]); + localMoment[2] = coordSys_.localVector(moment_[2]); if (binCumulative_) { @@ -277,8 +327,10 @@ void Foam::forces::writeBins() const { localForce[0][i] += localForce[0][i-1]; localForce[1][i] += localForce[1][i-1]; + localForce[2][i] += localForce[2][i-1]; localMoment[0][i] += localMoment[0][i-1]; localMoment[1][i] += localMoment[1][i-1]; + localMoment[2][i] += localMoment[2][i-1]; } } @@ -306,8 +358,8 @@ Foam::forces::forces obr_(obr), active_(true), log_(false), - force_(2), - moment_(2), + force_(3), + moment_(3), patchSet_(), pName_(word::null), UName_(word::null), @@ -318,6 +370,7 @@ Foam::forces::forces pRef_(0), coordSys_(), localSystem_(false), + porosity_(false), nBin_(1), binDir_(vector::zero), binDx_(0.0), @@ -365,8 +418,8 @@ Foam::forces::forces obr_(obr), active_(true), log_(false), - force_(2), - moment_(2), + force_(3), + moment_(3), patchSet_(patchSet), pName_(pName), UName_(UName), @@ -377,6 +430,7 @@ Foam::forces::forces pRef_(pRef), coordSys_(coordSys), localSystem_(false), + porosity_(false), nBin_(1), binDir_(vector::zero), binDx_(0.0), @@ -475,6 +529,18 @@ void Foam::forces::read(const dictionary& dict) localSystem_ = true; } + if (dict.readIfPresent("porosity", porosity_)) + { + if (porosity_) + { + Info<< "Including porosity effects" << endl; + } + else + { + Info<< "Not including porosity effects" << endl; + } + } + if (dict.found("binData")) { const dictionary& binDict(dict.subDict("binData")); @@ -491,10 +557,11 @@ void Foam::forces::read(const dictionary& dict) else if ((nBin_ == 0) || (nBin_ == 1)) { nBin_ = 1; - force_[0].setSize(1); - force_[1].setSize(1); - moment_[0].setSize(1); - moment_[1].setSize(1); + forAll(force_, i) + { + force_[i].setSize(1); + moment_[i].setSize(1); + } } if (nBin_ > 1) @@ -546,8 +613,10 @@ void Foam::forces::read(const dictionary& dict) // allocate storage for forces and moments force_[0].setSize(1); force_[1].setSize(1); + force_[2].setSize(1); moment_[0].setSize(1); moment_[1].setSize(1); + moment_[2].setSize(1); } } } @@ -581,28 +650,50 @@ void Foam::forces::write() if (log_) { Info<< type() << " output:" << nl - << " forces(pressure,viscous)" - << "(" << sum(force_[0]) << "," << sum(force_[1]) << ")" << nl - << " moment(pressure,viscous)" - << "(" << sum(moment_[0]) << "," << sum(moment_[1]) << ")" + << " forces(pressure,viscous,porous) = (" + << sum(force_[0]) << "," + << sum(force_[1]) << "," + << sum(force_[2]) << ")" << nl + << " moment(pressure,viscous,porous) = (" + << sum(moment_[0]) << "," + << sum(moment_[1]) << "," + << sum(moment_[2]) << ")" << nl; } file() << obr_.time().value() << tab - << "(" << sum(force_[0]) << "," << sum(force_[1]) << ") " - << "(" << sum(moment_[0]) << "," << sum(moment_[1]) << ")" + << "(" + << sum(force_[0]) << "," + << sum(force_[1]) << "," + << sum(force_[2]) + << ") " + << "(" + << sum(moment_[0]) << "," + << sum(moment_[1]) << "," + << sum(moment_[2]) + << ")" << endl; if (localSystem_) { - vectorField localForceP(coordSys_.localVector(force_[0])); - vectorField localForceV(coordSys_.localVector(force_[1])); - vectorField localMomentP(coordSys_.localVector(moment_[0])); - vectorField localMomentV(coordSys_.localVector(moment_[1])); + vectorField localForceN(coordSys_.localVector(force_[0])); + vectorField localForceT(coordSys_.localVector(force_[1])); + vectorField localForceP(coordSys_.localVector(force_[2])); + vectorField localMomentN(coordSys_.localVector(moment_[0])); + vectorField localMomentT(coordSys_.localVector(moment_[1])); + vectorField localMomentP(coordSys_.localVector(moment_[2])); file() << obr_.time().value() << tab - << "(" << sum(localForceP) << "," << sum(localForceV) << ") " - << "(" << sum(localMomentP) << "," << sum(localMomentV) << ")" + << "(" + << sum(localForceN) << "," + << sum(localForceT) << "," + << sum(localForceP) + << ") " + << "(" + << sum(localMomentN) << "," + << sum(localMomentT) << "," + << sum(localMomentP) + << ")" << endl; } @@ -620,9 +711,11 @@ void Foam::forces::calcForcesMoment() { force_[0] = vector::zero; force_[1] = vector::zero; + force_[2] = vector::zero; moment_[0] = vector::zero; moment_[1] = vector::zero; + moment_[2] = vector::zero; if (directForceDensity_) { @@ -656,7 +749,10 @@ void Foam::forces::calcForcesMoment() // Tangential force (total force minus normal fN) vectorField fT(sA*fD.boundaryField()[patchI] - fN); - applyBins(patchI, fN, Md, fT); + //- Porous force + vectorField fP(Md.size(), vector::zero); + + applyBins(Md, fN, fT, fP, mesh.C().boundaryField()[patchI]); } } else @@ -685,14 +781,51 @@ void Foam::forces::calcForcesMoment() mesh.C().boundaryField()[patchI] - coordSys_.origin() ); - vectorField pf + vectorField fN ( rho(p)*Sfb[patchI]*(p.boundaryField()[patchI] - pRef) ); - vectorField vf(Sfb[patchI] & devRhoReffb[patchI]); + vectorField fT(Sfb[patchI] & devRhoReffb[patchI]); - applyBins(patchI, pf, Md, vf); + vectorField fP(Md.size(), vector::zero); + + applyBins(Md, fN, fT, fP, mesh.C().boundaryField()[patchI]); + } + } + + if (porosity_) + { + const volVectorField& U = obr_.lookupObject(UName_); + const volScalarField rho(this->rho()); + const volScalarField mu(this->mu()); + + const fvMesh& mesh = U.mesh(); + + const HashTable models = + obr_.lookupClass(); + + forAllConstIter(HashTable, models, iter) + { + const porosityModel& pm = *iter(); + + vectorField fPTot(pm.force(U, rho, mu)); + + const labelList& cellZoneIDs = pm.cellZoneIDs(); + + forAll(cellZoneIDs, i) + { + label zoneI = cellZoneIDs[i]; + const cellZone& cZone = mesh.cellZones()[zoneI]; + + const vectorField d(mesh.C(), cZone); + const vectorField fP(fPTot, cZone); + const vectorField Md(d - coordSys_.origin()); + + const vectorField fDummy(Md.size(), vector::zero); + + applyBins(Md, fDummy, fDummy, fP, d); + } } } @@ -703,13 +836,13 @@ void Foam::forces::calcForcesMoment() Foam::vector Foam::forces::forceEff() const { - return sum(force_[0]) + sum(force_[1]); + return sum(force_[0]) + sum(force_[1]) + sum(force_[2]); } Foam::vector Foam::forces::momentEff() const { - return sum(moment_[0]) + sum(moment_[1]); + return sum(moment_[0]) + sum(moment_[1]) + sum(moment_[2]); } diff --git a/src/postProcessing/functionObjects/forces/forces/forces.H b/src/postProcessing/functionObjects/forces/forces/forces.H index 1f8ff06cb4..dd87184c2f 100644 --- a/src/postProcessing/functionObjects/forces/forces/forces.H +++ b/src/postProcessing/functionObjects/forces/forces/forces.H @@ -149,10 +149,10 @@ protected: //- Switch to send output to Info as well as to file Switch log_; - //- Pressure and viscous force per bin + //- Pressure, viscous and porous force per bin List > force_; - //- Pressure and viscous pressure per bin + //- Pressure, viscous and porous moment per bin List > moment_; @@ -188,6 +188,9 @@ protected: //- Flag to indicate whether we are using a local co-ordinate sys bool localSystem_; + //- Flag to include porosity effects + bool porosity_; + // Bin information @@ -221,6 +224,9 @@ protected: //- Return the effective viscous stress (laminar + turbulent). tmp devRhoReff() const; + //- Dynamic viscosity field + tmp mu() const; + //- Return rho if rhoName is specified otherwise rhoRef tmp rho() const; @@ -231,10 +237,11 @@ protected: //- Accumulate bin data void applyBins ( - const label patchI, - const vectorField fN, - const vectorField Md, - const vectorField fT + const vectorField& Md, + const vectorField& fN, + const vectorField& fT, + const vectorField& fP, + const vectorField& d ); //- Helper function to write bin data