surfaceFieldValue.C 27.6 KB
Newer Older
Henry's avatar
Henry committed
1 2 3 4
/*---------------------------------------------------------------------------*\
  =========                 |
  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   \\    /   O peration     |
5
    \\  /    A nd           | Copyright (C) 2011-2017 OpenFOAM Foundation
6
     \\/     M anipulation  | Copyright (C) 2017-2018 OpenCFD Ltd.
Henry's avatar
Henry committed
7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25
-------------------------------------------------------------------------------
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 3 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, see <http://www.gnu.org/licenses/>.

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

26
#include "surfaceFieldValue.H"
Henry's avatar
Henry committed
27 28 29 30 31 32 33
#include "fvMesh.H"
#include "emptyPolyPatch.H"
#include "coupledPolyPatch.H"
#include "sampledSurface.H"
#include "mergePoints.H"
#include "indirectPrimitivePatch.H"
#include "PatchTools.H"
34
#include "meshedSurfRef.H"
Henry's avatar
Henry committed
35 36 37 38 39 40
#include "addToRunTimeSelectionTable.H"

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

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


53
const Foam::Enum
54
<
55 56 57
    Foam::functionObjects::fieldValues::surfaceFieldValue::regionTypes
>
Foam::functionObjects::fieldValues::surfaceFieldValue::regionTypeNames_
Mark Olesen's avatar
Mark Olesen committed
58
({
59 60 61 62
    { regionTypes::stFaceZone, "faceZone" },
    { regionTypes::stPatch, "patch" },
    { regionTypes::stSurface, "surface" },
    { regionTypes::stSampledSurface, "sampledSurface" },
Mark Olesen's avatar
Mark Olesen committed
63
});
64

65 66

const Foam::Enum
67
<
68 69 70
    Foam::functionObjects::fieldValues::surfaceFieldValue::operationType
>
Foam::functionObjects::fieldValues::surfaceFieldValue::operationTypeNames_
Mark Olesen's avatar
Mark Olesen committed
71
({
72
    // Normal operations
73
    { operationType::opNone, "none" },
74 75
    { operationType::opMin, "min" },
    { operationType::opMax, "max" },
76 77 78 79 80 81 82 83 84 85
    { 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" },
86
    { operationType::opUniformity, "uniformity" },
87 88 89 90 91 92

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

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

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

Henry's avatar
Henry committed
113

Henry Weller's avatar
Henry Weller committed
114
// * * * * * * * * * * * * Protected Member Functions  * * * * * * * * * * * //
Henry's avatar
Henry committed
115

116 117 118 119 120 121 122
const Foam::objectRegistry&
Foam::functionObjects::fieldValues::surfaceFieldValue::obr() const
{
    if (regionType_ == stSurface)
    {
        return mesh_.lookupObject<objectRegistry>(regionName_);
    }
123 124

    return mesh_;
125 126 127
}


128
void Foam::functionObjects::fieldValues::surfaceFieldValue::setFaceZoneFaces()
Henry's avatar
Henry committed
129
{
130
    const label zoneId = mesh_.faceZones().findZoneID(regionName_);
Henry's avatar
Henry committed
131 132 133

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

142
    const faceZone& fZone = mesh_.faceZones()[zoneId];
Henry's avatar
Henry committed
143 144 145

    DynamicList<label> faceIds(fZone.size());
    DynamicList<label> facePatchIds(fZone.size());
146
    DynamicList<bool> faceFlip(fZone.size());
Henry's avatar
Henry committed
147 148 149

    forAll(fZone, i)
    {
150
        const label facei = fZone[i];
Henry's avatar
Henry committed
151 152 153

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

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

    faceId_.transfer(faceIds);
    facePatchId_.transfer(facePatchIds);
195
    faceFlip_.transfer(faceFlip);
Henry's avatar
Henry committed
196 197 198 199 200 201 202 203 204 205
    nFaces_ = returnReduce(faceId_.size(), sumOp<label>());

    if (debug)
    {
        Pout<< "Original face zone size = " << fZone.size()
            << ", new size = " << faceId_.size() << endl;
    }
}


206
void Foam::functionObjects::fieldValues::surfaceFieldValue::setPatchFaces()
Henry's avatar
Henry committed
207
{
208
    const label patchid = mesh_.boundaryMesh().findPatchID(regionName_);
Henry's avatar
Henry committed
209

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

221
    const polyPatch& pp = mesh_.boundaryMesh()[patchid];
Henry's avatar
Henry committed
222 223 224 225 226 227 228 229 230

    label nFaces = pp.size();
    if (isA<emptyPolyPatch>(pp))
    {
        nFaces = 0;
    }

    faceId_.setSize(nFaces);
    facePatchId_.setSize(nFaces);
231
    faceFlip_.setSize(nFaces);
Henry's avatar
Henry committed
232 233
    nFaces_ = returnReduce(faceId_.size(), sumOp<label>());

234
    forAll(faceId_, facei)
Henry's avatar
Henry committed
235
    {
236 237
        faceId_[facei] = facei;
        facePatchId_[facei] = patchid;
238
        faceFlip_[facei] = false;
Henry's avatar
Henry committed
239 240 241 242
    }
}


243
void Foam::functionObjects::fieldValues::surfaceFieldValue::combineMeshGeometry
Henry's avatar
Henry committed
244 245 246 247 248 249 250 251 252 253 254 255 256
(
    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)
        {
257
            const label patchi = facePatchId_[i];
258
            globalFacesIs[i] += mesh_.boundaryMesh()[patchi].start();
Henry's avatar
Henry committed
259 260 261 262 263 264
        }
    }

    // Add local faces and points to the all* lists
    indirectPrimitivePatch pp
    (
265 266
        IndirectList<face>(mesh_.faces(), globalFacesIs),
        mesh_.points()
Henry's avatar
Henry committed
267 268 269 270 271 272 273 274 275 276
    );
    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;
277
    forAll(allFaces, proci)
Henry's avatar
Henry committed
278
    {
279 280
        nFaces += allFaces[proci].size();
        nPoints += allPoints[proci].size();
Henry's avatar
Henry committed
281 282 283 284 285 286 287 288 289 290 291
    }

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

    nFaces = 0;
    nPoints = 0;

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

        const pointField& pts = allPoints[Pstream::myProcNo()];
303
        for (const point& pt : pts)
Henry's avatar
Henry committed
304
        {
305
            points[nPoints++] = pt;
Henry's avatar
Henry committed
306 307 308 309
        }
    }

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

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

    // 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);
354
        for (face& f : faces)
Henry's avatar
Henry committed
355
        {
356
            inplaceRenumber(oldToNew, f);
Henry's avatar
Henry committed
357 358 359 360 361
        }
    }
}


362 363
void Foam::functionObjects::fieldValues::surfaceFieldValue::
combineSurfaceGeometry
Henry's avatar
Henry committed
364 365 366 367 368
(
    faceList& faces,
    pointField& points
) const
{
369 370 371 372 373 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
    if (regionType_ == stSurface)
    {
        const surfMesh& s = dynamicCast<const surfMesh>(obr());

        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();
        }
    }
    else if (surfacePtr_.valid())
Henry's avatar
Henry committed
400 401 402 403 404
    {
        const sampledSurface& s = surfacePtr_();

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

            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();
        }
    }
}


432
Foam::scalar
433
Foam::functionObjects::fieldValues::surfaceFieldValue::totalArea() const
434 435 436
{
    scalar totalArea;

437 438 439 440 441 442 443
    if (regionType_ == stSurface)
    {
        const surfMesh& s = dynamicCast<const surfMesh>(obr());

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

    return totalArea;
}


Henry's avatar
Henry committed
456 457
// * * * * * * * * * * * * Protected Member Functions  * * * * * * * * * * * //

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


480
void Foam::functionObjects::fieldValues::surfaceFieldValue::initialise
481 482 483
(
    const dictionary& dict
)
Henry's avatar
Henry committed
484
{
485
    dict.readEntry("name", regionName_);
486

487
    switch (regionType_)
Henry's avatar
Henry committed
488 489 490 491
    {
        case stFaceZone:
        {
            setFaceZoneFaces();
492
            surfacePtr_.clear();
Henry's avatar
Henry committed
493 494 495 496 497
            break;
        }
        case stPatch:
        {
            setPatchFaces();
498 499 500 501 502 503 504 505 506 507 508 509
            surfacePtr_.clear();
            break;
        }
        case stSurface:
        {
            const surfMesh& s = dynamicCast<const surfMesh>(obr());
            nFaces_ = returnReduce(s.size(), sumOp<label>());

            faceId_.clear();
            facePatchId_.clear();
            faceFlip_.clear();
            surfacePtr_.clear();
Henry's avatar
Henry committed
510 511 512 513
            break;
        }
        case stSampledSurface:
        {
514 515 516 517 518 519 520 521 522 523 524 525 526 527
            faceId_.clear();
            facePatchId_.clear();
            faceFlip_.clear();

            surfacePtr_ = sampledSurface::New
            (
                name(),
                mesh_,
                dict.subDict("sampledSurfaceDict")
            );
            surfacePtr_().update();
            nFaces_ =
                returnReduce(surfacePtr_().faces().size(), sumOp<label>());

Henry's avatar
Henry committed
528 529 530 531
            break;
        }
        default:
        {
532
            FatalErrorInFunction
533
                << type() << " " << name() << ": "
Mark Olesen's avatar
Mark Olesen committed
534
                << int(regionType_) << "(" << regionName_ << "):"
535
                << nl << "    Unknown region type. Valid region types are:"
Mark Olesen's avatar
Mark Olesen committed
536 537
                << regionTypeNames_ << nl
                << exit(FatalError);
Henry's avatar
Henry committed
538 539 540 541 542
        }
    }

    if (nFaces_ == 0)
    {
543
        FatalErrorInFunction
544
            << type() << " " << name() << ": "
545
            << regionTypeNames_[regionType_] << "(" << regionName_ << "):" << nl
546
            << "    Region has no faces" << exit(FatalError);
Henry's avatar
Henry committed
547 548 549 550 551 552 553
    }

    if (surfacePtr_.valid())
    {
        surfacePtr_().update();
    }

554 555
    totalArea_ = totalArea();

556
    Info<< type() << " " << name() << ":" << nl
557
        << "    operation     = ";
Henry's avatar
Henry committed
558

559
    if (postOperation_ != postOpNone)
Henry's avatar
Henry committed
560
    {
561 562 563 564 565 566
        Info<< postOperationTypeNames_[postOperation_] << '('
            << operationTypeNames_[operation_] << ')'  << nl;
    }
    else
    {
        Info<< operationTypeNames_[operation_] << nl;
Henry's avatar
Henry committed
567
    }
568 569 570
    Info<< "    total faces   = " << nFaces_ << nl
        << "    total area    = " << totalArea_ << nl;

Henry's avatar
Henry committed
571

572
    weightFieldName_ = "none";
573
    if (usesWeight())
Henry's avatar
Henry committed
574
    {
575
        if (regionType_ == stSampledSurface)
Henry's avatar
Henry committed
576
        {
577
            FatalIOErrorInFunction(dict)
578 579 580
                << "Cannot use weighted operation '"
                << operationTypeNames_[operation_]
                << "' for sampledSurface"
581
                << exit(FatalIOError);
Henry's avatar
Henry committed
582
        }
583

584 585 586 587 588 589 590 591 592 593 594 595 596 597 598
        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);
        }
Henry's avatar
Henry committed
599 600
    }

601
    // Backwards compatibility for v1612 and older
Henry's avatar
Henry committed
602 603 604
    List<word> orientedFields;
    if (dict.readIfPresent("orientedFields", orientedFields))
    {
605 606 607 608 609
        WarningInFunction
            << "The 'orientedFields' option is deprecated.  These fields can "
            << "and have been added to the standard 'fields' list."
            << endl;

Henry's avatar
Henry committed
610 611 612
        fields_.append(orientedFields);
    }

613

614
    surfaceWriterPtr_.clear();
615
    if (writeFields_)
Henry's avatar
Henry committed
616
    {
617
        const word surfaceFormat(dict.get<word>("surfaceFormat"));
Henry's avatar
Henry committed
618

619 620 621
        if (surfaceFormat != "none")
        {
            surfaceWriterPtr_.reset
Henry's avatar
Henry committed
622
            (
623 624 625 626 627 628 629 630 631 632
                surfaceWriter::New
                (
                    surfaceFormat,
                    dict.subOrEmptyDict("formatOptions").
                        subOrEmptyDict(surfaceFormat)
                ).ptr()
            );

            Info<< "    surfaceFormat = " << surfaceFormat << nl;
        }
Henry's avatar
Henry committed
633
    }
634 635

    Info<< nl << endl;
Henry's avatar
Henry committed
636 637 638
}


639 640 641 642
void Foam::functionObjects::fieldValues::surfaceFieldValue::writeFileHeader
(
    Ostream& os
) const
Henry's avatar
Henry committed
643
{
644
    if (operation_ != opNone)
645
    {
646 647
        writeCommented(os, "Region type : ");
        os << regionTypeNames_[regionType_] << " " << regionName_ << endl;
Henry's avatar
Henry committed
648

649 650 651 652
        writeHeaderValue(os, "Faces", nFaces_);
        writeHeaderValue(os, "Area", totalArea_);
        writeHeaderValue(os, "Scale factor", scaleFactor_);

653 654 655 656 657
        if (weightFieldName_ != "none")
        {
            writeHeaderValue(os, "Weight field", weightFieldName_);
        }

658
        writeCommented(os, "Time");
659 660
        if (writeArea_)
        {
661
            os  << tab << "Area";
662
        }
Henry's avatar
Henry committed
663

664
        for (const word& fieldName : fields_)
665
        {
666
            os  << tab << operationTypeNames_[operation_]
667
                << "(" << fieldName << ")";
668
        }
Henry's avatar
Henry committed
669

670
        os  << endl;
671
    }
Henry's avatar
Henry committed
672 673 674 675
}


template<>
676 677
Foam::scalar
Foam::functionObjects::fieldValues::surfaceFieldValue::processValues
Henry's avatar
Henry committed
678 679 680 681 682 683 684 685 686 687
(
    const Field<scalar>& values,
    const vectorField& Sf,
    const scalarField& weightField
) const
{
    switch (operation_)
    {
        case opSumDirection:
        {
688
            const vector n(dict_.get<vector>("direction"));
Henry Weller's avatar
Henry Weller committed
689
            return gSum(pos0(values*(Sf & n))*mag(values));
Henry's avatar
Henry committed
690 691 692
        }
        case opSumDirectionBalance:
        {
693
            const vector n(dict_.get<vector>("direction"));
Henry's avatar
Henry committed
694 695
            const scalarField nv(values*(Sf & n));

Henry Weller's avatar
Henry Weller committed
696
            return gSum(pos0(nv)*mag(values) - neg(nv)*mag(values));
Henry's avatar
Henry committed
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 723 724 725 726 727 728 729 730 731 732 733 734 735 736

        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);
        }

Henry's avatar
Henry committed
737 738 739 740 741 742 743 744 745 746
        default:
        {
            // Fall through to other operations
            return processSameTypeValues(values, Sf, weightField);
        }
    }
}


template<>
747 748
Foam::vector
Foam::functionObjects::fieldValues::surfaceFieldValue::processValues
Henry's avatar
Henry committed
749 750 751 752 753 754 755 756 757 758
(
    const Field<vector>& values,
    const vectorField& Sf,
    const scalarField& weightField
) const
{
    switch (operation_)
    {
        case opSumDirection:
        {
759
            const vector n(dict_.get<vector>("direction").normalise());
Henry's avatar
Henry committed
760

761
            const scalarField nv(n & values);
Henry Weller's avatar
Henry Weller committed
762
            return gSum(pos0(nv)*n*(nv));
Henry's avatar
Henry committed
763 764 765
        }
        case opSumDirectionBalance:
        {
766
            const vector n(dict_.get<vector>("direction").normalise());
Henry's avatar
Henry committed
767

768
            const scalarField nv(n & values);
Henry Weller's avatar
Henry Weller committed
769
            return gSum(pos0(nv)*n*(nv));
Henry's avatar
Henry committed
770 771 772
        }
        case opAreaNormalAverage:
        {
773 774
            const scalar val = gSum(values & Sf)/gSum(mag(Sf));
            return vector(val, 0, 0);
Henry's avatar
Henry committed
775 776 777
        }
        case opAreaNormalIntegrate:
        {
778 779
            const scalar val = gSum(values & Sf);
            return vector(val, 0, 0);
Henry's avatar
Henry committed
780
        }
781 782 783 784 785 786 787 788 789 790 791 792 793 794 795 796 797 798 799 800 801 802 803 804 805 806 807 808 809 810 811 812 813 814 815 816 817 818 819

        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);
        }

Henry's avatar
Henry committed
820 821 822 823 824 825 826 827 828
        default:
        {
            // Fall through to other operations
            return processSameTypeValues(values, Sf, weightField);
        }
    }
}


829 830 831 832 833
template<>
Foam::tmp<Foam::scalarField>
Foam::functionObjects::fieldValues::surfaceFieldValue::weightingFactor
(
    const Field<scalar>& weightField
834
) const
835
{
836 837 838 839
    if (usesMag())
    {
        return mag(weightField);
    }
840 841 842

    // pass through
    return weightField;
843 844 845 846 847 848 849 850 851
}


template<>
Foam::tmp<Foam::scalarField>
Foam::functionObjects::fieldValues::surfaceFieldValue::weightingFactor
(
    const Field<scalar>& weightField,
    const vectorField& Sf
852
) const
853 854 855 856
{
    // scalar * Area
    if (returnReduce(weightField.empty(), andOp<bool>()))
    {
857
        // No weight field - revert to unweighted form
858 859
        return mag(Sf);
    }
860 861 862 863
    else if (usesMag())
    {
        return mag(weightField * mag(Sf));
    }
864 865

    return (weightField * mag(Sf));
866 867 868 869 870 871 872 873 874
}


template<>
Foam::tmp<Foam::scalarField>
Foam::functionObjects::fieldValues::surfaceFieldValue::weightingFactor
(
    const Field<vector>& weightField,
    const vectorField& Sf
875
) const
876 877 878 879
{
    // vector (dot) Area
    if (returnReduce(weightField.empty(), andOp<bool>()))
    {
880
        // No weight field - revert to unweighted form
881 882
        return mag(Sf);
    }
883 884 885 886
    else if (usesMag())
    {
        return mag(weightField & Sf);
    }
887 888

    return (weightField & Sf);
889 890 891
}


Henry's avatar
Henry committed
892 893
// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //

894
Foam::functionObjects::fieldValues::surfaceFieldValue::surfaceFieldValue
895 896 897 898 899 900 901
(
    const word& name,
    const Time& runTime,
    const dictionary& dict
)
:
    fieldValue(name, runTime, dict, typeName),
Mark Olesen's avatar
Mark Olesen committed
902 903
    regionType_(regionTypeNames_.get("regionType", dict)),
    operation_(operationTypeNames_.get("operation", dict)),
904 905
    postOperation_
    (
906 907 908 909
        postOperationTypeNames_.lookupOrDefault
        (
            "postOperation",
            dict,
Mark Olesen's avatar
Mark Olesen committed
910 911
            postOperationType::postOpNone,
            true  // Failsafe behaviour
912
        )
913
    ),
914 915 916 917 918
    weightFieldName_("none"),
    writeArea_(dict.lookupOrDefault("writeArea", false)),
    nFaces_(0),
    faceId_(),
    facePatchId_(),
919
    faceFlip_()
920 921
{
    read(dict);
922
    writeFileHeader(file());
923 924
}

925

926
Foam::functionObjects::fieldValues::surfaceFieldValue::surfaceFieldValue
Henry's avatar
Henry committed
927 928 929
(
    const word& name,
    const objectRegistry& obr,
930
    const dictionary& dict
Henry's avatar
Henry committed
931 932
)
:
933
    fieldValue(name, obr, dict, typeName),
Mark Olesen's avatar
Mark Olesen committed
934 935
    regionType_(regionTypeNames_.get("regionType", dict)),
    operation_(operationTypeNames_.get("operation", dict)),
936 937
    postOperation_
    (
938 939 940 941
        postOperationTypeNames_.lookupOrDefault
        (
            "postOperation",
            dict,
Mark Olesen's avatar
Mark Olesen committed
942 943
            postOperationType::postOpNone,
            true  // Failsafe behaviour
944
        )
945
    ),
Henry's avatar
Henry committed
946
    weightFieldName_("none"),
947
    writeArea_(dict.lookupOrDefault("writeArea", false)),
Henry's avatar
Henry committed
948 949 950
    nFaces_(0),
    faceId_(),
    facePatchId_(),
951
    faceFlip_()
Henry's avatar
Henry committed
952
{
953
    read(dict);
954
    writeFileHeader(file());
Henry's avatar
Henry committed
955 956 957 958 959
}


// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //

960
Foam::functionObjects::fieldValues::surfaceFieldValue::~surfaceFieldValue()
Henry's avatar
Henry committed
961 962 963 964 965
{}


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

966
bool Foam::functionObjects::fieldValues::surfaceFieldValue::read
967 968 969
(
    const dictionary& dict
)
Henry's avatar
Henry committed
970 971
{
    fieldValue::read(dict);
972
    initialise(dict);
973

974
    return true;
Henry's avatar
Henry committed
975 976 977
}


978
bool Foam::functionObjects::fieldValues::surfaceFieldValue::write()
Henry's avatar
Henry committed
979
{
980
    if (surfacePtr_.valid())
Henry's avatar
Henry committed
981
    {
982 983
        surfacePtr_().update();
    }
984

985
    if (operation_ != opNone)
986
    {
987 988 989 990 991 992
        fieldValue::write();

        if (Pstream::master())
        {
            writeTime(file());
        }
993
    }
Henry's avatar
Henry committed
994

995 996 997