From 485523eab5345e153c7ed2d21d2d244b907b8d49 Mon Sep 17 00:00:00 2001
From: Andrew Heather <>
Date: Tue, 18 Dec 2018 10:35:34 +0000
Subject: [PATCH] ENH: Added new columnAverage function object

    Averages columns of cells for layered meshes.

    For each patch face, calculates the average value of all cells attached in
    the patch face normal direction, and then pushes the average value back
    to all cells in the column.

    Useful for channel-like cases where we want to average fields in the
    spanwise direction.

    Example of function object specification:
    columnAverage1
    {
        type        columnAverage;
        libs        ("libfieldFunctionObjects.so");
        ...
        patches     (front side);
        fields      (U p);
    }

    Where the entries comprise:
    \table
        Property     | Description               | Required    | Default value
        type         | type name: fieldMinMax    | yes         |
        patches      | list of patches to collapse onto | yes  |
        fields       | list of fields to process | yes         |
    \endtable
---
 src/functionObjects/field/Make/files          |   2 +
 .../field/columnAverage/columnAverage.C       | 185 ++++++++++++++++++
 .../field/columnAverage/columnAverage.H       | 171 ++++++++++++++++
 .../columnAverage/columnAverageTemplates.C    | 114 +++++++++++
 4 files changed, 472 insertions(+)
 create mode 100644 src/functionObjects/field/columnAverage/columnAverage.C
 create mode 100644 src/functionObjects/field/columnAverage/columnAverage.H
 create mode 100644 src/functionObjects/field/columnAverage/columnAverageTemplates.C

diff --git a/src/functionObjects/field/Make/files b/src/functionObjects/field/Make/files
index fbbac57cb6d..1ca07b37660 100644
--- a/src/functionObjects/field/Make/files
+++ b/src/functionObjects/field/Make/files
@@ -1,5 +1,7 @@
 AMIWeights/AMIWeights.C
 
+columnAverage/columnAverage.C
+
 fieldAverage/fieldAverage.C
 fieldAverage/fieldAverageItem/fieldAverageItem.C
 fieldAverage/fieldAverageItem/fieldAverageItemIO.C
diff --git a/src/functionObjects/field/columnAverage/columnAverage.C b/src/functionObjects/field/columnAverage/columnAverage.C
new file mode 100644
index 00000000000..c0b00cfdab3
--- /dev/null
+++ b/src/functionObjects/field/columnAverage/columnAverage.C
@@ -0,0 +1,185 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2018 OpenCFD Ltd.
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 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/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "columnAverage.H"
+#include "volFields.H"
+#include "addToRunTimeSelectionTable.H"
+#include "meshStructure.H"
+#include "globalIndex.H"
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace functionObjects
+{
+    defineTypeNameAndDebug(columnAverage, 0);
+    addToRunTimeSelectionTable(functionObject, columnAverage, dictionary);
+}
+}
+
+
+// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
+
+const Foam::meshStructure&
+Foam::functionObjects::columnAverage::meshAddressing(const polyMesh& mesh) const
+{
+    if (!meshStructurePtr_.valid())
+    {
+        const polyBoundaryMesh& pbm = mesh.boundaryMesh();
+        const labelList patchIDs(patchSet_.sortedToc());
+
+        // Count
+        label sz = 0;
+        for (const label patchi : patchIDs)
+        {
+            sz += pbm[patchi].size();
+        }
+
+        // Fill
+        labelList meshFaces(sz);
+        sz = 0;
+        for (const label patchi : patchIDs)
+        {
+            label start = pbm[patchi].start();
+            label size = pbm[patchi].size();
+            for (label i = 0; i < size; ++i)
+            {
+                meshFaces[sz++] = start+i;
+            }
+        }
+
+        if (sz == 0)
+        {
+            // TODO: If source patch is a cyclic it may have have been
+            // converted to a processorCyclic for parallel runs
+
+            WarningInFunction
+                << "Requested patches have zero faces"
+                << endl;
+        }
+
+        uindirectPrimitivePatch uip
+        (
+            UIndirectList<face>(mesh.faces(), meshFaces),
+            mesh.points()
+        );
+
+        globalFaces_.set(new globalIndex(uip.size()));
+        globalEdges_.set(new globalIndex(uip.nEdges()));
+        globalPoints_.set(new globalIndex(uip.nPoints()));
+        meshStructurePtr_.set
+        (
+            new meshStructure
+            (
+                mesh,
+                uip,
+                globalFaces_(),
+                globalEdges_(),
+                globalPoints_()
+            )
+        );
+    }
+    return meshStructurePtr_();
+}
+
+
+const Foam::word Foam::functionObjects::columnAverage::averageName
+(
+    const word& fieldName
+) const
+{
+    return name() + ":columnAverage(" + fieldName + ")"; 
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+Foam::functionObjects::columnAverage::columnAverage
+(
+    const word& name,
+    const Time& runTime,
+    const dictionary& dict
+)
+:
+    fvMeshFunctionObject(name, runTime, dict),
+    patchSet_(),
+    fieldSet_(mesh_)
+{
+    read(dict);
+}
+
+
+// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+
+bool Foam::functionObjects::columnAverage::read(const dictionary& dict)
+{
+    fvMeshFunctionObject::read(dict);
+
+    patchSet_ = mesh_.boundaryMesh().patchSet(dict.get<wordRes>("patches"));
+
+    fieldSet_.read(dict);
+
+    return true;
+}
+
+
+bool Foam::functionObjects::columnAverage::execute()
+{
+    // Make fields up to date with current selection
+    fieldSet_.updateSelection();
+
+    for (const word& fieldName : fieldSet_.selection())
+    {
+        columnAverageField<scalar>(fieldName);
+        columnAverageField<vector>(fieldName);
+        columnAverageField<sphericalTensor>(fieldName);
+        columnAverageField<symmTensor>(fieldName);
+        columnAverageField<tensor>(fieldName);
+    }
+
+    return true;
+}
+
+
+bool Foam::functionObjects::columnAverage::write()
+{
+    for (const word& fieldName : fieldSet_.selection())
+    {
+        const word resultName("columnAverage(" + fieldName + ")");
+        const regIOobject* obj =
+            obr_.lookupObjectPtr<regIOobject>(averageName(fieldName));
+
+        if (obj)
+        {
+            obj->write();
+        }
+    }
+
+    return true;
+}
+
+
+// ************************************************************************* //
diff --git a/src/functionObjects/field/columnAverage/columnAverage.H b/src/functionObjects/field/columnAverage/columnAverage.H
new file mode 100644
index 00000000000..00c28ea9400
--- /dev/null
+++ b/src/functionObjects/field/columnAverage/columnAverage.H
@@ -0,0 +1,171 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2018 OpenCFD Ltd.
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 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/>.
+
+Class
+    Foam::functionObjects::columnAverage
+
+Group
+    grpFieldFunctionObjects
+
+Description
+    Averages columns of cells for layered meshes.
+
+    For each patch face, calculates the average value of all cells attached in
+    the patch face normal direction, and then pushes the average value back
+    to all cells in the column.
+
+    Useful for channel-like cases where we want to average fields in the
+    spanwise direction.
+
+Usage
+    Example of function object specification:
+    \verbatim
+    columnAverage1
+    {
+        type        columnAverage;
+        libs        ("libfieldFunctionObjects.so");
+        ...
+        patches     (front side);
+        fields      (U p);
+    }
+    \endverbatim
+
+    Where the entries comprise:
+    \table
+        Property     | Description               | Required    | Default value
+        type         | type name: fieldMinMax    | yes         |
+        patches      | list of patches to collapse onto | yes  |
+        fields       | list of fields to process | yes         |
+    \endtable
+
+See also
+    Foam::functionObjects::fvMeshFunctionObject
+
+SourceFiles
+    columnAverage.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef functionObjects_columnAverage_H
+#define functionObjects_columnAverage_H
+
+#include "volFieldsFwd.H"
+#include "fvMeshFunctionObject.H"
+#include "volFieldSelection.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+class globalIndex;
+class meshStructure;
+
+namespace functionObjects
+{
+
+/*---------------------------------------------------------------------------*\
+                        Class columnAverage Declaration
+\*---------------------------------------------------------------------------*/
+
+class columnAverage
+:
+    public fvMeshFunctionObject
+{
+    // Private data
+
+        //- Patches on which to collapse the fields
+        labelHashSet patchSet_;
+
+        //- Fields to collapse
+        volFieldSelection fieldSet_;
+
+        mutable autoPtr<globalIndex> globalFaces_;
+        mutable autoPtr<globalIndex> globalEdges_;
+        mutable autoPtr<globalIndex> globalPoints_;
+        mutable autoPtr<meshStructure> meshStructurePtr_;
+
+
+    // Private Member Functions
+
+        //- Create the column average field name
+        const word averageName(const word& fieldName) const;
+
+        //- Return the column-based addressing
+        const meshStructure& meshAddressing(const polyMesh&) const;
+
+        //- Calculate the averaged field and return true if successful
+        template<class Type>
+        bool columnAverageField(const word& fieldName);
+
+
+public:
+
+    //- Runtime type information
+    TypeName("columnAverage");
+
+
+    // Constructors
+
+        //- Construct from Time and dictionary
+        columnAverage
+        (
+            const word& name,
+            const Time& runTime,
+            const dictionary&
+        );
+
+
+    //- Destructor
+    virtual ~columnAverage() = default;
+
+
+    // Member Functions
+
+        //- Read the settings
+        virtual bool read(const dictionary&);
+
+        //- Execute, currently does nothing
+        virtual bool execute();
+
+        //- Write the results
+        virtual bool write();
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace functionObjects
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#ifdef NoRepository
+    #include "columnAverageTemplates.C"
+#endif
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/src/functionObjects/field/columnAverage/columnAverageTemplates.C b/src/functionObjects/field/columnAverage/columnAverageTemplates.C
new file mode 100644
index 00000000000..345307d0511
--- /dev/null
+++ b/src/functionObjects/field/columnAverage/columnAverageTemplates.C
@@ -0,0 +1,114 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2018 OpenCFD Ltd.
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 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/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "volFields.H"
+#include "meshStructure.H"
+#include "globalIndex.H"
+
+// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
+
+template<class Type>
+bool Foam::functionObjects::columnAverage::columnAverageField
+(
+    const word& fieldName
+)
+{
+    typedef GeometricField<Type, fvPatchField, volMesh> fieldType;
+
+    const fieldType* fldPtr = lookupObjectPtr<fieldType>(fieldName);
+
+    if (fldPtr)
+    {
+        const fieldType& fld = *fldPtr;
+
+        const word resultName(averageName(fieldName));
+
+        if (!obr_.foundObject<fieldType>(resultName))
+        {
+            fieldType* ptr = new fieldType
+            (
+                IOobject
+                (
+                    resultName,
+                    fld.mesh().time().timeName(),
+                    fld.mesh(),
+                    IOobject::NO_READ,
+                    IOobject::NO_WRITE
+                ),
+                fld
+            );
+            obr_.objectRegistry::store(ptr);
+        }
+
+        fieldType& res = obr_.lookupObjectRef<fieldType>(resultName);
+
+        const meshStructure& ms = meshAddressing(fld.mesh());
+        if (globalFaces_().empty())
+        {
+            return false;
+        }
+
+        const labelList& cellToPatchFace = ms.cellToPatchFaceAddressing();
+
+        // Brute force: collect per-global-patchface on all processors
+        Field<Type> regionField(globalFaces_().size(), Zero);
+        labelList regionCount(globalFaces_().size(), 0);
+
+        forAll(cellToPatchFace, celli)
+        {
+            const label regioni = cellToPatchFace[celli];
+            regionField[regioni] += fld[celli];
+            regionCount[regioni]++;
+        }
+
+        // Global sum
+        Pstream::listCombineGather(regionField, plusEqOp<Type>());
+        Pstream::listCombineScatter(regionField);
+        Pstream::listCombineGather(regionCount, plusEqOp<label>());
+        Pstream::listCombineScatter(regionCount);
+
+        forAll(regionField, regioni)
+        {
+            regionField[regioni] /= regionCount[regioni];
+        }
+
+        // And send result back
+        forAll(cellToPatchFace, celli)
+        {
+            const label regioni = cellToPatchFace[celli];
+            res[celli] = regionField[regioni];
+        }
+        res.correctBoundaryConditions();
+
+        return true;
+    }
+    else
+    {
+        return false;
+    }
+}
+
+
+// ************************************************************************* //
-- 
GitLab