Forward compatibility: flex
[foam-extend-3.2.git] / src / finiteVolume / fvMesh / fvPatches / constraint / regionCouple / regionCoupleFvPatch.C
blob28aaa6d8671efbb9acfa4c544ea1fae8ec408692
1 /*---------------------------------------------------------------------------*\
2   =========                 |
3   \\      /  F ield         | foam-extend: Open Source CFD
4    \\    /   O peration     | Version:     3.2
5     \\  /    A nd           | Web:         http://www.foam-extend.org
6      \\/     M anipulation  | For copyright notice see file Copyright
7 -------------------------------------------------------------------------------
8 License
9     This file is part of foam-extend.
11     foam-extend is free software: you can redistribute it and/or modify it
12     under the terms of the GNU General Public License as published by the
13     Free Software Foundation, either version 3 of the License, or (at your
14     option) any later version.
16     foam-extend is distributed in the hope that it will be useful, but
17     WITHOUT ANY WARRANTY; without even the implied warranty of
18     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
19     General Public License for more details.
21     You should have received a copy of the GNU General Public License
22     along with foam-extend.  If not, see <http://www.gnu.org/licenses/>.
24 Description
25     Region couple patch
27 Author
28     Hrvoje Jasak, Wikki Ltd.  All rights reserved
30 \*---------------------------------------------------------------------------*/
32 #include "regionCoupleFvPatch.H"
33 #include "addToRunTimeSelectionTable.H"
34 #include "fvMesh.H"
35 #include "fvBoundaryMesh.H"
36 #include "foamTime.H"
38 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
40 namespace Foam
42     defineTypeNameAndDebug(regionCoupleFvPatch, 0);
43     addToRunTimeSelectionTable(fvPatch, regionCoupleFvPatch, polyPatch);
47 // * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
49 Foam::regionCoupleFvPatch::~regionCoupleFvPatch()
53 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
55 // Make patch weighting factors
56 void Foam::regionCoupleFvPatch::makeWeights(scalarField& w) const
58     if (rcPolyPatch_.coupled())
59     {
60         if (rcPolyPatch_.master())
61         {
62             vectorField n = nf();
64             // Note: mag in the dot-product.
65             // For all valid meshes, the non-orthogonality will be less than
66             // 90 deg and the dot-product will be positive.  For invalid
67             // meshes (d & s <= 0), this will stabilise the calculation
68             // but the result will be poor.  HJ, 24/Aug/2011
69             scalarField nfc =
70                 mag(n & (rcPolyPatch_.reconFaceCellCentres() - Cf()));
72             w = nfc/(mag(n & (Cf() - Cn())) + nfc);
74             if (bridgeOverlap())
75             {
76                 // Set overlap weights to 0.5 and use mirrored neighbour field
77                 // for interpolation.  HJ, 21/Jan/2009
78                 bridge(scalarField(size(), 0.5), w);
79             }
80         }
81         else
82         {
83             // Pick up weights from the master side
84             scalarField masterWeights(shadow().size());
85             shadow().makeWeights(masterWeights);
87             scalarField oneMinusW = 1 - masterWeights;
89             w = interpolate(oneMinusW);
91             if (bridgeOverlap())
92             {
93                 // Set overlap weights to 0.5 and use mirrored neighbour field
94                 // for interpolation.  HJ, 21/Jan/2009
95                 bridge(scalarField(size(), 0.5), w);
96             }
97         }
98     }
99     else
100     {
101         fvPatch::makeWeights(w);
102     }
106 // Make patch face - neighbour cell distances
107 void Foam::regionCoupleFvPatch::makeDeltaCoeffs(scalarField& dc) const
109     if (rcPolyPatch_.coupled())
110     {
111         if (rcPolyPatch_.master())
112         {
113             // Stabilised form for bad meshes.  HJ, 24/Aug/2011
114             vectorField d = delta();
116             dc = 1.0/max(nf() & d, 0.05*mag(d));
118             if (bridgeOverlap())
119             {
120                 scalarField bridgeDeltas = nf() & fvPatch::delta();
122                 bridge(bridgeDeltas, dc);
123             }
124         }
125         else
126         {
127             scalarField masterDeltas(shadow().size());
128             shadow().makeDeltaCoeffs(masterDeltas);
129             dc = interpolate(masterDeltas);
131             if (bridgeOverlap())
132             {
133                 scalarField bridgeDeltas = nf() & fvPatch::delta();
135                 bridge(bridgeDeltas, dc);
136             }
137         }
138     }
139     else
140     {
141         fvPatch::makeDeltaCoeffs(dc);
142     }
146 // Make patch face non-orthogonality correction vectors
147 void Foam::regionCoupleFvPatch::makeCorrVecs(vectorField& cv) const
149     if (rcPolyPatch_.coupled())
150     {
151         // Non-orthogonality correction in attached state identical to ggi
152         // interface
154         // Calculate correction vectors on coupled patches
155         const scalarField& patchDeltaCoeffs = deltaCoeffs();
157         vectorField patchDeltas = delta();
158         vectorField n = nf();
160         // If non-orthogonality is over 90 deg, kill correction vector
161         // HJ, 6/Jan/2011
162         cv = pos(patchDeltas & n)*(n - patchDeltas*patchDeltaCoeffs);
163     }
164     else
165     {
166         // No correction in detached state.  HJ, 26/Jul/2011
167         cv = vector::zero;
168     }
172 // Return delta (P to N) vectors across coupled patch
173 Foam::tmp<Foam::vectorField> Foam::regionCoupleFvPatch::delta() const
175     if (rcPolyPatch_.coupled())
176     {
177         if (rcPolyPatch_.master())
178         {
179             tmp<vectorField> tDelta =
180                 rcPolyPatch_.reconFaceCellCentres() - Cn();
182             if (bridgeOverlap())
183             {
184                 vectorField bridgeDeltas = Cf() - Cn();
186                 bridge(bridgeDeltas, tDelta());
187             }
189             return tDelta;
190         }
191         else
192         {
193             tmp<vectorField> tDelta = interpolate
194             (
195                 shadow().Cn() - rcPolyPatch_.shadow().reconFaceCellCentres()
196             );
198             if (bridgeOverlap())
199             {
200                 vectorField bridgeDeltas = Cf() - Cn();
202                 bridge(bridgeDeltas, tDelta());
203             }
205             return tDelta;
206         }
207     }
208     else
209     {
210         return fvPatch::delta();
211     }
215 bool Foam::regionCoupleFvPatch::coupled() const
217     return rcPolyPatch_.coupled();
221 const Foam::fvMesh& Foam::regionCoupleFvPatch::shadowRegion() const
223     return
224         boundaryMesh().mesh().objectRegistry::parent().lookupObject<fvMesh>
225         (
226             rcPolyPatch_.shadowRegionName()
227         );
231 const Foam::regionCoupleFvPatch& Foam::regionCoupleFvPatch::shadow() const
233     const fvPatch& p =
234         shadowRegion().boundary()[rcPolyPatch_.shadowIndex()];
236     return refCast<const regionCoupleFvPatch>(p);
240 bool Foam::regionCoupleFvPatch::master() const
242     return rcPolyPatch_.master();
246 bool Foam::regionCoupleFvPatch::fineLevel() const
248     return true;
252 Foam::label Foam::regionCoupleFvPatch::shadowIndex() const
254     return rcPolyPatch_.shadowIndex();
258 const Foam::ggiLduInterface&
259 Foam::regionCoupleFvPatch::shadowInterface() const
261     const fvPatch& p =
262         shadowRegion().boundary()[rcPolyPatch_.shadowIndex()];
264     return refCast<const ggiLduInterface>(p);
268 Foam::label Foam::regionCoupleFvPatch::zoneSize() const
270     return rcPolyPatch_.zone().size();
274 const Foam::labelList& Foam::regionCoupleFvPatch::zoneAddressing() const
276     return rcPolyPatch_.zoneAddressing();
280 const Foam::labelListList& Foam::regionCoupleFvPatch::addressing() const
282     if (rcPolyPatch_.master())
283     {
284         return rcPolyPatch_.patchToPatch().masterAddr();
285     }
286     else
287     {
288         return rcPolyPatch_.patchToPatch().slaveAddr();
289     }
293 bool Foam::regionCoupleFvPatch::localParallel() const
295     return rcPolyPatch_.localParallel();
299 const Foam::scalarListList& Foam::regionCoupleFvPatch::weights() const
301     if (rcPolyPatch_.master())
302     {
303         return rcPolyPatch_.patchToPatch().masterWeights();
304     }
305     else
306     {
307         return rcPolyPatch_.patchToPatch().slaveWeights();
308     }
312 Foam::tmp<Foam::labelField> Foam::regionCoupleFvPatch::interfaceInternalField
314     const unallocLabelList& internalData
315 ) const
317     return patchInternalField(internalData);
321 void Foam::regionCoupleFvPatch::initTransfer
323     const Pstream::commsTypes commsType,
324     const unallocLabelList& interfaceData
325 ) const
327     labelTransferBuffer_ = interfaceData;
331 Foam::tmp<Foam::labelField> Foam::regionCoupleFvPatch::transfer
333     const Pstream::commsTypes commsType,
334     const unallocLabelList& interfaceData
335 ) const
338     return shadow().labelTransferBuffer();
342 void Foam::regionCoupleFvPatch::initInternalFieldTransfer
344     const Pstream::commsTypes commsType,
345     const unallocLabelList& iF
346 ) const
348     labelTransferBuffer_ = patchInternalField(iF);
352 Foam::tmp<Foam::labelField> Foam::regionCoupleFvPatch::internalFieldTransfer
354     const Pstream::commsTypes commsType,
355     const unallocLabelList& iF
356 ) const
358     return shadow().labelTransferBuffer();
362 // ************************************************************************* //