patchInjectionBase: ensure areaFraction is the same on all processors

Resolves bug-report http://www.openfoam.org/mantisbt/view.php?id=1556
This commit is contained in:
Henry
2015-03-05 17:20:34 +00:00
parent 6c0cc89cf5
commit 68c14d6bfe

View File

@ -2,7 +2,7 @@
========= | ========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | \\ / O peration |
\\ / A nd | Copyright (C) 2013-2014 OpenFOAM Foundation \\ / A nd | Copyright (C) 2013-2015 OpenFOAM Foundation
\\/ M anipulation | \\/ M anipulation |
------------------------------------------------------------------------------- -------------------------------------------------------------------------------
License License
@ -51,14 +51,14 @@ Foam::patchInjectionBase::patchInjectionBase
{ {
FatalErrorIn FatalErrorIn
( (
"Foam::patchInjectionBase::patchInjectionBase" "patchInjectionBase::patchInjectionBase"
"(" "("
"const polyMesh&, " "const polyMesh& mesh, "
"const word&" "const word& patchName"
")" ")"
) << "Requested patch " << patchName_ << " not found" << nl ) << "Requested patch " << patchName_ << " not found" << nl
<< "Available patches are: " << mesh.boundaryMesh().names() << "Available patches are: " << mesh.boundaryMesh().names() << nl
<< nl << exit(FatalError); << exit(FatalError);
} }
updateMesh(mesh); updateMesh(mesh);
@ -95,13 +95,13 @@ void Foam::patchInjectionBase::updateMesh(const polyMesh& mesh)
cellOwners_ = patch.faceCells(); cellOwners_ = patch.faceCells();
// triangulate the patch faces and create addressing // Triangulate the patch faces and create addressing
DynamicList<label> triToFace(2*patch.size()); DynamicList<label> triToFace(2*patch.size());
DynamicList<scalar> triMagSf(2*patch.size()); DynamicList<scalar> triMagSf(2*patch.size());
DynamicList<face> triFace(2*patch.size()); DynamicList<face> triFace(2*patch.size());
DynamicList<face> tris(5); DynamicList<face> tris(5);
// set zero value at the start of the tri area list // Set zero value at the start of the tri area list
triMagSf.append(0.0); triMagSf.append(0.0);
forAll(patch, faceI) forAll(patch, faceI)
@ -134,12 +134,12 @@ void Foam::patchInjectionBase::updateMesh(const polyMesh& mesh)
triMagSf[i] += triMagSf[i-1]; triMagSf[i] += triMagSf[i-1];
} }
// transfer to persistent storage // Transfer to persistent storage
triFace_.transfer(triFace); triFace_.transfer(triFace);
triToFace_.transfer(triToFace); triToFace_.transfer(triToFace);
triCumulativeMagSf_.transfer(triMagSf); triCumulativeMagSf_.transfer(triMagSf);
// convert sumTriMagSf_ into cumulative sum of areas per proc // Convert sumTriMagSf_ into cumulative sum of areas per proc
for (label i = 1; i < sumTriMagSf_.size(); i++) for (label i = 1; i < sumTriMagSf_.size(); i++)
{ {
sumTriMagSf_[i] += sumTriMagSf_[i-1]; sumTriMagSf_[i] += sumTriMagSf_[i-1];
@ -162,15 +162,22 @@ void Foam::patchInjectionBase::setPositionAndCell
label& tetPtI label& tetPtI
) )
{ {
scalar areaFraction = 0;
if (Pstream::master())
{
areaFraction = rnd.position<scalar>(0, patchArea_);
}
Pstream::scatter(areaFraction);
if (cellOwners_.size() > 0) if (cellOwners_.size() > 0)
{ {
// determine which processor to inject from // Determine which processor to inject from
scalar areaFraction = rnd.position<scalar>(0, patchArea_);
label procI = 0; label procI = 0;
forAllReverse(sumTriMagSf_, i) forAllReverse(sumTriMagSf_, i)
{ {
if (areaFraction > sumTriMagSf_[i]) if (areaFraction >= sumTriMagSf_[i])
{ {
procI = i; procI = i;
break; break;
@ -179,7 +186,7 @@ void Foam::patchInjectionBase::setPositionAndCell
if (Pstream::myProcNo() == procI) if (Pstream::myProcNo() == procI)
{ {
// find corresponding decomposed face triangle // Find corresponding decomposed face triangle
label triI = 0; label triI = 0;
scalar offset = sumTriMagSf_[procI]; scalar offset = sumTriMagSf_[procI];
forAllReverse(triCumulativeMagSf_, i) forAllReverse(triCumulativeMagSf_, i)
@ -191,18 +198,18 @@ void Foam::patchInjectionBase::setPositionAndCell
} }
} }
// set cellOwner // Set cellOwner
label faceI = triToFace_[triI]; label faceI = triToFace_[triI];
cellOwner = cellOwners_[faceI]; cellOwner = cellOwners_[faceI];
// find random point in triangle // Find random point in triangle
const polyPatch& patch = mesh.boundaryMesh()[patchId_]; const polyPatch& patch = mesh.boundaryMesh()[patchId_];
const pointField& points = patch.points(); const pointField& points = patch.points();
const face& tf = triFace_[triI]; const face& tf = triFace_[triI];
const triPointRef tri(points[tf[0]], points[tf[1]], points[tf[2]]); const triPointRef tri(points[tf[0]], points[tf[1]], points[tf[2]]);
const point pf(tri.randomPoint(rnd)); const point pf(tri.randomPoint(rnd));
// position perturbed away from face (into domain) // Position perturbed away from face (into domain)
const scalar a = rnd.position(scalar(0.1), scalar(0.5)); const scalar a = rnd.position(scalar(0.1), scalar(0.5));
const vector& pc = mesh.cellCentres()[cellOwner]; const vector& pc = mesh.cellCentres()[cellOwner];
const vector d = mag(pf - pc)*patchNormal_[faceI]; const vector d = mag(pf - pc)*patchNormal_[faceI];
@ -223,7 +230,7 @@ void Foam::patchInjectionBase::setPositionAndCell
tetFaceI = -1; tetFaceI = -1;
tetPtI = -1; tetPtI = -1;
// dummy position // Dummy position
position = pTraits<vector>::max; position = pTraits<vector>::max;
} }
} }
@ -233,7 +240,7 @@ void Foam::patchInjectionBase::setPositionAndCell
tetFaceI = -1; tetFaceI = -1;
tetPtI = -1; tetPtI = -1;
// dummy position // Dummy position
position = pTraits<vector>::max; position = pTraits<vector>::max;
} }
} }