surfaceFieldValue.C 27.6 KB
Newer Older
1
2
3
4
/*---------------------------------------------------------------------------*\
  =========                 |
  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   \\    /   O peration     |
OpenFOAM bot's avatar
OpenFOAM bot committed
5
    \\  /    A nd           | www.openfoam.com
OpenFOAM bot's avatar
OpenFOAM bot committed
6
7
     \\/     M anipulation  |
-------------------------------------------------------------------------------
OpenFOAM bot's avatar
OpenFOAM bot committed
8
9
    Copyright (C) 2011-2017 OpenFOAM Foundation
    Copyright (C) 2017-2019 OpenCFD Ltd.
10
11
12
13
-------------------------------------------------------------------------------
License
    This file is part of OpenFOAM.

14
15
16
17
    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.
18
19
20
21
22
23
24

    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
25
    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
26
27
28

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

29
#include "surfaceFieldValue.H"
30
31
#include "fvMesh.H"
#include "emptyPolyPatch.H"
32
#include "coupledPolyPatch.H"
33
#include "sampledSurface.H"
34
35
#include "mergePoints.H"
#include "indirectPrimitivePatch.H"
36
#include "PatchTools.H"
37
#include "addToRunTimeSelectionTable.H"
38
39
40

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

41
namespace Foam
42
{
43
44
45
46
namespace functionObjects
{
namespace fieldValues
{
47
48
49
    defineTypeNameAndDebug(surfaceFieldValue, 0);
    addToRunTimeSelectionTable(fieldValue, surfaceFieldValue, dictionary);
    addToRunTimeSelectionTable(functionObject, surfaceFieldValue, dictionary);
50
51
}
}
52
53
54
}


55
const Foam::Enum
56
<
57
58
59
    Foam::functionObjects::fieldValues::surfaceFieldValue::regionTypes
>
Foam::functionObjects::fieldValues::surfaceFieldValue::regionTypeNames_
Mark OLESEN's avatar
Mark OLESEN committed
60
({
61
62
    { regionTypes::stFaceZone, "faceZone" },
    { regionTypes::stPatch, "patch" },
63
    { regionTypes::stObject, "functionObjectSurface" },
64
    { regionTypes::stSampled, "sampledSurface" },
Mark OLESEN's avatar
Mark OLESEN committed
65
});
66

67
68

const Foam::Enum
69
<
70
71
72
    Foam::functionObjects::fieldValues::surfaceFieldValue::operationType
>
Foam::functionObjects::fieldValues::surfaceFieldValue::operationTypeNames_
Mark OLESEN's avatar
Mark OLESEN committed
73
({
74
    // Normal operations
75
    { operationType::opNone, "none" },
76
77
    { operationType::opMin, "min" },
    { operationType::opMax, "max" },
78
79
80
81
82
83
84
85
86
87
    { operationType::opSum, "sum" },
    { operationType::opSumMag, "sumMag" },
    { operationType::opSumDirection, "sumDirection" },
    { operationType::opSumDirectionBalance, "sumDirectionBalance" },
    { operationType::opAverage, "average" },
    { operationType::opAreaAverage, "areaAverage" },
    { operationType::opAreaIntegrate, "areaIntegrate" },
    { operationType::opCoV, "CoV" },
    { operationType::opAreaNormalAverage, "areaNormalAverage" },
    { operationType::opAreaNormalIntegrate, "areaNormalIntegrate" },
88
    { operationType::opUniformity, "uniformity" },
89
90
91
92
93
94

    // Using weighting
    { operationType::opWeightedSum, "weightedSum" },
    { operationType::opWeightedAverage, "weightedAverage" },
    { operationType::opWeightedAreaAverage, "weightedAreaAverage" },
    { operationType::opWeightedAreaIntegrate, "weightedAreaIntegrate" },
95
    { operationType::opWeightedUniformity, "weightedUniformity" },
96
97
98
99
100
101

    // Using absolute weighting
    { operationType::opAbsWeightedSum, "absWeightedSum" },
    { operationType::opAbsWeightedAverage, "absWeightedAverage" },
    { operationType::opAbsWeightedAreaAverage, "absWeightedAreaAverage" },
    { operationType::opAbsWeightedAreaIntegrate, "absWeightedAreaIntegrate" },
102
    { operationType::opAbsWeightedUniformity, "absWeightedUniformity" },
Mark OLESEN's avatar
Mark OLESEN committed
103
});
104

105
const Foam::Enum
106
<
107
    Foam::functionObjects::fieldValues::surfaceFieldValue::postOperationType
108
>
109
Foam::functionObjects::fieldValues::surfaceFieldValue::postOperationTypeNames_
Mark OLESEN's avatar
Mark OLESEN committed
110
({
111
112
    { postOperationType::postOpNone, "none" },
    { postOperationType::postOpSqrt, "sqrt" },
Mark OLESEN's avatar
Mark OLESEN committed
113
});
114

115

Henry Weller's avatar
Henry Weller committed
116
// * * * * * * * * * * * * Protected Member Functions  * * * * * * * * * * * //
117

118
119
120
const Foam::objectRegistry&
Foam::functionObjects::fieldValues::surfaceFieldValue::obr() const
{
121
    if (stObject == regionType_)
122
    {
123
        return storedObjects().lookupObject<polySurface>(regionName_);
124
    }
125
126

    return mesh_;
127
128
129
}


130
void Foam::functionObjects::fieldValues::surfaceFieldValue::setFaceZoneFaces()
131
{
132
    const label zoneId = mesh_.faceZones().findZoneID(regionName_);
133
134
135

    if (zoneId < 0)
    {
136
        FatalErrorInFunction
137
            << type() << " " << name() << ": "
138
            << regionTypeNames_[regionType_] << '(' << regionName_ << "):" << nl
139
            << "    Unknown face zone name: " << regionName_
140
            << ". Valid face zones are: " << mesh_.faceZones().names()
141
142
143
            << nl << exit(FatalError);
    }

144
    const faceZone& fZone = mesh_.faceZones()[zoneId];
145

146
147
    DynamicList<label> faceIds(fZone.size());
    DynamicList<label> facePatchIds(fZone.size());
148
    DynamicList<bool> faceFlip(fZone.size());
149
150
151

    forAll(fZone, i)
    {
152
        const label facei = fZone[i];
153
154
155

        label faceId = -1;
        label facePatchId = -1;
156
        if (mesh_.isInternalFace(facei))
157
        {
158
            faceId = facei;
159
160
161
162
            facePatchId = -1;
        }
        else
        {
163
164
            facePatchId = mesh_.boundaryMesh().whichPatch(facei);
            const polyPatch& pp = mesh_.boundaryMesh()[facePatchId];
165
            if (isA<coupledPolyPatch>(pp))
166
            {
167
                if (refCast<const coupledPolyPatch>(pp).owner())
168
                {
169
                    faceId = pp.whichFace(facei);
170
171
172
173
174
175
176
177
                }
                else
                {
                    faceId = -1;
                }
            }
            else if (!isA<emptyPolyPatch>(pp))
            {
178
                faceId = facei - pp.start();
179
180
181
182
183
184
185
186
187
188
            }
            else
            {
                faceId = -1;
                facePatchId = -1;
            }
        }

        if (faceId >= 0)
        {
189
190
            faceIds.append(faceId);
            facePatchIds.append(facePatchId);
191
            faceFlip.append(fZone.flipMap()[i] ? true : false);
192
193
194
        }
    }

195
196
    faceId_.transfer(faceIds);
    facePatchId_.transfer(facePatchIds);
197
    faceFlip_.transfer(faceFlip);
198
    nFaces_ = returnReduce(faceId_.size(), sumOp<label>());
199
200
201

    if (debug)
    {
202
        Pout<< "Original face zone size = " << fZone.size()
203
            << ", new size = " << faceId_.size() << endl;
204
205
206
207
    }
}


208
void Foam::functionObjects::fieldValues::surfaceFieldValue::setPatchFaces()
209
{
210
    const label patchid = mesh_.boundaryMesh().findPatchID(regionName_);
211

212
    if (patchid < 0)
213
    {
214
        FatalErrorInFunction
215
            << type() << " " << name() << ": "
216
            << regionTypeNames_[regionType_] << '(' << regionName_ << "):" << nl
217
            << "    Unknown patch name: " << regionName_
218
            << ". Valid patch names are: "
219
            << mesh_.boundaryMesh().names() << nl
220
221
222
            << exit(FatalError);
    }

223
    const polyPatch& pp = mesh_.boundaryMesh()[patchid];
224
225

    label nFaces = pp.size();
226
    if (isA<emptyPolyPatch>(pp))
227
228
229
230
231
232
    {
        nFaces = 0;
    }

    faceId_.setSize(nFaces);
    facePatchId_.setSize(nFaces);
233
    faceFlip_.setSize(nFaces);
234
    nFaces_ = returnReduce(faceId_.size(), sumOp<label>());
235

236
    forAll(faceId_, facei)
237
    {
238
239
        faceId_[facei] = facei;
        facePatchId_[facei] = patchid;
240
        faceFlip_[facei] = false;
241
242
243
244
    }
}


245
void Foam::functionObjects::fieldValues::surfaceFieldValue::combineMeshGeometry
246
247
248
249
250
251
252
253
254
255
256
257
258
(
    faceList& faces,
    pointField& points
) const
{
    List<faceList> allFaces(Pstream::nProcs());
    List<pointField> allPoints(Pstream::nProcs());

    labelList globalFacesIs(faceId_);
    forAll(globalFacesIs, i)
    {
        if (facePatchId_[i] != -1)
        {
259
            const label patchi = facePatchId_[i];
260
            globalFacesIs[i] += mesh_.boundaryMesh()[patchi].start();
261
262
263
264
265
266
        }
    }

    // Add local faces and points to the all* lists
    indirectPrimitivePatch pp
    (
267
268
        IndirectList<face>(mesh_.faces(), globalFacesIs),
        mesh_.points()
269
270
271
272
273
274
275
276
277
278
    );
    allFaces[Pstream::myProcNo()] = pp.localFaces();
    allPoints[Pstream::myProcNo()] = pp.localPoints();

    Pstream::gatherList(allFaces);
    Pstream::gatherList(allPoints);

    // Renumber and flatten
    label nFaces = 0;
    label nPoints = 0;
279
    forAll(allFaces, proci)
280
    {
281
282
        nFaces += allFaces[proci].size();
        nPoints += allPoints[proci].size();
283
284
285
286
287
288
289
290
291
292
293
    }

    faces.setSize(nFaces);
    points.setSize(nPoints);

    nFaces = 0;
    nPoints = 0;

    // My own data first
    {
        const faceList& fcs = allFaces[Pstream::myProcNo()];
294
        for (const face& f : fcs)
295
296
297
298
299
300
301
302
303
304
        {
            face& newF = faces[nFaces++];
            newF.setSize(f.size());
            forAll(f, fp)
            {
                newF[fp] = f[fp] + nPoints;
            }
        }

        const pointField& pts = allPoints[Pstream::myProcNo()];
305
        for (const point& pt : pts)
306
        {
307
            points[nPoints++] = pt;
308
309
310
311
        }
    }

    // Other proc data follows
312
    forAll(allFaces, proci)
313
    {
314
        if (proci != Pstream::myProcNo())
315
        {
316
            const faceList& fcs = allFaces[proci];
317
            for (const face& f : fcs)
318
319
320
321
322
323
324
325
326
            {
                face& newF = faces[nFaces++];
                newF.setSize(f.size());
                forAll(f, fp)
                {
                    newF[fp] = f[fp] + nPoints;
                }
            }

327
            const pointField& pts = allPoints[proci];
328
            for (const point& pt : pts)
329
            {
330
                points[nPoints++] = pt;
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
            }
        }
    }

    // Merge
    labelList oldToNew;
    pointField newPoints;
    bool hasMerged = mergePoints
    (
        points,
        SMALL,
        false,
        oldToNew,
        newPoints
    );

    if (hasMerged)
    {
        if (debug)
        {
            Pout<< "Merged from " << points.size()
                << " down to " << newPoints.size() << " points" << endl;
        }

        points.transfer(newPoints);
356
        for (face& f : faces)
357
        {
358
            inplaceRenumber(oldToNew, f);
359
360
361
362
363
        }
    }
}


364
365
void Foam::functionObjects::fieldValues::surfaceFieldValue::
combineSurfaceGeometry
366
367
368
369
370
(
    faceList& faces,
    pointField& points
) const
{
371
    if (stObject == regionType_)
372
    {
373
        const polySurface& s = dynamicCast<const polySurface>(obr());
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400

        if (Pstream::parRun())
        {
            // Dimension as fraction of surface
            const scalar mergeDim = 1e-10*boundBox(s.points(), true).mag();

            labelList pointsMap;

            PatchTools::gatherAndMerge
            (
                mergeDim,
                primitivePatch
                (
                    SubList<face>(s.faces(), s.faces().size()),
                    s.points()
                ),
                points,
                faces,
                pointsMap
            );
        }
        else
        {
            faces = s.faces();
            points = s.points();
        }
    }
401
    else if (sampledPtr_.valid())
402
    {
403
        const sampledSurface& s = sampledPtr_();
404
405
406

        if (Pstream::parRun())
        {
407
            // Dimension as fraction of mesh bounding box
408
            const scalar mergeDim = 1e-10*mesh_.bounds().mag();
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433

            labelList pointsMap;

            PatchTools::gatherAndMerge
            (
                mergeDim,
                primitivePatch
                (
                    SubList<face>(s.faces(), s.faces().size()),
                    s.points()
                ),
                points,
                faces,
                pointsMap
            );
        }
        else
        {
            faces = s.faces();
            points = s.points();
        }
    }
}


434
Foam::scalar
435
Foam::functionObjects::fieldValues::surfaceFieldValue::totalArea() const
436
{
437
    scalar totalArea = 0;
438

439
    if (stObject == regionType_)
440
    {
441
        const polySurface& s = dynamicCast<const polySurface>(obr());
442
443
444

        totalArea = gSum(s.magSf());
    }
445
    else if (sampledPtr_.valid())
446
    {
447
        totalArea = gSum(sampledPtr_().magSf());
448
449
450
    }
    else
    {
451
        totalArea = gSum(filterField(mesh_.magSf()));
452
453
454
455
456
457
    }

    return totalArea;
}


458
459
// * * * * * * * * * * * * Protected Member Functions  * * * * * * * * * * * //

460
bool Foam::functionObjects::fieldValues::surfaceFieldValue::usesSf() const
461
{
462
    // Only a few operations do not require the Sf field
463
464
465
    switch (operation_)
    {
        case opNone:
466
467
        case opMin:
        case opMax:
468
469
470
        case opSum:
        case opSumMag:
        case opAverage:
471
        {
472
            return false;
473
        }
474
        default:
475
        {
476
            return true;
477
        }
478
479
480
481
    }
}


482
bool Foam::functionObjects::fieldValues::surfaceFieldValue::update()
483
{
484
485
486
487
    if (sampledPtr_.valid())
    {
        sampledPtr_->update();
    }
488

489
490
491
492
    if (!needsUpdate_)
    {
        return false;
    }
493

494
    switch (regionType_)
495
496
497
498
499
500
501
502
503
    {
        case stFaceZone:
        {
            setFaceZoneFaces();
            break;
        }
        case stPatch:
        {
            setPatchFaces();
504
505
            break;
        }
506
        case stObject:
507
        {
508
            const polySurface& s = dynamicCast<const polySurface>(obr());
509
            nFaces_ = returnReduce(s.size(), sumOp<label>());
510
511
            break;
        }
512
        case stSampled:
513
        {
514
            nFaces_ = returnReduce(sampledPtr_->faces().size(), sumOp<label>());
515
516
            break;
        }
517
518

        // Compiler warning if we forgot an enumeration
519
520
    }

521
522
    if (nFaces_ == 0)
    {
523
        FatalErrorInFunction
524
            << type() << " " << name() << ": "
525
            << regionTypeNames_[regionType_] << '(' << regionName_ << "):" << nl
526
            << "    Region has no faces" << exit(FatalError);
527
528
    }

529
530
    totalArea_ = totalArea();

Andrew Heather's avatar
Andrew Heather committed
531
532
533
    Log << "    total faces   = " << nFaces_ << nl
        << "    total area    = " << totalArea_ << nl
        << endl;
534

535
    writeFileHeader(file());
536

537
538
    needsUpdate_ = false;
    return true;
539
540
541
}


542
543
544
545
void Foam::functionObjects::fieldValues::surfaceFieldValue::writeFileHeader
(
    Ostream& os
) const
546
{
547
    if (operation_ != opNone)
548
    {
549
550
        writeCommented(os, "Region type : ");
        os << regionTypeNames_[regionType_] << " " << regionName_ << endl;
551

552
553
554
555
        writeHeaderValue(os, "Faces", nFaces_);
        writeHeaderValue(os, "Area", totalArea_);
        writeHeaderValue(os, "Scale factor", scaleFactor_);

556
557
558
559
560
        if (weightFieldName_ != "none")
        {
            writeHeaderValue(os, "Weight field", weightFieldName_);
        }

561
        writeCommented(os, "Time");
562
563
        if (writeArea_)
        {
564
            os  << tab << "Area";
565
        }
566

567
        for (const word& fieldName : fields_)
568
        {
569
            os  << tab << operationTypeNames_[operation_]
570
                << "(" << fieldName << ")";
571
        }
572

573
        os  << endl;
574
    }
575
576
577
}


578
template<>
579
580
Foam::scalar
Foam::functionObjects::fieldValues::surfaceFieldValue::processValues
581
582
583
584
585
586
587
588
589
590
(
    const Field<scalar>& values,
    const vectorField& Sf,
    const scalarField& weightField
) const
{
    switch (operation_)
    {
        case opSumDirection:
        {
591
            const vector n(dict_.get<vector>("direction"));
Henry Weller's avatar
Henry Weller committed
592
            return gSum(pos0(values*(Sf & n))*mag(values));
593
594
595
        }
        case opSumDirectionBalance:
        {
596
            const vector n(dict_.get<vector>("direction"));
597
            const scalarField nv(values*(Sf & n));
598

Henry Weller's avatar
Henry Weller committed
599
            return gSum(pos0(nv)*mag(values) - neg(nv)*mag(values));
600
        }
601
602
603
604
605
606
607
608
609
610
611
612
613
614
615
616
617
618
619
620
621
622
623
624
625
626
627
628
629
630
631
632
633
634
635
636
637
638
639

        case opUniformity:
        case opWeightedUniformity:
        case opAbsWeightedUniformity:
        {
            const scalar areaTotal = gSum(mag(Sf));
            tmp<scalarField> areaVal(values * mag(Sf));

            scalar mean, numer;

            if (canWeight(weightField))
            {
                // Weighted quantity = (Weight * phi * dA)

                tmp<scalarField> weight(weightingFactor(weightField));

                // Mean weighted value (area-averaged)
                mean = gSum(weight()*areaVal()) / areaTotal;

                // Abs. deviation from weighted mean value
                numer = gSum(mag(weight*areaVal - (mean * mag(Sf))));
            }
            else
            {
                // Unweighted quantity = (1 * phi * dA)

                // Mean value (area-averaged)
                mean = gSum(areaVal()) / areaTotal;

                // Abs. deviation from mean value
                numer = gSum(mag(areaVal - (mean * mag(Sf))));
            }

            // Uniformity index
            const scalar ui = 1 - numer/(2*mag(mean*areaTotal) + ROOTVSMALL);

            return min(max(ui, 0), 1);
        }

640
641
642
643
644
645
646
647
648
        default:
        {
            // Fall through to other operations
            return processSameTypeValues(values, Sf, weightField);
        }
    }
}


649
template<>
650
651
Foam::vector
Foam::functionObjects::fieldValues::surfaceFieldValue::processValues
652
653
654
655
656
657
658
659
(
    const Field<vector>& values,
    const vectorField& Sf,
    const scalarField& weightField
) const
{
    switch (operation_)
    {
660
661
        case opSumDirection:
        {
662
            const vector n(dict_.get<vector>("direction").normalise());
663

664
            const scalarField nv(n & values);
Henry Weller's avatar
Henry Weller committed
665
            return gSum(pos0(nv)*n*(nv));
666
667
668
        }
        case opSumDirectionBalance:
        {
669
            const vector n(dict_.get<vector>("direction").normalise());
670

671
            const scalarField nv(n & values);
Henry Weller's avatar
Henry Weller committed
672
            return gSum(pos0(nv)*n*(nv));
673
        }
674
675
        case opAreaNormalAverage:
        {
676
677
            const scalar val = gSum(values & Sf)/gSum(mag(Sf));
            return vector(val, 0, 0);
678
679
680
        }
        case opAreaNormalIntegrate:
        {
681
682
            const scalar val = gSum(values & Sf);
            return vector(val, 0, 0);
683
        }
684
685
686
687
688
689
690
691
692
693
694
695
696
697
698
699
700
701
702
703
704
705
706
707
708
709
710
711
712
713
714
715
716
717
718
719
720
721
722

        case opUniformity:
        case opWeightedUniformity:
        case opAbsWeightedUniformity:
        {
            const scalar areaTotal = gSum(mag(Sf));
            tmp<scalarField> areaVal(values & Sf);

            scalar mean, numer;

            if (canWeight(weightField))
            {
                // Weighted quantity = (Weight * phi . dA)

                tmp<scalarField> weight(weightingFactor(weightField));

                // Mean weighted value (area-averaged)
                mean = gSum(weight()*areaVal()) / areaTotal;

                // Abs. deviation from weighted mean value
                numer = gSum(mag(weight*areaVal - (mean * mag(Sf))));
            }
            else
            {
                // Unweighted quantity = (1 * phi . dA)

                // Mean value (area-averaged)
                mean = gSum(areaVal()) / areaTotal;

                // Abs. deviation from mean value
                numer = gSum(mag(areaVal - (mean * mag(Sf))));
            }

            // Uniformity index
            const scalar ui = 1 - numer/(2*mag(mean*areaTotal) + ROOTVSMALL);

            return vector(min(max(ui, 0), 1), 0, 0);
        }

723
724
725
726
727
728
729
730
731
        default:
        {
            // Fall through to other operations
            return processSameTypeValues(values, Sf, weightField);
        }
    }
}


732
733
734
735
736
template<>
Foam::tmp<Foam::scalarField>
Foam::functionObjects::fieldValues::surfaceFieldValue::weightingFactor
(
    const Field<scalar>& weightField
737
) const
738
{
739
740
741
742
    if (usesMag())
    {
        return mag(weightField);
    }
743
744
745

    // pass through
    return weightField;
746
747
748
749
750
751
752
753
754
}


template<>
Foam::tmp<Foam::scalarField>
Foam::functionObjects::fieldValues::surfaceFieldValue::weightingFactor
(
    const Field<scalar>& weightField,
    const vectorField& Sf
755
) const
756
757
758
759
{
    // scalar * Area
    if (returnReduce(weightField.empty(), andOp<bool>()))
    {
760
        // No weight field - revert to unweighted form
761
762
        return mag(Sf);
    }
763
764
765
766
    else if (usesMag())
    {
        return mag(weightField * mag(Sf));
    }
767
768

    return (weightField * mag(Sf));
769
770
771
772
773
774
775
776
777
}


template<>
Foam::tmp<Foam::scalarField>
Foam::functionObjects::fieldValues::surfaceFieldValue::weightingFactor
(
    const Field<vector>& weightField,
    const vectorField& Sf
778
) const
779
780
781
782
{
    // vector (dot) Area
    if (returnReduce(weightField.empty(), andOp<bool>()))
    {
783
        // No weight field - revert to unweighted form
784
785
        return mag(Sf);
    }
786
787
788
789
    else if (usesMag())
    {
        return mag(weightField & Sf);
    }
790
791

    return (weightField & Sf);
792
793
794
}


795
796
// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //

797
Foam::functionObjects::fieldValues::surfaceFieldValue::surfaceFieldValue
798
799
800
801
802
803
804
(
    const word& name,
    const Time& runTime,
    const dictionary& dict
)
:
    fieldValue(name, runTime, dict, typeName),
Mark OLESEN's avatar
Mark OLESEN committed
805
806
    regionType_(regionTypeNames_.get("regionType", dict)),
    operation_(operationTypeNames_.get("operation", dict)),
807
808
    postOperation_
    (
809
810
811
812
        postOperationTypeNames_.lookupOrDefault
        (
            "postOperation",
            dict,
Mark OLESEN's avatar
Mark OLESEN committed
813
814
            postOperationType::postOpNone,
            true  // Failsafe behaviour
815
        )
816
    ),
817
    weightFieldName_("none"),
818
819
    needsUpdate_(true),
    writeArea_(false),
820
    totalArea_(0),
821
822
823
    nFaces_(0),
    faceId_(),
    facePatchId_(),
824
    faceFlip_()
825
826
827
828
{
    read(dict);
}

829

830
Foam::functionObjects::fieldValues::surfaceFieldValue::surfaceFieldValue
831
832
833
(
    const word& name,
    const objectRegistry& obr,
834
    const dictionary& dict
835
836
)
:
837
    fieldValue(name, obr, dict, typeName),
Mark OLESEN's avatar
Mark OLESEN committed
838
839
    regionType_(regionTypeNames_.get("regionType", dict)),
    operation_(operationTypeNames_.get("operation", dict)),
840
841
    postOperation_
    (
842
843
844
845
        postOperationTypeNames_.lookupOrDefault
        (
            "postOperation",
            dict,
Mark OLESEN's avatar
Mark OLESEN committed
846
847
            postOperationType::postOpNone,
            true  // Failsafe behaviour
848
        )
849
    ),
850
    weightFieldName_("none"),
851
852
    needsUpdate_(true),
    writeArea_(false),
853
    totalArea_(0),
854
    nFaces_(0),
855
856
    faceId_(),
    facePatchId_(),
857
    faceFlip_()
858
{
859
    read(dict);
860
861
862
863
864
}


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

865
bool Foam::functionObjects::fieldValues::surfaceFieldValue::read
866
867
868
(
    const dictionary& dict
)
869
{
Andrew Heather's avatar
Andrew Heather committed
870
    fieldValue::read(dict);
871
872
873
874
875
876
877
878
879
880
881
882
883
884
885
886
887
888
889
890
891
892
893
894
895
896
897
898
899
900
901
902
903
904
905
906
907
908
909
910
911
912
913
914
915
916
917
918
919
920
921
922
923
924
925
926
927
928
929
930
931
932
933
934
935
936
937
938
939
940
941
942
943
944
945
946
947
948
949
950
951
952
953

    weightFieldName_ = "none";
    needsUpdate_ = true;
    writeArea_ = dict.lookupOrDefault("writeArea", false);
    totalArea_ = 0;
    nFaces_ = 0;
    faceId_.clear();
    facePatchId_.clear();
    faceFlip_.clear();
    sampledPtr_.clear();
    surfaceWriterPtr_.clear();

    dict.readEntry("name", regionName_);

    // Create sampled surface, but leave 'expired' (ie, no update) since it
    // may depend on fields or data that do not yet exist
    if (stSampled == regionType_)
    {
        sampledPtr_ = sampledSurface::New
        (
            name(),
            mesh_,
            dict.subDict("sampledSurfaceDict")
        );
    }

    Info<< type() << " " << name() << ":" << nl
        << "    operation     = ";

    if (postOperation_ != postOpNone)
    {
        Info<< postOperationTypeNames_[postOperation_] << '('
            << operationTypeNames_[operation_] << ')'  << nl;
    }
    else
    {
        Info<< operationTypeNames_[operation_] << nl;
    }

    if (usesWeight())
    {
        if (stSampled == regionType_)
        {
            FatalIOErrorInFunction(dict)
                << "Cannot use weighted operation '"
                << operationTypeNames_[operation_]
                << "' for sampledSurface"
                << exit(FatalIOError);
        }

        if (dict.readIfPresent("weightField", weightFieldName_))
        {
            Info<< "    weight field  = " << weightFieldName_ << nl;
        }
        else
        {
            // Suggest possible alternative unweighted operation?
            FatalIOErrorInFunction(dict)
                << "The '" << operationTypeNames_[operation_]
                << "' operation is missing a weightField." << nl
                << "Either provide the weightField, "
                << "use weightField 'none' to suppress weighting," << nl
                << "or use a different operation."
                << exit(FatalIOError);
        }
    }

    // Backwards compatibility for v1612 and older
    List<word> orientedFields;
    if (dict.readIfPresent("orientedFields", orientedFields))
    {
        fields_.append(orientedFields);

        WarningInFunction
            << "The 'orientedFields' option is deprecated.  These fields can "
            << "and have been added to the standard 'fields' list."
            << endl;
    }

    if (writeFields_)
    {
        const word formatName(dict.get<word>("surfaceFormat"));

954
955
956
        surfaceWriterPtr_.reset
        (
            surfaceWriter::New
957
            (
958
959
960
961
                formatName,
                dict.subOrEmptyDict("formatOptions").subOrEmptyDict(formatName)
            )
        );
962

963
964
965
966
967
        if (debug)
        {
            surfaceWriterPtr_->verbose() = true;
        }

968
969
        if (surfaceWriterPtr_->enabled())
        {
970
971
            Info<< "    surfaceFormat = " << formatName << nl;
        }
972
973
974
975
        else
        {
            surfaceWriterPtr_->clear();
        }
976
977
978
    }

    Info<< nl << endl;
979

980
    return true;
981
982
983
}


984
bool Foam::functionObjects::fieldValues::surfaceFieldValue::write()
985
{
Andrew Heather's avatar
Andrew Heather committed
986
987
988
989
990
    if (needsUpdate_ || operation_ != opNone)
    {
        fieldValue::write();
    }

991
    update();
992

993
    if (operation_ != opNone)
994
    {
Andrew Heather's avatar
Andrew Heather committed
995
        writeTime(file());
996
    }
997

998
999
1000
    if (writeArea_)
    {
        totalArea_ = totalArea();
1001
1002
        Log << "    total area = " << totalArea_ << endl;

1003
        if (operation_ != opNone && Pstream::master())
1004
        {
1005
            file() << tab << totalArea_;
1006
        }
1007
    }
1008

1009
1010
    // Many operations use the Sf field
    vectorField Sf;
1011
    if (usesSf())
1012
    {
1013
        if (stObject == regionType_)
1014
        {
1015
            const polySurface& s = dynamicCast<const polySurface>(obr());
1016
1017
            Sf = s.Sf();
        }
1018
        else if (sampledPtr_.valid())
1019
        {
1020
            Sf = sampledPtr_().Sf();
1021
1022
1023
        }
        else
        {
1024
            Sf = filterField(mesh_.Sf());
1025
1026
1027
1028
1029
1030
1031
1032
1033
        }
    }

    // Faces and points for surface format (if specified)
    faceList faces;
    pointField points;

    if (surfaceWriterPtr_.valid())
    {
1034
        if (withTopologicalMerge())
1035
        {
1036
            combineMeshGeometry(faces, points);
1037
1038
1039
        }
        else
        {
1040
            combineSurfaceGeometry(faces, points);
1041
1042
1043
1044
        }
    }

    // Only a few weight types (scalar, vector)
1045
1046
    if (weightFieldName_ != "none")
    {
1047
1048
        if (validField<scalar>(weightFieldName_))
        {
1049
            scalarField weightField
1050
            (
1051
                getFieldValues<scalar>(weightFieldName_, true)
1052
            );
1053

1054
            // Process the fields
1055
            writeAll(Sf, weightField, points, faces);
1056
1057
1058
        }
        else if (validField<vector>(weightFieldName_))
        {
1059
            vectorField weightField
1060
            (
1061
                getFieldValues<vector>(weightFieldName_, true)
1062
            );
1063

1064
            // Process the fields
1065
            writeAll(Sf, weightField, points, faces);
1066
1067
        }
        else
1068
        {
1069
1070
1071
1072
            FatalErrorInFunction
                << "weightField " << weightFieldName_
                << " not found or an unsupported ty