structuredRenumber: Corrected to run in parallel

Patch contributed by Mattijs Janssens
This commit is contained in:
Henry Weller
2016-09-09 12:20:25 +01:00
parent 8e0981d9ed
commit cb1a012b66
5 changed files with 651 additions and 48 deletions

View File

@ -0,0 +1,331 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2016 OpenFOAM Foundation
\\/ 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 3 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, see <http://www.gnu.org/licenses/>.
\*---------------------------------------------------------------------------*/
#include "OppositeFaceCellWave.H"
#include "polyMesh.H"
// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
template<class Type, class TrackingData>
void Foam::OppositeFaceCellWave<Type, TrackingData>::opposingFaceLabels
(
const label celli,
const label masterFaceLabel,
DynamicList<label>& oppositeFaceLabels
) const
{
// Variant of cell::opposingFaceLabel
// Algorithm:
// Go through all the faces of the cell and find the one which
// does not share a single vertex with the master face. If there
// are two or more such faces, return the first one and issue a
// warning; if there is no opposite face, return -1;
const face& masterFace = this->mesh_.faces()[masterFaceLabel];
const labelList& curFaceLabels = this->mesh_.cells()[celli];
oppositeFaceLabels.clear();
forAll(curFaceLabels, facei)
{
// Compare the face with the master
const face& curFace = this->mesh_.faces()[curFaceLabels[facei]];
// Skip the master face
if (curFaceLabels[facei] != masterFaceLabel)
{
bool sharedPoint = false;
// Compare every vertex of the current face against the
// vertices of the master face
forAll(curFace, pointi)
{
const label l = curFace[pointi];
forAll(masterFace, masterPointi)
{
if (masterFace[masterPointi] == l)
{
sharedPoint = true;
break;
}
}
if (sharedPoint) break;
}
// If no points are shared, this is the opposite face
if (!sharedPoint)
{
// Found opposite face
oppositeFaceLabels.append(curFaceLabels[facei]);
}
}
}
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
// Iterate, propagating changedFacesInfo across mesh, until no change (or
// maxIter reached). Initial cell values specified.
template<class Type, class TrackingData>
Foam::OppositeFaceCellWave<Type, TrackingData>::OppositeFaceCellWave
(
const polyMesh& mesh,
const labelList& changedFaces,
const List<Type>& changedFacesInfo,
UList<Type>& allFaceInfo,
UList<Type>& allCellInfo,
const label maxIter,
TrackingData& td
)
:
FaceCellWave<Type, TrackingData>
(
mesh,
changedFaces,
changedFacesInfo,
allFaceInfo,
allCellInfo,
0, //maxIter,
td
),
changedOppositeFaces_(this->mesh_.nCells())
{
// Iterate until nothing changes
label iter = this->iterate(maxIter);
if ((maxIter > 0) && (iter >= maxIter))
{
FatalErrorInFunction
<< "Maximum number of iterations reached. Increase maxIter."
<< endl
<< " maxIter:" << maxIter << endl
<< " nChangedCells:" << this->changedCells_.size() << endl
<< " nChangedFaces:" << this->changedFaces_.size() << endl
<< exit(FatalError);
}
}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
template<class Type, class TrackingData>
Foam::label Foam::OppositeFaceCellWave<Type, TrackingData>::faceToCell()
{
const labelList& owner = this->mesh_.faceOwner();
const labelList& neighbour = this->mesh_.faceNeighbour();
label nInternalFaces = this->mesh_.nInternalFaces();
DynamicList<label> oppositeFaceLabels;
forAll(this->changedFaces_, changedFacei)
{
label facei = this->changedFaces_[changedFacei];
if (!this->changedFace_[facei])
{
FatalErrorInFunction
<< "Face " << facei
<< " not marked as having been changed"
<< abort(FatalError);
}
const Type& neighbourWallInfo = this->allFaceInfo_[facei];
// Evaluate all connected cells
// Owner
{
label celli = owner[facei];
Type& currentWallInfo = this->allCellInfo_[celli];
if (!currentWallInfo.equal(neighbourWallInfo, this->td_))
{
// Check if cell is prismatic w.r.t facei
opposingFaceLabels(celli, facei, oppositeFaceLabels);
if (oppositeFaceLabels.size())
{
label sz = this->changedCells_.size();
this->updateCell
(
celli,
facei,
neighbourWallInfo,
this->propagationTol_,
currentWallInfo
);
if (this->changedCells_.size() > sz)
{
label oppFacei = -1;
if (oppositeFaceLabels.size() == 1)
{
oppFacei = oppositeFaceLabels[0];
}
changedOppositeFaces_.append(oppFacei);
}
}
}
}
// Neighbour.
if (facei < nInternalFaces)
{
label celli = neighbour[facei];
Type& currentWallInfo2 = this->allCellInfo_[celli];
if (!currentWallInfo2.equal(neighbourWallInfo, this->td_))
{
// Check if cell is prismatic w.r.t facei
opposingFaceLabels(celli, facei, oppositeFaceLabels);
if (oppositeFaceLabels.size())
{
label sz = this->changedCells_.size();
this->updateCell
(
celli,
facei,
neighbourWallInfo,
this->propagationTol_,
currentWallInfo2
);
if (this->changedCells_.size() > sz)
{
label oppFacei = -1;
if (oppositeFaceLabels.size() == 1)
{
oppFacei = oppositeFaceLabels[0];
}
changedOppositeFaces_.append(oppFacei);
}
}
}
}
// Reset status of face
this->changedFace_[facei] = false;
}
// Handled all changed faces by now
this->changedFaces_.clear();
if (debug & 2)
{
Pout<< " Changed cells : " << this->changedCells_.size()
<< endl;
}
// Sum changedCells over all procs
label totNChanged = this->changedCells_.size();
reduce(totNChanged, sumOp<label>());
return totNChanged;
}
template<class Type, class TrackingData>
Foam::label Foam::OppositeFaceCellWave<Type, TrackingData>::cellToFace()
{
forAll(this->changedCells_, changedCelli)
{
label celli = this->changedCells_[changedCelli];
label facei = changedOppositeFaces_[changedCelli];
if (!this->changedCell_[celli])
{
FatalErrorInFunction
<< "Cell " << celli << " not marked as having been changed"
<< abort(FatalError);
}
if (facei != -1)
{
const Type& neighbourWallInfo = this->allCellInfo_[celli];
// Evaluate facei
Type& currentWallInfo = this->allFaceInfo_[facei];
if (!currentWallInfo.equal(neighbourWallInfo, this->td_))
{
this->updateFace
(
facei,
celli,
neighbourWallInfo,
this->propagationTol_,
currentWallInfo
);
}
}
// Reset status of cell
this->changedCell_[celli] = false;
}
// Handled all changed cells by now
this->changedCells_.clear();
changedOppositeFaces_.clear();
if (this->hasCyclicPatches_)
{
// Transfer changed faces across cyclic halves
this->handleCyclicPatches();
}
if (this->hasCyclicAMIPatches_)
{
this->handleAMICyclicPatches();
}
if (Pstream::parRun())
{
// Transfer changed faces from neighbouring processors.
this->handleProcPatches();
}
if (debug & 2)
{
Pout<< " Changed faces : " << this->changedFaces_.size()
<< endl;
}
// Sum nChangedFaces over all procs
label totNChanged = this->changedFaces_.size();
reduce(totNChanged, sumOp<label>());
return totNChanged;
}
// ************************************************************************* //

View File

@ -0,0 +1,143 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2016 OpenFOAM Foundation
\\/ 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 3 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, see <http://www.gnu.org/licenses/>.
Class
Foam::OppositeFaceCellWave
Description
Version of FaceCellWave that walks through prismatic cells only.
Used to determine mesh structure. In the front walking routines
(faceToCell and faceToCell) it
- walks across prismatic cells only
- and only to a single opposite face
Notes:
A cell with a split faces will be marked but not walked through (since
there is no single opposite face.
SourceFiles
OppositeFaceCellWave.C
\*---------------------------------------------------------------------------*/
#ifndef OppositeFaceCellWave_H
#define OppositeFaceCellWave_H
#include "FaceCellWave.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
/*---------------------------------------------------------------------------*\
Class OppositeFaceCellWaveName Declaration
\*---------------------------------------------------------------------------*/
TemplateName(OppositeFaceCellWave);
/*---------------------------------------------------------------------------*\
Class OppositeFaceCellWave Declaration
\*---------------------------------------------------------------------------*/
template<class Type, class TrackingData = int>
class OppositeFaceCellWave
:
public FaceCellWave<Type, TrackingData>,
public OppositeFaceCellWaveName
{
protected:
// Protected data
//- For every entry in changedCells (i.e. the cell front) gives
// the face that it needs to transfer to
DynamicList<label> changedOppositeFaces_;
// Protected Member Functions
//- Determine 'opposite' faces (= faces not sharing a vertex) on cell
void opposingFaceLabels
(
const label celli,
const label facei,
DynamicList<label>&
) const;
public:
// Constructors
//- Construct from mesh and list of changed faces with the Type
// for these faces. Iterates until nothing changes or maxIter reached.
// (maxIter can be 0)
OppositeFaceCellWave
(
const polyMesh&,
const labelList& initialChangedFaces,
const List<Type>& changedFacesInfo,
UList<Type>& allFaceInfo,
UList<Type>& allCellInfo,
const label maxIter,
TrackingData& td = FaceCellWave<Type, TrackingData>::dummyTrackData_
);
//- Destructor
virtual ~OppositeFaceCellWave()
{};
// Member Functions
//- Propagate from face to cell. Returns total number of cells
// (over all processors) changed.
virtual label faceToCell();
//- Propagate from cell to face. Returns total number of faces
// (over all processors) changed. (Faces on processorpatches are
// counted double)
virtual label cellToFace();
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#ifdef NoRepository
#include "OppositeFaceCellWave.C"
#endif
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //

View File

@ -0,0 +1,36 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2016 OpenFOAM Foundation
\\/ 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 3 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, see <http://www.gnu.org/licenses/>.
\*---------------------------------------------------------------------------*/
#include "OppositeFaceCellWave.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
namespace Foam
{
defineTypeNameAndDebug(OppositeFaceCellWaveName, 0);
}
// ************************************************************************* //

View File

@ -27,7 +27,7 @@ License
#include "addToRunTimeSelectionTable.H"
#include "topoDistanceData.H"
#include "fvMeshSubset.H"
#include "FaceCellWave.H"
#include "OppositeFaceCellWave.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@ -54,7 +54,7 @@ Foam::structuredRenumber::structuredRenumber
renumberMethod(renumberDict),
methodDict_(renumberDict.subDict(typeName + "Coeffs")),
patches_(methodDict_.lookup("patches")),
//nLayers_(readLabel(methodDict_.lookup("nLayers"))),
nLayers_(methodDict_.lookupOrDefault<label>("nLayers", labelMax)),
depthFirst_(methodDict_.lookup("depthFirst")),
method_(renumberMethod::New(methodDict_)),
reverse_(methodDict_.lookup("reverse"))
@ -63,6 +63,75 @@ Foam::structuredRenumber::structuredRenumber
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
bool Foam::structuredRenumber::layerLess::operator()
(
const label a,
const label b
)
{
const topoDistanceData& ta = distance_[a];
const topoDistanceData& tb = distance_[b];
int dummy;
if (ta.valid(dummy))
{
if (tb.valid(dummy))
{
if (depthFirst_)
{
if (ta.data() < tb.data())
{
// Sort column first
return true;
}
else if (ta.data() > tb.data())
{
return false;
}
else
{
// Same column. Sort according to layer
return ta.distance() < tb.distance();
}
}
else
{
if (ta.distance() < tb.distance())
{
return true;
}
else if (ta.distance() > tb.distance())
{
return false;
}
else
{
// Same layer; sort according to current values
return ta.data() < tb.data();
}
}
}
else
{
return true;
}
}
else
{
if (tb.valid(dummy))
{
return false;
}
else
{
// Both not valid; fall back to cell indices for sorting
return order_[a] < order_[b];
}
}
}
Foam::labelList Foam::structuredRenumber::renumber
(
const polyMesh& mesh,
@ -104,18 +173,13 @@ Foam::labelList Foam::structuredRenumber::renumber
const label nLayers = nTotalCells/nTotalSeeds;
Info<< type() << " : seeding " << nTotalSeeds
<< " cells on " << nLayers << " layers" << nl
<< " cells on (estimated) " << nLayers << " layers" << nl
<< endl;
// Avoid subsetMesh, FaceCellWave going through proc boundaries
bool oldParRun = Pstream::parRun();
Pstream::parRun() = false;
// Work array. Used here to temporarily store the original-to-ordered
// index. Later on used to store the ordered-to-original.
labelList orderedToOld(points.size(), -1);
labelList orderedToOld(mesh.nCells(), -1);
// Subset the layer of cells next to the patch
{
@ -125,20 +189,23 @@ Foam::labelList Foam::structuredRenumber::renumber
pointField subPoints(points, subsetter.cellMap());
// Decompose the layer of cells
// Locally renumber the layer of cells
labelList subOrder(method_().renumber(subMesh, subPoints));
labelList subOrigToOrdered(invert(subOrder.size(), subOrder));
// Transfer to final decomposition
globalIndex globalSubCells(subOrder.size());
// Transfer to final decomposition and convert into global numbering
forAll(subOrder, i)
{
orderedToOld[subsetter.cellMap()[i]] = subOrigToOrdered[i];
orderedToOld[subsetter.cellMap()[i]] =
globalSubCells.toGlobal(subOrigToOrdered[i]);
}
}
// Walk out.
// Walk sub-ordering (=column index) out.
labelList patchFaces(nFaces);
List<topoDistanceData> patchData(nFaces);
nFaces = 0;
@ -151,7 +218,7 @@ Foam::labelList Foam::structuredRenumber::renumber
patchFaces[nFaces] = pp.start()+i;
patchData[nFaces] = topoDistanceData
(
orderedToOld[fc[i]],// passive data: order of originating face
orderedToOld[fc[i]],// passive data: global column
0 // distance: layer
);
nFaces++;
@ -163,50 +230,43 @@ Foam::labelList Foam::structuredRenumber::renumber
List<topoDistanceData> faceData(mesh.nFaces());
// Propagate information inwards
FaceCellWave<topoDistanceData> deltaCalc
OppositeFaceCellWave<topoDistanceData> deltaCalc
(
mesh,
patchFaces,
patchData,
faceData,
cellData,
nTotalCells+1
0
);
deltaCalc.iterate(nLayers_);
Pstream::parRun() = oldParRun;
Info<< type() << " : did not visit "
<< deltaCalc.getUnsetCells()
<< " cells out of " << nTotalCells
<< "; using " << method_().type() << " renumbering for these" << endl;
// Get cell order using the method(). These values will get overwitten
// by any visited cell so are used only if the number of nLayers is limited.
labelList oldToOrdered
(
invert
(
mesh.nCells(),
method_().renumber(mesh, points)
)
);
// And extract.
// Note that distance is distance from face so starts at 1.
bool haveWarned = false;
forAll(orderedToOld, celli)
{
if (!cellData[celli].valid(deltaCalc.data()))
{
if (!haveWarned)
{
WarningInFunction
<< "Did not visit some cells, e.g. cell " << celli
<< " at " << mesh.cellCentres()[celli] << endl
<< "Assigning these cells to domain 0." << endl;
haveWarned = true;
}
orderedToOld[celli] = 0;
}
else
{
label layerI = cellData[celli].distance();
if (depthFirst_)
{
orderedToOld[nLayers*cellData[celli].data()+layerI] = celli;
}
else
{
orderedToOld[cellData[celli].data()+nLayers*layerI] = celli;
}
}
}
// Use specialised sorting to sorted either layers or columns first
// Done so that at no point we need to combine both into a single
// index and we might run out of label size.
sortedOrder
(
cellData,
orderedToOld,
layerLess(depthFirst_, oldToOrdered, cellData)
);
// Return furthest away cell first
if (reverse_)

View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2012-2015 OpenFOAM Foundation
\\ / A nd | Copyright (C) 2012-2016 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -41,6 +41,7 @@ SourceFiles
#define structuredRenumber_H
#include "renumberMethod.H"
#include "topoDistanceData.H"
#include "Switch.H"
namespace Foam
@ -54,12 +55,44 @@ class structuredRenumber
:
public renumberMethod
{
public:
// Public classes
//- Less function class that can be used for sorting according to
// column and layer
class layerLess
{
const Switch depthFirst_;
const labelList& order_;
const List<topoDistanceData>& distance_;
public:
layerLess
(
const Switch depthFirst,
const labelList& order,
const List<topoDistanceData>& distance
)
:
depthFirst_(depthFirst),
order_(order),
distance_(distance)
{}
bool operator()(const label a, const label b);
};
// Private data
const dictionary methodDict_;
const wordReList patches_;
const label nLayers_;
const Switch depthFirst_;
const autoPtr<renumberMethod> method_;