/*---------------------------------------------------------------------------*\ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | \\ / A nd | Copyright (C) 1991-2007 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 "PointEdgeWave.H" #include "pointMesh.H" #include "processorPolyPatch.H" #include "cyclicPolyPatch.H" #include "OPstream.H" #include "IPstream.H" #include "PstreamCombineReduceOps.H" #include "debug.H" #include "typeInfo.H" #include "ListOps.H" // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // template Foam::scalar Foam::PointEdgeWave::propagationTol_ = 0.01; // Offset labelList. Used for transferring from one cyclic half to the other. template void Foam::PointEdgeWave::offset(const label val, labelList& elems) { forAll(elems, i) { elems[i] += val; } } // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * // // Gets point-point correspondence. Is // - list of halfA points (in cyclic patch points) // - list of halfB points (can overlap with A!) // - for every patchPoint its corresponding point template void Foam::PointEdgeWave::calcCyclicAddressing() { const polyMesh& mesh = pMesh_(); label cycHalf = 0; forAll(mesh.boundaryMesh(), patchI) { const polyPatch& patch = mesh.boundaryMesh()[patchI]; if (isA(patch)) { label halfSize = patch.size()/2; SubList halfAFaces ( mesh.faces(), halfSize, patch.start() ); cycHalves_.set ( cycHalf++, new primitivePatch(halfAFaces, mesh.points()) ); SubList halfBFaces ( mesh.faces(), halfSize, patch.start() + halfSize ); cycHalves_.set ( cycHalf++, new primitivePatch(halfBFaces, mesh.points()) ); } } } // Handle leaving domain. Implementation referred to Type template void Foam::PointEdgeWave::leaveDomain ( const polyPatch& meshPatch, const primitivePatch& patch, const labelList& patchPointLabels, List& pointInfo ) const { const labelList& meshPoints = patch.meshPoints(); forAll(patchPointLabels, i) { label patchPointI = patchPointLabels[i]; const point& pt = patch.points()[meshPoints[patchPointI]]; pointInfo[i].leaveDomain(meshPatch, patchPointI, pt); } } // Handle entering domain. Implementation referred to Type template void Foam::PointEdgeWave::enterDomain ( const polyPatch& meshPatch, const primitivePatch& patch, const labelList& patchPointLabels, List& pointInfo ) const { const labelList& meshPoints = patch.meshPoints(); forAll(patchPointLabels, i) { label patchPointI = patchPointLabels[i]; const point& pt = patch.points()[meshPoints[patchPointI]]; pointInfo[i].enterDomain(meshPatch, patchPointI, pt); } } // Transform. Implementation referred to Type template void Foam::PointEdgeWave::transform ( const tensorField& rotTensor, List& pointInfo ) const { if (rotTensor.size() == 1) { const tensor& T = rotTensor[0]; forAll(pointInfo, i) { pointInfo[i].transform(T); } } else { FatalErrorIn ( "PointEdgeWave::transform(const tensorField&, List&)" ) << "Parallel cyclics not supported" << abort(FatalError); forAll(pointInfo, i) { pointInfo[i].transform(rotTensor[i]); } } } // Update info for pointI, at position pt, with information from // neighbouring edge. // Updates: // - changedPoint_, changedPoints_, nChangedPoints_, // - statistics: nEvals_, nUnvisitedPoints_ template bool Foam::PointEdgeWave::updatePoint ( const label pointI, const label neighbourEdgeI, const Type& neighbourInfo, const scalar tol, Type& pointInfo ) { nEvals_++; bool wasValid = pointInfo.valid(); bool propagate = pointInfo.updatePoint ( pMesh_(), pointI, neighbourEdgeI, neighbourInfo, tol ); if (propagate) { if (!changedPoint_[pointI]) { changedPoint_[pointI] = true; changedPoints_[nChangedPoints_++] = pointI; } } if (!wasValid && pointInfo.valid()) { --nUnvisitedPoints_; } return propagate; } // Update info for pointI, at position pt, with information from // same point. // Updates: // - changedPoint_, changedPoints_, nChangedPoints_, // - statistics: nEvals_, nUnvisitedPoints_ template bool Foam::PointEdgeWave::updatePoint ( const label pointI, const Type& neighbourInfo, const scalar tol, Type& pointInfo ) { nEvals_++; bool wasValid = pointInfo.valid(); bool propagate = pointInfo.updatePoint ( pMesh_(), pointI, neighbourInfo, tol ); if (propagate) { if (!changedPoint_[pointI]) { changedPoint_[pointI] = true; changedPoints_[nChangedPoints_++] = pointI; } } if (!wasValid && pointInfo.valid()) { --nUnvisitedPoints_; } return propagate; } // Update info for edgeI, at position pt, with information from // neighbouring point. // Updates: // - changedEdge_, changedEdges_, nChangedEdges_, // - statistics: nEvals_, nUnvisitedEdge_ template bool Foam::PointEdgeWave::updateEdge ( const label edgeI, const label neighbourPointI, const Type& neighbourInfo, const scalar tol, Type& edgeInfo ) { nEvals_++; bool wasValid = edgeInfo.valid(); bool propagate = edgeInfo.updateEdge ( pMesh_(), edgeI, neighbourPointI, neighbourInfo, tol ); if (propagate) { if (!changedEdge_[edgeI]) { changedEdge_[edgeI] = true; changedEdges_[nChangedEdges_++] = edgeI; } } if (!wasValid && edgeInfo.valid()) { --nUnvisitedEdges_; } return propagate; } // Check if patches of given type name are present template template Foam::label Foam::PointEdgeWave::countPatchType() const { label nPatches = 0; forAll(pMesh_().boundaryMesh(), patchI) { if (isA(pMesh_().boundaryMesh()[patchI])) { nPatches++; } } return nPatches; } // Collect changed patch points template void Foam::PointEdgeWave::getChangedPatchPoints ( const primitivePatch& patch, DynamicList& patchInfo, DynamicList