exaFOAM/foam-extend4.1-coherent-io
Archived
1
1
Fork 0
Code Issues Pull requests Projects Releases Packages Wiki Activity

Fixed parallel problem: doubly cut edges for face decomposition

Browse source
This commit is contained in:
Hrvoje Jasak 2011-02-25 11:51:12 +00:00
parent 66ff6d52ca
commit 3882171409
16 changed files with 138 additions and 55 deletions
Show stats Download patch file Download diff file Expand all files Collapse all files

2
src/tetDecompositionFiniteElement/tetFemMatrix/tetFemMatrix.H
View file

@ -117,7 +117,7 @@ class tetFemMatrix
// parallel patch and protruding into the neighbouring domain.
// This is treated by doing a "local" multiplication on the
// neighbouring side and then adding the result.
//
// HJ, 13/Nov/2001
void addCouplingCoeffs();
void eliminateCouplingCoeffs();

7
src/tetDecompositionFiniteElement/tetFemMatrix/tetFemMatrixCheck.C
View file

@ -33,7 +33,7 @@ Description
namespace Foam
{
// * * * * * * * * * * * * * * Member functions * * * * * * * * * * * * * * * //
// * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
template<class Type>
void tetFemMatrix<Type>::check()
@ -128,7 +128,7 @@ void tetFemMatrix<Type>::check()
*this,
coupledBouCoeffs[interfaceI],
0,
true
Pstream::defaultCommsType
);
}
}
@ -144,7 +144,8 @@ void tetFemMatrix<Type>::check()
matrixSumOffDiag,
*this,
coupledBouCoeffs[interfaceI],
0
0,
Pstream::defaultCommsType
);
}
}

17
src/tetDecompositionFiniteElement/tetFemMatrix/tetFemMatrixSolve.C
View file

@ -47,6 +47,12 @@ lduSolverPerformance tetFemMatrix<Type>::solve
<< endl;
}
// Check the matrix
if (debug > 1)
{
this->check();
}
lduSolverPerformance solverPerfVec
(
"tetFemMatrix<Type>::solve",
@ -88,8 +94,15 @@ lduSolverPerformance tetFemMatrix<Type>::solve
addCouplingSource(sourceCmpt);
// Prepare for coupled interface update
FieldField<Field, scalar> coupledBouCoeffs(psi_.boundaryField().size());
FieldField<Field, scalar> coupledIntCoeffs(psi_.boundaryField().size());
FieldField<Field, scalar> coupledBouCoeffs
(
psi_.boundaryField().size()
);
FieldField<Field, scalar> coupledIntCoeffs
(
psi_.boundaryField().size()
);
forAll(psi_.boundaryField(), patchI)
{

2
src/tetDecompositionFiniteElement/tetFemMatrix/tetFemMatrixTools.C
View file

@ -60,7 +60,7 @@ void tetFemMatrix<Type>::storeBoundaryCoeffs() const
// requested. Once the list of requested fixes is complete, collapse the
// list. Note the use of list of pointers for the first part of operation
// and creation of collapsed list by copying the pointers into a pointer
// list.
// list. HJ, 28/Feb/2001
if (!boundaryConditionsSet_)
{

6
src/tetDecompositionFiniteElement/tetPolyMeshCellDecomp/calcTetPolyMeshCellDecompParPointData.C
View file

@ -56,8 +56,7 @@ void Foam::tetPolyMeshCellDecomp::calcParPointData() const
{
if
(
typeid(mesh_.boundaryMesh()[patchI])
== typeid(processorPolyPatch)
isA<processorPolyPatch>(mesh_.boundaryMesh()[patchI])
)
{
const labelList& p = mesh_.boundaryMesh()[patchI].meshPoints();
@ -93,8 +92,7 @@ void Foam::tetPolyMeshCellDecomp::calcParPointData() const
{
if
(
typeid(mesh_.boundaryMesh()[patchI])
== typeid(processorPolyPatch)
isA<processorPolyPatch>(mesh_.boundaryMesh()[patchI])
)
{
const labelList& p = mesh_.boundaryMesh()[patchI].meshPoints();

4
src/tetDecompositionFiniteElement/tetPolyMeshCellDecomp/tetPolyPatches/basic/coupledFace/coupledFaceTetPolyPatchCellDecomp.C
View file

@ -89,7 +89,7 @@ const Foam::pointField&
Foam::coupledFaceTetPolyPatchCellDecomp::localPoints() const
{
notImplemented("coupledFaceTetPolyPatchCellDecomp::localPoints() const");
return Field<point>::null();
return pointField::null();
}
@ -97,7 +97,7 @@ const Foam::vectorField&
Foam::coupledFaceTetPolyPatchCellDecomp::pointNormals() const
{
notImplemented("coupledFaceTetPolyPatchCellDecomp::pointNormals() const");
return Field<vector>::null();
return vectorField::null();
}

4
src/tetDecompositionFiniteElement/tetPolyMeshCellDecomp/tetPolyPatches/constraint/global/globalTetPolyPatchCellDecomp.C
View file

@ -111,14 +111,14 @@ void globalTetPolyPatchCellDecomp::clearCutEdgeAddressing() const
const pointField& globalTetPolyPatchCellDecomp::localPoints() const
{
notImplemented("globalTetPolyPatchCellDecomp::localPoints() const");
return Field<point>::null();
return pointField::null();
}
const vectorField& globalTetPolyPatchCellDecomp::pointNormals() const
{
notImplemented("globalTetPolyPatchCellDecomp::pointNormals() const");
return Field<point>::null();
return pointField::null();
}

5
src/tetDecompositionFiniteElement/tetPolyMeshCellDecomp/tetPolyPatches/constraint/processor/calcProcessorTetPolyPatchCellDecompAddr.C
View file

@ -388,8 +388,7 @@ void processorTetPolyPatchCellDecomp::calcOwnNeiDoubleMask() const
{
if
(
typeid(boundaryMesh()[patchI])
== typeid(processorTetPolyPatchCellDecomp)
isA<processorTetPolyPatchCellDecomp>(boundaryMesh()[patchI])
)
{
// Get the local points and mark up the list
@ -469,7 +468,7 @@ void processorTetPolyPatchCellDecomp::calcOwnNeiDoubleMask() const
// Doubly cut coefficients
// ~~~~~~~~~~~~~~~~~~~~~~~
// Not needed.
// Not needed. HJ, 18/Apr/2002
}

2
src/tetDecompositionFiniteElement/tetPolyMeshCellDecomp/tetPolyPatches/constraint/processor/processorTetPolyPatchCellDecomp.H
View file

@ -85,7 +85,7 @@ class processorTetPolyPatchCellDecomp
// so on; the breaks in the list are stored in the "start"
// mechanism, described in lduAddressing. This will be used to
// calculate the cross-processor contribution of the vector-matrix
// multiply.
// multiply. HJ, 13/Nov/2001
//- Cut edges owner list
mutable labelList* cutEdgeOwnerIndicesPtr_;

4
src/tetDecompositionFiniteElement/tetPolyMeshCellDecomp/tetPolyPatches/constraint/processor/processorTetPolyPatchCellDecompM.C
View file

@ -104,14 +104,14 @@ void processorTetPolyPatchCellDecomp::clearCutEdgeAddressing() const
const pointField& processorTetPolyPatchCellDecomp::localPoints() const
{
notImplemented("processorTetPolyPatchCellDecomp::localPoints() const");
return Field<point>::null();
return pointField::null();
}
const vectorField& processorTetPolyPatchCellDecomp::pointNormals() const
{
notImplemented("processorTetPolyPatchCellDecomp::pointNormals() const");
return Field<vector>::null();
return vectorField::null();
}

6
src/tetDecompositionFiniteElement/tetPolyMeshFaceDecomp/calcTetPolyMeshFaceDecompParPointData.C
View file

@ -56,8 +56,7 @@ void Foam::tetPolyMeshFaceDecomp::calcParPointData() const
{
if
(
typeid(mesh_.boundaryMesh()[patchI])
== typeid(processorPolyPatch)
isA<processorPolyPatch>(mesh_.boundaryMesh()[patchI])
)
{
const labelList& p = mesh_.boundaryMesh()[patchI].meshPoints();
@ -93,8 +92,7 @@ void Foam::tetPolyMeshFaceDecomp::calcParPointData() const
{
if
(
typeid(mesh_.boundaryMesh()[patchI])
== typeid(processorPolyPatch)
isA<processorPolyPatch>(mesh_.boundaryMesh()[patchI])
)
{
const labelList& p = mesh_.boundaryMesh()[patchI].meshPoints();

4
src/tetDecompositionFiniteElement/tetPolyMeshFaceDecomp/tetPolyPatches/basic/coupledFace/coupledFaceTetPolyPatchFaceDecomp.C
View file

@ -89,7 +89,7 @@ const Foam::pointField&
Foam::coupledFaceTetPolyPatchFaceDecomp::localPoints() const
{
notImplemented("coupledFaceTetPolyPatchFaceDecomp::localPoints() const");
return Field<point>::null();
return pointField::null();
}
@ -97,7 +97,7 @@ const Foam::vectorField&
Foam::coupledFaceTetPolyPatchFaceDecomp::pointNormals() const
{
notImplemented("coupledFaceTetPolyPatchFaceDecomp::pointNormals() const");
return Field<vector>::null();
return vectorField::null();
}

4
src/tetDecompositionFiniteElement/tetPolyMeshFaceDecomp/tetPolyPatches/constraint/global/globalTetPolyPatchFaceDecomp.C
View file

@ -111,14 +111,14 @@ void globalTetPolyPatchFaceDecomp::clearCutEdgeAddressing() const
const pointField& globalTetPolyPatchFaceDecomp::localPoints() const
{
notImplemented("globalTetPolyPatchFaceDecomp::localPoints() const");
return Field<point>::null();
return pointField::null();
}
const vectorField& globalTetPolyPatchFaceDecomp::pointNormals() const
{
notImplemented("globalTetPolyPatchFaceDecomp::pointNormals() const");
return Field<point>::null();
return pointField::null();
}

105
src/tetDecompositionFiniteElement/tetPolyMeshFaceDecomp/tetPolyPatches/constraint/processor/calcProcessorTetPolyPatchFaceDecompAddr.C
View file

@ -147,6 +147,9 @@ void processorTetPolyPatchFaceDecomp::calcCutEdgeAddressing() const
|| cutEdgeOwnerStartPtr_
|| cutEdgeNeighbourIndicesPtr_
|| cutEdgeNeighbourStartPtr_
|| doubleCutEdgeIndicesPtr_
|| doubleCutOwnerPtr_
|| doubleCutNeighbourPtr_
)
{
FatalErrorIn
@ -160,14 +163,20 @@ void processorTetPolyPatchFaceDecomp::calcCutEdgeAddressing() const
// Make a list over all edges in the mesh. Mark the ones that are local
// to the patch and then collect the rest
// For doubly cut edges, mark up the local points
const tetPolyMeshFaceDecomp& mesh = boundaryMesh().mesh();
boolList isLocal(mesh.nEdges(), false);
// get reference to local edge indices
// Get reference to local edge indices
const labelList& localEdges = localEdgeIndices();
// Get reference to mesh point addressing for the patch
const labelList& mp = meshPoints();
// Mark up the local edges
boolList isLocal(mesh.nEdges(), false);
forAll (localEdges, edgeI)
{
isLocal[localEdges[edgeI]] = true;
@ -190,13 +199,16 @@ void processorTetPolyPatchFaceDecomp::calcCutEdgeAddressing() const
isLocal[sharedCutEdges[sharedCutI]] = true;
}
labelList localPointLabel(mesh.nPoints(), -1);
forAll (mp, pointI)
{
localPointLabel[mp[pointI]] = pointI;
}
// Count the maximum number of edges coming from the patch
label maxEdgesOnPatch = 0;
// Get reference to mesh point addressing for the patch
const labelList& mp = meshPoints();
forAll (mp, pointI)
{
maxEdgesOnPatch += mesh.nEdgesForPoint(mp[pointI]);
@ -205,6 +217,21 @@ void processorTetPolyPatchFaceDecomp::calcCutEdgeAddressing() const
// Get reference to addressing
const lduAddressing& ldu = mesh.lduAddr();
// Allocate doubly cut arrays
doubleCutEdgeIndicesPtr_ = new labelList(maxEdgesOnPatch, -1);
labelList& doubleCutEdges = *doubleCutEdgeIndicesPtr_;
doubleCutOwnerPtr_ = new labelList(maxEdgesOnPatch, -1);
labelList& doubleCutOwn = *doubleCutOwnerPtr_;
doubleCutNeighbourPtr_ = new labelList(maxEdgesOnPatch, -1);
labelList& doubleCutNei = *doubleCutNeighbourPtr_;
label nDoubleCut = 0;
const labelList& globalOwner = ldu.lowerAddr();
const labelList& globalNeighbour = ldu.upperAddr();
// Owner side
//~~~~~~~~~~~
@ -238,9 +265,24 @@ void processorTetPolyPatchFaceDecomp::calcCutEdgeAddressing() const
{
if (!isLocal[edgeLabel])
{
own[nOwn] = edgeLabel;
isLocal[edgeLabel] = true;
nOwn++;
if (localPointLabel[globalNeighbour[edgeLabel]] == -1)
{
// Singly cut edge
own[nOwn] = edgeLabel;
nOwn++;
}
else
{
// Doubly cut edge
doubleCutEdges[nDoubleCut] = edgeLabel;
doubleCutOwn[nDoubleCut] = pointI;
doubleCutNei[nDoubleCut] =
localPointLabel[globalNeighbour[edgeLabel]];
nDoubleCut++;
}
}
}
}
@ -287,9 +329,24 @@ void processorTetPolyPatchFaceDecomp::calcCutEdgeAddressing() const
{
if (!isLocal[losort[edgeLabel]])
{
nei[nNei] = losort[edgeLabel];
isLocal[losort[edgeLabel]] = true;
nNei++;
if (localPointLabel[globalOwner[losort[edgeLabel]]] == -1)
{
// Singly cut edge
nei[nNei] = losort[edgeLabel];
nNei++;
}
else
{
// Doubly cut edge
doubleCutEdges[nDoubleCut] = losort[edgeLabel];
doubleCutOwn[nDoubleCut] =
localPointLabel[globalOwner[losort[edgeLabel]]];
doubleCutNei[nDoubleCut] = pointI;
nDoubleCut++;
}
}
}
}
@ -299,6 +356,11 @@ void processorTetPolyPatchFaceDecomp::calcCutEdgeAddressing() const
// set the last start label by hand
neiStart[meshPoints().size()] = nNei;
// Reset the size of double cut edge data
doubleCutEdges.setSize(nDoubleCut);
doubleCutOwn.setSize(nDoubleCut);
doubleCutNei.setSize(nDoubleCut);
}
@ -328,8 +390,7 @@ void processorTetPolyPatchFaceDecomp::calcOwnNeiDoubleMask() const
{
if
(
typeid(boundaryMesh()[patchI])
== typeid(processorTetPolyPatchFaceDecomp)
isA<processorTetPolyPatchFaceDecomp>(boundaryMesh()[patchI])
)
{
// Get the local points and mark up the list
@ -409,7 +470,7 @@ void processorTetPolyPatchFaceDecomp::calcOwnNeiDoubleMask() const
// Doubly cut coefficients
// ~~~~~~~~~~~~~~~~~~~~~~~
// Not needed.
// Not needed. HJ, 18/Apr/2002
}
@ -512,21 +573,35 @@ const labelList& processorTetPolyPatchFaceDecomp::cutEdgeNeighbourStart() const
}
// For face decompositon, doubly cut coefficients cannot occur.
const labelList& processorTetPolyPatchFaceDecomp::doubleCutEdgeIndices() const
{
if (!doubleCutEdgeIndicesPtr_)
{
calcCutEdgeAddressing();
}
return *doubleCutEdgeIndicesPtr_;
}
const labelList& processorTetPolyPatchFaceDecomp::doubleCutOwner() const
{
if (!doubleCutEdgeIndicesPtr_)
{
calcCutEdgeAddressing();
}
return *doubleCutOwnerPtr_;
}
const labelList& processorTetPolyPatchFaceDecomp::doubleCutNeighbour() const
{
if (!doubleCutEdgeIndicesPtr_)
{
calcCutEdgeAddressing();
}
return *doubleCutNeighbourPtr_;
}

2
src/tetDecompositionFiniteElement/tetPolyMeshFaceDecomp/tetPolyPatches/constraint/processor/processorTetPolyPatchFaceDecomp.H
View file

@ -89,7 +89,7 @@ class processorTetPolyPatchFaceDecomp
// so on; the breaks in the list are stored in the "start"
// mechanism, described in lduAddressing. This will be used to
// calculate the cross-processor contribution of the vector-matrix
// multiply.
// multiply. HJ, 13/Nov/2001
//- Cut edges owner list
mutable labelList* cutEdgeOwnerIndicesPtr_;

19
src/tetDecompositionFiniteElement/tetPolyMeshFaceDecomp/tetPolyPatches/constraint/processor/processorTetPolyPatchFaceDecompM.C
View file

@ -66,9 +66,9 @@ processorTetPolyPatchFaceDecomp::processorTetPolyPatchFaceDecomp
cutEdgeOwnerStartPtr_(NULL),
cutEdgeNeighbourIndicesPtr_(NULL),
cutEdgeNeighbourStartPtr_(NULL),
doubleCutEdgeIndicesPtr_(new labelList(0)),
doubleCutOwnerPtr_(new labelList(0)),
doubleCutNeighbourPtr_(new labelList(0)),
doubleCutEdgeIndicesPtr_(NULL),
doubleCutOwnerPtr_(NULL),
doubleCutNeighbourPtr_(NULL),
ownNeiDoubleMaskPtr_(NULL)
{}
@ -80,11 +80,6 @@ processorTetPolyPatchFaceDecomp::~processorTetPolyPatchFaceDecomp()
deleteDemandDrivenData(localEdgeIndicesPtr_);
clearCutEdgeAddressing();
// Clear storage for nonexistent things
deleteDemandDrivenData(doubleCutEdgeIndicesPtr_);
deleteDemandDrivenData(doubleCutOwnerPtr_);
deleteDemandDrivenData(doubleCutNeighbourPtr_);
}
@ -96,6 +91,10 @@ void processorTetPolyPatchFaceDecomp::clearCutEdgeAddressing() const
deleteDemandDrivenData(cutEdgeNeighbourIndicesPtr_);
deleteDemandDrivenData(cutEdgeNeighbourStartPtr_);
deleteDemandDrivenData(doubleCutEdgeIndicesPtr_);
deleteDemandDrivenData(doubleCutOwnerPtr_);
deleteDemandDrivenData(doubleCutNeighbourPtr_);
deleteDemandDrivenData(ownNeiDoubleMaskPtr_);
}
@ -105,14 +104,14 @@ void processorTetPolyPatchFaceDecomp::clearCutEdgeAddressing() const
const pointField& processorTetPolyPatchFaceDecomp::localPoints() const
{
notImplemented("processorTetPolyPatchFaceDecomp::localPoints() const");
return Field<point>::null();
return pointField::null();
}
const vectorField& processorTetPolyPatchFaceDecomp::pointNormals() const
{
notImplemented("processorTetPolyPatchFaceDecomp::pointNormals() const");
return Field<vector>::null();
return vectorField::null();
}