diff --git a/src/OpenFOAM/meshes/polyMesh/polyPatches/constraint/regionCouple/regionCouplePolyPatch.C b/src/OpenFOAM/meshes/polyMesh/polyPatches/constraint/regionCouple/regionCouplePolyPatch.C index 7ffcbdfca..42485878e 100644 --- a/src/OpenFOAM/meshes/polyMesh/polyPatches/constraint/regionCouple/regionCouplePolyPatch.C +++ b/src/OpenFOAM/meshes/polyMesh/polyPatches/constraint/regionCouple/regionCouplePolyPatch.C @@ -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(); diff --git a/src/finiteVolume/fields/fvPatchFields/constraint/regionCoupling/regionCouplingFvPatchField.C b/src/finiteVolume/fields/fvPatchFields/constraint/regionCoupling/regionCouplingFvPatchField.C index 80e0636da..e18f9ec1b 100644 --- a/src/finiteVolume/fields/fvPatchFields/constraint/regionCoupling/regionCouplingFvPatchField.C +++ b/src/finiteVolume/fields/fvPatchFields/constraint/regionCoupling/regionCouplingFvPatchField.C @@ -48,6 +48,15 @@ regionCouplingFvPatchField::weights const Field& 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 tweights(new scalarField(fOwn.size(), 0.5)); scalarField& weights = tweights(); @@ -125,7 +134,8 @@ regionCouplingFvPatchField::weights } else { - weights[faceI] = 0.5; + // Use 0.5 weights + //weights[faceI] = 0.5; } } @@ -135,6 +145,23 @@ regionCouplingFvPatchField::weights } +template +const Foam::Field& regionCouplingFvPatchField::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 @@ -288,8 +315,8 @@ tmp > regionCouplingFvPatchField::patchNeighbourField() const { Field sField = shadowPatchField().patchInternalField(); - tmp > tpnf - ( + tmp > tpnf + ( regionCouplePatch_.interpolate ( shadowPatchField().patchInternalField() @@ -341,20 +368,24 @@ void regionCouplingFvPatchField::initEvaluate const Pstream::commsTypes commsType ) { - if (!this->updated()) + if(debug) { - this->updateCoeffs(); + Info << "In regionCouplingFvPatchField::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& fOwn = this->originalPatchField(); const Field fNei = regionCouplePatch_.interpolate @@ -397,6 +428,14 @@ void regionCouplingFvPatchField::evaluate template void regionCouplingFvPatchField::updateCoeffs() { + if(debug) + { + Info << "In regionCouplingFvPatchField::updateCoeffs() on " + << this->dimensionedInternalField().name() + << " in " << this->patch().boundaryMesh().mesh().name() + << " " << this->updated() << endl; + } + if (this->updated()) { return; @@ -422,14 +461,6 @@ void regionCouplingFvPatchField::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::initInterfaceMatrixUpdate FatalErrorIn ( "regionCouplingFvPatchField::initInterfaceMatrixUpdate" - ) << "init matrix update calld in detached state" + ) << "init matrix update called in detached state" << abort(FatalError); } @@ -530,7 +561,7 @@ void regionCouplingFvPatchField::updateInterfaceMatrix FatalErrorIn ( "regionCouplingFvPatchField::updateInterfaceMatrix" - ) << "Matrix update calld in detached state" + ) << "Matrix update called in detached state" << abort(FatalError); } diff --git a/src/finiteVolume/fields/fvPatchFields/constraint/regionCoupling/regionCouplingFvPatchField.H b/src/finiteVolume/fields/fvPatchFields/constraint/regionCoupling/regionCouplingFvPatchField.H index 2d78f7b8d..6f1356c99 100644 --- a/src/finiteVolume/fields/fvPatchFields/constraint/regionCoupling/regionCouplingFvPatchField.H +++ b/src/finiteVolume/fields/fvPatchFields/constraint/regionCoupling/regionCouplingFvPatchField.H @@ -71,26 +71,16 @@ class regionCouplingFvPatchField //- Original patch field. Required for correct evaluation // in harmonic averaging - Field originalPatchField_; + mutable Field originalPatchField_; //- Current time index used to store originalPatchField_ - label curTimeIndex_; + mutable label curTimeIndex_; protected: //- Return original patch field - const Field& originalPatchField() const - { - if (curTimeIndex_ == this->db().time().timeIndex()) - { - return originalPatchField_; - } - else - { - return *this; - } - } + const Field& 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();