CONFIG: inject -no-recursion into the argument list (#3198)

- sourcing a file with '-no-recursion "$@"' does not work with dash.
  Need to modify the argument list directly.

Allwmake-modules

@@ -1,13 +1,15 @@
 #!/bin/sh
+cd "${0%/*}" || exit            # Run from this directory
+set -- -no-recursion "$@"       # Parse arguments only
+
 # Run from OPENFOAM top-level directory only
-cd "${0%/*}" || exit
 wmake -check-dir "$WM_PROJECT_DIR" 2>/dev/null || {
     echo "Error (${0##*/}) : not located in \$WM_PROJECT_DIR"
     echo "    Check your OpenFOAM environment and installation"
     exit 1
 }
 if [ -f "$WM_PROJECT_DIR"/wmake/scripts/AllwmakeParseArguments ]
-then  . "$WM_PROJECT_DIR"/wmake/scripts/AllwmakeParseArguments -no-recursion "$@" || \
+then  . "$WM_PROJECT_DIR"/wmake/scripts/AllwmakeParseArguments || \
     echo "Argument parse error"
 else
     echo "Error (${0##*/}) : WM_PROJECT_DIR appears to be incorrect"

Allwmake-plugins

@@ -1,13 +1,15 @@
 #!/bin/sh
+cd "${0%/*}" || exit            # Run from this directory
+set -- -no-recursion "$@"       # Parse arguments only
+
 # Run from OPENFOAM top-level directory only
-cd "${0%/*}" || exit
 wmake -check-dir "$WM_PROJECT_DIR" 2>/dev/null || {
     echo "Error (${0##*/}) : not located in \$WM_PROJECT_DIR"
     echo "    Check your OpenFOAM environment and installation"
     exit 1
 }
 if [ -f "$WM_PROJECT_DIR"/wmake/scripts/AllwmakeParseArguments ]
-then  . "$WM_PROJECT_DIR"/wmake/scripts/AllwmakeParseArguments -no-recursion "$@" || \
+then  . "$WM_PROJECT_DIR"/wmake/scripts/AllwmakeParseArguments || \
     echo "Argument parse error"
 else
     echo "Error (${0##*/}) : WM_PROJECT_DIR appears to be incorrect"

META-INFO/api-info

@@ -1,2 +1,2 @@
 api=2406
-patch=0
+patch=241212

applications/solvers/heatTransfer/chtMultiRegionFoam/chtMultiRegionTwoPhaseEulerFoam/fluid/setRegionFluidFields.H

@@ -8,13 +8,11 @@
 
     volVectorField& U1 = phase1.URef();
     surfaceScalarField& phi1 = phase1.phiRef();
-    const tmp<surfaceScalarField> talphaPhi1 = phase1.alphaPhi();
-    const auto& alphaPhi1 = talphaPhi1();
+    const surfaceScalarField& alphaPhi1 = phase1.alphaPhi();
 
     volVectorField& U2 = phase2.URef();
     surfaceScalarField& phi2 = phase2.phiRef();
-    const tmp<surfaceScalarField> talphaPhi2 = phase2.alphaPhi();
-    const auto& alphaPhi2 = talphaPhi2();
+    const surfaceScalarField& alphaPhi2 = phase2.alphaPhi();
 
     surfaceScalarField& phi = fluid.phi();
 

applications/solvers/multiphase/compressibleInterFoam/VoFphaseCompressibleTurbulenceModels/compressibleInterPhaseTransportModel.C

@@ -74,11 +74,8 @@ Foam::compressibleInterPhaseTransportModel::compressibleInterPhaseTransportModel
         const volScalarField& alpha1(mixture_.alpha1());
         const volScalarField& alpha2(mixture_.alpha2());
 
-        const tmp<volScalarField> trho1 = mixture_.thermo1().rho();
-        const tmp<volScalarField> trho2 = mixture_.thermo2().rho();
-
-        const auto& rho1 = trho1();
-        const auto& rho2 = trho2();
+        const volScalarField& rho1 = mixture_.thermo1().rho();
+        const volScalarField& rho2 = mixture_.thermo2().rho();
 
         alphaRhoPhi1_ =
         (
@@ -188,8 +185,8 @@ void Foam::compressibleInterPhaseTransportModel::correctPhasePhi()
 {
     if (twoPhaseTransport_)
     {
-        const tmp<volScalarField> rho1 = mixture_.thermo1().rho();
-        const tmp<volScalarField> rho2 = mixture_.thermo2().rho();
+        const volScalarField& rho1 = mixture_.thermo1().rho();
+        const volScalarField& rho2 = mixture_.thermo2().rho();
 
         alphaRhoPhi1_.ref() = fvc::interpolate(rho1)*alphaPhi10_;
         alphaRhoPhi2_.ref() = fvc::interpolate(rho2)*(phi_ - alphaPhi10_);

applications/solvers/multiphase/reactingMultiphaseEulerFoam/EEqns.H

@@ -12,11 +12,8 @@ for (int Ecorr=0; Ecorr<nEnergyCorrectors; Ecorr++)
         phaseModel& phase = fluid.anisothermalPhases()[anisothermalPhasei];
 
         const volScalarField& alpha = phase;
-        const tmp<volScalarField> trho = phase.rho();
-        const tmp<volVectorField> tU = phase.U();
-
-        const auto& rho = trho();
-        const auto& U = tU();
+        const volScalarField& rho = phase.rho();
+        const volVectorField& U = phase.U();
 
         fvScalarMatrix EEqn
         (

applications/solvers/multiphase/reactingMultiphaseEulerFoam/YEqns.H

@@ -11,9 +11,7 @@
 
         UPtrList<volScalarField>& Y = phase.YActiveRef();
         const volScalarField& alpha = phase;
-        const tmp<volScalarField> trho = phase.rho();
-
-        const auto& rho = trho();
+        const volScalarField& rho = phase.rho();
 
         forAll(Y, i)
         {

applications/solvers/multiphase/reactingMultiphaseEulerFoam/pU/UEqns.H

@@ -14,11 +14,9 @@ PtrList<fvVectorMatrix> UEqns(phases.size());
         phaseModel& phase = fluid.movingPhases()[movingPhasei];
 
         const volScalarField& alpha = phase;
-        const tmp<volScalarField> trho = phase.rho();
+        const volScalarField& rho = phase.rho();
         volVectorField& U = phase.URef();
 
-        const auto& rho = trho();
-
         UEqns.set
         (
             phase.index(),

applications/solvers/multiphase/reactingMultiphaseEulerFoam/pUf/UEqns.H

@@ -17,11 +17,9 @@ PtrList<fvVectorMatrix> UEqns(phases.size());
         phaseModel& phase = fluid.movingPhases()[movingPhasei];
 
         const volScalarField& alpha = phase;
-        const tmp<volScalarField> trho = phase.rho();
+        const volScalarField& rho = phase.rho();
         volVectorField& U = phase.URef();
 
-        const auto& rho = trho();
-
         UEqns.set
         (
             phase.index(),

applications/solvers/multiphase/reactingTwoPhaseEulerFoam/createFieldRefs.H

@@ -6,13 +6,11 @@ const volScalarField& alpha2 = phase2;
 
 volVectorField& U1 = phase1.URef();
 surfaceScalarField& phi1 = phase1.phiRef();
-const tmp<surfaceScalarField> talphaPhi1 = phase1.alphaPhi();
-const auto& alphaPhi1 = talphaPhi1();
+const surfaceScalarField& alphaPhi1 = phase1.alphaPhi();
 
 volVectorField& U2 = phase2.URef();
 surfaceScalarField& phi2 = phase2.phiRef();
-const tmp<surfaceScalarField> talphaPhi2 = phase2.alphaPhi();
-const auto& alphaPhi2 = talphaPhi2();
+const surfaceScalarField& alphaPhi2 = phase2.alphaPhi();
 
 surfaceScalarField& phi = fluid.phi();
 

applications/test/00-dummy/library/Allwmake-mpi

@@ -1,6 +1,7 @@
 #!/bin/sh
-cd "${0%/*}" || exit                                # Run from this directory
-. ${WM_PROJECT_DIR:?}/wmake/scripts/AllwmakeParseArguments -no-recursion "$@"
+cd "${0%/*}" || exit            # Run from this directory
+set -- -no-recursion "$@"       # Parse arguments only
+. ${WM_PROJECT_DIR:?}/wmake/scripts/AllwmakeParseArguments
 . ${WM_PROJECT_DIR:?}/wmake/scripts/wmakeFunctions  # Require wmake functions
 
 #------------------------------------------------------------------------------

applications/utilities/mesh/conversion/vtkUnstructuredToFoam/Make/options

@@ -1,7 +1,5 @@
 EXE_INC = \
-    -I$(LIB_SRC)/finiteVolume/lnInclude \
     -I$(LIB_SRC)/fileFormats/lnInclude
 
 EXE_LIBS = \
-    -lfiniteVolume \
     -lfileFormats

applications/utilities/mesh/conversion/vtkUnstructuredToFoam/vtkUnstructuredToFoam.C

@@ -34,84 +34,23 @@ Description
     Convert legacy VTK file (ascii) containing an unstructured grid
     to an OpenFOAM mesh without boundary information.
 
-Usage
-    \b vtkUnstructuredToFoam \<XXX.vtk\>
-
-    Options:
-      - \par -no-fields
-        Do not attempt to recreate volFields
-
 Note
     The .vtk format does not contain any boundary information.
-    It is purely a description of the internal mesh. This also limits the
-    usefulness of reconstructing the volFields.
-
+    It is purely a description of the internal mesh.
     Not extensively tested.
 
 \*---------------------------------------------------------------------------*/
 
 #include "argList.H"
 #include "Time.H"
-#include "fvMesh.H"
+#include "polyMesh.H"
 #include "IFstream.H"
 #include "vtkUnstructuredReader.H"
 
-#include "columnFvMesh.H"
-#include "scalarIOField.H"
-#include "vectorIOField.H"
-#include "volFields.H"
-
 using namespace Foam;
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
-template<class Type>
-void constructVolFields(fvMesh& mesh, const vtkUnstructuredReader& reader)
-{
-    const auto fields(reader.cellData().csorted<IOField<Type>>());
-    for (const auto& field : fields)
-    {
-        Info<< "Constructing volField " << field.name() << endl;
-
-        // field is
-        //  - cell data followed by
-        //  - boundary face data
-
-
-        auto tfld = GeometricField<Type, fvPatchField, volMesh>::New
-        (
-            field.name(),
-            mesh,
-            dimless
-        );
-        auto& fld = tfld.ref();
-        fld.instance() = mesh.time().timeName();
-        fld.writeOpt() = IOobject::AUTO_WRITE;
-
-        // Fill cell values
-        fld.internalFieldRef().field() =
-            UIndirectList<Type>(field, reader.cellMap());
-
-        // Fill boundary values
-        const auto& map = reader.faceMap();
-        if (map.size())
-        {
-            for (auto& pfld : fld.boundaryFieldRef())
-            {
-                const auto& pp = pfld.patch();
-
-                forAll(pfld, i)
-                {
-                    const label bFacei = pp.patch().offset()+i;
-                    pfld[i] = field[map[bFacei]];
-                }
-            }
-        }
-
-        regIOobject::store(std::move(tfld));
-    }
-}
-
 
 int main(int argc, char *argv[])
 {
@@ -122,14 +61,11 @@ int main(int argc, char *argv[])
     );
 
     argList::noParallel();
-    argList::addOptionCompat("no-fields", {"noFields", 2106});
     argList::addArgument("vtk-file", "The input legacy ascii vtk file");
 
     #include "setRootCase.H"
     #include "createTime.H"
 
-    const bool doFields = !args.found("no-fields");
-
     IFstream mshStream(args.get<fileName>(1));
 
     vtkUnstructuredReader reader(runTime, mshStream);
@@ -160,24 +96,6 @@ int main(int argc, char *argv[])
     mesh.removeFiles();
     mesh.write();
 
-
-    if (doFields)
-    {
-        // Re-read mesh as fvMesh so we can have fields
-        Info<< "Re-reading mesh ..." << endl;
-        #include "createMesh.H"
-
-        constructVolFields<scalar>(mesh, reader);
-        constructVolFields<vector>(mesh, reader);
-        constructVolFields<sphericalTensor>(mesh, reader);
-        constructVolFields<symmTensor>(mesh, reader);
-        constructVolFields<tensor>(mesh, reader);
-
-        // No need to write the mesh, only fields
-        mesh.thisDb().write();
-    }
-
-
     Info<< "End\n" << endl;
 
     return 0;

applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/cellShapeControl/cellShapeControlMesh/cellShapeControlMesh.C

@@ -743,12 +743,12 @@ Foam::label Foam::cellShapeControlMesh::estimateCellCount
                 cit->vertex(3)->point()
             );
 
-            const auto tetCentre = CGAL::centroid(tet);
+            pointFromPoint centre = topoint(CGAL::centroid(tet));
 
             if
             (
-                UPstream::parRun()
-             && !decomposition().positionOnThisProcessor(topoint(tetCentre))
+                Pstream::parRun()
+             && !decomposition().positionOnThisProcessor(centre)
             )
             {
                 continue;

applications/utilities/mesh/manipulation/transformPoints/transformPoints.C

@@ -83,10 +83,6 @@ Note
     pitch (rotation about y)
     yaw (rotation about z)
 
-    - with -rotate and two exactly opposing vectors it will actually mirror
-    the geometry. Use any of the other rotation options instead or use
-    two steps, each 90 degrees.
-
 \*---------------------------------------------------------------------------*/
 
 #include "argList.H"

applications/utilities/parallelProcessing/decomposePar/decomposePar.C

@@ -176,7 +176,7 @@ autoPtr<labelIOList> procAddressing
     const objectRegistry& procRegistry,
     const word& name,
     const word& instance,
-    const word& local = polyMesh::meshSubDir
+    const word& local = fvMesh::meshSubDir
 )
 {
     return autoPtr<labelIOList>::New
@@ -218,22 +218,17 @@ const labelIOList& procAddressing
     PtrList<labelIOList>& procAddressingList
 )
 {
-    const auto& procMesh = procMeshList[proci];
+    const fvMesh& procMesh = procMeshList[proci];
 
-    return procAddressingList.try_emplace
-    (
-        proci,
-        IOobject
+    if (!procAddressingList.set(proci))
+    {
+        procAddressingList.set
         (
-            name,
-            procMesh.facesInstance(),
-            polyMesh::meshSubDir,
-            procMesh,
-            IOobject::MUST_READ,
-            IOobject::NO_WRITE,
-            IOobject::NO_REGISTER
-        )
-    );
+            proci,
+            procAddressing(procMesh, name, procMesh.facesInstance())
+        );
+    }
+    return procAddressingList[proci];
 }
 
 

src/OpenFOAM/db/IOstreams/memory/ISpanStream.H

@@ -119,10 +119,10 @@ public:
         #if __cplusplus >= 201703L
         //- Construct (shallow copy) from std::string_view content
         explicit ispanstream(std::string_view s)
-        {
+        :
             buffer_type(const_cast<char*>(s.data()), s.size()),
-            stream_type(static_cast<buffer_type*>(this));
-        }
+            stream_type(static_cast<buffer_type*>(this))
+        {}
         #endif
 
         //- Construct (shallow copy) from span character content

src/OpenFOAM/db/error/error.C

@@ -157,7 +157,10 @@ Foam::error::error(const error& err)
     throwing_(err.throwing_),
     messageStreamPtr_(nullptr)
 {
-    if (err.messageStreamPtr_ && (err.messageStreamPtr_->count() > 0))
+    // FIXME: OStringStream copy construct does not adjust tellp and
+    // thus the count is wrong! (#3281)
+    // // if (err.messageStreamPtr_ && (err.messageStreamPtr_->count() > 0))
+    if (err.messageStreamPtr_)
     {
         messageStreamPtr_.reset(new OStringStream(*err.messageStreamPtr_));
     }

src/OpenFOAM/fields/GeometricFields/GeometricField/GeometricBoundaryField.C

@@ -888,7 +888,7 @@ boundaryInternalField() const
         *this
     );
 
-    auto& result = tresult;
+    auto& result = tresult.ref();
 
     forAll(result, patchi)
     {

src/OpenFOAM/interpolations/interpolateSplineXY/interpolateSplineXY.C

@@ -6,7 +6,6 @@
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
     Copyright (C) 2011 OpenFOAM Foundation
-    Copyright (C) 2024 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -61,11 +60,11 @@ Type Foam::interpolateSplineXY
     label n = xOld.size();
 
     // early exit if out of bounds or only one value
-    if (n == 1 || x <= xOld[0])
+    if (n == 1 || x < xOld[0])
     {
         return yOld[0];
     }
-    if (x >= xOld[n - 1])
+    if (x > xOld[n - 1])
     {
         return yOld[n - 1];
     }

src/Pstream/Allwmake-mpi

@@ -1,6 +1,7 @@
 #!/bin/sh
-cd "${0%/*}" || exit                                # Run from this directory
-. ${WM_PROJECT_DIR:?}/wmake/scripts/AllwmakeParseArguments -no-recursion "$@"
+cd "${0%/*}" || exit            # Run from this directory
+set -- -no-recursion "$@"       # Parse arguments only
+. ${WM_PROJECT_DIR:?}/wmake/scripts/AllwmakeParseArguments
 . ${WM_PROJECT_DIR:?}/wmake/scripts/wmakeFunctions  # Require wmake functions
 
 #------------------------------------------------------------------------------

src/dynamicMesh/polyTopoChange/polyTopoChange/hexRef8/hexRef8.C

@@ -2443,45 +2443,26 @@ Foam::labelList Foam::hexRef8::consistentSlowRefinement
     // Seed all boundary faces with owner value. This is to make sure that
     // they are visited (probably only important for coupled faces since
     // these need to be visited from both sides)
-    List<refinementData> nbrCellInfo;
-    syncTools::swapBoundaryCellList(mesh_, allCellInfo, nbrCellInfo);
-
     for (label facei = mesh_.nInternalFaces(); facei < mesh_.nFaces(); facei++)
     {
         // Check if face already handled in loop above
         if (!allFaceInfo[facei].valid(dummyTrackData))
         {
-            const label own = faceOwner[facei];
-            const auto& nbrInfo = nbrCellInfo[facei-mesh_.nInternalFaces()];
+            label own = faceOwner[facei];
 
-            if (allCellInfo[own].count() > nbrInfo.count())
-            {
-                allFaceInfo[facei].updateFace
-                (
-                    mesh_,
-                    facei,
-                    own,
-                    allCellInfo[own],
-                    FaceCellWave<refinementData, int>::propagationTol(),
-                    dummyTrackData
-                );
-                seedFaces.append(facei);
-                seedFacesInfo.append(allFaceInfo[facei]);
-            }
-            else if (allCellInfo[own].count() < nbrInfo.count())
-            {
-                allFaceInfo[facei].updateFace
-                (
-                    mesh_,
-                    facei,
-                    -1,         // Lucky! neighbCelli not used!
-                    nbrInfo,
-                    FaceCellWave<refinementData, int>::propagationTol(),
-                    dummyTrackData
-                );
-                seedFaces.append(facei);
-                seedFacesInfo.append(allFaceInfo[facei]);
-            }
+            // Seed face with transported data from owner.
+            refinementData faceData;
+            faceData.updateFace
+            (
+                mesh_,
+                facei,
+                own,
+                allCellInfo[own],
+                FaceCellWave<refinementData, int>::propagationTol(),
+                dummyTrackData
+            );
+            seedFaces.append(facei);
+            seedFacesInfo.append(faceData);
         }
     }
 

src/dynamicMesh/polyTopoChange/polyTopoChange/refinementData.C

@@ -6,7 +6,7 @@
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
     Copyright (C) 2011 OpenFOAM Foundation
-    Copyright (C) 2019-2022,2024 OpenCFD Ltd.
+    Copyright (C) 2019-2022 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -28,13 +28,6 @@ License
 
 #include "refinementData.H"
 
-// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
-
-Foam::refinementData Foam::transform(const tensor&, const refinementData val)
-{
-    return val;
-}
-
 // * * * * * * * * * * * * * * * IOstream Operators  * * * * * * * * * * * * //
 
 Foam::Ostream& Foam::operator<<

src/dynamicMesh/polyTopoChange/polyTopoChange/refinementData.H

@@ -6,7 +6,7 @@
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
     Copyright (C) 2011-2016 OpenFOAM Foundation
-    Copyright (C) 2019-2020,2024 OpenCFD Ltd.
+    Copyright (C) 2019-2020 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -216,9 +216,6 @@ public:
 
 // * * * * * * * * * * * * * * * * * Traits  * * * * * * * * * * * * * * * * //
 
-//- No-op rotational transform for base types
-refinementData transform(const tensor&, const refinementData val);
-
 //- Contiguous data for refinementData
 template<> struct is_contiguous<refinementData> : std::true_type {};
 

src/fileFormats/nastran/NASCore.C

@@ -5,7 +5,7 @@
     \\  /    A nd           | www.openfoam.com
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
-    Copyright (C) 2017-2023 OpenCFD Ltd.
+    Copyright (C) 2017-2024 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -136,7 +136,8 @@ std::string Foam::fileFormats::NASCore::nextNasField
 (
     const std::string& str,
     std::string::size_type& pos,
-    std::string::size_type len
+    const std::string::size_type width,
+    const bool free_format
 )
 {
     const auto beg = pos;
@@ -144,15 +145,23 @@ std::string Foam::fileFormats::NASCore::nextNasField
 
     if (end == std::string::npos)
     {
-        pos = beg + len;    // Continue after field width
+        if (free_format)
+        {
+            // Nothing left
+            pos = str.size();
+            return str.substr(beg);
+        }
+
+        // Fixed format - continue after field width
+        pos = beg + width;
+        return str.substr(beg, width);
     }
     else
     {
-        len = (end - beg);  // Efffective width
-        pos = end + 1;      // Continue after comma
+        // Free format - continue after comma
+        pos = end + 1;
+        return str.substr(beg, (end - beg));
     }
-
-    return str.substr(beg, len);
 }
 
 
@@ -249,8 +258,8 @@ void Foam::fileFormats::NASCore::writeCoord
     // 2 ID   : point ID - requires starting index of 1
     // 3 CP   : coordinate system ID                (blank)
     // 4 X1   : point x coordinate
-    // 5 X2   : point x coordinate
-    // 6 X3   : point x coordinate
+    // 5 X2   : point y coordinate
+    // 6 X3   : point z coordinate
     // 7 CD   : coordinate system for displacements (blank)
     // 8 PS   : single point constraints            (blank)
     // 9 SEID : super-element ID

src/fileFormats/nastran/NASCore.H

@@ -5,7 +5,7 @@
     \\  /    A nd           | www.openfoam.com
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
-    Copyright (C) 2017-2023 OpenCFD Ltd.
+    Copyright (C) 2017-2024 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -47,7 +47,6 @@ SourceFiles
 
 namespace Foam
 {
-
 namespace fileFormats
 {
 
@@ -79,18 +78,18 @@ public:
         //- Output load format
         enum loadFormat
         {
-            PLOAD2,
-            PLOAD4
+            PLOAD2,     //!< Face load (eg, pressure)
+            PLOAD4      //!< Vertex load
         };
 
-        //- Selection names for the NASTRAN file field formats
+        //- Selection names for the NASTRAN load formats
         static const Enum<loadFormat> loadFormatNames;
 
 
     // Constructors
 
         //- Default construct
-        NASCore() = default;
+        NASCore() noexcept = default;
 
 
     // Public Static Member Functions
@@ -98,18 +97,20 @@ public:
         //- Extract numbers from things like "-2.358-8" (same as "-2.358e-8")
         static scalar readNasScalar(const std::string& str);
 
-        //- A string::substr() to handle fixed-format and free-format NASTRAN.
-        //  Returns the substr to the next comma (if found) or the given length
-        //
-        //  \param str The string to extract from
-        //  \param pos On input, the position of the first character of the
-        //      substring. On output, advances to the next position to use.
-        //  \param len The fixed-format length to use if a comma is not found.
+        //- A std::string::substr() variant to handle fixed-format and
+        //- free-format NASTRAN.
+        //  Returns the substr until the next comma (if found)
+        //  or the given fixed width
         static std::string nextNasField
         (
+            //! The string to extract from
             const std::string& str,
+            //! [in,out] The parse position within \p str
             std::string::size_type& pos,
-            std::string::size_type len
+            //! The fixed-format width to use (if comma is not found)
+            const std::string::size_type width,
+            //! The input is known to be free-format
+            const bool free_format = false
         );
 
 

src/finiteVolume/fvMesh/wallDist/patchDistMethods/advectionDiffusion/advectionDiffusionPatchDistMethod.C

@@ -79,7 +79,7 @@ Foam::patchDistMethods::advectionDiffusion::advectionDiffusion
     tolerance_(coeffs_.getOrDefault<scalar>("tolerance", 1e-3)),
     maxIter_(coeffs_.getOrDefault<int>("maxIter", 10)),
     predicted_(false),
-    checkAndWriteMesh_(coeffs_.getOrDefault("checkAndWriteMesh", false))
+    checkAndWriteMesh_(coeffs_.getOrDefault("checkAndWriteMesh", true))
 {}
 
 

src/functionObjects/field/Curle/Curle.C

@@ -67,8 +67,7 @@ Foam::functionObjects::Curle::Curle
     c0_(0),
     rawFilePtrs_(),
     inputSurface_(),
-    surfaceWriterPtr_(nullptr),
-    warnOnNoPatch_(true)
+    surfaceWriterPtr_(nullptr)
 {
     read(dict);
 }
@@ -86,10 +85,8 @@ bool Foam::functionObjects::Curle::read(const dictionary& dict)
     }
 
     dict.readIfPresent("p", pName_);
-    dict.readIfPresent("warnOnNoPatch", warnOnNoPatch_);
 
-    patchIDs_ =
-        pbm.patchSet(dict.get<wordRes>("patches"), warnOnNoPatch_).sortedToc();
+    patchIDs_ = pbm.patchSet(dict.get<wordRes>("patches")).sortedToc();
 
     if (patchIDs_.empty())
     {

src/functionObjects/field/Curle/Curle.H

@@ -181,9 +181,6 @@ class Curle
         //- Ouput surface when sampling a surface
         autoPtr<surfaceWriter> surfaceWriterPtr_;
 
-        //- Flag to suppress warnings for missing patches
-        bool warnOnNoPatch_;
-
 
 public:
 

src/functionObjects/field/binField/binModels/binModel/binModel.C

@@ -122,7 +122,6 @@ Foam::binModel::binModel
     mesh_(mesh),
     decomposePatchValues_(false),
     cumulative_(false),
-    warnOnNoPatch_(true),
     coordSysPtr_(nullptr),
     nBin_(1)
 {}
@@ -139,11 +138,8 @@ bool Foam::binModel::read(const dictionary& dict)
         return false;
     }
 
-    dict.readIfPresent("warnOnNoPatch", warnOnNoPatch_);
-
     // Can also use pbm.indices(), but no warnings...
-    patchIDs_ =
-        pbm.patchSet(dict.get<wordRes>("patches"), warnOnNoPatch_).sortedToc();
+    patchIDs_ = pbm.patchSet(dict.get<wordRes>("patches")).sortedToc();
     fieldNames_ = dict.get<wordHashSet>("fields").sortedToc();
 
     wordRes zoneNames;

src/functionObjects/field/binField/binModels/binModel/binModel.H

@@ -77,9 +77,6 @@ protected:
         //- increasing distance in binning direction
         bool cumulative_;
 
-        //- Flag to suppress warnings for missing patches
-        bool warnOnNoPatch_;
-
         //- Local coordinate system of bins
         autoPtr<coordinateSystem> coordSysPtr_;
 

src/functionObjects/field/columnAverage/columnAverage.C

@@ -128,8 +128,7 @@ Foam::functionObjects::columnAverage::columnAverage
 :
     fvMeshFunctionObject(name, runTime, dict),
     patchIDs_(),
-    fieldSet_(mesh_),
-    warnOnNoPatch_(true)
+    fieldSet_(mesh_)
 {
     read(dict);
 }
@@ -141,13 +140,10 @@ bool Foam::functionObjects::columnAverage::read(const dictionary& dict)
 {
     fvMeshFunctionObject::read(dict);
 
-    dict.readIfPresent("warnOnNoPatch", warnOnNoPatch_);
-
     patchIDs_ =
         mesh_.boundaryMesh().patchSet
         (
-            dict.get<wordRes>("patches"),
-            warnOnNoPatch_
+            dict.get<wordRes>("patches")
         ).sortedToc();
 
     fieldSet_.read(dict);

src/functionObjects/field/columnAverage/columnAverage.H

@@ -132,9 +132,6 @@ class columnAverage
         mutable autoPtr<globalIndex> globalPoints_;
         mutable autoPtr<meshStructure> meshStructurePtr_;
 
-        //- Flag to suppress warnings for missing patches
-        bool warnOnNoPatch_;
-
 
     // Private Member Functions
 

src/functionObjects/field/heatTransferCoeff/heatTransferCoeffModels/heatTransferCoeffModel/heatTransferCoeffModel.C

@@ -51,8 +51,7 @@ Foam::heatTransferCoeffModel::heatTransferCoeffModel
 :
     mesh_(mesh),
     TName_(TName),
-    qrName_("qr"),
-    warnOnNoPatch_(true)
+    qrName_("qr")
 {}
 
 
@@ -151,10 +150,8 @@ bool Foam::heatTransferCoeffModel::read(const dictionary& dict)
     const polyBoundaryMesh& pbm = mesh_.boundaryMesh();
 
     dict.readIfPresent("qr", qrName_);
-    dict.readIfPresent("warnOnNoPatch", warnOnNoPatch_);
 
-    patchIDs_ =
-        pbm.patchSet(dict.get<wordRes>("patches"), warnOnNoPatch_).sortedToc();
+    patchIDs_ = pbm.patchSet(dict.get<wordRes>("patches")).sortedToc();
 
     return true;
 }

src/functionObjects/field/heatTransferCoeff/heatTransferCoeffModels/heatTransferCoeffModel/heatTransferCoeffModel.H

@@ -76,9 +76,6 @@ protected:
         //- Name of radiative heat flux field
         word qrName_;
 
-        //- Flag to suppress warnings for missing patches
-        bool warnOnNoPatch_;
-
 
     // Protected Member Functions
 

src/functionObjects/field/nearWallFields/nearWallFields.C

@@ -245,8 +245,7 @@ Foam::functionObjects::nearWallFields::nearWallFields
 )
 :
     fvMeshFunctionObject(name, runTime, dict),
-    fieldSet_(),
-    warnOnNoPatch_(true)
+    fieldSet_()
 {
     read(dict);
 }
@@ -262,11 +261,9 @@ bool Foam::functionObjects::nearWallFields::read(const dictionary& dict)
 
     dict.readEntry("fields", fieldSet_);
     dict.readEntry("distance", distance_);
-    dict.readIfPresent("warnOnNoPatch", warnOnNoPatch_);
 
     // Can also use pbm.indices(), but no warnings...
-    patchIDs_ =
-        pbm.patchSet(dict.get<wordRes>("patches"), warnOnNoPatch_).sortedToc();
+    patchIDs_ = pbm.patchSet(dict.get<wordRes>("patches")).sortedToc();
 
 
     // Clear out any previously loaded fields

src/functionObjects/field/nearWallFields/nearWallFields.H

@@ -134,9 +134,6 @@ protected:
             //- From resulting back to original field
             HashTable<word> reverseFieldMap_;
 
-            //- Flag to suppress warnings for missing patches
-            bool warnOnNoPatch_;
-
 
         // Calculated addressing
 

src/functionObjects/field/regionSizeDistribution/regionSizeDistribution.C

@@ -206,10 +206,7 @@ Foam::functionObjects::regionSizeDistribution::findPatchRegions
 
     labelHashSet patchRegions(2*regions.nRegions());
 
-    labelHashSet patchSet
-    (
-        mesh_.boundaryMesh().patchSet(patchNames_, warnOnNoPatch_)
-    );
+    labelHashSet patchSet(mesh_.boundaryMesh().patchSet(patchNames_));
 
     for (const label patchi : patchSet)
     {
@@ -372,8 +369,7 @@ Foam::functionObjects::regionSizeDistribution::regionSizeDistribution
     writeFile(obr_, name),
     alphaName_(dict.get<word>("field")),
     patchNames_(dict.get<wordRes>("patches")),
-    isoPlanes_(dict.getOrDefault("isoPlanes", false)),
-    warnOnNoPatch_(true)
+    isoPlanes_(dict.getOrDefault("isoPlanes", false))
 {
     read(dict);
 }
@@ -420,8 +416,6 @@ bool Foam::functionObjects::regionSizeDistribution::read(const dictionary& dict)
          direction_.normalise();
     }
 
-    dict.readIfPresent("warnOnNoPatch", warnOnNoPatch_);
-
     return true;
 }
 

src/functionObjects/field/regionSizeDistribution/regionSizeDistribution.H

@@ -249,9 +249,6 @@ class regionSizeDistribution
             //- Switch to enable iso-planes sampling
             bool isoPlanes_;
 
-        //- Flag to suppress warnings for missing patches
-        bool warnOnNoPatch_;
-
 
     // Private Member Functions
 

src/functionObjects/field/wallHeatFlux/wallHeatFlux.C

@@ -106,8 +106,7 @@ Foam::functionObjects::wallHeatFlux::wallHeatFlux
 :
     fvMeshFunctionObject(name, runTime, dict),
     writeFile(obr_, name, typeName, dict),
-    qrName_("qr"),
-    warnOnNoPatch_(true)
+    qrName_("qr")
 {
     read(dict);
 
@@ -145,13 +144,12 @@ bool Foam::functionObjects::wallHeatFlux::read(const dictionary& dict)
     writeFile::read(dict);
 
     dict.readIfPresent("qr", qrName_);
-    dict.readIfPresent("warnOnNoPatch", warnOnNoPatch_);
 
     wordRes patchNames;
     labelHashSet patchSet;
     if (dict.readIfPresent("patches", patchNames) && !patchNames.empty())
     {
-        patchSet = pbm.patchSet(patchNames, warnOnNoPatch_);
+        patchSet = pbm.patchSet(patchNames);
     }
 
     labelHashSet allWalls(pbm.findPatchIDs<wallPolyPatch>());

src/functionObjects/field/wallHeatFlux/wallHeatFlux.H

@@ -123,9 +123,6 @@ protected:
         //- Name of radiative heat flux name
         word qrName_;
 
-        //- Flag to suppress warnings for missing patches
-        bool warnOnNoPatch_;
-
 
     // Protected Member Functions
 

src/functionObjects/field/wallShearStress/wallShearStress.C

@@ -90,8 +90,7 @@ Foam::functionObjects::wallShearStress::wallShearStress
 :
     fvMeshFunctionObject(name, runTime, dict),
     writeFile(mesh_, name, typeName, dict),
-    writeFields_(true),  // May change in the future
-    warnOnNoPatch_(true)
+    writeFields_(true)  // May change in the future
 {
     read(dict);
 
@@ -128,7 +127,6 @@ bool Foam::functionObjects::wallShearStress::read(const dictionary& dict)
 
     writeFields_ = true;   // May change in the future
     dict.readIfPresent("writeFields", writeFields_);
-    dict.readIfPresent("warnOnNoPatch", warnOnNoPatch_);
 
     const polyBoundaryMesh& pbm = mesh_.boundaryMesh();
 
@@ -136,7 +134,7 @@ bool Foam::functionObjects::wallShearStress::read(const dictionary& dict)
     labelHashSet patchSet;
     if (dict.readIfPresent("patches", patchNames) && !patchNames.empty())
     {
-        patchSet = pbm.patchSet(patchNames, warnOnNoPatch_);
+        patchSet = pbm.patchSet(patchNames);
     }
 
     labelHashSet allWalls(pbm.findPatchIDs<wallPolyPatch>());

src/functionObjects/field/wallShearStress/wallShearStress.H

@@ -133,9 +133,6 @@ class wallShearStress
         //- Write the shear stress field ?
         bool writeFields_;
 
-        //- Flag to suppress warnings for missing patches
-        bool warnOnNoPatch_;
-
 
 protected:
 

src/functionObjects/forces/forces/forces.C

@@ -550,8 +550,7 @@ Foam::functionObjects::forces::forces
     directForceDensity_(false),
     porosity_(false),
     writeFields_(false),
-    initialised_(false),
-    warnOnNoPatch_(true)
+    initialised_(false)
 {
     if (readFields)
     {
@@ -590,8 +589,7 @@ Foam::functionObjects::forces::forces
     directForceDensity_(false),
     porosity_(false),
     writeFields_(false),
-    initialised_(false),
-    warnOnNoPatch_(true)
+    initialised_(false)
 {
     if (readFields)
     {
@@ -615,13 +613,10 @@ bool Foam::functionObjects::forces::read(const dictionary& dict)
 
     initialised_ = false;
 
-    dict.readIfPresent("warnOnNoPatch", warnOnNoPatch_);
-
     Info<< type() << ' ' << name() << ':' << endl;
 
     // Can also use pbm.indices(), but no warnings...
-    patchIDs_ =
-        pbm.patchSet(dict.get<wordRes>("patches"), warnOnNoPatch_).sortedToc();
+    patchIDs_ = pbm.patchSet(dict.get<wordRes>("patches")).sortedToc();
 
     dict.readIfPresent("directForceDensity", directForceDensity_);
     if (directForceDensity_)

src/functionObjects/forces/forces/forces.H

@@ -264,9 +264,6 @@ protected:
             //- Flag of initialisation (internal)
             bool initialised_;
 
-            //- Flag to suppress warnings for missing patches
-            bool warnOnNoPatch_;
-
 
     // Protected Member Functions
 

src/functionObjects/solvers/electricPotential/electricPotential.C

@@ -208,7 +208,6 @@ Foam::functionObjects::electricPotential::electricPotential
         )
     ),
     fvOptions_(mesh_),
-    tol_(1),
     nCorr_(1),
     writeDerivedFields_(false),
     electricField_(false)
@@ -259,7 +258,6 @@ bool Foam::functionObjects::electricPotential::read(const dictionary& dict)
     dict.readIfPresent("sigma", sigma_);
     dict.readIfPresent("epsilonr", epsilonr_);
     dict.readIfPresent("nCorr", nCorr_);
-    dict.readIfPresent("tolerance", tol_);
     dict.readIfPresent("writeDerivedFields", writeDerivedFields_);
     dict.readIfPresent("electricField", electricField_);
 
@@ -348,10 +346,6 @@ bool Foam::functionObjects::electricPotential::execute()
 {
     Log << type() << " execute: " << name() << endl;
 
-    // Convergence monitor parameters
-    bool converged = false;
-    label iter = 0;
-
     tmp<volScalarField> tsigma = this->sigma();
     const auto& sigma = tsigma();
 
@@ -368,9 +362,7 @@ bool Foam::functionObjects::electricPotential::execute()
 
         fvOptions_.constrain(eVEqn);
 
-        ++iter;
-        converged = (eVEqn.solve().initialResidual() < tol_);
-        if (converged) break;
+        eVEqn.solve();
     }
 
     if (electricField_)
@@ -379,14 +371,6 @@ bool Foam::functionObjects::electricPotential::execute()
         E == -fvc::grad(eV);
     }
 
-    if (converged)
-    {
-        Log << type() << ": " << name() << ": "
-            << eV.name() << " is converged." << nl
-            << tab << "initial-residual tolerance: " << tol_ << nl
-            << tab << "outer iteration: " << iter << nl;
-    }
-
     Log << endl;
 
     return true;

src/functionObjects/solvers/electricPotential/electricPotential.H

@@ -124,7 +124,6 @@ Usage
         electricField         <bool>;
         E                     <word>;
         fvOptions             <dict>;
-        tolerance             <scalar>;
 
         // Inherited entries
         ...
@@ -144,7 +143,6 @@ Usage
       electricField | Flag to calculate electric field  | bool | no | false
       E    | Name of electric field           | word | no | electricPotential:E
       fvOptions | List of finite-volume options       | dict | no   | -
-      tolerance | Outer-loop initial-residual tolerance | scalar | no | 1
     \endtable
 
     The inherited entries are elaborated in:
@@ -220,9 +218,6 @@ class electricPotential
         //- Run-time selectable finite volume options
         fv::optionList fvOptions_;
 
-        //- Outer-loop initial-residual tolerance
-        scalar tol_;
-
         //- Number of corrector iterations
         int nCorr_;
 

src/functionObjects/solvers/energyTransport/energyTransport.C

@@ -26,6 +26,11 @@ License
 \*---------------------------------------------------------------------------*/
 
 #include "energyTransport.H"
+#include "surfaceFields.H"
+#include "fvmDdt.H"
+#include "fvmDiv.H"
+#include "fvmLaplacian.H"
+#include "fvmSup.H"
 #include "turbulentTransportModel.H"
 #include "turbulentFluidThermoModel.H"
 #include "addToRunTimeSelectionTable.H"
@@ -183,6 +188,23 @@ Foam::functionObjects::energyTransport::energyTransport
 )
 :
     fvMeshFunctionObject(name, runTime, dict),
+    fieldName_(dict.getOrDefault<word>("field", "T")),
+    phiName_(dict.getOrDefault<word>("phi", "phi")),
+    rhoName_(dict.getOrDefault<word>("rho", "rho")),
+    nCorr_(0),
+    schemesField_("unknown-schemesField"),
+    fvOptions_(mesh_),
+    multiphaseThermo_(dict.subOrEmptyDict("phaseThermos")),
+    Cp_("Cp", dimEnergy/dimMass/dimTemperature, 0, dict),
+    kappa_
+    (
+        "kappa",
+        dimEnergy/dimTime/dimLength/dimTemperature,
+        0,
+        dict
+    ),
+    rho_("rhoInf", dimDensity, 0, dict),
+    Prt_("Prt", dimless, 1, dict),
     rhoCp_
     (
         IOobject
@@ -196,25 +218,7 @@ Foam::functionObjects::energyTransport::energyTransport
         ),
         mesh_,
         dimensionedScalar(dimEnergy/dimTemperature/dimVolume, Zero)
-    ),
-    fvOptions_(mesh_),
-    multiphaseThermo_(dict.subOrEmptyDict("phaseThermos")),
-    Cp_("Cp", dimEnergy/dimMass/dimTemperature, 0, dict),
-    kappa_
-    (
-        "kappa",
-        dimEnergy/dimTime/dimLength/dimTemperature,
-        0,
-        dict
-    ),
-    rho_("rhoInf", dimDensity, 0, dict),
-    Prt_("Prt", dimless, 1, dict),
-    fieldName_(dict.getOrDefault<word>("field", "T")),
-    schemesField_("unknown-schemesField"),
-    phiName_(dict.getOrDefault<word>("phi", "phi")),
-    rhoName_(dict.getOrDefault<word>("rho", "rho")),
-    tol_(1),
-    nCorr_(0)
+    )
 {
     read(dict);
 
@@ -301,14 +305,17 @@ Foam::functionObjects::energyTransport::energyTransport
 }
 
 
+// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
+
+Foam::functionObjects::energyTransport::~energyTransport()
+{}
+
+
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
 bool Foam::functionObjects::energyTransport::read(const dictionary& dict)
 {
-    if (!fvMeshFunctionObject::read(dict))
-    {
-        return false;
-    }
+    fvMeshFunctionObject::read(dict);
 
     dict.readIfPresent("phi", phiName_);
     dict.readIfPresent("rho", rhoName_);
@@ -316,7 +323,6 @@ bool Foam::functionObjects::energyTransport::read(const dictionary& dict)
     schemesField_ = dict.getOrDefault("schemesField", fieldName_);
 
     dict.readIfPresent("nCorr", nCorr_);
-    dict.readIfPresent("tolerance", tol_);
 
     if (dict.found("fvOptions"))
     {
@@ -349,16 +355,12 @@ bool Foam::functionObjects::energyTransport::execute()
     scalar relaxCoeff = 0;
     mesh_.relaxEquation(schemesField_, relaxCoeff);
 
-    // Convergence monitor parameters
-    bool converged = false;
-    int iter = 0;
-
     if (phi.dimensions() == dimMass/dimTime)
     {
         rhoCp_ = rho()*Cp();
         const surfaceScalarField rhoCpPhi(fvc::interpolate(Cp())*phi);
 
-        for (int i = 0; i <= nCorr_; ++i)
+        for (label i = 0; i <= nCorr_; i++)
         {
             fvScalarMatrix sEqn
             (
@@ -374,9 +376,7 @@ bool Foam::functionObjects::energyTransport::execute()
 
             fvOptions_.constrain(sEqn);
 
-            ++iter;
-            converged = (sEqn.solve(schemesField_).initialResidual() < tol_);
-            if (converged) break;
+            sEqn.solve(schemesField_);
         }
     }
     else if (phi.dimensions() == dimVolume/dimTime)
@@ -393,7 +393,7 @@ bool Foam::functionObjects::energyTransport::execute()
             rhoCp
         );
 
-        for (int i = 0; i <= nCorr_; ++i)
+        for (label i = 0; i <= nCorr_; i++)
         {
             fvScalarMatrix sEqn
             (
@@ -408,9 +408,7 @@ bool Foam::functionObjects::energyTransport::execute()
 
             fvOptions_.constrain(sEqn);
 
-            ++iter;
-            converged = (sEqn.solve(schemesField_).initialResidual() < tol_);
-            if (converged) break;
+            sEqn.solve(schemesField_);
         }
     }
     else
@@ -421,14 +419,6 @@ bool Foam::functionObjects::energyTransport::execute()
             << dimVolume/dimTime << exit(FatalError);
     }
 
-    if (converged)
-    {
-        Log << type() << ": " << name() << ": "
-            << s.name() << " is converged." << nl
-            << tab << "initial-residual tolerance: " << tol_ << nl
-            << tab << "outer iteration: " << iter << nl;
-    }
-
     Log << endl;
 
     return true;

src/functionObjects/solvers/energyTransport/energyTransport.H

@@ -30,180 +30,151 @@ Group
     grpSolversFunctionObjects
 
 Description
-    Computes the simplified energy transport equation in single-phase or
-    two-phase flow, considering incompressible cases:
+    Evolves a simplified energy transport equation for incompressible flows.
+    It takes into account the inertia, conduction and convection terms plus
+    a source.
 
-    \f[
-        \frac{\partial \rho \, C_p \, T}{\partial t}
-      + \nabla \cdot \left(\rho \, C_p \, \phi \, T \right)
-      - \nabla \cdot \left(\rho \, C_p \, \phi \right) \, T
-      - \nabla \cdot \left(\kappa_{eff} \, \nabla T \right)
-      = S_T
-    \f]
+    - The field name must be temperature and its BC's specified in the time
+      directory.
+    - The turbulence model should be incompressible
+    - In order to use in a incompressible multi phase a list of thermal
+      properties are needed. See below
 
-    where:
-    \vartable
-      T       | Scalar field
-      \rho    | (Generic) Fluid density which is unity when not specified
-      C_p     | Specific heat capacity at constant pressure
-      \phi    | (Generic) Flux field
-      \kappa_{eff}  | Effective thermal conductivity
-      S_T     | Scalar field source term
-    \endvartable
 
 Usage
-    Minimal example in \c system/controlDict.functions:
+    Example of function object specification to solve a energy transport
+    equation for a single phase flow plus a source term
     \verbatim
-    energyTransport1
+    functions
     {
-        // Mandatory entries
-        type            energyTransport;
-        libs            (solverFunctionObjects);
+        energy
+        {
+            type            energyTransport;
+            libs            (energyTransportFunctionObjects);
 
-        // Optional entries
-        field           <word>;
-        phi             <word>;
-        rho             <word>;
-        Cp              <scalar>;
-        kappa           <scalar>;
-        rhoInf          <scalar>;
-        Prt             <scalar>;
-        schemesField    <word>;
-        tolerance       <scalar>;
-        nCorr           <int>;
-        fvOptions       <dict>;
-        phaseThermos    <dict>;
+            enabled         true;
+            writeControl    writeTime;
+            writeInterval   1;
 
-        // Inherited entries
-        ...
+            field           T;
+
+            // volumetric Flux
+            phi             phi;
+
+            // Thermal properties
+            Cp              Cp    [J/kg/K]  1e3;
+            kappa           kappa [W/m/K]   0.0257;
+            rhoInf          rho   [kg/m^3]  1.2;
+
+            write           true;
+
+            fvOptions
+            {
+                viscousDissipation
+                {
+                    type            viscousDissipation;
+                    enabled         true;
+
+                    viscousDissipationCoeffs
+                    {
+                        fields          (T);
+                        rhoInf          $....rhoInf;
+                    }
+                }
+            }
+        }
     }
     \endverbatim
 
-    where:
+    Example of function object specification to solve a energy transport
+    equation for a multiphase phase flow plus a source term
+
+    equation:
+    \verbatim
+    functions
+    {
+        energy
+        {
+            type            energyTransport;
+            libs            (energyTransportFunctionObjects);
+
+            enabled         true;
+            writeControl    writeTime;
+            writeInterval   1;
+
+            field           T;
+
+            // rho field name
+            rho             rho;
+            // mass flux for multiphase
+            phi             rhoPhi;
+
+            write           true;
+
+            // Thermal properties of the phases
+            phaseThermos
+            {
+                alpha.air
+                {
+                    Cp          1e3;
+                    kappa       0.0243;
+                }
+                alpha.mercury
+                {
+                    Cp          140;
+                    kappa       8.2;
+                }
+                alpha.oil
+                {
+                    Cp          2e3;
+                    kappa       0.2;
+                }
+                alpha.water
+                {
+                    Cp          4e3;
+                    kappa       0.6;
+                }
+            }
+
+
+            fvOptions
+            {
+                viscousDissipation
+                {
+                    type            viscousDissipation;
+                    enabled         true;
+
+                    viscousDissipationCoeffs
+                    {
+                        fields          (T);
+                        rho             rho; //rho Field
+                    }
+                }
+            }
+        }
+    }
+    \endverbatim
+
+    Where the entries comprise:
     \table
-      Property  | Description                         | Type | Reqd | Deflt
-      type      | Type name: energyTransport          | word | yes  | -
-      libs      | Library name: solverFunctionObjects | word | yes  | -
-      field     | Name of the passive-scalar field    | word | no   | s
-      phi       | Name of flux field                  | word | no   | phi
-      rho       | Name of density field               | word | no   | rho
-      Cp        | Specific heat capacity at constant pressure | scalar | no | 0
-      kappa     | Thermal conductivity                | scalar | no | 0
-      rhoInf    | Fluid density                       | scalar | no | 0
-      Prt       | Turbulent Prandtl number            | scalar | no | 1
-      schemesField | Name of field to specify schemes | word   | no | field
-      tolerance | Outer-loop initial-residual tolerance | scalar | no | 1
-      nCorr     | Number of outer-loop correctors     | int    | no | 0
-      fvOptions | List of finite-volume options       | dict   | no | -
-      phaseThermos | Dictionary for multi-phase thermo | dict  | no | null
+        Property     | Description             | Required    | Default value
+        type         | Type name: energyTransport | yes      |
+        field        | Name of the scalar field | no         | T
+        phi          | Name of flux field      | no          | phi
+        rho          | Name of density field   | no          | rho
+        nCorr        | Number of correctors    | no          | 0
+        schemesField | Name of field to specify schemes | no | field name
+        fvOptions    | List of scalar sources  | no          |
+        Cp           | Heat capacity for single phase | no   | 0
+        rhoInf       | Density for single phase | no         | 0
+        kappa        | Thermal conductivity for single phase | no   | 0
+        Prt          | Turbulent Prandlt number | no         | 1.0
+        phaseThermos | Dictionary for multi-phase thermo |no | null
+        fvOptions    | Opotional extra sources | no          | null
     \endtable
 
-    The inherited entries are elaborated in:
-      - \link fvMeshFunctionObject.H \endlink
-      - \link fvOption.H \endlink
-
-    An example of function object specification to solve a energy transport
-    equation for a single phase flow plus a source term:
-    \verbatim
-    energyTransport1
-    {
-        // Mandatory entries
-        type            energyTransport;
-        libs            (solverFunctionObjects);
-
-        // Optional entries
-        field           T;
-        phi             phi;
-        Cp              Cp    [J/kg/K]  1e3;
-        kappa           kappa [W/m/K]   0.0257;
-        rhoInf          rho   [kg/m^3]  1.2;
-        fvOptions
-        {
-            viscousDissipation
-            {
-                type            viscousDissipation;
-                enabled         true;
-
-                viscousDissipationCoeffs
-                {
-                    fields          (T);
-                    rhoInf          $....rhoInf;
-                }
-            }
-        }
-
-        // Inherited entries
-        enabled         true;
-        writeControl    writeTime;
-        writeInterval   1;
-    }
-    \endverbatim
-
-    An example of function object specification to solve a energy transport
-    equation for a multiphase phase flow plus a source term:
-    \verbatim
-    energyTransport1
-    {
-        // Mandatory entries
-        type            energyTransport;
-        libs            (solverFunctionObjects);
-
-        // Optional entries
-        field           T;
-        rho             rho;
-        phi             rhoPhi;
-
-        // Thermal properties of the phases
-        phaseThermos
-        {
-            alpha.air
-            {
-                Cp          1e3;
-                kappa       0.0243;
-            }
-            alpha.mercury
-            {
-                Cp          140;
-                kappa       8.2;
-            }
-            alpha.oil
-            {
-                Cp          2e3;
-                kappa       0.2;
-            }
-            alpha.water
-            {
-                Cp          4e3;
-                kappa       0.6;
-            }
-        }
-
-        fvOptions
-        {
-            viscousDissipation
-            {
-                type            viscousDissipation;
-                enabled         true;
-
-                viscousDissipationCoeffs
-                {
-                    fields          (T);
-                    rho             rho;
-                }
-            }
-        }
-
-        // Inherited entries
-        enabled         true;
-        writeControl    writeTime;
-        writeInterval   1;
-    }
-    \endverbatim
-
-Note
-  - The field name must be temperature and its boundary conditions
-    specified in the time directory.
-  - The turbulence model should be incompressible.
+See also
+    Foam::functionObjects::fvMeshFunctionObject
 
 SourceFiles
     energyTransport.C
@@ -232,10 +203,22 @@ class energyTransport
 :
     public fvMeshFunctionObject
 {
-    // Private Data
+    // Private data
 
-        //- Volumetric heat capacity field [J/m^3/K]
-        volScalarField rhoCp_;
+        //- Name of the transport field.
+        word fieldName_;
+
+        //- Name of flux field
+        word phiName_;
+
+        //- Name of density field
+        word rhoName_;
+
+        //- Number of corrector iterations (optional)
+        label nCorr_;
+
+        //- Name of field whose schemes are used (optional)
+        word schemesField_;
 
         //- Run-time selectable finite volume options, e.g. sources, constraints
         fv::optionList fvOptions_;
@@ -246,7 +229,7 @@ class energyTransport
         //- List of phase names
         wordList phaseNames_;
 
-        //- List of phase specific heat capacities at constant pressure
+        //- List of phase heat capacities
         PtrList<dimensionedScalar> Cps_;
 
         //- List of phase thermal diffusivity for temperature [J/m/s/K]
@@ -255,7 +238,7 @@ class energyTransport
         //- Unallocated phase list
         UPtrList<volScalarField> phases_;
 
-        //- Specific heat capacity at constant pressure for single phase flows
+        //- Heat capacity for single phase flows
         dimensionedScalar Cp_;
 
         //- Thermal diffusivity for temperature for single phase flows
@@ -267,23 +250,8 @@ class energyTransport
         //- Turbulent Prandt number
         dimensionedScalar Prt_;
 
-        //- Name of the transport field
-        word fieldName_;
-
-        //- Name of field whose schemes are used
-        word schemesField_;
-
-        //- Name of flux field
-        word phiName_;
-
-        //- Name of density field
-        word rhoName_;
-
-        //- Outer-loop initial-residual tolerance
-        scalar tol_;
-
-        //- Number of corrector iterations
-        int nCorr_;
+        //- rhoCp
+        volScalarField rhoCp_;
 
 
     // Private Member Functions
@@ -294,15 +262,21 @@ class energyTransport
         //- Return the diffusivity field
         tmp<volScalarField> kappaEff() const;
 
-        //- Return the density field, rho
+        //- Return rho field
         tmp<volScalarField> rho() const;
 
-        //- Return the specific heat capacity at constant pressure field, Cp
+        //- Return Cp
         tmp<volScalarField> Cp() const;
 
-        //- Return the thermal diffusivity field
+        //- Return kappa
         tmp<volScalarField> kappa() const;
 
+        //- No copy construct
+        energyTransport(const energyTransport&) = delete;
+
+        //- No copy assignment
+        void operator=(const energyTransport&) = delete;
+
 
 public:
 
@@ -322,7 +296,7 @@ public:
 
 
     //- Destructor
-    virtual ~energyTransport() = default;
+    virtual ~energyTransport();
 
 
     // Member Functions

src/functionObjects/solvers/scalarTransport/scalarTransport.C

@@ -27,6 +27,11 @@ License
 \*---------------------------------------------------------------------------*/
 
 #include "scalarTransport.H"
+#include "surfaceFields.H"
+#include "fvmDdt.H"
+#include "fvmDiv.H"
+#include "fvmLaplacian.H"
+#include "fvmSup.H"
 #include "CMULES.H"
 #include "turbulentTransportModel.H"
 #include "turbulentFluidThermoModel.H"
@@ -167,9 +172,7 @@ Foam::functionObjects::scalarTransport::scalarTransport
 )
 :
     fvMeshFunctionObject(name, runTime, dict),
-    fvOptions_(mesh_),
     fieldName_(dict.getOrDefault<word>("field", "s")),
-    schemesField_("unknown-schemesField"),
     phiName_(dict.getOrDefault<word>("phi", "phi")),
     rhoName_(dict.getOrDefault<word>("rho", "rho")),
     nutName_(dict.getOrDefault<word>("nut", "none")),
@@ -179,12 +182,11 @@ Foam::functionObjects::scalarTransport::scalarTransport
         dict.getOrDefault<word>("phasePhiCompressed", "alphaPhiUn")
     ),
     D_(0),
-    alphaD_(1),
-    alphaDt_(1),
-    tol_(1),
+    constantD_(false),
     nCorr_(0),
     resetOnStartUp_(false),
-    constantD_(false),
+    schemesField_("unknown-schemesField"),
+    fvOptions_(mesh_),
     bounded01_(dict.getOrDefault("bounded01", true))
 {
     read(dict);
@@ -200,30 +202,31 @@ Foam::functionObjects::scalarTransport::scalarTransport
 }
 
 
+// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
+
+Foam::functionObjects::scalarTransport::~scalarTransport()
+{}
+
+
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
 bool Foam::functionObjects::scalarTransport::read(const dictionary& dict)
 {
-    if (!fvMeshFunctionObject::read(dict))
-    {
-        return false;
-    }
+    fvMeshFunctionObject::read(dict);
 
     dict.readIfPresent("phi", phiName_);
     dict.readIfPresent("rho", rhoName_);
     dict.readIfPresent("nut", nutName_);
     dict.readIfPresent("phase", phaseName_);
-    dict.readIfPresent("phasePhiCompressed", phasePhiCompressedName_);
+    dict.readIfPresent("bounded01", bounded01_);
 
     schemesField_ = dict.getOrDefault("schemesField", fieldName_);
+    constantD_ = dict.readIfPresent("D", D_);
+    alphaD_ = dict.getOrDefault<scalar>("alphaD", 1);
+    alphaDt_ = dict.getOrDefault<scalar>("alphaDt", 1);
 
-    dict.readIfPresent("alphaD", alphaD_);
-    dict.readIfPresent("alphaDt", alphaDt_);
-    dict.readIfPresent("tolerance", tol_);
     dict.readIfPresent("nCorr", nCorr_);
     dict.readIfPresent("resetOnStartUp", resetOnStartUp_);
-    constantD_ = dict.readIfPresent("D", D_);
-    dict.readIfPresent("bounded01", bounded01_);
 
     if (dict.found("fvOptions"))
     {
@@ -253,10 +256,6 @@ bool Foam::functionObjects::scalarTransport::execute()
     scalar relaxCoeff = 0;
     mesh_.relaxEquation(schemesField_, relaxCoeff);
 
-    // Convergence monitor parameters
-    bool converged = false;
-    int iter = 0;
-
     // Two phase scalar transport
     if (phaseName_ != "none")
     {
@@ -273,7 +272,7 @@ bool Foam::functionObjects::scalarTransport::execute()
 
         // Solve
         tmp<surfaceScalarField> tTPhiUD;
-        for (int i = 0; i <= nCorr_; ++i)
+        for (label i = 0; i <= nCorr_; i++)
         {
             fvScalarMatrix sEqn
             (
@@ -286,13 +285,9 @@ bool Foam::functionObjects::scalarTransport::execute()
 
             sEqn.relax(relaxCoeff);
             fvOptions_.constrain(sEqn);
-
-            ++iter;
-            converged = (sEqn.solve(schemesField_).initialResidual() < tol_);
+            sEqn.solve(schemesField_);
 
             tTPhiUD = sEqn.flux();
-
-            if (converged) break;
         }
 
         if (bounded01_)
@@ -312,8 +307,9 @@ bool Foam::functionObjects::scalarTransport::execute()
     {
         const volScalarField& rho = lookupObject<volScalarField>(rhoName_);
 
-        for (int i = 0; i <= nCorr_; ++i)
+        for (label i = 0; i <= nCorr_; i++)
         {
+
             fvScalarMatrix sEqn
             (
                 fvm::ddt(rho, s)
@@ -327,14 +323,12 @@ bool Foam::functionObjects::scalarTransport::execute()
 
             fvOptions_.constrain(sEqn);
 
-            ++iter;
-            converged = (sEqn.solve(schemesField_).initialResidual() < tol_);
-            if (converged) break;
+            sEqn.solve(schemesField_);
         }
     }
     else if (phi.dimensions() == dimVolume/dimTime)
     {
-        for (int i = 0; i <= nCorr_; ++i)
+        for (label i = 0; i <= nCorr_; i++)
         {
             fvScalarMatrix sEqn
             (
@@ -349,9 +343,7 @@ bool Foam::functionObjects::scalarTransport::execute()
 
             fvOptions_.constrain(sEqn);
 
-            ++iter;
-            converged = (sEqn.solve(schemesField_).initialResidual() < tol_);
-            if (converged) break;
+            sEqn.solve(schemesField_);
         }
     }
     else
@@ -362,14 +354,6 @@ bool Foam::functionObjects::scalarTransport::execute()
             << dimVolume/dimTime << exit(FatalError);
     }
 
-    if (converged)
-    {
-        Log << type() << ": " << name() << ": "
-            << s.name() << " is converged." << nl
-            << tab << "initial-residual tolerance: " << tol_ << nl
-            << tab << "outer iteration: " << iter << nl;
-    }
-
     Log << endl;
 
     return true;

src/functionObjects/solvers/scalarTransport/scalarTransport.H

@@ -31,134 +31,110 @@ Group
     grpSolversFunctionObjects
 
 Description
-    Computes the transport equation for a passive scalar in single-phase or
-    two-phase flow, considering both incompressible and compressible cases:
+    Evolves a passive scalar transport equation.
 
-    \f[
-        \frac{\partial \rho \, T}{\partial t}
-      + \nabla \cdot \left( \phi_\alpha \, T \right)
-      - \nabla \cdot (D_T \, \nabla T)
-      = \alpha \, S_T
-    \f]
-
-    where:
-    \vartable
-      T       | Passive scalar field
-      \rho    | (Generic) Fluid density which is unity when not specified
-      \phi_\alpha | (Generic) Flux field
-      \alpha  | Phase fraction which is unity for single-phase flows
-      D_T     | Diffusivity representing the diffusive transport of T
-      S_T     | Passive-scalar field source term
-    \endvartable
+    - To specify the field name set the \c field entry
+    - To employ the same numerical schemes as another field set
+      the \c schemesField entry,
+    - The diffusivity can be set manually using the 'D' entry, retrieved
+      from the turbulence model or specified nut
+    - Alternatively if a turbulence model is available a turbulent diffusivity
+      may be constructed from the laminar and turbulent viscosities using the
+      optional diffusivity coefficients \c alphaD and \c alphaDt (which default
+      to 1):
+      \verbatim
+          D = alphaD*nu + alphaDt*nut
+      \endverbatim
+    - To specify a transport quantity within a phase enter phase.
+    - bounded01 bounds the transported scalar within 0 and 1.
 
 Usage
-    Minimal example in \c system/controlDict.functions:
+    Example of function object specification to solve a scalar transport
+    equation:
     \verbatim
-    scalarTransport1
+    functions
     {
-        // Mandatory entries
-        type            scalarTransport;
-        libs            (solverFunctionObjects);
-
-        // Optional entries
-        field           <word>;
-        phi             <word>;
-        rho             <word>;
-        nut             <word>;
-        phase           <word>;
-        phasePhiCompressed  <word>;
-        schemesField    <word>;
-        bounded01       <bool>;
-        D               <scalar>;
-        alphaD          <scalar>;
-        alphaDt         <scalar>;
-        tolerance       <scalar>;
-        nCorr           <int>;
-        resetOnStartUp  <bool>;
-        fvOptions       <dict>;
-
-        // Inherited entries
-        ...
-    }
-
-    where:
-    \table
-      Property  | Description                         | Type | Reqd | Deflt
-      type      | Type name: scalarTransport          | word | yes  | -
-      libs      | Library name: solverFunctionObjects | word | yes  | -
-      field     | Name of the passive-scalar field    | word | no   | s
-      phi       | Name of flux field                  | word | no   | phi
-      rho       | Name of density field               | word | no   | rho
-      nut       | Name of the turbulence viscosity    | word | no   | none
-      phase     | Name of the phase                   | word | no   | none
-      phasePhiCompressed | Name of compressed VOF flux | word | no | alphaPhiUn
-      schemesField | Name of field to specify schemes | word  | no | field
-      bounded01 | Bounds scalar between 0-1 for multiphase  | bool | no | true
-      D         | Diffusion coefficient               | scalar | no | -
-      alphaD    | Laminar diffusivity coefficient     | scalar | no | 1
-      alphaDt   | Turbulent diffusivity coefficient   | scalar | no | 1
-      tolerance | Outer-loop initial-residual tolerance | scalar | no | 1
-      nCorr     | Number of outer-loop correctors     | int    | no | 0
-      resetOnStartUp | Flag to reset field to zero on start-up | bool | no | no
-      fvOptions | List of finite-volume options       | dict   | no | -
-    \endtable
-
-    The inherited entries are elaborated in:
-      - \link fvMeshFunctionObject.H \endlink
-      - \link fvOption.H \endlink
-
-    An example of function object specification to solve a residence time
-    in a two-phase flow:
-    \verbatim
-    scalarTransport1
-    {
-        // Mandatory entries
-        type            scalarTransport;
-        libs            (solverFunctionObjects);
-
-        // Optional entries
-        field           s;
-        bounded01       false;
-        phase           alpha.water;
-        tolerance       1e-5;
-        resetOnStartUp  false;
-        fvOptions
+        scalar1
         {
-            unitySource
+            type            scalarTransport;
+            libs            (solverFunctionObjects);
+
+            resetOnStartUp  no;
+            region          cabin;
+            field           H2O;
+
+
+            fvOptions
             {
-                type            scalarSemiImplicitSource;
-                enabled         true;
-
-                selectionMode   all;
-                volumeMode      specific;
-
-                sources
-                {
-                    s           (1 0);
-                }
+                ...
             }
         }
-
-        // Inherited entries
-        enabled         true;
-        writeControl    writeTime;
-        writeInterval   1;
     }
     \endverbatim
 
-Note
-  - To use the same numerical schemes as another field,
-  set the \c schemesField entry.
-  - The diffusivity can be set manually using the \c D entry, obtained
-  from the turbulence model or specified as `nut`.
-  - Alternatively, if a turbulence model is available, turbulent diffusivity
-  can be constructed from the laminar and turbulent viscosities using the
-  optional diffusivity coefficients \c alphaD and \c alphaDt
-  (which default to 1):
+    Example of function object specification to solve a residence time
+    in a two phase flow:
+    equation:
+    \verbatim
+    functions
+    {
+        sTransport
+        {
+            type            scalarTransport;
+            libs            (solverFunctionObjects);
 
-    \f[
-        D = \alpha_D \, \nu + \alpha_{Dt} \, \nu_t
-    \f]
+            enabled         true;
+            writeControl    writeTime;
+            writeInterval   1;
+
+            field           s;
+            bounded01       false;
+            phase           alpha.water;
+
+            write           true;
+
+            fvOptions
+            {
+                unitySource
+                {
+                    type            scalarSemiImplicitSource;
+                    enabled         true;
+
+                    selectionMode   all;
+                    volumeMode      specific;
+
+                    sources
+                    {
+                        s           (1 0);
+                    }
+                }
+            }
+
+            resetOnStartUp  false;
+        }
+    }
+    \endverbatim
+
+    Where the entries comprise:
+    \table
+        Property     | Description             | Required    | Default value
+        type         | Type name: scalarTransport | yes      |
+        field        | Name of the scalar field | no         | s
+        phi          | Name of flux field      | no          | phi
+        rho          | Name of density field   | no          | rho
+        phase        | Name of the phase       | no          | none
+        nut          | Name of the turbulence viscosity | no | none
+        D            | Diffusion coefficient   | no          | auto generated
+        nCorr        | Number of correctors    | no          | 0
+        resetOnStartUp | Reset scalar to zero on start-up | no | no
+        schemesField | Name of field to specify schemes | no | field name
+        fvOptions    | List of scalar sources  | no          |
+        bounded01    | Bounds scalar between 0-1 for multiphase | no | true
+        phasePhiCompressed | Compressed flux for VOF | no | alphaPhiUn
+    \endtable
+
+See also
+    Foam::functionObjects::fvMeshFunctionObject
 
 SourceFiles
     scalarTransport.C
@@ -187,52 +163,49 @@ class scalarTransport
 :
     public fvMeshFunctionObject
 {
-    // Private Data
-
-        //- Run-time selectable finite volume options, e.g. sources, constraints
-        fv::optionList fvOptions_;
+    // Private data
 
         //- Name of the transport field.
         word fieldName_;
 
-        //- Name of field whose schemes are used
-        word schemesField_;
-
-        //- Name of flux field
+        //- Name of flux field (optional)
         word phiName_;
 
-        //- Name of density field
+        //- Name of density field (optional)
         word rhoName_;
 
-        //- Name of turbulent viscosity field
+        //- Name of turbulent viscosity field (optional)
         word nutName_;
 
-        //- Name of phase field
+        //- Name of phase field (optional)
         word phaseName_;
 
-        //- Name of phase field compressed flux
+        //- Name of phase field compressed flux (optional)
         word phasePhiCompressedName_;
 
-        //- Diffusion coefficient
+        //- Diffusion coefficient (optional)
         scalar D_;
 
-        //- Laminar diffusion coefficient
+        //- Flag to indicate whether a constant, uniform D_ is specified
+        bool constantD_;
+
+        //- Laminar diffusion coefficient (optional)
         scalar alphaD_;
 
-        //- Turbulent diffusion coefficient
+        //- Turbulent diffusion coefficient (optional)
         scalar alphaDt_;
 
-        //- Outer-loop initial-residual tolerance
-        scalar tol_;
-
-        //- Number of corrector iterations
-        int nCorr_;
+        //- Number of corrector iterations (optional)
+        label nCorr_;
 
         //- Flag to reset the scalar to zero on start-up
         bool resetOnStartUp_;
 
-        //- Flag to indicate whether a constant, uniform D_ is specified
-        bool constantD_;
+        //- Name of field whose schemes are used (optional)
+        word schemesField_;
+
+        //- Run-time selectable finite volume options, e.g. sources, constraints
+        fv::optionList fvOptions_;
 
         //- Bound scalar between 0-1 using MULES for multiphase case
         bool bounded01_;
@@ -250,6 +223,12 @@ class scalarTransport
             const surfaceScalarField& phi
         ) const;
 
+        //- No copy construct
+        scalarTransport(const scalarTransport&) = delete;
+
+        //- No copy assignment
+        void operator=(const scalarTransport&) = delete;
+
 
 public:
 
@@ -269,7 +248,7 @@ public:
 
 
     //- Destructor
-    virtual ~scalarTransport() = default;
+    virtual ~scalarTransport();
 
 
     // Member Functions

src/mesh/snappyHexMesh/meshRefinement/meshRefinement.C

@@ -2054,6 +2054,7 @@ Foam::labelList Foam::meshRefinement::intersectedPoints() const
 
     // Mark all points on faces that will become baffles
     bitSet isBoundaryPoint(mesh_.nPoints());
+    label nBoundaryPoints = 0;
 
     const labelList& surfIndex = surfaceIndex();
 
@@ -2061,7 +2062,15 @@ Foam::labelList Foam::meshRefinement::intersectedPoints() const
     {
         if (surfIndex[facei] != -1)
         {
-            isBoundaryPoint.set(faces[facei]);
+            const face& f = faces[facei];
+
+            forAll(f, fp)
+            {
+                if (isBoundaryPoint.set(f[fp]))
+                {
+                    nBoundaryPoints++;
+                }
+            }
         }
     }
 
@@ -2074,17 +2083,37 @@ Foam::labelList Foam::meshRefinement::intersectedPoints() const
     //    if (patchi != -1)
     //    {
     //        const polyPatch& pp = mesh_.boundaryMesh()[patchi];
+    //
+    //        label facei = pp.start();
+    //
     //        forAll(pp, i)
     //        {
-    //            isBoundaryPoint.set(faces[pp.start()+i]);
+    //            const face& f = faces[facei];
+    //
+    //            forAll(f, fp)
+    //            {
+    //                if (isBoundaryPoint.set(f[fp]))
+    //                    nBoundaryPoints++;
+    //                }
+    //            }
+    //            facei++;
     //        }
     //    }
     //}
 
-    // Make sure all processors have the same data
-    syncTools::syncPointList(mesh_, isBoundaryPoint, orEqOp<unsigned int>(), 0);
 
-    return isBoundaryPoint.sortedToc();
+    // Pack
+    labelList boundaryPoints(nBoundaryPoints);
+    nBoundaryPoints = 0;
+    forAll(isBoundaryPoint, pointi)
+    {
+        if (isBoundaryPoint.test(pointi))
+        {
+            boundaryPoints[nBoundaryPoints++] = pointi;
+        }
+    }
+
+    return boundaryPoints;
 }
 
 

src/meshTools/coordSet/writers/nastran/nastranCoordSetWriter.C

@@ -5,7 +5,7 @@
     \\  /    A nd           | www.openfoam.com
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
-    Copyright (C) 2018-2024 OpenCFD Ltd.
+    Copyright (C) 2018-2022 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -250,12 +250,12 @@ Foam::fileName Foam::coordSetWriters::nastranWriter::writeTemplate
             Info<< "Writing nastran geometry to " << outputFile << endl;
         }
 
-        if (!Foam::isDir(outputFile.path()))
+        if (!isDir(outputFile.path()))
         {
-            Foam::mkDir(outputFile.path());
+            mkDir(outputFile.path());
         }
 
-        OFstream os(IOstreamOption::ATOMIC, outputFile);
+        OFstream os(outputFile);
         fileFormats::NASCore::setPrecision(os, writeFormat_);
 
         os  << "TITLE=OpenFOAM " << outputFile.stem()
@@ -293,12 +293,12 @@ Foam::fileName Foam::coordSetWriters::nastranWriter::writeTemplate
             Info<< "Writing nastran geometry to " << outputFile << endl;
         }
 
-        if (!Foam::isDir(outputFile.path()))
+        if (!isDir(outputFile.path()))
         {
-            Foam::mkDir(outputFile.path());
+            mkDir(outputFile.path());
         }
 
-        OFstream os(IOstreamOption::ATOMIC, outputFile);
+        OFstream os(outputFile);
         fileFormats::NASCore::setPrecision(os, writeFormat_);
 
         os  << "TITLE=OpenFOAM " << outputFile.stem()

src/meshTools/edgeMesh/edgeFormats/nastran/NASedgeFormat.C

@@ -62,13 +62,17 @@ bool Foam::fileFormats::NASedgeFormat::read
 
     while (is.good())
     {
-        string::size_type linei = 0;  // parsing position within current line
         string line;
         is.getLine(line);
 
-        if (line.empty() || line[0] == '$')
+        if (line.empty())
         {
-            continue; // Skip empty or comment
+            continue;  // Ignore empty
+        }
+        else if (line[0] == '$')
+        {
+            // Ignore comment
+            continue;
         }
 
         // Check if character 72 is continuation
@@ -94,41 +98,69 @@ bool Foam::fileFormats::NASedgeFormat::read
         }
 
 
+        // Parsing position within current line
+        std::string::size_type linei = 0;
+
+        // Is free format if line contains a comma
+        const bool freeFormat = line.contains(',');
+
         // First word (column 0-8)
         const word cmd(word::validate(nextNasField(line, linei, 8)));
 
         if (cmd == "CBEAM" || cmd == "CROD")
         {
-            // discard elementId (8-16)
-            (void) nextNasField(line, linei, 8); // 8-16
-            // discard groupId (16-24)
-            (void) nextNasField(line, linei, 8); // 16-24
+            // Fixed format:
+            //  8-16 : element id
+            // 16-24 : group id
+            // 24-32 : vertex
+            // 32-40 : vertex
 
-            label a = readLabel(nextNasField(line, linei, 8)); // 24-32
-            label b = readLabel(nextNasField(line, linei, 8)); // 32-40
+            // discard elementId
+            (void) nextNasField(line, linei, 8, freeFormat);
+            // discard groupId
+            (void) nextNasField(line, linei, 8, freeFormat);
 
-            dynEdges.append(edge(a,b));
+            label a = readLabel(nextNasField(line, linei, 8, freeFormat));
+            label b = readLabel(nextNasField(line, linei, 8, freeFormat));
+
+            dynEdges.emplace_back(a,b);
         }
         else if (cmd == "PLOTEL")
         {
+            // Fixed format:
+            //  8-16 : element id
+            // 16-24 : vertex
+            // 24-32 : vertex
+            // 32-40 : vertex
+
             // discard elementId (8-16)
-            (void) nextNasField(line, linei, 8); // 8-16
+            (void) nextNasField(line, linei, 8, freeFormat);
 
-            label a = readLabel(nextNasField(line, linei, 8)); // 16-24
-            label b = readLabel(nextNasField(line, linei, 8)); // 24-32
+            label a = readLabel(nextNasField(line, linei, 8, freeFormat));
+            label b = readLabel(nextNasField(line, linei, 8, freeFormat));
 
-            dynEdges.append(edge(a,b));
+            dynEdges.emplace_back(a,b);
         }
         else if (cmd == "GRID")
         {
-            label index = readLabel(nextNasField(line, linei, 8)); // 8-16
-            (void) nextNasField(line, linei, 8); // 16-24
-            scalar x = readNasScalar(nextNasField(line, linei, 8)); // 24-32
-            scalar y = readNasScalar(nextNasField(line, linei, 8)); // 32-40
-            scalar z = readNasScalar(nextNasField(line, linei, 8)); // 40-48
+            // Fixed (short) format:
+            //  8-16 : point id
+            // 16-24 : coordinate system (unsupported)
+            // 24-32 : point x coordinate
+            // 32-40 : point y coordinate
+            // 40-48 : point z coordinate
+            // 48-56 : displacement coordinate system (optional, unsupported)
+            // 56-64 : single point constraints (optional, unsupported)
+            // 64-70 : super-element id (optional, unsupported)
 
-            pointId.append(index);
-            dynPoints.append(point(x, y, z));
+            label index = readLabel(nextNasField(line, linei, 8, freeFormat));
+            (void) nextNasField(line, linei, 8, freeFormat);
+            scalar x = readNasScalar(nextNasField(line, linei, 8, freeFormat));
+            scalar y = readNasScalar(nextNasField(line, linei, 8, freeFormat));
+            scalar z = readNasScalar(nextNasField(line, linei, 8, freeFormat));
+
+            pointId.push_back(index);
+            dynPoints.emplace_back(x, y, z);
         }
         else if (cmd == "GRID*")
         {
@@ -138,6 +170,8 @@ bool Foam::fileFormats::NASedgeFormat::read
             // GRID*      126   0 -5.55999875E+02 -5.68730474E+02
             // *         2.14897901E+02
 
+            // Cannot be long format and free format at the same time!
+
             label index = readLabel(nextNasField(line, linei, 16)); // 8-24
             (void) nextNasField(line, linei, 16); // 24-40
             scalar x = readNasScalar(nextNasField(line, linei, 16)); // 40-56
@@ -157,8 +191,8 @@ bool Foam::fileFormats::NASedgeFormat::read
             (void) nextNasField(line, linei, 8); // 0-8
             scalar z = readNasScalar(nextNasField(line, linei, 16)); // 8-16
 
-            pointId.append(index);
-            dynPoints.append(point(x, y, z));
+            pointId.push_back(index);
+            dynPoints.emplace_back(x, y, z);
         }
     }
 

src/optimisation/adjointOptimisation/adjoint/optimisation/optimisationManager/optimisationManager/designVariablesUpdate/designVariablesUpdate.C

@@ -136,82 +136,6 @@ void Foam::designVariablesUpdate::writeToCostFile(bool zeroAdjointSolns)
 }
 
 
-void Foam::designVariablesUpdate::checkConvergence
-(
-    const scalarField& oldCorrection
-)
-{
-    bool converged(false);
-    // Design variables convergence check
-    if (designVarsThreshold_)
-    {
-        const labelList& activeVarIDs =
-            designVars_->activeDesignVariables();
-        const scalarField correction(oldCorrection, activeVarIDs);
-        const scalarField activeVars(designVars_->getVars(), activeVarIDs);
-        const scalar scaledCorrection =
-            gMax(mag(correction)/(mag(activeVars) + SMALL));
-        DebugInfo
-            << "Current correction " << correction << nl
-            << "Active vars " << activeVars << endl;
-        Info<< "Max. scaled correction of the design variables = "
-            << scaledCorrection
-            << endl;
-        if (scaledCorrection < designVarsThreshold_())
-        {
-            Info<< tab << "Design variables have converged " << endl;
-            converged = true;
-        }
-    }
-    // Objective convergence check
-    if (objectiveThreshold_)
-    {
-        const scalar newObjective = computeObjectiveValue();
-        const scalar oldObjective = updateMethod_->getObjectiveValueOld();
-        const scalar relativeUpdate =
-            mag(newObjective - oldObjective)/(mag(oldObjective) + SMALL);
-        Info<< "Relative change of the objective value = "
-            << relativeUpdate
-            << endl;
-        if (relativeUpdate < objectiveThreshold_())
-        {
-            Info<< tab << "Objective function has converged " << endl;
-            converged = true;
-        }
-    }
-    // Feasibility check
-    const scalarField& constraints = updateMethod_->getConstraintValues();
-    const scalar feasibility = sum(pos(constraints)*constraints);
-    Info<< "Feasibility = " << feasibility << endl;
-    if (converged && feasibility < feasibilityThreshold_)
-    {
-        Info<< "Stopping criteria met and all constraints satisfied." << nl
-            << "Optimisation has converged, stopping ..." << nl << nl
-            << "End" << nl
-            << endl;
-        // Force writing of all objective and constraint functions, to get
-        // the converged results to files
-        for (adjointSolverManager& am : adjointSolvManagers_)
-        {
-            for (adjointSolver& as : am.adjointSolvers())
-            {
-                // Use dummy weighted objective
-                as.getObjectiveManager().writeObjectives();
-            }
-        }
-        writeToCostFile(true);
-        if (UPstream::parRun())
-        {
-            UPstream::exit(0);
-        }
-        else
-        {
-            std::exit(0);
-        }
-    }
-}
-
-
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 
 Foam::designVariablesUpdate::designVariablesUpdate
@@ -570,10 +494,86 @@ void Foam::designVariablesUpdate::postUpdate(const scalarField& oldCorrection)
                 }
             }
         }
+        checkConvergence();
     }
-    if (convergenceCriteriaDefined_)
+}
+
+
+void Foam::designVariablesUpdate::checkConvergence()
+{
+    if (!convergenceCriteriaDefined_)
     {
-        checkConvergence(oldCorrection);
+        return;
+    }
+
+    bool converged(false);
+    const scalarField& oldCorrection = updateMethod_->returnCorrection();
+    // Design variables convergence check
+    if (designVarsThreshold_)
+    {
+        const labelList& activeVarIDs =
+            designVars_->activeDesignVariables();
+        const scalarField correction(oldCorrection, activeVarIDs);
+        const scalarField activeVars(designVars_->getVars(), activeVarIDs);
+        const scalar scaledCorrection =
+            gMax(mag(correction)/(mag(activeVars) + SMALL));
+        DebugInfo
+            << "Current correction " << correction << nl
+            << "Active vars " << activeVars << endl;
+        Info<< "Max. scaled correction of the design variables = "
+            << scaledCorrection
+            << endl;
+        if (scaledCorrection < designVarsThreshold_())
+        {
+            Info<< tab << "Design variables have converged " << endl;
+            converged = true;
+        }
+    }
+    // Objective convergence check
+    if (objectiveThreshold_ && updateMethod_->getObjectiveValueOld())
+    {
+        const scalar newObjective = computeObjectiveValue();
+        const scalar oldObjective = updateMethod_->getObjectiveValueOld()();
+        const scalar relativeUpdate =
+            mag(newObjective - oldObjective)/(mag(oldObjective) + SMALL);
+        Info<< "Relative change of the objective value = "
+            << relativeUpdate
+            << endl;
+        if (relativeUpdate < objectiveThreshold_())
+        {
+            Info<< tab << "Objective function has converged " << endl;
+            converged = true;
+        }
+    }
+    // Feasibility check
+    const scalarField& constraints = updateMethod_->getConstraintValues();
+    const scalar feasibility = sum(pos(constraints)*constraints);
+    Info<< "Feasibility = " << feasibility << endl;
+    if (converged && feasibility < feasibilityThreshold_)
+    {
+        Info<< "Stopping criteria met and all constraints satisfied." << nl
+            << "Optimisation has converged, stopping ..." << nl << nl
+            << "End" << nl
+            << endl;
+        // Force writing of all objective and constraint functions, to get
+        // the converged results to files
+        for (adjointSolverManager& am : adjointSolvManagers_)
+        {
+            for (adjointSolver& as : am.adjointSolvers())
+            {
+                // Use dummy weighted objective
+                as.getObjectiveManager().writeObjectives();
+            }
+        }
+        writeToCostFile(true);
+        if (UPstream::parRun())
+        {
+            UPstream::exit(0);
+        }
+        else
+        {
+            std::exit(0);
+        }
     }
 }
 

src/optimisation/adjointOptimisation/adjoint/optimisation/optimisationManager/optimisationManager/designVariablesUpdate/designVariablesUpdate.H

@@ -138,11 +138,6 @@ protected:
         //- Write to cost file
         void writeToCostFile(bool zeroAdjointSolns = false);
 
-        //- Check if the optimisation loop has converged based on the provided
-        //- criteria
-        //  May terminate the program
-        void checkConvergence(const scalarField& oldCorrection);
-
 
 private:
 
@@ -216,6 +211,11 @@ public:
         //- in the fixedStep approach
         void postUpdate(const scalarField& oldCorrection);
 
+        //- Check if the optimisation loop has converged based on the provided
+        //- criteria
+        //  May terminate the program
+        void checkConvergence();
+
         //- Get access to design variables
         inline autoPtr<designVariables>& getDesignVariables()
         {

src/optimisation/adjointOptimisation/adjoint/optimisation/optimisationManager/optimisationManager/optimisationManager.C

@@ -172,6 +172,9 @@ void Foam::optimisationManager::fixedStepUpdate()
     // Solve primal equations
     solvePrimalEquations();
 
+    // Check the convergence criteria of the optimisation loop
+    dvUpdate_->checkConvergence();
+
     // Reset adjoint sensitivities in all adjoint solver managers
     clearSensitivities();
 

src/optimisation/adjointOptimisation/adjoint/optimisation/updateMethod/updateMethod/updateMethod.C

@@ -237,7 +237,7 @@ Foam::updateMethod::updateMethod
     objectiveDerivatives_(designVars().getVars().size(), Zero),
     constraintDerivatives_(0),
     objectiveValue_(0),
-    objectiveValueOld_(0),
+    objectiveValueOld_(nullptr),
     cValues_(0),
     correction_(readOrZeroField("correction", designVars().getVars().size())),
     cumulativeCorrection_(0),
@@ -334,7 +334,11 @@ void Foam::updateMethod::setObjectiveValue(const scalar value)
 
 void Foam::updateMethod::setObjectiveValueOld(const scalar value)
 {
-    objectiveValueOld_ = value;
+    if (!objectiveValueOld_)
+    {
+        objectiveValueOld_.reset(new scalar(Zero));
+    }
+    objectiveValueOld_.ref() = value;
 }
 
 
@@ -350,7 +354,8 @@ Foam::scalar Foam::updateMethod::getObjectiveValue() const
 }
 
 
-Foam::scalar Foam::updateMethod::getObjectiveValueOld() const
+const Foam::autoPtr<Foam::scalar>&
+Foam::updateMethod::getObjectiveValueOld() const
 {
     return objectiveValueOld_;
 }

src/optimisation/adjointOptimisation/adjoint/optimisation/updateMethod/updateMethod/updateMethod.H

@@ -84,7 +84,7 @@ protected:
 
         //- Old objective value
         //  Used for convergence checking
-        scalar objectiveValueOld_;
+        autoPtr<scalar> objectiveValueOld_;
 
         //- Constraint values
         scalarField cValues_;
@@ -239,7 +239,7 @@ public:
         scalar getObjectiveValue() const;
 
         //- Get old objective value
-        scalar getObjectiveValueOld() const;
+        const autoPtr<scalar>& getObjectiveValueOld() const;
 
         //- Get values of constraints
         const scalarField& getConstraintValues() const;

src/overset/cellCellStencil/cellCellStencil/cellCellStencilObject.C

@@ -5,7 +5,7 @@
     \\  /    A nd           | www.openfoam.com
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
-    Copyright (C) 2017,2024 OpenCFD Ltd.
+    Copyright (C) 2017 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -26,7 +26,7 @@ License
 \*---------------------------------------------------------------------------*/
 
 #include "cellCellStencilObject.H"
-#include "dynamicOversetFvMesh.H"
+
 
 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
 
@@ -36,25 +36,4 @@ namespace Foam
 }
 
 
-// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
-
-bool Foam::cellCellStencilObject::movePoints()
-{
-    // Update (if needed) underlying stencil
-    const bool changed = stencilPtr_().update();
-
-    if (changed)
-    {
-        const auto* oversetMeshPtr = isA<oversetFvMeshBase>(mesh());
-        if (oversetMeshPtr)
-        {
-            // Clear out any additional (ldu)addressing
-            const_cast<oversetFvMeshBase&>(*oversetMeshPtr).clearOut();
-        }
-    }
-
-    return changed;
-}
-
-
 // ************************************************************************* //

src/overset/cellCellStencil/cellCellStencil/cellCellStencilObject.H

@@ -5,7 +5,7 @@
     \\  /    A nd           | www.openfoam.com
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
-    Copyright (C) 2017-2019,2024 OpenCFD Ltd.
+    Copyright (C) 2017-2019 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -54,9 +54,19 @@ typedef MeshObject<fvMesh, MoveableMeshObject, cellCellStencilObject> Stencil;
 
 class cellCellStencilObject
 :
-    public Stencil,
+    public MeshObject<fvMesh, MoveableMeshObject, cellCellStencilObject>,
     public cellCellStencil
 {
+    // Private Typedefs
+
+        typedef MeshObject
+        <
+            fvMesh,
+            MoveableMeshObject,
+            cellCellStencilObject
+        > MeshObject_type;
+
+
     // Private Data
 
         autoPtr<cellCellStencil> stencilPtr_;
@@ -76,7 +86,8 @@ public:
             const bool update = true
         )
         :
-            Stencil(mesh),
+            MeshObject_type(mesh),
+
             cellCellStencil(mesh),
             stencilPtr_
             (
@@ -100,7 +111,10 @@ public:
     // Member Functions
 
         //- Callback for geometry motion
-        virtual bool movePoints();
+        virtual bool movePoints()
+        {
+            return stencilPtr_().update();
+        }
 
         //- Update stencils. Return false if nothing changed.
         virtual bool update()

src/overset/oversetFvMesh/dynamicOversetFvMesh/dynamicOversetFvMesh.C

@@ -5,7 +5,7 @@
     \\  /    A nd           | www.openfoam.com
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
-    Copyright (C) 2022,2024 OpenCFD Ltd.
+    Copyright (C) 2022 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -65,10 +65,6 @@ bool Foam::dynamicOversetFvMesh::update()
         return false;
     }
 
-    //Note: too late to do oversetFvMeshBase::clearOut() to get it
-    //      consistent with any new cell-cell stencil since 
-    //      dynamicMotionSolverListFvMesh already triggers
-    //      meshObject::movePoints on cellCellStencilObject
     oversetFvMeshBase::update();
 
     return true;

src/overset/oversetFvMesh/oversetFvMeshBase.C

@@ -5,7 +5,7 @@
     \\  /    A nd           | www.openfoam.com
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
-    Copyright (C) 2014-2022,2024 OpenCFD Ltd.
+    Copyright (C) 2014-2022 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -578,14 +578,6 @@ Foam::oversetFvMeshBase::primitiveLduAddr() const
 }
 
 
-void Foam::oversetFvMeshBase::clearOut()
-{
-    // Cell-cell stencil is already mesh object. Clear out local
-    // addressing to force rebuilding addressing
-    lduPtr_.clear();
-}
-
-
 bool Foam::oversetFvMeshBase::update()
 {
     {

src/overset/oversetFvMesh/oversetFvMeshBase.H

@@ -5,7 +5,7 @@
     \\  /    A nd           | www.openfoam.com
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
-    Copyright (C) 2015-2022,2024 OpenCFD Ltd.
+    Copyright (C) 2015-2022 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -37,6 +37,7 @@ SourceFiles
 #ifndef oversetFvMeshBase_H
 #define oversetFvMeshBase_H
 
+#include "fvMeshPrimitiveLduAddressing.H"
 #include "fvMeshPrimitiveLduAddressing.H"
 #include "lduInterfaceFieldPtrsList.H"
 #include "volFieldsFwd.H"
@@ -237,9 +238,6 @@ public:
             ) const;
 
 
-        //- Clear out local storage
-        void clearOut();
-
         //- Update the mesh for both mesh motion and topology change
         virtual bool update();
 

src/parallel/decompose/Allwmake-mpi

@@ -1,6 +1,7 @@
 #!/bin/sh
-cd "${0%/*}" || exit                                # Run from this directory
-. ${WM_PROJECT_DIR:?}/wmake/scripts/AllwmakeParseArguments -no-recursion "$@"
+cd "${0%/*}" || exit            # Run from this directory
+set -- -no-recursion "$@"       # Parse arguments only
+. ${WM_PROJECT_DIR:?}/wmake/scripts/AllwmakeParseArguments
 . ${WM_PROJECT_DIR:?}/wmake/scripts/wmakeFunctions  # Require wmake functions
 . ${WM_PROJECT_DIR:?}/wmake/scripts/have_scotch
 

src/rigidBodyDynamics/rigidBodySolvers/symplectic/symplectic.C

@@ -75,7 +75,6 @@ void Foam::RBD::rigidBodySolvers::symplectic::solve
     qDot() = qDot0() + 0.5*deltaT0()*qDdot();
     q() = q0() + deltaT()*qDot();
 
-    // Update joints from new locations (q, q0)
     correctQuaternionJoints();
 
     // Update the body-state prior to the evaluation of the restraints

src/surfMesh/surfaceFormats/nas/NASsurfaceFormat.C

@@ -183,8 +183,6 @@ bool Foam::fileFormats::NASsurfaceFormat<Face>::read
 
     while (is.good())
     {
-        // Parsing position within current line
-        std::string::size_type linei = 0;
         is.getLine(line);
 
         if (NASCore::debug > 1) Info<< "Process: " << line << nl;
@@ -320,16 +318,30 @@ bool Foam::fileFormats::NASsurfaceFormat<Face>::read
             }
         }
 
+
+        // Parsing position within current line
+        std::string::size_type linei = 0;
+
+        // Is free format if line contains a comma
+        const bool freeFormat = line.contains(',');
+
         // First word (column 0-8)
         const word cmd(word::validate(nextNasField(line, linei, 8)));
 
         if (cmd == "CTRIA3")
         {
-            label elemId = readLabel(nextNasField(line, linei, 8)); // 8-16
-            label groupId = readLabel(nextNasField(line, linei, 8)); // 16-24
-            const auto a = readLabel(nextNasField(line, linei, 8)); // 24-32
-            const auto b = readLabel(nextNasField(line, linei, 8)); // 32-40
-            const auto c = readLabel(nextNasField(line, linei, 8)); // 40-48
+            // Fixed format:
+            //  8-16 : element id
+            // 16-24 : group id
+            // 24-32 : vertex
+            // 32-40 : vertex
+            // 40-48 : vertex
+
+            label elemId = readLabel(nextNasField(line, linei, 8, freeFormat));
+            label groupId = readLabel(nextNasField(line, linei, 8, freeFormat));
+            const auto a = readLabel(nextNasField(line, linei, 8, freeFormat));
+            const auto b = readLabel(nextNasField(line, linei, 8, freeFormat));
+            const auto c = readLabel(nextNasField(line, linei, 8, freeFormat));
 
             // Convert groupId into zoneId
             const auto iterZone = zoneLookup.cfind(groupId);
@@ -358,12 +370,20 @@ bool Foam::fileFormats::NASsurfaceFormat<Face>::read
         }
         else if (cmd == "CQUAD4")
         {
-            label elemId = readLabel(nextNasField(line, linei, 8)); // 8-16
-            label groupId = readLabel(nextNasField(line, linei, 8)); // 16-24
-            const auto a = readLabel(nextNasField(line, linei, 8)); // 24-32
-            const auto b = readLabel(nextNasField(line, linei, 8)); // 32-40
-            const auto c = readLabel(nextNasField(line, linei, 8)); // 40-48
-            const auto d = readLabel(nextNasField(line, linei, 8)); // 48-56
+            // Fixed format:
+            //  8-16 : element id
+            // 16-24 : group id
+            // 24-32 : vertex
+            // 32-40 : vertex
+            // 40-48 : vertex
+            // 48-56 : vertex
+
+            label elemId = readLabel(nextNasField(line, linei, 8, freeFormat));
+            label groupId = readLabel(nextNasField(line, linei, 8, freeFormat));
+            const auto a = readLabel(nextNasField(line, linei, 8, freeFormat));
+            const auto b = readLabel(nextNasField(line, linei, 8, freeFormat));
+            const auto c = readLabel(nextNasField(line, linei, 8, freeFormat));
+            const auto d = readLabel(nextNasField(line, linei, 8, freeFormat));
 
             // Convert groupId into zoneId
             const auto iterZone = zoneLookup.cfind(groupId);
@@ -407,11 +427,21 @@ bool Foam::fileFormats::NASsurfaceFormat<Face>::read
         }
         else if (cmd == "GRID")
         {
-            label index = readLabel(nextNasField(line, linei, 8)); // 8-16
-            (void) nextNasField(line, linei, 8); // 16-24
-            scalar x = readNasScalar(nextNasField(line, linei, 8)); // 24-32
-            scalar y = readNasScalar(nextNasField(line, linei, 8)); // 32-40
-            scalar z = readNasScalar(nextNasField(line, linei, 8)); // 40-48
+            // Fixed (short) format:
+            //  8-16 : point id
+            // 16-24 : coordinate system (not supported)
+            // 24-32 : point x coordinate
+            // 32-40 : point y coordinate
+            // 40-48 : point z coordinate
+            // 48-56 : displacement coordinate system (optional, unsupported)
+            // 56-64 : single point constraints (optional, unsupported)
+            // 64-70 : super-element id (optional, unsupported)
+
+            label index = readLabel(nextNasField(line, linei, 8, freeFormat));
+            (void) nextNasField(line, linei, 8, freeFormat);
+            scalar x = readNasScalar(nextNasField(line, linei, 8, freeFormat));
+            scalar y = readNasScalar(nextNasField(line, linei, 8, freeFormat));
+            scalar z = readNasScalar(nextNasField(line, linei, 8, freeFormat));
 
             pointId.push_back(index);
             dynPoints.emplace_back(x, y, z);
@@ -424,6 +454,8 @@ bool Foam::fileFormats::NASsurfaceFormat<Face>::read
             // GRID*      126   0 -5.55999875E+02 -5.68730474E+02
             // *         2.14897901E+02
 
+            // Cannot be long format and free format at the same time!
+
             label index = readLabel(nextNasField(line, linei, 16)); // 8-24
             (void) nextNasField(line, linei, 16); // 24-40
             scalar x = readNasScalar(nextNasField(line, linei, 16)); // 40-56
@@ -453,7 +485,10 @@ bool Foam::fileFormats::NASsurfaceFormat<Face>::read
             // have the 'weird' format where the immediately preceeding
             // comment contains the information.
 
-            label groupId = readLabel(nextNasField(line, linei, 8)); // 8-16
+            // Fixed format:
+            //  8-16 : pshell id
+
+            label groupId = readLabel(nextNasField(line, linei, 8, freeFormat));
 
             if (lastComment.size() > 1 && !nameLookup.contains(groupId))
             {

src/surfMesh/writers/nastran/nastranSurfaceWriter.C

@@ -335,30 +335,17 @@ Foam::surfaceWriters::nastranWriter::nastranWriter
         separator_ = ",";
     }
 
-    // Explicit PLOAD2, PLOAD4 field specification
-    options.readIfPresent("PLOAD2", pload2_);
-    options.readIfPresent("PLOAD4", pload4_);
-
-    // Compatibility:
-    // Optional pairs of (field name => load format)
-    // Could make conditional on (pload2_.empty() && pload4_.empty())
-
     List<Pair<word>> fieldPairs;
-    options.readIfPresent("fields", fieldPairs);
+    options.readEntry("fields", fieldPairs);
 
     for (const Pair<word>& item : fieldPairs)
     {
-        const loadFormat format =
-            fileFormats::NASCore::loadFormatNames[item.second()];
-
-        if (format == loadFormat::PLOAD2)
-        {
-            pload2_.push_back(item.first());
-        }
-        else if (format == loadFormat::PLOAD4)
-        {
-            pload4_.push_back(item.first());
-        }
+        // (field name => load format)
+        fieldMap_.insert
+        (
+            item.first(),
+            fileFormats::NASCore::loadFormatNames[item.second()]
+        );
     }
 }
 
@@ -419,12 +406,12 @@ Foam::fileName Foam::surfaceWriters::nastranWriter::write()
 
     if (UPstream::master() || !parallel_)
     {
-        if (!Foam::isDir(outputFile.path()))
+        if (!isDir(outputFile.path()))
         {
-            Foam::mkDir(outputFile.path());
+            mkDir(outputFile.path());
         }
 
-        OFstream os(IOstreamOption::ATOMIC, outputFile);
+        OFstream os(outputFile);
         fileFormats::NASCore::setPrecision(os, writeFormat_);
 
         os  << "TITLE=OpenFOAM " << outputPath_.name() << " geometry" << nl

src/surfMesh/writers/nastran/nastranSurfaceWriter.H

@@ -39,9 +39,7 @@ Description
         fieldLevel  | Subtract field level before scaling   | no  | empty dict
         fieldScale  | Output field scaling                  | no  | empty dict
         commonGeometry | use separate geometry files        | no  | false
-        PLOAD2      | Field selection (words/regex) for PLOAD2 | no  |
-        PLOAD4      | Field selection (words/regex) for PLOAD4 | no  |
-        fields      | Compat: Field pairs for PLOAD2/PLOAD4 | no  |
+        fields      | Field pairs for PLOAD2/PLOAD4         | yes |
     \endtable
 
     For example,
@@ -58,11 +56,7 @@ Description
                "p.*"   0.01;  // [Pa] -> [mbar]
             }
 
-            // Specific NASTRAN load types
-            PLOAD2 ( pMean );
-            PLOAD4 ( "p.*" );
-
-            // old style specification
+            // OpenFOAM field name to NASTRAN load types
             fields
             (
                 (pMean  PLOAD2)
@@ -72,12 +66,6 @@ Description
     }
     \endverbatim
 
-    Unless otherwise specified, all fields will be treated as PLOAD2
-    output. Can optionally specify PLOAD4 output using a combination
-    of \c PLOAD4 (accept) and \c PLOAD2 (deny) entries. The older \c fields
-    specification is also accepted and will be transcribed to
-    corresponding PLOAD4, PLOAD2 entries.
-
     \section Output file locations
 
     The \c rootdir normally corresponds to something like
@@ -117,7 +105,7 @@ SourceFiles
 
 #include "surfaceWriter.H"
 #include "NASCore.H"
-#include "wordRes.H"
+#include "HashTable.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
@@ -150,18 +138,15 @@ private:
         //- Field format (width and separator)
         fieldFormat writeFormat_;
 
+        //- Mapping from field name to data format enumeration
+        HashTable<loadFormat> fieldMap_;
+
         //- Use common geometry file
         bool commonGeometry_;
 
         //- Separator (used for free format)
         word separator_;
 
-        //- Explicit selection for PLOAD2 output (deselects for PLOAD4)
-        wordRes pload2_;
-
-        //- Explicit selection for PLOAD4 output
-        wordRes pload4_;
-
 
     // Private Member Functions
 

src/surfMesh/writers/nastran/nastranSurfaceWriterImpl.C

@@ -103,13 +103,13 @@ Foam::Ostream& Foam::surfaceWriters::nastranWriter::writeFaceValue
     {
         case loadFormat::PLOAD2 :
         {
-            if (pTraits<Type>::nComponents > 1)
+            if (pTraits<Type>::nComponents == 1)
             {
-                writeValue(os, Foam::mag(value));
+                writeValue(os, value);
             }
             else
             {
-                writeValue(os, value);
+                writeValue(os, mag(value));
             }
 
             os << separator_;
@@ -121,9 +121,6 @@ Foam::Ostream& Foam::surfaceWriters::nastranWriter::writeFaceValue
         {
             writeValue(os, elemId);
 
-            // NOTE: these should actually be vertex values,
-            // but misused here to provide vector quantities!
-
             for (direction d = 0; d < pTraits<Type>::nComponents; ++d)
             {
                 os  << separator_;
@@ -159,23 +156,43 @@ Foam::fileName Foam::surfaceWriters::nastranWriter::writeTemplate
 
     checkOpen();
 
-    // Default is PLOAD2
-    loadFormat format = loadFormat::PLOAD2;
+    const loadFormat format
+    (
+        fieldMap_.lookup
+        (
+            fieldName,
+            // Default format
+            (
+                pTraits<Type>::nComponents == 1
+              ? loadFormat::PLOAD2
+              : loadFormat::PLOAD4
+            )
+        )
+    );
 
     if
     (
         !std::is_integral<Type>::value  // Handle 'Ids' etc silently
-     && !pload4_.empty()
+     && !fieldMap_.empty()
+     && !fieldMap_.found(fieldName)
     )
     {
-        const wordRes::filter matcher(pload4_, pload2_);
-
-        if (matcher(fieldName))
-        {
-            format = loadFormat::PLOAD4;
-        }
+        WarningInFunction
+            << "No mapping found between field " << fieldName
+            << " and corresponding Nastran field.  Available types:"
+            << fieldMap_ << nl;
     }
 
+    // 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;
+    }
+
+
     // Common geometry
     // Field:  rootdir/<TIME>/<field>_surfaceName.bdf
 
@@ -222,9 +239,9 @@ Foam::fileName Foam::surfaceWriters::nastranWriter::writeTemplate
     {
         const auto& values = tfield();
 
-        if (!Foam::isDir(outputFile.path()))
+        if (!isDir(outputFile.path()))
         {
-            Foam::mkDir(outputFile.path());
+            mkDir(outputFile.path());
         }
 
         const scalar timeValue(0);
@@ -234,7 +251,7 @@ Foam::fileName Foam::surfaceWriters::nastranWriter::writeTemplate
         DynamicList<face> decompFaces;
 
 
-        OFstream os(IOstreamOption::ATOMIC, outputFile);
+        OFstream os(outputFile);
         fileFormats::NASCore::setPrecision(os, writeFormat_);
 
         os  << "TITLE=OpenFOAM " << outputFile.name()

src/topoChangerFvMesh/rawTopoChangerFvMesh/rawTopoChangerFvMesh.C

@@ -65,6 +65,12 @@ Foam::rawTopoChangerFvMesh::rawTopoChangerFvMesh
 {}
 
 
+// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
+
+Foam::rawTopoChangerFvMesh::~rawTopoChangerFvMesh()
+{}
+
+
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
 bool Foam::rawTopoChangerFvMesh::update()
@@ -115,24 +121,30 @@ bool Foam::rawTopoChangerFvMesh::update()
             }
         }
 
-        for (const auto& map : topoChangeMap().facesFromFacesMap())
+        const List<objectMap>& fromFaces = topoChangeMap().facesFromFacesMap();
+
+        forAll(fromFaces, i)
         {
-            mappedFace.set(map.index());
+            mappedFace.set(fromFaces[i].index());
         }
 
-        for (const auto& map : topoChangeMap().facesFromEdgesMap())
+        const List<objectMap>& fromEdges = topoChangeMap().facesFromEdgesMap();
+
+        forAll(fromEdges, i)
         {
-            mappedFace.set(map.index());
+            mappedFace.set(fromEdges[i].index());
         }
 
-        for (const auto& map : topoChangeMap().facesFromPointsMap())
+        const List<objectMap>& fromPts = topoChangeMap().facesFromPointsMap();
+
+        forAll(fromPts, i)
         {
-            mappedFace.set(map.index());
+            mappedFace.set(fromPts[i].index());
         }
 
         // Set unmapped faces to zero
-        Info<< "rawTopoChangerFvMesh : zeroing unmapped boundary values." << nl;
-
+        Info<< "rawTopoChangerFvMesh : zeroing unmapped boundary values."
+            << endl;
         zeroUnmappedValues<scalar, fvPatchField, volMesh>(mappedFace);
         zeroUnmappedValues<vector, fvPatchField, volMesh>(mappedFace);
         zeroUnmappedValues<sphericalTensor, fvPatchField, volMesh>(mappedFace);
@@ -143,8 +155,8 @@ bool Foam::rawTopoChangerFvMesh::update()
         Info<< "rawTopoChangerFvMesh :"
             << " recreating phi for unmapped boundary values." << endl;
 
-        const auto& U = lookupObject<volVectorField>("U");
-        auto& phi = lookupObjectRef<surfaceScalarField>("phi");
+        const volVectorField& U = lookupObject<volVectorField>("U");
+        surfaceScalarField& phi = lookupObjectRef<surfaceScalarField>("phi");
 
         setUnmappedValues
         (

src/topoChangerFvMesh/rawTopoChangerFvMesh/rawTopoChangerFvMesh.H

@@ -38,8 +38,8 @@ SourceFiles
 
 \*---------------------------------------------------------------------------*/
 
-#ifndef Foam_rawTopoChangerFvMesh_H
-#define Foam_rawTopoChangerFvMesh_H
+#ifndef rawTopoChangerFvMesh_H
+#define rawTopoChangerFvMesh_H
 
 #include "topoChangerFvMesh.H"
 #include "bitSet.H"
@@ -49,6 +49,8 @@ SourceFiles
 namespace Foam
 {
 
+// Forward declaration of classes
+
 /*---------------------------------------------------------------------------*\
                       Class rawTopoChangerFvMesh Declaration
 \*---------------------------------------------------------------------------*/
@@ -92,9 +94,8 @@ public:
             const bool doInit=true
         );
 
-
     //- Destructor
-    virtual ~rawTopoChangerFvMesh() = default;
+    virtual ~rawTopoChangerFvMesh();
 
 
     // Member Functions

src/topoChangerFvMesh/rawTopoChangerFvMesh/rawTopoChangerFvMeshTemplates.C

@@ -6,7 +6,6 @@
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
     Copyright (C) 2011-2016 OpenFOAM Foundation
-    Copyright (C) 2024 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -26,6 +25,9 @@ License
 
 \*---------------------------------------------------------------------------*/
 
+#include "rawTopoChangerFvMesh.H"
+#include "Time.H"
+
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
 template<class Type, template<class> class PatchField, class GeoMesh>
@@ -40,7 +42,10 @@ void Foam::rawTopoChangerFvMesh::setUnmappedValues
 
     forAll(fld.boundaryField(), patchi)
     {
-        auto& fvp = const_cast<PatchField<Type>&>(fld.boundaryField()[patchi]);
+        PatchField<Type>& fvp = const_cast<PatchField<Type>&>
+        (
+            fld.boundaryField()[patchi]
+        );
 
         const label start = fvp.patch().start();
         forAll(fvp, i)
@@ -66,27 +71,33 @@ void Foam::rawTopoChangerFvMesh::zeroUnmappedValues
 {
     typedef GeometricField<Type, PatchField, GeoMesh> FieldType;
 
-    std::unique_ptr<FieldType> zeroFieldPtr;
+    const wordList fldNames(names(FieldType::typeName));
 
-    for (const word& fldName : names<FieldType>())
+    forAll(fldNames, i)
     {
-        FieldType& fld = lookupObjectRef<FieldType>(fldName);
-        //Pout<< "Checking field " << fld.name() << endl;
+        //Pout<< "Checking field " << fldNames[i] << endl;
 
-        if (!zeroFieldPtr)
-        {
-            zeroFieldPtr = std::make_unique<FieldType>
+        FieldType& fld = lookupObjectRef<FieldType>(fldNames[i]);
+
+        setUnmappedValues
+        (
+            fld,
+            mappedFace,
+            FieldType
             (
-                this->newIOobject("zero"),
+                IOobject
+                (
+                    "zero",
+                    time().timeName(),
+                    *this,
+                    IOobject::NO_READ,
+                    IOobject::NO_WRITE,
+                    IOobject::NO_REGISTER
+                ),
                 *this,
-                Foam::zero{},
-                dimless
-            );
-        }
-
-        zeroFieldPtr->dimensions().reset(fld.dimensions());
-
-        setUnmappedValues(fld, mappedFace, *zeroFieldPtr);
+                dimensioned<Type>(fld.dimensions(), Zero)
+            )
+        );
     }
 }
 

tutorials/mesh/blockMesh/sphere/system/blockMeshDict

@@ -35,7 +35,7 @@ geometry
 }
 
 // Box size
-vo   #eval{ $outerRadius/sqrt(3) };
+vo   #eval{ sqrt($outerRadius/3) };
 
 vertices
 (

tutorials/mesh/blockMesh/sphere7/system/blockMeshDict

@@ -38,7 +38,7 @@ geometry
 }
 
 // Outer box size
-vo   #eval{ $outerRadius/sqrt(3) };
+vo   #eval{ sqrt($outerRadius/3) };
 
 // Inner box size - % of overall dimension
 vi   #eval{ $vo * $innerRatio };

tutorials/mesh/blockMesh/sphere7ProjectedEdges/system/blockMeshDict

@@ -45,7 +45,7 @@ geometry
 
 
 // Outer box size (approximate)
-vo   #eval{ $outerRadius/sqrt(3) };
+vo   #eval{ sqrt($outerRadius/3) };
 
 // Inner box size - % of overall dimension
 vi   #eval{ $vo * $innerRatio };

wmake/rules/linux64Icx/c++

@@ -1,7 +1,7 @@
 #------------------------------------------------------------------------------
 include $(GENERAL_RULES)/Icx/c++
 
-c++ARCH    := -fp-model precise
+c++ARCH    := -pthread -fp-model precise
 
 ifneq (,$(strip $(WM_COMPILE_OPTION)))
     sinclude $(DEFAULT_RULES)/c++$(WM_COMPILE_OPTION)