Skip to content
Snippets Groups Projects
Commit dbcb3812 authored by laurence's avatar laurence
Browse files

Merge branch 'master' of /home/dm4/OpenFOAM/OpenFOAM-dev

parents 8473c4b2 9426345a
Branches
Tags
No related merge requests found
......@@ -35,116 +35,6 @@ Description
#include "meshTools.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
void checkConnectedAgglomeration
(
const lduMesh& mesh,
const labelUList& restrict,
const label nCoarse
)
{
if (mesh.lduAddr().size() != restrict.size())
{
FatalErrorIn
(
"checkConnectedAgglomeration(const lduMesh&, const labelList&)"
) << "nCells:" << mesh.lduAddr().size()
<< " agglom:" << restrict.size()
<< abort(FatalError);
}
// Seed (master) for every region
labelList regionToMaster(nCoarse, -1);
labelList master(mesh.lduAddr().size(), -1);
forAll(restrict, cellI)
{
label region = restrict[cellI];
if (regionToMaster[region] == -1)
{
// Set cell to be master for region
//Pout<< "For region " << region
// << " allocating local master " << cellI
// << endl;
regionToMaster[region] = cellI;
master[cellI] = cellI;
}
}
// Now loop and transport master through region
const labelUList& lower = mesh.lduAddr().lowerAddr();
const labelUList& upper = mesh.lduAddr().upperAddr();
while (true)
{
label nChanged = 0;
forAll(lower, faceI)
{
label own = lower[faceI];
label nei = upper[faceI];
if (restrict[own] == restrict[nei])
{
// Region-internal face
if (master[own] != -1)
{
if (master[nei] == -1)
{
master[nei] = master[own];
nChanged++;
}
else if (master[nei] != master[own])
{
FatalErrorIn("checkConnectedAgglomeration(..)")
<< "problem" << abort(FatalError);
}
}
else if (master[nei] != -1)
{
master[own] = master[nei];
nChanged++;
}
}
}
reduce(nChanged, sumOp<label>());
if (nChanged == 0)
{
break;
}
}
// Check that master is set for all cells
boolList singleRegion(nCoarse, true);
label nSet = nCoarse;
forAll(master, cellI)
{
if (master[cellI] == -1)
{
label region = restrict[cellI];
if (singleRegion[region] == true)
{
singleRegion[region] = false;
nSet--;
}
}
}
label totalNCoarse = returnReduce(nCoarse, sumOp<label>());
label totalNVisited = returnReduce(nSet, sumOp<label>());
if (totalNVisited < totalNCoarse)
{
WarningIn("checkConnectedAgglomeration(..)")
<< "out of " << totalNCoarse
<< " agglomerated cells have " << totalNCoarse-totalNVisited
<< " cells that are not a single connected region" << endl;
}
}
// Main program:
int main(int argc, char *argv[])
......@@ -227,12 +117,28 @@ int main(int argc, char *argv[])
<< " agglomerated size : "
<< returnReduce(coarseSize, sumOp<label>()) << endl;
checkConnectedAgglomeration
labelList newAddr;
label newCoarseSize = 0;
bool ok = GAMGAgglomeration::checkRestriction
(
agglom.meshLevel(level),
newAddr,
newCoarseSize,
agglom.meshLevel(level).lduAddr(),
addr,
coarseSize
);
if (!ok)
{
WarningIn(args.executable())
<< "At level " << level
<< " there are " << coarseSize
<< " agglomerated cells but " << newCoarseSize
<< " disconnected regions" << endl
<< " This means that some agglomerations (coarse cells)"
<< " consist of multiple disconnected regions."
<< endl;
}
forAll(addr, fineI)
......
......@@ -535,4 +535,123 @@ const Foam::labelListListList& Foam::GAMGAgglomeration::boundaryFaceMap
}
bool Foam::GAMGAgglomeration::checkRestriction
(
labelList& newRestrict,
label& nNewCoarse,
const lduAddressing& fineAddressing,
const labelUList& restrict,
const label nCoarse
)
{
if (fineAddressing.size() != restrict.size())
{
FatalErrorIn
(
"checkRestriction(..)"
) << "nCells:" << fineAddressing.size()
<< " agglom:" << restrict.size()
<< abort(FatalError);
}
// Seed (master) for every region
labelList master(identity(fineAddressing.size()));
// Now loop and transport master through region
const labelUList& lower = fineAddressing.lowerAddr();
const labelUList& upper = fineAddressing.upperAddr();
while (true)
{
label nChanged = 0;
forAll(lower, faceI)
{
label own = lower[faceI];
label nei = upper[faceI];
if (restrict[own] == restrict[nei])
{
// coarse-mesh-internal face
if (master[own] < master[nei])
{
master[nei] = master[own];
nChanged++;
}
else if (master[own] > master[nei])
{
master[own] = master[nei];
nChanged++;
}
}
}
reduce(nChanged, sumOp<label>());
if (nChanged == 0)
{
break;
}
}
// Count number of regions/masters per coarse cell
labelListList coarseToMasters(nCoarse);
nNewCoarse = 0;
forAll(restrict, cellI)
{
labelList& masters = coarseToMasters[restrict[cellI]];
if (findIndex(masters, master[cellI]) == -1)
{
masters.append(master[cellI]);
nNewCoarse++;
}
}
if (nNewCoarse > nCoarse)
{
//WarningIn("GAMGAgglomeration::checkRestriction(..)")
// << "Have " << nCoarse
// << " agglomerated cells but " << nNewCoarse
// << " disconnected regions" << endl;
// Keep coarseToMasters[0] the original coarse, allocate new ones
// for the others
labelListList coarseToNewCoarse(coarseToMasters.size());
nNewCoarse = nCoarse;
forAll(coarseToMasters, coarseI)
{
const labelList& masters = coarseToMasters[coarseI];
labelList& newCoarse = coarseToNewCoarse[coarseI];
newCoarse.setSize(masters.size());
newCoarse[0] = coarseI;
for (label i = 1; i < newCoarse.size(); i++)
{
newCoarse[i] = nNewCoarse++;
}
}
newRestrict.setSize(fineAddressing.size());
forAll(restrict, cellI)
{
label coarseI = restrict[cellI];
label index = findIndex(coarseToMasters[coarseI], master[cellI]);
newRestrict[cellI] = coarseToNewCoarse[coarseI][index];
}
return false;
}
else
{
return true;
}
}
// ************************************************************************* //
......@@ -474,6 +474,16 @@ public:
const labelListListList& boundaryFaceMap(const label fineLeveli)
const;
//- Given restriction determines if coarse cells are connected.
// Return ok is so, otherwise creates new restriction that is
static bool checkRestriction
(
labelList& newRestrict,
label& nNewCoarse,
const lduAddressing& fineAddressing,
const labelUList& restrict,
const label nCoarse
);
};
......
......@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2011 OpenFOAM Foundation
\\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
......@@ -156,6 +156,26 @@ Foam::tmp<Foam::labelField> Foam::MGridGenGAMGAgglomeration::agglomerate
finalAgglom.begin()
);
{
label nNewCoarseCells = 0;
labelList newRestrictAddr;
bool ok = checkRestriction
(
newRestrictAddr,
nNewCoarseCells
,
fineAddressing,
finalAgglom,
nCoarseCells
);
if (!ok)
{
nCoarseCells = nNewCoarseCells;
finalAgglom.transfer(newRestrictAddr);
}
}
return tmp<labelField>(new labelField(finalAgglom));
}
......
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment