From ee057765e9d6d84dbf8c693f9ccf10380e82e3cd Mon Sep 17 00:00:00 2001 From: Hrvoje Jasak Date: Fri, 8 Jul 2011 12:06:08 +0100 Subject: [PATCH 001/154] Minor clean-up --- .../polyMesh/globalMeshData/globalMeshData.C | 41 ++++++++++++------- .../polyMesh/globalMeshData/globalMeshData.H | 6 +-- .../laplaceTetDecompositionMotionSolver.H | 2 +- .../addParallelPointPatchCellDecomp.C | 4 +- .../tetPolyBoundaryMeshCellDecomp.C | 5 ++- 5 files changed, 36 insertions(+), 22 deletions(-) diff --git a/src/OpenFOAM/meshes/polyMesh/globalMeshData/globalMeshData.C b/src/OpenFOAM/meshes/polyMesh/globalMeshData/globalMeshData.C index 9da720900..fb34d12c5 100644 --- a/src/OpenFOAM/meshes/polyMesh/globalMeshData/globalMeshData.C +++ b/src/OpenFOAM/meshes/polyMesh/globalMeshData/globalMeshData.C @@ -43,7 +43,10 @@ License defineTypeNameAndDebug(Foam::globalMeshData, 0); // Geometric matching tolerance. Factor of mesh bounding box. -const Foam::scalar Foam::globalMeshData::matchTol_ = 1E-8; +const Foam::scalar Foam::globalMeshData::matchTol_ +( + debug::tolerances("globalMeshDataMatchTol", 1e-8) +); // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * // @@ -73,7 +76,6 @@ void Foam::globalMeshData::initProcAddr() } processorPatches_.setSize(nNeighbours); - if (Pstream::parRun()) { // Send indices of my processor patches to my neighbours @@ -96,7 +98,7 @@ void Foam::globalMeshData::initProcAddr() forAll(processorPatches_, i) { label patchi = processorPatches_[i]; - + IPstream fromNeighbour ( Pstream::blocking, @@ -105,7 +107,7 @@ void Foam::globalMeshData::initProcAddr() mesh_.boundaryMesh()[patchi] ).neighbProcNo() ); - + fromNeighbour >> processorPatchNeighbours_[patchi]; } } @@ -210,8 +212,8 @@ void Foam::globalMeshData::calcSharedEdges() const { // Found edge which uses shared points. Probably shared. - // Construct the edge in shared points (or rather global indices - // of the shared points) + // Construct the edge in shared points (or rather global + // indices of the shared points) edge sharedEdge ( sharedPtAddr[e0Fnd()], @@ -628,7 +630,7 @@ Foam::pointField Foam::globalMeshData::sharedPoints() const OPstream toMaster(Pstream::blocking, Pstream::masterNo()); toMaster - << sharedPointAddr_ + << sharedPointAddr_ << UIndirectList(mesh_.points(), sharedPointLabels_)(); } @@ -728,7 +730,7 @@ void Foam::globalMeshData::updateMesh() // Note: boundBox does reduce bb_ = boundBox(mesh_.points()); - scalar tolDim = matchTol_ * bb_.mag(); + scalar tolDim = matchTol_*bb_.mag(); if (debug) { @@ -747,6 +749,7 @@ void Foam::globalMeshData::updateMesh() sharedPointLabels_ = parallelPoints.sharedPointLabels(); sharedPointAddr_ = parallelPoints.sharedPointAddr(); } + //// Option 2. Geometric //{ // // Calculate all shared points. This does all the hard work. @@ -837,22 +840,32 @@ void Foam::globalMeshData::updateMesh() label patchI = processorPatches_[i]; const processorPolyPatch& procPatch = - refCast(mesh_.boundaryMesh()[patchI]); + refCast + ( + mesh_.boundaryMesh()[patchI] + ); OPstream toNeighbour(Pstream::blocking, procPatch.neighbProcNo()); toNeighbour << procPatch.localPoints(); } - // Receive patch local points and uncount if coincident (and not shared) + // Receive patch local points and uncount if coincident and not shared forAll(processorPatches_, i) { label patchI = processorPatches_[i]; const processorPolyPatch& procPatch = - refCast(mesh_.boundaryMesh()[patchI]); + refCast + ( + mesh_.boundaryMesh()[patchI] + ); - IPstream fromNeighbour(Pstream::blocking, procPatch.neighbProcNo()); + IPstream fromNeighbour + ( + Pstream::blocking, + procPatch.neighbProcNo() + ); pointField nbrPoints(fromNeighbour); @@ -876,8 +889,8 @@ void Foam::globalMeshData::updateMesh() if (stat == UNSET) { - // Mark point as visited so if point is on multiple proc - // patches it only gets uncounted once. + // Mark point as visited so if point is on multipl + // e processor patches it only gets uncounted once pointStatus.set(meshPointI, VISITED); if (pMap[patchPointI] != -1) diff --git a/src/OpenFOAM/meshes/polyMesh/globalMeshData/globalMeshData.H b/src/OpenFOAM/meshes/polyMesh/globalMeshData/globalMeshData.H index 0eeadf1a1..595fa5342 100644 --- a/src/OpenFOAM/meshes/polyMesh/globalMeshData/globalMeshData.H +++ b/src/OpenFOAM/meshes/polyMesh/globalMeshData/globalMeshData.H @@ -64,10 +64,8 @@ Description your first call to one of the access functions synchronous amongst all processors! - SourceFiles globalMeshData.C - globalMeshDataMorph.C \*---------------------------------------------------------------------------*/ @@ -297,7 +295,7 @@ public: // Processor patch addressing (be careful when not running in parallel) - //- Return list of processor patch labels + //- Return list of processor patch labels // (size of list = number of processor patches) const labelList& processorPatches() const { @@ -360,7 +358,7 @@ public: pointField sharedPoints() const; //- Like sharedPoints but keeps cyclic points separate. - // (does geometric merging; uses matchTol_*bb as merging tolerance) + // (does geometric merging; uses matchTol_*bb as merge tolerance) // Use sharedPoints() instead. pointField geometricSharedPoints() const; diff --git a/src/dynamicMesh/meshMotion/tetDecompositionMotionSolver/tetDecompositionMotionSolver/laplace/laplaceTetDecompositionMotionSolver.H b/src/dynamicMesh/meshMotion/tetDecompositionMotionSolver/tetDecompositionMotionSolver/laplace/laplaceTetDecompositionMotionSolver.H index 29d8aabbb..f79df053b 100644 --- a/src/dynamicMesh/meshMotion/tetDecompositionMotionSolver/tetDecompositionMotionSolver/laplace/laplaceTetDecompositionMotionSolver.H +++ b/src/dynamicMesh/meshMotion/tetDecompositionMotionSolver/tetDecompositionMotionSolver/laplace/laplaceTetDecompositionMotionSolver.H @@ -54,7 +54,7 @@ namespace Foam class motionDiff; /*---------------------------------------------------------------------------*\ - Class laplaceTetDecompositionMotionSolver Declaration + Class laplaceTetDecompositionMotionSolver Declaration \*---------------------------------------------------------------------------*/ class laplaceTetDecompositionMotionSolver diff --git a/src/tetDecompositionFiniteElement/tetPolyMeshCellDecomp/addParallelPointPatchCellDecomp.C b/src/tetDecompositionFiniteElement/tetPolyMeshCellDecomp/addParallelPointPatchCellDecomp.C index e7e5cad5e..c5ee6336b 100644 --- a/src/tetDecompositionFiniteElement/tetPolyMeshCellDecomp/addParallelPointPatchCellDecomp.C +++ b/src/tetDecompositionFiniteElement/tetPolyMeshCellDecomp/addParallelPointPatchCellDecomp.C @@ -56,7 +56,7 @@ public: for ( - typename SLList::const_iterator myObjectsIter = + typename SLList::const_iterator myObjectsIter = myObjects.begin(); myObjectsIter != myObjects.end(); ++myObjectsIter @@ -68,7 +68,7 @@ public: for ( - typename SLList::iterator globalObjectsIter = + typename SLList::iterator globalObjectsIter = globalObjects.begin(); globalObjectsIter != globalObjects.end(); ++globalObjectsIter diff --git a/src/tetDecompositionFiniteElement/tetPolyMeshCellDecomp/tetPolyBoundaryMesh/tetPolyBoundaryMeshCellDecomp.C b/src/tetDecompositionFiniteElement/tetPolyMeshCellDecomp/tetPolyBoundaryMesh/tetPolyBoundaryMeshCellDecomp.C index 59ec3e0ff..3b0e0ff96 100644 --- a/src/tetDecompositionFiniteElement/tetPolyMeshCellDecomp/tetPolyBoundaryMesh/tetPolyBoundaryMeshCellDecomp.C +++ b/src/tetDecompositionFiniteElement/tetPolyMeshCellDecomp/tetPolyBoundaryMesh/tetPolyBoundaryMeshCellDecomp.C @@ -91,7 +91,10 @@ tetPolyBoundaryMeshCellDecomp::globalPatch() const { if (isType(patches[patchI])) { - return refCast(patches[patchI]); + return refCast + ( + patches[patchI] + ); } } From f7c17e0efe4d8082c162638ceb52ad1e9bfc7915 Mon Sep 17 00:00:00 2001 From: Hrvoje Jasak Date: Wed, 3 Aug 2011 09:14:17 +0100 Subject: [PATCH 002/154] Protecting class data --- src/OpenFOAM/meshes/MeshObject/MeshObject.H | 3 +-- 1 file changed, 1 insertion(+), 2 deletions(-) diff --git a/src/OpenFOAM/meshes/MeshObject/MeshObject.H b/src/OpenFOAM/meshes/MeshObject/MeshObject.H index e0ea58e82..0b03a4fe0 100644 --- a/src/OpenFOAM/meshes/MeshObject/MeshObject.H +++ b/src/OpenFOAM/meshes/MeshObject/MeshObject.H @@ -56,8 +56,7 @@ class MeshObject public meshObjectBase, public regIOobject { - -protected: + // Private data // Reference to Mesh const Mesh& mesh_; From 3920c00f71c57f022408615efa09049c5c2ba473 Mon Sep 17 00:00:00 2001 From: Hrvoje Jasak Date: Wed, 3 Aug 2011 09:16:14 +0100 Subject: [PATCH 003/154] Formatting --- .../OutputFilterFunctionObject/OutputFilterFunctionObject.C | 2 ++ .../OutputFilterFunctionObject/OutputFilterFunctionObject.H | 2 +- .../turbulentIntensityKineticEnergyInletFvPatchScalarField.C | 3 +++ 3 files changed, 6 insertions(+), 1 deletion(-) diff --git a/src/OpenFOAM/db/functionObjects/OutputFilterFunctionObject/OutputFilterFunctionObject.C b/src/OpenFOAM/db/functionObjects/OutputFilterFunctionObject/OutputFilterFunctionObject.C index 072bba60a..18b547a9b 100644 --- a/src/OpenFOAM/db/functionObjects/OutputFilterFunctionObject/OutputFilterFunctionObject.C +++ b/src/OpenFOAM/db/functionObjects/OutputFilterFunctionObject/OutputFilterFunctionObject.C @@ -40,6 +40,7 @@ void Foam::OutputFilterFunctionObject::readDict() dict_.readIfPresent("storeFilter", storeFilter_); } + template void Foam::OutputFilterFunctionObject::allocateFilter() { @@ -69,6 +70,7 @@ void Foam::OutputFilterFunctionObject::allocateFilter() } } + template void Foam::OutputFilterFunctionObject::destroyFilter() { diff --git a/src/OpenFOAM/db/functionObjects/OutputFilterFunctionObject/OutputFilterFunctionObject.H b/src/OpenFOAM/db/functionObjects/OutputFilterFunctionObject/OutputFilterFunctionObject.H index 87e9b308e..b5ebdf774 100644 --- a/src/OpenFOAM/db/functionObjects/OutputFilterFunctionObject/OutputFilterFunctionObject.H +++ b/src/OpenFOAM/db/functionObjects/OutputFilterFunctionObject/OutputFilterFunctionObject.H @@ -93,7 +93,7 @@ class OutputFilterFunctionObject //- Read relevant dictionary entries void readDict(); - + //- Creates most of the data associated with this object. void allocateFilter(); diff --git a/src/finiteVolume/fields/fvPatchFields/derived/turbulentIntensityKineticEnergyInlet/turbulentIntensityKineticEnergyInletFvPatchScalarField.C b/src/finiteVolume/fields/fvPatchFields/derived/turbulentIntensityKineticEnergyInlet/turbulentIntensityKineticEnergyInletFvPatchScalarField.C index e84c58a5d..7c5bc9eff 100644 --- a/src/finiteVolume/fields/fvPatchFields/derived/turbulentIntensityKineticEnergyInlet/turbulentIntensityKineticEnergyInletFvPatchScalarField.C +++ b/src/finiteVolume/fields/fvPatchFields/derived/turbulentIntensityKineticEnergyInlet/turbulentIntensityKineticEnergyInletFvPatchScalarField.C @@ -44,6 +44,7 @@ turbulentIntensityKineticEnergyInletFvPatchScalarField intensity_(0.05) {} + Foam::turbulentIntensityKineticEnergyInletFvPatchScalarField:: turbulentIntensityKineticEnergyInletFvPatchScalarField ( @@ -58,6 +59,7 @@ turbulentIntensityKineticEnergyInletFvPatchScalarField intensity_(ptf.intensity_) {} + Foam::turbulentIntensityKineticEnergyInletFvPatchScalarField:: turbulentIntensityKineticEnergyInletFvPatchScalarField ( @@ -88,6 +90,7 @@ turbulentIntensityKineticEnergyInletFvPatchScalarField } } + Foam::turbulentIntensityKineticEnergyInletFvPatchScalarField:: turbulentIntensityKineticEnergyInletFvPatchScalarField ( From 4a48eea9a2959b86ef2518be88a14e8a5b608290 Mon Sep 17 00:00:00 2001 From: Hrvoje Jasak Date: Wed, 3 Aug 2011 09:16:19 +0100 Subject: [PATCH 004/154] Extended functionality of fixedInternalValue --- .../fixedInternalValueFvPatchField.C | 106 ++++++++++++++---- .../fixedInternalValueFvPatchField.H | 47 +++++++- .../epsilonWallFunctionFvPatchScalarField.C | 16 ++- .../omegaWallFunctionFvPatchScalarField.C | 15 ++- .../epsilonWallFunctionFvPatchScalarField.C | 12 +- .../omegaWallFunctionFvPatchScalarField.C | 15 ++- 6 files changed, 171 insertions(+), 40 deletions(-) diff --git a/src/finiteVolume/fields/fvPatchFields/derived/fixedInternalValueFvPatchField/fixedInternalValueFvPatchField.C b/src/finiteVolume/fields/fvPatchFields/derived/fixedInternalValueFvPatchField/fixedInternalValueFvPatchField.C index 5171c3b1d..d15614c16 100644 --- a/src/finiteVolume/fields/fvPatchFields/derived/fixedInternalValueFvPatchField/fixedInternalValueFvPatchField.C +++ b/src/finiteVolume/fields/fvPatchFields/derived/fixedInternalValueFvPatchField/fixedInternalValueFvPatchField.C @@ -37,20 +37,8 @@ Foam::fixedInternalValueFvPatchField::fixedInternalValueFvPatchField const DimensionedField& iF ) : - zeroGradientFvPatchField(p, iF) -{} - - -template -Foam::fixedInternalValueFvPatchField::fixedInternalValueFvPatchField -( - const fixedInternalValueFvPatchField& ptf, - const fvPatch& p, - const DimensionedField& iF, - const fvPatchFieldMapper& mapper -) -: - zeroGradientFvPatchField(ptf, p, iF, mapper) + zeroGradientFvPatchField(p, iF), + refValue_(p.size()) {} @@ -62,7 +50,31 @@ Foam::fixedInternalValueFvPatchField::fixedInternalValueFvPatchField const dictionary& dict ) : - zeroGradientFvPatchField(p, iF, dict) + zeroGradientFvPatchField(p, iF, dict), + refValue_(p.size()) +{ + if (dict.found("refValue")) + { + refValue_ = Field("refValue", dict, p.size()); + } + else + { + refValue_ = this->patchInternalField(); + } +} + + +template +Foam::fixedInternalValueFvPatchField::fixedInternalValueFvPatchField +( + const fixedInternalValueFvPatchField& ptf, + const fvPatch& p, + const DimensionedField& iF, + const fvPatchFieldMapper& mapper +) +: + zeroGradientFvPatchField(ptf, p, iF, mapper), + refValue_(ptf.refValue_, mapper) {} @@ -72,7 +84,8 @@ Foam::fixedInternalValueFvPatchField::fixedInternalValueFvPatchField const fixedInternalValueFvPatchField& fivpf ) : - zeroGradientFvPatchField(fivpf) + zeroGradientFvPatchField(fivpf), + refValue_(fivpf.refValue_) {} @@ -83,20 +96,75 @@ Foam::fixedInternalValueFvPatchField::fixedInternalValueFvPatchField const DimensionedField& iF ) : - zeroGradientFvPatchField(fivpf, iF) + zeroGradientFvPatchField(fivpf, iF), + refValue_(fivpf.refValue_) {} // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // +template +void Foam::fixedInternalValueFvPatchField::autoMap +( + const fvPatchFieldMapper& m +) +{ + fvPatchField::autoMap(m); + refValue_.autoMap(m); +} + + +template +void Foam::fixedInternalValueFvPatchField::rmap +( + const fvPatchField& ptf, + const labelList& addr +) +{ + fvPatchField::rmap(ptf, addr); + + const fixedInternalValueFvPatchField& mptf = + refCast >(ptf); + + refValue_.rmap(mptf.refValue_, addr); +} + + +template +void Foam::fixedInternalValueFvPatchField::updateCoeffs() +{ + Field& vf = const_cast& > + ( + this->internalField() + ); + + const labelList& faceCells = this->patch().faceCells(); + + // Apply the refValue into the cells next to the boundary + forAll (faceCells, faceI) + { + vf[faceCells[faceI]] = refValue_[faceI]; + } +} + + template void Foam::fixedInternalValueFvPatchField::manipulateMatrix ( fvMatrix& matrix ) { - // Apply the patch internal field as a constraint in the matrix - matrix.setValues(this->patch().faceCells(), this->patchInternalField()); + // Apply the refValue as a constraint in the matrix + matrix.setValues(this->patch().faceCells(), refValue_); +} + + +template +void Foam::fixedInternalValueFvPatchField::write(Ostream& os) const +{ + fvPatchField::write(os); + refValue_.writeEntry("refValue", os); + this->writeEntry("value", os); } diff --git a/src/finiteVolume/fields/fvPatchFields/derived/fixedInternalValueFvPatchField/fixedInternalValueFvPatchField.H b/src/finiteVolume/fields/fvPatchFields/derived/fixedInternalValueFvPatchField/fixedInternalValueFvPatchField.H index 378ffeed2..e28d98736 100644 --- a/src/finiteVolume/fields/fvPatchFields/derived/fixedInternalValueFvPatchField/fixedInternalValueFvPatchField.H +++ b/src/finiteVolume/fields/fvPatchFields/derived/fixedInternalValueFvPatchField/fixedInternalValueFvPatchField.H @@ -46,7 +46,7 @@ namespace Foam { /*---------------------------------------------------------------------------*\ - Class fixedInternalValueFvPatchField Declaration + Class fixedInternalValueFvPatchField Declaration \*---------------------------------------------------------------------------*/ template @@ -54,6 +54,11 @@ class fixedInternalValueFvPatchField : public zeroGradientFvPatchField { + // Private data + + //- Value to be set in internal field + Field refValue_; + public: @@ -125,10 +130,48 @@ public: // Member functions + // Access + + //- Return reference field value + const Field& refValue() const + { + return refValue_; + } + + //- Return access to reference field value + Field& refValue() + { + return refValue_; + } + + + // Mapping functions + + //- Map (and resize as needed) from self given a mapping object + virtual void autoMap + ( + const fvPatchFieldMapper& + ); + + //- Reverse map the given fvPatchField onto this fvPatchField + virtual void rmap + ( + const fvPatchField&, + const labelList& + ); + + // Evaluation functions - //-Manipulate a matrix + //- Update the coefficients associated with the patch field + virtual void updateCoeffs(); + + //- Manipulate a matrix virtual void manipulateMatrix(fvMatrix& matrix); + + + //- Write + virtual void write(Ostream&) const; }; diff --git a/src/turbulenceModels/compressible/RAS/derivedFvPatchFields/wallFunctions/epsilonWallFunctions/epsilonWallFunction/epsilonWallFunctionFvPatchScalarField.C b/src/turbulenceModels/compressible/RAS/derivedFvPatchFields/wallFunctions/epsilonWallFunctions/epsilonWallFunction/epsilonWallFunctionFvPatchScalarField.C index 89f79aa55..e7c273e2c 100644 --- a/src/turbulenceModels/compressible/RAS/derivedFvPatchFields/wallFunctions/epsilonWallFunctions/epsilonWallFunction/epsilonWallFunctionFvPatchScalarField.C +++ b/src/turbulenceModels/compressible/RAS/derivedFvPatchFields/wallFunctions/epsilonWallFunctions/epsilonWallFunction/epsilonWallFunctionFvPatchScalarField.C @@ -193,8 +193,10 @@ void epsilonWallFunctionFvPatchScalarField::updateCoeffs() volScalarField& G = const_cast (db().lookupObject(GName_)); - volScalarField& epsilon = const_cast - (db().lookupObject(dimensionedInternalField().name())); + // Note: epsilon is now a refValue and set in + // fixedInternalValueFvPatchField + // HJ, 3/Aug/2011 + scalarField& epsilon = refValue(); const volScalarField& k = db().lookupObject(kName_); @@ -212,16 +214,21 @@ void epsilonWallFunctionFvPatchScalarField::updateCoeffs() const scalarField magGradUw = mag(Uw.snGrad()); + const labelList& faceCells = patch().faceCells(); + // Set epsilon and G forAll(mutw, faceI) { - label faceCellI = patch().faceCells()[faceI]; + label faceCellI = faceCells[faceI]; scalar yPlus = Cmu25*y[faceI]*sqrt(k[faceCellI]) /(muw[faceI]/rhow[faceI]); - epsilon[faceCellI] = Cmu75*pow(k[faceCellI], 1.5)/(kappa_*y[faceI]); + // Note: epsilon is now a refValue and set in + // fixedInternalValueFvPatchField + // HJ, 3/Aug/2011 + epsilon[faceI] = Cmu75*pow(k[faceCellI], 1.5)/(kappa_*y[faceI]); if (yPlus > yPlusLam) { @@ -264,7 +271,6 @@ void epsilonWallFunctionFvPatchScalarField::write(Ostream& os) const os.writeKeyword("Cmu") << Cmu_ << token::END_STATEMENT << nl; os.writeKeyword("kappa") << kappa_ << token::END_STATEMENT << nl; os.writeKeyword("E") << E_ << token::END_STATEMENT << nl; - writeEntry("value", os); } diff --git a/src/turbulenceModels/compressible/RAS/derivedFvPatchFields/wallFunctions/omegaWallFunctions/omegaWallFunction/omegaWallFunctionFvPatchScalarField.C b/src/turbulenceModels/compressible/RAS/derivedFvPatchFields/wallFunctions/omegaWallFunctions/omegaWallFunction/omegaWallFunctionFvPatchScalarField.C index 51bef4764..f8f639a84 100644 --- a/src/turbulenceModels/compressible/RAS/derivedFvPatchFields/wallFunctions/omegaWallFunctions/omegaWallFunction/omegaWallFunctionFvPatchScalarField.C +++ b/src/turbulenceModels/compressible/RAS/derivedFvPatchFields/wallFunctions/omegaWallFunctions/omegaWallFunction/omegaWallFunctionFvPatchScalarField.C @@ -196,8 +196,9 @@ void omegaWallFunctionFvPatchScalarField::updateCoeffs() volScalarField& G = const_cast (db().lookupObject(GName_)); - volScalarField& omega = const_cast - (db().lookupObject(dimensionedInternalField().name())); + // Note: omega is now a refValue and set in fixedInternalValueFvPatchField + // HJ, 3/Aug/2011 + scalarField& omega = refValue(); const scalarField& k = db().lookupObject(kName_); @@ -215,10 +216,12 @@ void omegaWallFunctionFvPatchScalarField::updateCoeffs() const scalarField magGradUw = mag(Uw.snGrad()); + const labelList& faceCells = patch().faceCells(); + // Set omega and G forAll(mutw, faceI) { - label faceCellI = patch().faceCells()[faceI]; + label faceCellI = faceCells[faceI]; scalar yPlus = Cmu25*y[faceI]*sqrt(k[faceCellI]) @@ -228,7 +231,10 @@ void omegaWallFunctionFvPatchScalarField::updateCoeffs() scalar omegaLog = sqrt(k[faceCellI])/(Cmu25*kappa_*y[faceI]); - omega[faceCellI] = sqrt(sqr(omegaVis) + sqr(omegaLog)); + // Note: omega is now a refValue and set in + // fixedInternalValueFvPatchField + // HJ, 3/Aug/2011 + omega[faceI] = sqrt(sqr(omegaVis) + sqr(omegaLog)); if (yPlus > yPlusLam) { @@ -263,7 +269,6 @@ void omegaWallFunctionFvPatchScalarField::write(Ostream& os) const os.writeKeyword("kappa") << kappa_ << token::END_STATEMENT << nl; os.writeKeyword("E") << E_ << token::END_STATEMENT << nl; os.writeKeyword("beta1") << beta1_ << token::END_STATEMENT << nl; - writeEntry("value", os); } diff --git a/src/turbulenceModels/incompressible/RAS/derivedFvPatchFields/wallFunctions/epsilonWallFunctions/epsilonWallFunction/epsilonWallFunctionFvPatchScalarField.C b/src/turbulenceModels/incompressible/RAS/derivedFvPatchFields/wallFunctions/epsilonWallFunctions/epsilonWallFunction/epsilonWallFunctionFvPatchScalarField.C index 4de3a7105..8aa6e9c2d 100644 --- a/src/turbulenceModels/incompressible/RAS/derivedFvPatchFields/wallFunctions/epsilonWallFunctions/epsilonWallFunction/epsilonWallFunctionFvPatchScalarField.C +++ b/src/turbulenceModels/incompressible/RAS/derivedFvPatchFields/wallFunctions/epsilonWallFunctions/epsilonWallFunction/epsilonWallFunctionFvPatchScalarField.C @@ -187,8 +187,10 @@ void epsilonWallFunctionFvPatchScalarField::updateCoeffs() volScalarField& G = const_cast (db().lookupObject(GName_)); - volScalarField& epsilon = const_cast - (db().lookupObject(dimensionedInternalField().name())); + // Note: epsilon is now a refValue and set in + // fixedInternalValueFvPatchField + // HJ, 3/Aug/2011 + scalarField& epsilon = refValue(); const volScalarField& k = db().lookupObject(kName_); @@ -210,7 +212,10 @@ void epsilonWallFunctionFvPatchScalarField::updateCoeffs() scalar yPlus = Cmu25*y[faceI]*sqrt(k[faceCellI])/nuw[faceI]; - epsilon[faceCellI] = Cmu75*pow(k[faceCellI], 1.5)/(kappa_*y[faceI]); + // Note: epsilon is now a refValue and set in + // fixedInternalValueFvPatchField + // HJ, 3/Aug/2011 + epsilon[faceI] = Cmu75*pow(k[faceCellI], 1.5)/(kappa_*y[faceI]); if (yPlus > yPlusLam) { @@ -252,7 +257,6 @@ void epsilonWallFunctionFvPatchScalarField::write(Ostream& os) const os.writeKeyword("Cmu") << Cmu_ << token::END_STATEMENT << nl; os.writeKeyword("kappa") << kappa_ << token::END_STATEMENT << nl; os.writeKeyword("E") << E_ << token::END_STATEMENT << nl; - writeEntry("value", os); } diff --git a/src/turbulenceModels/incompressible/RAS/derivedFvPatchFields/wallFunctions/omegaWallFunctions/omegaWallFunction/omegaWallFunctionFvPatchScalarField.C b/src/turbulenceModels/incompressible/RAS/derivedFvPatchFields/wallFunctions/omegaWallFunctions/omegaWallFunction/omegaWallFunctionFvPatchScalarField.C index 0c6f547da..766dcb6d8 100644 --- a/src/turbulenceModels/incompressible/RAS/derivedFvPatchFields/wallFunctions/omegaWallFunctions/omegaWallFunction/omegaWallFunctionFvPatchScalarField.C +++ b/src/turbulenceModels/incompressible/RAS/derivedFvPatchFields/wallFunctions/omegaWallFunctions/omegaWallFunction/omegaWallFunctionFvPatchScalarField.C @@ -191,8 +191,9 @@ void omegaWallFunctionFvPatchScalarField::updateCoeffs() volScalarField& G = const_cast (db().lookupObject(GName_)); - volScalarField& omega = const_cast - (db().lookupObject(dimensionedInternalField().name())); + // Note: omega is now a refValue and set in fixedInternalValueFvPatchField + // HJ, 3/Aug/2011 + scalarField& omega = refValue(); const scalarField& k = db().lookupObject(kName_); @@ -209,10 +210,12 @@ void omegaWallFunctionFvPatchScalarField::updateCoeffs() const scalarField magGradUw = mag(Uw.snGrad()); + const labelList& faceCells = patch().faceCells(); + // Set omega and G forAll(nutw, faceI) { - label faceCellI = patch().faceCells()[faceI]; + label faceCellI = faceCells[faceI]; scalar yPlus = Cmu25*y[faceI]*sqrt(k[faceCellI])/nuw[faceI]; @@ -220,7 +223,10 @@ void omegaWallFunctionFvPatchScalarField::updateCoeffs() scalar omegaLog = sqrt(k[faceCellI])/(Cmu25*kappa_*y[faceI]); - omega[faceCellI] = sqrt(sqr(omegaVis) + sqr(omegaLog)); + // Note: omega is now a refValue and set in + // fixedInternalValueFvPatchField + // HJ, 3/Aug/2011 + omega[faceI] = sqrt(sqr(omegaVis) + sqr(omegaLog)); if (yPlus > yPlusLam) { @@ -254,7 +260,6 @@ void omegaWallFunctionFvPatchScalarField::write(Ostream& os) const os.writeKeyword("kappa") << kappa_ << token::END_STATEMENT << nl; os.writeKeyword("E") << E_ << token::END_STATEMENT << nl; os.writeKeyword("beta1") << beta1_ << token::END_STATEMENT << nl; - writeEntry("value", os); } From a2a3c47d72eebd3b6e5bb12d4571e8182fe2d4ae Mon Sep 17 00:00:00 2001 From: Hrvoje Jasak Date: Wed, 3 Aug 2011 09:17:24 +0100 Subject: [PATCH 005/154] Add access to motion time index --- src/OpenFOAM/meshes/polyMesh/polyMesh.H | 6 ++++++ 1 file changed, 6 insertions(+) diff --git a/src/OpenFOAM/meshes/polyMesh/polyMesh.H b/src/OpenFOAM/meshes/polyMesh/polyMesh.H index 4b4226aac..baded20f9 100644 --- a/src/OpenFOAM/meshes/polyMesh/polyMesh.H +++ b/src/OpenFOAM/meshes/polyMesh/polyMesh.H @@ -394,6 +394,12 @@ public: return moving_; } + //- Return current motion time index + label curMotionTimeIndex() const + { + return moving_; + } + //- Set the mesh to be moving bool moving(const bool m) { From e75b0ee525eb05670c6d0932d40b05afd86e742a Mon Sep 17 00:00:00 2001 From: Hrvoje Jasak Date: Wed, 3 Aug 2011 09:18:21 +0100 Subject: [PATCH 006/154] Formatting --- src/finiteVolume/fvMesh/fvMeshGeometry.C | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/finiteVolume/fvMesh/fvMeshGeometry.C b/src/finiteVolume/fvMesh/fvMeshGeometry.C index 1a13302d1..233504647 100644 --- a/src/finiteVolume/fvMesh/fvMeshGeometry.C +++ b/src/finiteVolume/fvMesh/fvMeshGeometry.C @@ -227,7 +227,7 @@ void fvMesh::makePhi() const << abort(FatalError); } - // Reading old time mesh motion flux if it exists and + // Reading old time mesh motion flux if it exists and // creating zero current time mesh motion flux scalar t0 = this->time().value() - this->time().deltaT().value(); From 52ddb89ac7ac176c82866fb26d720517546d3ca7 Mon Sep 17 00:00:00 2001 From: Hrvoje Jasak Date: Wed, 3 Aug 2011 09:19:40 +0100 Subject: [PATCH 007/154] Resolve naming clash; full implementation of linear upwind on processor boundaries --- .../linearUpwind/linearUpwind.C | 1 - .../linearUpwind/linearUpwindV.C | 231 ++++++++++++------ 2 files changed, 157 insertions(+), 75 deletions(-) diff --git a/src/finiteVolume/interpolation/surfaceInterpolation/limitedSchemes/linearUpwind/linearUpwind.C b/src/finiteVolume/interpolation/surfaceInterpolation/limitedSchemes/linearUpwind/linearUpwind.C index cb9f7d986..2564a6670 100644 --- a/src/finiteVolume/interpolation/surfaceInterpolation/limitedSchemes/linearUpwind/linearUpwind.C +++ b/src/finiteVolume/interpolation/surfaceInterpolation/limitedSchemes/linearUpwind/linearUpwind.C @@ -85,7 +85,6 @@ Foam::linearUpwind::correction } } - typename GeometricField:: GeometricBoundaryField& bSfCorr = sfCorr.boundaryField(); diff --git a/src/finiteVolume/interpolation/surfaceInterpolation/limitedSchemes/linearUpwind/linearUpwindV.C b/src/finiteVolume/interpolation/surfaceInterpolation/limitedSchemes/linearUpwind/linearUpwindV.C index 53b42bea4..2b0915962 100644 --- a/src/finiteVolume/interpolation/surfaceInterpolation/limitedSchemes/linearUpwind/linearUpwindV.C +++ b/src/finiteVolume/interpolation/surfaceInterpolation/limitedSchemes/linearUpwind/linearUpwindV.C @@ -38,93 +38,176 @@ Foam::linearUpwindV::correction const GeometricField& vf ) const { - const fvMesh& mesh = this->mesh(); + const fvMesh& mesh = this->mesh(); - tmp > tsfCorr - ( - new GeometricField - ( - IOobject - ( - vf.name(), - mesh.time().timeName(), - mesh - ), - mesh, - dimensioned - ( - vf.name(), - vf.dimensions(), - pTraits::zero - ) - ) - ); + tmp > tsfCorr + ( + new GeometricField + ( + IOobject + ( + "linearUpwindCorrection(" + vf.name() + ')', + mesh.time().timeName(), + mesh, + IOobject::NO_READ, + IOobject::NO_WRITE, + false + ), + mesh, + dimensioned(vf.name(), vf.dimensions(), pTraits::zero) + ) + ); - GeometricField& sfCorr = tsfCorr(); + GeometricField& sfCorr = tsfCorr(); - const surfaceScalarField& faceFlux = this->faceFlux_; - const surfaceScalarField& w = mesh.weights(); + const surfaceScalarField& faceFlux = this->faceFlux_; + const surfaceScalarField& w = mesh.weights(); - const labelList& own = mesh.owner(); - const labelList& nei = mesh.neighbour(); + const labelList& own = mesh.owner(); + const labelList& nei = mesh.neighbour(); - const vectorField& C = mesh.C(); - const vectorField& Cf = mesh.Cf(); + const volVectorField& C = mesh.C(); + const surfaceVectorField& Cf = mesh.Cf(); - GeometricField - ::type, fvPatchField, volMesh> - gradVf = gradScheme_().grad(vf); + GeometricField + ::type, fvPatchField, volMesh> + gradVf = gradScheme_().grad(vf); - forAll(faceFlux, facei) - { - vector maxCorr; + forAll(faceFlux, facei) + { + vector maxCorr; - if (faceFlux[facei] > 0.0) - { - maxCorr = - (1.0 - w[facei]) - *(vf[nei[facei]] - vf[own[facei]]); + if (faceFlux[facei] > 0) + { + maxCorr = + (1 - w[facei])*(vf[nei[facei]] - vf[own[facei]]); - sfCorr[facei] = - (Cf[facei] - C[own[facei]]) & gradVf[own[facei]]; - } - else - { - maxCorr = - w[facei]*(vf[own[facei]] - vf[nei[facei]]); + sfCorr[facei] = + (Cf[facei] - C[own[facei]]) & gradVf[own[facei]]; + } + else + { + maxCorr = + w[facei]*(vf[own[facei]] - vf[nei[facei]]); - sfCorr[facei] = - (Cf[facei] - C[nei[facei]]) & gradVf[nei[facei]]; - } + sfCorr[facei] = + (Cf[facei] - C[nei[facei]]) & gradVf[nei[facei]]; + } - scalar sfCorrs = magSqr(sfCorr[facei]); - scalar maxCorrs = sfCorr[facei] & maxCorr; + scalar sfCorrs = magSqr(sfCorr[facei]); + scalar maxCorrs = sfCorr[facei] & maxCorr; - if (sfCorrs > 0) - { - if (maxCorrs < 0) - { - sfCorr[facei] = vector::zero; - } - else if (sfCorrs > maxCorrs) - { - sfCorr[facei] *= maxCorrs/(sfCorrs + VSMALL); - } - } - else if (sfCorrs < 0) - { - if (maxCorrs > 0) - { - sfCorr[facei] = vector::zero; - } - else if (sfCorrs < maxCorrs) - { - sfCorr[facei] *= maxCorrs/(sfCorrs - VSMALL); - } - } - } + if (sfCorrs > 0) + { + if (maxCorrs < 0) + { + sfCorr[facei] = vector::zero; + } + else if (sfCorrs > maxCorrs) + { + sfCorr[facei] *= maxCorrs/(sfCorrs + VSMALL); + } + } + else if (sfCorrs < 0) + { + if (maxCorrs > 0) + { + sfCorr[facei] = vector::zero; + } + else if (sfCorrs < maxCorrs) + { + sfCorr[facei] *= maxCorrs/(sfCorrs - VSMALL); + } + } + } - return tsfCorr; + // Added missing treatment of coupled boundaries. HJ, 27/Jul/2011 + + typename GeometricField:: + GeometricBoundaryField& bSfCorr = sfCorr.boundaryField(); + + forAll(bSfCorr, patchi) + { + fvsPatchField& pSfCorr = bSfCorr[patchi]; + + if (pSfCorr.coupled()) + { + const fvPatch& p = mesh.boundary()[patchi]; + + const unallocLabelList& pOwner = p.faceCells(); + + const vectorField& pCf = Cf.boundaryField()[patchi]; + + const scalarField& pFaceFlux = faceFlux.boundaryField()[patchi]; + + const scalarField& pW = w.boundaryField()[patchi]; + + Field vfOwn = + vf.boundaryField()[patchi].patchInternalField(); + + Field vfNei = + vf.boundaryField()[patchi].patchNeighbourField(); + + + Field::type> pGradVfNei = + gradVf.boundaryField()[patchi].patchNeighbourField(); + + // Build the d-vectors + // Better version of d-vectors: Zeljko Tukovic, 25/Apr/2010 + vectorField pd = p.delta(); + + forAll(pOwner, facei) + { + vector maxCorr; + + label own = pOwner[facei]; + + if (pFaceFlux[facei] > 0) + { + maxCorr = + (1.0 - pW[facei])*(vfNei[facei] - vfOwn[facei]); + + pSfCorr[facei] = (pCf[facei] - C[own]) & gradVf[own]; + } + else + { + maxCorr = + pW[facei]*(vfOwn[facei] - vfNei[facei]); + + pSfCorr[facei] = + (pCf[facei] - pd[facei] - C[own]) & pGradVfNei[facei]; + } + + scalar pSfCorrs = magSqr(pSfCorr[facei]); + scalar maxCorrs = pSfCorr[facei] & maxCorr; + + if (pSfCorrs > 0) + { + if (maxCorrs < 0) + { + pSfCorr[facei] = vector::zero; + } + else if (pSfCorrs > maxCorrs) + { + pSfCorr[facei] *= maxCorrs/(pSfCorrs + VSMALL); + } + } + else if (pSfCorrs < 0) + { + if (maxCorrs > 0) + { + pSfCorr[facei] = vector::zero; + } + else if (pSfCorrs < maxCorrs) + { + pSfCorr[facei] *= maxCorrs/(pSfCorrs - VSMALL); + } + } + } + } + } + + return tsfCorr; } From 74cdf1f9c03b3f800be5cf5b36a87d78e913f6bc Mon Sep 17 00:00:00 2001 From: Hrvoje Jasak Date: Wed, 3 Aug 2011 09:20:09 +0100 Subject: [PATCH 008/154] Formatting --- .../finiteVolume/ddtSchemes/EulerDdtScheme/EulerDdtScheme.C | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/finiteVolume/finiteVolume/ddtSchemes/EulerDdtScheme/EulerDdtScheme.C b/src/finiteVolume/finiteVolume/ddtSchemes/EulerDdtScheme/EulerDdtScheme.C index 73cc88d29..c1e2e0579 100644 --- a/src/finiteVolume/finiteVolume/ddtSchemes/EulerDdtScheme/EulerDdtScheme.C +++ b/src/finiteVolume/finiteVolume/ddtSchemes/EulerDdtScheme/EulerDdtScheme.C @@ -1,4 +1,4 @@ -/*---------------------------------------------------------------------------* \ +/*---------------------------------------------------------------------------*\ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | From df400b9a7fe2c6d0c8ac446fd5a570517cc345dd Mon Sep 17 00:00:00 2001 From: Hrvoje Jasak Date: Wed, 3 Aug 2011 09:20:33 +0100 Subject: [PATCH 009/154] Formatting --- .../basic/mixed/mixedFvPatchField.C | 32 +++++++++---------- 1 file changed, 16 insertions(+), 16 deletions(-) diff --git a/src/finiteVolume/fields/fvPatchFields/basic/mixed/mixedFvPatchField.C b/src/finiteVolume/fields/fvPatchFields/basic/mixed/mixedFvPatchField.C index f9b38edf0..3c036ef30 100644 --- a/src/finiteVolume/fields/fvPatchFields/basic/mixed/mixedFvPatchField.C +++ b/src/finiteVolume/fields/fvPatchFields/basic/mixed/mixedFvPatchField.C @@ -47,22 +47,6 @@ mixedFvPatchField::mixedFvPatchField {} -template -mixedFvPatchField::mixedFvPatchField -( - const mixedFvPatchField& ptf, - const fvPatch& p, - const DimensionedField& iF, - const fvPatchFieldMapper& mapper -) -: - fvPatchField(ptf, p, iF, mapper), - refValue_(ptf.refValue_, mapper), - refGrad_(ptf.refGrad_, mapper), - valueFraction_(ptf.valueFraction_, mapper) -{} - - template mixedFvPatchField::mixedFvPatchField ( @@ -80,6 +64,22 @@ mixedFvPatchField::mixedFvPatchField } +template +mixedFvPatchField::mixedFvPatchField +( + const mixedFvPatchField& ptf, + const fvPatch& p, + const DimensionedField& iF, + const fvPatchFieldMapper& mapper +) +: + fvPatchField(ptf, p, iF, mapper), + refValue_(ptf.refValue_, mapper), + refGrad_(ptf.refGrad_, mapper), + valueFraction_(ptf.valueFraction_, mapper) +{} + + template mixedFvPatchField::mixedFvPatchField ( From 011b3bad4157a8a005cbb245596971b50f85e96f Mon Sep 17 00:00:00 2001 From: Hrvoje Jasak Date: Wed, 3 Aug 2011 09:20:46 +0100 Subject: [PATCH 010/154] Mesh object updates: mesh access --- .../GAMGAgglomeration/GAMGAgglomeration.C | 2 +- .../extendedLeastSquaresVectors.C | 42 +++++++++---------- .../leastSquaresGrad/leastSquaresVectors.C | 36 ++++++++-------- .../skewCorrected/skewCorrectionVectors.C | 16 +++---- 4 files changed, 47 insertions(+), 49 deletions(-) diff --git a/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGAgglomerations/GAMGAgglomeration/GAMGAgglomeration.C b/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGAgglomerations/GAMGAgglomeration/GAMGAgglomeration.C index baafa4ad4..0edde0170 100644 --- a/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGAgglomerations/GAMGAgglomeration/GAMGAgglomeration.C +++ b/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGAgglomerations/GAMGAgglomeration/GAMGAgglomeration.C @@ -223,7 +223,7 @@ const Foam::lduMesh& Foam::GAMGAgglomeration::meshLevel { if (i == 0) { - return mesh_; + return mesh(); } else { diff --git a/src/finiteVolume/finiteVolume/gradSchemes/extendedLeastSquaresGrad/extendedLeastSquaresVectors.C b/src/finiteVolume/finiteVolume/gradSchemes/extendedLeastSquaresGrad/extendedLeastSquaresVectors.C index 3ddf08095..631872126 100644 --- a/src/finiteVolume/finiteVolume/gradSchemes/extendedLeastSquaresGrad/extendedLeastSquaresVectors.C +++ b/src/finiteVolume/finiteVolume/gradSchemes/extendedLeastSquaresGrad/extendedLeastSquaresVectors.C @@ -77,20 +77,18 @@ void Foam::extendedLeastSquaresVectors::makeLeastSquaresVectors() const << endl; } - const fvMesh& mesh = mesh_; - pVectorsPtr_ = new surfaceVectorField ( IOobject ( "LeastSquaresP", - mesh_.pointsInstance(), - mesh_, + mesh().pointsInstance(), + mesh(), IOobject::NO_READ, IOobject::NO_WRITE, false ), - mesh_, + mesh(), dimensionedVector("zero", dimless/dimLength, vector::zero) ); surfaceVectorField& lsP = *pVectorsPtr_; @@ -100,24 +98,24 @@ void Foam::extendedLeastSquaresVectors::makeLeastSquaresVectors() const IOobject ( "LeastSquaresN", - mesh_.pointsInstance(), - mesh_, + mesh().pointsInstance(), + mesh(), IOobject::NO_READ, IOobject::NO_WRITE, false ), - mesh_, + mesh(), dimensionedVector("zero", dimless/dimLength, vector::zero) ); surfaceVectorField& lsN = *nVectorsPtr_; // Set local references to mesh data - const unallocLabelList& owner = mesh_.owner(); - const unallocLabelList& neighbour = mesh_.neighbour(); - const volVectorField& C = mesh_.C(); + const unallocLabelList& owner = mesh().owner(); + const unallocLabelList& neighbour = mesh().neighbour(); + const volVectorField& C = mesh().C(); // Set up temporary storage for the dd tensor (before inversion) - symmTensorField dd(mesh_.nCells(), symmTensor::zero); + symmTensorField dd(mesh().nCells(), symmTensor::zero); forAll(owner, faceI) { @@ -137,7 +135,7 @@ void Foam::extendedLeastSquaresVectors::makeLeastSquaresVectors() const // in the construction of d vectors. forAll(lsP.boundaryField(), patchI) { - const fvPatch& p = mesh_.boundary()[patchI]; + const fvPatch& p = mesh().boundary()[patchI]; const unallocLabelList& faceCells = p.faceCells(); // Better version of d-vectors: Zeljko Tukovic, 25/Apr/2010 @@ -177,7 +175,7 @@ void Foam::extendedLeastSquaresVectors::makeLeastSquaresVectors() const << exit(FatalError); } - labelList pointLabels = mesh.cells()[i].labels(mesh.faces()); + labelList pointLabels = mesh().cells()[i].labels(mesh().faces()); scalar maxDetddij = 0.0; @@ -185,15 +183,15 @@ void Foam::extendedLeastSquaresVectors::makeLeastSquaresVectors() const forAll(pointLabels, j) { - forAll(mesh.pointCells()[pointLabels[j]], k) + forAll(mesh().pointCells()[pointLabels[j]], k) { - label cellj = mesh.pointCells()[pointLabels[j]][k]; + label cellj = mesh().pointCells()[pointLabels[j]][k]; if (cellj != i) { if (cellj != -1) { - vector dCij = (mesh.C()[cellj] - mesh.C()[i]); + vector dCij = (C[cellj] - C[i]); symmTensor ddij = dd[i] + (1.0/magSqr(dCij))*sqr(dCij); @@ -214,7 +212,7 @@ void Foam::extendedLeastSquaresVectors::makeLeastSquaresVectors() const { additionalCells_[nAddCells][0] = i; additionalCells_[nAddCells++][1] = addCell; - vector dCij = mesh.C()[addCell] - mesh.C()[i]; + vector dCij = C[addCell] - C[i]; dd[i] += (1.0/magSqr(dCij))*sqr(dCij); detdd[i] = det(dd[i]); } @@ -226,7 +224,7 @@ void Foam::extendedLeastSquaresVectors::makeLeastSquaresVectors() const Info<< "max(detdd) = " << max(detdd) << endl; Info<< "min(detdd) = " << min(detdd) << endl; Info<< "average(detdd) = " << average(detdd) << endl; - Info<< "nAddCells/nCells = " << scalar(nAddCells)/mesh.nCells() << endl; + Info<< "nAddCells/nCells = " << scalar(nAddCells)/mesh().nCells() << endl; // Invert the dd tensor symmTensorField invDd = inv(dd); @@ -248,7 +246,7 @@ void Foam::extendedLeastSquaresVectors::makeLeastSquaresVectors() const forAll(lsP.boundaryField(), patchI) { - vectorField pd = mesh_.boundary()[patchI].delta(); + vectorField pd = mesh().boundary()[patchI].delta(); fvsPatchVectorField& patchLsP = lsP.boundaryField()[patchI]; @@ -269,8 +267,8 @@ void Foam::extendedLeastSquaresVectors::makeLeastSquaresVectors() const forAll(additionalCells_, i) { - vector dCij = mesh.C()[additionalCells_[i][1]] - - mesh.C()[additionalCells_[i][0]]; + vector dCij = C[additionalCells_[i][1]] + - C[additionalCells_[i][0]]; additionalVectors_[i] = (1.0/magSqr(dCij))*(invDd[additionalCells_[i][0]] & dCij); diff --git a/src/finiteVolume/finiteVolume/gradSchemes/leastSquaresGrad/leastSquaresVectors.C b/src/finiteVolume/finiteVolume/gradSchemes/leastSquaresGrad/leastSquaresVectors.C index 3138098d8..6b7b3ac68 100644 --- a/src/finiteVolume/finiteVolume/gradSchemes/leastSquaresGrad/leastSquaresVectors.C +++ b/src/finiteVolume/finiteVolume/gradSchemes/leastSquaresGrad/leastSquaresVectors.C @@ -72,13 +72,13 @@ void Foam::leastSquaresVectors::makeLeastSquaresVectors() const IOobject ( "LeastSquaresP", - mesh_.pointsInstance(), - mesh_, + mesh().pointsInstance(), + mesh(), IOobject::NO_READ, IOobject::NO_WRITE, false ), - mesh_, + mesh(), dimensionedVector("zero", dimless/dimLength, vector::zero) ); surfaceVectorField& lsP = *pVectorsPtr_; @@ -88,28 +88,28 @@ void Foam::leastSquaresVectors::makeLeastSquaresVectors() const IOobject ( "LeastSquaresN", - mesh_.pointsInstance(), - mesh_, + mesh().pointsInstance(), + mesh(), IOobject::NO_READ, IOobject::NO_WRITE, false ), - mesh_, + mesh(), dimensionedVector("zero", dimless/dimLength, vector::zero) ); surfaceVectorField& lsN = *nVectorsPtr_; // Set local references to mesh data - const unallocLabelList& owner = mesh_.owner(); - const unallocLabelList& neighbour = mesh_.neighbour(); + const unallocLabelList& owner = mesh().owner(); + const unallocLabelList& neighbour = mesh().neighbour(); - const volVectorField& C = mesh_.C(); - const surfaceScalarField& w = mesh_.weights(); -// const surfaceScalarField& magSf = mesh_.magSf(); + const volVectorField& C = mesh().C(); + const surfaceScalarField& w = mesh().weights(); +// const surfaceScalarField& magSf = mesh().magSf(); // Set up temporary storage for the dd tensor (before inversion) - symmTensorField dd(mesh_.nCells(), symmTensor::zero); + symmTensorField dd(mesh().nCells(), symmTensor::zero); forAll(owner, facei) { @@ -138,16 +138,16 @@ void Foam::leastSquaresVectors::makeLeastSquaresVectors() const // Original version: closest distance to boundary vectorField pd = - mesh_.Sf().boundaryField()[patchi] + mesh().Sf().boundaryField()[patchi] /( - mesh_.magSf().boundaryField()[patchi] - *mesh_.deltaCoeffs().boundaryField()[patchi] + mesh().magSf().boundaryField()[patchi] + *mesh().deltaCoeffs().boundaryField()[patchi] ); - if (!mesh_.orthogonal()) + if (!mesh().orthogonal()) { - pd -= mesh_.correctionVectors().boundaryField()[patchi] - /mesh_.deltaCoeffs().boundaryField()[patchi]; + pd -= mesh().correctionVectors().boundaryField()[patchi] + /mesh().deltaCoeffs().boundaryField()[patchi]; } // Better version of d-vectors: Zeljko Tukovic, 25/Apr/2010 diff --git a/src/finiteVolume/interpolation/surfaceInterpolation/schemes/skewCorrected/skewCorrectionVectors.C b/src/finiteVolume/interpolation/surfaceInterpolation/schemes/skewCorrected/skewCorrectionVectors.C index d3a78c805..a9ac84b94 100644 --- a/src/finiteVolume/interpolation/surfaceInterpolation/schemes/skewCorrected/skewCorrectionVectors.C +++ b/src/finiteVolume/interpolation/surfaceInterpolation/schemes/skewCorrected/skewCorrectionVectors.C @@ -67,24 +67,24 @@ void Foam::skewCorrectionVectors::makeSkewCorrectionVectors() const IOobject ( "skewCorrectionVectors", - mesh_.pointsInstance(), - mesh_, + mesh().pointsInstance(), + mesh(), IOobject::NO_READ, IOobject::NO_WRITE, false ), - mesh_, + mesh(), dimLength ); surfaceVectorField& SkewCorrVecs = *skewCorrectionVectors_; // Set local references to mesh data - const volVectorField& C = mesh_.C(); - const surfaceVectorField& Cf = mesh_.Cf(); - const surfaceVectorField& Sf = mesh_.Sf(); + const volVectorField& C = mesh().C(); + const surfaceVectorField& Cf = mesh().Cf(); + const surfaceVectorField& Sf = mesh().Sf(); - const unallocLabelList& owner = mesh_.owner(); - const unallocLabelList& neighbour = mesh_.neighbour(); + const unallocLabelList& owner = mesh().owner(); + const unallocLabelList& neighbour = mesh().neighbour(); // Build the d-vectors. Changed to exact vectors. HJ, 24/Apr/2010 From 97baced83c5fabf5b9b63767b7a6e3b3ec226fe1 Mon Sep 17 00:00:00 2001 From: Hrvoje Jasak Date: Wed, 3 Aug 2011 09:22:45 +0100 Subject: [PATCH 011/154] Bug fix for rotor motion in parallel --- .../multiMixerFvMesh/mixerRotor.C | 76 ++----------------- 1 file changed, 6 insertions(+), 70 deletions(-) diff --git a/src/dynamicMesh/topoChangerFvMesh/multiMixerFvMesh/mixerRotor.C b/src/dynamicMesh/topoChangerFvMesh/multiMixerFvMesh/mixerRotor.C index 503be205d..a3197a2c4 100644 --- a/src/dynamicMesh/topoChangerFvMesh/multiMixerFvMesh/mixerRotor.C +++ b/src/dynamicMesh/topoChangerFvMesh/multiMixerFvMesh/mixerRotor.C @@ -28,7 +28,6 @@ License #include "regionSplit.H" #include "polyTopoChanger.H" #include "slidingInterface.H" -#include "ggiPolyPatch.H" #include "Time.H" // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * // @@ -222,17 +221,8 @@ void Foam::mixerRotor::calcMovingMask() const const cellList& c = mesh_.cells(); const faceList& f = mesh_.allFaces(); - label movingCellZoneID = - mesh_.cellZones().findZoneID("movingCellsZone" + name_); - - if (movingCellZoneID < 0) - { - FatalErrorIn("void mixerRotor::calcMovingMask() const") - << "Cannot find moving cell zone ID" - << abort(FatalError); - } - - const labelList& cellAddr = mesh_.cellZones()[movingCellZoneID]; + const labelList& cellAddr = mesh_.cellZones() + [mesh_.cellZones().findZoneID("movingCellsZone" + name_)]; forAll (cellAddr, cellI) { @@ -251,35 +241,8 @@ void Foam::mixerRotor::calcMovingMask() const } // Attempt to enforce motion on sliders if zones exist - - // Master side - label msI = -1; - - if (useTopoSliding_) - { - // For topological changes, find the zone - msI = mesh_.faceZones().findZoneID(movingSliderName_ + "Zone" + name_); - } - else - { - // For GGI, take face zone from ggi interpolator - label movingSliderIndex = - mesh_.boundaryMesh().findPatchID(movingSliderName_); - - if (movingSliderIndex > -1) - { - if (isA(mesh_.boundaryMesh()[movingSliderIndex])) - { - const ggiPolyPatch& movingSliderGgi = - refCast - ( - mesh_.boundaryMesh()[movingSliderIndex] - ); - - msI = mesh_.faceZones().findZoneID(movingSliderGgi.zoneName()); - } - } - } + const label msI = + mesh_.faceZones().findZoneID(movingSliderName_ + "Zone" + name_); if (msI > -1) { @@ -296,35 +259,8 @@ void Foam::mixerRotor::calcMovingMask() const } } - // Slave side - label ssI = -1; - - if (useTopoSliding_) - { - // For topological changes, find the zone - ssI = mesh_.faceZones().findZoneID(staticSliderName_ + "Zone" + name_); - } - else - { - // For GGI, take face zone from ggi interpolator - label staticSliderIndex = - mesh_.boundaryMesh().findPatchID(staticSliderName_); - - if (staticSliderIndex > -1) - { - if (isA(mesh_.boundaryMesh()[staticSliderIndex])) - { - const ggiPolyPatch& staticSliderGgi = - refCast - ( - mesh_.boundaryMesh()[staticSliderIndex] - ); - - ssI = mesh_.faceZones().findZoneID(staticSliderGgi.zoneName()); - } - } - } - + const label ssI = + mesh_.faceZones().findZoneID(staticSliderName_ + "Zone" + name_); if (ssI > -1) { From 9d030bef9154dac7fe3cb2215f7882a325488315 Mon Sep 17 00:00:00 2001 From: Hrvoje Jasak Date: Wed, 3 Aug 2011 09:23:02 +0100 Subject: [PATCH 012/154] New form of moving body topo motion --- .../mixerFvMesh/mixerFvMesh.C | 310 +++++++------ .../movingBodyTopoFvMesh.C | 306 +++++++++++++ .../movingBodyTopoFvMesh.H} | 71 ++- .../movingConeTopoFvMesh.C | 429 ------------------ .../icoDyMFoam/movingConeTopo/Allclean | 2 +- .../icoDyMFoam/movingConeTopo/Allrun | 6 +- .../icoDyMFoam/movingConeTopo/setBatch | 5 + .../movingConeTopo/system/controlDict | 4 +- .../movingConeTopo/system/decomposeParDict | 42 ++ 9 files changed, 576 insertions(+), 599 deletions(-) create mode 100644 src/dynamicMesh/topoChangerFvMesh/movingBodyTopoFvMesh/movingBodyTopoFvMesh.C rename src/dynamicMesh/topoChangerFvMesh/{movingConeTopoFvMesh/movingConeTopoFvMesh.H => movingBodyTopoFvMesh/movingBodyTopoFvMesh.H} (63%) delete mode 100644 src/dynamicMesh/topoChangerFvMesh/movingConeTopoFvMesh/movingConeTopoFvMesh.C create mode 100644 tutorials/incompressible/icoDyMFoam/movingConeTopo/system/decomposeParDict diff --git a/src/dynamicMesh/topoChangerFvMesh/mixerFvMesh/mixerFvMesh.C b/src/dynamicMesh/topoChangerFvMesh/mixerFvMesh/mixerFvMesh.C index eda0674b7..6c7b60b3a 100644 --- a/src/dynamicMesh/topoChangerFvMesh/mixerFvMesh/mixerFvMesh.C +++ b/src/dynamicMesh/topoChangerFvMesh/mixerFvMesh/mixerFvMesh.C @@ -46,50 +46,35 @@ namespace Foam void Foam::mixerFvMesh::addZonesAndModifiers() { // Add zones and modifiers for motion action - - if - ( - pointZones().size() > 0 - || faceZones().size() > 0 - || cellZones().size() > 0 - ) + if (cellZones().size() > 0) { Info<< "void mixerFvMesh::addZonesAndModifiers() : " << "Zones and modifiers already present. Skipping." << endl; - if (topoChanger_.size() == 0) + // Check definition of the modifier + if + ( + pointZones().size() > 0 + || faceZones().size() > 0 + ) { - FatalErrorIn - ( - "void mixerFvMesh::addZonesAndModifiers()" - ) << "Mesh modifiers not read properly" - << abort(FatalError); + if (topoChanger_.size() == 0) + { + FatalErrorIn + ( + "void mixerFvMesh::addZonesAndModifiers()" + ) << "Mesh modifiers not read properly. " + << "pointZones = " << pointZones().size() + << " faceZones = " << faceZones().size() + << abort(FatalError); + } } return; } - Info<< "Time = " << time().timeName() << endl - << "Adding zones and modifiers to the mesh" << endl; - - // Add zones - List pz(1); - - // Add an empty zone for cut points - - pz[0] = new pointZone - ( - "cutPointZone", - labelList(0), - 0, - pointZones() - ); - - - // Do face zones for slider - - List fz(3); + // Find patches and check sizes // Moving slider const word movingSliderName(dict_.subDict("slider").lookup("moving")); @@ -114,51 +99,84 @@ void Foam::mixerFvMesh::addZonesAndModifiers() } const polyPatch& movingSlider = boundaryMesh()[movingSliderIndex]; - - labelList isf(movingSlider.size()); - - forAll (isf, i) - { - isf[i] = movingSlider.start() + i; - } - - fz[0] = new faceZone - ( - movingSliderName + "Zone", - isf, - boolList(movingSlider.size(), false), - 0, - faceZones() - ); - - // Static slider const polyPatch& staticSlider = boundaryMesh()[staticSliderIndex]; - labelList osf(staticSlider.size()); + List pz; + List fz; - forAll (osf, i) + bool addSlider = false; + + if (movingSlider.size() > 0 && staticSlider.size() > 0) { - osf[i] = staticSlider.start() + i; + addSlider = true; + + pz.setSize(1); + fz.setSize(3); + + Info<< "Time = " << time().timeName() << endl + << "Adding zones and modifiers to the mesh" << endl; + + // Add zones + + // Add an empty zone for cut points + + pz[0] = new pointZone + ( + "cutPointZone", + labelList(0), + 0, + pointZones() + ); + + + // Do face zones for slider + + // Moving slider + labelList isf(movingSlider.size()); + + forAll (isf, i) + { + isf[i] = movingSlider.start() + i; + } + + fz[0] = new faceZone + ( + movingSliderName + "Zone", + isf, + boolList(movingSlider.size(), false), + 0, + faceZones() + ); + + // Static slider + labelList osf(staticSlider.size()); + + forAll (osf, i) + { + osf[i] = staticSlider.start() + i; + } + + fz[1] = new faceZone + ( + staticSliderName + "Zone", + osf, + boolList(staticSlider.size(), false), + 1, + faceZones() + ); + + // Add empty zone for cut faces + fz[2] = new faceZone + ( + "cutFaceZone", + labelList(0), + boolList(0, false), + 2, + faceZones() + ); } - fz[1] = new faceZone - ( - staticSliderName + "Zone", - osf, - boolList(staticSlider.size(), false), - 1, - faceZones() - ); - - // Add empty zone for cut faces - fz[2] = new faceZone - ( - "cutFaceZone", - labelList(0), - boolList(0, false), - 2, - faceZones() - ); + // Add cell zone even if slider does not exist List cz(1); @@ -194,28 +212,32 @@ void Foam::mixerFvMesh::addZonesAndModifiers() Info << "Adding point, face and cell zones" << endl; addZones(pz, fz, cz); - // Add a topology modifier - Info << "Adding topology modifiers" << endl; - topoChanger_.setSize(1); - topoChanger_.set - ( - 0, - new slidingInterface + if (addSlider) + { + // Add a topology modifier + Info << "Adding topology modifiers" << endl; + topoChanger_.setSize(1); + + topoChanger_.set ( - "mixerSlider", 0, - topoChanger_, - staticSliderName + "Zone", - movingSliderName + "Zone", - "cutPointZone", - "cutFaceZone", - staticSliderName, - movingSliderName, - slidingInterface::INTEGRAL, // Edge matching algorithm - attachDetach_, // Attach-detach action - intersection::VISIBLE // Projection algorithm - ) - ); + new slidingInterface + ( + "mixerSlider", + 0, + topoChanger_, + staticSliderName + "Zone", + movingSliderName + "Zone", + "cutPointZone", + "cutFaceZone", + staticSliderName, + movingSliderName, + slidingInterface::INTEGRAL, // Edge matching algorithm + attachDetach_, // Attach-detach action + intersection::VISIBLE // Projection algorithm + ) + ); + } // Write mesh and modifiers topoChanger_.writeOpt() = IOobject::AUTO_WRITE; @@ -226,7 +248,16 @@ void Foam::mixerFvMesh::addZonesAndModifiers() bool Foam::mixerFvMesh::attached() const { - return refCast(topoChanger_[0]).attached(); + bool result = false; + + if (topoChanger_.size() > 0) + { + result = refCast(topoChanger_[0]).attached(); + } + + reduce(result, orOp()); + + return result; } @@ -272,41 +303,43 @@ void Foam::mixerFvMesh::calcMovingMask() const } } - const word movingSliderZoneName - ( - word(dict_.subDict("slider").lookup("moving")) - + "Zone" - ); - - const labelList& movingSliderAddr = - faceZones()[faceZones().findZoneID(movingSliderZoneName)]; - - forAll (movingSliderAddr, faceI) + if (topoChanger_.size() > 0) { - const face& curFace = f[movingSliderAddr[faceI]]; + // Topo changer present. Use zones for motion + const word movingSliderZoneName + ( + word(dict_.subDict("slider").lookup("moving")) + "Zone" + ); - forAll (curFace, pointI) + const labelList& movingSliderAddr = + faceZones()[faceZones().findZoneID(movingSliderZoneName)]; + + forAll (movingSliderAddr, faceI) { - movingPointsMask[curFace[pointI]] = 1; + const face& curFace = f[movingSliderAddr[faceI]]; + + forAll (curFace, pointI) + { + movingPointsMask[curFace[pointI]] = 1; + } } - } - const word staticSliderZoneName - ( - word(dict_.subDict("slider").lookup("static")) - + "Zone" - ); + const word staticSliderZoneName + ( + word(dict_.subDict("slider").lookup("static")) + "Zone" + ); - const labelList& staticSliderAddr = - faceZones()[faceZones().findZoneID(staticSliderZoneName)]; + const labelList& staticSliderAddr = + faceZones()[faceZones().findZoneID(staticSliderZoneName)]; - forAll (staticSliderAddr, faceI) - { - const face& curFace = f[staticSliderAddr[faceI]]; - - forAll (curFace, pointI) + forAll (staticSliderAddr, faceI) { - movingPointsMask[curFace[pointI]] = 0; + const face& curFace = f[staticSliderAddr[faceI]]; + + forAll (curFace, pointI) + { + movingPointsMask[curFace[pointI]] = 0; + } } } } @@ -415,22 +448,43 @@ bool Foam::mixerFvMesh::update() autoPtr topoChangeMap = topoChanger_.changeMesh(); - if (topoChangeMap->morphing()) + bool localMorphing = topoChangeMap->morphing(); + bool globalMorphing = localMorphing; + + reduce(globalMorphing, orOp()); + + if (globalMorphing) { Info << "Attaching rotors" << endl; deleteDemandDrivenData(movingPointsMaskPtr_); // Move the sliding interface points to correct position - pointField mappedOldPointsNew(allPoints().size()); - mappedOldPointsNew.map(oldPointsNew, topoChangeMap->pointMap()); + if (localMorphing) + { + pointField mappedOldPointsNew(allPoints().size()); + mappedOldPointsNew.map(oldPointsNew, topoChangeMap->pointMap()); - movePoints(mappedOldPointsNew); - resetMotion(); - setV0(); + movePoints(mappedOldPointsNew); + + resetMotion(); + setV0(); + + // Move the sliding interface points to correct position + movePoints(topoChangeMap->preMotionPoints()); + } + else + { + pointField newPoints = allPoints(); + movePoints(oldPointsNew); + + resetMotion(); + setV0(); + + // Move the sliding interface points to correct position + movePoints(newPoints); + } - // Move the sliding interface points to correct position - movePoints(topoChangeMap->preMotionPoints()); } return topoChangeMap->morphing(); diff --git a/src/dynamicMesh/topoChangerFvMesh/movingBodyTopoFvMesh/movingBodyTopoFvMesh.C b/src/dynamicMesh/topoChangerFvMesh/movingBodyTopoFvMesh/movingBodyTopoFvMesh.C new file mode 100644 index 000000000..b9f4cdb58 --- /dev/null +++ b/src/dynamicMesh/topoChangerFvMesh/movingBodyTopoFvMesh/movingBodyTopoFvMesh.C @@ -0,0 +1,306 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright held by original author + \\/ 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 "movingBodyTopoFvMesh.H" +#include "Time.H" +#include "mapPolyMesh.H" +#include "layerAdditionRemoval.H" +#include "volMesh.H" +#include "transformField.H" +#include "addToRunTimeSelectionTable.H" + +// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // + +namespace Foam +{ + defineTypeNameAndDebug(movingBodyTopoFvMesh, 0); + + addToRunTimeSelectionTable + ( + topoChangerFvMesh, + movingBodyTopoFvMesh, + IOobject + ); +} + + +// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * // + +Foam::tmp +Foam::movingBodyTopoFvMesh::calcMotionMask() const +{ + Info<< "Updating vertex markup" << endl; + + tmp tvertexMarkup(new scalarField(allPoints().size(), 0)); + scalarField& vertexMarkup = tvertexMarkup(); + + cellZoneID movingCellsID(movingCellsName_, cellZones()); + + // In order to do a correct update on a mask on processor boundaries, + // Detection of moving cells should use patchNeighbourField for + // processor (not coupled!) boundaries. This is done by expanding + // a moving cell set into a field and making sure that processor patch + // points move in sync. Not done at the moment, probably best to do + // using parallel update of pointFields. HJ, 19/Feb/2011 + + // If moving cells are found, perform mark-up + if (movingCellsID.active()) + { + // Get cell-point addressing + const labelListList& cp = cellPoints(); + + // Get labels of all moving cells + const labelList& movingCells = cellZones()[movingCellsID.index()]; + + forAll (movingCells, cellI) + { + const labelList& curCp = cp[movingCells[cellI]]; + + forAll (curCp, pointI) + { + vertexMarkup[curCp[pointI]] = 1; + } + } + } + + faceZoneID frontFacesID(frontFacesName_, faceZones()); + + if (frontFacesID.active()) + { + const faceZone& frontFaces = faceZones()[frontFacesID.index()]; + + const labelList& mp = frontFaces().meshPoints(); + + forAll (mp, mpI) + { + vertexMarkup[mp[mpI]] = 1; + } + } + + faceZoneID backFacesID(backFacesName_, faceZones()); + + if (backFacesID.active()) + { + const faceZone& backFaces = faceZones()[backFacesID.index()]; + + const labelList& mp = backFaces().meshPoints(); + + forAll (mp, mpI) + { + vertexMarkup[mp[mpI]] = 1; + } + } + + return tvertexMarkup; +} + + +void Foam::movingBodyTopoFvMesh::addZonesAndModifiers() +{ + // Add zones and modifiers for motion action + + if (topoChanger_.size() > 0) + { + Info<< "void movingBodyTopoFvMesh::addZonesAndModifiers() : " + << "Zones and modifiers already present. Skipping." + << endl; + + return; + } + + // Add layer addition/removal interfaces + topoChanger_.setSize(2); + label nMods = 0; + + + faceZoneID frontFacesID(frontFacesName_, faceZones()); + faceZoneID backFacesID(backFacesName_, faceZones()); + + if (frontFacesID.active()) + { + const faceZone& frontFaces = faceZones()[frontFacesID.index()]; + + if (!frontFaces.empty()) + { + topoChanger_.set + ( + nMods, + new layerAdditionRemoval + ( + frontFacesName_ + "Layer", + nMods, + topoChanger_, + frontFacesName_, + readScalar + ( + dict_.subDict("front").lookup("minThickness") + ), + readScalar + ( + dict_.subDict("front").lookup("maxThickness") + ) + ) + ); + + nMods++; + } + } + + if (backFacesID.active()) + { + const faceZone& backFaces = faceZones()[backFacesID.index()]; + + if (!backFaces.empty()) + { + topoChanger_.set + ( + nMods, + new layerAdditionRemoval + ( + backFacesName_ + "Layer", + nMods, + topoChanger_, + backFacesName_, + readScalar + ( + dict_.subDict("back").lookup("minThickness") + ), + readScalar + ( + dict_.subDict("back").lookup("maxThickness") + ) + ) + ); + + nMods++; + } + } + + topoChanger_.setSize(nMods); + + reduce(nMods, sumOp