FaceCellWave.C 24.5 KB
Newer Older
1
2
3
4
/*---------------------------------------------------------------------------*\
  =========                 |
  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   \\    /   O peration     |
5
    \\  /    A nd           | Copyright (C) 1991-2010 OpenCFD Ltd.
6
7
8
9
10
     \\/     M anipulation  |
-------------------------------------------------------------------------------
License
    This file is part of OpenFOAM.

11
12
13
14
    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 3 of the License, or
    (at your option) any later version.
15
16
17
18
19
20
21

    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
22
    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37

\*---------------------------------------------------------------------------*/

#include "FaceCellWave.H"
#include "polyMesh.H"
#include "processorPolyPatch.H"
#include "cyclicPolyPatch.H"
#include "OPstream.H"
#include "IPstream.H"
#include "PstreamReduceOps.H"
#include "debug.H"
#include "typeInfo.H"

// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //

38
39
40
template <class Type, class TrackingData>
const Foam::scalar
Foam::FaceCellWave<Type, TrackingData>::geomTol_ = 1e-6;
41

42
43
44
template <class Type, class TrackingData>
const Foam::scalar
Foam::FaceCellWave<Type, TrackingData>::propagationTol_ = 0.01;
45

46
47
template <class Type, class TrackingData>
Foam::label Foam::FaceCellWave<Type, TrackingData>::dummyTrackData_ = 12345;
48
49
50
51
52
53
54
55
56


// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //

// Update info for cellI, at position pt, with information from
// neighbouring face/cell.
// Updates:
//      - changedCell_, changedCells_, nChangedCells_,
//      - statistics: nEvals_, nUnvisitedCells_
57
58
template <class Type, class TrackingData>
bool Foam::FaceCellWave<Type, TrackingData>::updateCell
59
60
61
62
63
64
65
66
67
68
(
    const label cellI,
    const label neighbourFaceI,
    const Type& neighbourInfo,
    const scalar tol,
    Type& cellInfo
)
{
    nEvals_++;

69
    bool wasValid = cellInfo.valid(td_);
70
71
72
73
74
75
76
77

    bool propagate =
        cellInfo.updateCell
        (
            mesh_,
            cellI,
            neighbourFaceI,
            neighbourInfo,
78
79
            tol,
            td_
80
81
82
83
84
85
86
87
88
89
90
        );

    if (propagate)
    {
        if (!changedCell_[cellI])
        {
            changedCell_[cellI] = true;
            changedCells_[nChangedCells_++] = cellI;
        }
    }

91
    if (!wasValid && cellInfo.valid(td_))
92
93
94
95
96
97
98
99
100
101
102
103
104
    {
        --nUnvisitedCells_;
    }

    return propagate;
}


// Update info for faceI, at position pt, with information from
// neighbouring face/cell.
// Updates:
//      - changedFace_, changedFaces_, nChangedFaces_,
//      - statistics: nEvals_, nUnvisitedFaces_
105
106
template <class Type, class TrackingData>
bool Foam::FaceCellWave<Type, TrackingData>::updateFace
107
108
109
110
111
112
113
114
115
116
(
    const label faceI,
    const label neighbourCellI,
    const Type& neighbourInfo,
    const scalar tol,
    Type& faceInfo
)
{
    nEvals_++;

117
    bool wasValid = faceInfo.valid(td_);
118
119
120
121
122
123
124
125

    bool propagate =
        faceInfo.updateFace
        (
            mesh_,
            faceI,
            neighbourCellI,
            neighbourInfo,
126
127
            tol,
            td_
128
129
130
131
132
133
134
135
136
137
138
        );

    if (propagate)
    {
        if (!changedFace_[faceI])
        {
            changedFace_[faceI] = true;
            changedFaces_[nChangedFaces_++] = faceI;
        }
    }

139
    if (!wasValid && faceInfo.valid(td_))
140
141
142
143
144
145
146
147
148
149
150
151
152
    {
        --nUnvisitedFaces_;
    }

    return propagate;
}


// Update info for faceI, at position pt, with information from
// same face.
// Updates:
//      - changedFace_, changedFaces_, nChangedFaces_,
//      - statistics: nEvals_, nUnvisitedFaces_
153
154
template <class Type, class TrackingData>
bool Foam::FaceCellWave<Type, TrackingData>::updateFace
155
156
157
158
159
160
161
162
163
(
    const label faceI,
    const Type& neighbourInfo,
    const scalar tol,
    Type& faceInfo
)
{
    nEvals_++;

164
    bool wasValid = faceInfo.valid(td_);
165
166
167
168
169
170
171

    bool propagate =
        faceInfo.updateFace
        (
            mesh_,
            faceI,
            neighbourInfo,
172
173
            tol,
            td_
174
175
176
177
178
179
180
181
182
183
184
        );

    if (propagate)
    {
        if (!changedFace_[faceI])
        {
            changedFace_[faceI] = true;
            changedFaces_[nChangedFaces_++] = faceI;
        }
    }

185
    if (!wasValid && faceInfo.valid(td_))
186
187
188
189
190
191
192
193
194
    {
        --nUnvisitedFaces_;
    }

    return propagate;
}


// For debugging: check status on both sides of cyclic
195
196
197
198
199
template <class Type, class TrackingData>
void Foam::FaceCellWave<Type, TrackingData>::checkCyclic
(
    const polyPatch& patch
) const
200
{
201
202
    const cyclicPolyPatch& nbrPatch =
        refCast<const cyclicPolyPatch>(patch).neighbPatch();
203

204
    forAll(patch, patchFaceI)
205
206
    {
        label i1 = patch.start() + patchFaceI;
207
        label i2 = nbrPatch.start() + patchFaceI;
208

209
210
211
212
213
214
215
216
217
218
        if
        (
           !allFaceInfo_[i1].sameGeometry
            (
                mesh_,
                allFaceInfo_[i2],
                geomTol_,
                td_
            )
        )
219
        {
220
221
222
223
224
            FatalErrorIn
            (
                "FaceCellWave<Type, TrackingData>"
                "::checkCyclic(const polyPatch&)"
            )   << "problem: i:" << i1 << "  otheri:" << i2
225
226
227
228
229
230
231
                << "   faceInfo:" << allFaceInfo_[i1]
                << "   otherfaceInfo:" << allFaceInfo_[i2]
                << abort(FatalError);
        }

        if (changedFace_[i1] != changedFace_[i2])
        {
232
233
234
235
236
            FatalErrorIn
            (
                "FaceCellWave<Type, TrackingData>"
                "::checkCyclic(const polyPatch&)"
            )   << " problem: i:" << i1 << "  otheri:" << i2
237
238
239
240
241
242
243
244
245
246
                << "   faceInfo:" << allFaceInfo_[i1]
                << "   otherfaceInfo:" << allFaceInfo_[i2]
                << "   changedFace:" << changedFace_[i1]
                << "   otherchangedFace:" << changedFace_[i2]
                << abort(FatalError);
        }
    }
}


247
// Check if has cyclic patches
248
249
template <class Type, class TrackingData>
bool Foam::FaceCellWave<Type, TrackingData>::hasCyclicPatch() const
250
251
252
{
    forAll(mesh_.boundaryMesh(), patchI)
    {
253
        if (isA<cyclicPolyPatch>(mesh_.boundaryMesh()[patchI]))
254
255
256
257
258
259
260
261
262
        {
            return true;
        }
    }
    return false;
}


// Copy face information into member data
263
264
template <class Type, class TrackingData>
void Foam::FaceCellWave<Type, TrackingData>::setFaceInfo
265
266
267
268
269
270
271
272
273
(
    const labelList& changedFaces,
    const List<Type>& changedFacesInfo
)
{
    forAll(changedFaces, changedFaceI)
    {
        label faceI = changedFaces[changedFaceI];

274
        bool wasValid = allFaceInfo_[faceI].valid(td_);
275
276
277
278
279

        // Copy info for faceI
        allFaceInfo_[faceI] = changedFacesInfo[changedFaceI];

        // Maintain count of unset faces
280
        if (!wasValid && allFaceInfo_[faceI].valid(td_))
281
282
283
284
285
286
287
288
289
290
291
292
293
        {
            --nUnvisitedFaces_;
        }

        // Mark faceI as changed, both on list and on face itself.

        changedFace_[faceI] = true;
        changedFaces_[nChangedFaces_++] = faceI;
    }
}


// Merge face information into member data
294
295
template <class Type, class TrackingData>
void Foam::FaceCellWave<Type, TrackingData>::mergeFaceInfo
296
297
298
299
(
    const polyPatch& patch,
    const label nFaces,
    const labelList& changedFaces,
300
    const List<Type>& changedFacesInfo
301
302
)
{
303
    for (label changedFaceI = 0; changedFaceI < nFaces; changedFaceI++)
304
305
306
307
308
309
310
311
    {
        const Type& neighbourWallInfo = changedFacesInfo[changedFaceI];
        label patchFaceI = changedFaces[changedFaceI];

        label meshFaceI = patch.start() + patchFaceI;

        Type& currentWallInfo = allFaceInfo_[meshFaceI];

312
        if (!currentWallInfo.equal(neighbourWallInfo, td_))
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
        {
            updateFace
            (
                meshFaceI,
                neighbourWallInfo,
                propagationTol_,
                currentWallInfo
            );
        }
    }
}


// Construct compact patchFace change arrays for a (slice of a) single patch.
// changedPatchFaces in local patch numbering.
// Return length of arrays.
329
330
template <class Type, class TrackingData>
Foam::label Foam::FaceCellWave<Type, TrackingData>::getChangedPatchFaces
331
332
333
334
335
336
337
338
339
340
(
    const polyPatch& patch,
    const label startFaceI,
    const label nFaces,
    labelList& changedPatchFaces,
    List<Type>& changedPatchFacesInfo
) const
{
    label nChangedPatchFaces = 0;

341
    for (label i = 0; i < nFaces; i++)
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
    {
        label patchFaceI = i + startFaceI;

        label meshFaceI = patch.start() + patchFaceI;

        if (changedFace_[meshFaceI])
        {
            changedPatchFaces[nChangedPatchFaces] = patchFaceI;
            changedPatchFacesInfo[nChangedPatchFaces] = allFaceInfo_[meshFaceI];
            nChangedPatchFaces++;
        }
    }
    return nChangedPatchFaces;
}


// Handle leaving domain. Implementation referred to Type
359
360
template <class Type, class TrackingData>
void Foam::FaceCellWave<Type, TrackingData>::leaveDomain
361
362
363
364
365
366
367
368
(
    const polyPatch& patch,
    const label nFaces,
    const labelList& faceLabels,
    List<Type>& faceInfo
) const
{
    const vectorField& fc = mesh_.faceCentres();
369

370
    for (label i = 0; i < nFaces; i++)
371
372
373
374
    {
        label patchFaceI = faceLabels[i];

        label meshFaceI = patch.start() + patchFaceI;
375
        faceInfo[i].leaveDomain(mesh_, patch, patchFaceI, fc[meshFaceI], td_);
376
377
378
379
380
    }
}


// Handle entering domain. Implementation referred to Type
381
382
template <class Type, class TrackingData>
void Foam::FaceCellWave<Type, TrackingData>::enterDomain
383
384
385
386
387
388
389
390
(
    const polyPatch& patch,
    const label nFaces,
    const labelList& faceLabels,
    List<Type>& faceInfo
) const
{
    const vectorField& fc = mesh_.faceCentres();
391

392
    for (label i = 0; i < nFaces; i++)
393
394
395
396
    {
        label patchFaceI = faceLabels[i];

        label meshFaceI = patch.start() + patchFaceI;
397
        faceInfo[i].enterDomain(mesh_, patch, patchFaceI, fc[meshFaceI], td_);
398
399
400
401
402
    }
}


// Transform. Implementation referred to Type
403
404
template <class Type, class TrackingData>
void Foam::FaceCellWave<Type, TrackingData>::transform
405
406
407
408
409
410
411
412
413
414
(
    const tensorField& rotTensor,
    const label nFaces,
    List<Type>& faceInfo
)
{
    if (rotTensor.size() == 1)
    {
        const tensor& T = rotTensor[0];

415
        for (label faceI = 0; faceI < nFaces; faceI++)
416
        {
417
            faceInfo[faceI].transform(mesh_, T, td_);
418
419
420
421
        }
    }
    else
    {
422
        for (label faceI = 0; faceI < nFaces; faceI++)
423
        {
424
            faceInfo[faceI].transform(mesh_, rotTensor[faceI], td_);
425
426
427
428
429
430
        }
    }
}


// Offset mesh face. Used for transferring from one cyclic half to the other.
431
432
template <class Type, class TrackingData>
void Foam::FaceCellWave<Type, TrackingData>::offset
433
434
435
436
437
438
439
(
    const polyPatch&,
    const label cycOffset,
    const label nFaces,
    labelList& faces
)
{
440
    for (label faceI = 0; faceI < nFaces; faceI++)
441
442
443
444
445
446
447
    {
        faces[faceI] += cycOffset;
    }
}


// Tranfer all the information to/from neighbouring processors
448
449
template <class Type, class TrackingData>
void Foam::FaceCellWave<Type, TrackingData>::handleProcPatches()
450
451
452
{
    // Send all

453
454
    PstreamBuffers pBufs(Pstream::nonBlocking);

455
456
457
458
    forAll(mesh_.boundaryMesh(), patchI)
    {
        const polyPatch& patch = mesh_.boundaryMesh()[patchI];

mattijs's avatar
mattijs committed
459
        if (isA<processorPolyPatch>(patch))
460
461
462
463
464
465
466
467
468
469
470
471
472
473
474
475
476
477
478
479
480
481
482
483
484
485
486
487
488
489
490
491
492
493
494
495
        {
            // Allocate buffers
            label nSendFaces;
            labelList sendFaces(patch.size());
            List<Type> sendFacesInfo(patch.size());

            // Determine which faces changed on current patch
            nSendFaces = getChangedPatchFaces
            (
                patch,
                0,
                patch.size(),
                sendFaces,
                sendFacesInfo
            );

            // Adapt wallInfo for leaving domain
            leaveDomain
            (
                patch,
                nSendFaces,
                sendFaces,
                sendFacesInfo
            );

            const processorPolyPatch& procPatch =
                refCast<const processorPolyPatch>(patch);

            if (debug)
            {
                Pout<< " Processor patch " << patchI << ' ' << patch.name()
                    << " communicating with " << procPatch.neighbProcNo()
                    << "  Sending:" << nSendFaces
                    << endl;
            }

496
            UOPstream toNeighbour(procPatch.neighbProcNo(), pBufs);
497
498
499
500
501
            //writeFaces(nSendFaces, sendFaces, sendFacesInfo, toNeighbour);
            toNeighbour
                << SubList<label>(sendFaces, nSendFaces)
                << SubList<Type>(sendFacesInfo, nSendFaces);

502
503
504
        }
    }

505
506
    pBufs.finishedSends();

507
508
509
510
511
512
    // Receive all

    forAll(mesh_.boundaryMesh(), patchI)
    {
        const polyPatch& patch = mesh_.boundaryMesh()[patchI];

mattijs's avatar
mattijs committed
513
        if (isA<processorPolyPatch>(patch))
514
515
516
517
518
519
520
521
        {
            const processorPolyPatch& procPatch =
                refCast<const processorPolyPatch>(patch);

            // Allocate buffers
            labelList receiveFaces(patch.size());
            List<Type> receiveFacesInfo(patch.size());

522
523
            {
                UIPstream fromNeighbour(procPatch.neighbProcNo(), pBufs);
524
                fromNeighbour >> receiveFaces >> receiveFacesInfo;
525
            }
526
527
528
529
530

            if (debug)
            {
                Pout<< " Processor patch " << patchI << ' ' << patch.name()
                    << " communicating with " << procPatch.neighbProcNo()
531
                    << "  Receiving:" << receiveFaces.size()
532
533
534
535
536
537
538
539
540
                    << endl;
            }

            // Apply transform to received data for non-parallel planes
            if (!procPatch.parallel())
            {
                transform
                (
                    procPatch.reverseT(),
541
                    receiveFaces.size(),
542
543
544
545
546
547
548
549
                    receiveFacesInfo
                );
            }

            // Adapt wallInfo for entering domain
            enterDomain
            (
                patch,
550
                receiveFaces.size(),
551
552
553
554
555
556
557
558
                receiveFaces,
                receiveFacesInfo
            );

            // Merge received info
            mergeFaceInfo
            (
                patch,
559
                receiveFaces.size(),
560
                receiveFaces,
561
                receiveFacesInfo
562
563
564
565
566
567
568
            );
        }
    }
}


// Transfer information across cyclic halves.
569
570
template <class Type, class TrackingData>
void Foam::FaceCellWave<Type, TrackingData>::handleCyclicPatches()
571
572
573
574
575
576
577
{
    forAll(mesh_.boundaryMesh(), patchI)
    {
        const polyPatch& patch = mesh_.boundaryMesh()[patchI];

        if (isA<cyclicPolyPatch>(patch))
        {
578
579
            const cyclicPolyPatch& nbrPatch =
                refCast<const cyclicPolyPatch>(patch).neighbPatch();
580
581
582

            // Allocate buffers
            label nReceiveFaces;
583
584
            labelList receiveFaces(patch.size());
            List<Type> receiveFacesInfo(patch.size());
585

586
            // Determine which faces changed
587
588
            nReceiveFaces = getChangedPatchFaces
            (
589
590
591
                nbrPatch,
                0,
                nbrPatch.size(),
592
593
594
595
                receiveFaces,
                receiveFacesInfo
            );

596
            // Adapt wallInfo for leaving domain
597
598
            leaveDomain
            (
599
                nbrPatch,
600
601
602
603
604
605
606
607
608
609
                nReceiveFaces,
                receiveFaces,
                receiveFacesInfo
            );

            const cyclicPolyPatch& cycPatch =
                refCast<const cyclicPolyPatch>(patch);

            if (!cycPatch.parallel())
            {
610
                // received data from other half
611
612
613
614
615
616
617
618
619
620
                transform
                (
                    cycPatch.forwardT(),
                    nReceiveFaces,
                    receiveFacesInfo
                );
            }

            if (debug)
            {
621
622
                Pout<< " Cyclic patch " << patchI << ' ' << cycPatch.name()
                    << "  Changed : " << nReceiveFaces
623
624
625
626
627
628
                    << endl;
            }

            // Half2: Adapt wallInfo for entering domain
            enterDomain
            (
629
                cycPatch,
630
631
632
633
634
635
636
637
                nReceiveFaces,
                receiveFaces,
                receiveFacesInfo
            );

            // Merge into global storage
            mergeFaceInfo
            (
638
                cycPatch,
639
640
                nReceiveFaces,
                receiveFaces,
641
                receiveFacesInfo
642
643
644
645
            );

            if (debug)
            {
646
                checkCyclic(cycPatch);
647
648
649
650
651
652
653
654
655
            }
        }
    }
}


// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //

// Set up only. Use setFaceInfo and iterate() to do actual calculation.
656
657
template <class Type, class TrackingData>
Foam::FaceCellWave<Type, TrackingData>::FaceCellWave
658
659
660
(
    const polyMesh& mesh,
    UList<Type>& allFaceInfo,
661
662
    UList<Type>& allCellInfo,
    TrackingData& td
663
664
665
666
667
)
:
    mesh_(mesh),
    allFaceInfo_(allFaceInfo),
    allCellInfo_(allCellInfo),
668
    td_(td),
669
670
671
672
673
674
    changedFace_(mesh_.nFaces(), false),
    changedFaces_(mesh_.nFaces()),
    nChangedFaces_(0),
    changedCell_(mesh_.nCells(), false),
    changedCells_(mesh_.nCells()),
    nChangedCells_(0),
675
    hasCyclicPatches_(hasCyclicPatch()),
676
677
678
679
680
681
682
    nEvals_(0),
    nUnvisitedCells_(mesh_.nCells()),
    nUnvisitedFaces_(mesh_.nFaces()),
    iter_(0)
{}


683
// Iterate, propagating changedFacesInfo across mesh, until no change (or
684
// maxIter reached). Initial cell values specified.
685
686
template <class Type, class TrackingData>
Foam::FaceCellWave<Type, TrackingData>::FaceCellWave
687
688
689
690
691
692
(
    const polyMesh& mesh,
    const labelList& changedFaces,
    const List<Type>& changedFacesInfo,
    UList<Type>& allFaceInfo,
    UList<Type>& allCellInfo,
693
694
    const label maxIter,
    TrackingData& td
695
696
697
698
699
)
:
    mesh_(mesh),
    allFaceInfo_(allFaceInfo),
    allCellInfo_(allCellInfo),
700
    td_(td),
701
702
703
704
705
706
    changedFace_(mesh_.nFaces(), false),
    changedFaces_(mesh_.nFaces()),
    nChangedFaces_(0),
    changedCell_(mesh_.nCells(), false),
    changedCells_(mesh_.nCells()),
    nChangedCells_(0),
707
    hasCyclicPatches_(hasCyclicPatch()),
708
709
710
711
712
713
714
715
716
717
718
719
720
721
722
    nEvals_(0),
    nUnvisitedCells_(mesh_.nCells()),
    nUnvisitedFaces_(mesh_.nFaces()),
    iter_(0)
{
    // Copy initial changed faces data
    setFaceInfo(changedFaces, changedFacesInfo);

    // Iterate until nothing changes
    iterate(maxIter);

    if ((maxIter > 0) && (iter_ >= maxIter))
    {
        FatalErrorIn
        (
723
            "FaceCellWave<Type, TrackingData>::FaceCellWave"
724
725
726
727
728
729
730
731
732
733
734
735
736
737
738
            "(const polyMesh&, const labelList&, const List<Type>,"
            " UList<Type>&, UList<Type>&, const label maxIter)"
        )
            << "Maximum number of iterations reached. Increase maxIter." << endl
            << "    maxIter:" << maxIter << endl
            << "    nChangedCells:" << nChangedCells_ << endl
            << "    nChangedFaces:" << nChangedFaces_ << endl
            << exit(FatalError);
    }
}


// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //


739
740
template <class Type, class TrackingData>
Foam::label Foam::FaceCellWave<Type, TrackingData>::getUnsetCells() const
741
742
743
744
745
{
    return nUnvisitedCells_;
}


746
747
template <class Type, class TrackingData>
Foam::label Foam::FaceCellWave<Type, TrackingData>::getUnsetFaces() const
748
749
750
751
752
753
754
{
    return nUnvisitedFaces_;
}



// Propagate cell to face
755
756
template <class Type, class TrackingData>
Foam::label Foam::FaceCellWave<Type, TrackingData>::faceToCell()
757
758
759
760
761
762
763
764
765
766
767
768
769
770
771
{
    const labelList& owner = mesh_.faceOwner();
    const labelList& neighbour = mesh_.faceNeighbour();
    label nInternalFaces = mesh_.nInternalFaces();

    for
    (
        label changedFaceI = 0;
        changedFaceI < nChangedFaces_;
        changedFaceI++
    )
    {
        label faceI = changedFaces_[changedFaceI];
        if (!changedFace_[faceI])
        {
772
            FatalErrorIn("FaceCellWave<Type, TrackingData>::faceToCell()")
773
774
775
776
777
778
779
780
781
782
783
784
785
786
                << "Face " << faceI
                << " not marked as having been changed"
                << abort(FatalError);
        }


        const Type& neighbourWallInfo = allFaceInfo_[faceI];

        // Evaluate all connected cells

        // Owner
        label cellI = owner[faceI];
        Type& currentWallInfo = allCellInfo_[cellI];

787
        if (!currentWallInfo.equal(neighbourWallInfo, td_))
788
789
790
791
792
793
794
795
796
797
798
799
800
801
802
803
804
        {
            updateCell
            (
                cellI,
                faceI,
                neighbourWallInfo,
                propagationTol_,
                currentWallInfo
            );
        }

        // Neighbour. Hack for check if face has neighbour.
        if (faceI < nInternalFaces)
        {
            cellI = neighbour[faceI];
            Type& currentWallInfo2 = allCellInfo_[cellI];

805
            if (!currentWallInfo2.equal(neighbourWallInfo, td_))
806
807
808
809
810
811
812
813
814
815
816
817
818
819
820
821
822
823
824
825
826
827
828
829
830
831
832
833
834
835
836
837
838
839
            {
                updateCell
                (
                    cellI,
                    faceI,
                    neighbourWallInfo,
                    propagationTol_,
                    currentWallInfo2
                );
            }
        }

        // Reset status of face
        changedFace_[faceI] = false;
    }

    // Handled all changed faces by now
    nChangedFaces_ = 0;

    if (debug)
    {
        Pout<< " Changed cells            : " << nChangedCells_ << endl;
    }

    // Sum nChangedCells over all procs
    label totNChanged = nChangedCells_;

    reduce(totNChanged, sumOp<label>());

    return totNChanged;
}


// Propagate cell to face
840
841
template <class Type, class TrackingData>
Foam::label Foam::FaceCellWave<Type, TrackingData>::cellToFace()
842
843
844
845
846
847
848
849
850
851
852
853
854
{
    const cellList& cells = mesh_.cells();

    for
    (
        label changedCellI = 0;
        changedCellI < nChangedCells_;
        changedCellI++
    )
    {
        label cellI = changedCells_[changedCellI];
        if (!changedCell_[cellI])
        {
mattijs's avatar
mattijs committed
855
856
            FatalErrorIn("FaceCellWave<Type, TrackingData>::cellToFace()")
                << "Cell " << cellI << " not marked as having been changed"
857
858
859
860
861
862
863
864
865
866
867
868
869
                << abort(FatalError);
        }

        const Type& neighbourWallInfo = allCellInfo_[cellI];

        // Evaluate all connected faces

        const labelList& faceLabels = cells[cellI];
        forAll(faceLabels, faceLabelI)
        {
            label faceI = faceLabels[faceLabelI];
            Type& currentWallInfo = allFaceInfo_[faceI];

870
            if (!currentWallInfo.equal(neighbourWallInfo, td_))
871
872
873
874
875
876
877
878
879
880
881
882
883
884
885
886
887
888
889
890
891
892
893
894
            {
                updateFace
                (
                    faceI,
                    cellI,
                    neighbourWallInfo,
                    propagationTol_,
                    currentWallInfo
                );
            }
        }

        // Reset status of cell
        changedCell_[cellI] = false;
    }

    // Handled all changed cells by now
    nChangedCells_ = 0;

    if (hasCyclicPatches_)
    {
        // Transfer changed faces across cyclic halves
        handleCyclicPatches();
    }
mattijs's avatar
mattijs committed
895
    if (Pstream::parRun())
896
897
898
899
900
901
902
903
904
905
906
907
908
909
910
911
912
913
914
915
    {
        // Transfer changed faces from neighbouring processors.
        handleProcPatches();
    }

    if (debug)
    {
        Pout<< " Changed faces            : " << nChangedFaces_ << endl;
    }

    // Sum nChangedFaces over all procs
    label totNChanged = nChangedFaces_;

    reduce(totNChanged, sumOp<label>());

    return totNChanged;
}


// Iterate
916
917
template <class Type, class TrackingData>
Foam::label Foam::FaceCellWave<Type, TrackingData>::iterate(const label maxIter)
918
919
920
921
922
923
{
    if (hasCyclicPatches_)
    {
        // Transfer changed faces across cyclic halves
        handleCyclicPatches();
    }
mattijs's avatar
mattijs committed
924
    if (Pstream::parRun())
925
926
927
928
929
    {
        // Transfer changed faces from neighbouring processors.
        handleProcPatches();
    }

930
    while (iter_ < maxIter)
931
932
933
934
935
936
937
938
939
940
941
942
943
944
945
946
947
948
949
950
951
952
953
954
955
956
957
958
959
960
961
962
963
964
965
966
967
968
969
970
971
972
    {
        if (debug)
        {
            Pout<< " Iteration " << iter_ << endl;
        }

        nEvals_ = 0;

        label nCells = faceToCell();

        if (debug)
        {
            Pout<< " Total changed cells      : " << nCells << endl;
        }

        if (nCells == 0)
        {
            break;
        }

        label nFaces = cellToFace();

        if (debug)
        {
            Pout<< " Total changed faces      : " << nFaces << nl
                << " Total evaluations        : " << nEvals_ << nl
                << " Remaining unvisited cells: " << nUnvisitedCells_ << nl
                << " Remaining unvisited faces: " << nUnvisitedFaces_ << endl;
        }

        if (nFaces == 0)
        {
            break;
        }

        ++iter_;
    }

    return nUnvisitedCells_;
}

// ************************************************************************* //