diff --git a/src/OpenFOAM/db/IOstreams/memory/ISpanStream.H b/src/OpenFOAM/db/IOstreams/memory/ISpanStream.H
index 9b03bf7bbda5fc58e058b13fd46f2ed3d00a2e26..d1d259b1938c397b5df2e9499e31549f53714d48 100644
--- a/src/OpenFOAM/db/IOstreams/memory/ISpanStream.H
+++ b/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
diff --git a/src/OpenFOAM/db/error/error.C b/src/OpenFOAM/db/error/error.C
index c4f327c95d42915ee2343902330f92ebd7a52ae4..22b81ca642fae6e7b1c33acc5d2b15227f935a34 100644
--- a/src/OpenFOAM/db/error/error.C
+++ b/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_));
     }
diff --git a/src/OpenFOAM/fields/GeometricFields/GeometricField/GeometricBoundaryField.C b/src/OpenFOAM/fields/GeometricFields/GeometricField/GeometricBoundaryField.C
index aff9204eefe4ae764c2de34454293b1841182d52..b6cbc8ae6ed78dcd42c413aa1022e6cd4913557a 100644
--- a/src/OpenFOAM/fields/GeometricFields/GeometricField/GeometricBoundaryField.C
+++ b/src/OpenFOAM/fields/GeometricFields/GeometricField/GeometricBoundaryField.C
@@ -888,7 +888,7 @@ boundaryInternalField() const
         *this
     );
 
-    auto& result = tresult;
+    auto& result = tresult.ref();
 
     forAll(result, patchi)
     {
diff --git a/src/fileFormats/nastran/NASCore.C b/src/fileFormats/nastran/NASCore.C
index 80d8117bca7cbb0883692c54943e1e9847146759..8713f243bd3fad6e6b96a7414ab74d46af734340 100644
--- a/src/fileFormats/nastran/NASCore.C
+++ b/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
diff --git a/src/fileFormats/nastran/NASCore.H b/src/fileFormats/nastran/NASCore.H
index 422127fefe7f08a6ba7bb5b1623288b242879d3f..fc72ede6fa0a5fb2b8310b33c6c88866b06cdb10 100644
--- a/src/fileFormats/nastran/NASCore.H
+++ b/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
         );
 
 
diff --git a/src/meshTools/edgeMesh/edgeFormats/nastran/NASedgeFormat.C b/src/meshTools/edgeMesh/edgeFormats/nastran/NASedgeFormat.C
index 916d898cd65dbc23aff2556e19184d2c3a721068..62e1fe1b2b06d2f1824d292ecf122ba4732ab11e 100644
--- a/src/meshTools/edgeMesh/edgeFormats/nastran/NASedgeFormat.C
+++ b/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
+
+            // discard elementId
+            (void) nextNasField(line, linei, 8, freeFormat);
+            // discard groupId
+            (void) nextNasField(line, linei, 8, freeFormat);
 
-            label a = readLabel(nextNasField(line, linei, 8)); // 24-32
-            label b = readLabel(nextNasField(line, linei, 8)); // 32-40
+            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 == "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
-
-            pointId.append(index);
-            dynPoints.append(point(x, y, z));
+            // 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)
+
+            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);
         }
     }
 
diff --git a/src/optimisation/adjointOptimisation/adjoint/optimisation/optimisationManager/optimisationManager/designVariablesUpdate/designVariablesUpdate.C b/src/optimisation/adjointOptimisation/adjoint/optimisation/optimisationManager/optimisationManager/designVariablesUpdate/designVariablesUpdate.C
index 690b2dde1560e0ee732bbb60ce2345950580bbce..0247ae205e30f4f3cc86772e33b73bf8f1c7f26d 100644
--- a/src/optimisation/adjointOptimisation/adjoint/optimisation/optimisationManager/optimisationManager/designVariablesUpdate/designVariablesUpdate.C
+++ b/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);
+        }
     }
 }
 
diff --git a/src/optimisation/adjointOptimisation/adjoint/optimisation/optimisationManager/optimisationManager/designVariablesUpdate/designVariablesUpdate.H b/src/optimisation/adjointOptimisation/adjoint/optimisation/optimisationManager/optimisationManager/designVariablesUpdate/designVariablesUpdate.H
index 2ed9edddf2b0f9268c6bef087663840fc15595b9..69fd7c39f225eabe85c7123f4bd232b0585ace50 100644
--- a/src/optimisation/adjointOptimisation/adjoint/optimisation/optimisationManager/optimisationManager/designVariablesUpdate/designVariablesUpdate.H
+++ b/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()
         {
diff --git a/src/optimisation/adjointOptimisation/adjoint/optimisation/optimisationManager/optimisationManager/optimisationManager.C b/src/optimisation/adjointOptimisation/adjoint/optimisation/optimisationManager/optimisationManager/optimisationManager.C
index 796a373eaa423c39d396ab8f4a6bc0cc92c59fca..20671a895c5a09b4b4cc1351264087df6547e706 100644
--- a/src/optimisation/adjointOptimisation/adjoint/optimisation/optimisationManager/optimisationManager/optimisationManager.C
+++ b/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();
 
diff --git a/src/optimisation/adjointOptimisation/adjoint/optimisation/updateMethod/updateMethod/updateMethod.C b/src/optimisation/adjointOptimisation/adjoint/optimisation/updateMethod/updateMethod/updateMethod.C
index b08a2ffe2c6595f4d5e7bc44d2462d498dd36221..d7d2a78962657520c6679531d813a9bf42f0667b 100644
--- a/src/optimisation/adjointOptimisation/adjoint/optimisation/updateMethod/updateMethod/updateMethod.C
+++ b/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_;
 }
diff --git a/src/optimisation/adjointOptimisation/adjoint/optimisation/updateMethod/updateMethod/updateMethod.H b/src/optimisation/adjointOptimisation/adjoint/optimisation/updateMethod/updateMethod/updateMethod.H
index 49ee8c61c622e2fbd33c813c0a06f9de2d926065..759512346f5e5dde38939e78acf161e0b43627cd 100644
--- a/src/optimisation/adjointOptimisation/adjoint/optimisation/updateMethod/updateMethod/updateMethod.H
+++ b/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;
diff --git a/src/surfMesh/surfaceFormats/nas/NASsurfaceFormat.C b/src/surfMesh/surfaceFormats/nas/NASsurfaceFormat.C
index a836708018c58e5d643cc40f502f06363df0dc15..0e19dbc55a22f6cfa4733928eadc7772279389c5 100644
--- a/src/surfMesh/surfaceFormats/nas/NASsurfaceFormat.C
+++ b/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))
             {
diff --git a/wmake/rules/linux64Icx/c++ b/wmake/rules/linux64Icx/c++
index 315ee1b50916ad8c1495992cb53b3ad0b16b3f9c..98237d5904cb276bdf810ed7979897738020adfa 100644
--- a/wmake/rules/linux64Icx/c++
+++ b/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)