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

Improvements in conjugate heat solver: non-matching interface, parallelisation

Browse source
This commit is contained in:
Hrvoje Jasak 2011-09-29 02:12:36 +01:00
parent 2cabd52ff0
commit da49cac84d
25 changed files with 2765 additions and 755 deletions
Show stats Download patch file Download diff file Expand all files Collapse all files

9
src/OpenFOAM/meshes/polyMesh/polyPatches/constraint/ggi/ggiPolyPatch.C
View file

@ -72,12 +72,6 @@ void Foam::ggiPolyPatch::calcZoneAddressing() const
<< abort(FatalError);
}
if (debug)
{
Pout<< "ggiPolyPatch::calcZoneAddressing() const for patch "
<< index() << endl;
}
// Calculate patch-to-zone addressing
zoneAddressingPtr_ = new labelList(size());
labelList& zAddr = *zoneAddressingPtr_;
@ -195,7 +189,7 @@ void Foam::ggiPolyPatch::calcPatchToPatch() const
new ggiZoneInterpolation
(
zone()(), // This zone reference
shadow().zone()(), // This shadow zone reference
shadow().zone()(), // Shadow zone reference
forwardT(),
reverseT(),
shadow().separation(), // Slave-to-master separation. Bug fix
@ -640,6 +634,7 @@ Foam::label Foam::ggiPolyPatch::zoneIndex() const
{
FatalErrorIn("label ggiPolyPatch::zoneIndex() const")
<< "Face zone name " << zoneName_
<< " for GGI patch " << name()
<< " not found. Please check your GGI interface definition."
<< abort(FatalError);
}

188
src/OpenFOAM/meshes/polyMesh/polyPatches/constraint/ggi/ggiPolyPatch.H
View file

@ -279,116 +279,128 @@ public:
// Member functions
//- Return shadow patch name
const word& shadowName() const
{
return shadowName_;
}
// Basic info
//- Return name of interpolation face zone
const word& zoneName() const
{
return zoneName_;
}
//- Return shadow patch name
const word& shadowName() const
{
return shadowName_;
}
//- Return shadow patch index
label shadowIndex() const;
//- Return name of interpolation face zone
const word& zoneName() const
{
return zoneName_;
}
//- Return zone patch index
label zoneIndex() const;
//- Return shadow patch index
label shadowIndex() const;
//- Return shadow patch
const ggiPolyPatch& shadow() const;
//- Return zone patch index
label zoneIndex() const;
//- Return interpolation face zone
const faceZone& zone() const;
//- Return shadow patch
const ggiPolyPatch& shadow() const;
//- Is this the master side?
bool master() const
{
return index() < shadowIndex();
}
//- Return interpolation face zone
const faceZone& zone() const;
//- Is this the slave side?
bool slave() const
{
return !master();
}
// Interpolation data
//- Use bridging to fix overlap error in interpolation
bool bridgeOverlap() const
{
return bridgeOverlap_;
}
//- Is this the master side?
bool master() const
{
return index() < shadowIndex();
}
//- Return patch-to-zone addressing
const labelList& zoneAddressing() const;
//- Is this the slave side?
bool slave() const
{
return !master();
}
//- Return remote patch-to-zone addressing
const labelList& remoteZoneAddressing() const;
//- Use bridging to fix overlap error in interpolation
bool bridgeOverlap() const
{
return bridgeOverlap_;
}
//- Is the patch localised on a single processor
bool localParallel() const;
//- Return patch-to-zone addressing
const labelList& zoneAddressing() const;
//- Return reference to patch-to-patch interpolation
// Used only for addressing
const ggiZoneInterpolation& patchToPatch() const;
//- Return remote patch-to-zone addressing
const labelList& remoteZoneAddressing() const;
//- Expand face field to full zone size, including reduction
// New. HJ, 12/Jun/2011
template<class Type>
tmp<Field<Type> > fastExpand(const Field<Type>& pf) const;
//- Is the patch localised on a single processor
bool localParallel() const;
//- Expand face field to full zone size, including reduction
// Obsolete. HJ, 12/Jun/2011
template<class Type>
tmp<Field<Type> > expand(const Field<Type>& pf) const;
//- Filter zone field to patch size
// Obsolete. HJ, 12/Jun/2011
template<class Type>
tmp<Field<Type> > filter(const Field<Type>& zf) const;
//- Return reference to patch-to-patch interpolation
// Used only for addressing
const ggiZoneInterpolation& patchToPatch() const;
//- Interpolate face field: given field on a the shadow side,
// create an interpolated field on this side
template<class Type>
tmp<Field<Type> > interpolate(const Field<Type>& pf) const;
// Interpolation functions
template<class Type>
tmp<Field<Type> > interpolate(const tmp<Field<Type> >& tpf) const;
//- Expand face field to full zone size, including reduction
// New. HJ, 12/Jun/2011
template<class Type>
tmp<Field<Type> > fastExpand(const Field<Type>& pf) const;
//- Bridge interpolated face field for uncovered faces
template<class Type>
void bridge
(
const Field<Type>& bridgeField,
Field<Type>& ff
) const;
//- Expand face field to full zone size, including reduction
// Obsolete. HJ, 12/Jun/2011
template<class Type>
tmp<Field<Type> > expand(const Field<Type>& pf) const;
//- Return reconstructed cell centres
const vectorField& reconFaceCellCentres() const;
//- Filter zone field to patch size
// Obsolete. HJ, 12/Jun/2011
template<class Type>
tmp<Field<Type> > filter(const Field<Type>& zf) const;
//- Initialize ordering for primitivePatch. Does not
// refer to *this (except for name() and type() etc.)
virtual void initOrder(const primitivePatch&) const;
//- Interpolate face field: given field on a the shadow side,
// create an interpolated field on this side
template<class Type>
tmp<Field<Type> > interpolate(const Field<Type>& pf) const;
//- Return new ordering for primitivePatch.
// Ordering is -faceMap: for every face
// index of the new face -rotation: for every new face the clockwise
// shift of the original face. Return false if nothing changes
// (faceMap is identity, rotation is 0), true otherwise.
virtual bool order
(
const primitivePatch&,
labelList& faceMap,
labelList& rotation
) const;
template<class Type>
tmp<Field<Type> > interpolate(const tmp<Field<Type> >& tpf) const;
//- Synchronise communications of ordering for primitivePatch
// Used in cases when no topological change happens locally,
// but is happening on other processors
virtual void syncOrder() const;
//- Bridge interpolated face field for uncovered faces
template<class Type>
void bridge
(
const Field<Type>& bridgeField,
Field<Type>& ff
) const;
// Geometric data
//- Return reconstructed cell centres
const vectorField& reconFaceCellCentres() const;
// Patch ordering
//- Initialize ordering for primitivePatch. Does not
// refer to *this (except for name() and type() etc.)
virtual void initOrder(const primitivePatch&) const;
//- Return new ordering for primitivePatch.
// Ordering is -faceMap: for every face
// index of the new face -rotation: for every new face the
// clockwise shift of the original face. Return false if nothing
// changes (faceMap is identity, rotation is 0), true otherwise.
virtual bool order
(
const primitivePatch&,
labelList& faceMap,
labelList& rotation
) const;
//- Synchronise communications of ordering for primitivePatch
// Used in cases when no topological change happens locally,
// but is happening on other processors
virtual void syncOrder() const;
//- Write

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

File diff suppressed because it is too large Load diff

307
src/OpenFOAM/meshes/polyMesh/polyPatches/constraint/regionCouple/regionCouplePolyPatch.H
View file

@ -40,7 +40,8 @@ SourceFiles
#define regionCouplePolyPatch_H
#include "coupledPolyPatch.H"
#include "patchToPatchInterpolation.H"
#include "ggiInterpolation.H"
#include "faceZone.H"
#include "Switch.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@ -66,33 +67,108 @@ class regionCouplePolyPatch
//- Shadow patch name
const word shadowPatchName_;
//- Interpolation zone name
const word zoneName_;
//- Are the regions attached
mutable Switch attached_;
//- Master side. Defined by user, checked for consistency
Switch master_;
//- Are the region attached walls
Switch isWall_;
//- Use bridging to fix overlap error in interpolation
Switch bridgeOverlap_;
//- Shadow patch index. Delayed evaluation for construction
mutable label shadowIndex_;
//- Interpolation zone index. Delayed evaluation for construction
mutable label zoneIndex_;
//- Patch-to-patch interpolation
mutable patchToPatchInterpolation* patchToPatchPtr_;
mutable ggiZoneInterpolation* patchToPatchPtr_;
//- Patch-to-zone addressing
mutable labelList* zoneAddressingPtr_;
//- Remote zone addressing: data needed for interpolation
mutable labelList* remoteZoneAddressingPtr_;
//- Reconstructed patch neighbour cell centres
mutable vectorField* reconFaceCellCentresPtr_;
// Parallel communication optimisation, stored on master processor
//- Is the patch localised on a single processor
// (single processor in a parallel run)?
// Used for parallel optimisation
mutable bool* localParallelPtr_;
//- List of zone faces indices received from each processor
mutable labelListList* receiveAddrPtr_;
//- List of zone faces indices to send to each processor
mutable labelListList* sendAddrPtr_;
// Private member functions
//- Calculate patch-to-zone addressing
virtual void calcZoneAddressing() const;
//- Calculate remote patch-to-zone addressing
virtual void calcRemoteZoneAddressing() const;
//- Calculate interpolation
void calcInterpolation() const;
virtual void calcPatchToPatch() const;
// Geometry
//- Calculate reconstructed cell centres
void calcReconFaceCellCentres() const;
//- Force calculation of transformation tensors
virtual void calcTransforms();
// Parallel communication optimisation, stored on master processor
//- Calculate local parallel switch
void calcLocalParallel() const;
//- Calculate send and receive addressing
void calcSendReceive() const;
//- Return receive addressing
const labelListList& receiveAddr() const;
//- Return send addressing
const labelListList& sendAddr() const;
// Memory management
//- Clear geometry
void clearGeom() const;
//- Clear out
void clearOut();
protected:
// Protected Member functions
//- Is the region couple active? (zone and shadow present)
bool active() const;
//- Initialise the calculation of the patch addressing
virtual void initAddressing();
@ -121,18 +197,6 @@ protected:
//- Return non-constant access to shadow patch
regionCouplePolyPatch& shadow();
//- Return reference to patch-to-patch interpolation
const patchToPatchInterpolation& patchToPatch() const;
//- Calculate reconstructed cell centres
void calcReconFaceCellCentres() const;
//- Clear geometry
void clearGeom() const;
//- Clear out
void clearOut() const;
public:
@ -142,6 +206,16 @@ public:
// Constructors
//- Construct from components
regionCouplePolyPatch
(
const word& name,
const label size,
const label start,
const label index,
const polyBoundaryMesh& bm
);
//- Construct from components
regionCouplePolyPatch
(
@ -152,8 +226,11 @@ public:
const polyBoundaryMesh& bm,
const word& shadowRegionName,
const word& shadowPatchName,
const word& zoneName,
const bool attached,
const bool attachedWall
const bool master,
const bool isWall,
const bool bridgeOverlap
);
//- Construct from dictionary
@ -220,93 +297,143 @@ public:
// Member functions
//- Return true if coupled
virtual bool coupled() const;
// Basic info
//- Return shadow region name
const word& shadowRegionName() const
{
return shadowRegionName_;
}
//- Return true if coupled
virtual bool coupled() const;
//- Return shadow patch name
const word& shadowPatchName() const
{
return shadowPatchName_;
}
//- Return shadow region name
const word& shadowRegionName() const
{
return shadowRegionName_;
}
//- Return shadow region
const polyMesh& shadowRegion() const;
//- Return shadow patch name
const word& shadowPatchName() const
{
return shadowPatchName_;
}
//- Return shadow patch index
int shadowIndex() const;
//- Return name of interpolation face zone
const word& zoneName() const
{
return zoneName_;
}
//- Return true if regions are attached
bool attached() const
{
return attached_;
}
//- Return shadow region
const polyMesh& shadowRegion() const;
// Return true if patch is considered a wall
bool isWall() const
{
return isWall_;
}
//- Return shadow patch index
label shadowIndex() const;
//- Attach regions
void attach() const;
//- Return zone patch index
label zoneIndex() const;
//- Attach regions
void detach() const;
// Return true if patch is considered a wall
bool isWall() const
{
return isWall_;
}
//- Return shadow patch
const regionCouplePolyPatch& shadow() const;
//- Return shadow patch
const regionCouplePolyPatch& shadow() const;
//- Interpolate face field
template<class Type>
tmp<Field<Type> > interpolate(const Field<Type>& pf) const
{
return patchToPatch().faceInterpolate(pf);
}
//- Return interpolation face zone
const faceZone& zone() const;
template<class Type>
tmp<Field<Type> > interpolate(const tmp<Field<Type> >& tpf) const
{
return patchToPatch().faceInterpolate(tpf);
}
//- Return reconstructed cell centres
const vectorField& reconFaceCellCentres() const;
// Edit state
//- Force calculation of transformation tensors
void calcTransformTensors
(
const vectorField& Cf,
const vectorField& Cr,
const vectorField& nf,
const vectorField& nr
) const;
//- Attach regions
void attach() const;
//- Initialize ordering for primitivePatch. Does not
// refer to *this (except for name() and type() etc.)
virtual void initOrder(const primitivePatch&) const;
//- Attach regions
void detach() const;
//- Return new ordering for primitivePatch.
// Ordering is -faceMap: for every face
// index of the new face -rotation: for every new face the clockwise
// shift of the original face. Return false if nothing changes
// (faceMap is identity, rotation is 0), true otherwise.
virtual bool order
(
const primitivePatch&,
labelList& faceMap,
labelList& rotation
) const;
//- Synchronise communications of ordering for primitivePatch
// Used in cases when no topological change happens locally,
// but is happening on other processors
virtual void syncOrder() const;
// Interpolation data
//- Is this the master side?
bool master() const
{
return master_;
}
//- Is this the slave side?
bool slave() const
{
return !master();
}
//- Use bridging to fix overlap error in interpolation
bool bridgeOverlap() const
{
return bridgeOverlap_;
}
//- Return patch-to-zone addressing
const labelList& zoneAddressing() const;
//- Return remote patch-to-zone addressing
const labelList& remoteZoneAddressing() const;
//- Is the patch localised on a single processor
bool localParallel() const;
//- Return reference to patch-to-patch interpolation
const ggiZoneInterpolation& patchToPatch() const;
// Interpolation functions
//- Expand face field to full zone size, including reduction
template<class Type>
tmp<Field<Type> > fastExpand(const Field<Type>& pf) const;
//- Interpolate face field
template<class Type>
tmp<Field<Type> > interpolate(const Field<Type>& pf) const;
template<class Type>
tmp<Field<Type> > interpolate(const tmp<Field<Type> >& tpf) const;
//- Bridge interpolated face field for uncovered faces
template<class Type>
void bridge
(
const Field<Type>& bridgeField,
Field<Type>& ff
) const;
// Geometric data
//- Return reconstructed cell centres
const vectorField& reconFaceCellCentres() const;
// Patch ordering
//- Initialize ordering for primitivePatch. Does not
// refer to *this (except for name() and type() etc.)
virtual void initOrder(const primitivePatch&) const;
//- Return new ordering for primitivePatch.
// Ordering is -faceMap: for every face
// index of the new face -rotation: for every new face the
// clockwise shift of the original face. Return false if nothing
// changes (faceMap is identity, rotation is 0), true otherwise.
virtual bool order
(
const primitivePatch&,
labelList& faceMap,
labelList& rotation
) const;
//- Synchronise communications of ordering for primitivePatch
// Used in cases when no topological change happens locally,
// but is happening on other processors
virtual void syncOrder() const;
//- Write
@ -320,6 +447,12 @@ public:
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#ifdef NoRepository
# include "regionCouplePolyPatchTemplates.C"
#endif
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //

350
src/OpenFOAM/meshes/polyMesh/polyPatches/constraint/regionCouple/regionCouplePolyPatchTemplates.C Normal file
View file

@ -0,0 +1,350 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright held by original author
\\/ M anipulation |
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
OpenFOAM is free software; you can redistribute it and/or modify it
under the terms of the GNU General Public License as published by the
Free Software Foundation; either version 2 of the License, or (at your
option) any later version.
OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
for more details.
You should have received a copy of the GNU General Public License
along with OpenFOAM; if not, write to the Free Software Foundation,
Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
Author
Hrvoje Jasak, Wikki Ltd. All rights reserved.
\*---------------------------------------------------------------------------*/
#include "regionCouplePolyPatch.H"
#include "OSspecific.H"
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
template<class Type>
Foam::tmp<Foam::Field<Type> > Foam::regionCouplePolyPatch::fastExpand
(
const Field<Type>& ff
) const
{
// Check and expand the field from patch size to zone size
// with communication
// Algorithm:
// 1) Master processor holds maps of all zone addressing (data provided)
// and all remote zone addressing (data required)
// 2) Each processor will send the locally active data to the master
// 3) Master assembles all the data
// 4) Master sends to all processors the data they need to receive
//
// Notes:
// A) If the size of zone addressing is zero, data is not sent
// B) Communicated data on each processor has the size of live faces
// C) Expanded data will be equal to actual data fronm other processors
// only for the faces marked in remote; for other faces, it will be
// equal to zero
// D) On processor zero, complete data is available
// HJ, 4/Jun/2011
if (ff.size() != size())
{
FatalErrorIn
(
"tmp<Field<Type> > regionCouplePolyPatch::fastExpand\n"
"(\n"
" const Field<Type>& ff\n"
") const"
) << "Incorrect patch field size. Field size: "
<< ff.size() << " patch size: " << size()
<< abort(FatalError);
}
if (localParallel())
{
FatalErrorIn
(
"tmp<Field<Type> > regionCouplePolyPatch::fastExpand"
"("
" const Field<Type>& ff"
") const"
) << "Requested expand on local parallel. This is not allowed"
<< abort(FatalError);
}
// Expand the field to zone size
tmp<Field<Type> > texpandField
(
new Field<Type>(zone().size(), pTraits<Type>::zero)
);
Field<Type>& expandField = texpandField();
if (Pstream::master())
{
// Insert master processor
const labelList& za = zoneAddressing();
forAll (za, i)
{
expandField[za[i]] = ff[i];
}
// Master receives and inserts data from all processors for which
// receiveAddr contains entries
for (label procI = 1; procI < Pstream::nProcs(); procI++)
{
const labelList& curRAddr = receiveAddr()[procI];
if (!curRAddr.empty())
{
Field<Type> receiveBuf(curRAddr.size());
// Opt: reconsider mode of communication
IPstream::read
(
Pstream::blocking,
procI,
reinterpret_cast<char*>(receiveBuf.begin()),
receiveBuf.byteSize()
);
// Insert received information
forAll (curRAddr, i)
{
expandField[curRAddr[i]] = receiveBuf[i];
}
}
}
// Expanded field complete, send required data to other processors
for (label procI = 1; procI < Pstream::nProcs(); procI++)
{
const labelList& curSAddr = shadow().sendAddr()[procI];
if (!curSAddr.empty())
{
Field<Type> sendBuf(curSAddr.size());
forAll (curSAddr, i)
{
sendBuf[i] = expandField[curSAddr[i]];
}
// Opt: reconsider mode of communication
OPstream::write
(
Pstream::blocking,
procI,
reinterpret_cast<const char*>(sendBuf.begin()),
sendBuf.byteSize()
);
}
}
}
else
{
// Send local data to master and receive remote data
// If patch is empty, communication is avoided
// HJ, 4/Jun/2011
if (size())
{
// Opt: reconsider mode of communication
OPstream::write
(
Pstream::blocking,
Pstream::masterNo(),
reinterpret_cast<const char*>(ff.begin()),
ff.byteSize()
);
}
// Prepare to receive remote data
const labelList& rza = shadow().remoteZoneAddressing();
if (!rza.empty())
{
Field<Type> receiveBuf(rza.size());
// Opt: reconsider mode of communication
IPstream::read
(
Pstream::blocking,
Pstream::masterNo(),
reinterpret_cast<char*>(receiveBuf.begin()),
receiveBuf.byteSize()
);
// Insert the data into expanded field
forAll (rza, i)
{
expandField[rza[i]] = receiveBuf[i];
}
}
}
return texpandField;
}
template<class Type>
Foam::tmp<Foam::Field<Type> > Foam::regionCouplePolyPatch::interpolate
(
const Field<Type>& ff
) const
{
// Check and expand the field from patch size to zone size
if (ff.size() != shadow().size())
{
FatalErrorIn
(
"tmp<Field<Type> > regionCouplePolyPatch::interpolate\n"
"(\n"
" const Field<Type>& ff\n"
") const"
) << "Incorrect slave patch field size. Field size: "
<< ff.size() << " patch size: " << shadow().size()
<< abort(FatalError);
}
if (localParallel())
{
// No expansion or filtering needed. HJ, 4/Jun/2011
if (empty())
{
// Patch empty, no interpolation
return tmp<Field<Type> >(new Field<Type>());
}
// Interpolate field
if (master())
{
return patchToPatch().slaveToMaster(ff);
}
else
{
return patchToPatch().masterToSlave(ff);
}
}
else
{
// Expand shadow
Field<Type> expandField = shadow().fastExpand(ff);
tmp<Field<Type> > tresult(new Field<Type>(size()));
Field<Type>& result = tresult();
if (master())
{
patchToPatch().maskedSlaveToMaster
(
expandField,
result,
zoneAddressing()
);
}
else
{
patchToPatch().maskedMasterToSlave
(
expandField,
result,
zoneAddressing()
);
}
return tresult;
}
}
template<class Type>
Foam::tmp<Foam::Field<Type> > Foam::regionCouplePolyPatch::interpolate
(
const tmp<Field<Type> >& tff
) const
{
tmp<Field<Type> > tint = interpolate(tff());
tff.clear();
return tint;
}
template<class Type>
void Foam::regionCouplePolyPatch::bridge
(
const Field<Type>& bridgeField,
Field<Type>& ff
) const
{
// Check and expand the field from patch size to zone size
if (ff.size() != size())
{
FatalErrorIn
(
"tmp<Field<Type> > regionCouplePolyPatch::bridge\n"
"(\n"
" const Field<Type>& ff,\n"
" Field<Type>& ff\n"
") const"
) << "Incorrect patch field size for bridge. Field size: "
<< ff.size() << " patch size: " << size()
<< abort(FatalError);
}
if (bridgeOverlap())
{
if (empty())
{
// Patch empty, no bridging
return;
}
if (localParallel())
{
if (master())
{
patchToPatch().bridgeMaster(bridgeField, ff);
}
else
{
patchToPatch().bridgeSlave(bridgeField, ff);
}
}
else
{
// Note: since bridging is only a local operation
if (master())
{
patchToPatch().maskedBridgeMaster
(
bridgeField,
ff,
zoneAddressing()
);
}
else
{
patchToPatch().maskedBridgeSlave
(
bridgeField,
ff,
zoneAddressing()
);
}
}
}
}
// ************************************************************************* //

5
src/finiteVolume/Make/files
View file

@ -110,7 +110,7 @@ $(constraintFvPatchFields)/wedge/wedgeFvPatchScalarField.C
$(constraintFvPatchFields)/ggi/ggiFvPatchFields.C
$(constraintFvPatchFields)/cyclicGgi/cyclicGgiFvPatchFields.C
$(constraintFvPatchFields)/overlapGgi/overlapGgiFvPatchFields.C
$(constraintFvPatchFields)/regionCouple/regionCoupleFvPatchFields.C
$(constraintFvPatchFields)/regionCoupling/regionCouplingFvPatchFields.C
derivedFvPatchFields = $(fvPatchFields)/derived
$(derivedFvPatchFields)/activeBaffleVelocity/activeBaffleVelocityFvPatchVectorField.C
@ -185,7 +185,7 @@ $(constraintFvsPatchFields)/wedge/wedgeFvsPatchFields.C
$(constraintFvsPatchFields)/ggi/ggiFvsPatchFields.C
$(constraintFvsPatchFields)/cyclicGgi/cyclicGgiFvsPatchFields.C
$(constraintFvsPatchFields)/overlapGgi/overlapGgiFvsPatchFields.C
$(constraintFvsPatchFields)/regionCouple/regionCoupleFvsPatchFields.C
$(constraintFvsPatchFields)/regionCoupling/regionCouplingFvsPatchFields.C
fields/volFields/volFields.C
@ -227,6 +227,7 @@ $(schemes)/skewCorrected/skewCorrected.C
$(schemes)/outletStabilised/outletStabilised.C
$(schemes)/reverseLinear/reverseLinear.C
$(schemes)/clippedLinear/clippedLinear.C
$(schemes)/harmonic/magLongDelta.C
$(schemes)/harmonic/harmonic.C
$(schemes)/fixedBlended/fixedBlended.C
$(schemes)/localBlended/localBlended.C

12
src/finiteVolume/fields/fvPatchFields/constraint/ggi/ggiFvPatchField.C
View file

@ -140,6 +140,15 @@ ggiFvPatchField<Type>::ggiFvPatchField
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
// Return shadow field
template<class Type>
const Foam::ggiFvPatchField<Type>&
ggiFvPatchField<Type>::shadowPatchField() const
{
return this->boundaryField()[ggiPatch_.shadowIndex()];
}
// Return neighbour field
template<class Type>
tmp<Field<Type> > ggiFvPatchField<Type>::patchNeighbourField() const
@ -147,6 +156,9 @@ tmp<Field<Type> > ggiFvPatchField<Type>::patchNeighbourField() const
const Field<Type>& iField = this->internalField();
// Get shadow face-cells and assemble shadow field
// This is a patchInternalField of neighbour but access is inconvenient.
// Assemble by hand.
// HJ, 27/Sep/2011
const unallocLabelList& sfc = ggiPatch_.shadow().faceCells();
Field<Type> sField(sfc.size());

6
src/finiteVolume/fields/fvPatchFields/constraint/ggi/ggiFvPatchField.H
View file

@ -127,6 +127,12 @@ public:
// Member functions
// Access
//- Return shadow patch field
const ggiFvPatchField<Type>& shadowPatchField() const;
// Evaluation functions
//- Return neighbour field given internal cell data

275
src/finiteVolume/fields/fvPatchFields/constraint/regionCouple/regionCoupleFvPatchField.C
View file

@ -1,275 +0,0 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright held by original author
\\/ M anipulation |
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
OpenFOAM is free software; you can redistribute it and/or modify it
under the terms of the GNU General Public License as published by the
Free Software Foundation; either version 2 of the License, or (at your
option) any later version.
OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
for more details.
You should have received a copy of the GNU General Public License
along with OpenFOAM; if not, write to the Free Software Foundation,
Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
Author
Hrvoje Jasak, Wikki Ltd. All rights reserved
\*---------------------------------------------------------------------------*/
#include "regionCoupleFvPatchField.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
template<class Type>
regionCoupleFvPatchField<Type>::regionCoupleFvPatchField
(
const fvPatch& p,
const DimensionedField<Type, volMesh>& iF
)
:
coupledFvPatchField<Type>(p, iF),
regionCouplePatch_(refCast<const regionCoupleFvPatch>(p)),
remoteFieldName_(iF.name()),
matrixUpdateBuffer_()
{}
template<class Type>
regionCoupleFvPatchField<Type>::regionCoupleFvPatchField
(
const fvPatch& p,
const DimensionedField<Type, volMesh>& iF,
const dictionary& dict
)
:
coupledFvPatchField<Type>(p, iF, dict),
regionCouplePatch_(refCast<const regionCoupleFvPatch>(p)),
remoteFieldName_(dict.lookup("remoteField")),
matrixUpdateBuffer_()
{
if (!isType<regionCoupleFvPatch>(p))
{
FatalIOErrorIn
(
"regionCoupleFvPatchField<Type>::regionCoupleFvPatchField\n"
"(\n"
" const fvPatch& p,\n"
" const DimensionedField<Type, volMesh>& iF,\n"
" const dictionary& dict\n"
")\n",
dict
) << "patch " << this->patch().index() << " not regionCouple type. "
<< "Patch type = " << p.type()
<< exit(FatalIOError);
}
}
template<class Type>
regionCoupleFvPatchField<Type>::regionCoupleFvPatchField
(
const regionCoupleFvPatchField<Type>& ptf,
const fvPatch& p,
const DimensionedField<Type, volMesh>& iF,
const fvPatchFieldMapper& mapper
)
:
coupledFvPatchField<Type>(ptf, p, iF, mapper),
regionCouplePatch_(refCast<const regionCoupleFvPatch>(p)),
remoteFieldName_(ptf.remoteFieldName_),
matrixUpdateBuffer_()
{
if (!isType<regionCoupleFvPatch>(this->patch()))
{
FatalErrorIn
(
"regionCoupleFvPatchField<Type>::regionCoupleFvPatchField\n"
"(\n"
" const regionCoupleFvPatchField<Type>& ptf,\n"
" const fvPatch& p,\n"
" const DimensionedField<Type, volMesh>& iF,\n"
" const fvPatchFieldMapper& mapper\n"
")\n"
) << "Field type does not correspond to patch type for patch "
<< this->patch().index() << "." << endl
<< "Field type: " << typeName << endl
<< "Patch type: " << this->patch().type()
<< exit(FatalError);
}
}
template<class Type>
regionCoupleFvPatchField<Type>::regionCoupleFvPatchField
(
const regionCoupleFvPatchField<Type>& ptf,
const DimensionedField<Type, volMesh>& iF
)
:
regionCoupleLduInterfaceField(),
coupledFvPatchField<Type>(ptf, iF),
regionCouplePatch_(refCast<const regionCoupleFvPatch>(ptf.patch())),
remoteFieldName_(ptf.remoteFieldName_),
matrixUpdateBuffer_()
{}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
// Return neighbour field
template<class Type>
const regionCoupleFvPatchField<Type>&
regionCoupleFvPatchField<Type>::shadowPatchField() const
{
// Lookup neighbour field
typedef GeometricField<Type, fvPatchField, volMesh> GeoField;
const GeoField& coupleField =
regionCouplePatch_.shadowRegion().
objectRegistry::lookupObject<GeoField>(remoteFieldName_);
return refCast<const regionCoupleFvPatchField<Type> >
(
coupleField.boundaryField()[regionCouplePatch_.shadowIndex()]
);
}
// Return neighbour field
template<class Type>
tmp<Field<Type> > regionCoupleFvPatchField<Type>::patchNeighbourField() const
{
return regionCouplePatch_.interpolate
(
shadowPatchField().patchInternalField()
);
}
template<class Type>
void regionCoupleFvPatchField<Type>::evaluate
(
const Pstream::commsTypes
)
{
// 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
Field<Type> fOwn = this->patchInternalField();
Field<Type> fNei = this->patchNeighbourField();
scalarField magFOwn = mag(fOwn);
scalarField magFNei = mag(fNei);
// Calculate internal weights using field magnitude
scalarField weights(fOwn.size());
forAll (weights, faceI)
{
scalar mOwn = magFOwn[faceI];
scalar mNei = magFNei[faceI];
scalar den1 = mOwn + mNei;
scalar den2 = mOwn - mNei;
// Strange optimisation check: floating point failure if only
// checking for den2. den1 check should not be needed
// HJ, 24/Aug/2011
if (mag(den1) > SMALL && mag(den2) > SMALL)
{
scalar mean = 2.0*mOwn*mNei/den1;
weights[faceI] = (mean - mNei)/den2;
}
else
{
weights[faceI] = 0.5;
}
}
// Do interpolation
Field<Type>::operator=(weights*fOwn + (1.0 - weights)*fNei);
}
// Initialise neighbour processor internal cell data
template<class Type>
void regionCoupleFvPatchField<Type>::initInterfaceMatrixUpdate
(
const scalarField& psiInternal,
scalarField&,
const lduMatrix&,
const scalarField&,
const direction,
const Pstream::commsTypes
) const
{
matrixUpdateBuffer_ = this->patch().patchInternalField(psiInternal);
}
// Return matrix product for coupled boundary
template<class Type>
void regionCoupleFvPatchField<Type>::updateInterfaceMatrix
(
const scalarField& psiInternal,
scalarField& result,
const lduMatrix&,
const scalarField& coeffs,
const direction cmpt,
const Pstream::commsTypes
) const
{
scalarField pnf =
regionCouplePatch_.interpolate
(
this->shadowPatchField().matrixUpdateBuffer()
);
// Multiply the field by coefficients and add into the result
const unallocLabelList& fc = regionCouplePatch_.faceCells();
forAll(fc, elemI)
{
result[fc[elemI]] -= coeffs[elemI]*pnf[elemI];
}
}
// Write
template<class Type>
void regionCoupleFvPatchField<Type>::write(Ostream& os) const
{
fvPatchField<Type>::write(os);
os.writeKeyword("remoteField")
<< remoteFieldName_ << token::END_STATEMENT << nl;
this->writeEntry("value", os);
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam
// ************************************************************************* //

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

@ -0,0 +1,468 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright held by original author
\\/ M anipulation |
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
OpenFOAM is free software; you can redistribute it and/or modify it
under the terms of the GNU General Public License as published by the
Free Software Foundation; either version 2 of the License, or (at your
option) any later version.
OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
for more details.
You should have received a copy of the GNU General Public License
along with OpenFOAM; if not, write to the Free Software Foundation,
Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
Author
Hrvoje Jasak, Wikki Ltd. All rights reserved
\*---------------------------------------------------------------------------*/
#include "regionCouplingFvPatchField.H"
#include "symmTransformField.H"
#include "magLongDelta.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
template<class Type>
regionCouplingFvPatchField<Type>::regionCouplingFvPatchField
(
const fvPatch& p,
const DimensionedField<Type, volMesh>& iF
)
:
coupledFvPatchField<Type>(p, iF),
regionCouplePatch_(refCast<const regionCoupleFvPatch>(p)),
remoteFieldName_(iF.name()),
matrixUpdateBuffer_(),
originalPatchField_(),
curTimeIndex_(-1)
{}
template<class Type>
regionCouplingFvPatchField<Type>::regionCouplingFvPatchField
(
const fvPatch& p,
const DimensionedField<Type, volMesh>& iF,
const dictionary& dict
)
:
coupledFvPatchField<Type>(p, iF, dict),
regionCouplePatch_(refCast<const regionCoupleFvPatch>(p)),
remoteFieldName_(dict.lookup("remoteField")),
matrixUpdateBuffer_(),
originalPatchField_(),
curTimeIndex_(-1)
{
if (!isType<regionCoupleFvPatch>(p))
{
FatalIOErrorIn
(
"regionCouplingFvPatchField<Type>::regionCouplingFvPatchField\n"
"(\n"
" const fvPatch& p,\n"
" const DimensionedField<Type, volMesh>& iF,\n"
" const dictionary& dict\n"
")\n",
dict
) << "patch " << this->patch().index() << " not regionCouple type. "
<< "Patch type = " << p.type()
<< exit(FatalIOError);
}
if (!dict.found("value"))
{
// Grab the internal value for initialisation. (?) HJ, 27/Feb/2009
fvPatchField<Type>::operator=(this->patchInternalField()());
}
}
template<class Type>
regionCouplingFvPatchField<Type>::regionCouplingFvPatchField
(
const regionCouplingFvPatchField<Type>& ptf,
const fvPatch& p,
const DimensionedField<Type, volMesh>& iF,
const fvPatchFieldMapper& mapper
)
:
coupledFvPatchField<Type>(ptf, p, iF, mapper),
regionCouplePatch_(refCast<const regionCoupleFvPatch>(p)),
remoteFieldName_(ptf.remoteFieldName_),
matrixUpdateBuffer_(),
originalPatchField_(),
curTimeIndex_(-1)
{
if (!isType<regionCoupleFvPatch>(this->patch()))
{
FatalErrorIn
(
"regionCouplingFvPatchField<Type>::regionCouplingFvPatchField\n"
"(\n"
" const regionCouplingFvPatchField<Type>& ptf,\n"
" const fvPatch& p,\n"
" const DimensionedField<Type, volMesh>& iF,\n"
" const fvPatchFieldMapper& mapper\n"
")\n"
) << "Field type does not correspond to patch type for patch "
<< this->patch().index() << "." << endl
<< "Field type: " << typeName << endl
<< "Patch type: " << this->patch().type()
<< exit(FatalError);
}
}
template<class Type>
regionCouplingFvPatchField<Type>::regionCouplingFvPatchField
(
const regionCouplingFvPatchField<Type>& ptf,
const DimensionedField<Type, volMesh>& iF
)
:
ggiLduInterfaceField(),
coupledFvPatchField<Type>(ptf, iF),
regionCouplePatch_(refCast<const regionCoupleFvPatch>(ptf.patch())),
remoteFieldName_(ptf.remoteFieldName_),
matrixUpdateBuffer_(),
originalPatchField_(),
curTimeIndex_(-1)
{}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
// Return a named shadow patch field
template<class Type>
template<class GeometricField, class Type2>
const typename GeometricField::PatchFieldType&
regionCouplingFvPatchField<Type>::lookupShadowPatchField
(
const word& name,
const GeometricField*,
const Type2*
) const
{
// Lookup neighbour field
const GeometricField& shadowField =
regionCouplePatch_.shadowRegion().
objectRegistry::lookupObject<GeometricField>(name);
return shadowField.boundaryField()[regionCouplePatch_.shadowIndex()];
}
// Return shadow patch field
template<class Type>
const regionCouplingFvPatchField<Type>&
regionCouplingFvPatchField<Type>::shadowPatchField() const
{
// Lookup neighbour field
typedef GeometricField<Type, fvPatchField, volMesh> GeoField;
return refCast<const regionCouplingFvPatchField<Type> >
(
lookupShadowPatchField<GeoField, Type>(remoteFieldName_)
);
}
// Return neighbour field
template<class Type>
tmp<Field<Type> > regionCouplingFvPatchField<Type>::patchNeighbourField() const
{
Field<Type> sField = shadowPatchField().patchInternalField();
tmp<Field<Type> > tpnf
(
regionCouplePatch_.interpolate
(
shadowPatchField().patchInternalField()
)
);
Field<Type>& pnf = tpnf();
if (regionCouplePatch_.bridgeOverlap())
{
// Symmetry treatment used for overlap
vectorField nHat = this->patch().nf();
// Use mirrored neighbour field for interpolation
// HJ, 21/Jan/2009
Field<Type> bridgeField =
transform(I - 2.0*sqr(nHat), this->patchInternalField());
regionCouplePatch_.bridge(bridgeField, pnf);
}
return tpnf;
}
// Return neighbour field given internal cell data
template<class Type>
tmp<Field<Type> > regionCouplingFvPatchField<Type>::patchNeighbourField
(
const word& name
) const
{
// Lookup neighbour field
typedef GeometricField<Type, fvPatchField, volMesh> GeoField;
return regionCouplePatch_.interpolate
(
lookupShadowPatchField<GeoField, Type>(name).patchInternalField()
);
// Note: this field is not bridged because local data does not exist
// for named field. HJ, 27/Sep/2011
}
template<class Type>
void regionCouplingFvPatchField<Type>::initEvaluate
(
const Pstream::commsTypes commsType
)
{
if (!this->updated())
{
this->updateCoeffs();
}
// 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
const Field<Type>& fOwn = this->originalPatchField();
Field<Type> fNei = regionCouplePatch_.interpolate
(
this->shadowPatchField().originalPatchField()
);
// Larger small for complex arithmetic accuracy
const scalar kSmall = 1000*SMALL;
# if 0
// Hrv's treatment
scalarField mOwn = mag(fOwn);
scalarField mNei = mag(fNei);
scalarField mean = 2*(mOwn*mNei)/(mOwn + mNei);
scalarField weights(fOwn.size(), 0.5);
scalar den;
forAll (weights, faceI)
{
den = (mNei[faceI] - mOwn[faceI]);
// Note: complex arithmetic requires extra accuracy
// This is a division of two close subtractions
// HJ, 28/Sep/2011
if (mag(den) > kSmall)
{
weights[faceI] = (mNei[faceI] - mean[faceI])/den;
}
else
{
// Use 0.5 weights
}
}
# else
// Henrik's treatment
const fvPatch& p = this->patch();
// Note: for interpolation, work with face fields, to allow wall-corrected
// diffusivity (eg wall functions) to operate correctly.
// HJ, 28/Sep/2011
Field<Type> f = *this;
const regionCouplingFvPatchField<Type>& spf = shadowPatchField();
const fvPatch& sp = spf.patch();
// Mag long deltas are identical on both sides. HJ, 28/Sep/2011
const magLongDelta& mld = magLongDelta::New(p.boundaryMesh().mesh());
scalarField magPhiOwn = mag(fOwn);
scalarField magPhiNei = mag(fNei);
const scalarField& pWeights = p.weights();
const scalarField& pDeltaCoeffs = p.deltaCoeffs();
const scalarField& pLongDelta = mld.magDelta(p.index());
// Calculate internal weights using field magnitude
scalarField weights(fOwn.size());
forAll (weights, faceI)
{
scalar mOwn = magPhiOwn[faceI]/(1 - pWeights[faceI]);
scalar mNei = magPhiNei[faceI]/pWeights[faceI];
scalar den = magPhiNei[faceI] - magPhiOwn[faceI];
// Note: complex arithmetic requires extra accuracy
// This is a division of two close subtractions
// HJ, 28/Sep/2011
if (mag(den) > kSmall)
{
scalar mean = mOwn*mNei/
(
(mOwn + mNei)*
pLongDelta[faceI]*
pDeltaCoeffs[faceI]
);
weights[faceI] = (magPhiNei[faceI] - mean)/den;
}
else
{
weights[faceI] = 0.5;
}
}
#endif
// Do interpolation
Field<Type>::operator=(weights*fOwn + (1.0 - weights)*fNei);
}
template<class Type>
void regionCouplingFvPatchField<Type>::evaluate
(
const Pstream::commsTypes
)
{
// No interpolation allowed
fvPatchField<Type>::evaluate();
}
template<class Type>
void regionCouplingFvPatchField<Type>::updateCoeffs()
{
if (this->updated())
{
return;
}
// Store original field for symmetric evaluation
// Henrik Rusche, Aug/2011
if (curTimeIndex_ != this->db().time().timeIndex())
{
originalPatchField_ = *this;
curTimeIndex_ = this->db().time().timeIndex();
}
}
template<class Type>
tmp<Field<Type> > regionCouplingFvPatchField<Type>::snGrad() const
{
if (regionCouplePatch_.coupled())
{
return coupledFvPatchField<Type>::snGrad();
}
else
{
return fvPatchField<Type>::snGrad();
}
}
// Initialise neighbour processor internal cell data
template<class Type>
void regionCouplingFvPatchField<Type>::initInterfaceMatrixUpdate
(
const scalarField& psiInternal,
scalarField& result,
const lduMatrix&,
const scalarField& coeffs,
const direction cmpt,
const Pstream::commsTypes
) const
{
// Prepare local matrix update buffer for the remote side.
// Note that only remote side will have access to its psiInternal
// as they are on different regions
// Since interpolation needs to happen on the shadow, and within the
// init, prepare interpolation for the other side.
matrixUpdateBuffer_ =
this->shadowPatchField().regionCouplePatch().interpolate
(
this->patch().patchInternalField(psiInternal)
);
}
// Return matrix product for coupled boundary
template<class Type>
void regionCouplingFvPatchField<Type>::updateInterfaceMatrix
(
const scalarField& psiInternal,
scalarField& result,
const lduMatrix&,
const scalarField& coeffs,
const direction ,
const Pstream::commsTypes
) const
{
// Note: interpolation involves parallel communications and needs to
// happen during init. This changes the use of matrix update buffer
// compared to earlier versions
// HJ, 28/Sep/2011
scalarField pnf = this->shadowPatchField().matrixUpdateBuffer();
// Multiply the field by coefficients and add into the result
const unallocLabelList& fc = regionCouplePatch_.faceCells();
forAll(fc, elemI)
{
result[fc[elemI]] -= coeffs[elemI]*pnf[elemI];
}
}
// Write
template<class Type>
void regionCouplingFvPatchField<Type>::write(Ostream& os) const
{
fvPatchField<Type>::write(os);
os.writeKeyword("remoteField")
<< remoteFieldName_ << token::END_STATEMENT << nl;
this->writeEntry("value", os);
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam
// ************************************************************************* //

127
src/finiteVolume/fields/fvPatchFields/constraint/regionCouple/regionCoupleFvPatchField.H → src/finiteVolume/fields/fvPatchFields/constraint/regionCoupling/regionCouplingFvPatchField.H
View file

@ -23,7 +23,7 @@ License
Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
Class
regionCoupleFvPatchField
regionCouplingFvPatchField
Description
Region couple patch field
@ -32,15 +32,15 @@ Author
Hrvoje Jasak, Wikki Ltd. All rights reserved
SourceFiles
regionCoupleFvPatchField.C
regionCouplingFvPatchField.C
\*---------------------------------------------------------------------------*/
#ifndef regionCoupleFvPatchField_H
#define regionCoupleFvPatchField_H
#ifndef regionCouplingFvPatchField_H
#define regionCouplingFvPatchField_H
#include "coupledFvPatchField.H"
#include "regionCoupleLduInterfaceField.H"
#include "ggiLduInterfaceField.H"
#include "regionCoupleFvPatch.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@ -49,13 +49,13 @@ namespace Foam
{
/*---------------------------------------------------------------------------*\
Class regionCoupleFvPatchField Declaration
Class regionCouplingFvPatchField Declaration
\*---------------------------------------------------------------------------*/
template<class Type>
class regionCoupleFvPatchField
class regionCouplingFvPatchField
:
public regionCoupleLduInterfaceField,
public ggiLduInterfaceField,
public coupledFvPatchField<Type>
{
// Private data
@ -69,9 +69,29 @@ class regionCoupleFvPatchField
//- Matrix update buffer
mutable scalarField matrixUpdateBuffer_;
//- Original patch field. Required for correct evaluation
// in harmonic averaging
Field<Type> originalPatchField_;
//- Current time index used to store originalPatchField_
label curTimeIndex_;
protected:
//- Return original patch field
const Field<Type>& originalPatchField() const
{
if (curTimeIndex_ == this->db().time().timeIndex())
{
return originalPatchField_;
}
else
{
return *this;
}
}
//- Return contents of a matrix update buffer
const scalarField& matrixUpdateBuffer() const
{
@ -82,40 +102,40 @@ protected:
public:
//- Runtime type information
// TypeName(regionCoupleFvPatch::typeName_());
TypeName("regionCoupling");
// Constructors
//- Construct from patch and internal field
regionCoupleFvPatchField
regionCouplingFvPatchField
(
const fvPatch&,
const DimensionedField<Type, volMesh>&
);
//- Construct from patch, internal field and dictionary
regionCoupleFvPatchField
regionCouplingFvPatchField
(
const fvPatch&,
const DimensionedField<Type, volMesh>&,
const dictionary&
);
//- Construct by mapping given regionCoupleFvPatchField onto a new patch
regionCoupleFvPatchField
//- Construct by mapping given regionCouplingFvPatchField
// onto a new patch
regionCouplingFvPatchField
(
const regionCoupleFvPatchField<Type>&,
const regionCouplingFvPatchField<Type>&,
const fvPatch&,
const DimensionedField<Type, volMesh>&,
const fvPatchFieldMapper&
);
//- Construct as copy setting internal field reference
regionCoupleFvPatchField
regionCouplingFvPatchField
(
const regionCoupleFvPatchField<Type>&,
const regionCouplingFvPatchField<Type>&,
const DimensionedField<Type, volMesh>&
);
@ -124,7 +144,7 @@ public:
{
return tmp<fvPatchField<Type> >
(
new regionCoupleFvPatchField<Type>(*this)
new regionCouplingFvPatchField<Type>(*this)
);
}
@ -136,14 +156,14 @@ public:
{
return tmp<fvPatchField<Type> >
(
new regionCoupleFvPatchField<Type>(*this, iF)
new regionCouplingFvPatchField<Type>(*this, iF)
);
}
// Member functions
// Evaluation functions
// Access functions
//- Return remote field name
const word& remoteFieldName() const
@ -151,18 +171,52 @@ public:
return remoteFieldName_;
}
//- Return a named shadow patch field
template<class GeometricField, class Type2>
const typename GeometricField::PatchFieldType&
lookupShadowPatchField
(
const word& name,
const GeometricField* = NULL,
const Type2* = NULL
) const;
const regionCoupleFvPatch& regionCouplePatch() const
{
return regionCouplePatch_;
}
//- Return shadow patch field
const regionCoupleFvPatchField<Type>& shadowPatchField() const;
const regionCouplingFvPatchField<Type>& shadowPatchField() const;
// Evaluation functions
//- Return neighbour field given internal cell data
virtual tmp<Field<Type> > patchNeighbourField() const;
//- Return neighbour field given internal cell data
virtual tmp<Field<Type> > patchNeighbourField
(
const word& name
) const;
//- Initialise the evaluation of the patch field
virtual void initEvaluate(const Pstream::commsTypes commsType);
//- Evaluate the patch field
virtual void evaluate
(
const Pstream::commsTypes commsType
);
//- Update the coefficients associated with the patch field
virtual void updateCoeffs();
//- Return patch-normal gradient
virtual tmp<Field<Type> > snGrad() const;
//- Initialise neighbour matrix update
virtual void initInterfaceMatrixUpdate
(
@ -186,6 +240,37 @@ public:
) const;
// GGI coupled interface functions
//- Does the patch field perform the transfromation
virtual bool doTransform() const
{
return
!(
regionCouplePatch_.parallel()
|| pTraits<Type>::rank == 0
);
}
//- Return face transformation tensor
virtual const tensorField& forwardT() const
{
return regionCouplePatch_.forwardT();
}
//- Return neighbour-cell transformation tensor
virtual const tensorField& reverseT() const
{
return regionCouplePatch_.reverseT();
}
//- Return rank of component for transform
virtual int rank() const
{
return pTraits<Type>::rank;
}
//- Write
virtual void write(Ostream&) const;
@ -199,7 +284,7 @@ public:
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#ifdef NoRepository
# include "regionCoupleFvPatchField.C"
# include "regionCouplingFvPatchField.C"
#endif
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

4
src/finiteVolume/fields/fvPatchFields/constraint/regionCouple/regionCoupleFvPatchFields.C → src/finiteVolume/fields/fvPatchFields/constraint/regionCoupling/regionCouplingFvPatchFields.C
View file

@ -30,7 +30,7 @@ Author
\*---------------------------------------------------------------------------*/
#include "regionCoupleFvPatchFields.H"
#include "regionCouplingFvPatchFields.H"
#include "addToRunTimeSelectionTable.H"
#include "volFields.H"
@ -41,7 +41,7 @@ namespace Foam
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
makePatchFields(regionCouple);
makePatchFields(regionCoupling);
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

12
src/finiteVolume/fields/fvsPatchFields/constraint/regionCouple/regionCoupleFvsPatchFields.H → src/finiteVolume/fields/fvPatchFields/constraint/regionCoupling/regionCouplingFvPatchFields.H
View file

@ -23,7 +23,7 @@ License
Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
Class
regionCouplewFvsPatchFields
regionCouplingFvPatchFields
Description
Region couple patch field
@ -32,14 +32,14 @@ Author
Hrvoje Jasak, Wikki Ltd. All rights reserved
SourceFiles
regionCoupleFvsPatchFields.C
regionCouplingFvPatchFields.C
\*---------------------------------------------------------------------------*/
#ifndef regionCoupleFvsPatchFields_H
#define regionCoupleFvsPatchFields_H
#ifndef regionCouplingFvPatchFields_H
#define regionCouplingFvPatchFields_H
#include "regionCoupleFvsPatchField.H"
#include "regionCouplingFvPatchField.H"
#include "fieldTypes.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@ -49,7 +49,7 @@ namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
makeFvsPatchTypeFieldTypedefs(regionCouple)
makePatchTypeFieldTypedefs(regionCoupling)
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

10
src/finiteVolume/fields/fvPatchFields/constraint/regionCouple/regionCoupleFvPatchFieldsFwd.H → src/finiteVolume/fields/fvPatchFields/constraint/regionCoupling/regionCouplingFvPatchFieldsFwd.H
View file

@ -23,7 +23,7 @@ License
Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
Class
regionCoupleFvPatchField
regionCouplingFvPatchField
Description
Region couple patch field
@ -33,8 +33,8 @@ Author
\*---------------------------------------------------------------------------*/
#ifndef regionCoupleFvPatchFieldsFwd_H
#define regionCoupleFvPatchFieldsFwd_H
#ifndef regionCouplingFvPatchFieldsFwd_H
#define regionCouplingFvPatchFieldsFwd_H
#include "fvPatchField.H"
#include "fieldTypes.H"
@ -46,9 +46,9 @@ namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
template<class Type> class regionCoupleFvPatchField;
template<class Type> class regionCouplingFvPatchField;
makePatchTypeFieldTypedefs(regionCouple)
makePatchTypeFieldTypedefs(regionCoupling)
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

28
src/finiteVolume/fields/fvsPatchFields/constraint/regionCouple/regionCoupleFvsPatchField.C → src/finiteVolume/fields/fvsPatchFields/constraint/regionCoupling/regionCouplingFvsPatchField.C
View file

@ -27,7 +27,7 @@ Author
\*---------------------------------------------------------------------------*/
#include "regionCoupleFvsPatchField.H"
#include "regionCouplingFvsPatchField.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@ -37,7 +37,7 @@ namespace Foam
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
template<class Type>
regionCoupleFvsPatchField<Type>::regionCoupleFvsPatchField
regionCouplingFvsPatchField<Type>::regionCouplingFvsPatchField
(
const fvPatch& p,
const DimensionedField<Type, surfaceMesh>& iF
@ -51,7 +51,7 @@ regionCoupleFvsPatchField<Type>::regionCoupleFvsPatchField
template<class Type>
regionCoupleFvsPatchField<Type>::regionCoupleFvsPatchField
regionCouplingFvsPatchField<Type>::regionCouplingFvsPatchField
(
const fvPatch& p,
const DimensionedField<Type, surfaceMesh>& iF,
@ -67,7 +67,7 @@ regionCoupleFvsPatchField<Type>::regionCoupleFvsPatchField
{
FatalIOErrorIn
(
"regionCoupleFvsPatchField<Type>::regionCoupleFvsPatchField\n"
"regionCouplingFvsPatchField<Type>::regionCouplingFvsPatchField\n"
"(\n"
" const fvPatch& p,\n"
" const DimensionedField<Type, surfaceMesh>& iF,\n"
@ -82,9 +82,9 @@ regionCoupleFvsPatchField<Type>::regionCoupleFvsPatchField
template<class Type>
regionCoupleFvsPatchField<Type>::regionCoupleFvsPatchField
regionCouplingFvsPatchField<Type>::regionCouplingFvsPatchField
(
const regionCoupleFvsPatchField<Type>& ptf,
const regionCouplingFvsPatchField<Type>& ptf,
const fvPatch& p,
const DimensionedField<Type, surfaceMesh>& iF,
const fvPatchFieldMapper& mapper
@ -99,9 +99,9 @@ regionCoupleFvsPatchField<Type>::regionCoupleFvsPatchField
{
FatalErrorIn
(
"regionCoupleFvsPatchField<Type>::regionCoupleFvsPatchField\n"
"regionCouplingFvsPatchField<Type>::regionCouplingFvsPatchField\n"
"(\n"
" const regionCoupleFvsPatchField<Type>& ptf,\n"
" const regionCouplingFvsPatchField<Type>& ptf,\n"
" const fvPatch& p,\n"
" const DimensionedField<Type, surfaceMesh>& iF,\n"
" const fvPatchFieldMapper& mapper\n"
@ -116,9 +116,9 @@ regionCoupleFvsPatchField<Type>::regionCoupleFvsPatchField
template<class Type>
regionCoupleFvsPatchField<Type>::regionCoupleFvsPatchField
regionCouplingFvsPatchField<Type>::regionCouplingFvsPatchField
(
const regionCoupleFvsPatchField<Type>& ptf,
const regionCouplingFvsPatchField<Type>& ptf,
const DimensionedField<Type, surfaceMesh>& iF
)
:
@ -133,8 +133,8 @@ regionCoupleFvsPatchField<Type>::regionCoupleFvsPatchField
// Return neighbour field
template<class Type>
const regionCoupleFvsPatchField<Type>&
regionCoupleFvsPatchField<Type>::shadowPatchField() const
const regionCouplingFvsPatchField<Type>&
regionCouplingFvsPatchField<Type>::shadowPatchField() const
{
// Lookup neighbour field
typedef GeometricField<Type, fvsPatchField, surfaceMesh> GeoField;
@ -143,7 +143,7 @@ regionCoupleFvsPatchField<Type>::shadowPatchField() const
regionCouplePatch_.shadowRegion().
objectRegistry::lookupObject<GeoField>(remoteFieldName_);
return refCast<const regionCoupleFvsPatchField<Type> >
return refCast<const regionCouplingFvsPatchField<Type> >
(
coupleField.boundaryField()[regionCouplePatch_.shadowIndex()]
);
@ -152,7 +152,7 @@ regionCoupleFvsPatchField<Type>::shadowPatchField() const
// Write
template<class Type>
void regionCoupleFvsPatchField<Type>::write(Ostream& os) const
void regionCouplingFvsPatchField<Type>::write(Ostream& os) const
{
fvsPatchField<Type>::write(os);
os.writeKeyword("remoteField")

41
src/finiteVolume/fields/fvsPatchFields/constraint/regionCouple/regionCoupleFvsPatchField.H → src/finiteVolume/fields/fvsPatchFields/constraint/regionCoupling/regionCouplingFvsPatchField.H
View file

@ -23,7 +23,7 @@ License
Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
Class
regionCoupleFvsPatchField
regionCouplingFvsPatchField
Description
Region couple patch field
@ -31,13 +31,17 @@ Description
Author
Hrvoje Jasak, Wikki Ltd. All rights reserved
Note
The name of the class needs to be different from the name of the patch
(regionCouple) to avoid constraint patch behaviour
SourceFiles
regionCoupleFvsPatchField.C
regionCouplingFvsPatchField.C
\*---------------------------------------------------------------------------*/
#ifndef regionCoupleFvsPatchField_H
#define regionCoupleFvsPatchField_H
#ifndef regionCouplingFvsPatchField_H
#define regionCouplingFvsPatchField_H
#include "coupledFvsPatchField.H"
#include "regionCoupleLduInterfaceField.H"
@ -49,17 +53,17 @@ namespace Foam
{
/*---------------------------------------------------------------------------*\
Class regionCoupleFvsPatchField Declaration
Class regionCouplingFvsPatchField Declaration
\*---------------------------------------------------------------------------*/
template<class Type>
class regionCoupleFvsPatchField
class regionCouplingFvsPatchField
:
public coupledFvsPatchField<Type>
{
// Private data
//- Local reference cast into the regionCouple patch
//- Local reference cast into the regionCoupling patch
const regionCoupleFvPatch& regionCouplePatch_;
//- Name of remote field to couple to
@ -81,41 +85,40 @@ protected:
public:
//- Runtime type information
// TypeName(regionCoupleFvPatch::typeName_());
TypeName("regionCoupling");
// Constructors
//- Construct from patch and internal field
regionCoupleFvsPatchField
regionCouplingFvsPatchField
(
const fvPatch&,
const DimensionedField<Type, surfaceMesh>&
);
//- Construct from patch, internal field and dictionary
regionCoupleFvsPatchField
regionCouplingFvsPatchField
(
const fvPatch&,
const DimensionedField<Type, surfaceMesh>&,
const dictionary&
);
//- Construct by mapping given regionCoupleFvsPatchField
//- Construct by mapping given regionCouplingFvsPatchField
// onto a new patch
regionCoupleFvsPatchField
regionCouplingFvsPatchField
(
const regionCoupleFvsPatchField<Type>&,
const regionCouplingFvsPatchField<Type>&,
const fvPatch&,
const DimensionedField<Type, surfaceMesh>&,
const fvPatchFieldMapper&
);
//- Construct as copy setting internal field reference
regionCoupleFvsPatchField
regionCouplingFvsPatchField
(
const regionCoupleFvsPatchField<Type>&,
const regionCouplingFvsPatchField<Type>&,
const DimensionedField<Type, surfaceMesh>&
);
@ -124,7 +127,7 @@ public:
{
return tmp<fvsPatchField<Type> >
(
new regionCoupleFvsPatchField<Type>(*this)
new regionCouplingFvsPatchField<Type>(*this)
);
}
@ -136,7 +139,7 @@ public:
{
return tmp<fvsPatchField<Type> >
(
new regionCoupleFvsPatchField<Type>(*this, iF)
new regionCouplingFvsPatchField<Type>(*this, iF)
);
}
@ -152,7 +155,7 @@ public:
}
//- Return shadow patch field
const regionCoupleFvsPatchField<Type>& shadowPatchField() const;
const regionCouplingFvsPatchField<Type>& shadowPatchField() const;
//- Write
@ -168,7 +171,7 @@ public:
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#ifdef NoRepository
# include "regionCoupleFvsPatchField.C"
# include "regionCouplingFvsPatchField.C"
#endif
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

4
src/finiteVolume/fields/fvsPatchFields/constraint/regionCouple/regionCoupleFvsPatchFields.C → src/finiteVolume/fields/fvsPatchFields/constraint/regionCoupling/regionCouplingFvsPatchFields.C
View file

@ -30,7 +30,7 @@ Author
\*---------------------------------------------------------------------------*/
#include "regionCoupleFvsPatchFields.H"
#include "regionCouplingFvsPatchFields.H"
#include "fvsPatchFields.H"
#include "addToRunTimeSelectionTable.H"
@ -41,7 +41,7 @@ namespace Foam
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
makeFvsPatchFields(regionCouple);
makeFvsPatchFields(regionCoupling);
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

12
src/finiteVolume/fields/fvPatchFields/constraint/regionCouple/regionCoupleFvPatchFields.H → src/finiteVolume/fields/fvsPatchFields/constraint/regionCoupling/regionCouplingFvsPatchFields.H
View file

@ -23,7 +23,7 @@ License
Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
Class
regionCouplewFvPatchFields
regionCouplingFvsPatchFields
Description
Region couple patch field
@ -32,14 +32,14 @@ Author
Hrvoje Jasak, Wikki Ltd. All rights reserved
SourceFiles
regionCoupleFvPatchFields.C
regionCouplingFvsPatchFields.C
\*---------------------------------------------------------------------------*/
#ifndef regionCoupleFvPatchFields_H
#define regionCoupleFvPatchFields_H
#ifndef regionCouplingFvsPatchFields_H
#define regionCouplingFvsPatchFields_H
#include "regionCoupleFvPatchField.H"
#include "regionCouplingFvsPatchField.H"
#include "fieldTypes.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@ -49,7 +49,7 @@ namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
makePatchTypeFieldTypedefs(regionCouple)
makeFvsPatchTypeFieldTypedefs(regionCoupling)
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

10
src/finiteVolume/fields/fvsPatchFields/constraint/regionCouple/regionCoupleFvsPatchFieldsFwd.H → src/finiteVolume/fields/fvsPatchFields/constraint/regionCoupling/regionCouplingFvsPatchFieldsFwd.H
View file

@ -23,7 +23,7 @@ License
Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
Class
regionCoupleFvsPatchField
regionCouplingFvsPatchField
Description
Region couple patch field
@ -33,8 +33,8 @@ Author
\*---------------------------------------------------------------------------*/
#ifndef regionCoupleFvsPatchFieldsFwd_H
#define regionCoupleFvsPatchFieldsFwd_H
#ifndef regionCouplingFvsPatchFieldsFwd_H
#define regionCouplingFvsPatchFieldsFwd_H
#include "fvPatchField.H"
#include "fieldTypes.H"
@ -46,9 +46,9 @@ namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
template<class Type> class regionCoupleFvsPatchField;
template<class Type> class regionCouplingFvsPatchField;
makeFvsPatchTypeFieldTypedefs(regionCouple)
makeFvsPatchTypeFieldTypedefs(regionCoupling)
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

17
src/finiteVolume/fvMesh/fvPatches/constraint/ggi/ggiFvPatch.C
View file

@ -135,6 +135,7 @@ void Foam::ggiFvPatch::makeDeltaCoeffs(scalarField& dc) const
}
// Make patch face non-orthogonality correction vectors
void Foam::ggiFvPatch::makeCorrVecs(vectorField& cv) const
{
// Non-orthogonality correction on a ggi interface
@ -149,14 +150,6 @@ void Foam::ggiFvPatch::makeCorrVecs(vectorField& cv) const
}
const Foam::ggiFvPatch& Foam::ggiFvPatch::shadow() const
{
const fvPatch& p = this->boundaryMesh()[ggiPolyPatch_.shadowIndex()];
return refCast<const ggiFvPatch>(p);
}
// Return delta (P to N) vectors across coupled patch
Foam::tmp<Foam::vectorField> Foam::ggiFvPatch::delta() const
{
@ -192,6 +185,14 @@ Foam::tmp<Foam::vectorField> Foam::ggiFvPatch::delta() const
}
const Foam::ggiFvPatch& Foam::ggiFvPatch::shadow() const
{
const fvPatch& p = this->boundaryMesh()[ggiPolyPatch_.shadowIndex()];
return refCast<const ggiFvPatch>(p);
}
bool Foam::ggiFvPatch::master() const
{
return ggiPolyPatch_.master();

208
src/finiteVolume/fvMesh/fvPatches/constraint/regionCouple/regionCoupleFvPatch.C
View file

@ -45,24 +45,57 @@ namespace Foam
}
// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
Foam::regionCoupleFvPatch::~regionCoupleFvPatch()
{}
// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
// Make patch weighting factors
void Foam::regionCoupleFvPatch::makeWeights(scalarField& w) const
{
if (rcPolyPatch_.attached())
if (rcPolyPatch_.coupled())
{
vectorField n = nf();
if (rcPolyPatch_.master())
{
vectorField n = nf();
// Note: mag in the dot-product.
// For all valid meshes, the non-orthogonality will be less that
// 90 deg and the dot-product will be positive. For invalid
// meshes (d & s <= 0), this will stabilise the calculation
// but the result will be poor. HJ, 24/Aug/2011
scalarField nfc =
mag(n & (rcPolyPatch_.reconFaceCellCentres() - Cf()));
// Note: mag in the dot-product.
// For all valid meshes, the non-orthogonality will be less that
// 90 deg and the dot-product will be positive. For invalid
// meshes (d & s <= 0), this will stabilise the calculation
// but the result will be poor. HJ, 24/Aug/2011
scalarField nfc =
mag(n & (rcPolyPatch_.reconFaceCellCentres() - Cf()));
w = nfc/(mag(n & (Cf() - Cn())) + nfc);
w = nfc/(mag(n & (Cf() - Cn())) + nfc);
if (bridgeOverlap())
{
// Set overlap weights to 0.5 and use mirrored neighbour field
// for interpolation. HJ, 21/Jan/2009
bridge(scalarField(size(), 0.5), w);
}
}
else
{
// Pick up weights from the master side
scalarField masterWeights(shadow().size());
shadow().makeWeights(masterWeights);
scalarField oneMinusW = 1 - masterWeights;
w = interpolate(oneMinusW);
if (bridgeOverlap())
{
// Set overlap weights to 0.5 and use mirrored neighbour field
// for interpolation. HJ, 21/Jan/2009
bridge(scalarField(size(), 0.5), w);
}
}
}
else
{
@ -74,12 +107,35 @@ void Foam::regionCoupleFvPatch::makeWeights(scalarField& w) const
// Make patch face - neighbour cell distances
void Foam::regionCoupleFvPatch::makeDeltaCoeffs(scalarField& dc) const
{
if (rcPolyPatch_.attached())
if (rcPolyPatch_.coupled())
{
// Stabilised form for bad meshes. HJ, 24/Aug/2011
vectorField d = delta();
if (rcPolyPatch_.master())
{
// Stabilised form for bad meshes. HJ, 24/Aug/2011
vectorField d = delta();
dc = 1.0/max(nf() & d, 0.05*mag(d));
dc = 1.0/max(nf() & d, 0.05*mag(d));
if (bridgeOverlap())
{
scalarField bridgeDeltas = nf() & fvPatch::delta();
bridge(bridgeDeltas, dc);
}
}
else
{
scalarField masterDeltas(shadow().size());
shadow().makeDeltaCoeffs(masterDeltas);
dc = interpolate(masterDeltas);
if (bridgeOverlap())
{
scalarField bridgeDeltas = nf() & fvPatch::delta();
bridge(bridgeDeltas, dc);
}
}
}
else
{
@ -91,7 +147,7 @@ void Foam::regionCoupleFvPatch::makeDeltaCoeffs(scalarField& dc) const
// Make patch face non-orthogonality correction vectors
void Foam::regionCoupleFvPatch::makeCorrVecs(vectorField& cv) const
{
if (rcPolyPatch_.attached())
if (rcPolyPatch_.coupled())
{
// Non-orthogonality correction in attached state identical to ggi
// interface
@ -111,7 +167,48 @@ void Foam::regionCoupleFvPatch::makeCorrVecs(vectorField& cv) const
}
// * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
// Return delta (P to N) vectors across coupled patch
Foam::tmp<Foam::vectorField> Foam::regionCoupleFvPatch::delta() const
{
if (rcPolyPatch_.coupled())
{
if (rcPolyPatch_.master())
{
tmp<vectorField> tDelta =
rcPolyPatch_.reconFaceCellCentres() - Cn();
if (bridgeOverlap())
{
vectorField bridgeDeltas = Cf() - Cn();
bridge(bridgeDeltas, tDelta());
}
return tDelta;
}
else
{
tmp<vectorField> tDelta = interpolate
(
shadow().Cn() - rcPolyPatch_.shadow().reconFaceCellCentres()
);
if (bridgeOverlap())
{
vectorField bridgeDeltas = Cf() - Cn();
bridge(bridgeDeltas, tDelta());
}
return tDelta;
}
}
else
{
return fvPatch::delta();
}
}
bool Foam::regionCoupleFvPatch::coupled() const
{
@ -138,16 +235,74 @@ const Foam::regionCoupleFvPatch& Foam::regionCoupleFvPatch::shadow() const
}
// Return delta (P to N) vectors across coupled patch
Foam::tmp<Foam::vectorField> Foam::regionCoupleFvPatch::delta() const
bool Foam::regionCoupleFvPatch::master() const
{
if (rcPolyPatch_.attached())
return rcPolyPatch_.master();
}
bool Foam::regionCoupleFvPatch::fineLevel() const
{
return true;
}
Foam::label Foam::regionCoupleFvPatch::shadowIndex() const
{
return rcPolyPatch_.shadowIndex();
}
const Foam::ggiLduInterface&
Foam::regionCoupleFvPatch::shadowInterface() const
{
const fvPatch& p =
shadowRegion().boundary()[rcPolyPatch_.shadowIndex()];
return refCast<const ggiLduInterface>(p);
}
Foam::label Foam::regionCoupleFvPatch::zoneSize() const
{
return rcPolyPatch_.zone().size();
}
const Foam::labelList& Foam::regionCoupleFvPatch::zoneAddressing() const
{
return rcPolyPatch_.zoneAddressing();
}
const Foam::labelListList& Foam::regionCoupleFvPatch::addressing() const
{
if (rcPolyPatch_.master())
{
return rcPolyPatch_.reconFaceCellCentres() - Cn();
return rcPolyPatch_.patchToPatch().masterAddr();
}
else
{
return fvPatch::delta();
return rcPolyPatch_.patchToPatch().slaveAddr();
}
}
bool Foam::regionCoupleFvPatch::localParallel() const
{
return rcPolyPatch_.localParallel();
}
const Foam::scalarListList& Foam::regionCoupleFvPatch::weights() const
{
if (rcPolyPatch_.master())
{
return rcPolyPatch_.patchToPatch().masterWeights();
}
else
{
return rcPolyPatch_.patchToPatch().slaveWeights();
}
}
@ -167,7 +322,7 @@ void Foam::regionCoupleFvPatch::initTransfer
const unallocLabelList& interfaceData
) const
{
transferBuffer_ = interfaceData;
labelTransferBuffer_ = interfaceData;
}
@ -177,8 +332,8 @@ Foam::tmp<Foam::labelField> Foam::regionCoupleFvPatch::transfer
const unallocLabelList& interfaceData
) const
{
//HJ Should this be mapped? 22/Jun/2007
return shadow().transferBuffer();
return shadow().labelTransferBuffer();
}
@ -188,7 +343,7 @@ void Foam::regionCoupleFvPatch::initInternalFieldTransfer
const unallocLabelList& iF
) const
{
transferBuffer_ = patchInternalField(iF);
labelTransferBuffer_ = patchInternalField(iF);
}
@ -198,8 +353,7 @@ Foam::tmp<Foam::labelField> Foam::regionCoupleFvPatch::internalFieldTransfer
const unallocLabelList& iF
) const
{
//HJ Should this be mapped? 22/Jun/2007
return shadow().transferBuffer();
return shadow().labelTransferBuffer();
}

124
src/finiteVolume/fvMesh/fvPatches/constraint/regionCouple/regionCoupleFvPatch.H
View file

@ -41,7 +41,7 @@ SourceFiles
#define regionCoupleFvPatch_H
#include "coupledFvPatch.H"
#include "regionCoupleLduInterface.H"
#include "ggiLduInterface.H"
#include "regionCouplePolyPatch.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@ -58,16 +58,13 @@ class fvMesh;
class regionCoupleFvPatch
:
public coupledFvPatch,
public regionCoupleLduInterface
public ggiLduInterface
{
// Private Data
//- Reference to polyPatch
const regionCouplePolyPatch& rcPolyPatch_;
//- Transfer buffer
mutable labelField transferBuffer_;
protected:
@ -82,12 +79,6 @@ protected:
//- Make patch face non-orthogonality correction vectors
virtual void makeCorrVecs(vectorField&) const;
//- Return contents of transfer buffer
const labelField& transferBuffer() const
{
return transferBuffer_;
}
public:
@ -101,53 +92,108 @@ public:
regionCoupleFvPatch(const polyPatch& patch, const fvBoundaryMesh& bm)
:
coupledFvPatch(patch, bm),
rcPolyPatch_(refCast<const regionCouplePolyPatch>(patch)),
transferBuffer_()
rcPolyPatch_(refCast<const regionCouplePolyPatch>(patch))
{}
// Destructor
virtual ~regionCoupleFvPatch()
{}
virtual ~regionCoupleFvPatch();
// Member functions
//- Return true if coupled
virtual bool coupled() const;
// Access
//- Return shadow patch index
int shadowIndex() const
{
return rcPolyPatch_.shadowIndex();
}
//- Return true if coupled
virtual bool coupled() const;
//- Return shadow region
const fvMesh& shadowRegion() const;
//- Return shadow patch
const fvMesh& shadowRegion() const;
//- Return shadow patch
const regionCoupleFvPatch& shadow() const;
//- Return shadow patch
const regionCoupleFvPatch& shadow() const;
//- Interpolate face field
template<class Type>
tmp<Field<Type> > interpolate(const Field<Type>& pf) const
{
return rcPolyPatch_.interpolate(pf);
}
//- Use bridging to fix overlap error in interpolation
bool bridgeOverlap() const
{
return rcPolyPatch_.bridgeOverlap();
}
template<class Type>
tmp<Field<Type> > interpolate(const tmp<Field<Type> >& tpf) const
{
return rcPolyPatch_.interpolate(tpf);
}
//- Return delta (P to N) vectors across coupled patch
virtual tmp<vectorField> delta() const;
//- Return delta (P to N) vectors across coupled patch
virtual tmp<vectorField> delta() const;
// Interpolation
//- Interpolate face field
template<class Type>
tmp<Field<Type> > interpolate(const Field<Type>& pf) const
{
return rcPolyPatch_.interpolate(pf);
}
template<class Type>
tmp<Field<Type> > interpolate(const tmp<Field<Type> >& tpf) const
{
return rcPolyPatch_.interpolate(tpf);
}
//- Bridge interpolated face field for uncovered faces
template<class Type>
void bridge
(
const Field<Type>& bridgeField,
Field<Type>& ff
) const
{
return rcPolyPatch_.bridge(bridgeField, ff);
}
// Interface transfer functions
//- Is this the master side?
virtual bool master() const;
//- Is this the fine level?
virtual bool fineLevel() const;
//- Return shadow patch index
virtual label shadowIndex() const;
//- Return shadow interface
virtual const ggiLduInterface& shadowInterface() const;
//- Return zone size
virtual label zoneSize() const;
//- Return zone addressing
virtual const labelList& zoneAddressing() const;
//- Return addressing. Master side returns own addressing and
// slave side returns addressing from master
virtual const labelListList& addressing() const;
//- Is the patch localised on a single processor
virtual bool localParallel() const;
//- Return weights. Master side returns own weights and
// slave side returns weights from master
virtual const scalarListList& weights() const;
//- Return face transformation tensor
virtual const tensorField& forwardT() const
{
return coupledFvPatch::forwardT();
}
//- Return neighbour-cell transformation tensor
virtual const tensorField& reverseT() const
{
return coupledFvPatch::reverseT();
}
//- Return the values of the given internal data adjacent to
// the interface as a field
virtual tmp<labelField> interfaceInternalField

71
src/finiteVolume/interpolation/surfaceInterpolation/schemes/harmonic/harmonic.H
View file

@ -42,6 +42,7 @@ SourceFiles
#include "surfaceInterpolationScheme.H"
#include "volFields.H"
#include "surfaceFields.H"
#include "magLongDelta.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@ -124,6 +125,12 @@ public:
)
);
const surfaceScalarField& deltaCoeffs = this->mesh().deltaCoeffs();
const surfaceScalarField& weights = this->mesh().weights();
const magLongDelta& mld = magLongDelta::New(this->mesh());
const surfaceScalarField& longDelta = mld.magDelta();
surfaceScalarField& w = tw();
const unallocLabelList& owner = this->mesh().owner();
@ -133,18 +140,30 @@ public:
scalarField& wIn = w.internalField();
// Larger small for complex arithmetic accuracy
const scalar kSmall = 1000*SMALL;
// Calculate internal weights using field magnitude
forAll (owner, faceI)
{
scalar mOwn = magPhi[owner[faceI]];
scalar mNei = magPhi[neighbour[faceI]];
label own = owner[faceI];
label nei = neighbour[faceI];
scalar den = mOwn - mNei;
scalar mOwn = magPhi[own]/(1 - weights[faceI]);
scalar mNei = magPhi[nei]/weights[faceI];
if (mag(den) > SMALL)
scalar den = magPhi[nei] - magPhi[own];
scalar mean = mOwn*mNei/
((mOwn + mNei)*longDelta[faceI]*deltaCoeffs[faceI]);
// Note: complex arithmetic requires extra accuracy
// This is a division of two close subtractions
// HJ, 28/Sep/2011
if (mag(den) > kSmall)
{
scalar mean = 2.0*mOwn*mNei/(mOwn + mNei);
wIn[faceI] = (mean - mNei)/den;
// mOwn > mNei
wIn[faceI] = (magPhi[nei] - mean)/den;
}
else
{
@ -156,29 +175,45 @@ public:
{
fvsPatchScalarField& wp = w.boundaryField()[pi];
const fvPatchField<Type>& patchPhi = phi.boundaryField()[pi];
const fvPatchField<Type>& pPhi = phi.boundaryField()[pi];
if (patchPhi.coupled())
if (pPhi.coupled())
{
scalarField magPhiOwn = mag(patchPhi.patchInternalField());
scalarField magPhiNei = mag(patchPhi.patchNeighbourField());
const scalarField& pWeights = weights.boundaryField()[pi];
scalarField magPhiOwn = mag(pPhi.patchInternalField());
scalarField magPhiNei = mag(pPhi.patchNeighbourField());
const scalarField& pDeltaCoeffs =
deltaCoeffs.boundaryField()[pi];
const scalarField& pLongDelta = mld.magDelta(pi);
// Calculate internal weights using field magnitude
forAll (patchPhi, faceI)
forAll (pPhi, faceI)
{
scalar mOwn = magPhiOwn[faceI];
scalar mNei = magPhiNei[faceI];
scalar mOwn = magPhiOwn[faceI]/(1 - pWeights[faceI]);
scalar mNei = magPhiNei[faceI]/pWeights[faceI];
scalar den = mOwn - mNei;
scalar den = magPhiNei[faceI] - magPhiOwn[faceI];
if (mag(den) > SMALL)
// Note: complex arithmetic requires extra accuracy
// This is a division of two close subtractions
// HJ, 28/Sep/2011
if (mag(den) > kSmall)
{
scalar mean = 2.0*mOwn*mNei/(mOwn + mNei);
wp[faceI] = (mean - mNei)/den;
scalar mean = mOwn*mNei/
(
(mOwn + mNei)*
pLongDelta[faceI]*
pDeltaCoeffs[faceI]
);
wp[faceI] = (magPhiNei[faceI] - mean)/den;
}
else
{
wp[faceI] = 0.5;
wp[faceI] = 0.5;
}
}
}

246
src/finiteVolume/interpolation/surfaceInterpolation/schemes/harmonic/magLongDelta.C Normal file
View file

@ -0,0 +1,246 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright held by original author
\\/ M anipulation |
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
OpenFOAM is free software; you can redistribute it and/or modify it
under the terms of the GNU General Public License as published by the
Free Software Foundation; either version 2 of the License, or (at your
option) any later version.
OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
for more details.
You should have received a copy of the GNU General Public License
along with OpenFOAM; if not, write to the Free Software Foundation,
Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
\*---------------------------------------------------------------------------*/
#include "magLongDelta.H"
#include "surfaceFields.H"
#include "volFields.H"
#include "mapPolyMesh.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
namespace Foam
{
defineTypeNameAndDebug(magLongDelta, 0);
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
Foam::magLongDelta::magLongDelta(const fvMesh& mesh)
:
MeshObject<fvMesh, magLongDelta>(mesh),
magLongDeltaPtr_(NULL),
magLongDeltaBnd_(mesh.boundary().size(), NULL)
{}
// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
Foam::magLongDelta::~magLongDelta()
{
clearData();
}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
void Foam::magLongDelta::clearData() const
{
if (magLongDeltaPtr_)
{
delete magLongDeltaPtr_;
}
else
{
forAll (magLongDeltaBnd_, i)
{
if (magLongDeltaBnd_[i])
{
delete magLongDeltaBnd_[i];
}
}
}
}
void Foam::magLongDelta::makeMagLongDistance() const
{
if (debug)
{
Info<< "magLongDelta::makeMagLongDistance() :"
<< "Constructing magnitude of long cell distance"
<< endl;
}
magLongDeltaPtr_ = new surfaceScalarField
(
IOobject
(
"magLongDelta",
mesh().pointsInstance(),
mesh(),
IOobject::NO_READ,
IOobject::NO_WRITE,
false
),
mesh(),
dimensionedScalar("zero", dimLength, 0.0)
);
surfaceScalarField& mldp = *magLongDeltaPtr_;
// Set local references to mesh data
const unallocLabelList& owner = mesh().owner();
const unallocLabelList& neighbour = mesh().neighbour();
const vectorField& Cf = mesh().faceCentres();
const vectorField& C = mesh().cellCentres();
const vectorField& Sf = mesh().faceAreas();
const scalarField& magSf = mesh().magSf();
forAll (owner, facei)
{
// This must be the same as in surfaceInterpolation.C
scalar SfdOwn = mag(Sf[facei] & (Cf[facei] - C[owner[facei]]));
scalar SfdNei = mag(Sf[facei] & (C[neighbour[facei]] - Cf[facei]));
mldp[facei] = (SfdOwn + SfdNei)/magSf[facei];
}
forAll (mldp.boundaryField(), patchi)
{
const fvPatch& p = mesh().boundary()[patchi];
if (p.coupled())
{
if (!magLongDeltaBnd_[patchi])
{
makeMagLongDistance(patchi);
}
mldp.boundaryField()[patchi] = *magLongDeltaBnd_[patchi];
delete magLongDeltaBnd_[patchi];
magLongDeltaBnd_[patchi] = &mldp.boundaryField()[patchi];
}
else
{
delete magLongDeltaBnd_[patchi];
}
}
if (debug)
{
Info<< "magLongDelta::makeMagLongDistance() :"
<< "Finished magnitude of long cell distance"
<< endl;
}
}
void Foam::magLongDelta::makeMagLongDistance(label patchi) const
{
if (debug)
{
Info<< "magLongDelta::makeMagLongDistanceBnd(label patchi) :"
<< "Constructing magnitude of long cell distance"
<< endl;
}
const fvPatch& p = mesh().boundary()[patchi];
vectorField d = p.fvPatch::delta();
magLongDeltaBnd_[patchi] =
new scalarField
(
(mag(p.Sf() & d) + mag(p.Sf() & (p.delta() - d)))/p.magSf()
);
if (debug)
{
Info<< "magLongDelta::makeMagLongDistanceBnd(label patchi) :"
<< "Finished magnitude of long cell distance"
<< endl;
}
}
const Foam::surfaceScalarField& Foam::magLongDelta::magDelta() const
{
if (!magLongDeltaPtr_)
{
makeMagLongDistance();
}
return *magLongDeltaPtr_;
}
const Foam::scalarField& Foam::magLongDelta::magDelta
(
const label patchi
) const
{
if (!mesh().boundary()[patchi].coupled())
{
FatalErrorIn
(
"const Foam::scalarField& Foam::magLongDelta::magDelta("
"const label patchi) const"
) << "patch is not a coupled. Cannot calculate long distance"
<< abort(FatalError);
}
if (!magLongDeltaBnd_[patchi])
{
makeMagLongDistance(patchi);
}
return *magLongDeltaBnd_[patchi];
}
bool Foam::magLongDelta::movePoints() const
{
if (debug)
{
InfoIn("bool magLongDelta::movePoints() const")
<< "Clearing long cell distance data" << endl;
}
clearData();
magLongDeltaBnd_.setSize(mesh().boundary().size(), NULL);
return true;
}
bool Foam::magLongDelta::updateMesh(const mapPolyMesh& mpm) const
{
if (debug)
{
InfoIn("bool magLongDelta::updateMesh(const mapPolyMesh&) const")
<< "Clearing long cell distance data" << endl;
}
clearData();
magLongDeltaBnd_.setSize(mesh().boundary().size(), NULL);
return true;
}
// ************************************************************************* //

124
src/finiteVolume/interpolation/surfaceInterpolation/schemes/harmonic/magLongDelta.H Normal file
View file

@ -0,0 +1,124 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright held by original author
\\/ M anipulation |
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
OpenFOAM is free software; you can redistribute it and/or modify it
under the terms of the GNU General Public License as published by the
Free Software Foundation; either version 2 of the License, or (at your
option) any later version.
OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
for more details.
You should have received a copy of the GNU General Public License
along with OpenFOAM; if not, write to the Free Software Foundation,
Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
Class
Foam::magLongDelta
Description
Long distance to the neighbour via face centre
SourceFiles
magLongDelta.C
\*---------------------------------------------------------------------------*/
#ifndef magLongDelta_H
#define magLongDelta_H
#include "MeshObject.H"
#include "fvMesh.H"
#include "primitiveFieldsFwd.H"
#include "surfaceFieldsFwd.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
class mapPolyMesh;
/*---------------------------------------------------------------------------*\
Class magLongDelta Declaration
\*---------------------------------------------------------------------------*/
class magLongDelta
:
public MeshObject<fvMesh, magLongDelta>
{
// Private data
//- Magnitude of the normal distances to the face
mutable surfaceScalarField* magLongDeltaPtr_;
//- Magnitude of the normal distances on coupled patches
//- Data in magLongDeltaPtr_ is reused, if it exists
mutable List<scalarField*> magLongDeltaBnd_;
// Private member functions
//- Construct Long distance
void makeMagLongDistance() const;
//- Construct Long distance
void makeMagLongDistance(label patchi) const;
//- Clear Data
void clearData() const;
public:
// Declare name of the class and its debug switch
TypeName("magLongDelta");
// Constructors
//- Construct given an fvMesh
explicit magLongDelta(const fvMesh&);
// Destructor
virtual ~magLongDelta();
// Member functions
//- Return reference to the magnitude of the long distance
const surfaceScalarField& magDelta() const;
//- Return reference to the magnitude of the long distance
//- for patch i
const scalarField& magDelta(const label patchi) const;
//- Update after mesh motion:
// Delete data when the mesh moves
virtual bool movePoints() const;
//- Update after topo change:
// Delete data when mesh changes
virtual bool updateMesh(const mapPolyMesh&) const;
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //