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

BUG: regionCouple in parallel & originalField storage

Browse source
This commit is contained in:
Henrik Rusche 2012-06-08 16:08:42 +02:00
parent 1ea9d81566
commit f2ef2a57db
3 changed files with 91 additions and 51 deletions
Show stats Download patch file Download diff file Expand all files Collapse all files

44
src/OpenFOAM/meshes/polyMesh/polyPatches/constraint/regionCouple/regionCouplePolyPatch.C
View file

@ -800,7 +800,10 @@ void Foam::regionCouplePolyPatch::attach() const
// Patch-to-patch interpolation does not need to be cleared,
// only face/cell centres and interpolation factors
// HJ, 6/Jun/2011
clearGeom();
// Breaks simple parallel CHT cases
// HR, 8/Jun/2012
//clearGeom();
}
}
@ -815,7 +818,10 @@ void Foam::regionCouplePolyPatch::detach() const
// Patch-to-patch interpolation does not need to be cleared,
// only face/cell centres and interpolation factors
// HJ, 6/Jun/2011
clearGeom();
// Breaks simple parallel CHT cases
// HR, 8/Jun/2012
//clearGeom();
}
}
@ -883,15 +889,17 @@ Foam::regionCouplePolyPatch::patchToPatch() const
const Foam::vectorField&
Foam::regionCouplePolyPatch::reconFaceCellCentres() const
{
if (!attached_)
{
FatalErrorIn
(
"const vectorField& "
"regionCouplePolyPatch::reconFaceCellCentres() const"
) << "Requesting reconFaceCellCentres in detached state"
<< abort(FatalError);
}
// Breaks simple parallel CHT cases
// HR, 8/Jun/2012
// if (!attached_)
// {
// FatalErrorIn
// (
// "const vectorField& "
// "regionCouplePolyPatch::reconFaceCellCentres() const"
// ) << "Requesting reconFaceCellCentres in detached state"
// << abort(FatalError);
// }
if (!reconFaceCellCentresPtr_)
{
@ -919,6 +927,10 @@ void Foam::regionCouplePolyPatch::initAddressing()
{
// Calculate send addressing
sendAddr();
// Deferred execution on startup
// HR, 8/Jun/2012
shadow().sendAddr();
}
}
@ -943,6 +955,12 @@ void Foam::regionCouplePolyPatch::initGeometry()
{
reconFaceCellCentres();
}
else
{
// Deferred execution on startup
// HR, 8/Jun/2012
shadow().reconFaceCellCentres();
}
}
polyPatch::initGeometry();
@ -982,10 +1000,14 @@ void Foam::regionCouplePolyPatch::initMovePoints(const pointField& p)
if (Pstream::parRun() && !localParallel())
{
// Calculate send addressing
sendAddr();
}
}
// Reconsider: Two communication in one function call may not work in parallel
// HR, 8/Jun/2012
if (active() && master())
{
reconFaceCellCentres();

77
src/finiteVolume/fields/fvPatchFields/constraint/regionCoupling/regionCouplingFvPatchField.C
View file

@ -48,6 +48,15 @@ regionCouplingFvPatchField<Type>::weights
const Field<Type>& fNei
) const
{
// Implement weights-based stabilised harmonic interpolation using
// magnitude of type
// Algorithm:
// 1) calculate magnitude of internal field and neighbour field
// 2) calculate harmonic mean magnitude
// 3) express harmonic mean magnitude as: mean = w*mOwn + (1 - w)*mNei
// 4) Based on above, calculate w = (mean - mNei)/(mOwn - mNei)
// 5) Use weights to interpolate values
tmp<scalarField> tweights(new scalarField(fOwn.size(), 0.5));
scalarField& weights = tweights();
@ -125,7 +134,8 @@ regionCouplingFvPatchField<Type>::weights
}
else
{
weights[faceI] = 0.5;
// Use 0.5 weights
//weights[faceI] = 0.5;
}
}
@ -135,6 +145,23 @@ regionCouplingFvPatchField<Type>::weights
}
template<class Type>
const Foam::Field<Type>& regionCouplingFvPatchField<Type>::originalPatchField() const
{
if (curTimeIndex_ != this->db().time().timeIndex())
{
// Store original field for symmetric evaluation
// Henrik Rusche, Aug/2011
Info << "store original field" << endl;
originalPatchField_ = *this;
curTimeIndex_ = this->db().time().timeIndex();
}
return originalPatchField_;
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
template<class Type>
@ -288,8 +315,8 @@ tmp<Field<Type> > regionCouplingFvPatchField<Type>::patchNeighbourField() const
{
Field<Type> sField = shadowPatchField().patchInternalField();
tmp<Field<Type> > tpnf
(
tmp<Field<Type> > tpnf
(
regionCouplePatch_.interpolate
(
shadowPatchField().patchInternalField()
@ -341,20 +368,24 @@ void regionCouplingFvPatchField<Type>::initEvaluate
const Pstream::commsTypes commsType
)
{
if (!this->updated())
if(debug)
{
this->updateCoeffs();
Info << "In regionCouplingFvPatchField<Type>::initEvaluate() on "
<< this->dimensionedInternalField().name()
<< " in " << this->patch().boundaryMesh().mesh().name()
<< " " << this->updated() << endl;
}
// Interpolation must happen at init
// Implement weights-based stabilised harmonic interpolation using
// magnitude of type
// Algorithm:
// 1) calculate magnitude of internal field and neighbour field
// 2) calculate harmonic mean magnitude
// 3) express harmonic mean magnitude as: mean = w*mOwn + (1 - w)*mNei
// 4) Based on above, calculate w = (mean - mNei)/(mOwn - mNei)
// 5) Use weights to interpolate values
// Note: If used with interpolation - either on explicitly or called by the
// laplacian operator, the values set here are overridden by the interpolation
// scheme. In order to get the same diffusivities on both sides an identical
// interpolation scheme must be used.
// Note^2: Even if harmonic used, the interpolation is still wrong for most
// CHT cases since (cell values vs. face values)
// Note^3: None of this is intuitiv - fix requires low-level changes!
// HR, 8/Jun/2012
const Field<Type>& fOwn = this->originalPatchField();
const Field<Type> fNei = regionCouplePatch_.interpolate
@ -397,6 +428,14 @@ void regionCouplingFvPatchField<Type>::evaluate
template<class Type>
void regionCouplingFvPatchField<Type>::updateCoeffs()
{
if(debug)
{
Info << "In regionCouplingFvPatchField<Type>::updateCoeffs() on "
<< this->dimensionedInternalField().name()
<< " in " << this->patch().boundaryMesh().mesh().name()
<< " " << this->updated() << endl;
}
if (this->updated())
{
return;
@ -422,14 +461,6 @@ void regionCouplingFvPatchField<Type>::updateCoeffs()
regionCouplePatch_.bridge(bridgeField, *this);
}
// Store original field for symmetric evaluation
// Henrik Rusche, Aug/2011
if (curTimeIndex_ != this->db().time().timeIndex())
{
originalPatchField_ = *this;
curTimeIndex_ = this->db().time().timeIndex();
}
}
@ -479,7 +510,7 @@ void regionCouplingFvPatchField<Type>::initInterfaceMatrixUpdate
FatalErrorIn
(
"regionCouplingFvPatchField<Type>::initInterfaceMatrixUpdate"
) << "init matrix update calld in detached state"
) << "init matrix update called in detached state"
<< abort(FatalError);
}
@ -530,7 +561,7 @@ void regionCouplingFvPatchField<Type>::updateInterfaceMatrix
FatalErrorIn
(
"regionCouplingFvPatchField<Type>::updateInterfaceMatrix"
) << "Matrix update calld in detached state"
) << "Matrix update called in detached state"
<< abort(FatalError);
}

21
src/finiteVolume/fields/fvPatchFields/constraint/regionCoupling/regionCouplingFvPatchField.H
View file

@ -71,26 +71,16 @@ class regionCouplingFvPatchField
//- Original patch field. Required for correct evaluation
// in harmonic averaging
Field<Type> originalPatchField_;
mutable Field<Type> originalPatchField_;
//- Current time index used to store originalPatchField_
label curTimeIndex_;
mutable label curTimeIndex_;
protected:
//- Return original patch field
const Field<Type>& originalPatchField() const
{
if (curTimeIndex_ == this->db().time().timeIndex())
{
return originalPatchField_;
}
else
{
return *this;
}
}
const Field<Type>& originalPatchField() const;
//- Return contents of a matrix update buffer
const scalarField& matrixUpdateBuffer() const
@ -213,10 +203,7 @@ public:
virtual void initEvaluate(const Pstream::commsTypes commsType);
//- Evaluate the patch field
virtual void evaluate
(
const Pstream::commsTypes commsType
);
virtual void evaluate(const Pstream::commsTypes commsType);
//- Update the coefficients associated with the patch field
virtual void updateCoeffs();