diff --git a/applications/utilities/mesh/generation/extrudeMesh/Make/options b/applications/utilities/mesh/generation/extrudeMesh/Make/options index ce0e27f401..1cdc68a2c7 100644 --- a/applications/utilities/mesh/generation/extrudeMesh/Make/options +++ b/applications/utilities/mesh/generation/extrudeMesh/Make/options @@ -1,10 +1,14 @@ EXE_INC = \ -IextrudedMesh \ -IextrudeModel/lnInclude \ + -I$(LIB_SRC)/finiteVolume/lnInclude \ + -I$(LIB_SRC)/surfMesh/lnInclude \ -I$(LIB_SRC)/meshTools/lnInclude \ -I$(LIB_SRC)/dynamicMesh/lnInclude EXE_LIBS = \ + -lfiniteVolume \ + -lsurfMesh \ -lmeshTools \ -ldynamicMesh \ -lextrudeModel diff --git a/applications/utilities/mesh/generation/extrudeMesh/extrudeMesh.C b/applications/utilities/mesh/generation/extrudeMesh/extrudeMesh.C index d149f85923..4d5bca0e82 100644 --- a/applications/utilities/mesh/generation/extrudeMesh/extrudeMesh.C +++ b/applications/utilities/mesh/generation/extrudeMesh/extrudeMesh.C @@ -36,13 +36,15 @@ Description #include "Time.H" #include "dimensionedTypes.H" #include "IFstream.H" -#include "faceMesh.H" #include "polyTopoChange.H" #include "polyTopoChanger.H" #include "edgeCollapser.H" #include "mathematicalConstants.H" #include "globalMeshData.H" #include "perfectInterface.H" +#include "addPatchCellLayer.H" +#include "fvMesh.H" +#include "MeshedSurfaces.H" #include "extrudedMesh.H" #include "extrudeModel.H" @@ -50,14 +52,148 @@ Description using namespace Foam; // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +enum ExtrudeMode +{ + MESH, + PATCH, + SURFACE +}; + +template<> +const char* NamedEnum::names[] = +{ + "mesh", + "patch", + "surface" +}; +static const NamedEnum ExtrudeModeNames; + + +void createDummyFvMeshFiles(const polyMesh& mesh, const word& regionName) +{ + // Create dummy system/fv* + { + IOobject io + ( + "fvSchemes", + mesh.time().system(), + regionName, + mesh, + IOobject::NO_READ, + IOobject::NO_WRITE, + false + ); + + Info<< "Testing:" << io.objectPath() << endl; + + if (!io.headerOk()) + { + Info<< "Writing dummy " << regionName/io.name() << endl; + dictionary dummyDict; + dictionary divDict; + dummyDict.add("divSchemes", divDict); + dictionary gradDict; + dummyDict.add("gradSchemes", gradDict); + dictionary laplDict; + dummyDict.add("laplacianSchemes", laplDict); + + IOdictionary(io, dummyDict).regIOobject::write(); + } + } + { + IOobject io + ( + "fvSolution", + mesh.time().system(), + regionName, + mesh, + IOobject::NO_READ, + IOobject::NO_WRITE, + false + ); + + if (!io.headerOk()) + { + Info<< "Writing dummy " << regionName/io.name() << endl; + dictionary dummyDict; + IOdictionary(io, dummyDict).regIOobject::write(); + } + } +} + + +label findPatchID(const polyBoundaryMesh& patches, const word& name) +{ + label patchID = patches.findPatchID(name); + + if (patchID == -1) + { + FatalErrorIn("findPatchID(const polyBoundaryMesh&, const word&)") + << "Cannot find patch " << name + << " in the source mesh.\n" + << "Valid patch names are " << patches.names() + << exit(FatalError); + } + return patchID; +} + + +labelList patchFaces(const polyBoundaryMesh& patches, const word& name) +{ + label patchID = findPatchID(patches, name); + const polyPatch& pp = patches[patchID]; + + return identity(pp.size()) + pp.start(); +} + + +void updateFaceLabels(const mapPolyMesh& map, labelList& faceLabels) +{ + const labelList& reverseMap = map.reverseFaceMap(); + + labelList newFaceLabels(faceLabels.size()); + label newI = 0; + + forAll(faceLabels, i) + { + label oldFaceI = faceLabels[i]; + + if (reverseMap[oldFaceI] >= 0) + { + newFaceLabels[newI++] = reverseMap[oldFaceI]; + } + } + newFaceLabels.setSize(newI); + faceLabels.transfer(newFaceLabels); +} + + + // Main program: int main(int argc, char *argv[]) { + #include "addRegionOption.H" #include "setRootCase.H" #include "createTimeExtruded.H" - autoPtr meshPtr(NULL); + // Get optional regionName + word regionName; + word regionDir; + if (args.optionReadIfPresent("region", regionName)) + { + regionDir = regionName; + Info<< "Create mesh " << regionName << " for time = " + << runTimeExtruded.timeName() << nl << endl; + } + else + { + regionName = fvMesh::defaultRegion; + Info<< "Create mesh for time = " + << runTimeExtruded.timeName() << nl << endl; + } + IOdictionary dict ( @@ -65,26 +201,44 @@ int main(int argc, char *argv[]) ( "extrudeProperties", runTimeExtruded.constant(), + regionDir, runTimeExtruded, IOobject::MUST_READ ) ); + // Point generator autoPtr model(extrudeModel::New(dict)); - const word sourceType(dict.lookup("constructFrom")); + // Whether to flip normals + const Switch flipNormals(dict.lookup("flipNormals")); - autoPtr fMesh; + // What to extrude + const ExtrudeMode mode = ExtrudeModeNames.read + ( + dict.lookup("constructFrom") + ); - if (sourceType == "patch") + + // Generated mesh (one of either) + autoPtr meshFromMesh; + autoPtr meshFromSurface; + + // Faces on front and back for stitching (in case of mergeFaces) + word frontPatchName; + labelList frontPatchFaces; + word backPatchName; + labelList backPatchFaces; + + if (mode == PATCH || mode == MESH) { fileName sourceCasePath(dict.lookup("sourceCase")); sourceCasePath.expand(); fileName sourceRootDir = sourceCasePath.path(); fileName sourceCaseDir = sourceCasePath.name(); - word patchName(dict.lookup("sourcePatch")); + dict.lookup("sourcePatch") >> frontPatchName; - Info<< "Extruding patch " << patchName + Info<< "Extruding patch " << frontPatchName << " on mesh " << sourceCasePath << nl << endl; @@ -94,31 +248,183 @@ int main(int argc, char *argv[]) sourceRootDir, sourceCaseDir ); - #include "createPolyMesh.H" + #include "createMesh.H" - label patchID = mesh.boundaryMesh().findPatchID(patchName); + const polyBoundaryMesh& patches = mesh.boundaryMesh(); - if (patchID == -1) + + // Topo change container. Either copy an existing mesh or start + // with empty storage (number of patches only needed for checking) + autoPtr meshMod + ( + ( + mode == MESH + ? new polyTopoChange(mesh) + : new polyTopoChange(patches.size()) + ) + ); + + // Extrusion engine. Either adding to existing mesh or + // creating separate mesh. + addPatchCellLayer layerExtrude(mesh, (mode == MESH)); + + indirectPrimitivePatch extrudePatch + ( + IndirectList + ( + mesh.faces(), + patchFaces(patches, frontPatchName) + ), + mesh.points() + ); + + // Only used for addPatchCellLayer into new mesh + labelList exposedPatchIDs; + if (mode == PATCH) { - FatalErrorIn(args.executable()) - << "Cannot find patch " << patchName - << " in the source mesh.\n" - << "Valid patch names are " << mesh.boundaryMesh().names() - << exit(FatalError); + dict.lookup("exposedPatchName") >> backPatchName; + exposedPatchIDs.setSize + ( + extrudePatch.size(), + findPatchID(patches, backPatchName) + ); } - const polyPatch& pp = mesh.boundaryMesh()[patchID]; - fMesh.reset(new faceMesh(pp.localFaces(), pp.localPoints())); + pointField layer0Points(extrudePatch.nPoints()); + pointField displacement(extrudePatch.nPoints()); + forAll(displacement, pointI) { - fileName surfName(runTime.path()/patchName + ".sMesh"); - Info<< "Writing patch as surfaceMesh to " - << surfName << nl << endl; - OFstream os(surfName); - os << fMesh() << nl; + const vector& patchNormal = extrudePatch.pointNormals()[pointI]; + + // layer0 point + layer0Points[pointI] = model() + ( + extrudePatch.localPoints()[pointI], + patchNormal, + 0 + ); + // layerN point + point extrudePt = model() + ( + extrudePatch.localPoints()[pointI], + patchNormal, + model().nLayers() + ); + displacement[pointI] = extrudePt - layer0Points[pointI]; } + + if (flipNormals) + { + Info<< "Flipping faces." << nl << endl; + displacement = -displacement; + } + + // Check if wedge (has layer0 different from original patch points) + // If so move the mesh to starting position. + if (gMax(mag(layer0Points-extrudePatch.localPoints())) > SMALL) + { + Info<< "Moving mesh to layer0 points since differ from original" + << " points - this can happen for wedge extrusions." << nl + << endl; + + pointField newPoints(mesh.points()); + forAll(extrudePatch.meshPoints(), i) + { + newPoints[extrudePatch.meshPoints()[i]] = layer0Points[i]; + } + mesh.movePoints(newPoints); + } + + + // Layers per face + labelList nFaceLayers(extrudePatch.size(), model().nLayers()); + // Layers per point + labelList nPointLayers(extrudePatch.nPoints(), model().nLayers()); + // Displacement for first layer + vectorField firstLayerDisp = displacement*model().sumThickness(1); + // Expansion ratio not used. + scalarField ratio(extrudePatch.nPoints(), 1.0); + + layerExtrude.setRefinement + ( + ratio, // expansion ratio + extrudePatch, // patch faces to extrude + exposedPatchIDs, // if new mesh: patches for exposed faces + nFaceLayers, + nPointLayers, + firstLayerDisp, + meshMod() + ); + + // Reset points according to extrusion model + forAll(layerExtrude.addedPoints(), pointI) + { + const labelList& pPoints = layerExtrude.addedPoints()[pointI]; + forAll(pPoints, pPointI) + { + label meshPointI = pPoints[pPointI]; + + point modelPt + ( + model() + ( + extrudePatch.localPoints()[pointI], + extrudePatch.pointNormals()[pointI], + pPointI+1 // layer + ) + ); + + const_cast&> + ( + meshMod().points() + )[meshPointI] = modelPt; + } + } + + // Store faces on front and exposed patch (if mode=patch there are + // only added faces so cannot used map to old faces) + const labelListList& layerFaces = layerExtrude.layerFaces(); + backPatchFaces.setSize(layerFaces.size()); + frontPatchFaces.setSize(layerFaces.size()); + forAll(backPatchFaces, i) + { + backPatchFaces[i] = layerFaces[i][0]; + frontPatchFaces[i] = layerFaces[i][layerFaces[i].size()-1]; + } + + // Create dummy fvSchemes, fvSolution + createDummyFvMeshFiles(mesh, regionDir); + + // Create actual mesh from polyTopoChange container + autoPtr map = meshMod().makeMesh + ( + meshFromMesh, + IOobject + ( + regionName, + runTimeExtruded.constant(), + runTimeExtruded, + IOobject::NO_READ, + IOobject::AUTO_WRITE, + false + ), + mesh + ); + + // Calculate face labels for front and back. + frontPatchFaces = renumber + ( + map().reverseFaceMap(), + frontPatchFaces + ); + backPatchFaces = renumber + ( + map().reverseFaceMap(), + backPatchFaces + ); } - else if (sourceType == "surface") + else { // Read from surface fileName surfName(dict.lookup("surface")); @@ -126,52 +432,62 @@ int main(int argc, char *argv[]) Info<< "Extruding surfaceMesh read from file " << surfName << nl << endl; - IFstream is(surfName); + MeshedSurface fMesh(surfName); - fMesh.reset(new faceMesh(is)); - - Info<< "Read patch from file " << surfName << nl - << endl; - } - else - { - FatalErrorIn(args.executable()) - << "Illegal 'constructFrom' specification. Should either be " - << "patch or surface." << exit(FatalError); - } - - Switch flipNormals(dict.lookup("flipNormals")); - - if (flipNormals) - { - Info<< "Flipping faces." << nl << endl; - - faceList faces(fMesh().size()); - forAll(faces, i) + if (flipNormals) { - faces[i] = fMesh()[i].reverseFace(); + Info<< "Flipping faces." << nl << endl; + faceList& faces = const_cast(fMesh.faces()); + forAll(faces, i) + { + faces[i] = fMesh[i].reverseFace(); + } } - fMesh.reset(new faceMesh(faces, fMesh().localPoints())); + + Info<< "Extruding surface with :" << nl + << " points : " << fMesh.points().size() << nl + << " faces : " << fMesh.size() << nl + << " normals[0] : " << fMesh.faceNormals()[0] + << nl + << endl; + + meshFromSurface.reset + ( + new extrudedMesh + ( + IOobject + ( + extrudedMesh::defaultRegion, + runTimeExtruded.constant(), + runTimeExtruded + ), + fMesh, + model() + ) + ); + + + // Get the faces on front and back + frontPatchName = "originalPatch"; + frontPatchFaces = patchFaces + ( + meshFromSurface().boundaryMesh(), + frontPatchName + ); + backPatchName = "otherSide"; + backPatchFaces = patchFaces + ( + meshFromSurface().boundaryMesh(), + backPatchName + ); } - Info<< "Extruding patch with :" << nl - << " points : " << fMesh().points().size() << nl - << " faces : " << fMesh().size() << nl - << " normals[0] : " << fMesh().faceNormals()[0] - << nl - << endl; - - extrudedMesh mesh + polyMesh& mesh = ( - IOobject - ( - extrudedMesh::defaultRegion, - runTimeExtruded.constant(), - runTimeExtruded - ), - fMesh(), - model() + meshFromMesh.valid() + ? meshFromMesh() + : meshFromSurface() ); @@ -184,17 +500,6 @@ int main(int argc, char *argv[]) << "Merge distance : " << mergeDim << nl << endl; - const polyBoundaryMesh& patches = mesh.boundaryMesh(); - - const label origPatchID = patches.findPatchID("originalPatch"); - const label otherPatchID = patches.findPatchID("otherSide"); - - if (origPatchID == -1 || otherPatchID == -1) - { - FatalErrorIn(args.executable()) - << "Cannot find patch originalPatch or otherSide." << nl - << "Valid patches are " << patches.names() << exit(FatalError); - } // Collapse edges // ~~~~~~~~~~~~~~ @@ -237,6 +542,10 @@ int main(int argc, char *argv[]) // Update fields mesh.updateMesh(map); + // Update stored data + updateFaceLabels(map(), frontPatchFaces); + updateFaceLabels(map(), backPatchFaces); + // Move mesh (if inflation used) if (map().hasMotionPoints()) { @@ -252,22 +561,33 @@ int main(int argc, char *argv[]) Switch mergeFaces(dict.lookup("mergeFaces")); if (mergeFaces) { + if (mode == MESH) + { + FatalErrorIn(args.executable()) + << "Cannot stitch front and back of extrusion since" + << " in 'mesh' mode (extrusion appended to mesh)." + << exit(FatalError); + } + Info<< "Assuming full 360 degree axisymmetric case;" << " stitching faces on patches " - << patches[origPatchID].name() << " and " - << patches[otherPatchID].name() << " together ..." << nl << endl; + << frontPatchName << " and " + << backPatchName << " together ..." << nl << endl; + + if (frontPatchFaces.size() != backPatchFaces.size()) + { + FatalErrorIn(args.executable()) + << "Differing number of faces on front (" + << frontPatchFaces.size() << ") and back (" + << backPatchFaces.size() << ")" + << exit(FatalError); + } + + polyTopoChanger stitcher(mesh); stitcher.setSize(1); - // Make list of masterPatch faces - labelList isf(patches[origPatchID].size()); - - forAll (isf, i) - { - isf[i] = patches[origPatchID].start() + i; - } - const word cutZoneName("originalCutFaceZone"); List fz @@ -276,8 +596,8 @@ int main(int argc, char *argv[]) new faceZone ( cutZoneName, - isf, - boolList(isf.size(), false), + frontPatchFaces, + boolList(frontPatchFaces.size(), false), 0, mesh.faceZones() ) @@ -286,26 +606,58 @@ int main(int argc, char *argv[]) mesh.addZones(List(0), fz, List(0)); // Add the perfect interface mesh modifier - stitcher.set + perfectInterface perfectStitcher ( + "couple", 0, - new perfectInterface - ( - "couple", - 0, - stitcher, - cutZoneName, - patches[origPatchID].name(), - patches[otherPatchID].name() - ) + stitcher, + cutZoneName, + word::null, // dummy patch name + word::null // dummy patch name ); - // Execute all polyMeshModifiers - autoPtr morphMap = stitcher.changeMesh(true); + // Topo change container + polyTopoChange meshMod(mesh); - mesh.movePoints(morphMap->preMotionPoints()); + perfectStitcher.setRefinement + ( + indirectPrimitivePatch + ( + IndirectList + ( + mesh.faces(), + frontPatchFaces + ), + mesh.points() + ), + indirectPrimitivePatch + ( + IndirectList + ( + mesh.faces(), + backPatchFaces + ), + mesh.points() + ), + meshMod + ); + + // Construct new mesh from polyTopoChange. + autoPtr map = meshMod.changeMesh(mesh, false); + + // Update fields + mesh.updateMesh(map); + + // Move mesh (if inflation used) + if (map().hasMotionPoints()) + { + mesh.movePoints(map().preMotionPoints()); + } } + mesh.setInstance(runTimeExtruded.constant()); + Info<< "Writing mesh to " << mesh.objectPath() << nl << endl; + if (!mesh.write()) { FatalErrorIn(args.executable()) << "Failed writing mesh" diff --git a/applications/utilities/mesh/generation/extrudeMesh/extrudeModel/extrudeModel/extrudeModel.C b/applications/utilities/mesh/generation/extrudeMesh/extrudeModel/extrudeModel/extrudeModel.C index 7305431b96..5550c56b15 100644 --- a/applications/utilities/mesh/generation/extrudeMesh/extrudeModel/extrudeModel/extrudeModel.C +++ b/applications/utilities/mesh/generation/extrudeMesh/extrudeModel/extrudeModel/extrudeModel.C @@ -43,6 +43,7 @@ Foam::extrudeModel::extrudeModel ) : nLayers_(readLabel(dict.lookup("nLayers"))), + expansionRatio_(readScalar(dict.lookup("expansionRatio"))), dict_(dict), coeffDict_(dict.subDict(modelType + "Coeffs")) {} @@ -62,5 +63,28 @@ Foam::label Foam::extrudeModel::nLayers() const } +Foam::scalar Foam::extrudeModel::expansionRatio() const +{ + return expansionRatio_; +} + + +Foam::scalar Foam::extrudeModel::sumThickness(const label layer) const +{ + // 1+r+r^2+ .. +r^(n-1) = (1-r^n)/(1-r) + + if (mag(1.0-expansionRatio_) < SMALL) + { + return scalar(layer)/nLayers_; + } + else + { + return + (1.0-pow(expansionRatio_, layer)) + / (1.0-pow(expansionRatio_, nLayers_)); + } +} + + // ************************************************************************* // diff --git a/applications/utilities/mesh/generation/extrudeMesh/extrudeModel/extrudeModel/extrudeModel.H b/applications/utilities/mesh/generation/extrudeMesh/extrudeModel/extrudeModel/extrudeModel.H index d8f3699fd5..ac96bbeea6 100644 --- a/applications/utilities/mesh/generation/extrudeMesh/extrudeModel/extrudeModel/extrudeModel.H +++ b/applications/utilities/mesh/generation/extrudeMesh/extrudeModel/extrudeModel/extrudeModel.H @@ -58,6 +58,8 @@ protected: const label nLayers_; + const scalar expansionRatio_; + const dictionary& dict_; const dictionary& coeffDict_; @@ -113,9 +115,15 @@ public: label nLayers() const; + scalar expansionRatio() const; + // Member Operators + //- Helper: calculate cumulative relative thickness for layer. + // (layer=0 -> 0; layer=nLayers -> 1) + scalar sumThickness(const label layer) const; + virtual point operator() ( const point& surfacePoint, diff --git a/applications/utilities/mesh/generation/extrudeMesh/extrudeModel/linearNormal/linearNormal.C b/applications/utilities/mesh/generation/extrudeMesh/extrudeModel/linearNormal/linearNormal.C index c36098ab1e..cd42e91816 100644 --- a/applications/utilities/mesh/generation/extrudeMesh/extrudeModel/linearNormal/linearNormal.C +++ b/applications/utilities/mesh/generation/extrudeMesh/extrudeModel/linearNormal/linearNormal.C @@ -72,7 +72,8 @@ point linearNormal::operator() const label layer ) const { - scalar d = thickness_*layer/nLayers_; + //scalar d = thickness_*layer/nLayers_; + scalar d = thickness_*sumThickness(layer); return surfacePoint + d*surfaceNormal; } diff --git a/applications/utilities/mesh/generation/extrudeMesh/extrudeModel/linearNormal/linearNormal.H b/applications/utilities/mesh/generation/extrudeMesh/extrudeModel/linearNormal/linearNormal.H index 949547311b..2dbcdccb16 100644 --- a/applications/utilities/mesh/generation/extrudeMesh/extrudeModel/linearNormal/linearNormal.H +++ b/applications/utilities/mesh/generation/extrudeMesh/extrudeModel/linearNormal/linearNormal.H @@ -64,7 +64,7 @@ public: // Constructors - //- Construct from components + //- Construct from dictionary linearNormal(const dictionary& dict); diff --git a/applications/utilities/mesh/generation/extrudeMesh/extrudeModel/linearRadial/linearRadial.C b/applications/utilities/mesh/generation/extrudeMesh/extrudeModel/linearRadial/linearRadial.C index e8d406b110..7d4fcff9ba 100644 --- a/applications/utilities/mesh/generation/extrudeMesh/extrudeModel/linearRadial/linearRadial.C +++ b/applications/utilities/mesh/generation/extrudeMesh/extrudeModel/linearRadial/linearRadial.C @@ -67,8 +67,7 @@ point linearRadial::operator() scalar rs = mag(surfacePoint); vector rsHat = surfacePoint/rs; - scalar delta = (R_ - rs)/nLayers_; - scalar r = rs + layer*delta; + scalar r = rs + (R_ - rs)*sumThickness(layer); return r*rsHat; } diff --git a/applications/utilities/mesh/generation/extrudeMesh/extrudeModel/linearRadial/linearRadial.H b/applications/utilities/mesh/generation/extrudeMesh/extrudeModel/linearRadial/linearRadial.H index 0fda12e1b1..66758c1b12 100644 --- a/applications/utilities/mesh/generation/extrudeMesh/extrudeModel/linearRadial/linearRadial.H +++ b/applications/utilities/mesh/generation/extrudeMesh/extrudeModel/linearRadial/linearRadial.H @@ -61,7 +61,7 @@ public: // Constructors - //- Construct from components + //- Construct from dictionary linearRadial(const dictionary& dict); diff --git a/applications/utilities/mesh/generation/extrudeMesh/extrudeModel/sigmaRadial/sigmaRadial.C b/applications/utilities/mesh/generation/extrudeMesh/extrudeModel/sigmaRadial/sigmaRadial.C index fc6f040305..d35ad4e539 100644 --- a/applications/utilities/mesh/generation/extrudeMesh/extrudeModel/sigmaRadial/sigmaRadial.C +++ b/applications/utilities/mesh/generation/extrudeMesh/extrudeModel/sigmaRadial/sigmaRadial.C @@ -49,7 +49,13 @@ sigmaRadial::sigmaRadial(const dictionary& dict) RTbyg_(readScalar(coeffDict_.lookup("RTbyg"))), pRef_(readScalar(coeffDict_.lookup("pRef"))), pStrat_(readScalar(coeffDict_.lookup("pStrat"))) -{} +{ + if (mag(expansionRatio() - 1.0) > SMALL) + { + WarningIn("sigmaRadial::sigmaRadial(const dictionary&)") + << "Ignoring expansionRatio setting." << endl; + } +} // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * // diff --git a/applications/utilities/mesh/generation/extrudeMesh/extrudeModel/sigmaRadial/sigmaRadial.H b/applications/utilities/mesh/generation/extrudeMesh/extrudeModel/sigmaRadial/sigmaRadial.H index 525e89d11c..c74ad02457 100644 --- a/applications/utilities/mesh/generation/extrudeMesh/extrudeModel/sigmaRadial/sigmaRadial.H +++ b/applications/utilities/mesh/generation/extrudeMesh/extrudeModel/sigmaRadial/sigmaRadial.H @@ -63,7 +63,7 @@ public: // Constructors - //- Construct from components + //- Construct from dictionary sigmaRadial(const dictionary& dict); diff --git a/applications/utilities/mesh/generation/extrudeMesh/extrudeModel/wedge/wedge.C b/applications/utilities/mesh/generation/extrudeMesh/extrudeModel/wedge/wedge.C index 3d2c883ea1..5d0a3621b4 100644 --- a/applications/utilities/mesh/generation/extrudeMesh/extrudeModel/wedge/wedge.C +++ b/applications/utilities/mesh/generation/extrudeMesh/extrudeModel/wedge/wedge.C @@ -88,8 +88,8 @@ point wedge::operator() } else { - //sliceAngle = angle_*(layer + 1)/nLayers_; - sliceAngle = angle_*layer/nLayers_; + //sliceAngle = angle_*layer/nLayers_; + sliceAngle = angle_*sumThickness(layer); } // Find projection onto axis (or rather decompose surfacePoint diff --git a/applications/utilities/mesh/generation/extrudeMesh/extrudeModel/wedge/wedge.H b/applications/utilities/mesh/generation/extrudeMesh/extrudeModel/wedge/wedge.H index 39c4b4a2d3..18707c32d1 100644 --- a/applications/utilities/mesh/generation/extrudeMesh/extrudeModel/wedge/wedge.H +++ b/applications/utilities/mesh/generation/extrudeMesh/extrudeModel/wedge/wedge.H @@ -77,7 +77,7 @@ public: // Constructors - //- Construct from components + //- Construct from dictionary wedge(const dictionary& dict); diff --git a/applications/utilities/mesh/generation/extrudeMesh/extrudeProperties b/applications/utilities/mesh/generation/extrudeMesh/extrudeProperties index 45618034de..e1283ab8b5 100644 --- a/applications/utilities/mesh/generation/extrudeMesh/extrudeProperties +++ b/applications/utilities/mesh/generation/extrudeMesh/extrudeProperties @@ -14,24 +14,28 @@ FoamFile } // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // -// Where to get surface from: either from surface ('surface') or -// from (flipped) patch of existing case ('patch') -constructFrom patch; //surface; +// What to extrude: +// patch : from patch of another case ('sourceCase') +// mesh : as above but with original case included +// surface : from externally read surface -// If construct from (flipped) patch -sourceCase "$FOAM_RUN/icoFoam/cavity"; +//constructFrom mesh; +constructFrom patch; +//constructFrom surface; + +// If construct from patch/mesh: +sourceCase "../cavity"; sourcePatch movingWall; +// If construct from patch: patch to use for back (can be same as sourcePatch) +exposedPatchName movingWall; +// If construct from surface: +surface "movingWall.stl"; + + // Flip surface normals before usage. flipNormals false; -// If construct from surface -surface "movingWall.sMesh"; - - -// Do front and back need to be merged? Usually only makes sense for 360 -// degree wedges. -mergeFaces true; //- Linear extrusion in point-normal direction @@ -46,11 +50,13 @@ extrudeModel wedge; //- Extrudes into sphere with grading according to pressure (atmospherics) //extrudeModel sigmaRadial; -nLayers 20; +nLayers 10; + +expansionRatio 1.0; //0.9; wedgeCoeffs { - axisPt (0 0.1 0); + axisPt (0 0.1 -0.05); axis (-1 0 0); angle 360; // For nLayers=1 assume symmetry so angle/2 on each side } @@ -73,5 +79,10 @@ sigmaRadialCoeffs } +// Do front and back need to be merged? Usually only makes sense for 360 +// degree wedges. +mergeFaces true; + + // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // diff --git a/applications/utilities/mesh/manipulation/createPatch/createPatch.C b/applications/utilities/mesh/manipulation/createPatch/createPatch.C index 7d74e7208c..40ba0efade 100644 --- a/applications/utilities/mesh/manipulation/createPatch/createPatch.C +++ b/applications/utilities/mesh/manipulation/createPatch/createPatch.C @@ -208,8 +208,8 @@ void dumpCyclicMatch(const fileName& prefix, const polyMesh& mesh) // Dump halves { OFstream str(prefix+cycPatch.name()+"_half0.obj"); - Pout<< "Dumping cycPatch.name() half0 faces to " << str.name() - << endl; + Pout<< "Dumping " << cycPatch.name() + << " half0 faces to " << str.name() << endl; meshTools::writeOBJ ( str, @@ -226,8 +226,8 @@ void dumpCyclicMatch(const fileName& prefix, const polyMesh& mesh) } { OFstream str(prefix+cycPatch.name()+"_half1.obj"); - Pout<< "Dumping cycPatch.name() half1 faces to " << str.name() - << endl; + Pout<< "Dumping " << cycPatch.name() + << " half1 faces to " << str.name() << endl; meshTools::writeOBJ ( str, @@ -563,7 +563,7 @@ int main(int argc, char *argv[]) dumpCyclicMatch("initial_", mesh); // Read patch construct info from dictionary - PtrList patchSources(dict.lookup("patches")); + PtrList patchSources(dict.lookup("patchInfo")); diff --git a/applications/utilities/mesh/manipulation/createPatch/createPatchDict b/applications/utilities/mesh/manipulation/createPatch/createPatchDict index 9704392a10..b8676f134d 100644 --- a/applications/utilities/mesh/manipulation/createPatch/createPatchDict +++ b/applications/utilities/mesh/manipulation/createPatch/createPatchDict @@ -46,7 +46,7 @@ matchTolerance 1E-3; pointSync true; // Patches to create. -patches +patchInfo ( { // Name of new patch diff --git a/applications/utilities/mesh/manipulation/stitchMesh/stitchMesh.C b/applications/utilities/mesh/manipulation/stitchMesh/stitchMesh.C index 859ddf20a2..a5adac74e8 100644 --- a/applications/utilities/mesh/manipulation/stitchMesh/stitchMesh.C +++ b/applications/utilities/mesh/manipulation/stitchMesh/stitchMesh.C @@ -105,6 +105,8 @@ int main(int argc, char *argv[]) Foam::argList::validOptions.insert("overwrite", ""); + Foam::argList::validOptions.insert("toleranceDict", "file with tolerances"); + # include "setRootCase.H" # include "createTime.H" runTime.functionObjects().off(); @@ -168,6 +170,22 @@ int main(int argc, char *argv[]) << "If this is not the case use the -partial option" << nl << endl; } + // set up the tolerances for the sliding mesh + dictionary slidingTolerances; + if (args.options().found("toleranceDict")) + { + IOdictionary toleranceFile( + IOobject( + args.options()["toleranceDict"], + runTime.constant(), + mesh, + IOobject::MUST_READ, + IOobject::NO_WRITE + ) + ); + slidingTolerances += toleranceFile; + } + // Check for non-empty master and slave patches checkPatch(mesh.boundaryMesh(), masterPatchName); checkPatch(mesh.boundaryMesh(), slavePatchName); @@ -320,6 +338,11 @@ int main(int argc, char *argv[]) true // couple/decouple mode ) ); + static_cast(stitcher[0]).setTolerances + ( + slidingTolerances, + true + ); } diff --git a/applications/utilities/postProcessing/dataConversion/foamToEnsight/Make/options b/applications/utilities/postProcessing/dataConversion/foamToEnsight/Make/options index 823dae1ba7..01bd39cf76 100644 --- a/applications/utilities/postProcessing/dataConversion/foamToEnsight/Make/options +++ b/applications/utilities/postProcessing/dataConversion/foamToEnsight/Make/options @@ -5,5 +5,6 @@ EXE_INC = \ EXE_LIBS = \ -lfiniteVolume \ + -lgenericPatchFields \ -llagrangian diff --git a/applications/utilities/postProcessing/dataConversion/foamToEnsightParts/Make/options b/applications/utilities/postProcessing/dataConversion/foamToEnsightParts/Make/options index 518c2a0099..9243ffdb03 100644 --- a/applications/utilities/postProcessing/dataConversion/foamToEnsightParts/Make/options +++ b/applications/utilities/postProcessing/dataConversion/foamToEnsightParts/Make/options @@ -6,4 +6,5 @@ EXE_INC = \ EXE_LIBS = \ -lfiniteVolume \ -llagrangian \ + -lgenericPatchFields \ -lconversion diff --git a/applications/utilities/postProcessing/dataConversion/foamToFieldview9/Make/options b/applications/utilities/postProcessing/dataConversion/foamToFieldview9/Make/options index 1adeeefa12..a61c912ba0 100644 --- a/applications/utilities/postProcessing/dataConversion/foamToFieldview9/Make/options +++ b/applications/utilities/postProcessing/dataConversion/foamToFieldview9/Make/options @@ -4,4 +4,5 @@ EXE_INC = \ EXE_LIBS = \ -lfiniteVolume \ + -lgenericPatchFields \ -llagrangian diff --git a/applications/utilities/postProcessing/dataConversion/foamToGMV/Make/options b/applications/utilities/postProcessing/dataConversion/foamToGMV/Make/options index 3aab32c268..db1117a9b6 100644 --- a/applications/utilities/postProcessing/dataConversion/foamToGMV/Make/options +++ b/applications/utilities/postProcessing/dataConversion/foamToGMV/Make/options @@ -14,4 +14,5 @@ EXE_INC = \ EXE_LIBS = \ -lfiniteVolume \ + -lgenericPatchFields \ -llagrangian diff --git a/applications/utilities/postProcessing/dataConversion/foamToVTK/Make/options b/applications/utilities/postProcessing/dataConversion/foamToVTK/Make/options index 7d09e7e83f..2b32e21ed0 100644 --- a/applications/utilities/postProcessing/dataConversion/foamToVTK/Make/options +++ b/applications/utilities/postProcessing/dataConversion/foamToVTK/Make/options @@ -6,5 +6,6 @@ EXE_INC = \ EXE_LIBS = \ -lfiniteVolume \ -llagrangian \ + -lgenericPatchFields \ -lmeshTools diff --git a/applications/utilities/postProcessing/dataConversion/smapToFoam/Make/options b/applications/utilities/postProcessing/dataConversion/smapToFoam/Make/options index fa15f12452..44268cc287 100644 --- a/applications/utilities/postProcessing/dataConversion/smapToFoam/Make/options +++ b/applications/utilities/postProcessing/dataConversion/smapToFoam/Make/options @@ -2,4 +2,5 @@ EXE_INC = \ -I$(LIB_SRC)/finiteVolume/lnInclude EXE_LIBS = \ + -lgenericPatchFields \ -lfiniteVolume diff --git a/applications/utilities/postProcessing/sampling/sample/sampleDict b/applications/utilities/postProcessing/sampling/sample/sampleDict index f3659e728a..9c288d7cbf 100644 --- a/applications/utilities/postProcessing/sampling/sample/sampleDict +++ b/applications/utilities/postProcessing/sampling/sample/sampleDict @@ -22,6 +22,7 @@ FoamFile // jplot // gnuplot // raw +// vtk setFormat raw; // Surface output format. Choice of @@ -62,6 +63,7 @@ fields // curve specified points, not nessecary on line, uses // tracking // cloud specified points, uses findCell +// triSurfaceMeshPointSet points of triSurface // // axis: how to write point coordinate. Choice of // - x/y/z: x/y/z coordinate only diff --git a/bin/foamClearPolyMesh b/bin/foamClearPolyMesh index 1af677459d..ce3e563d31 100755 --- a/bin/foamClearPolyMesh +++ b/bin/foamClearPolyMesh @@ -130,6 +130,22 @@ for i in \ pointLevel \ refinementHistory \ surfaceIndex \ + points.gz \ + faces.gz \ + owner.gz \ + neighbour.gz \ + cells.gz \ + boundary.gz \ + pointZones.gz \ + faceZones.gz \ + cellZones.gz \ + meshModifiers.gz \ + parallelData.gz \ + sets.gz \ + cellLevel.gz \ + pointLevel.gz \ + refinementHistory.gz \ + surfaceIndex.gz \ ; do rm -rf $meshDir/$i diff --git a/bin/tools/CleanFunctions b/bin/tools/CleanFunctions index 4a2015febf..5e0ea64c5e 100644 --- a/bin/tools/CleanFunctions +++ b/bin/tools/CleanFunctions @@ -72,6 +72,7 @@ cleanCase () constant/{cellToRegion,cellLevel*,pointLevel*} \ constant/polyMesh/sets/ \ VTK \ + sets/streamLines \ > /dev/null 2>&1 for f in `find . -name "*Dict"` diff --git a/src/OpenFOAM/db/IOstreams/Pstreams/IPstream.H b/src/OpenFOAM/db/IOstreams/Pstreams/IPstream.H index 09db5931b8..db40ccfd4f 100644 --- a/src/OpenFOAM/db/IOstreams/Pstreams/IPstream.H +++ b/src/OpenFOAM/db/IOstreams/Pstreams/IPstream.H @@ -117,12 +117,6 @@ public: const std::streamsize bufSize ); - //- Non-blocking receives: wait until all have finished. - static void waitRequests(); - - //- Non-blocking receives: has request i finished? - static bool finishedRequest(const label i); - //- Return next token from stream Istream& read(token&); diff --git a/src/OpenFOAM/db/IOstreams/Pstreams/OPstream.H b/src/OpenFOAM/db/IOstreams/Pstreams/OPstream.H index 6db9a827b1..bccb5b08ef 100644 --- a/src/OpenFOAM/db/IOstreams/Pstreams/OPstream.H +++ b/src/OpenFOAM/db/IOstreams/Pstreams/OPstream.H @@ -115,12 +115,6 @@ public: const std::streamsize bufSize ); - //- Non-blocking writes: wait until all have finished. - static void waitRequests(); - - //- Non-blocking writes: has request i finished? - static bool finishedRequest(const label i); - //- Write next token to stream Ostream& write(const token&); diff --git a/src/OpenFOAM/db/IOstreams/Pstreams/Pstream.H b/src/OpenFOAM/db/IOstreams/Pstreams/Pstream.H index 5e8e5cbbb7..f54f02b271 100644 --- a/src/OpenFOAM/db/IOstreams/Pstreams/Pstream.H +++ b/src/OpenFOAM/db/IOstreams/Pstreams/Pstream.H @@ -264,6 +264,12 @@ public: // Spawns slave processes and initialises inter-communication static bool init(int& argc, char**& argv); + //- Non-blocking comms: wait until all have finished. + static void waitRequests(); + + //- Non-blocking comms: has request i finished? + static bool finishedRequest(const label i); + //- Is this a parallel run? static bool parRun() { diff --git a/src/OpenFOAM/fields/GeometricFields/GeometricField/GeometricBoundaryField.C b/src/OpenFOAM/fields/GeometricFields/GeometricField/GeometricBoundaryField.C index b3b6bb1d72..b2aa0d93e3 100644 --- a/src/OpenFOAM/fields/GeometricFields/GeometricField/GeometricBoundaryField.C +++ b/src/OpenFOAM/fields/GeometricFields/GeometricField/GeometricBoundaryField.C @@ -300,8 +300,7 @@ evaluate() // Block for any outstanding requests if (Pstream::defaultCommsType == Pstream::nonBlocking) { - IPstream::waitRequests(); - OPstream::waitRequests(); + Pstream::waitRequests(); } forAll(*this, patchi) diff --git a/src/OpenFOAM/meshes/polyMesh/mapPolyMesh/mapDistribute/mapDistribute.C b/src/OpenFOAM/meshes/polyMesh/mapPolyMesh/mapDistribute/mapDistribute.C index e2c2dbf773..0c524920af 100644 --- a/src/OpenFOAM/meshes/polyMesh/mapPolyMesh/mapDistribute/mapDistribute.C +++ b/src/OpenFOAM/meshes/polyMesh/mapPolyMesh/mapDistribute/mapDistribute.C @@ -339,8 +339,7 @@ void Foam::mapDistribute::compact(const boolList& elemIsUsed) // Wait for all to finish - OPstream::waitRequests(); - IPstream::waitRequests(); + Pstream::waitRequests(); // Compact out all submap entries that are referring to unused elements diff --git a/src/OpenFOAM/meshes/polyMesh/mapPolyMesh/mapDistribute/mapDistributeTemplates.C b/src/OpenFOAM/meshes/polyMesh/mapPolyMesh/mapDistribute/mapDistributeTemplates.C index da1ed19d1e..bd3e6b744e 100644 --- a/src/OpenFOAM/meshes/polyMesh/mapPolyMesh/mapDistribute/mapDistributeTemplates.C +++ b/src/OpenFOAM/meshes/polyMesh/mapPolyMesh/mapDistribute/mapDistributeTemplates.C @@ -287,8 +287,7 @@ void Foam::mapDistribute::distribute // Wait till all finished - IPstream::waitRequests(); - OPstream::waitRequests(); + Pstream::waitRequests(); // Consume for (label domain = 0; domain < Pstream::nProcs(); domain++) @@ -413,8 +412,7 @@ void Foam::mapDistribute::distribute // Wait for all to finish - OPstream::waitRequests(); - IPstream::waitRequests(); + Pstream::waitRequests(); // Collect neighbour fields @@ -717,8 +715,7 @@ void Foam::mapDistribute::distribute // Wait till all finished - IPstream::waitRequests(); - OPstream::waitRequests(); + Pstream::waitRequests(); // Consume for (label domain = 0; domain < Pstream::nProcs(); domain++) @@ -842,8 +839,7 @@ void Foam::mapDistribute::distribute // Wait for all to finish - OPstream::waitRequests(); - IPstream::waitRequests(); + Pstream::waitRequests(); // Collect neighbour fields diff --git a/src/OpenFOAM/meshes/polyMesh/polyPatches/constraint/cyclic/cyclicPolyPatch.C b/src/OpenFOAM/meshes/polyMesh/polyPatches/constraint/cyclic/cyclicPolyPatch.C index 0b9c08b39e..265b61cc4b 100644 --- a/src/OpenFOAM/meshes/polyMesh/polyPatches/constraint/cyclic/cyclicPolyPatch.C +++ b/src/OpenFOAM/meshes/polyMesh/polyPatches/constraint/cyclic/cyclicPolyPatch.C @@ -1305,7 +1305,7 @@ bool Foam::cyclicPolyPatch::order { label baffleI = 0; - forAll(*this, faceI) + forAll(pp, faceI) { const face& f = pp.localFaces()[faceI]; const labelList& pFaces = pp.pointFaces()[f[0]]; diff --git a/src/Pstream/dummy/IPread.C b/src/Pstream/dummy/IPread.C index 6ab24bc152..8e3d64982b 100644 --- a/src/Pstream/dummy/IPread.C +++ b/src/Pstream/dummy/IPread.C @@ -86,17 +86,4 @@ int Foam::IPstream::read } -void Foam::IPstream::waitRequests() -{} - - -bool Foam::IPstream::finishedRequest(const label) -{ - notImplemented("IPstream::finishedRequest()"); - return false; -} - - -// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // - // ************************************************************************* // diff --git a/src/Pstream/dummy/OPwrite.C b/src/Pstream/dummy/OPwrite.C index 1b7326bdd8..ad12d738fd 100644 --- a/src/Pstream/dummy/OPwrite.C +++ b/src/Pstream/dummy/OPwrite.C @@ -65,17 +65,4 @@ bool Foam::OPstream::write } -void Foam::OPstream::waitRequests() -{} - - -bool Foam::OPstream::finishedRequest(const label) -{ - notImplemented("OPstream::finishedRequest()"); - return false; -} - - -// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // - // ************************************************************************* // diff --git a/src/Pstream/dummy/Pstream.C b/src/Pstream/dummy/Pstream.C index 2c3825eac4..bd22701cfe 100644 --- a/src/Pstream/dummy/Pstream.C +++ b/src/Pstream/dummy/Pstream.C @@ -61,6 +61,17 @@ void Foam::Pstream::abort() void Foam::reduce(scalar&, const sumOp&) {} -// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + + +void Foam::Pstream::waitRequests() +{} + + +bool Foam::Pstream::finishedRequest(const label i) +{ + notImplemented("Pstream::finishedRequest()"); + return false; +} + // ************************************************************************* // diff --git a/src/Pstream/mpi/IPread.C b/src/Pstream/mpi/IPread.C index c08c6b56ab..966ae7ddc9 100644 --- a/src/Pstream/mpi/IPread.C +++ b/src/Pstream/mpi/IPread.C @@ -30,6 +30,7 @@ Description #include "mpi.h" #include "IPstream.H" +#include "PstreamGlobals.H" // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // @@ -37,7 +38,7 @@ Description // Outstanding non-blocking operations. //! @cond fileScope -Foam::DynamicList IPstream_outstandingRequests_; +//Foam::DynamicList IPstream_outstandingRequests_; //! @endcond fileScope // * * * * * * * * * * * * * * * * Constructor * * * * * * * * * * * * * * * // @@ -185,7 +186,7 @@ Foam::label Foam::IPstream::read return 0; } - IPstream_outstandingRequests_.append(request); + PstreamGlobals::outstandingRequests_.append(request); // Assume the message is completely received. return bufSize; @@ -204,52 +205,6 @@ Foam::label Foam::IPstream::read } -void Foam::IPstream::waitRequests() -{ - if (IPstream_outstandingRequests_.size()) - { - if - ( - MPI_Waitall - ( - IPstream_outstandingRequests_.size(), - IPstream_outstandingRequests_.begin(), - MPI_STATUSES_IGNORE - ) - ) - { - FatalErrorIn - ( - "IPstream::waitRequests()" - ) << "MPI_Waitall returned with error" << endl; - } - - IPstream_outstandingRequests_.clear(); - } -} - - -bool Foam::IPstream::finishedRequest(const label i) -{ - if (i >= IPstream_outstandingRequests_.size()) - { - FatalErrorIn - ( - "IPstream::finishedRequest(const label)" - ) << "There are " << IPstream_outstandingRequests_.size() - << " outstanding send requests and you are asking for i=" << i - << nl - << "Maybe you are mixing blocking/non-blocking comms?" - << Foam::abort(FatalError); - } - - int flag; - MPI_Test(&IPstream_outstandingRequests_[i], &flag, MPI_STATUS_IGNORE); - - return flag != 0; -} - - // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // // ************************************************************************* // diff --git a/src/Pstream/mpi/Make/files b/src/Pstream/mpi/Make/files index 65016b83f8..6b1d3f4c57 100644 --- a/src/Pstream/mpi/Make/files +++ b/src/Pstream/mpi/Make/files @@ -1,5 +1,6 @@ OPwrite.C IPread.C Pstream.C +PstreamGlobals.C LIB = $(FOAM_MPI_LIBBIN)/libPstream diff --git a/src/Pstream/mpi/OPwrite.C b/src/Pstream/mpi/OPwrite.C index 4836b92c8b..1e963ec13f 100644 --- a/src/Pstream/mpi/OPwrite.C +++ b/src/Pstream/mpi/OPwrite.C @@ -30,13 +30,7 @@ Description #include "mpi.h" #include "OPstream.H" - -// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // - -// Outstanding non-blocking operations. -//! @cond fileScope -Foam::DynamicList OPstream_outstandingRequests_; -//! @endcond fileScope +#include "PstreamGlobals.H" // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * // @@ -126,7 +120,7 @@ bool Foam::OPstream::write &request ); - OPstream_outstandingRequests_.append(request); + PstreamGlobals::outstandingRequests_.append(request); } else { @@ -142,52 +136,6 @@ bool Foam::OPstream::write } -void Foam::OPstream::waitRequests() -{ - if (OPstream_outstandingRequests_.size()) - { - if - ( - MPI_Waitall - ( - OPstream_outstandingRequests_.size(), - OPstream_outstandingRequests_.begin(), - MPI_STATUSES_IGNORE - ) - ) - { - FatalErrorIn - ( - "OPstream::waitRequests()" - ) << "MPI_Waitall returned with error" << Foam::endl; - } - - OPstream_outstandingRequests_.clear(); - } -} - - -bool Foam::OPstream::finishedRequest(const label i) -{ - if (i >= OPstream_outstandingRequests_.size()) - { - FatalErrorIn - ( - "OPstream::finishedRequest(const label)" - ) << "There are " << OPstream_outstandingRequests_.size() - << " outstanding send requests and you are asking for i=" << i - << nl - << "Maybe you are mixing blocking/non-blocking comms?" - << Foam::abort(FatalError); - } - - int flag; - MPI_Test(&OPstream_outstandingRequests_[i], &flag, MPI_STATUS_IGNORE); - - return flag != 0; -} - - // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // // ************************************************************************* // diff --git a/src/Pstream/mpi/Pstream.C b/src/Pstream/mpi/Pstream.C index ba4629d754..a535f4804e 100644 --- a/src/Pstream/mpi/Pstream.C +++ b/src/Pstream/mpi/Pstream.C @@ -29,6 +29,7 @@ License #include "Pstream.H" #include "PstreamReduceOps.H" #include "OSspecific.H" +#include "PstreamGlobals.H" #include #include @@ -130,6 +131,19 @@ void Foam::Pstream::exit(int errnum) delete[] buff; # endif + if (PstreamGlobals::outstandingRequests_.size()) + { + label n = PstreamGlobals::outstandingRequests_.size(); + PstreamGlobals::outstandingRequests_.clear(); + + WarningIn("Pstream::exit(int)") + << "There are still " << n << " outstanding MPI_Requests." << endl + << "This means that your code exited before doing a" + << " Pstream::waitRequests()." << endl + << "This should not happen for a normal code exit." + << endl; + } + if (errnum == 0) { MPI_Finalize(); @@ -422,6 +436,57 @@ void Foam::reduce(scalar& Value, const sumOp& bop) } +void Foam::Pstream::waitRequests() +{ + if (PstreamGlobals::outstandingRequests_.size()) + { + if + ( + MPI_Waitall + ( + PstreamGlobals::outstandingRequests_.size(), + PstreamGlobals::outstandingRequests_.begin(), + MPI_STATUSES_IGNORE + ) + ) + { + FatalErrorIn + ( + "Pstream::waitRequests()" + ) << "MPI_Waitall returned with error" << Foam::endl; + } + + PstreamGlobals::outstandingRequests_.clear(); + } +} + + +bool Foam::Pstream::finishedRequest(const label i) +{ + if (i >= PstreamGlobals::outstandingRequests_.size()) + { + FatalErrorIn + ( + "Pstream::finishedRequest(const label)" + ) << "There are " << PstreamGlobals::outstandingRequests_.size() + << " outstanding send requests and you are asking for i=" << i + << nl + << "Maybe you are mixing blocking/non-blocking comms?" + << Foam::abort(FatalError); + } + + int flag; + MPI_Test + ( + &PstreamGlobals::outstandingRequests_[i], + &flag, + MPI_STATUS_IGNORE + ); + + return flag != 0; +} + + // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // // ************************************************************************* // diff --git a/src/Pstream/mpi/PstreamGlobals.C b/src/Pstream/mpi/PstreamGlobals.C new file mode 100644 index 0000000000..f9271e5fe4 --- /dev/null +++ b/src/Pstream/mpi/PstreamGlobals.C @@ -0,0 +1,45 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 1991-2009 OpenCFD Ltd. + \\/ M anipulation | +------------------------------------------------------------------------------- +License + This file is part of OpenFOAM. + + OpenFOAM is free software; you can redistribute it and/or modify it + under the terms of the GNU General Public License as published by the + Free Software Foundation; either version 2 of the License, or (at your + option) any later version. + + OpenFOAM is distributed in the hope that it will be useful, but WITHOUT + ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + for more details. + + You should have received a copy of the GNU General Public License + along with OpenFOAM; if not, write to the Free Software Foundation, + Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA + +\*---------------------------------------------------------------------------*/ + +#include "PstreamGlobals.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ + +// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // + +// Outstanding non-blocking operations. +//! @cond fileScope +DynamicList PstreamGlobals::outstandingRequests_; +//! @endcond fileScope + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace Foam + +// ************************************************************************* // diff --git a/src/Pstream/mpi/PstreamGlobals.H b/src/Pstream/mpi/PstreamGlobals.H new file mode 100644 index 0000000000..07aa35a48f --- /dev/null +++ b/src/Pstream/mpi/PstreamGlobals.H @@ -0,0 +1,69 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 1991-2009 OpenCFD Ltd. + \\/ M anipulation | +------------------------------------------------------------------------------- +License + This file is part of OpenFOAM. + + OpenFOAM is free software; you can redistribute it and/or modify it + under the terms of the GNU General Public License as published by the + Free Software Foundation; either version 2 of the License, or (at your + option) any later version. + + OpenFOAM is distributed in the hope that it will be useful, but WITHOUT + ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + for more details. + + You should have received a copy of the GNU General Public License + along with OpenFOAM; if not, write to the Free Software Foundation, + Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA + +Namespace + Foam::PstreamGlobals + +Description + Global functions and variables for working with parallel streams, + but principally for gamma/mpi + +SourceFiles + PstreamGlobals.C + +\*---------------------------------------------------------------------------*/ + +#ifndef PstreamGlobals_H +#define PstreamGlobals_H + +#include "mpi.h" + +#include "DynamicList.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ + +/*---------------------------------------------------------------------------*\ + Class PstreamGlobals Declaration +\*---------------------------------------------------------------------------*/ + +namespace PstreamGlobals +{ + +extern DynamicList outstandingRequests_; + +}; + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace Foam + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#endif + +// ************************************************************************* // diff --git a/src/autoMesh/autoHexMesh/autoHexMeshDriver/autoLayerDriver.C b/src/autoMesh/autoHexMesh/autoHexMeshDriver/autoLayerDriver.C index fea25e35eb..716bdc6ddf 100644 --- a/src/autoMesh/autoHexMesh/autoHexMeshDriver/autoLayerDriver.C +++ b/src/autoMesh/autoHexMesh/autoHexMeshDriver/autoLayerDriver.C @@ -2987,6 +2987,7 @@ void Foam::autoLayerDriver::addLayers ( invExpansionRatio, pp, + labelList(0), // exposed patchIDs, not used for adding layers nPatchFaceLayers, // layers per face nPatchPointLayers, // layers per point firstDisp, // thickness of layer nearest internal mesh diff --git a/src/dynamicMesh/perfectInterface/perfectInterface.C b/src/dynamicMesh/perfectInterface/perfectInterface.C index 80c9ed6583..ee892b7e87 100644 --- a/src/dynamicMesh/perfectInterface/perfectInterface.C +++ b/src/dynamicMesh/perfectInterface/perfectInterface.C @@ -38,6 +38,7 @@ Description #include "polyModifyFace.H" #include "polyRemovePoint.H" #include "polyRemoveFace.H" +#include "indirectPrimitivePatch.H" // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // @@ -61,7 +62,7 @@ const Foam::scalar Foam::perfectInterface::tol_ = 1E-3; Foam::pointField Foam::perfectInterface::calcFaceCentres ( - const primitivePatch& pp + const indirectPrimitivePatch& pp ) { const pointField& points = pp.points(); @@ -155,6 +156,295 @@ bool Foam::perfectInterface::changeTopology() const } +void Foam::perfectInterface::setRefinement +( + const indirectPrimitivePatch& pp0, + const indirectPrimitivePatch& pp1, + polyTopoChange& ref +) const +{ + const polyMesh& mesh = topoChanger().mesh(); + + const polyBoundaryMesh& patches = mesh.boundaryMesh(); + + // Some aliases + const edgeList& edges0 = pp0.edges(); + const pointField& pts0 = pp0.localPoints(); + const pointField& pts1 = pp1.localPoints(); + const labelList& meshPts0 = pp0.meshPoints(); + const labelList& meshPts1 = pp1.meshPoints(); + + + // Get local dimension as fraction of minimum edge length + + scalar minLen = GREAT; + + forAll(edges0, edgeI) + { + minLen = min(minLen, edges0[edgeI].mag(pts0)); + } + scalar typDim = tol_*minLen; + + if (debug) + { + Pout<< "typDim:" << typDim << " edges0:" << edges0.size() + << " pts0:" << pts0.size() << " pts1:" << pts1.size() + << " pp0:" << pp0.size() << " pp1:" << pp1.size() << endl; + } + + + // Determine pointMapping in mesh point labels. Uses geometric + // comparison to find correspondence between patch points. + + labelList renumberPoints(mesh.points().size()); + forAll(renumberPoints, i) + { + renumberPoints[i] = i; + } + { + labelList from1To0Points(pts1.size()); + + bool matchOk = matchPoints + ( + pts1, + pts0, + scalarField(pts1.size(), typDim), // tolerance + true, // verbose + from1To0Points + ); + + if (!matchOk) + { + FatalErrorIn + ( + "perfectInterface::setRefinement(polyTopoChange& ref) const" + ) << "Points on patch sides do not match to within tolerance " + << typDim << exit(FatalError); + } + + forAll(pts1, i) + { + renumberPoints[meshPts1[i]] = meshPts0[from1To0Points[i]]; + } + } + + + + // Calculate correspondence between patch faces + + labelList from0To1Faces(pp1.size()); + + bool matchOk = matchPoints + ( + calcFaceCentres(pp0), + calcFaceCentres(pp1), + scalarField(pp0.size(), typDim), // tolerance + true, // verbose + from0To1Faces + ); + + if (!matchOk) + { + FatalErrorIn + ( + "perfectInterface::setRefinement(polyTopoChange& ref) const" + ) << "Face centres of patch sides do not match to within tolerance " + << typDim << exit(FatalError); + } + + + + // Now + // - renumber faces using pts1 (except patch1 faces) + // - remove patch1 faces. Remember cell label on owner side. + // - modify patch0 faces to be internal. + + // 1. Get faces to be renumbered + labelHashSet affectedFaces(2*pp1.size()); + forAll(meshPts1, i) + { + label meshPointI = meshPts1[i]; + + if (meshPointI != renumberPoints[meshPointI]) + { + const labelList& pFaces = mesh.pointFaces()[meshPointI]; + + forAll(pFaces, pFaceI) + { + affectedFaces.insert(pFaces[pFaceI]); + } + } + } + forAll(pp1, i) + { + affectedFaces.erase(pp1.addressing()[i]); + } + // Remove patch0 from renumbered faces. Should not be nessecary since + // patch0 and 1 should not share any point (if created by mergeMeshing) + // so affectedFaces should not contain any patch0 faces but you can + // never be sure what the user is doing. + forAll(pp0, i) + { + label faceI = pp0.addressing()[i]; + + if (affectedFaces.erase(faceI)) + { + WarningIn + ( + "perfectInterface::setRefinement(polyTopoChange&) const" + ) << "Found face " << faceI << " vertices " + << mesh.faces()[faceI] << " whose points are" + << " used both by master patch and slave patch" << endl; + } + } + + + // 2. Renumber (non patch0/1) faces. + for + ( + labelHashSet::const_iterator iter = affectedFaces.begin(); + iter != affectedFaces.end(); + ++iter + ) + { + label faceI = iter.key(); + + const face& f = mesh.faces()[faceI]; + + face newFace(f.size()); + + forAll(newFace, fp) + { + newFace[fp] = renumberPoints[f[fp]]; + } + + label nbr = -1; + + label patchI = -1; + + if (mesh.isInternalFace(faceI)) + { + nbr = mesh.faceNeighbour()[faceI]; + } + else + { + patchI = patches.whichPatch(faceI); + } + + label zoneID = mesh.faceZones().whichZone(faceI); + + bool zoneFlip = false; + + if (zoneID >= 0) + { + const faceZone& fZone = mesh.faceZones()[zoneID]; + + zoneFlip = fZone.flipMap()[fZone.whichFace(faceI)]; + } + + ref.setAction + ( + polyModifyFace + ( + newFace, // modified face + faceI, // label of face being modified + mesh.faceOwner()[faceI], // owner + nbr, // neighbour + false, // face flip + patchI, // patch for face + false, // remove from zone + zoneID, // zone for face + zoneFlip // face flip in zone + ) + ); + } + + + // 3. Remove patch1 points + forAll(meshPts1, i) + { + label meshPointI = meshPts1[i]; + + if (meshPointI != renumberPoints[meshPointI]) + { + ref.setAction(polyRemovePoint(meshPointI)); + } + } + + + // 4. Remove patch1 faces + forAll(pp1, i) + { + label faceI = pp1.addressing()[i]; + ref.setAction(polyRemoveFace(faceI)); + } + + + // 5. Modify patch0 faces for new points (not really nessecary; see + // comment above about patch1 and patch0 never sharing points) and + // becoming internal. + const boolList& mfFlip = + mesh.faceZones()[faceZoneID_.index()].flipMap(); + + forAll(pp0, i) + { + label faceI = pp0.addressing()[i]; + + const face& f = mesh.faces()[faceI]; + + face newFace(f.size()); + + forAll(newFace, fp) + { + newFace[fp] = renumberPoints[f[fp]]; + } + + label own = mesh.faceOwner()[faceI]; + + label pp1FaceI = pp1.addressing()[from0To1Faces[i]]; + + label nbr = mesh.faceOwner()[pp1FaceI]; + + if (own < nbr) + { + ref.setAction + ( + polyModifyFace + ( + newFace, // modified face + faceI, // label of face being modified + own, // owner + nbr, // neighbour + false, // face flip + -1, // patch for face + false, // remove from zone + faceZoneID_.index(), // zone for face + mfFlip[i] // face flip in zone + ) + ); + } + else + { + ref.setAction + ( + polyModifyFace + ( + newFace.reverseFace(), // modified face + faceI, // label of face being modified + nbr, // owner + own, // neighbour + true, // face flip + -1, // patch for face + false, // remove from zone + faceZoneID_.index(), // zone for face + !mfFlip[i] // face flip in zone + ) + ); + } + } +} + + void Foam::perfectInterface::setRefinement(polyTopoChange& ref) const { if (debug) @@ -176,286 +466,25 @@ void Foam::perfectInterface::setRefinement(polyTopoChange& ref) const const polyMesh& mesh = topoChanger().mesh(); const polyBoundaryMesh& patches = mesh.boundaryMesh(); - - const polyPatch& pp0 = patches[masterPatchID_.index()]; - const polyPatch& pp1 = patches[slavePatchID_.index()]; - - // Some aliases - const edgeList& edges0 = pp0.edges(); - const pointField& pts0 = pp0.localPoints(); - const pointField& pts1 = pp1.localPoints(); - const labelList& meshPts0 = pp0.meshPoints(); - const labelList& meshPts1 = pp1.meshPoints(); - - - // Get local dimension as fraction of minimum edge length - - scalar minLen = GREAT; - - forAll(edges0, edgeI) - { - minLen = min(minLen, edges0[edgeI].mag(pts0)); - } - scalar typDim = tol_*minLen; - - if (debug) - { - Pout<< "typDim:" << typDim << " edges0:" << edges0.size() - << " pts0:" << pts0.size() << " pts1:" << pts1.size() - << " pp0:" << pp0.size() << " pp1:" << pp1.size() << endl; - } + const polyPatch& patch0 = patches[masterPatchID_.index()]; + const polyPatch& patch1 = patches[slavePatchID_.index()]; - // Determine pointMapping in mesh point labels. Uses geometric - // comparison to find correspondence between patch points. - - labelList renumberPoints(mesh.points().size()); - forAll(renumberPoints, i) - { - renumberPoints[i] = i; - } - { - labelList from1To0Points(pts1.size()); - - bool matchOk = matchPoints - ( - pts1, - pts0, - scalarField(pts1.size(), typDim), // tolerance - true, // verbose - from1To0Points - ); - - if (!matchOk) - { - FatalErrorIn - ( - "perfectInterface::setRefinement(polyTopoChange& ref) const" - ) << "Points on patches " << pp0.name() << " and " - << pp1.name() << " do not match to within tolerance " - << typDim << exit(FatalError); - } - - forAll(pts1, i) - { - renumberPoints[meshPts1[i]] = meshPts0[from1To0Points[i]]; - } - } - - - - // Calculate correspondence between patch faces - - labelList from0To1Faces(pp1.size()); - - bool matchOk = matchPoints + labelList pp0Labels(identity(patch0.size())+patch0.start()); + indirectPrimitivePatch pp0 ( - calcFaceCentres(pp0), - calcFaceCentres(pp1), - scalarField(pp0.size(), typDim), // tolerance - true, // verbose - from0To1Faces + IndirectList(mesh.faces(), pp0Labels), + mesh.points() ); - if (!matchOk) - { - FatalErrorIn - ( - "perfectInterface::setRefinement(polyTopoChange& ref) const" - ) << "Face centres of patches " << pp0.name() << " and " - << pp1.name() << " do not match to within tolerance " << typDim - << exit(FatalError); - } - - - - // Now - // - renumber faces using pts1 (except patch1 faces) - // - remove patch1 faces. Remember cell label on owner side. - // - modify patch0 faces to be internal. - - // 1. Get faces to be renumbered - labelHashSet affectedFaces(2*pp1.size()); - forAll(meshPts1, i) - { - label meshPointI = meshPts1[i]; - - if (meshPointI != renumberPoints[meshPointI]) - { - const labelList& pFaces = mesh.pointFaces()[meshPointI]; - - forAll(pFaces, pFaceI) - { - affectedFaces.insert(pFaces[pFaceI]); - } - } - } - forAll(pp1, i) - { - affectedFaces.erase(pp1.start() + i); - } - // Remove patch0 from renumbered faces. Should not be nessecary since - // patch0 and 1 should not share any point (if created by mergeMeshing) - // so affectedFaces should not contain any patch0 faces but you can - // never be sure what the user is doing. - forAll(pp0, i) - { - if (affectedFaces.erase(pp0.start() + i)) - { - WarningIn - ( - "perfectInterface::setRefinement(polyTopoChange&) const" - ) << "Found face " << pp0.start() + i << " vertices " - << mesh.faces()[pp0.start() + i] << " whose points are" - << " used both by master patch " << pp0.name() - << " and slave patch " << pp1.name() - << endl; - } - } - - - // 2. Renumber (non patch0/1) faces. - for + labelList pp1Labels(identity(patch1.size())+patch1.start()); + indirectPrimitivePatch pp1 ( - labelHashSet::const_iterator iter = affectedFaces.begin(); - iter != affectedFaces.end(); - ++iter - ) - { - label faceI = iter.key(); + IndirectList(mesh.faces(), pp1Labels), + mesh.points() + ); - const face& f = mesh.faces()[faceI]; - - face newFace(f.size()); - - forAll(newFace, fp) - { - newFace[fp] = renumberPoints[f[fp]]; - } - - label nbr = -1; - - label patchI = -1; - - if (mesh.isInternalFace(faceI)) - { - nbr = mesh.faceNeighbour()[faceI]; - } - else - { - patchI = patches.whichPatch(faceI); - } - - label zoneID = mesh.faceZones().whichZone(faceI); - - bool zoneFlip = false; - - if (zoneID >= 0) - { - const faceZone& fZone = mesh.faceZones()[zoneID]; - - zoneFlip = fZone.flipMap()[fZone.whichFace(faceI)]; - } - - ref.setAction - ( - polyModifyFace - ( - newFace, // modified face - faceI, // label of face being modified - mesh.faceOwner()[faceI], // owner - nbr, // neighbour - false, // face flip - patchI, // patch for face - false, // remove from zone - zoneID, // zone for face - zoneFlip // face flip in zone - ) - ); - } - - - // 3. Remove patch1 points - forAll(meshPts1, i) - { - label meshPointI = meshPts1[i]; - - if (meshPointI != renumberPoints[meshPointI]) - { - ref.setAction(polyRemovePoint(meshPointI)); - } - } - - - // 4. Remove patch1 faces - forAll(pp1, i) - { - ref.setAction(polyRemoveFace(pp1.start() + i)); - } - - - // 5. Modify patch0 faces for new points (not really nessecary; see - // comment above about patch1 and patch0 never sharing points) and - // becoming internal. - const boolList& mfFlip = - mesh.faceZones()[faceZoneID_.index()].flipMap(); - - forAll(pp0, i) - { - label faceI = pp0.start() + i; - - const face& f = mesh.faces()[faceI]; - - face newFace(f.size()); - - forAll(newFace, fp) - { - newFace[fp] = renumberPoints[f[fp]]; - } - - label own = mesh.faceOwner()[faceI]; - - label pp1FaceI = pp1.start() + from0To1Faces[i]; - - label nbr = mesh.faceOwner()[pp1FaceI]; - - if (own < nbr) - { - ref.setAction - ( - polyModifyFace - ( - newFace, // modified face - faceI, // label of face being modified - own, // owner - nbr, // neighbour - false, // face flip - -1, // patch for face - false, // remove from zone - faceZoneID_.index(), // zone for face - mfFlip[i] // face flip in zone - ) - ); - } - else - { - ref.setAction - ( - polyModifyFace - ( - newFace.reverseFace(), // modified face - faceI, // label of face being modified - nbr, // owner - own, // neighbour - false, // face flip - -1, // patch for face - false, // remove from zone - faceZoneID_.index(), // zone for face - !mfFlip[i] // face flip in zone - ) - ); - } - } + setRefinement(pp0, pp1, ref); } } diff --git a/src/dynamicMesh/perfectInterface/perfectInterface.H b/src/dynamicMesh/perfectInterface/perfectInterface.H index 3916bb3e74..dec4f6a48d 100644 --- a/src/dynamicMesh/perfectInterface/perfectInterface.H +++ b/src/dynamicMesh/perfectInterface/perfectInterface.H @@ -40,6 +40,7 @@ SourceFiles #include "polyMeshModifier.H" #include "polyPatchID.H" #include "ZoneIDs.H" +#include "indirectPrimitivePatch.H" // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // @@ -75,7 +76,7 @@ class perfectInterface // Private Member Functions //- Calculate face centres on patch - static pointField calcFaceCentres(const primitivePatch&); + static pointField calcFaceCentres(const indirectPrimitivePatch&); //- Disallow default bitwise copy construct @@ -128,6 +129,17 @@ public: // into the topological change virtual void setRefinement(polyTopoChange&) const; + //- Insert the layer addition/removal instructions + // into the topological change. Uses only mesh, not any of the + // patch and zone indices. Bit of a workaround - used in extruding + // a mesh. + virtual void setRefinement + ( + const indirectPrimitivePatch& pp0, + const indirectPrimitivePatch& pp1, + polyTopoChange& + ) const; + //- Modify motion points to comply with the topological change virtual void modifyMotionPoints(pointField& motionPoints) const; diff --git a/src/dynamicMesh/polyTopoChange/polyTopoChange/addPatchCellLayer.C b/src/dynamicMesh/polyTopoChange/polyTopoChange/addPatchCellLayer.C index 2f842e6d66..1a725f0ebf 100644 --- a/src/dynamicMesh/polyTopoChange/polyTopoChange/addPatchCellLayer.C +++ b/src/dynamicMesh/polyTopoChange/polyTopoChange/addPatchCellLayer.C @@ -335,17 +335,19 @@ Foam::label Foam::addPatchCellLayer::addSideFace label inflateEdgeI = -1; // Check mesh faces using edge - forAll(meshFaces, i) + if (addToMesh_) { - if (mesh_.isInternalFace(meshFaces[i])) + forAll(meshFaces, i) { - // meshEdge uses internal faces so ok to inflate from it - inflateEdgeI = meshEdgeI; - break; + if (mesh_.isInternalFace(meshFaces[i])) + { + // meshEdge uses internal faces so ok to inflate from it + inflateEdgeI = meshEdgeI; + break; + } } } - // Get my mesh face and its zone. label meshFaceI = pp.addressing()[ownFaceI]; label zoneI = mesh_.faceZones().whichZone(meshFaceI); @@ -504,9 +506,14 @@ Foam::label Foam::addPatchCellLayer::addSideFace // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // // Construct from mesh -Foam::addPatchCellLayer::addPatchCellLayer(const polyMesh& mesh) +Foam::addPatchCellLayer::addPatchCellLayer +( + const polyMesh& mesh, + const bool addToMesh +) : mesh_(mesh), + addToMesh_(addToMesh), addedPoints_(0), layerFaces_(0) {} @@ -551,6 +558,7 @@ void Foam::addPatchCellLayer::setRefinement ( const scalarField& expansionRatio, const indirectPrimitivePatch& pp, + const labelList& exposedPatchID, const labelList& nFaceLayers, const labelList& nPointLayers, const vectorField& firstLayerDisp, @@ -827,6 +835,19 @@ void Foam::addPatchCellLayer::setRefinement } + const polyBoundaryMesh& patches = mesh_.boundaryMesh(); + + // Precalculated patchID for each patch face + labelList patchID(pp.size()); + + forAll(pp, patchFaceI) + { + label meshFaceI = pp.addressing()[patchFaceI]; + + patchID[patchFaceI] = patches.whichPatch(meshFaceI); + } + + // From master point (in patch point label) to added points (in mesh point // label) addedPoints_.setSize(pp.nPoints()); @@ -857,6 +878,33 @@ void Foam::addPatchCellLayer::setRefinement // Create new points // + // If creating new mesh: copy existing patch points + labelList copiedPatchPoints; + if (!addToMesh_) + { + copiedPatchPoints.setSize(firstLayerDisp.size()); + forAll(firstLayerDisp, patchPointI) + { + if (addedPoints_[patchPointI].size()) + { + label meshPointI = meshPoints[patchPointI]; + label zoneI = mesh_.pointZones().whichZone(meshPointI); + copiedPatchPoints[patchPointI] = meshMod.setAction + ( + polyAddPoint + ( + mesh_.points()[meshPointI], // point + -1, // master point + zoneI, // zone for point + true // supports a cell + ) + ); + } + } + } + + + // Create points for additional layers forAll(firstLayerDisp, patchPointI) { if (addedPoints_[patchPointI].size()) @@ -878,7 +926,7 @@ void Foam::addPatchCellLayer::setRefinement polyAddPoint ( pt, // point - meshPointI, // master point + (addToMesh_ ? meshPointI : -1), // master point zoneI, // zone for point true // supports a cell ) @@ -922,34 +970,15 @@ void Foam::addPatchCellLayer::setRefinement -1, // master point -1, // master edge -1, // master face - mesh_.faceOwner()[meshFaceI], // master cell id + (addToMesh_ ? mesh_.faceOwner()[meshFaceI] : -1),//master ownZoneI // zone for cell ) ); - - //Pout<< "For patchFace:" << patchFaceI - // << " meshFace:" << pp.addressing()[patchFaceI] - // << " layer:" << i << " added cell:" - // << addedCells[patchFaceI][i] - // << endl; } } } - const polyBoundaryMesh& patches = mesh_.boundaryMesh(); - - // Precalculated patchID for each patch face - labelList patchID(pp.size()); - - forAll(pp, patchFaceI) - { - label meshFaceI = pp.addressing()[patchFaceI]; - - patchID[patchFaceI] = patches.whichPatch(meshFaceI); - } - - // Create faces on top of the original patch faces. // These faces are created from original patch faces outwards so in order @@ -965,7 +994,6 @@ void Foam::addPatchCellLayer::setRefinement if (addedCells[patchFaceI].size()) { layerFaces_[patchFaceI].setSize(addedCells[patchFaceI].size() + 1); - layerFaces_[patchFaceI][0] = meshFaceI; label zoneI = mesh_.faceZones().whichZone(meshFaceI); @@ -981,7 +1009,12 @@ void Foam::addPatchCellLayer::setRefinement if (addedPoints_[f[fp]].empty()) { // Keep original point - newFace[fp] = meshPoints[f[fp]]; + newFace[fp] = + ( + addToMesh_ + ? meshPoints[f[fp]] + : copiedPatchPoints[f[fp]] + ); } else { @@ -1021,18 +1054,13 @@ void Foam::addPatchCellLayer::setRefinement nei, // neighbour -1, // master point -1, // master edge - meshFaceI, // master face for addition + (addToMesh_ ? meshFaceI : -1), // master face false, // flux flip patchI, // patch for face zoneI, // zone for face false // face zone flip ) ); - //Pout<< "Added inbetween face " << newFace - // << " own:" << addedCells[patchFaceI][i] - // << " nei:" << nei - // << " patch:" << patchI - // << endl; } } } @@ -1040,35 +1068,81 @@ void Foam::addPatchCellLayer::setRefinement // // Modify old patch faces to be on the inside // - forAll(pp, patchFaceI) + + if (addToMesh_) { - if (addedCells[patchFaceI].size()) + forAll(pp, patchFaceI) { - label meshFaceI = pp.addressing()[patchFaceI]; + if (addedCells[patchFaceI].size()) + { + label meshFaceI = pp.addressing()[patchFaceI]; - label zoneI = mesh_.faceZones().whichZone(meshFaceI); + layerFaces_[patchFaceI][0] = meshFaceI; - meshMod.setAction - ( - polyModifyFace + label zoneI = mesh_.faceZones().whichZone(meshFaceI); + + meshMod.setAction ( - pp[patchFaceI], // modified face - meshFaceI, // label of face - mesh_.faceOwner()[meshFaceI], // owner - addedCells[patchFaceI][0], // neighbour - false, // face flip - -1, // patch for face - false, // remove from zone - zoneI, // zone for face - false // face flip in zone - ) - ); - //Pout<< "Modified old patch face " << meshFaceI - // << " own:" << mesh_.faceOwner()[meshFaceI] - // << " nei:" << addedCells[patchFaceI][0] - // << endl; + polyModifyFace + ( + pp[patchFaceI], // modified face + meshFaceI, // label of face + mesh_.faceOwner()[meshFaceI], // owner + addedCells[patchFaceI][0], // neighbour + false, // face flip + -1, // patch for face + false, // remove from zone + zoneI, // zone for face + false // face flip in zone + ) + ); + } } } + else + { + // If creating new mesh: reverse original faces and put them + // in the exposed patch ID. + forAll(pp, patchFaceI) + { + if (nFaceLayers[patchFaceI] > 0) + { + label meshFaceI = pp.addressing()[patchFaceI]; + label zoneI = mesh_.faceZones().whichZone(meshFaceI); + bool zoneFlip = false; + if (zoneI != -1) + { + const faceZone& fz = mesh_.faceZones()[zoneI]; + zoneFlip = !fz.flipMap()[fz.whichFace(meshFaceI)]; + } + + // Reverse and renumber old patch face. + face f(pp.localFaces()[patchFaceI].reverseFace()); + forAll(f, fp) + { + f[fp] = copiedPatchPoints[f[fp]]; + } + + layerFaces_[patchFaceI][0] = meshMod.setAction + ( + polyAddFace + ( + f, // modified face + addedCells[patchFaceI][0], // owner + -1, // neighbour + -1, // masterPoint + -1, // masterEdge + -1, // masterFace + true, // face flip + exposedPatchID[patchFaceI], // patch for face + zoneI, // zone for face + zoneFlip // face flip in zone + ) + ); + } + } + } + // @@ -1255,7 +1329,16 @@ void Foam::addPatchCellLayer::setRefinement forAll(stringedVerts, stringedI) { label v = stringedVerts[stringedI]; - addVertex(meshPoints[v], newFace, newFp); + addVertex + ( + ( + addToMesh_ + ? meshPoints[v] + : copiedPatchPoints[v] + ), + newFace, + newFp + ); } } else @@ -1276,7 +1359,16 @@ void Foam::addPatchCellLayer::setRefinement } else { - addVertex(meshPoints[v], newFace, newFp); + addVertex + ( + ( + addToMesh_ + ? meshPoints[v] + : copiedPatchPoints[v] + ), + newFace, + newFp + ); } } } @@ -1316,7 +1408,16 @@ void Foam::addPatchCellLayer::setRefinement } else { - addVertex(meshPoints[v], newFace, newFp); + addVertex + ( + ( + addToMesh_ + ? meshPoints[v] + : copiedPatchPoints[v] + ), + newFace, + newFp + ); } } diff --git a/src/dynamicMesh/polyTopoChange/polyTopoChange/addPatchCellLayer.H b/src/dynamicMesh/polyTopoChange/polyTopoChange/addPatchCellLayer.H index f03d23e214..f4b5bb2e79 100644 --- a/src/dynamicMesh/polyTopoChange/polyTopoChange/addPatchCellLayer.H +++ b/src/dynamicMesh/polyTopoChange/polyTopoChange/addPatchCellLayer.H @@ -26,7 +26,8 @@ Class Foam::addPatchCellLayer Description - Adds layers of cells to outside of polyPatch. + Adds layers of cells to outside of polyPatch. Can optionally create + stand-alone extruded mesh (addToMesh=false). Call setRefinement with offset vector for every patch point and number of layers per patch face and number of layers per patch point. @@ -164,6 +165,9 @@ class addPatchCellLayer //- Reference to mesh const polyMesh& mesh_; + //- Add layers to existing mesh or create new mesh + const bool addToMesh_; + //- For all patchpoints: list of added points (size 0 or nLayers) // First point in list is one nearest to original point in patch, // last one is the new point on the surface. @@ -253,7 +257,7 @@ public: // Constructors //- Construct from mesh. - addPatchCellLayer(const polyMesh& mesh); + addPatchCellLayer(const polyMesh& mesh, const bool addToMesh = true); // Member Functions @@ -291,8 +295,10 @@ public: //- Play commands into polyTopoChange to create layers on top // of indirectPrimitivePatch (have to be outside faces). // Gets displacement per patch point. - // - nPointLayers : number of layers per (patch)point - // - nFaceLayers : number of layers per (patch) face + // - exposedPatchID : only used if creating a new mesh (addToMesh=false) + // gives per pp face the patch the exposed face should get. + // - nPointLayers : number of layers per (patch)point. + // - nFaceLayers : number of layers per (patch) face. // - firstDisplacement : displacement per point for first // layer of points (i.e. nearest to original mesh). If zero // do not add point. @@ -305,14 +311,15 @@ public: // get a cell should firstDisplacement be <> 0 // Note: cells get added from owner cells of patch faces // (instead of e.g. from patch faces) - void setRefinement - ( - const scalarField& expansionRatio, - const indirectPrimitivePatch& pp, - const labelList& nFaceLayers, - const labelList& nPointLayers, - const vectorField& firstLayerDisp, - polyTopoChange& meshMod + void setRefinement + ( + const scalarField& expansionRatio, + const indirectPrimitivePatch& pp, + const labelList& exposedPatchID, + const labelList& nFaceLayers, + const labelList& nPointLayers, + const vectorField& firstLayerDisp, + polyTopoChange& meshMod ); @@ -329,6 +336,7 @@ public: ( scalarField(pp.nPoints(), 1.0), // expansion ration pp, + labelList(0), labelList(pp.size(), nLayers), labelList(pp.nPoints(), nLayers), overallDisplacement / nLayers, diff --git a/src/dynamicMesh/polyTopoChange/polyTopoChange/polyTopoChange.C b/src/dynamicMesh/polyTopoChange/polyTopoChange/polyTopoChange.C index f6c397640f..6b1874dc4b 100644 --- a/src/dynamicMesh/polyTopoChange/polyTopoChange/polyTopoChange.C +++ b/src/dynamicMesh/polyTopoChange/polyTopoChange/polyTopoChange.C @@ -749,6 +749,11 @@ void Foam::polyTopoChange::getFaceOrder " const" ) << "Did not determine new position" << " for face " << faceI + << " owner " << faceOwner_[faceI] + << " neighbour " << faceNeighbour_[faceI] + << " region " << region_[faceI] << endl + << "This is usually caused by not specifying a patch for" + << " a boundary face." << abort(FatalError); } } @@ -3176,7 +3181,7 @@ Foam::autoPtr Foam::polyTopoChange::makeMesh ( autoPtr& newMeshPtr, const IOobject& io, - const fvMesh& mesh, + const polyMesh& mesh, const bool syncParallel, const bool orderCells, const bool orderPoints diff --git a/src/dynamicMesh/polyTopoChange/polyTopoChange/polyTopoChange.H b/src/dynamicMesh/polyTopoChange/polyTopoChange/polyTopoChange.H index 77801e21f0..e7a0ba8fb3 100644 --- a/src/dynamicMesh/polyTopoChange/polyTopoChange/polyTopoChange.H +++ b/src/dynamicMesh/polyTopoChange/polyTopoChange/polyTopoChange.H @@ -585,7 +585,7 @@ public: ( autoPtr& newMesh, const IOobject& io, - const fvMesh& mesh, + const polyMesh& mesh, const bool syncParallel = true, const bool orderCells = false, const bool orderPoints = false diff --git a/src/dynamicMesh/slidingInterface/coupleSlidingInterface.C b/src/dynamicMesh/slidingInterface/coupleSlidingInterface.C index 18d1170c00..f90e9a036b 100644 --- a/src/dynamicMesh/slidingInterface/coupleSlidingInterface.C +++ b/src/dynamicMesh/slidingInterface/coupleSlidingInterface.C @@ -43,7 +43,7 @@ License // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // -const Foam::scalar Foam::slidingInterface::edgeCoPlanarTol_ = 0.8; +const Foam::scalar Foam::slidingInterface::edgeCoPlanarTolDefault_ = 0.8; // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * // diff --git a/src/dynamicMesh/slidingInterface/slidingInterface.C b/src/dynamicMesh/slidingInterface/slidingInterface.C index 80d68749e8..316dd389b6 100644 --- a/src/dynamicMesh/slidingInterface/slidingInterface.C +++ b/src/dynamicMesh/slidingInterface/slidingInterface.C @@ -27,10 +27,8 @@ License #include "slidingInterface.H" #include "polyTopoChanger.H" #include "polyMesh.H" -#include "primitiveMesh.H" #include "polyTopoChange.H" #include "addToRunTimeSelectionTable.H" -#include "triPointRef.H" #include "plane.H" // Index of debug signs: @@ -173,6 +171,14 @@ Foam::slidingInterface::slidingInterface attached_(false), projectionAlgo_(algo), trigger_(false), + pointMergeTol_(pointMergeTolDefault_), + edgeMergeTol_(edgeMergeTolDefault_), + nFacesPerSlaveEdge_(nFacesPerSlaveEdgeDefault_), + edgeFaceEscapeLimit_(edgeFaceEscapeLimitDefault_), + integralAdjTol_(integralAdjTolDefault_), + edgeMasterCatchFraction_(edgeMasterCatchFractionDefault_), + edgeCoPlanarTol_(edgeCoPlanarTolDefault_), + edgeEndCutoffTol_(edgeEndCutoffTolDefault_), cutFaceMasterPtr_(NULL), cutFaceSlavePtr_(NULL), masterFaceCellsPtr_(NULL), @@ -280,6 +286,9 @@ Foam::slidingInterface::slidingInterface masterPointEdgeHitsPtr_(NULL), projectedSlavePointsPtr_(NULL) { + // Optionally default tolerances from dictionary + setTolerances(dict); + checkDefinition(); // If the interface is attached, the master and slave face zone addressing @@ -686,6 +695,63 @@ const Foam::pointField& Foam::slidingInterface::pointProjection() const return *projectedSlavePointsPtr_; } +void Foam::slidingInterface::setTolerances(const dictionary&dict, bool report) +{ + pointMergeTol_ = dict.lookupOrDefault + ( + "pointMergeTol", + pointMergeTol_ + ); + edgeMergeTol_ = dict.lookupOrDefault + ( + "edgeMergeTol", + edgeMergeTol_ + ); + nFacesPerSlaveEdge_ = dict.lookupOrDefault + ( + "nFacesPerSlaveEdge", + nFacesPerSlaveEdge_ + ); + edgeFaceEscapeLimit_ = dict.lookupOrDefault + ( + "edgeFaceEscapeLimit", + edgeFaceEscapeLimit_ + ); + integralAdjTol_ = dict.lookupOrDefault + ( + "integralAdjTol", + integralAdjTol_ + ); + edgeMasterCatchFraction_ = dict.lookupOrDefault + ( + "edgeMasterCatchFraction", + edgeMasterCatchFraction_ + ); + edgeCoPlanarTol_ = dict.lookupOrDefault + ( + "edgeCoPlanarTol", + edgeCoPlanarTol_ + ); + edgeEndCutoffTol_ = dict.lookupOrDefault + ( + "edgeEndCutoffTol", + edgeEndCutoffTol_ + ); + + if (report) + { + Info<< "Sliding interface parameters:" << nl + << "pointMergeTol : " << pointMergeTol_ << nl + << "edgeMergeTol : " << edgeMergeTol_ << nl + << "nFacesPerSlaveEdge : " << nFacesPerSlaveEdge_ << nl + << "edgeFaceEscapeLimit : " << edgeFaceEscapeLimit_ << nl + << "integralAdjTol : " << integralAdjTol_ << nl + << "edgeMasterCatchFraction : " << edgeMasterCatchFraction_ << nl + << "edgeCoPlanarTol : " << edgeCoPlanarTol_ << nl + << "edgeEndCutoffTol : " << edgeEndCutoffTol_ << endl; + } +} + void Foam::slidingInterface::write(Ostream& os) const { @@ -703,6 +769,14 @@ void Foam::slidingInterface::write(Ostream& os) const } +// To write out all those tolerances +#define WRITE_NON_DEFAULT(name) \ + if( name ## _ != name ## Default_ )\ + { \ + os << " " #name " " << name ## _ << token::END_STATEMENT << nl; \ + } + + void Foam::slidingInterface::writeDict(Ostream& os) const { os << nl << name() << nl << token::BEGIN_BLOCK << nl @@ -743,6 +817,15 @@ void Foam::slidingInterface::writeDict(Ostream& os) const << token::END_STATEMENT << nl; } + WRITE_NON_DEFAULT(pointMergeTol) + WRITE_NON_DEFAULT(edgeMergeTol) + WRITE_NON_DEFAULT(nFacesPerSlaveEdge) + WRITE_NON_DEFAULT(edgeFaceEscapeLimit) + WRITE_NON_DEFAULT(integralAdjTol) + WRITE_NON_DEFAULT(edgeMasterCatchFraction) + WRITE_NON_DEFAULT(edgeCoPlanarTol) + WRITE_NON_DEFAULT(edgeEndCutoffTol) + os << token::END_BLOCK << endl; } diff --git a/src/dynamicMesh/slidingInterface/slidingInterface.H b/src/dynamicMesh/slidingInterface/slidingInterface.H index 376c969e38..f86e423c65 100644 --- a/src/dynamicMesh/slidingInterface/slidingInterface.H +++ b/src/dynamicMesh/slidingInterface/slidingInterface.H @@ -129,6 +129,32 @@ private: //- Trigger topological change mutable bool trigger_; + // Tolerances. Initialised to static ones below. + + //- Point merge tolerance + scalar pointMergeTol_; + + //- Edge merge tolerance + scalar edgeMergeTol_; + + //- Estimated number of faces an edge goes through + label nFacesPerSlaveEdge_; + + //- Edge-face interaction escape limit + label edgeFaceEscapeLimit_; + + //- Integral match point adjustment tolerance + scalar integralAdjTol_; + + //- Edge intersection master catch fraction + scalar edgeMasterCatchFraction_; + + //- Edge intersection co-planar tolerance + scalar edgeCoPlanarTol_; + + //- Edge end cut-off tolerance + scalar edgeEndCutoffTol_; + // Private addressing data. @@ -256,28 +282,28 @@ private: // Static data members //- Point merge tolerance - static const scalar pointMergeTol_; + static const scalar pointMergeTolDefault_; //- Edge merge tolerance - static const scalar edgeMergeTol_; + static const scalar edgeMergeTolDefault_; //- Estimated number of faces an edge goes through - static const label nFacesPerSlaveEdge_; + static const label nFacesPerSlaveEdgeDefault_; //- Edge-face interaction escape limit - static const label edgeFaceEscapeLimit_; + static const label edgeFaceEscapeLimitDefault_; //- Integral match point adjustment tolerance - static const scalar integralAdjTol_; + static const scalar integralAdjTolDefault_; //- Edge intersection master catch fraction - static const scalar edgeMasterCatchFraction_; + static const scalar edgeMasterCatchFractionDefault_; //- Edge intersection co-planar tolerance - static const scalar edgeCoPlanarTol_; + static const scalar edgeCoPlanarTolDefault_; //- Edge end cut-off tolerance - static const scalar edgeEndCutoffTol_; + static const scalar edgeEndCutoffTolDefault_; public: @@ -350,6 +376,8 @@ public: //- Return projected points for a slave patch const pointField& pointProjection() const; + //- Set the tolerances from the values in a dictionary + void setTolerances(const dictionary&, bool report=false); //- Write virtual void write(Ostream&) const; diff --git a/src/dynamicMesh/slidingInterface/slidingInterfaceProjectPoints.C b/src/dynamicMesh/slidingInterface/slidingInterfaceProjectPoints.C index a34fd0f7b9..cc4643bcde 100644 --- a/src/dynamicMesh/slidingInterface/slidingInterfaceProjectPoints.C +++ b/src/dynamicMesh/slidingInterface/slidingInterfaceProjectPoints.C @@ -26,22 +26,19 @@ License #include "slidingInterface.H" #include "polyMesh.H" -#include "primitiveMesh.H" #include "line.H" -#include "triPointRef.H" -#include "plane.H" #include "polyTopoChanger.H" // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // -const Foam::scalar Foam::slidingInterface::pointMergeTol_ = 0.05; -const Foam::scalar Foam::slidingInterface::edgeMergeTol_ = 0.01; -const Foam::label Foam::slidingInterface::nFacesPerSlaveEdge_ = 5; -const Foam::label Foam::slidingInterface::edgeFaceEscapeLimit_ = 10; +const Foam::scalar Foam::slidingInterface::pointMergeTolDefault_ = 0.05; +const Foam::scalar Foam::slidingInterface::edgeMergeTolDefault_ = 0.01; +const Foam::label Foam::slidingInterface::nFacesPerSlaveEdgeDefault_ = 5; +const Foam::label Foam::slidingInterface::edgeFaceEscapeLimitDefault_ = 10; -const Foam::scalar Foam::slidingInterface::integralAdjTol_ = 0.05; -const Foam::scalar Foam::slidingInterface::edgeMasterCatchFraction_ = 0.4; -const Foam::scalar Foam::slidingInterface::edgeEndCutoffTol_ = 0.0001; +const Foam::scalar Foam::slidingInterface::integralAdjTolDefault_ = 0.05; +const Foam::scalar Foam::slidingInterface::edgeMasterCatchFractionDefault_ = 0.4; +const Foam::scalar Foam::slidingInterface::edgeEndCutoffTolDefault_ = 0.0001; // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * // diff --git a/src/finiteVolume/interpolation/surfaceInterpolation/limitedSchemes/linearUpwind/linearUpwind.C b/src/finiteVolume/interpolation/surfaceInterpolation/limitedSchemes/linearUpwind/linearUpwind.C index 89843d397e..4cda19451b 100644 --- a/src/finiteVolume/interpolation/surfaceInterpolation/limitedSchemes/linearUpwind/linearUpwind.C +++ b/src/finiteVolume/interpolation/surfaceInterpolation/limitedSchemes/linearUpwind/linearUpwind.C @@ -45,7 +45,7 @@ Foam::linearUpwind::correction ( IOobject ( - "linearUpwindCorrection(" + vf.name() + ')', + "linearUpwind::correction(" + vf.name() + ')', mesh.time().timeName(), mesh, IOobject::NO_READ, diff --git a/src/finiteVolume/interpolation/surfaceInterpolation/limitedSchemes/linearUpwind/linearUpwindV.C b/src/finiteVolume/interpolation/surfaceInterpolation/limitedSchemes/linearUpwind/linearUpwindV.C index db4a31c551..80a137b396 100644 --- a/src/finiteVolume/interpolation/surfaceInterpolation/limitedSchemes/linearUpwind/linearUpwindV.C +++ b/src/finiteVolume/interpolation/surfaceInterpolation/limitedSchemes/linearUpwind/linearUpwindV.C @@ -38,93 +38,96 @@ 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 + ( + "linearUpwindV::correction(" + 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 vectorField& C = mesh.C(); + const vectorField& 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]) + if (faceFlux[facei] > 0.0) + { + maxCorr = + (1.0 - 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] = + 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; + return tsfCorr; } diff --git a/src/finiteVolume/interpolation/surfaceInterpolation/schemes/cubic/cubic.H b/src/finiteVolume/interpolation/surfaceInterpolation/schemes/cubic/cubic.H index 776313c552..f357602d5b 100644 --- a/src/finiteVolume/interpolation/surfaceInterpolation/schemes/cubic/cubic.H +++ b/src/finiteVolume/interpolation/surfaceInterpolation/schemes/cubic/cubic.H @@ -131,9 +131,12 @@ public: ( IOobject ( - vf.name(), + "cubic::correction(" + vf.name() +')', mesh.time().timeName(), - mesh + mesh, + IOobject::NO_READ, + IOobject::NO_WRITE, + false ), surfaceInterpolationScheme::interpolate(vf, kSc, -kSc) ) diff --git a/src/finiteVolume/interpolation/surfaceInterpolation/schemes/localMax/localMax.H b/src/finiteVolume/interpolation/surfaceInterpolation/schemes/localMax/localMax.H index 3fe3f21001..957f4a6920 100644 --- a/src/finiteVolume/interpolation/surfaceInterpolation/schemes/localMax/localMax.H +++ b/src/finiteVolume/interpolation/surfaceInterpolation/schemes/localMax/localMax.H @@ -133,7 +133,7 @@ public: ( IOobject ( - vf.name(), + "localMax::interpolate(" + vf.name() + ')', mesh.time().timeName(), mesh ), diff --git a/src/finiteVolume/interpolation/surfaceInterpolation/schemes/localMin/localMin.H b/src/finiteVolume/interpolation/surfaceInterpolation/schemes/localMin/localMin.H index 302e5bfd85..ef35ca64d3 100644 --- a/src/finiteVolume/interpolation/surfaceInterpolation/schemes/localMin/localMin.H +++ b/src/finiteVolume/interpolation/surfaceInterpolation/schemes/localMin/localMin.H @@ -133,7 +133,7 @@ public: ( IOobject ( - vf.name(), + "localMin::interpolate(" + vf.name() + ')', mesh.time().timeName(), mesh ), diff --git a/src/finiteVolume/interpolation/surfaceInterpolation/schemes/skewCorrected/skewCorrected.H b/src/finiteVolume/interpolation/surfaceInterpolation/schemes/skewCorrected/skewCorrected.H index 047f5c7df1..80122feb07 100644 --- a/src/finiteVolume/interpolation/surfaceInterpolation/schemes/skewCorrected/skewCorrected.H +++ b/src/finiteVolume/interpolation/surfaceInterpolation/schemes/skewCorrected/skewCorrected.H @@ -144,7 +144,7 @@ public: ( IOobject ( - vf.name(), + "skewCorrected::skewCorrection(" + vf.name() + ')', mesh.time().timeName(), mesh ), diff --git a/src/genericPatchFields/genericFvPatchField/genericFvPatchField.C b/src/genericPatchFields/genericFvPatchField/genericFvPatchField.C index 058021f612..ec4c323a91 100644 --- a/src/genericPatchFields/genericFvPatchField/genericFvPatchField.C +++ b/src/genericPatchFields/genericFvPatchField/genericFvPatchField.C @@ -79,7 +79,6 @@ Foam::genericFvPatchField::genericFvPatchField << " (Actual type " << actualTypeName_ << ")" << nl << "\n Please add the 'value' entry to the write function " "of the user-defined boundary-condition\n" - " or link the boundary-condition into libfoamUtil.so" << exit(FatalIOError); } diff --git a/src/meshTools/meshSearch/meshSearch.C b/src/meshTools/meshSearch/meshSearch.C index 2c803a2772..9b822b1151 100644 --- a/src/meshTools/meshSearch/meshSearch.C +++ b/src/meshTools/meshSearch/meshSearch.C @@ -428,7 +428,7 @@ Foam::meshSearch::meshSearch(const polyMesh& mesh, const bool faceDecomp) : mesh_(mesh), faceDecomp_(faceDecomp), - cloud_(mesh_, IDLList()), + cloud_(mesh_, "meshSearchCloud", IDLList()), boundaryTreePtr_(NULL), cellTreePtr_(NULL), cellCentreTreePtr_(NULL) diff --git a/src/meshTools/meshSearch/meshSearch.H b/src/meshTools/meshSearch/meshSearch.H index c8d44bdd48..2f33815756 100644 --- a/src/meshTools/meshSearch/meshSearch.H +++ b/src/meshTools/meshSearch/meshSearch.H @@ -38,8 +38,7 @@ SourceFiles #define meshSearch_H #include "pointIndexHit.H" -#include "Cloud.H" -#include "passiveParticle.H" +#include "passiveParticleCloud.H" // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // @@ -68,7 +67,7 @@ class meshSearch const bool faceDecomp_; //- Dummy cloud to put particles on for tracking. - Cloud cloud_; + passiveParticleCloud cloud_; //- demand driven octrees diff --git a/src/meshTools/triSurface/triSurfaceTools/triSurfaceTools.C b/src/meshTools/triSurface/triSurfaceTools/triSurfaceTools.C index b406fb0467..0e97a63316 100644 --- a/src/meshTools/triSurface/triSurfaceTools/triSurfaceTools.C +++ b/src/meshTools/triSurface/triSurfaceTools/triSurfaceTools.C @@ -2560,7 +2560,7 @@ Foam::triSurface Foam::triSurfaceTools::delaunay2D(const List& pts) // Triangulate int nTris = 0; - dtris2 + int err = dtris2 ( pts.size(), geompackVertices.begin(), @@ -2569,6 +2569,13 @@ Foam::triSurface Foam::triSurfaceTools::delaunay2D(const List& pts) triangle_neighbor.begin() ); + if (err != 0) + { + FatalErrorIn("triSurfaceTools::delaunay2D(const List&)") + << "Failed dtris2 with vertices:" << pts.size() + << abort(FatalError); + } + // Trim triangle_node.setSize(3*nTris); triangle_neighbor.setSize(3*nTris); diff --git a/src/postProcessing/functionObjects/field/Make/files b/src/postProcessing/functionObjects/field/Make/files index e300214731..e950e4d4d5 100644 --- a/src/postProcessing/functionObjects/field/Make/files +++ b/src/postProcessing/functionObjects/field/Make/files @@ -6,4 +6,16 @@ fieldAverage/fieldAverageFunctionObject/fieldAverageFunctionObject.C fieldMinMax/fieldMinMax.C fieldMinMax/fieldMinMaxFunctionObject.C +fieldValues/fieldValue/fieldValue.C +fieldValues/face/faceSource.C +fieldValues/face/faceSourceFunctionObject.C +fieldValues/cell/cellSource.C +fieldValues/cell/cellSourceFunctionObject.C + +streamLine/streamLine.C +streamLine/streamLineParticle.C +streamLine/streamLineParticleCloud.C +streamLine/streamLineFunctionObject.C + + LIB = $(FOAM_LIBBIN)/libfieldFunctionObjects diff --git a/src/postProcessing/functionObjects/field/Make/options b/src/postProcessing/functionObjects/field/Make/options index 5166bcc9e3..1e254fd4bb 100644 --- a/src/postProcessing/functionObjects/field/Make/options +++ b/src/postProcessing/functionObjects/field/Make/options @@ -1,9 +1,11 @@ EXE_INC = \ -I$(LIB_SRC)/finiteVolume/lnInclude \ -I$(LIB_SRC)/meshTools/lnInclude \ + -I$(LIB_SRC)/lagrangian/basic/lnInclude \ -I$(LIB_SRC)/sampling/lnInclude LIB_LIBS = \ -lfiniteVolume \ -lmeshTools \ + -llagrangian \ -lsampling diff --git a/src/postProcessing/functionObjects/field/fieldValues/cell/IOcellSource.H b/src/postProcessing/functionObjects/field/fieldValues/cell/IOcellSource.H new file mode 100644 index 0000000000..b116bb96fd --- /dev/null +++ b/src/postProcessing/functionObjects/field/fieldValues/cell/IOcellSource.H @@ -0,0 +1,51 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2009-2009 OpenCFD Ltd. + \\/ M anipulation | +------------------------------------------------------------------------------- +License + This file is part of OpenFOAM. + + OpenFOAM is free software; you can redistribute it and/or modify it + under the terms of the GNU General Public License as published by the + Free Software Foundation; either version 2 of the License, or (at your + option) any later version. + + OpenFOAM is distributed in the hope that it will be useful, but WITHOUT + ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + for more details. + + You should have received a copy of the GNU General Public License + along with OpenFOAM; if not, write to the Free Software Foundation, + Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA + +Typedef + Foam::IOcellSource + + +Description + Instance of the generic IOOutputFilter for cellSource. + +\*---------------------------------------------------------------------------*/ + +#ifndef IOcellSource_H +#define IOcellSource_H + +#include "cellSource.H" +#include "IOOutputFilter.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ + typedef IOOutputFilter IOcellSource; +} + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#endif + +// ************************************************************************* // diff --git a/src/postProcessing/functionObjects/field/fieldValues/cell/cellSource.C b/src/postProcessing/functionObjects/field/fieldValues/cell/cellSource.C new file mode 100644 index 0000000000..5e63d1287a --- /dev/null +++ b/src/postProcessing/functionObjects/field/fieldValues/cell/cellSource.C @@ -0,0 +1,257 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2009-2009 OpenCFD Ltd. + \\/ M anipulation | +------------------------------------------------------------------------------- +License + This file is part of OpenFOAM. + + OpenFOAM is free software; you can redistribute it and/or modify it + under the terms of the GNU General Public License as published by the + Free Software Foundation; either version 2 of the License, or (at your + option) any later version. + + OpenFOAM is distributed in the hope that it will be useful, but WITHOUT + ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + for more details. + + You should have received a copy of the GNU General Public License + along with OpenFOAM; if not, write to the Free Software Foundation, + Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA + +\*---------------------------------------------------------------------------*/ + +#include "cellSource.H" +#include "fvMesh.H" +#include "volFields.H" +#include "IOList.H" + +// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // + +namespace Foam +{ + namespace fieldValues + { + defineTypeNameAndDebug(cellSource, 0); + } + + template<> + const char* NamedEnum:: + names[] = {"cellZone"}; + + const NamedEnum + fieldValues::cellSource::sourceTypeNames_; + + template<> + const char* NamedEnum:: + names[] = {"none", "sum", "volAverage", "volIntegrate"}; + + const NamedEnum + fieldValues::cellSource::operationTypeNames_; + +} + + +// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * // + +void Foam::fieldValues::cellSource::setCellZoneCells() +{ + label zoneId = mesh().cellZones().findZoneID(sourceName_); + + if (zoneId < 0) + { + FatalErrorIn("cellSource::cellSource::setCellZoneCells()") + << "Unknown cell zone name: " << sourceName_ + << ". Valid cell zones are: " << mesh().cellZones().names() + << nl << exit(FatalError); + } + + const cellZone& cZone = mesh().cellZones()[zoneId]; + + cellId_.setSize(cZone.size()); + + label count = 0; + forAll(cZone, i) + { + label cellI = cZone[i]; + cellId_[count] = cellI; + count++; + } + + cellId_.setSize(count); + + if (debug) + { + Info<< "Original cell zone size = " << cZone.size() + << ", new size = " << count << endl; + } +} + + +// * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * // + +void Foam::fieldValues::cellSource::initialise() +{ + switch (source_) + { + case stCellZone: + { + setCellZoneCells(); + break; + } + default: + { + FatalErrorIn("cellSource::constructCellAddressing()") + << "Unknown source type. Valid source types are:" + << sourceTypeNames_ << nl << exit(FatalError); + } + } + + Info<< type() << " " << name_ << ":" << nl + << " total cells = " << cellId_.size() << nl + << " total volume = " << sum(filterField(mesh().V())) + << nl << endl; +} + + +void Foam::fieldValues::cellSource::makeFile() +{ + // Create the forces file if not already created + if (outputFilePtr_.empty()) + { + if (debug) + { + Info<< "Creating output file." << endl; + } + + // File update + if (Pstream::master()) + { + fileName outputDir; + if (Pstream::parRun()) + { + // Put in undecomposed case (Note: gives problems for + // distributed data running) + outputDir = + obr_.time().path()/".."/name_/obr_.time().timeName(); + } + else + { + outputDir = obr_.time().path()/name_/obr_.time().timeName(); + } + + // Create directory if does not exist + mkDir(outputDir); + + // Open new file at start up + outputFilePtr_.reset(new OFstream(outputDir/(type() + ".dat"))); + + // Add headers to output data + writeFileHeader(); + } + } +} + + +void Foam::fieldValues::cellSource::writeFileHeader() +{ + if (outputFilePtr_.valid()) + { + outputFilePtr_() + << "# Source : " << sourceTypeNames_[source_] << " " + << sourceName_ << nl << "# Cells : " << cellId_.size() << nl + << "# Time" << tab << "sum(V)"; + + forAll(fields_, i) + { + outputFilePtr_() + << tab << operationTypeNames_[operation_] + << "(" << fields_[i] << ")"; + } + + outputFilePtr_() << endl; + } +} + + +// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // + +Foam::fieldValues::cellSource::cellSource +( + const word& name, + const objectRegistry& obr, + const dictionary& dict, + const bool loadFromFiles +) +: + fieldValue(name, obr, dict, loadFromFiles), + source_(sourceTypeNames_.read(dict.lookup("source"))), + operation_(operationTypeNames_.read(dict.lookup("operation"))), + cellId_(), + outputFilePtr_(NULL) +{ + initialise(); + + if (active_) + { + // Create the output file if not already created + makeFile(); + } +} + + +// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * // + +Foam::fieldValues::cellSource::~cellSource() +{} + + +// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // + +void Foam::fieldValues::cellSource::read(const dictionary& dict) +{ + if (active_) + { + fieldValue::read(dict); + initialise(); + } +} + + +void Foam::fieldValues::cellSource::write() +{ + if (active_) + { + if (log_) + { + Info<< type() << " " << name_ << " output:" << nl; + } + + outputFilePtr_() + << obr_.time().value() << tab + << sum(filterField(mesh().V())); + + forAll(fields_, i) + { + writeValues(fields_[i]); + writeValues(fields_[i]); + writeValues(fields_[i]); + writeValues(fields_[i]); + writeValues(fields_[i]); + } + + outputFilePtr_()<< endl; + + if (log_) + { + Info<< endl; + } + } +} + + +// ************************************************************************* // + diff --git a/src/postProcessing/functionObjects/field/fieldValues/cell/cellSource.H b/src/postProcessing/functionObjects/field/fieldValues/cell/cellSource.H new file mode 100644 index 0000000000..d1d08cdbd9 --- /dev/null +++ b/src/postProcessing/functionObjects/field/fieldValues/cell/cellSource.H @@ -0,0 +1,227 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2009-2009 OpenCFD Ltd. + \\/ M anipulation | +------------------------------------------------------------------------------- +License + This file is part of OpenFOAM. + + OpenFOAM is free software; you can redistribute it and/or modify it + under the terms of the GNU General Public License as published by the + Free Software Foundation; either version 2 of the License, or (at your + option) any later version. + + OpenFOAM is distributed in the hope that it will be useful, but WITHOUT + ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + for more details. + + You should have received a copy of the GNU General Public License + along with OpenFOAM; if not, write to the Free Software Foundation, + Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA + +Class + Foam::fieldValues::cellSource + +Description + Cell source variant of field value function object. Values of user- + specified fields reported for collections of cells. + + cellObj1 // Name also used to identify output folder + { + type cellSource; + functionObjectLibs ("libfieldValueFunctionObjects.so"); + enabled true; + outputControl outputTime; + log true; // log to screen? + valueOutput true; // Write values at run-time output times? + source cellZone; // Type of cell source + sourceName c0; + operation volAverage; // none, sum, volAverage, volIntegrate + fields + ( + p + U + ); + } + +SourceFiles + cellSource.C + +\*---------------------------------------------------------------------------*/ + +#ifndef cellSource_H +#define cellSource_H + +#include "NamedEnum.H" +#include "fieldValue.H" +#include "labelList.H" +#include "OFstream.H" +#include "volFieldsFwd.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ +namespace fieldValues +{ + +/*---------------------------------------------------------------------------*\ + Class cellSource Declaration +\*---------------------------------------------------------------------------*/ + +class cellSource +: + public fieldValue +{ + +public: + + // Public data types + + //- Source type enumeration + enum sourceType + { + stCellZone + }; + + //- Source type names + static const NamedEnum sourceTypeNames_; + + + //- Operation type enumeration + enum operationType + { + opNone, + opSum, + opVolAverage, + opVolIntegrate + }; + + //- Operation type names + static const NamedEnum operationTypeNames_; + + +private: + + // Private member functions + + //- Set cells to evaluate based on a cell zone + void setCellZoneCells(); + + //- Set cells to evaluate based on a patch + void setPatchCells(); + + //- Create the output file if not already created + void makeFile(); + + +protected: + + // Protected data + + //- Source type + sourceType source_; + + //- Operation to apply to values + operationType operation_; + + //- Local list of cell IDs + labelList cellId_; + + //- Output file pointer + autoPtr outputFilePtr_; + + + // Protected member functions + + //- Initialise, e.g. cell addressing + void initialise(); + + //- Insert field values into values list + template + bool setFieldValues + ( + const word& fieldName, + List& values + ) const; + + //- Apply the 'operation' to the values + template + Type processValues(const List& values) const; + + //- Output file header information + virtual void writeFileHeader(); + + +public: + + //- Run-time type information + TypeName("cellSource"); + + + //- Construct from components + cellSource + ( + const word& name, + const objectRegistry& obr, + const dictionary& dict, + const bool loadFromFiles = false + ); + + + //- Destructor + virtual ~cellSource(); + + + // Public member functions + + // Access + + //- Return the source type + inline const sourceType& source() const; + + //- Return the local list of cell IDs + inline const labelList& cellId() const; + + + // Function object functions + + //- Read from dictionary + virtual void read(const dictionary&); + + //- Calculate and write + virtual void write(); + + //- Templated helper function to output field values + template + bool writeValues(const word& fieldName); + + //- Filter a field according to cellIds + template + tmp > filterField(const Field& field) const; +}; + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace fieldValues +} // End namespace Foam + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#include "cellSourceI.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#ifdef NoRepository + #include "cellSourceTemplates.C" +#endif + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#endif + +// ************************************************************************* // diff --git a/src/postProcessing/functionObjects/field/fieldValues/cell/cellSourceFunctionObject.C b/src/postProcessing/functionObjects/field/fieldValues/cell/cellSourceFunctionObject.C new file mode 100644 index 0000000000..24f4b7bbb8 --- /dev/null +++ b/src/postProcessing/functionObjects/field/fieldValues/cell/cellSourceFunctionObject.C @@ -0,0 +1,47 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2009-2009 OpenCFD Ltd. + \\/ M anipulation | +------------------------------------------------------------------------------- +License + This file is part of OpenFOAM. + + OpenFOAM is free software; you can redistribute it and/or modify it + under the terms of the GNU General Public License as published by the + Free Software Foundation; either version 2 of the License, or (at your + option) any later version. + + OpenFOAM is distributed in the hope that it will be useful, but WITHOUT + ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + for more details. + + You should have received a copy of the GNU General Public License + along with OpenFOAM; if not, write to the Free Software Foundation, + Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA + +\*---------------------------------------------------------------------------*/ + +#include "cellSourceFunctionObject.H" + +// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // + +namespace Foam +{ + defineNamedTemplateTypeNameAndDebug + ( + cellSourceFunctionObject, + 0 + ); + + addToRunTimeSelectionTable + ( + functionObject, + cellSourceFunctionObject, + dictionary + ); +} + +// ************************************************************************* // diff --git a/src/postProcessing/functionObjects/field/fieldValues/cell/cellSourceFunctionObject.H b/src/postProcessing/functionObjects/field/fieldValues/cell/cellSourceFunctionObject.H new file mode 100644 index 0000000000..eca172a58d --- /dev/null +++ b/src/postProcessing/functionObjects/field/fieldValues/cell/cellSourceFunctionObject.H @@ -0,0 +1,55 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2009-2009 OpenCFD Ltd. + \\/ M anipulation | +------------------------------------------------------------------------------- +License + This file is part of OpenFOAM. + + OpenFOAM is free software; you can redistribute it and/or modify it + under the terms of the GNU General Public License as published by the + Free Software Foundation; either version 2 of the License, or (at your + option) any later version. + + OpenFOAM is distributed in the hope that it will be useful, but WITHOUT + ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + for more details. + + You should have received a copy of the GNU General Public License + along with OpenFOAM; if not, write to the Free Software Foundation, + Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA + +Typedef + Foam::cellSourceFunctionObject + +Description + FunctionObject wrapper around cellSource to allow it to be + created via the functions list within controlDict. + +SourceFiles + cellSourceFunctionObject.C + +\*---------------------------------------------------------------------------*/ + +#ifndef cellSourceFunctionObject_H +#define cellSourceFunctionObject_H + +#include "cellSource.H" +#include "OutputFilterFunctionObject.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ + typedef OutputFilterFunctionObject + cellSourceFunctionObject; +} + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#endif + +// ************************************************************************* // diff --git a/src/postProcessing/functionObjects/field/fieldValues/cell/cellSourceI.H b/src/postProcessing/functionObjects/field/fieldValues/cell/cellSourceI.H new file mode 100644 index 0000000000..b86462d2cc --- /dev/null +++ b/src/postProcessing/functionObjects/field/fieldValues/cell/cellSourceI.H @@ -0,0 +1,46 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2009-2009 OpenCFD Ltd. + \\/ M anipulation | +------------------------------------------------------------------------------- +License + This file is part of OpenFOAM. + + OpenFOAM is free software; you can redistribute it and/or modify it + under the terms of the GNU General Public License as published by the + Free Software Foundation; either version 2 of the License, or (at your + option) any later version. + + OpenFOAM is distributed in the hope that it will be useful, but WITHOUT + ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + for more details. + + You should have received a copy of the GNU General Public License + along with OpenFOAM; if not, write to the Free Software Foundation, + Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA + +\*---------------------------------------------------------------------------*/ + +#include "cellSource.H" + +// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // + +inline const Foam::fieldValues::cellSource::sourceType& +Foam::fieldValues::cellSource::source() const +{ + return source_; +} + + +inline const Foam::labelList& +Foam::fieldValues::cellSource::cellId() const +{ + return cellId_; +} + + +// ************************************************************************* // + diff --git a/src/postProcessing/functionObjects/field/fieldValues/cell/cellSourceTemplates.C b/src/postProcessing/functionObjects/field/fieldValues/cell/cellSourceTemplates.C new file mode 100644 index 0000000000..fe495fcca6 --- /dev/null +++ b/src/postProcessing/functionObjects/field/fieldValues/cell/cellSourceTemplates.C @@ -0,0 +1,162 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2009-2009 OpenCFD Ltd. + \\/ M anipulation | +------------------------------------------------------------------------------- +License + This file is part of OpenFOAM. + + OpenFOAM is free software; you can redistribute it and/or modify it + under the terms of the GNU General Public License as published by the + Free Software Foundation; either version 2 of the License, or (at your + option) any later version. + + OpenFOAM is distributed in the hope that it will be useful, but WITHOUT + ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + for more details. + + You should have received a copy of the GNU General Public License + along with OpenFOAM; if not, write to the Free Software Foundation, + Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA + +\*---------------------------------------------------------------------------*/ + +#include "cellSource.H" +#include "volFields.H" +#include "IOList.H" +#include "ListListOps.H" + +// * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * // + +template +bool Foam::fieldValues::cellSource::setFieldValues +( + const word& fieldName, + List& values +) const +{ + values.setSize(cellId_.size(), pTraits::zero); + + typedef GeometricField vf; + + if (obr_.foundObject(fieldName)) + { + const vf& field = obr_.lookupObject(fieldName); + + values = UIndirectList(field, cellId_); + + return true; + } + + return false; +} + + +template +Type Foam::fieldValues::cellSource::processValues +( + const List& values +) const +{ + Type result = pTraits::zero; + switch (operation_) + { + case opSum: + { + result = sum(values); + break; + } + case opVolAverage: + { + tmp V = filterField(mesh().V()); + result = sum(values*V())/sum(V()); + break; + } + case opVolIntegrate: + { + result = sum(values*filterField(mesh().V())); + break; + } + default: + { + // Do nothing + } + } + + return result; +} + + +// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // + +template +bool Foam::fieldValues::cellSource::writeValues(const word& fieldName) +{ + List > allValues(Pstream::nProcs()); + + bool validField = + setFieldValues(fieldName, allValues[Pstream::myProcNo()]); + + if (validField) + { + Pstream::gatherList(allValues); + + if (Pstream::master()) + { + List values = + ListListOps::combine > + ( + allValues, + accessOp >() + ); + + Type result = processValues(values); + + if (valueOutput_) + { + IOList + ( + IOobject + ( + fieldName + "_" + sourceTypeNames_[source_] + "-" + + sourceName_, + obr_.time().timeName(), + obr_, + IOobject::NO_READ, + IOobject::NO_WRITE + ), + values + ).write(); + } + + + outputFilePtr_()<< tab << result; + + if (log_) + { + Info<< " " << operationTypeNames_[operation_] + << "(" << sourceName_ << ") for " << fieldName + << " = " << result << endl; + } + } + } + + return validField; +} + + +template +Foam::tmp > Foam::fieldValues::cellSource::filterField +( + const Field& field +) const +{ + return tmp >(new Field(field, cellId_)); +} + + +// ************************************************************************* // + diff --git a/src/postProcessing/functionObjects/field/fieldValues/controlDict b/src/postProcessing/functionObjects/field/fieldValues/controlDict new file mode 100644 index 0000000000..7d9ab1f656 --- /dev/null +++ b/src/postProcessing/functionObjects/field/fieldValues/controlDict @@ -0,0 +1,108 @@ +/*--------------------------------*- C++ -*----------------------------------*\ +| ========= | | +| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox | +| \\ / O peration | Version: 1.6 | +| \\ / A nd | Web: www.OpenFOAM.org | +| \\/ M anipulation | | +\*---------------------------------------------------------------------------*/ +FoamFile +{ + version 2.0; + format ascii; + class dictionary; + location "system"; + object controlDict; +} +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +application icoFoam; + +startFrom startTime; + +startTime 0; + +stopAt endTime; + +endTime 0.5; + +deltaT 0.005; + +writeControl timeStep; + +writeInterval 20; + +purgeWrite 0; + +writeFormat ascii; + +writePrecision 6; + +writeCompression uncompressed; + +timeFormat general; + +timePrecision 6; + +runTimeModifiable yes; + +functions +( + faceObj1 + { + type faceSource; + functionObjectLibs ("libfieldFunctionObjects.so"); + enabled true; + outputControl outputTime; + log true; + valueOutput true; + source patch; + sourceName movingWall; +// source faceZone; +// sourceName f0; + operation areaAverage; + fields + ( + p + phi + U + ); + } + + faceObj2 + { + type faceSource; + functionObjectLibs ("libfieldFunctionObjects.so"); + enabled true; + outputControl outputTime; + log true; + valueOutput true; + source faceZone; + sourceName f0; + operation sum; + fields + ( + phi + ); + } + + cellObj1 + { + type cellSource; + functionObjectLibs ("libfieldFunctionObjects.so"); + enabled true; + outputControl outputTime; + log true; + valueOutput true; + source cellZone; + sourceName c0; + operation volAverage; + fields + ( + p + U + ); + } +); + + +// ************************************************************************* // diff --git a/src/postProcessing/functionObjects/field/fieldValues/face/IOfaceSource.H b/src/postProcessing/functionObjects/field/fieldValues/face/IOfaceSource.H new file mode 100644 index 0000000000..caca9968c6 --- /dev/null +++ b/src/postProcessing/functionObjects/field/fieldValues/face/IOfaceSource.H @@ -0,0 +1,51 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2009-2009 OpenCFD Ltd. + \\/ M anipulation | +------------------------------------------------------------------------------- +License + This file is part of OpenFOAM. + + OpenFOAM is free software; you can redistribute it and/or modify it + under the terms of the GNU General Public License as published by the + Free Software Foundation; either version 2 of the License, or (at your + option) any later version. + + OpenFOAM is distributed in the hope that it will be useful, but WITHOUT + ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + for more details. + + You should have received a copy of the GNU General Public License + along with OpenFOAM; if not, write to the Free Software Foundation, + Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA + +Typedef + Foam::IOfaceSource + + +Description + Instance of the generic IOOutputFilter for faceSource. + +\*---------------------------------------------------------------------------*/ + +#ifndef IOfaceSource_H +#define IOfaceSource_H + +#include "faceSource.H" +#include "IOOutputFilter.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ + typedef IOOutputFilter IOfaceSource; +} + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#endif + +// ************************************************************************* // diff --git a/src/postProcessing/functionObjects/field/fieldValues/face/faceSource.C b/src/postProcessing/functionObjects/field/fieldValues/face/faceSource.C new file mode 100644 index 0000000000..76c1ac9515 --- /dev/null +++ b/src/postProcessing/functionObjects/field/fieldValues/face/faceSource.C @@ -0,0 +1,368 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2009-2009 OpenCFD Ltd. + \\/ M anipulation | +------------------------------------------------------------------------------- +License + This file is part of OpenFOAM. + + OpenFOAM is free software; you can redistribute it and/or modify it + under the terms of the GNU General Public License as published by the + Free Software Foundation; either version 2 of the License, or (at your + option) any later version. + + OpenFOAM is distributed in the hope that it will be useful, but WITHOUT + ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + for more details. + + You should have received a copy of the GNU General Public License + along with OpenFOAM; if not, write to the Free Software Foundation, + Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA + +\*---------------------------------------------------------------------------*/ + +#include "faceSource.H" +#include "fvMesh.H" +#include "cyclicPolyPatch.H" +#include "emptyPolyPatch.H" +#include "processorPolyPatch.H" +#include "surfaceFields.H" +#include "volFields.H" +#include "IOList.H" + +// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // + +namespace Foam +{ + namespace fieldValues + { + defineTypeNameAndDebug(faceSource, 0); + } + + template<> + const char* NamedEnum:: + names[] = {"faceZone", "patch"}; + + const NamedEnum + fieldValues::faceSource::sourceTypeNames_; + + template<> + const char* NamedEnum:: + names[] = {"none", "sum", "areaAverage", "areaIntegrate"}; + + const NamedEnum + fieldValues::faceSource::operationTypeNames_; + +} + + +// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * // + +void Foam::fieldValues::faceSource::setFaceZoneFaces() +{ + label zoneId = mesh().faceZones().findZoneID(sourceName_); + + if (zoneId < 0) + { + FatalErrorIn("faceSource::faceSource::setFaceZoneFaces()") + << "Unknown face zone name: " << sourceName_ + << ". Valid face zones are: " << mesh().faceZones().names() + << nl << exit(FatalError); + } + + const faceZone& fZone = mesh().faceZones()[zoneId]; + + faceId_.setSize(fZone.size()); + facePatchId_.setSize(fZone.size()); + flipMap_.setSize(fZone.size()); + + label count = 0; + forAll(fZone, i) + { + label faceI = fZone[i]; + + label faceId = -1; + label facePatchId = -1; + if (mesh().isInternalFace(faceI)) + { + faceId = faceI; + facePatchId = -1; + } + else + { + facePatchId = mesh().boundaryMesh().whichPatch(faceI); + const polyPatch& pp = mesh().boundaryMesh()[facePatchId]; + if (isA(pp)) + { + if (refCast(pp).owner()) + { + faceId = pp.whichFace(faceI); + } + else + { + faceId = -1; + } + } + else if (isA(pp)) + { + label patchFaceI = faceI - pp.start(); + if (patchFaceI < pp.size()/2) + { + faceId = patchFaceI; + } + else + { + faceId = -1; + } + } + else if (!isA(pp)) + { + faceId = faceI - pp.start(); + } + else + { + faceId = -1; + facePatchId = -1; + } + } + + if (faceId >= 0) + { + if (fZone.flipMap()[i]) + { + flipMap_[count] = -1; + } + else + { + flipMap_[count] = 1; + } + faceId_[count] = faceId; + facePatchId_[count] = facePatchId; + count++; + } + } + + faceId_.setSize(count); + facePatchId_.setSize(count); + flipMap_.setSize(count); + + if (debug) + { + Info<< "Original face zone size = " << fZone.size() + << ", new size = " << count << endl; + } +} + + +void Foam::fieldValues::faceSource::setPatchFaces() +{ + label patchId = mesh().boundaryMesh().findPatchID(sourceName_); + + if (patchId < 0) + { + FatalErrorIn("faceSource::constructFaceAddressing()") + << "Unknown patch name: " << sourceName_ + << ". Valid patch names are: " + << mesh().boundaryMesh().names() << nl + << exit(FatalError); + } + + const polyPatch& pp = mesh().boundaryMesh()[patchId]; + + label nFaces = pp.size(); + if (isA(pp)) + { + nFaces /= 2; + } + else if (isA(pp)) + { + nFaces = 0; + } + + faceId_.setSize(nFaces); + facePatchId_.setSize(nFaces); + flipMap_.setSize(nFaces); + + forAll(faceId_, faceI) + { + faceId_[faceI] = faceI; + facePatchId_[faceI] = patchId; + flipMap_[faceI] = 1; + } +} + + +// * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * // + +void Foam::fieldValues::faceSource::initialise() +{ + switch (source_) + { + case stFaceZone: + { + setFaceZoneFaces(); + break; + } + case stPatch: + { + setPatchFaces(); + break; + } + default: + { + FatalErrorIn("faceSource::constructFaceAddressing()") + << "Unknown source type. Valid source types are:" + << sourceTypeNames_ << nl << exit(FatalError); + } + } + + Info<< type() << " " << name_ << ":" << nl + << " total faces = " << faceId_.size() << nl + << " total area = " << sum(filterField(mesh().magSf())) + << nl << endl; +} + + +void Foam::fieldValues::faceSource::makeFile() +{ + // Create the forces file if not already created + if (outputFilePtr_.empty()) + { + if (debug) + { + Info<< "Creating output file." << endl; + } + + // File update + if (Pstream::master()) + { + fileName outputDir; + if (Pstream::parRun()) + { + // Put in undecomposed case (Note: gives problems for + // distributed data running) + outputDir = + obr_.time().path()/".."/name_/obr_.time().timeName(); + } + else + { + outputDir = obr_.time().path()/name_/obr_.time().timeName(); + } + + // Create directory if does not exist + mkDir(outputDir); + + // Open new file at start up + outputFilePtr_.reset(new OFstream(outputDir/(type() + ".dat"))); + + // Add headers to output data + writeFileHeader(); + } + } +} + + +void Foam::fieldValues::faceSource::writeFileHeader() +{ + if (outputFilePtr_.valid()) + { + outputFilePtr_() + << "# Source : " << sourceTypeNames_[source_] << " " + << sourceName_ << nl << "# Faces : " << faceId_.size() << nl + << "# Time" << tab << "sum(magSf)"; + + forAll(fields_, i) + { + outputFilePtr_() + << tab << operationTypeNames_[operation_] + << "(" << fields_[i] << ")"; + } + + outputFilePtr_() << endl; + } +} + + +// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // + +Foam::fieldValues::faceSource::faceSource +( + const word& name, + const objectRegistry& obr, + const dictionary& dict, + const bool loadFromFiles +) +: + fieldValue(name, obr, dict, loadFromFiles), + source_(sourceTypeNames_.read(dict.lookup("source"))), + operation_(operationTypeNames_.read(dict.lookup("operation"))), + faceId_(), + facePatchId_(), + flipMap_(), + outputFilePtr_(NULL) +{ + initialise(); + + if (active_) + { + // Create the output file if not already created + makeFile(); + } +} + + +// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * // + +Foam::fieldValues::faceSource::~faceSource() +{} + + +// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // + +void Foam::fieldValues::faceSource::read(const dictionary& dict) +{ + if (active_) + { + fieldValue::read(dict); + initialise(); + } +} + + +void Foam::fieldValues::faceSource::write() +{ + if (active_) + { + if (log_) + { + Info<< type() << " " << name_ << " output:" << nl; + } + + outputFilePtr_() + << obr_.time().value() << tab + << sum(filterField(mesh().magSf())); + + forAll(fields_, i) + { + writeValues(fields_[i]); + writeValues(fields_[i]); + writeValues(fields_[i]); + writeValues(fields_[i]); + writeValues(fields_[i]); + } + + outputFilePtr_()<< endl; + + if (log_) + { + Info<< endl; + } + } +} + + +// ************************************************************************* // + diff --git a/src/postProcessing/functionObjects/field/fieldValues/face/faceSource.H b/src/postProcessing/functionObjects/field/fieldValues/face/faceSource.H new file mode 100644 index 0000000000..fa89fcdd6c --- /dev/null +++ b/src/postProcessing/functionObjects/field/fieldValues/face/faceSource.H @@ -0,0 +1,252 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2009-2009 OpenCFD Ltd. + \\/ M anipulation | +------------------------------------------------------------------------------- +License + This file is part of OpenFOAM. + + OpenFOAM is free software; you can redistribute it and/or modify it + under the terms of the GNU General Public License as published by the + Free Software Foundation; either version 2 of the License, or (at your + option) any later version. + + OpenFOAM is distributed in the hope that it will be useful, but WITHOUT + ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + for more details. + + You should have received a copy of the GNU General Public License + along with OpenFOAM; if not, write to the Free Software Foundation, + Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA + +Class + Foam::fieldValues::faceSource + +Description + Face source variant of field value function object. Values of user- + specified fields reported for collections of faces. + + cellObj1 // Name also used to identify output folder + { + type cellSource; + functionObjectLibs ("libfieldValueFunctionObjects.so"); + enabled true; + outputControl outputTime; + log true; // log to screen? + valueOutput true; // Write values at run-time output times? + source faceZone; // Type of face source: faceZone, patch + sourceName f0; + operation sum; // none, sum, areaAverage, areaIntegrate + fields + ( + p + phi + U + ); + } + +SourceFiles + faceSource.C + +\*---------------------------------------------------------------------------*/ + +#ifndef faceSource_H +#define faceSource_H + +#include "NamedEnum.H" +#include "fieldValue.H" +#include "labelList.H" +#include "OFstream.H" +#include "surfaceFieldsFwd.H" +#include "volFieldsFwd.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ +namespace fieldValues +{ + +/*---------------------------------------------------------------------------*\ + Class faceSource Declaration +\*---------------------------------------------------------------------------*/ + +class faceSource +: + public fieldValue +{ + +public: + + // Public data types + + //- Source type enumeration + enum sourceType + { + stFaceZone, + stPatch + }; + + //- Source type names + static const NamedEnum sourceTypeNames_; + + + //- Operation type enumeration + enum operationType + { + opNone, + opSum, + opAreaAverage, + opAreaIntegrate + }; + + //- Operation type names + static const NamedEnum operationTypeNames_; + + +private: + + // Private member functions + + //- Set faces to evaluate based on a face zone + void setFaceZoneFaces(); + + //- Set faces to evaluate based on a patch + void setPatchFaces(); + + //- Create the output file if not already created + void makeFile(); + + +protected: + + // Protected data + + //- Source type + sourceType source_; + + //- Operation to apply to values + operationType operation_; + + //- Local list of face IDs + labelList faceId_; + + //- Local list of patch ID per face + labelList facePatchId_; + + //- List of +1/-1 representing face flip map + labelList flipMap_; + + //- Output file pointer + autoPtr outputFilePtr_; + + + // Protected member functions + + //- Initialise, e.g. face addressing + void initialise(); + + //- Insert field values into values list + template + bool setFieldValues + ( + const word& fieldName, + List& values + ) const; + + //- Apply the 'operation' to the values + template + Type processValues(const List& values) const; + + //- Output file header information + virtual void writeFileHeader(); + + +public: + + //- Run-time type information + TypeName("faceSource"); + + + //- Construct from components + faceSource + ( + const word& name, + const objectRegistry& obr, + const dictionary& dict, + const bool loadFromFiles = false + ); + + + //- Destructor + virtual ~faceSource(); + + + // Public member functions + + // Access + + //- Return the source type + inline const sourceType& source() const; + + //- Return the local list of face IDs + inline const labelList& faceId() const; + + //- Return the local list of patch ID per face + inline const labelList& facePatch() const; + + //- Return the list of +1/-1 representing face flip map + inline const labelList& flipMap() const; + + + // Function object functions + + //- Read from dictionary + virtual void read(const dictionary&); + + //- Calculate and write + virtual void write(); + + //- Templated helper function to output field values + template + bool writeValues(const word& fieldName); + + //- Filter a surface field according to faceIds + template + tmp > filterField + ( + const GeometricField& field + ) const; + + //- Filter a volume field according to faceIds + template + tmp > filterField + ( + const GeometricField& field + ) const; +}; + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace fieldValues +} // End namespace Foam + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#include "faceSourceI.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#ifdef NoRepository + #include "faceSourceTemplates.C" +#endif + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#endif + +// ************************************************************************* // diff --git a/src/postProcessing/functionObjects/field/fieldValues/face/faceSourceFunctionObject.C b/src/postProcessing/functionObjects/field/fieldValues/face/faceSourceFunctionObject.C new file mode 100644 index 0000000000..12019c86c2 --- /dev/null +++ b/src/postProcessing/functionObjects/field/fieldValues/face/faceSourceFunctionObject.C @@ -0,0 +1,47 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2009-2009 OpenCFD Ltd. + \\/ M anipulation | +------------------------------------------------------------------------------- +License + This file is part of OpenFOAM. + + OpenFOAM is free software; you can redistribute it and/or modify it + under the terms of the GNU General Public License as published by the + Free Software Foundation; either version 2 of the License, or (at your + option) any later version. + + OpenFOAM is distributed in the hope that it will be useful, but WITHOUT + ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + for more details. + + You should have received a copy of the GNU General Public License + along with OpenFOAM; if not, write to the Free Software Foundation, + Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA + +\*---------------------------------------------------------------------------*/ + +#include "faceSourceFunctionObject.H" + +// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // + +namespace Foam +{ + defineNamedTemplateTypeNameAndDebug + ( + faceSourceFunctionObject, + 0 + ); + + addToRunTimeSelectionTable + ( + functionObject, + faceSourceFunctionObject, + dictionary + ); +} + +// ************************************************************************* // diff --git a/src/postProcessing/functionObjects/field/fieldValues/face/faceSourceFunctionObject.H b/src/postProcessing/functionObjects/field/fieldValues/face/faceSourceFunctionObject.H new file mode 100644 index 0000000000..54d717f870 --- /dev/null +++ b/src/postProcessing/functionObjects/field/fieldValues/face/faceSourceFunctionObject.H @@ -0,0 +1,55 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2009-2009 OpenCFD Ltd. + \\/ M anipulation | +------------------------------------------------------------------------------- +License + This file is part of OpenFOAM. + + OpenFOAM is free software; you can redistribute it and/or modify it + under the terms of the GNU General Public License as published by the + Free Software Foundation; either version 2 of the License, or (at your + option) any later version. + + OpenFOAM is distributed in the hope that it will be useful, but WITHOUT + ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + for more details. + + You should have received a copy of the GNU General Public License + along with OpenFOAM; if not, write to the Free Software Foundation, + Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA + +Typedef + Foam::faceSourceFunctionObject + +Description + FunctionObject wrapper around faceSource to allow it to be + created via the functions list within controlDict. + +SourceFiles + faceSourceFunctionObject.C + +\*---------------------------------------------------------------------------*/ + +#ifndef faceSourceFunctionObject_H +#define faceSourceFunctionObject_H + +#include "faceSource.H" +#include "OutputFilterFunctionObject.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ + typedef OutputFilterFunctionObject + faceSourceFunctionObject; +} + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#endif + +// ************************************************************************* // diff --git a/src/postProcessing/functionObjects/field/fieldValues/face/faceSourceI.H b/src/postProcessing/functionObjects/field/fieldValues/face/faceSourceI.H new file mode 100644 index 0000000000..9d86255c36 --- /dev/null +++ b/src/postProcessing/functionObjects/field/fieldValues/face/faceSourceI.H @@ -0,0 +1,60 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2009-2009 OpenCFD Ltd. + \\/ M anipulation | +------------------------------------------------------------------------------- +License + This file is part of OpenFOAM. + + OpenFOAM is free software; you can redistribute it and/or modify it + under the terms of the GNU General Public License as published by the + Free Software Foundation; either version 2 of the License, or (at your + option) any later version. + + OpenFOAM is distributed in the hope that it will be useful, but WITHOUT + ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + for more details. + + You should have received a copy of the GNU General Public License + along with OpenFOAM; if not, write to the Free Software Foundation, + Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA + +\*---------------------------------------------------------------------------*/ + +#include "faceSource.H" + +// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // + +inline const Foam::fieldValues::faceSource::sourceType& +Foam::fieldValues::faceSource::source() const +{ + return source_; +} + + +inline const Foam::labelList& +Foam::fieldValues::faceSource::faceId() const +{ + return faceId_; +} + + +inline const Foam::labelList& +Foam::fieldValues::faceSource::facePatch() const +{ + return facePatchId_; +} + + +inline const Foam::labelList& +Foam::fieldValues::faceSource::flipMap() const +{ + return flipMap_; +} + + +// ************************************************************************* // + diff --git a/src/postProcessing/functionObjects/field/fieldValues/face/faceSourceTemplates.C b/src/postProcessing/functionObjects/field/fieldValues/face/faceSourceTemplates.C new file mode 100644 index 0000000000..ced2ee82ee --- /dev/null +++ b/src/postProcessing/functionObjects/field/fieldValues/face/faceSourceTemplates.C @@ -0,0 +1,269 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2009-2009 OpenCFD Ltd. + \\/ M anipulation | +------------------------------------------------------------------------------- +License + This file is part of OpenFOAM. + + OpenFOAM is free software; you can redistribute it and/or modify it + under the terms of the GNU General Public License as published by the + Free Software Foundation; either version 2 of the License, or (at your + option) any later version. + + OpenFOAM is distributed in the hope that it will be useful, but WITHOUT + ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + for more details. + + You should have received a copy of the GNU General Public License + along with OpenFOAM; if not, write to the Free Software Foundation, + Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA + +\*---------------------------------------------------------------------------*/ + +#include "faceSource.H" +#include "surfaceFields.H" +#include "volFields.H" +#include "IOList.H" +#include "ListListOps.H" + +// * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * // + +template +bool Foam::fieldValues::faceSource::setFieldValues +( + const word& fieldName, + List& values +) const +{ + values.setSize(faceId_.size(), pTraits::zero); + + typedef GeometricField sf; + typedef GeometricField vf; + + if (obr_.foundObject(fieldName)) + { + const sf& field = obr_.lookupObject(fieldName); + + forAll(values, i) + { + label faceI = faceId_[i]; + label patchI = facePatchId_[i]; + if (patchI >= 0) + { + values[i] = field.boundaryField()[patchI][faceI]; + } + else + { + values[i] = field[faceI]; + } + + values[i] *= flipMap_[i]; + } + + return true; + } + else if (obr_.foundObject(fieldName)) + { + const vf& field = obr_.lookupObject(fieldName); + + forAll(values, i) + { + label faceI = faceId_[i]; + label patchI = facePatchId_[i]; + if (patchI >= 0) + { + values[i] = field.boundaryField()[patchI][faceI]; + } + else + { + FatalErrorIn + ( + "fieldValues::faceSource::setFieldValues" + "(" + "const word&, " + "List&" + ") const" + ) << type() << " " << name_ << ": " + << sourceTypeNames_[source_] << "(" << sourceName_ << "):" + << nl + << " Unable to process internal faces for volume field " + << fieldName << nl << abort(FatalError); + } + + values[i] *= flipMap_[i]; + } + + return true; + } + + return false; +} + + +template +Type Foam::fieldValues::faceSource::processValues +( + const List& values +) const +{ + Type result = pTraits::zero; + switch (operation_) + { + case opSum: + { + result = sum(values); + break; + } + case opAreaAverage: + { + tmp magSf = filterField(mesh().magSf()); + result = sum(values*magSf())/sum(magSf()); + break; + } + case opAreaIntegrate: + { + result = sum(values*filterField(mesh().magSf())); + break; + } + default: + { + // Do nothing + } + } + + return result; +} + + +// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // + +template +bool Foam::fieldValues::faceSource::writeValues(const word& fieldName) +{ + List > allValues(Pstream::nProcs()); + + bool validField = + setFieldValues(fieldName, allValues[Pstream::myProcNo()]); + + if (validField) + { + Pstream::gatherList(allValues); + + if (Pstream::master()) + { + List values = + ListListOps::combine > + ( + allValues, + accessOp >() + ); + + Type result = processValues(values); + + if (valueOutput_) + { + IOList + ( + IOobject + ( + fieldName + "_" + sourceTypeNames_[source_] + "-" + + sourceName_, + obr_.time().timeName(), + obr_, + IOobject::NO_READ, + IOobject::NO_WRITE + ), + values + ).write(); + } + + + outputFilePtr_()<< tab << result; + + if (log_) + { + Info<< " " << operationTypeNames_[operation_] + << "(" << sourceName_ << ") for " << fieldName + << " = " << result << endl; + } + } + } + + return validField; +} + + +template +Foam::tmp > Foam::fieldValues::faceSource::filterField +( + const GeometricField& field +) const +{ + tmp > tvalues(new Field(faceId_.size())); + Field& values = tvalues(); + + forAll(values, i) + { + label faceI = faceId_[i]; + label patchI = facePatchId_[i]; + if (patchI >= 0) + { + values[i] = field.boundaryField()[patchI][faceI]; + } + else + { + FatalErrorIn + ( + "fieldValues::faceSource::filterField" + "(" + "const GeometricField&" + ") const" + ) << type() << " " << name_ << ": " + << sourceTypeNames_[source_] << "(" << sourceName_ << "):" + << nl + << " Unable to process internal faces for volume field " + << field.name() << nl << abort(FatalError); + } + + values[i] *= flipMap_[i]; + } + + return tvalues; +} + + +template +Foam::tmp > Foam::fieldValues::faceSource::filterField +( + const GeometricField& field +) const +{ + tmp > tvalues(new Field(faceId_.size())); + Field& values = tvalues(); + + forAll(values, i) + { + label faceI = faceId_[i]; + label patchI = facePatchId_[i]; + if (patchI >= 0) + { + values[i] = field.boundaryField()[patchI][faceI]; + } + else + { + values[i] = field[faceI]; + } + + values[i] *= flipMap_[i]; + } + + return tvalues; +} + + +// ************************************************************************* // + diff --git a/src/postProcessing/functionObjects/field/fieldValues/fieldValue/fieldValue.C b/src/postProcessing/functionObjects/field/fieldValues/fieldValue/fieldValue.C new file mode 100644 index 0000000000..36b88d5e81 --- /dev/null +++ b/src/postProcessing/functionObjects/field/fieldValues/fieldValue/fieldValue.C @@ -0,0 +1,177 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2009-2009 OpenCFD Ltd. + \\/ M anipulation | +------------------------------------------------------------------------------- +License + This file is part of OpenFOAM. + + OpenFOAM is free software; you can redistribute it and/or modify it + under the terms of the GNU General Public License as published by the + Free Software Foundation; either version 2 of the License, or (at your + option) any later version. + + OpenFOAM is distributed in the hope that it will be useful, but WITHOUT + ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + for more details. + + You should have received a copy of the GNU General Public License + along with OpenFOAM; if not, write to the Free Software Foundation, + Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA + +\*---------------------------------------------------------------------------*/ + +#include "fieldValue.H" +#include "fvMesh.H" + +// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // + +namespace Foam +{ + defineTypeNameAndDebug(fieldValue, 0); + + defineTemplateTypeNameAndDebug(IOList, 0); + defineTemplateTypeNameAndDebug(IOList, 0); + defineTemplateTypeNameAndDebug(IOList, 0); + defineTemplateTypeNameAndDebug(IOList, 0); +} + + +// * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * // + +void Foam::fieldValue::updateMesh(const mapPolyMesh&) +{ + // Do nothing +} + + +void Foam::fieldValue::movePoints(const Field&) +{ + // Do nothing +} + + +// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // + +Foam::fieldValue::fieldValue +( + const word& name, + const objectRegistry& obr, + const dictionary& dict, + const bool loadFromFiles +) +: + name_(name), + obr_(obr), + active_(true), + log_(false), + sourceName_(dict.lookup("sourceName")), + fields_(dict.lookup("fields")), + valueOutput_(dict.lookup("valueOutput")) +{ + // Only active if obr is an fvMesh + if (isA(obr_)) + { + read(dict); + } + else + { + WarningIn + ( + "fieldValue::fieldValue" + "(" + "const word&, " + "const objectRegistry&, " + "const dictionary&, " + "const bool" + ")" + ) << "No fvMesh available, deactivating." + << nl << endl; + active_ = false; + } +} + + +// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * // + +Foam::fieldValue::~fieldValue() +{} + + +// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // + +const Foam::word& Foam::fieldValue::name() const +{ + return name_; +} + + +const Foam::objectRegistry& Foam::fieldValue::obr() const +{ + return obr_; +} + + +bool Foam::fieldValue::active() const +{ + return active_; +} + + +const Foam::Switch& Foam::fieldValue::log() const +{ + return log_; +} + + +const Foam::word& Foam::fieldValue::sourceName() const +{ + return sourceName_; +} + + +const Foam::wordList& Foam::fieldValue::fields() const +{ + return fields_; +} + + +const Foam::Switch& Foam::fieldValue::valueOutput() const +{ + return valueOutput_; +} + + +const Foam::fvMesh& Foam::fieldValue::mesh() const +{ + return refCast(obr_); +} + + +void Foam::fieldValue::read(const dictionary& dict) +{ + if (active_) + { + log_ = dict.lookupOrDefault("log", false); + dict.lookup("fields") >> fields_; + dict.lookup("valueOutput") >> valueOutput_; + } +} + + +void Foam::fieldValue::execute() +{ + // Do nothing +} + + +void Foam::fieldValue::end() +{ + // Do nothing +} + + +// ************************************************************************* // diff --git a/src/postProcessing/functionObjects/field/fieldValues/fieldValue/fieldValue.H b/src/postProcessing/functionObjects/field/fieldValues/fieldValue/fieldValue.H new file mode 100644 index 0000000000..cdedaae131 --- /dev/null +++ b/src/postProcessing/functionObjects/field/fieldValues/fieldValue/fieldValue.H @@ -0,0 +1,165 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2009-2009 OpenCFD Ltd. + \\/ M anipulation | +------------------------------------------------------------------------------- +License + This file is part of OpenFOAM. + + OpenFOAM is free software; you can redistribute it and/or modify it + under the terms of the GNU General Public License as published by the + Free Software Foundation; either version 2 of the License, or (at your + option) any later version. + + OpenFOAM is distributed in the hope that it will be useful, but WITHOUT + ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + for more details. + + You should have received a copy of the GNU General Public License + along with OpenFOAM; if not, write to the Free Software Foundation, + Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA + +Class + Foam::fieldValue + +Description + Base class for field value -based function objects. + +SourceFiles + fieldValue.C + +\*---------------------------------------------------------------------------*/ + +#ifndef fieldValue_H +#define fieldValue_H + +#include "Switch.H" +#include "pointFieldFwd.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ + +// Forward declaration of classes +class dictionary; +class objectRegistry; +class fvMesh; +class mapPolyMesh; + +/*---------------------------------------------------------------------------*\ + Class fieldValue Declaration +\*---------------------------------------------------------------------------*/ + +class fieldValue +{ + +protected: + + // Protected data + + //- Name of this fieldValue object + word name_; + + //- Database this class is registered to + const objectRegistry& obr_; + + //- Active flag + bool active_; + + //- Switch to send output to Info as well as to file + Switch log_; + + //- Name of source object + word sourceName_; + + //- List of field names to operate on + wordList fields_; + + //- Output field values flag + Switch valueOutput_; + + + // Functions to be over-ridden from IOoutputFilter class + + //- Update mesh + virtual void updateMesh(const mapPolyMesh&); + + //- Move points + virtual void movePoints(const Field&); + + +public: + + //- Run-time type information + TypeName("fieldValue"); + + + //- Construct from components + fieldValue + ( + const word& name, + const objectRegistry& obr, + const dictionary& dict, + const bool loadFromFiles = false + ); + + + //- Destructor + virtual ~fieldValue(); + + + // Public member functions + + // Access + + //- Return the name of the geometric source + const word& name() const; + + //- Return the reference to the object registry + const objectRegistry& obr() const; + + //- Return the active flag + bool active() const; + + //- Return the switch to send output to Info as well as to file + const Switch& log() const; + + //- Return the source name + const word& sourceName() const; + + //- Return the list of field names + const wordList& fields() const; + + //- Return the output field values flag + const Switch& valueOutput() const; + + //- Helper funvction to return the reference to the mesh + const fvMesh& mesh() const; + + + // Function object functions + + //- Read from dictionary + virtual void read(const dictionary& dict); + + //- Execute + virtual void execute(); + + //- Execute the at the final time-loop, currently does nothing + virtual void end(); +}; + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace Foam + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#endif + +// ************************************************************************* // diff --git a/src/postProcessing/functionObjects/field/streamLine/streamLine.C b/src/postProcessing/functionObjects/field/streamLine/streamLine.C new file mode 100644 index 0000000000..f3aa202123 --- /dev/null +++ b/src/postProcessing/functionObjects/field/streamLine/streamLine.C @@ -0,0 +1,622 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 1991-2009 OpenCFD Ltd. + \\/ M anipulation | +------------------------------------------------------------------------------- +License + This file is part of OpenFOAM. + + OpenFOAM is free software; you can redistribute it and/or modify it + under the terms of the GNU General Public License as published by the + Free Software Foundation; either version 2 of the License, or (at your + option) any later version. + + OpenFOAM is distributed in the hope that it will be useful, but WITHOUT + ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + for more details. + + You should have received a copy of the GNU General Public License + along with OpenFOAM; if not, write to the Free Software Foundation, + Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA + +\*---------------------------------------------------------------------------*/ + +#include "Pstream.H" +#include "functionObjectList.H" +#include "streamLine.H" +#include "fvMesh.H" +#include "streamLineParticleCloud.H" +#include "ReadFields.H" +#include "meshSearch.H" +#include "sampledSet.H" +#include "globalIndex.H" +#include "mapDistribute.H" + +// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // + +defineTypeNameAndDebug(Foam::streamLine, 0); + +// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * // + +void Foam::streamLine::track() +{ + const Time& runTime = const_cast(obr_.time()); + const fvMesh& mesh = dynamic_cast(obr_); + + IDLList initialParticles; + streamLineParticleCloud particles + ( + mesh, + cloudName_, + initialParticles + ); + + //Pout<< "Seeding particles." << endl; + + const sampledSet& seedPoints = sampledSetPtr_(); + forAll(seedPoints, i) + { + //Pout<< "Seeded particle at " << seedPoints[i] + // << " at cell:" << seedPoints.cells()[i] + // << endl; + + particles.addParticle + ( + new streamLineParticle + ( + particles, + seedPoints[i], + seedPoints.cells()[i], + lifeTime_ // lifetime + ) + ); + } + label nSeeds = returnReduce(particles.size(), sumOp