ENH: Split src/parallel into decompse and reconstruct to remove cyclic build dependency

This commit is contained in:
andy
2010-03-18 12:23:01 +00:00
parent 0d4349bb7e
commit dc06f24a0c
55 changed files with 32 additions and 25 deletions

View File

@ -0,0 +1,7 @@
decompositionMethod/decompositionMethod.C
geomDecomp/geomDecomp.C
simpleGeomDecomp/simpleGeomDecomp.C
hierarchGeomDecomp/hierarchGeomDecomp.C
manualDecomp/manualDecomp.C
LIB = $(FOAM_LIBBIN)/libdecompositionMethods

View File

@ -0,0 +1,9 @@
EXE_INC =
LIB_LIBS = \
-L$(FOAM_LIBBIN)/dummy \
-L$(FOAM_MPI_LIBBIN) \
-lscotchDecomp \
-lmetisDecomp \
-lparMetisDecomp

View File

@ -0,0 +1,216 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / 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
InClass
decompositionMethod
\*---------------------------------------------------------------------------*/
#include "decompositionMethod.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
namespace Foam
{
defineTypeNameAndDebug(decompositionMethod, 0);
defineRunTimeSelectionTable(decompositionMethod, dictionary);
defineRunTimeSelectionTable(decompositionMethod, dictionaryMesh);
}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
Foam::autoPtr<Foam::decompositionMethod> Foam::decompositionMethod::New
(
const dictionary& decompositionDict
)
{
word decompositionMethodTypeName(decompositionDict.lookup("method"));
Info<< "Selecting decompositionMethod "
<< decompositionMethodTypeName << endl;
dictionaryConstructorTable::iterator cstrIter =
dictionaryConstructorTablePtr_->find(decompositionMethodTypeName);
if (cstrIter == dictionaryConstructorTablePtr_->end())
{
FatalErrorIn
(
"decompositionMethod::New"
"(const dictionary& decompositionDict)"
) << "Unknown decompositionMethod "
<< decompositionMethodTypeName << endl << endl
<< "Valid decompositionMethods are : " << endl
<< dictionaryConstructorTablePtr_->sortedToc()
<< exit(FatalError);
}
return autoPtr<decompositionMethod>(cstrIter()(decompositionDict));
}
Foam::autoPtr<Foam::decompositionMethod> Foam::decompositionMethod::New
(
const dictionary& decompositionDict,
const polyMesh& mesh
)
{
word decompositionMethodTypeName(decompositionDict.lookup("method"));
Info<< "Selecting decompositionMethod "
<< decompositionMethodTypeName << endl;
dictionaryMeshConstructorTable::iterator cstrIter =
dictionaryMeshConstructorTablePtr_->find(decompositionMethodTypeName);
if (cstrIter == dictionaryMeshConstructorTablePtr_->end())
{
FatalErrorIn
(
"decompositionMethod::New"
"(const dictionary& decompositionDict, "
"const polyMesh& mesh)"
) << "Unknown decompositionMethod "
<< decompositionMethodTypeName << endl << endl
<< "Valid decompositionMethods are : " << endl
<< dictionaryMeshConstructorTablePtr_->sortedToc()
<< exit(FatalError);
}
return autoPtr<decompositionMethod>(cstrIter()(decompositionDict, mesh));
}
Foam::labelList Foam::decompositionMethod::decompose
(
const pointField& points
)
{
scalarField weights(0);
return decompose(points, weights);
}
Foam::labelList Foam::decompositionMethod::decompose
(
const labelList& fineToCoarse,
const pointField& coarsePoints,
const scalarField& coarseWeights
)
{
// Decompose based on agglomerated points
labelList coarseDistribution(decompose(coarsePoints, coarseWeights));
// Rework back into decomposition for original mesh_
labelList fineDistribution(fineToCoarse.size());
forAll(fineDistribution, i)
{
fineDistribution[i] = coarseDistribution[fineToCoarse[i]];
}
return fineDistribution;
}
Foam::labelList Foam::decompositionMethod::decompose
(
const labelList& fineToCoarse,
const pointField& coarsePoints
)
{
// Decompose based on agglomerated points
labelList coarseDistribution(decompose(coarsePoints));
// Rework back into decomposition for original mesh_
labelList fineDistribution(fineToCoarse.size());
forAll(fineDistribution, i)
{
fineDistribution[i] = coarseDistribution[fineToCoarse[i]];
}
return fineDistribution;
}
void Foam::decompositionMethod::calcCellCells
(
const polyMesh& mesh,
const labelList& fineToCoarse,
const label nCoarse,
labelListList& cellCells
)
{
if (fineToCoarse.size() != mesh.nCells())
{
FatalErrorIn
(
"decompositionMethod::calcCellCells"
"(const labelList&, labelListList&) const"
) << "Only valid for mesh agglomeration." << exit(FatalError);
}
List<DynamicList<label> > dynCellCells(nCoarse);
forAll(mesh.faceNeighbour(), faceI)
{
label own = fineToCoarse[mesh.faceOwner()[faceI]];
label nei = fineToCoarse[mesh.faceNeighbour()[faceI]];
if (own != nei)
{
if (findIndex(dynCellCells[own], nei) == -1)
{
dynCellCells[own].append(nei);
}
if (findIndex(dynCellCells[nei], own) == -1)
{
dynCellCells[nei].append(own);
}
}
}
cellCells.setSize(dynCellCells.size());
forAll(dynCellCells, coarseI)
{
cellCells[coarseI].transfer(dynCellCells[coarseI]);
}
}
Foam::labelList Foam::decompositionMethod::decompose
(
const labelListList& globalCellCells,
const pointField& cc
)
{
scalarField cWeights(0);
return decompose(globalCellCells, cc, cWeights);
}
// ************************************************************************* //

View File

@ -0,0 +1,216 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / 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
Class
Foam::decompositionMethod
Description
Abstract base class for decomposition
SourceFiles
decompositionMethod.C
\*---------------------------------------------------------------------------*/
#ifndef decompositionMethod_H
#define decompositionMethod_H
#include "polyMesh.H"
#include "pointField.H"
namespace Foam
{
/*---------------------------------------------------------------------------*\
Class decompositionMethod Declaration
\*---------------------------------------------------------------------------*/
class decompositionMethod
{
protected:
// Protected data
const dictionary& decompositionDict_;
label nProcessors_;
//- Helper: determine (non-parallel) cellCells from mesh agglomeration.
static void calcCellCells
(
const polyMesh& mesh,
const labelList& agglom,
const label nCoarse,
labelListList& cellCells
);
private:
// Private Member Functions
//- Disallow default bitwise copy construct and assignment
decompositionMethod(const decompositionMethod&);
void operator=(const decompositionMethod&);
public:
//- Runtime type information
TypeName("decompositionMethod");
// Declare run-time constructor selection tables
declareRunTimeSelectionTable
(
autoPtr,
decompositionMethod,
dictionary,
(
const dictionary& decompositionDict
),
(decompositionDict)
);
declareRunTimeSelectionTable
(
autoPtr,
decompositionMethod,
dictionaryMesh,
(
const dictionary& decompositionDict,
const polyMesh& mesh
),
(decompositionDict, mesh)
);
// Selectors
//- Return a reference to the selected decomposition method
static autoPtr<decompositionMethod> New
(
const dictionary& decompositionDict
);
//- Return a reference to the selected decomposition method
static autoPtr<decompositionMethod> New
(
const dictionary& decompositionDict,
const polyMesh& mesh
);
// Constructors
//- Construct given the decomposition dictionary
decompositionMethod(const dictionary& decompositionDict)
:
decompositionDict_(decompositionDict),
nProcessors_
(
readLabel(decompositionDict.lookup("numberOfSubdomains"))
)
{}
// Destructor
virtual ~decompositionMethod()
{}
// Member Functions
//- Is method parallel aware (i.e. does it synchronize domains across
// proc boundaries)
virtual bool parallelAware() const = 0;
//- Return for every coordinate the wanted processor number. Use the
// mesh connectivity (if needed)
virtual labelList decompose
(
const pointField& points,
const scalarField& pointWeights
) = 0;
//- Like decompose but with uniform weights on the points
virtual labelList decompose(const pointField&);
//- Return for every coordinate the wanted processor number. Gets
// passed agglomeration map (from fine to coarse cells) and coarse cell
// location. Can be overridden by decomposers that provide this
// functionality natively. Coarse cells are local to the processor
// (if in parallel). If you want to have coarse cells spanning
// processors use the globalCellCells instead.
virtual labelList decompose
(
const labelList& cellToRegion,
const pointField& regionPoints,
const scalarField& regionWeights
);
//- Like decompose but with uniform weights on the regions
virtual labelList decompose
(
const labelList& cellToRegion,
const pointField& regionPoints
);
//- Return for every coordinate the wanted processor number. Explicitly
// provided connectivity - does not use mesh_.
// The connectivity is equal to mesh.cellCells() except for
// - in parallel the cell numbers are global cell numbers (starting
// from 0 at processor0 and then incrementing all through the
// processors)
// - the connections are across coupled patches
virtual labelList decompose
(
const labelListList& globalCellCells,
const pointField& cc,
const scalarField& cWeights
) = 0;
//- Like decompose but with uniform weights on the cells
virtual labelList decompose
(
const labelListList& globalCellCells,
const pointField& cc
);
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //

View File

@ -0,0 +1,72 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / 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 "geomDecomp.H"
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
Foam::geomDecomp::geomDecomp
(
const dictionary& decompositionDict,
const word& derivedType
)
:
decompositionMethod(decompositionDict),
geomDecomDict_(decompositionDict.subDict(derivedType + "Coeffs")),
n_(geomDecomDict_.lookup("n")),
delta_(readScalar(geomDecomDict_.lookup("delta"))),
rotDelta_(I)
{
// check that the case makes sense :
if (nProcessors_ != n_.x()*n_.y()*n_.z())
{
FatalErrorIn
(
"geomDecomp::geomDecomp"
"(const dictionary& decompositionDict)"
) << "Wrong number of processor divisions in geomDecomp:" << nl
<< "Number of domains : " << nProcessors_ << nl
<< "Wanted decomposition : " << n_
<< exit(FatalError);
}
scalar d = 1 - 0.5*delta_*delta_;
scalar d2 = sqr(d);
scalar a = delta_;
scalar a2 = sqr(a);
rotDelta_ = tensor
(
d2, -a*d, a,
a*d - a2*d, a*a2 + d2, -2*a*d,
a*d2 + a2, a*d - a2*d, d2 - a2
);
}
// ************************************************************************* //

View File

@ -0,0 +1,86 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / 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
Class
Foam::geomDecomp
Description
Geometrical domain decomposition
SourceFiles
geomDecomp.C
\*---------------------------------------------------------------------------*/
#ifndef geomDecomp_H
#define geomDecomp_H
#include "decompositionMethod.H"
#include "Vector.H"
namespace Foam
{
/*---------------------------------------------------------------------------*\
Class geomDecomp Declaration
\*---------------------------------------------------------------------------*/
class geomDecomp
:
public decompositionMethod
{
protected:
// Protected data
const dictionary& geomDecomDict_;
Vector<label> n_;
scalar delta_;
tensor rotDelta_;
public:
// Constructors
//- Construct given the decomposition dictionary
// and the derived type name
geomDecomp
(
const dictionary& decompositionDict,
const word& derivedType
);
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //

View File

@ -0,0 +1,801 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / 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 "hierarchGeomDecomp.H"
#include "addToRunTimeSelectionTable.H"
#include "PstreamReduceOps.H"
#include "SortableList.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
defineTypeNameAndDebug(hierarchGeomDecomp, 0);
addToRunTimeSelectionTable
(
decompositionMethod,
hierarchGeomDecomp,
dictionary
);
addToRunTimeSelectionTable
(
decompositionMethod,
hierarchGeomDecomp,
dictionaryMesh
);
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
void Foam::hierarchGeomDecomp::setDecompOrder()
{
word order(geomDecomDict_.lookup("order"));
if (order.size() != 3)
{
FatalIOErrorIn
(
"hierarchGeomDecomp::hierarchGeomDecomp"
"(const dictionary& decompositionDict)",
decompositionDict_
) << "number of characters in order (" << order << ") != 3"
<< exit(FatalIOError);
}
for (label i = 0; i < 3; i++)
{
if (order[i] == 'x')
{
decompOrder_[i] = 0;
}
else if (order[i] == 'y')
{
decompOrder_[i] = 1;
}
else if (order[i] == 'z')
{
decompOrder_[i] = 2;
}
else
{
FatalIOErrorIn
(
"hierarchGeomDecomp::hierarchGeomDecomp"
"(const dictionary& decompositionDict)",
decompositionDict_
) << "Illegal decomposition order " << order << endl
<< "It should only contain x, y or z" << exit(FatalError);
}
}
}
Foam::label Foam::hierarchGeomDecomp::findLower
(
const List<scalar>& l,
const scalar t,
const label initLow,
const label initHigh
)
{
if (initHigh <= initLow)
{
return initLow;
}
label low = initLow;
label high = initHigh;
while ((high - low) > 1)
{
label mid = (low + high)/2;
if (l[mid] < t)
{
low = mid;
}
else
{
high = mid;
}
}
// high and low can still differ by one. Choose best.
label tIndex = -1;
if (l[high-1] < t)
{
tIndex = high;
}
else
{
tIndex = low;
}
return tIndex;
}
// Create a mapping between the index and the weighted size.
// For convenience, sortedWeightedSize is one size bigger than current. This
// avoids extra tests.
void Foam::hierarchGeomDecomp::calculateSortedWeightedSizes
(
const labelList& current,
const labelList& indices,
const scalarField& weights,
const label globalCurrentSize,
scalarField& sortedWeightedSizes
)
{
// Evaluate cumulative weights.
sortedWeightedSizes[0] = 0;
forAll(current, i)
{
label pointI = current[indices[i]];
sortedWeightedSizes[i + 1] = sortedWeightedSizes[i] + weights[pointI];
}
// Non-dimensionalise and multiply by size.
scalar globalCurrentLength = returnReduce
(
sortedWeightedSizes[current.size()],
sumOp<scalar>()
);
// Normalise weights by global sum of weights and multiply through
// by global size.
sortedWeightedSizes *= (globalCurrentSize/globalCurrentLength);
}
// Find position in values so between minIndex and this position there
// are wantedSize elements.
void Foam::hierarchGeomDecomp::findBinary
(
const label sizeTol,
const List<scalar>& values,
const label minIndex, // index of previous value
const scalar minValue, // value at minIndex
const scalar maxValue, // global max of values
const scalar wantedSize, // wanted size
label& mid, // index where size of bin is
// wantedSize (to within sizeTol)
scalar& midValue // value at mid
)
{
label low = minIndex;
scalar lowValue = minValue;
scalar highValue = maxValue;
// (one beyond) index of highValue
label high = values.size();
// Safeguards to avoid infinite loop.
scalar midValuePrev = VGREAT;
while (true)
{
label size = returnReduce(mid-minIndex, sumOp<label>());
if (debug)
{
Pout<< " low:" << low << " lowValue:" << lowValue
<< " high:" << high << " highValue:" << highValue
<< " mid:" << mid << " midValue:" << midValue << endl
<< " globalSize:" << size << " wantedSize:" << wantedSize
<< " sizeTol:" << sizeTol << endl;
}
if (wantedSize < size - sizeTol)
{
high = mid;
highValue = midValue;
}
else if (wantedSize > size + sizeTol)
{
low = mid;
lowValue = midValue;
}
else
{
break;
}
// Update mid, midValue
midValue = 0.5*(lowValue+highValue);
mid = findLower(values, midValue, low, high);
// Safeguard if same as previous.
bool hasNotChanged = (mag(midValue-midValuePrev) < SMALL);
if (returnReduce(hasNotChanged, andOp<bool>()))
{
WarningIn("hierarchGeomDecomp::findBinary(..)")
<< "unable to find desired deomposition split, making do!"
<< endl;
break;
}
midValuePrev = midValue;
}
}
// Find position in values so between minIndex and this position there
// are wantedSize elements.
void Foam::hierarchGeomDecomp::findBinary
(
const label sizeTol,
const List<scalar>& sortedWeightedSizes,
const List<scalar>& values,
const label minIndex, // index of previous value
const scalar minValue, // value at minIndex
const scalar maxValue, // global max of values
const scalar wantedSize, // wanted size
label& mid, // index where size of bin is
// wantedSize (to within sizeTol)
scalar& midValue // value at mid
)
{
label low = minIndex;
scalar lowValue = minValue;
scalar highValue = maxValue;
// (one beyond) index of highValue
label high = values.size();
// Safeguards to avoid infinite loop.
scalar midValuePrev = VGREAT;
while (true)
{
label weightedSize = returnReduce
(
sortedWeightedSizes[mid] - sortedWeightedSizes[minIndex],
sumOp<label>()
);
if (debug)
{
Pout<< " low:" << low << " lowValue:" << lowValue
<< " high:" << high << " highValue:" << highValue
<< " mid:" << mid << " midValue:" << midValue << endl
<< " globalSize:" << weightedSize
<< " wantedSize:" << wantedSize
<< " sizeTol:" << sizeTol << endl;
}
if (wantedSize < weightedSize - sizeTol)
{
high = mid;
highValue = midValue;
}
else if (wantedSize > weightedSize + sizeTol)
{
low = mid;
lowValue = midValue;
}
else
{
break;
}
// Update mid, midValue
midValue = 0.5*(lowValue+highValue);
mid = findLower(values, midValue, low, high);
// Safeguard if same as previous.
bool hasNotChanged = (mag(midValue-midValuePrev) < SMALL);
if (returnReduce(hasNotChanged, andOp<bool>()))
{
WarningIn("hierarchGeomDecomp::findBinary(..)")
<< "unable to find desired deomposition split, making do!"
<< endl;
break;
}
midValuePrev = midValue;
}
}
// Sort points into bins according to one component. Recurses to next component.
void Foam::hierarchGeomDecomp::sortComponent
(
const label sizeTol,
const pointField& points,
const labelList& current, // slice of points to decompose
const direction componentIndex, // index in decompOrder_
const label mult, // multiplication factor for finalDecomp
labelList& finalDecomp
)
{
// Current component
label compI = decompOrder_[componentIndex];
if (debug)
{
Pout<< "sortComponent : Sorting slice of size " << current.size()
<< " in component " << compI << endl;
}
// Storage for sorted component compI
SortableList<scalar> sortedCoord(current.size());
forAll(current, i)
{
label pointI = current[i];
sortedCoord[i] = points[pointI][compI];
}
sortedCoord.sort();
label globalCurrentSize = returnReduce(current.size(), sumOp<label>());
scalar minCoord = returnReduce
(
(
sortedCoord.size()
? sortedCoord[0]
: GREAT
),
minOp<scalar>()
);
scalar maxCoord = returnReduce
(
(
sortedCoord.size()
? sortedCoord.last()
: -GREAT
),
maxOp<scalar>()
);
if (debug)
{
Pout<< "sortComponent : minCoord:" << minCoord
<< " maxCoord:" << maxCoord << endl;
}
// starting index (in sortedCoord) of bin (= local)
label leftIndex = 0;
// starting value of bin (= global since coordinate)
scalar leftCoord = minCoord;
// Sort bins of size n
for (label bin = 0; bin < n_[compI]; bin++)
{
// Now we need to determine the size of the bin (dx). This is
// determined by the 'pivot' values - everything to the left of this
// value goes in the current bin, everything to the right into the next
// bins.
// Local number of elements
label localSize = -1; // offset from leftOffset
// Value at right of bin (leftIndex+localSize-1)
scalar rightCoord = -GREAT;
if (bin == n_[compI]-1)
{
// Last bin. Copy all.
localSize = current.size()-leftIndex;
rightCoord = maxCoord; // note: not used anymore
}
else if (Pstream::nProcs() == 1)
{
// No need for binary searching of bin size
localSize = label(current.size()/n_[compI]);
rightCoord = sortedCoord[leftIndex+localSize];
}
else
{
// For the current bin (starting at leftCoord) we want a rightCoord
// such that the sum of all sizes are globalCurrentSize/n_[compI].
// We have to iterate to obtain this.
label rightIndex = current.size();
rightCoord = maxCoord;
// Calculate rightIndex/rightCoord to have wanted size
findBinary
(
sizeTol,
sortedCoord,
leftIndex,
leftCoord,
maxCoord,
globalCurrentSize/n_[compI], // wanted size
rightIndex,
rightCoord
);
localSize = rightIndex - leftIndex;
}
if (debug)
{
Pout<< "For component " << compI << ", bin " << bin
<< " copying" << endl
<< "from " << leftCoord << " at local index "
<< leftIndex << endl
<< "to " << rightCoord << " localSize:"
<< localSize << endl
<< endl;
}
// Copy localSize elements starting from leftIndex.
labelList slice(localSize);
forAll(slice, i)
{
label pointI = current[sortedCoord.indices()[leftIndex+i]];
// Mark point into correct bin
finalDecomp[pointI] += bin*mult;
// And collect for next sorting action
slice[i] = pointI;
}
// Sort slice in next component
if (componentIndex < 2)
{
string oldPrefix;
if (debug)
{
oldPrefix = Pout.prefix();
Pout.prefix() = " " + oldPrefix;
}
sortComponent
(
sizeTol,
points,
slice,
componentIndex+1,
mult*n_[compI], // Multiplier to apply to decomposition.
finalDecomp
);
if (debug)
{
Pout.prefix() = oldPrefix;
}
}
// Step to next bin.
leftIndex += localSize;
leftCoord = rightCoord;
}
}
// Sort points into bins according to one component. Recurses to next component.
void Foam::hierarchGeomDecomp::sortComponent
(
const label sizeTol,
const scalarField& weights,
const pointField& points,
const labelList& current, // slice of points to decompose
const direction componentIndex, // index in decompOrder_
const label mult, // multiplication factor for finalDecomp
labelList& finalDecomp
)
{
// Current component
label compI = decompOrder_[componentIndex];
if (debug)
{
Pout<< "sortComponent : Sorting slice of size " << current.size()
<< " in component " << compI << endl;
}
// Storage for sorted component compI
SortableList<scalar> sortedCoord(current.size());
forAll(current, i)
{
label pointI = current[i];
sortedCoord[i] = points[pointI][compI];
}
sortedCoord.sort();
label globalCurrentSize = returnReduce(current.size(), sumOp<label>());
// Now evaluate local cumulative weights, based on the sorting.
// Make one bigger than the nodes.
scalarField sortedWeightedSizes(current.size()+1, 0);
calculateSortedWeightedSizes
(
current,
sortedCoord.indices(),
weights,
globalCurrentSize,
sortedWeightedSizes
);
scalar minCoord = returnReduce
(
(
sortedCoord.size()
? sortedCoord[0]
: GREAT
),
minOp<scalar>()
);
scalar maxCoord = returnReduce
(
(
sortedCoord.size()
? sortedCoord.last()
: -GREAT
),
maxOp<scalar>()
);
if (debug)
{
Pout<< "sortComponent : minCoord:" << minCoord
<< " maxCoord:" << maxCoord << endl;
}
// starting index (in sortedCoord) of bin (= local)
label leftIndex = 0;
// starting value of bin (= global since coordinate)
scalar leftCoord = minCoord;
// Sort bins of size n
for (label bin = 0; bin < n_[compI]; bin++)
{
// Now we need to determine the size of the bin (dx). This is
// determined by the 'pivot' values - everything to the left of this
// value goes in the current bin, everything to the right into the next
// bins.
// Local number of elements
label localSize = -1; // offset from leftOffset
// Value at right of bin (leftIndex+localSize-1)
scalar rightCoord = -GREAT;
if (bin == n_[compI]-1)
{
// Last bin. Copy all.
localSize = current.size()-leftIndex;
rightCoord = maxCoord; // note: not used anymore
}
else
{
// For the current bin (starting at leftCoord) we want a rightCoord
// such that the sum of all weighted sizes are
// globalCurrentLength/n_[compI].
// We have to iterate to obtain this.
label rightIndex = current.size();
rightCoord = maxCoord;
// Calculate rightIndex/rightCoord to have wanted size
findBinary
(
sizeTol,
sortedWeightedSizes,
sortedCoord,
leftIndex,
leftCoord,
maxCoord,
globalCurrentSize/n_[compI], // wanted size
rightIndex,
rightCoord
);
localSize = rightIndex - leftIndex;
}
if (debug)
{
Pout<< "For component " << compI << ", bin " << bin
<< " copying" << endl
<< "from " << leftCoord << " at local index "
<< leftIndex << endl
<< "to " << rightCoord << " localSize:"
<< localSize << endl
<< endl;
}
// Copy localSize elements starting from leftIndex.
labelList slice(localSize);
forAll(slice, i)
{
label pointI = current[sortedCoord.indices()[leftIndex+i]];
// Mark point into correct bin
finalDecomp[pointI] += bin*mult;
// And collect for next sorting action
slice[i] = pointI;
}
// Sort slice in next component
if (componentIndex < 2)
{
string oldPrefix;
if (debug)
{
oldPrefix = Pout.prefix();
Pout.prefix() = " " + oldPrefix;
}
sortComponent
(
sizeTol,
weights,
points,
slice,
componentIndex+1,
mult*n_[compI], // Multiplier to apply to decomposition.
finalDecomp
);
if (debug)
{
Pout.prefix() = oldPrefix;
}
}
// Step to next bin.
leftIndex += localSize;
leftCoord = rightCoord;
}
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
Foam::hierarchGeomDecomp::hierarchGeomDecomp
(
const dictionary& decompositionDict
)
:
geomDecomp(decompositionDict, typeName),
decompOrder_()
{
setDecompOrder();
}
Foam::hierarchGeomDecomp::hierarchGeomDecomp
(
const dictionary& decompositionDict,
const polyMesh&
)
:
geomDecomp(decompositionDict, hierarchGeomDecomp::typeName),
decompOrder_()
{
setDecompOrder();
}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
Foam::labelList Foam::hierarchGeomDecomp::decompose
(
const pointField& points
)
{
// construct a list for the final result
labelList finalDecomp(points.size(), 0);
// Start off with every point sorted onto itself.
labelList slice(points.size());
forAll(slice, i)
{
slice[i] = i;
}
pointField rotatedPoints = rotDelta_ & points;
// Calculate tolerance of cell distribution. For large cases finding
// distibution to the cell exact would cause too many iterations so allow
// some slack.
label allSize = points.size();
reduce(allSize, sumOp<label>());
const label sizeTol = max(1, label(1E-3*allSize/nProcessors_));
// Sort recursive
sortComponent
(
sizeTol,
rotatedPoints,
slice,
0, // Sort first component in decompOrder.
1, // Offset for different x bins.
finalDecomp
);
return finalDecomp;
}
Foam::labelList Foam::hierarchGeomDecomp::decompose
(
const pointField& points,
const scalarField& weights
)
{
// construct a list for the final result
labelList finalDecomp(points.size(), 0);
// Start off with every point sorted onto itself.
labelList slice(points.size());
forAll(slice, i)
{
slice[i] = i;
}
pointField rotatedPoints = rotDelta_ & points;
// Calculate tolerance of cell distribution. For large cases finding
// distibution to the cell exact would cause too many iterations so allow
// some slack.
label allSize = points.size();
reduce(allSize, sumOp<label>());
const label sizeTol = max(1, label(1E-3*allSize/nProcessors_));
// Sort recursive
sortComponent
(
sizeTol,
weights,
rotatedPoints,
slice,
0, // Sort first component in decompOrder.
1, // Offset for different x bins.
finalDecomp
);
return finalDecomp;
}
// ************************************************************************* //

View File

@ -0,0 +1,245 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / 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
Class
Foam::hierarchGeomDecomp
Description
Does hierarchical decomposition of points. Works by first sorting the
points in x direction into equal sized bins, then in y direction and
finally in z direction.
Uses single array to hold decomposition which is indexed as if it is a
3 dimensional array:
finalDecomp[i,j,k] is indexed as
i*n[0]*n[1] + j*n[1] + k
E.g. if we're sorting 'xyz': the first sort (over the x-component)
determines in which x-domain the point goes. Then for each of the x-domains
the points are sorted in y direction and each individual x-domain gets
split into three y-domains. And similar for the z-direction.
Since the domains are of equal size the maximum difference in size is
n[0]*n[1] (or n[1]*n[2]?) (small anyway)
SourceFiles
hierarchGeomDecomp.C
\*---------------------------------------------------------------------------*/
#ifndef hierarchGeomDecomp_H
#define hierarchGeomDecomp_H
#include "geomDecomp.H"
#include "FixedList.H"
#include "direction.H"
namespace Foam
{
/*---------------------------------------------------------------------------*\
Class hierarchGeomDecomp Declaration
\*---------------------------------------------------------------------------*/
class hierarchGeomDecomp
:
public geomDecomp
{
// Private data
//- Decomposition order in terms of components.
FixedList<direction, 3> decompOrder_;
// Private Member Functions
//- Convert ordering string ("xyz") into list of components.
void setDecompOrder();
//- Evaluates the weighted sizes for each sorted point.
static void calculateSortedWeightedSizes
(
const labelList& current,
const labelList& indices,
const scalarField& weights,
const label globalCurrentSize,
scalarField& sortedWeightedSizes
);
//- Find index of t in list inbetween indices left and right
static label findLower
(
const List<scalar>&,
const scalar t,
const label left,
const label right
);
//- Find midValue (at local index mid) such that the number of
// elements between mid and leftIndex are (globally summed) the
// wantedSize. Binary search.
static void findBinary
(
const label sizeTol, // size difference considered acceptible
const List<scalar>&,
const label leftIndex, // index of previous value
const scalar leftValue, // value at leftIndex
const scalar maxValue, // global max of values
const scalar wantedSize, // wanted size
label& mid, // index where size of bin is wantedSize
scalar& midValue // value at mid
);
//- Find midValue (at local index mid) such that the number of
// elements between mid and leftIndex are (globally summed) the
// wantedSize. Binary search.
static void findBinary
(
const label sizeTol, // size difference considered acceptible
const List<scalar>& sortedWeightedSizes,
const List<scalar>&,
const label leftIndex, // index of previous value
const scalar leftValue, // value at leftIndex
const scalar maxValue, // global max of values
const scalar wantedSize, // wanted size
label& mid, // index where size of bin is wantedSize
scalar& midValue // value at mid
);
//- Recursively sort in x,y,z (or rather acc. to decompOrder_)
void sortComponent
(
const label sizeTol,
const pointField&,
const labelList& slice, // slice of points to decompose
const direction componentIndex, // index in decompOrder_
const label prevMult, // multiplication factor
labelList& finalDecomp // overall decomposition
);
//- Recursively sort in x,y,z (or rather acc. to decompOrder_)
//- using weighted points.
void sortComponent
(
const label sizeTol,
const scalarField& weights,
const pointField&,
const labelList& slice, // slice of points to decompose
const direction componentIndex, // index in decompOrder_
const label prevMult, // multiplication factor
labelList& finalDecomp // overall decomposition
);
//- Disallow default bitwise copy construct and assignment
void operator=(const hierarchGeomDecomp&);
hierarchGeomDecomp(const hierarchGeomDecomp&);
public:
//- Runtime type information
TypeName("hierarchical");
// Constructors
//- Construct given the decomposition dictionary
hierarchGeomDecomp(const dictionary& decompositionDict);
//- Construct given the decomposition dictionary and mesh
hierarchGeomDecomp
(
const dictionary& decompositionDict,
const polyMesh& mesh
);
// Destructor
virtual ~hierarchGeomDecomp()
{}
// Member Functions
//- hierarchgeom is aware of processor boundaries
virtual bool parallelAware() const
{
return true;
}
//- Return for every coordinate the wanted processor number. Use the
// mesh connectivity (if needed)
virtual labelList decompose
(
const pointField&,
const scalarField& weights
);
//- Without weights. Code for weighted decomposition is a bit complex
// so kept separate for now.
virtual labelList decompose(const pointField&);
//- Return for every coordinate the wanted processor number. Explicitly
// provided connectivity - does not use mesh_.
// The connectivity is equal to mesh.cellCells() except for
// - in parallel the cell numbers are global cell numbers (starting
// from 0 at processor0 and then incrementing all through the
// processors)
// - the connections are across coupled patches
virtual labelList decompose
(
const labelListList& globalCellCells,
const pointField& cc,
const scalarField& cWeights
)
{
return decompose(cc, cWeights);
}
virtual labelList decompose
(
const labelListList& globalCellCells,
const pointField& cc
)
{
return decompose(cc);
}
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //

View File

@ -0,0 +1,139 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / 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
Description
Decomposition given a cell-to-processor association in a file
\*---------------------------------------------------------------------------*/
#include "manualDecomp.H"
#include "addToRunTimeSelectionTable.H"
#include "IFstream.H"
#include "labelIOList.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
defineTypeNameAndDebug(manualDecomp, 0);
addToRunTimeSelectionTable
(
decompositionMethod,
manualDecomp,
dictionary
);
addToRunTimeSelectionTable
(
decompositionMethod,
manualDecomp,
dictionaryMesh
);
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
Foam::manualDecomp::manualDecomp(const dictionary& decompositionDict)
:
decompositionMethod(decompositionDict)
{
notImplemented("manualDecomp(const dictionary&)");
}
Foam::manualDecomp::manualDecomp
(
const dictionary& decompositionDict,
const polyMesh& mesh
)
:
decompositionMethod(decompositionDict),
meshPtr_(&mesh),
decompDataFile_
(
decompositionDict.subDict(word(decompositionDict.lookup("method"))
+ "Coeffs").lookup("dataFile")
)
{}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
Foam::labelList Foam::manualDecomp::decompose
(
const pointField& points,
const scalarField& pointWeights
)
{
const polyMesh& mesh = *meshPtr_;
labelIOList finalDecomp
(
IOobject
(
decompDataFile_,
mesh.facesInstance(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE,
false
)
);
// check if the final decomposition is OK
if (finalDecomp.size() != points.size())
{
FatalErrorIn
(
"manualDecomp::decompose(const pointField&, const scalarField&)"
) << "Size of decomposition list does not correspond "
<< "to the number of points. Size: "
<< finalDecomp.size() << " Number of points: "
<< points.size()
<< ".\n" << "Manual decomposition data read from file "
<< decompDataFile_ << "." << endl
<< exit(FatalError);
}
if (min(finalDecomp) < 0 || max(finalDecomp) > nProcessors_ - 1)
{
FatalErrorIn
(
"manualDecomp::decompose(const pointField&, const scalarField&)"
) << "According to the decomposition, cells assigned to "
<< "impossible processor numbers. Min processor = "
<< min(finalDecomp) << " Max processor = " << max(finalDecomp)
<< ".\n" << "Manual decomposition data read from file "
<< decompDataFile_ << "." << endl
<< exit(FatalError);
}
return finalDecomp;
}
// ************************************************************************* //

View File

@ -0,0 +1,135 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / 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
Class
Foam::manualDecomp
Description
Decomposition given a cell-to-processor association in a file
SourceFiles
manualDecomp.C
\*---------------------------------------------------------------------------*/
#ifndef manualDecomp_H
#define manualDecomp_H
#include "decompositionMethod.H"
namespace Foam
{
/*---------------------------------------------------------------------------*\
Class manualDecomp Declaration
\*---------------------------------------------------------------------------*/
class manualDecomp
:
public decompositionMethod
{
// Private data
const polyMesh* meshPtr_;
fileName decompDataFile_;
// Private Member Functions
//- Disallow default bitwise copy construct and assignment
void operator=(const manualDecomp&);
manualDecomp(const manualDecomp&);
public:
//- Runtime type information
TypeName("manual");
// Constructors
//- Construct given the decomposition dictionary
manualDecomp(const dictionary& decompositionDict);
//- Construct given the decomposition dictionary and mesh
manualDecomp
(
const dictionary& decompositionDict,
const polyMesh& mesh
);
// Destructor
virtual ~manualDecomp()
{}
// Member Functions
//- manual decompose does not care about proc boundaries - is all
// up to the user.
virtual bool parallelAware() const
{
return true;
}
//- Return for every coordinate the wanted processor number. Use the
// mesh connectivity (if needed)
virtual labelList decompose
(
const pointField& points,
const scalarField& pointWeights
);
//- Return for every coordinate the wanted processor number. Explicitly
// provided connectivity - does not use mesh_.
// The connectivity is equal to mesh.cellCells() except for
// - in parallel the cell numbers are global cell numbers (starting
// from 0 at processor0 and then incrementing all through the
// processors)
// - the connections are across coupled patches
virtual labelList decompose
(
const labelListList& globalCellCells,
const pointField& cc,
const scalarField& cWeights
)
{
return decompose(cc, cWeights);
}
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //

View File

@ -0,0 +1,319 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / 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 "simpleGeomDecomp.H"
#include "addToRunTimeSelectionTable.H"
#include "SortableList.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
defineTypeNameAndDebug(simpleGeomDecomp, 0);
addToRunTimeSelectionTable
(
decompositionMethod,
simpleGeomDecomp,
dictionary
);
addToRunTimeSelectionTable
(
decompositionMethod,
simpleGeomDecomp,
dictionaryMesh
);
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
// assignToProcessorGroup : given nCells cells and nProcGroup processor
// groups to share them, how do we share them out? Answer : each group
// gets nCells/nProcGroup cells, and the first few get one
// extra to make up the numbers. This should produce almost
// perfect load balancing
void Foam::simpleGeomDecomp::assignToProcessorGroup
(
labelList& processorGroup,
const label nProcGroup
)
{
label jump = processorGroup.size()/nProcGroup;
label jumpb = jump + 1;
label fstProcessorGroup = processorGroup.size() - jump*nProcGroup;
label ind = 0;
label j = 0;
// assign cells to the first few processor groups (those with
// one extra cell each
for (j=0; j<fstProcessorGroup; j++)
{
for (register label k=0; k<jumpb; k++)
{
processorGroup[ind++] = j;
}
}
// and now to the `normal' processor groups
for (; j<nProcGroup; j++)
{
for (register label k=0; k<jump; k++)
{
processorGroup[ind++] = j;
}
}
}
void Foam::simpleGeomDecomp::assignToProcessorGroup
(
labelList& processorGroup,
const label nProcGroup,
const labelList& indices,
const scalarField& weights,
const scalar summedWeights
)
{
// This routine gets the sorted points.
// Easiest to explain with an example.
// E.g. 400 points, sum of weights : 513.
// Now with number of divisions in this direction (nProcGroup) : 4
// gives the split at 513/4 = 128
// So summed weight from 0..128 goes into bin 0,
// ,, 128..256 goes into bin 1
// etc.
// Finally any remaining ones go into the last bin (3).
const scalar jump = summedWeights/nProcGroup;
const label nProcGroupM1 = nProcGroup - 1;
scalar sumWeights = 0;
label ind = 0;
label j = 0;
// assign cells to all except last group.
for (j=0; j<nProcGroupM1; j++)
{
const scalar limit = jump*scalar(j + 1);
while (sumWeights < limit)
{
sumWeights += weights[indices[ind]];
processorGroup[ind++] = j;
}
}
// Ensure last included.
while (ind < processorGroup.size())
{
processorGroup[ind++] = nProcGroupM1;
}
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
Foam::simpleGeomDecomp::simpleGeomDecomp(const dictionary& decompositionDict)
:
geomDecomp(decompositionDict, typeName)
{}
Foam::simpleGeomDecomp::simpleGeomDecomp
(
const dictionary& decompositionDict,
const polyMesh&
)
:
geomDecomp(decompositionDict, typeName)
{}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
Foam::labelList Foam::simpleGeomDecomp::decompose(const pointField& points)
{
// construct a list for the final result
labelList finalDecomp(points.size());
labelList processorGroups(points.size());
labelList pointIndices(points.size());
forAll(pointIndices, i)
{
pointIndices[i] = i;
}
pointField rotatedPoints = rotDelta_ & points;
// and one to take the processor group id's. For each direction.
// we assign the processors to groups of processors labelled
// 0..nX to give a banded structure on the mesh. Then we
// construct the actual processor number by treating this as
// the units part of the processor number.
sort
(
pointIndices,
UList<scalar>::less(rotatedPoints.component(vector::X))
);
assignToProcessorGroup(processorGroups, n_.x());
forAll(points, i)
{
finalDecomp[pointIndices[i]] = processorGroups[i];
}
// now do the same thing in the Y direction. These processor group
// numbers add multiples of nX to the proc. number (columns)
sort
(
pointIndices,
UList<scalar>::less(rotatedPoints.component(vector::Y))
);
assignToProcessorGroup(processorGroups, n_.y());
forAll(points, i)
{
finalDecomp[pointIndices[i]] += n_.x()*processorGroups[i];
}
// finally in the Z direction. Now we add multiples of nX*nY to give
// layers
sort
(
pointIndices,
UList<scalar>::less(rotatedPoints.component(vector::Z))
);
assignToProcessorGroup(processorGroups, n_.z());
forAll(points, i)
{
finalDecomp[pointIndices[i]] += n_.x()*n_.y()*processorGroups[i];
}
return finalDecomp;
}
Foam::labelList Foam::simpleGeomDecomp::decompose
(
const pointField& points,
const scalarField& weights
)
{
// construct a list for the final result
labelList finalDecomp(points.size());
labelList processorGroups(points.size());
labelList pointIndices(points.size());
forAll(pointIndices, i)
{
pointIndices[i] = i;
}
pointField rotatedPoints = rotDelta_ & points;
// and one to take the processor group id's. For each direction.
// we assign the processors to groups of processors labelled
// 0..nX to give a banded structure on the mesh. Then we
// construct the actual processor number by treating this as
// the units part of the processor number.
sort
(
pointIndices,
UList<scalar>::less(rotatedPoints.component(vector::X))
);
const scalar summedWeights = sum(weights);
assignToProcessorGroup
(
processorGroups,
n_.x(),
pointIndices,
weights,
summedWeights
);
forAll(points, i)
{
finalDecomp[pointIndices[i]] = processorGroups[i];
}
// now do the same thing in the Y direction. These processor group
// numbers add multiples of nX to the proc. number (columns)
sort
(
pointIndices,
UList<scalar>::less(rotatedPoints.component(vector::Y))
);
assignToProcessorGroup
(
processorGroups,
n_.y(),
pointIndices,
weights,
summedWeights
);
forAll(points, i)
{
finalDecomp[pointIndices[i]] += n_.x()*processorGroups[i];
}
// finally in the Z direction. Now we add multiples of nX*nY to give
// layers
sort
(
pointIndices,
UList<scalar>::less(rotatedPoints.component(vector::Z))
);
assignToProcessorGroup
(
processorGroups,
n_.z(),
pointIndices,
weights,
summedWeights
);
forAll(points, i)
{
finalDecomp[pointIndices[i]] += n_.x()*n_.y()*processorGroups[i];
}
return finalDecomp;
}
// ************************************************************************* //

View File

@ -0,0 +1,135 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / 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
Class
Foam::simpleGeomDecomp
Description
SourceFiles
simpleGeomDecomp.C
\*---------------------------------------------------------------------------*/
#ifndef simpleGeomDecomp_H
#define simpleGeomDecomp_H
#include "geomDecomp.H"
namespace Foam
{
/*---------------------------------------------------------------------------*\
Class simpleGeomDecomp Declaration
\*---------------------------------------------------------------------------*/
class simpleGeomDecomp
:
public geomDecomp
{
// Private Member Functions
void assignToProcessorGroup(labelList& processorGroup, const label);
void assignToProcessorGroup
(
labelList& processorGroup,
const label nProcGroup,
const labelList& indices,
const scalarField& weights,
const scalar summedWeights
);
//- Disallow default bitwise copy construct and assignment
void operator=(const simpleGeomDecomp&);
simpleGeomDecomp(const simpleGeomDecomp&);
public:
//- Runtime type information
TypeName("simple");
// Constructors
//- Construct given the decomposition dictionary
simpleGeomDecomp(const dictionary& decompositionDict);
//- Construct given the decomposition dictionary and mesh
simpleGeomDecomp
(
const dictionary& decompositionDict,
const polyMesh& mesh
);
// Destructor
virtual ~simpleGeomDecomp()
{}
// Member Functions
virtual bool parallelAware() const
{
// simpleDecomp does not attempt to do anything across proc
// boundaries
return false;
}
virtual labelList decompose
(
const pointField& points
);
virtual labelList decompose
(
const pointField& points,
const scalarField& pointWeights
);
//- Explicitly provided connectivity
virtual labelList decompose
(
const labelListList& globalCellCells,
const pointField& cc,
const scalarField& cWeights
)
{
return decompose(cc, cWeights);
}
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //