From f644d9d27719b7db3ed2640c1483fcdedc847707 Mon Sep 17 00:00:00 2001 From: mattijs Date: Thu, 26 Feb 2009 09:02:13 +0000 Subject: [PATCH 01/17] handling empty --- .../sampledSurface/isoSurface/isoSurface.C | 143 ++++++++---------- .../isoSurface/isoSurfaceTemplates.C | 50 +++++- 2 files changed, 108 insertions(+), 85 deletions(-) diff --git a/src/sampling/sampledSurface/isoSurface/isoSurface.C b/src/sampling/sampledSurface/isoSurface/isoSurface.C index 939d77cf1b..78081050ab 100644 --- a/src/sampling/sampledSurface/isoSurface/isoSurface.C +++ b/src/sampling/sampledSurface/isoSurface/isoSurface.C @@ -199,81 +199,47 @@ void Foam::isoSurface::calcCutTypes } } } + forAll(patches, patchI) { const polyPatch& pp = patches[patchI]; label faceI = pp.start(); - if (isA(pp)) + forAll(pp, i) { - // Still needs special treatment? + bool ownLower = (cVals[own[faceI]] < iso_); - forAll(pp, i) + scalar nbrValue; + point nbrPoint; + getNeighbour + ( + boundaryRegion, + cVals, + own[faceI], + faceI, + nbrValue, + nbrPoint + ); + + bool neiLower = (nbrValue < iso_); + + if (ownLower != neiLower) { - bool ownLower = (cVals[own[faceI]] < iso_); - - scalar nbrValue; - point nbrPoint; - getNeighbour - ( - boundaryRegion, - cVals, - own[faceI], - faceI, - nbrValue, - nbrPoint - ); - - bool neiLower = (nbrValue < iso_); - + faceCutType_[faceI] = CUT; + } + else + { + // Mesh edge. const face f = mesh_.faces()[faceI]; if (isEdgeOfFaceCut(pVals, f, ownLower, neiLower)) { faceCutType_[faceI] = CUT; } - - faceI++; } - } - else - { - forAll(pp, i) - { - bool ownLower = (cVals[own[faceI]] < iso_); - scalar nbrValue; - point nbrPoint; - getNeighbour - ( - boundaryRegion, - cVals, - own[faceI], - faceI, - nbrValue, - nbrPoint - ); - - bool neiLower = (nbrValue < iso_); - - if (ownLower != neiLower) - { - faceCutType_[faceI] = CUT; - } - else - { - // Mesh edge. - const face f = mesh_.faces()[faceI]; - - if (isEdgeOfFaceCut(pVals, f, ownLower, neiLower)) - { - faceCutType_[faceI] = CUT; - } - } - - faceI++; - } + faceI++; } } @@ -1589,26 +1555,6 @@ Foam::isoSurface::isoSurface const polyBoundaryMesh& patches = mesh_.boundaryMesh(); - // Check - forAll(patches, patchI) - { - if (isA(patches[patchI])) - { - FatalErrorIn - ( - "isoSurface::isoSurface\n" - "(\n" - " const volScalarField& cVals,\n" - " const scalarField& pVals,\n" - " const scalar iso,\n" - " const bool regularise,\n" - " const scalar mergeTol\n" - ")\n" - ) << "Iso surfaces not supported on case with empty patches." - << exit(FatalError); - } - } - // Pre-calculate patch-per-face to avoid whichPatch call. labelList boundaryRegion(mesh_.nFaces()-mesh_.nInternalFaces()); @@ -1724,7 +1670,7 @@ Foam::isoSurface::isoSurface // Generate field to interpolate. This is identical to the mesh.C() - // except on separated coupled patches. + // except on separated coupled patches and on empty patches. slicedVolVectorField meshC ( IOobject @@ -1746,10 +1692,12 @@ Foam::isoSurface::isoSurface const polyBoundaryMesh& patches = mesh_.boundaryMesh(); forAll(patches, patchI) { + const polyPatch& pp = patches[patchI]; + if ( - patches[patchI].coupled() - && refCast(patches[patchI]).separated() + pp.coupled() + && refCast(pp).separated() ) { fvPatchVectorField& pfld = const_cast @@ -1758,9 +1706,32 @@ Foam::isoSurface::isoSurface ); pfld.operator== ( - patches[patchI].patchSlice(mesh_.faceCentres()) + pp.patchSlice(mesh_.faceCentres()) ); } + else if (isA(pp)) + { + typedef slicedVolVectorField::GeometricBoundaryField bType; + + bType& bfld = const_cast(meshC.boundaryField()); + + // Clear old value. Cannot resize it since slice. + bfld.set(patchI, NULL); + + // Set new value we can change + bfld.set + ( + patchI, + new calculatedFvPatchField + ( + mesh_.boundary()[patchI], + meshC + ) + ); + + // Change to face centres + bfld[patchI] = pp.patchSlice(mesh_.faceCentres()); + } } } @@ -1885,6 +1856,14 @@ Foam::isoSurface::isoSurface orientSurface(*this, faceEdges, edgeFace0, edgeFace1, edgeFacesRest); //} + + + if (debug) + { + fileName stlFile = mesh_.time().path() + ".stl"; + Pout<< "Dumping surface to " << stlFile << endl; + triSurface::write(stlFile); + } } diff --git a/src/sampling/sampledSurface/isoSurface/isoSurfaceTemplates.C b/src/sampling/sampledSurface/isoSurface/isoSurfaceTemplates.C index 7e7206edcc..0760b06a52 100644 --- a/src/sampling/sampledSurface/isoSurface/isoSurfaceTemplates.C +++ b/src/sampling/sampledSurface/isoSurface/isoSurfaceTemplates.C @@ -514,9 +514,53 @@ void Foam::isoSurface::generateTriPoints } else if (isA(pp)) { - // Assume zero-gradient. But what about coordinates? + // Check if field is there (when generating geometry the + // empty patches have been rewritten to be the face centres), + // otherwise use zero-gradient. + label faceI = pp.start(); + + const fvPatchScalarField& fvp = cVals.boundaryField()[patchI]; + + // Owner value of cVals + scalarField internalVals; + if (fvp.size() == 0) + { + internalVals.setSize(pp.size()); + forAll(pp, i) + { + internalVals[i] = cVals[own[pp.start()+i]]; + } + } + const scalarField& bVals = + ( + fvp.size() > 0 + ? static_cast(fvp) + : internalVals + ); + + + const fvPatchField& pc = cCoords.boundaryField()[patchI]; + + // Owner value of cCoords + Field internalCoords; + if (pc.size() == 0) + { + internalCoords.setSize(pp.size()); + forAll(pp, i) + { + internalCoords[i] = cCoords[own[pp.start()+i]]; + } + } + const Field& bCoords = + ( + pc.size() > 0 + ? static_cast&>(pc) + : internalCoords + ); + + forAll(pp, i) { if (faceCutType_[faceI] != NOTCUT) @@ -534,8 +578,8 @@ void Foam::isoSurface::generateTriPoints snappedPoint, faceI, - cVals[own[faceI]], - cCoords.boundaryField()[patchI][i], + bVals[i], + bCoords[i], false, // fc not snapped pTraits::zero, From 0d899d924af19366b2badc817fcec8a30202bb4b Mon Sep 17 00:00:00 2001 From: mattijs Date: Thu, 26 Feb 2009 17:59:12 +0000 Subject: [PATCH 02/17] instance searching --- .../distributedTriSurfaceMesh.C | 79 +++++++++++++++---- .../distributedTriSurfaceMesh.H | 5 +- .../searchableSurface/searchableSphere.C | 9 ++- .../searchableSurface/triSurfaceMesh.C | 52 ++++++------ .../searchableSurface/triSurfaceMesh.H | 1 - 5 files changed, 100 insertions(+), 46 deletions(-) diff --git a/src/meshTools/searchableSurface/distributedTriSurfaceMesh.C b/src/meshTools/searchableSurface/distributedTriSurfaceMesh.C index 870f7e8324..501998cafc 100644 --- a/src/meshTools/searchableSurface/distributedTriSurfaceMesh.C +++ b/src/meshTools/searchableSurface/distributedTriSurfaceMesh.C @@ -1341,27 +1341,53 @@ Foam::distributedTriSurfaceMesh::distributedTriSurfaceMesh dict_(io, dict) { read(); + + if (debug) + { + Info<< "Constructed from triSurface:" << endl; + writeStats(Info); + } } Foam::distributedTriSurfaceMesh::distributedTriSurfaceMesh(const IOobject& io) : - triSurfaceMesh(io), + //triSurfaceMesh(io), + triSurfaceMesh + ( + IOobject + ( + io.name(), + io.time().findInstance(io.local(), word::null), + io.local(), + io.db(), + io.readOpt(), + io.writeOpt(), + io.registerObject() + ) + ), dict_ ( IOobject ( - searchableSurface::name() + "Dict", - searchableSurface::instance(), - searchableSurface::local(), - searchableSurface::db(), - searchableSurface::readOpt(), - searchableSurface::writeOpt(), - searchableSurface::registerObject() + triSurfaceMesh::name() + "Dict", + triSurfaceMesh::instance(), + triSurfaceMesh::local(), + triSurfaceMesh::db(), + triSurfaceMesh::readOpt(), + triSurfaceMesh::writeOpt(), + triSurfaceMesh::registerObject() ) ) { read(); + + if (debug) + { + Info<< "Read distributedTriSurface from " << io.objectPath() + << ':' << endl; + writeStats(Info); + } } @@ -1371,22 +1397,43 @@ Foam::distributedTriSurfaceMesh::distributedTriSurfaceMesh const dictionary& dict ) : - triSurfaceMesh(io, dict), + //triSurfaceMesh(io, dict), + triSurfaceMesh + ( + IOobject + ( + io.name(), + io.time().findInstance(io.local(), word::null), + io.local(), + io.db(), + io.readOpt(), + io.writeOpt(), + io.registerObject() + ), + dict + ), dict_ ( IOobject ( - searchableSurface::name() + "Dict", - searchableSurface::instance(), - searchableSurface::local(), - searchableSurface::db(), - searchableSurface::readOpt(), - searchableSurface::writeOpt(), - searchableSurface::registerObject() + triSurfaceMesh::name() + "Dict", + triSurfaceMesh::instance(), + triSurfaceMesh::local(), + triSurfaceMesh::db(), + triSurfaceMesh::readOpt(), + triSurfaceMesh::writeOpt(), + triSurfaceMesh::registerObject() ) ) { read(); + + if (debug) + { + Info<< "Read distributedTriSurface from " << io.objectPath() + << " and dictionary:" << endl; + writeStats(Info); + } } diff --git a/src/meshTools/searchableSurface/distributedTriSurfaceMesh.H b/src/meshTools/searchableSurface/distributedTriSurfaceMesh.H index 733cd24891..678f1199d4 100644 --- a/src/meshTools/searchableSurface/distributedTriSurfaceMesh.H +++ b/src/meshTools/searchableSurface/distributedTriSurfaceMesh.H @@ -317,10 +317,11 @@ public: const dictionary& dict ); - //- Construct read + //- Construct read. Does findInstance to find io.local(). distributedTriSurfaceMesh(const IOobject& io); - //- Construct from dictionary (used by searchableSurface) + //- Construct from dictionary (used by searchableSurface). + // Does read. Does findInstance to find io.local(). distributedTriSurfaceMesh ( const IOobject& io, diff --git a/src/meshTools/searchableSurface/searchableSphere.C b/src/meshTools/searchableSurface/searchableSphere.C index ba5a1888cf..287ae6743c 100644 --- a/src/meshTools/searchableSurface/searchableSphere.C +++ b/src/meshTools/searchableSurface/searchableSphere.C @@ -53,7 +53,14 @@ Foam::pointIndexHit Foam::searchableSphere::findNearest if (nearestDistSqr > sqr(magN-radius_)) { - info.rawPoint() = centre_ + n/magN*radius_; + if (magN < ROOTVSMALL) + { + info.rawPoint() = centre_ + vector(1,0,0)/magN*radius_; + } + else + { + info.rawPoint() = centre_ + n/magN*radius_; + } info.setHit(); info.setIndex(0); } diff --git a/src/meshTools/searchableSurface/triSurfaceMesh.C b/src/meshTools/searchableSurface/triSurfaceMesh.C index 121762d7c4..9958afd1b1 100644 --- a/src/meshTools/searchableSurface/triSurfaceMesh.C +++ b/src/meshTools/searchableSurface/triSurfaceMesh.C @@ -228,19 +228,19 @@ Foam::triSurfaceMesh::triSurfaceMesh(const IOobject& io, const triSurface& s) Foam::triSurfaceMesh::triSurfaceMesh(const IOobject& io) : // Find instance for triSurfaceMesh - searchableSurface - ( - IOobject - ( - io.name(), - io.time().findInstance(io.local(), word::null), - io.local(), - io.db(), - io.readOpt(), - io.writeOpt(), - io.registerObject() - ) - ), + searchableSurface(io), + //( + // IOobject + // ( + // io.name(), + // io.time().findInstance(io.local(), word::null), + // io.local(), + // io.db(), + // io.readOpt(), + // io.writeOpt(), + // io.registerObject() + // ) + //), // Reused found instance in objectRegistry objectRegistry ( @@ -273,19 +273,19 @@ Foam::triSurfaceMesh::triSurfaceMesh const dictionary& dict ) : - searchableSurface - ( - IOobject - ( - io.name(), - io.time().findInstance(io.local(), word::null), - io.local(), - io.db(), - io.readOpt(), - io.writeOpt(), - io.registerObject() - ) - ), + searchableSurface(io), + //( + // IOobject + // ( + // io.name(), + // io.time().findInstance(io.local(), word::null), + // io.local(), + // io.db(), + // io.readOpt(), + // io.writeOpt(), + // io.registerObject() + // ) + //), // Reused found instance in objectRegistry objectRegistry ( diff --git a/src/meshTools/searchableSurface/triSurfaceMesh.H b/src/meshTools/searchableSurface/triSurfaceMesh.H index d7aa33cfca..4f98f37064 100644 --- a/src/meshTools/searchableSurface/triSurfaceMesh.H +++ b/src/meshTools/searchableSurface/triSurfaceMesh.H @@ -117,7 +117,6 @@ public: triSurfaceMesh(const IOobject& io); //- Construct from IO and dictionary (used by searchableSurface). - // Does timeInstance search. // Dictionary may contain a 'scale' entry (eg, 0.001: mm -> m) triSurfaceMesh ( From 65baa26a804e46ca42ad93ab6a45a6c27f0d6110 Mon Sep 17 00:00:00 2001 From: mattijs Date: Thu, 26 Feb 2009 20:14:20 +0000 Subject: [PATCH 03/17] mergeDist never set. --- .../distributedTriSurfaceMesh.C | 75 ++++++++++++++----- .../distributedTriSurfaceMesh.H | 5 +- 2 files changed, 59 insertions(+), 21 deletions(-) diff --git a/src/meshTools/searchableSurface/distributedTriSurfaceMesh.C b/src/meshTools/searchableSurface/distributedTriSurfaceMesh.C index 501998cafc..7f86c0f658 100644 --- a/src/meshTools/searchableSurface/distributedTriSurfaceMesh.C +++ b/src/meshTools/searchableSurface/distributedTriSurfaceMesh.C @@ -76,10 +76,8 @@ bool Foam::distributedTriSurfaceMesh::read() // Distribution type distType_ = distributionTypeNames_.read(dict_.lookup("distributionType")); - if (dict_.found("mergeDistance")) - { - dict_.lookup("mergeDistance") >> mergeDist_; - } + // Merge distance + mergeDist_ = readScalar(dict_.lookup("mergeDistance")); return true; } @@ -1346,6 +1344,19 @@ Foam::distributedTriSurfaceMesh::distributedTriSurfaceMesh { Info<< "Constructed from triSurface:" << endl; writeStats(Info); + + labelList nTris(Pstream::nProcs()); + nTris[Pstream::myProcNo()] = triSurface::size(); + Pstream::gatherList(nTris); + Pstream::scatterList(nTris); + + Info<< endl<< "\tproc\ttris\tbb" << endl; + forAll(nTris, procI) + { + Info<< '\t' << procI << '\t' << nTris[procI] + << '\t' << procBb_[procI] << endl; + } + Info<< endl; } } @@ -1370,13 +1381,13 @@ Foam::distributedTriSurfaceMesh::distributedTriSurfaceMesh(const IOobject& io) ( IOobject ( - triSurfaceMesh::name() + "Dict", - triSurfaceMesh::instance(), - triSurfaceMesh::local(), - triSurfaceMesh::db(), - triSurfaceMesh::readOpt(), - triSurfaceMesh::writeOpt(), - triSurfaceMesh::registerObject() + searchableSurface::name() + "Dict", + searchableSurface::instance(), + searchableSurface::local(), + searchableSurface::db(), + searchableSurface::readOpt(), + searchableSurface::writeOpt(), + searchableSurface::registerObject() ) ) { @@ -1387,6 +1398,19 @@ Foam::distributedTriSurfaceMesh::distributedTriSurfaceMesh(const IOobject& io) Info<< "Read distributedTriSurface from " << io.objectPath() << ':' << endl; writeStats(Info); + + labelList nTris(Pstream::nProcs()); + nTris[Pstream::myProcNo()] = triSurface::size(); + Pstream::gatherList(nTris); + Pstream::scatterList(nTris); + + Info<< endl<< "\tproc\ttris\tbb" << endl; + forAll(nTris, procI) + { + Info<< '\t' << procI << '\t' << nTris[procI] + << '\t' << procBb_[procI] << endl; + } + Info<< endl; } } @@ -1416,13 +1440,13 @@ Foam::distributedTriSurfaceMesh::distributedTriSurfaceMesh ( IOobject ( - triSurfaceMesh::name() + "Dict", - triSurfaceMesh::instance(), - triSurfaceMesh::local(), - triSurfaceMesh::db(), - triSurfaceMesh::readOpt(), - triSurfaceMesh::writeOpt(), - triSurfaceMesh::registerObject() + searchableSurface::name() + "Dict", + searchableSurface::instance(), + searchableSurface::local(), + searchableSurface::db(), + searchableSurface::readOpt(), + searchableSurface::writeOpt(), + searchableSurface::registerObject() ) ) { @@ -1433,6 +1457,19 @@ Foam::distributedTriSurfaceMesh::distributedTriSurfaceMesh Info<< "Read distributedTriSurface from " << io.objectPath() << " and dictionary:" << endl; writeStats(Info); + + labelList nTris(Pstream::nProcs()); + nTris[Pstream::myProcNo()] = triSurface::size(); + Pstream::gatherList(nTris); + Pstream::scatterList(nTris); + + Info<< endl<< "\tproc\ttris\tbb" << endl; + forAll(nTris, procI) + { + Info<< '\t' << procI << '\t' << nTris[procI] + << '\t' << procBb_[procI] << endl; + } + Info<< endl; } } @@ -2436,6 +2473,8 @@ void Foam::distributedTriSurfaceMesh::writeStats(Ostream& os) const boundBox bb; label nPoints; calcBounds(bb, nPoints); + reduce(bb.min(), minOp()); + reduce(bb.max(), maxOp()); os << "Triangles : " << returnReduce(triSurface::size(), sumOp