ENH: added areaWrite function object (#1237)

- write finiteArea meshes and fields to standard surface output
  formats (Ensight, VTK, etc).

modules/avalanche

@@ -1 +1 @@
-Subproject commit fdf0f40d7a24c7b08c255ba2ee05ecbb847cc264
+Subproject commit 3e147ba27a59e6dc008c7d32cbebfd29f58b5a36

src/functionObjects/utilities/Make/files

@@ -2,6 +2,8 @@ abort/abort.C
 
 codedFunctionObject/codedFunctionObject.C
 
+areaWrite/areaWrite.C
+
 ensightWrite/ensightWrite.C
 ensightWrite/ensightWriteUpdate.C
 

src/functionObjects/utilities/Make/options

@@ -5,6 +5,7 @@ EXE_INC = \
     -I$(LIB_SRC)/meshTools/lnInclude \
     -I$(LIB_SRC)/conversion/lnInclude \
     -I$(LIB_SRC)/sampling/lnInclude \
+    -I$(LIB_SRC)/surfMesh/lnInclude \
     -I$(LIB_SRC)/dynamicMesh/lnInclude \
     -I$(LIB_SRC)/ODE/lnInclude \
     -I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \
@@ -16,6 +17,7 @@ LIB_LIBS = \
     -lmeshTools \
     -lconversion \
     -lsampling \
+    -lsurfMesh \
     -ldynamicMesh \
     -lfluidThermophysicalModels \
     -lcompressibleTransportModels \

src/functionObjects/utilities/areaWrite/areaWrite.C

@@ -0,0 +1,396 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2019 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 "areaWrite.H"
+#include "polySurface.H"
+
+#include "fvMesh.H"
+#include "mapPolyMesh.H"
+#include "areaFields.H"
+#include "HashOps.H"
+#include "Time.H"
+#include "UIndirectList.H"
+#include "addToRunTimeSelectionTable.H"
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+namespace Foam
+{
+    defineTypeNameAndDebug(areaWrite, 0);
+
+    addToRunTimeSelectionTable
+    (
+        functionObject,
+        areaWrite,
+        dictionary
+    );
+}
+
+Foam::scalar Foam::areaWrite::mergeTol_ = 1e-10;
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+Foam::areaWrite::areaWrite
+(
+    const word& name,
+    const Time& runTime,
+    const dictionary& dict
+)
+:
+    functionObjects::fvMeshFunctionObject(name, runTime, dict),
+    loadFromFiles_(false),
+    verbose_(false),
+    outputPath_
+    (
+        time_.globalPath()/functionObject::outputPrefix/name
+    ),
+    selectAreas_(),
+    fieldSelection_(),
+    meshes_(),
+    surfaces_(nullptr),
+    writers_()
+{
+    outputPath_.clean();  // Remove unneeded ".."
+
+    read(dict);
+}
+
+
+Foam::areaWrite::areaWrite
+(
+    const word& name,
+    const objectRegistry& obr,
+    const dictionary& dict,
+    const bool loadFromFiles
+)
+:
+    functionObjects::fvMeshFunctionObject(name, obr, dict),
+    loadFromFiles_(loadFromFiles),
+    verbose_(false),
+    outputPath_
+    (
+        time_.globalPath()/functionObject::outputPrefix/name
+    ),
+    selectAreas_(),
+    fieldSelection_(),
+    meshes_(),
+    surfaces_(nullptr),
+    writers_()
+{
+    outputPath_.clean();  // Remove unneeded ".."
+
+    read(dict);
+}
+
+
+// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+
+void Foam::areaWrite::verbose(const bool verbosity)
+{
+    verbose_ = verbosity;
+}
+
+
+bool Foam::areaWrite::read(const dictionary& dict)
+{
+    fvMeshFunctionObject::read(dict);
+
+    writers_.clear();
+    selectAreas_.clear();
+    fieldSelection_.clear();
+
+    surfaces_.reset
+    (
+        new objectRegistry
+        (
+            IOobject
+            (
+                "::areaWrite::",
+                obr_.time().constant(),
+                obr_,
+                IOobject::NO_READ,
+                IOobject::NO_WRITE,
+                false
+            )
+        )
+    );
+
+    verbose_ = dict.lookupOrDefault("verbose", false);
+
+    // All possible area meshes for the given fvMesh region
+    meshes_ = obr().lookupClass<faMesh>();
+
+    dict.readIfPresent("areas", selectAreas_);
+
+    if (selectAreas_.empty())
+    {
+        word areaName;
+        if (!dict.readIfPresent("area", areaName))
+        {
+            wordList available = obr().sortedNames<faMesh>();
+
+            if (available.size())
+            {
+                areaName = available.first();
+            }
+        }
+
+        if (!areaName.empty())
+        {
+            selectAreas_.resize(1);
+            selectAreas_.first() = areaName;
+        }
+    }
+
+    // Restrict to specified meshes
+    meshes_.filterKeys(selectAreas_);
+
+    dict.readEntry("fields", fieldSelection_);
+    fieldSelection_.uniq();
+
+
+    // Surface writer type and format options
+    const word writerType = dict.get<word>("surfaceFormat");
+
+    const dictionary writerOptions
+    (
+        dict.subOrEmptyDict("formatOptions").subOrEmptyDict(writerType)
+    );
+
+    for (const word& areaName : meshes_.keys())
+    {
+        // Define surface writer, but do NOT yet attach a surface
+
+        auto surfWriter = surfaceWriter::New(writerType, writerOptions);
+
+        // Use outputDir/TIME/surface-name
+        surfWriter->useTimeDir() = true;
+        surfWriter->verbose() = verbose_;
+
+        writers_.set(areaName, surfWriter);
+    }
+
+    // Ensure all surfaces and merge information are expired
+    expire();
+
+    return true;
+}
+
+
+bool Foam::areaWrite::execute()
+{
+    return true;
+}
+
+
+bool Foam::areaWrite::write()
+{
+    // Just needed for warnings
+    wordList allFields;
+    HashTable<wordHashSet> selected;
+    DynamicList<label> missed(fieldSelection_.size());
+
+
+    for (const word& areaName : meshes_.sortedToc())
+    {
+        const faMesh& areaMesh = *meshes_[areaName];
+
+        polySurface* surfptr = surfaces_->getObjectPtr<polySurface>(areaName);
+
+        if (!surfptr)
+        {
+            // Construct null and add to registry (owned by registry)
+            surfptr = new polySurface(areaName, *surfaces_, true);
+        }
+
+        pointField pts(areaMesh.patch().localPoints());
+        faceList fcs(areaMesh.patch().localFaces());
+
+        // Copy in geometry
+        surfptr->transfer(std::move(pts), std::move(fcs));
+
+        surfaceWriter& outWriter = *writers_[areaName];
+
+        if (outWriter.needsUpdate())
+        {
+            outWriter.setSurface(*surfptr);
+        }
+
+
+        // Determine the per-surface number of fields
+        // Only seems to be needed for VTK legacy
+
+        selected.clear();
+
+        const IOobjectList objects(areaMesh.thisDb(), obr_.time().timeName());
+
+        if (loadFromFiles_)
+        {
+            allFields = objects.names();
+            selected = objects.classes(fieldSelection_);
+        }
+        else
+        {
+            // Check currently available fields
+            allFields = areaMesh.thisDb().names();
+            selected = areaMesh.thisDb().classes(fieldSelection_);
+        }
+
+        if (Pstream::parRun())
+        {
+            Pstream::mapCombineGather(selected, HashSetOps::plusEqOp<word>());
+            Pstream::mapCombineScatter(selected);
+        }
+
+        missed.clear();
+
+        // Detect missing fields
+        forAll(fieldSelection_, i)
+        {
+            if (findStrings(fieldSelection_[i], allFields).empty())
+            {
+                missed.append(i);
+            }
+        }
+
+        if (missed.size())
+        {
+            WarningInFunction
+                << nl
+                << "Cannot find "
+                << (loadFromFiles_ ? "field file" : "registered field")
+                << " matching "
+                << UIndirectList<wordRe>(fieldSelection_, missed) << endl;
+        }
+
+
+        // Currently only support area field types
+        label nAreaFields = 0;
+
+        forAllConstIters(selected, iter)
+        {
+            const word& clsName = iter.key();
+            const label n = iter.val().size();
+
+            if (fieldTypes::area.found(clsName))
+            {
+                nAreaFields += n;
+            }
+        }
+
+
+        // Propagate field counts (per surface)
+        outWriter.nFields() = nAreaFields;
+
+
+        // Begin writing
+
+        outWriter.open(outputPath_/areaName);
+
+        outWriter.beginTime(obr_.time());
+
+        // Write fields
+
+        performAction<areaScalarField>(outWriter, areaMesh, objects);
+        performAction<areaVectorField>(outWriter, areaMesh, objects);
+        performAction<areaSphericalTensorField>(outWriter, areaMesh, objects);
+        performAction<areaSymmTensorField>(outWriter, areaMesh, objects);
+        performAction<areaTensorField>(outWriter, areaMesh, objects);
+
+        // Finish this time step
+
+        // Write geometry if no fields were written so that we still
+        // can have something to look at
+
+        if (!outWriter.wroteData())
+        {
+            outWriter.write();
+        }
+
+        outWriter.endTime();
+    }
+
+    return true;
+}
+
+
+void Foam::areaWrite::expire()
+{
+    surfaces_->clear();
+
+    // Dimension as fraction of mesh bounding box
+    const scalar mergeDim = mergeTol_ * mesh_.bounds().mag();
+
+    forAllIters(writers_, iter)
+    {
+        surfaceWriter& writer = *(iter.val());
+        writer.expire();
+        writer.mergeDim() = mergeDim;
+    }
+}
+
+
+void Foam::areaWrite::updateMesh(const mapPolyMesh& mpm)
+{
+    if (&mpm.mesh() == &mesh_)
+    {
+        expire();
+    }
+}
+
+
+void Foam::areaWrite::movePoints(const polyMesh& mesh)
+{
+    if (&mesh == &mesh_)
+    {
+        expire();
+    }
+}
+
+
+void Foam::areaWrite::readUpdate(const polyMesh::readUpdateState state)
+{
+    if (state != polyMesh::UNCHANGED)
+    {
+        expire();
+    }
+}
+
+
+Foam::scalar Foam::areaWrite::mergeTol()
+{
+    return mergeTol_;
+}
+
+
+Foam::scalar Foam::areaWrite::mergeTol(const scalar tol)
+{
+    const scalar prev(mergeTol_);
+    mergeTol_ = tol;
+    return prev;
+}
+
+
+// ************************************************************************* //

src/functionObjects/utilities/areaWrite/areaWrite.H

@@ -0,0 +1,249 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2019 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::areaWrite
+
+Description
+    Write finite area mesh/fields to standard output formats.
+
+    Example of function object specification:
+
+    \verbatim
+    surfaces
+    {
+        type    areaWrite;
+        libs    ("libutilityFunctionObjects.so");
+
+        // Write at same frequency as fields
+        writeControl    outputTime;
+        writeInterval   1;
+
+        // Fields to be sampled
+        fields          (p U);
+
+        // Output surface format
+        surfaceFormat   vtk;
+
+        formatOptions
+        {
+            vtk
+            {
+                precision  10;
+            }
+        }
+    }
+    \endverbatim
+
+    Entries:
+    \table
+        Property | Description                              | Required | Default
+        type    | Type-name: \c areaWrite                   | yes |
+        region  | name for a single region                  | no  | region0
+        area    | select a single area                      | no  |
+        areas   | wordRe list of multiple areas             | no  |
+        fields  | wordRe list of fields                     | yes |
+        surfaceFormat | output surface format               | yes |
+        formatOptions | dictionary of format options        | no  |
+    \endtable
+
+SourceFiles
+    areaWrite.C
+    areaWriteTemplates.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef areaWrite_H
+#define areaWrite_H
+
+#include "fvMeshFunctionObject.H"
+#include "polyMesh.H"
+#include "areaFieldsFwd.H"
+#include "surfaceWriter.H"
+#include "HashPtrTable.H"
+#include "IOobjectList.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+// Forward Declarations
+class faMesh;
+class polySurface;
+
+/*---------------------------------------------------------------------------*\
+                          Class areaWrite Declaration
+\*---------------------------------------------------------------------------*/
+
+class areaWrite
+:
+    public functionObjects::fvMeshFunctionObject
+{
+    // Static Data Members
+
+        //- Tolerance for merging points (fraction of mesh bounding box)
+        static scalar mergeTol_;
+
+
+    // Private Data
+
+        //- Load fields from files (not from objectRegistry)
+        const bool loadFromFiles_;
+
+        //- Output verbosity
+        bool verbose_;
+
+        //- Output path
+        fileName outputPath_;
+
+
+    // Read from dictionary
+
+        //- Names of areas to select
+        wordRes selectAreas_;
+
+        //- Names of fields to write
+        wordRes fieldSelection_;
+
+        //- Pointers to the requested mesh regions
+        HashTable<const faMesh*> meshes_;
+
+        //- Hold intermediate surfaces.
+        //  The faMesh has an indirect face list but we require real ones.
+        autoPtr<objectRegistry> surfaces_;
+
+
+    // Output control
+
+        //- Surface writers (one per surface)
+        HashPtrTable<surfaceWriter> writers_;
+
+
+    // Private Member Functions
+
+        //- Write fieldName on surface and on outputDir path
+        template<class Type>
+        void writeSurface
+        (
+            surfaceWriter& writer,
+            const Field<Type>& values,
+            const word& fieldName
+        );
+
+        //- Write all applicable fields
+        template<class GeoField>
+        void performAction
+        (
+            surfaceWriter& writer,
+            const faMesh& areaMesh,
+            const IOobjectList& objects
+        );
+
+        //- Mark intermediate surfaces and writers as needing an update.
+        void expire();
+
+        //- No copy construct
+        areaWrite(const areaWrite&) = delete;
+
+        //- No copy assignment
+        void operator=(const areaWrite&) = delete;
+
+
+public:
+
+    //- Runtime type information
+    TypeName("areaWrite");
+
+
+    // Constructors
+
+        //- Construct from Time and dictionary
+        areaWrite
+        (
+            const word& name,
+            const Time& runTime,
+            const dictionary& dict
+        );
+
+        //- Construct for given objectRegistry and dictionary
+        //  allow the possibility to load fields from files
+        areaWrite
+        (
+            const word& name,
+            const objectRegistry& obr,
+            const dictionary& dict,
+            const bool loadFromFiles = false
+        );
+
+
+    //- Destructor
+    virtual ~areaWrite() = default;
+
+
+    // Member Functions
+
+        //- Set verbosity level
+        void verbose(const bool verbosity = true);
+
+        //- Read the areaWrite dictionary
+        virtual bool read(const dictionary& dict);
+
+        //- Execute, currently does nothing
+        virtual bool execute();
+
+        //- Sample and write
+        virtual bool write();
+
+        //- Update for changes of mesh - expires the surfaces
+        virtual void updateMesh(const mapPolyMesh& mpm);
+
+        //- Update for mesh point-motion - expires the surfaces
+        virtual void movePoints(const polyMesh& mesh);
+
+        //- Update for changes of mesh due to readUpdate - expires the surfaces
+        virtual void readUpdate(const polyMesh::readUpdateState state);
+
+        //- Get merge tolerance
+        static scalar mergeTol();
+
+        //- Set merge tolerance and return old value
+        static scalar mergeTol(const scalar tol);
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#ifdef NoRepository
+    #include "areaWriteTemplates.C"
+#endif
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

src/functionObjects/utilities/areaWrite/areaWriteTemplates.C

@@ -0,0 +1,108 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2019 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 "areaWrite.H"
+#include "areaFields.H"
+#include "faMesh.H"
+
+// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
+
+template<class Type>
+void Foam::areaWrite::writeSurface
+(
+    surfaceWriter& writer,
+    const Field<Type>& values,
+    const word& fieldName
+)
+{
+    fileName outputName = writer.write(fieldName, values);
+
+    // Case-local file name with "<case>" to make relocatable
+
+    dictionary propsDict;
+    propsDict.add
+    (
+        "file",
+        time_.relativePath(outputName, true)
+    );
+    setProperty(fieldName, propsDict);
+}
+
+
+template<class GeoField>
+void Foam::areaWrite::performAction
+(
+    surfaceWriter& writer,
+    const faMesh& areaMesh,
+    const IOobjectList& objects
+)
+{
+    wordList fieldNames;
+    if (loadFromFiles_)
+    {
+        fieldNames = objects.sortedNames<GeoField>(fieldSelection_);
+    }
+    else
+    {
+        fieldNames = areaMesh.thisDb().sortedNames<GeoField>(fieldSelection_);
+    }
+
+    for (const word& fieldName : fieldNames)
+    {
+        if (verbose_)
+        {
+            Info<< "write: " << fieldName << endl;
+        }
+
+        if (loadFromFiles_)
+        {
+            const GeoField fld
+            (
+                IOobject
+                (
+                    fieldName,
+                    time_.timeName(),
+                    areaMesh.thisDb(),
+                    IOobject::MUST_READ
+                ),
+                areaMesh
+            );
+
+            writeSurface(writer, fld, fieldName);
+        }
+        else
+        {
+            writeSurface
+            (
+                writer,
+                areaMesh.thisDb().lookupObject<GeoField>(fieldName),
+                fieldName
+            );
+        }
+    }
+}
+
+
+// ************************************************************************* //

tutorials/finiteArea/surfactantFoam/planeTransport/system/areaWrite

@@ -0,0 +1,29 @@
+// -*- C++ -*-
+// Use the areaWrite function object
+
+areaWrite
+{
+    type    areaWrite;
+    libs    ("libutilityFunctionObjects.so");
+    log     true;
+
+    writeControl    writeTime;
+    writeInterval   1;
+
+    // Fields to output (words or regex)
+    fields  (".*");
+
+    surfaceFormat ensight;
+    // surfaceFormat vtk;
+
+    formatOptions
+    {
+        vtk
+        {
+            format  ascii;
+        }
+    }
+}
+
+
+// ************************************************************************* //

tutorials/finiteArea/surfactantFoam/planeTransport/system/controlDict

@@ -43,6 +43,9 @@ timePrecision       6;
 
 runTimeModifiable   yes;
 
+functions
+{
+    #include "areaWrite"
+}
 
 // ************************************************************************* //
-