ENH: improvements for nastran surface writer (#1571)

- avoid face copying.
  Maintain separate offsets/list for non tri/quad face decomposition,
  which eliminates copying for tri/quad types that represent the bulk
  of geometries

- report inappropriate use of PLOAD2 for higher-ranks only once per
  field instead of per face.  For this case, write its magnitude
  instead of 0.

- perform field output scaling prior to calling the write face
  function. This will make it easier to handle different per-field
  scaling in the future (#1612)

BUG: nastran quad written as "CTRIA3" instead of "CQUAD4"

src/fileFormats/nastran/NASCore.C

@@ -5,7 +5,7 @@
     \\  /    A nd           | www.openfoam.com
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
-    Copyright (C) 2017-2019 OpenCFD Ltd.
+    Copyright (C) 2017-2020 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -197,4 +197,38 @@ Foam::Ostream& Foam::fileFormats::NASCore::writeKeyword
 }
 
 
+Foam::label Foam::fileFormats::NASCore::faceDecomposition
+(
+    const UList<point>& points,
+    const UList<face>& faces,
+    labelList& decompOffsets,
+    DynamicList<face>& decompFaces
+)
+{
+    // On-demand face decomposition (triangulation)
+
+    decompOffsets.resize(faces.size()+1);
+    decompFaces.clear();
+
+    auto offsetIter = decompOffsets.begin();
+    *offsetIter = 0; // The first offset is always zero
+
+    for (const face& f : faces)
+    {
+        const label n = f.size();
+
+        if (n != 3 && n != 4)
+        {
+            // Decompose non-tri/quad into tris
+            f.triangles(points, decompFaces);
+        }
+
+        // The end offset, which is the next begin offset
+        *(++offsetIter) = decompFaces.size();
+    }
+
+    return decompFaces.size();
+}
+
+
 // ************************************************************************* //

src/fileFormats/nastran/NASCore.H

@@ -5,7 +5,7 @@
     \\  /    A nd           | www.openfoam.com
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
-    Copyright (C) 2017-2019 OpenCFD Ltd.
+    Copyright (C) 2017-2020 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -40,15 +40,15 @@ SourceFiles
 #include "scalar.H"
 #include "string.H"
 #include "Enum.H"
+#include "face.H"
+#include "point.H"
+#include "DynamicList.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
 namespace Foam
 {
 
-// Forward Declarations
-class Ostream;
-
 namespace fileFormats
 {
 
@@ -84,7 +84,7 @@ public:
 
     // Constructors
 
-        //- Construct null
+        //- Default construct
         NASCore() = default;
 
 
@@ -120,6 +120,22 @@ public:
             const word& keyword,
             const fieldFormat format
         );
+
+        //- Calculate face decomposition for non tri/quad faces
+        //
+        //  \param points the surface points
+        //  \param faces  the surface faces
+        //  \param decompOffsets begin/end offsets (size+1) into decompFaces
+        //  \param decompFaces  List of non-tri/quad decomposed into triangles
+        //
+        //  \return number of decomposed faces
+        static label faceDecomposition
+        (
+            const UList<point>& points,
+            const UList<face>& faces,
+            labelList& decompOffsets,
+            DynamicList<face>& decompFaces
+        );
 };
 
 

src/surfMesh/surfaceFormats/nas/NASsurfaceFormat.C

@@ -55,7 +55,7 @@ inline Foam::label Foam::fileFormats::NASsurfaceFormat<Face>::writeShell
     }
     else if (n == 4)
     {
-        os  << "CTRIA3" << ','
+        os  << "CQUAD4" << ','
             << ++elementId << ','
             << (groupId + 1) << ','
             << (f[0] + 1) << ','

src/surfMesh/writers/nastran/nastranSurfaceWriter.C

@@ -204,7 +204,8 @@ void Foam::surfaceWriters::nastranWriter::writeGeometry
 (
     Ostream& os,
     const meshedSurf& surf,
-    List<faceList>& decomposedFaces
+    labelList& decompOffsets,
+    DynamicList<face>& decompFaces
 ) const
 {
     const pointField& points = surf.points();
@@ -213,9 +214,9 @@ void Foam::surfaceWriters::nastranWriter::writeGeometry
 
     // Write points
 
-    os  << "$" << nl
+    os  << '$' << nl
         << "$ Points" << nl
-        << "$" << nl;
+        << '$' << nl;
 
     forAll(points, pointi)
     {
@@ -223,18 +224,19 @@ void Foam::surfaceWriters::nastranWriter::writeGeometry
     }
 
     // Write faces, with on-the-fly decomposition (triangulation)
-    decomposedFaces.clear();
+    decompOffsets.resize(faces.size()+1);
-    decomposedFaces.resize(faces.size());
+    decompFaces.clear();
 
-    os  << "$" << nl
+    decompOffsets[0] = 0; // The first offset is always zero
+
+    os  << '$' << nl
         << "$ Faces" << nl
-        << "$" << nl;
+        << '$' << nl;
 
-    label elemId = 0; // The element-id
+    label elemId = 0;  // The element-id
     forAll(faces, facei)
     {
         const face& f = faces[facei];
-        faceList& decomp = decomposedFaces[facei];
 
         // 1-offset for PID
         const label propId = 1 + (facei < zones.size() ? zones[facei] : 0);
@@ -242,28 +244,36 @@ void Foam::surfaceWriters::nastranWriter::writeGeometry
         if (f.size() == 3)
         {
             writeFace(os, "CTRIA3", f, ++elemId, propId);
-            decomp.resize(1);
-            decomp[0] = f;
         }
         else if (f.size() == 4)
         {
             writeFace(os, "CQUAD4", f, ++elemId, propId);
-            decomp.resize(1);
-            decomp[0] = f;
         }
         else
         {
-            // Decompose poly face into tris
+            // Decompose into tris
-            decomp.resize(f.nTriangles());
+            f.triangles(points, decompFaces);
 
-            label nTri = 0;
+            for
-            f.triangles(points, nTri, decomp);
+            (
-
+                label decompi = decompOffsets[facei];
-            for (const face& f2 : decomp)
+                decompi < decompFaces.size();
+                ++decompi
+            )
             {
-                writeFace(os, "CTRIA3", f2, ++elemId, propId);
+                writeFace
+                (
+                    os,
+                    "CTRIA3",
+                    decompFaces[decompi],
+                    ++elemId,
+                    propId
+                );
             }
         }
+
+        // The end offset, which is the next begin offset
+        decompOffsets[facei+1] = decompFaces.size();
     }
 }
 
@@ -323,7 +333,7 @@ Foam::surfaceWriters::nastranWriter::nastranWriter()
     surfaceWriter(),
     writeFormat_(fieldFormat::SHORT),
     fieldMap_(),
-    scale_(1.0),
+    scale_(1),
     separator_()
 {}
 
@@ -432,11 +442,13 @@ Foam::fileName Foam::surfaceWriters::nastranWriter::write()
 
         os  << "TITLE=OpenFOAM " << outputPath_.name()
             << " mesh" << nl
-            << "$" << nl
+            << '$' << nl
             << "BEGIN BULK" << nl;
 
-        List<faceList> decomposedFaces;
+        labelList decompOffsets;
-        writeGeometry(os, surf, decomposedFaces);
+        DynamicList<face> decompFaces;
+
+        writeGeometry(os, surf, decompOffsets, decompFaces);
 
         writeFooter(os, surf)
             << "ENDDATA" << nl;

src/surfMesh/writers/nastran/nastranSurfaceWriter.H

@@ -155,12 +155,16 @@ private:
             const label propId      //!< 1-based Property Id
         ) const;
 
-        //- Main driver to write the surface mesh geometry
+        //- Write the surface mesh geometry, tracking face decomposition
+        //
+        //  \param decompOffsets  begin/end offsets (size+1) into decompFaces
+        //  \param decompFaces    Non tri/quad decomposed into triangles
         void writeGeometry
         (
             Ostream& os,
             const meshedSurf& surf,
-            List<faceList>& decomposedFaces
+            labelList& decompOffsets,
+            DynamicList<face>& decompFaces
         ) const;
 
         //- Write the formatted keyword to the output stream

src/surfMesh/writers/nastran/nastranSurfaceWriterImpl.C

@@ -89,8 +89,6 @@ Foam::Ostream& Foam::surfaceWriters::nastranWriter::writeFaceValue
 
     const label setId = 1;
 
-    Type scaledValue = scale_*value;
-
     // Write keyword
     writeKeyword(os, fileFormats::NASCore::loadFormatNames[format])
         << separator_;
@@ -106,19 +104,14 @@ Foam::Ostream& Foam::surfaceWriters::nastranWriter::writeFaceValue
         {
             if (pTraits<Type>::nComponents == 1)
             {
-                writeValue(os, scaledValue) << separator_;
+                writeValue(os, value);
             }
             else
             {
-                WarningInFunction
+                writeValue(os, mag(value));
-                    << fileFormats::NASCore::loadFormatNames[format]
-                    << " requires scalar values"
-                    << " - it cannot be used for higher rank values"
-                    << endl;
-
-                writeValue(os, scalar(0)) << separator_;
             }
 
+            os << separator_;
             writeValue(os, elemId);
             break;
         }
@@ -130,7 +123,7 @@ Foam::Ostream& Foam::surfaceWriters::nastranWriter::writeFaceValue
             for (direction d = 0; d < pTraits<Type>::nComponents; ++d)
             {
                 os  << separator_;
-                writeValue(os, component(scaledValue, d));
+                writeValue(os, component(value, d));
             }
             break;
         }
@@ -164,7 +157,7 @@ Foam::fileName Foam::surfaceWriters::nastranWriter::writeTemplate
         return fileName::null;
     }
 
-    const loadFormat& format(fieldMap_[fieldName]);
+    const loadFormat format(fieldMap_[fieldName]);
 
     // Field:  rootdir/<TIME>/field/surfaceName.nas
 
@@ -177,9 +170,27 @@ Foam::fileName Foam::surfaceWriters::nastranWriter::writeTemplate
     outputFile /= fieldName / outputPath_.name();
     outputFile.ext("nas");
 
+
+    // Currently the same scaling for all variables
+    const scalar varScale = scale_;
+
     if (verbose_)
     {
-        Info<< "Writing field " << fieldName << " to " << outputFile << endl;
+        Info<< "Writing field " << fieldName;
+        if (!equal(varScale, 1))
+        {
+            Info<< " (scaling " << varScale << ')';
+        }
+        Info<< " to " << outputFile << endl;
+    }
+
+    // Emit any common warnings
+    if (format == loadFormat::PLOAD2 && pTraits<Type>::nComponents != 1)
+    {
+        WarningInFunction
+            << fileFormats::NASCore::loadFormatNames[format]
+            << " cannot be used for higher rank values"
+            << " - reverting to mag()" << endl;
     }
 
 
@@ -199,43 +210,84 @@ Foam::fileName Foam::surfaceWriters::nastranWriter::writeTemplate
 
         const scalar timeValue(0);
 
+        // Additional bookkeeping for decomposing non tri/quad
+        labelList decompOffsets;
+        DynamicList<face> decompFaces;
+
+
+        // Could handle separate geometry here
+
         OFstream os(outputFile);
         fileFormats::NASCore::setPrecision(os, writeFormat_);
 
-        if (verbose_)
+        os  << "TITLE=OpenFOAM " << outputFile.name()
+            << token::SPACE << fieldName << " data" << nl;
+
+        if (useTimeDir() && !timeName().empty())
         {
-            Info<< "Writing nastran file to " << os.name() << endl;
+            os  << '$' << nl
+                << "$ TIME " << timeName() << nl;
         }
 
-        os  << "TITLE=OpenFOAM " << outputFile.name()
+        os  << '$' << nl
-            << " " << fieldName << " data" << nl
-            << "$" << nl
             << "TIME " << timeValue << nl
-            << "$" << nl
+            << '$' << nl
             << "BEGIN BULK" << nl;
 
-        List<faceList> decomposedFaces;
-        writeGeometry(os, surf, decomposedFaces);
 
-        os  << "$" << nl
+        writeGeometry(os, surf, decompOffsets, decompFaces);
+
+
+        // Write field
+
+        os  << '$' << nl
             << "$ Field data" << nl
-            << "$" << nl;
+            << '$' << nl;
+
+
+        // Regular (undecomposed) faces
+        const faceList& faces = surf.faces();
 
         label elemId = 0;
 
         if (this->isPointData())
         {
-            for (const faceList& dFaces : decomposedFaces)
+            forAll(faces, facei)
             {
-                for (const face& f : dFaces)
+                const label beginElemId = elemId;
-                {
-                    Type v = Zero;
 
+                // Any face decomposition
+                for
+                (
+                    label decompi = decompOffsets[facei];
+                    decompi < decompOffsets[facei+1];
+                    ++decompi
+                )
+                {
+                    const face& f = decompFaces[decompi];
+
+                    Type v = Zero;
                     for (const label verti : f)
                     {
                         v += values[verti];
                     }
-                    v /= f.size();
+                    v *= (varScale / f.size());
+
+                    writeFaceValue(os, format, v, ++elemId);
+                }
+
+
+                // Face not decomposed
+                if (beginElemId == elemId)
+                {
+                    const face& f = faces[facei];
+
+                    Type v = Zero;
+                    for (const label verti : f)
+                    {
+                        v += values[verti];
+                    }
+                    v *= (varScale / f.size());
 
                     writeFaceValue(os, format, v, ++elemId);
                 }
@@ -245,13 +297,22 @@ Foam::fileName Foam::surfaceWriters::nastranWriter::writeTemplate
         {
             auto valIter = values.cbegin();
 
-            for (const faceList& dFaces : decomposedFaces)
+            forAll(faces, facei)
             {
-                forAll(dFaces, facei)
+                const Type v(varScale * *valIter);
-                {
-                    writeFaceValue(os, format, *valIter, ++elemId);
-                }
                 ++valIter;
+
+                label nValues =
+                    max
+                    (
+                        label(1),
+                        (decompOffsets[facei+1] - decompOffsets[facei])
+                    );
+
+                while (nValues--)
+                {
+                    writeFaceValue(os, format, v, ++elemId);
+                }
             }
         }