ENH: Updated forces and forceCoeffs so that p and U fields need not exist on construction

src/postProcessing/functionObjects/forces/forceCoeffs/forceCoeffs.C

@@ -156,80 +156,82 @@ void Foam::forceCoeffs::timeSet()
 
 void Foam::forceCoeffs::write()
 {
-    if (active_)
+    forces::calcForcesMoment();
+
+    if (!active_)
     {
-        forces::calcForcesMoment();
+        return;
+    }
 
-        if (Pstream::master())
+    if (Pstream::master())
+    {
+        functionObjectFile::write();
+
+        scalar pDyn = 0.5*rhoRef_*magUInf_*magUInf_;
+
+        Field<vector> totForce(force_[0] + force_[1] + force_[2]);
+        Field<vector> totMoment(moment_[0] + moment_[1] + moment_[2]);
+
+        List<Field<scalar> > coeffs(3);
+        coeffs[0].setSize(nBin_);
+        coeffs[1].setSize(nBin_);
+        coeffs[2].setSize(nBin_);
+
+        // lift, drag and moment
+        coeffs[0] = (totForce & liftDir_)/(Aref_*pDyn);
+        coeffs[1] = (totForce & dragDir_)/(Aref_*pDyn);
+        coeffs[2] = (totMoment & pitchAxis_)/(Aref_*lRef_*pDyn);
+
+        scalar Cl = sum(coeffs[0]);
+        scalar Cd = sum(coeffs[1]);
+        scalar Cm = sum(coeffs[2]);
+
+        scalar Clf = Cl/2.0 + Cm;
+        scalar Clr = Cl/2.0 - Cm;
+
+        file(0)
+            << obr_.time().value() << tab
+            << Cm << tab << Cd << tab << Cl << tab << Clf << tab << Clr
+            << endl;
+
+        if (log_)
         {
-            functionObjectFile::write();
+            Info<< type() << " output:" << nl
+                << "    Cm    = " << Cm << nl
+                << "    Cd    = " << Cd << nl
+                << "    Cl    = " << Cl << nl
+                << "    Cl(f) = " << Clf << nl
+                << "    Cl(r) = " << Clr << endl;
+        }
 
-            scalar pDyn = 0.5*rhoRef_*magUInf_*magUInf_;
+        if (nBin_ > 1)
-
+        {
-            Field<vector> totForce(force_[0] + force_[1] + force_[2]);
+            if (binCumulative_)
-            Field<vector> totMoment(moment_[0] + moment_[1] + moment_[2]);
-
-            List<Field<scalar> > coeffs(3);
-            coeffs[0].setSize(nBin_);
-            coeffs[1].setSize(nBin_);
-            coeffs[2].setSize(nBin_);
-
-            // lift, drag and moment
-            coeffs[0] = (totForce & liftDir_)/(Aref_*pDyn);
-            coeffs[1] = (totForce & dragDir_)/(Aref_*pDyn);
-            coeffs[2] = (totMoment & pitchAxis_)/(Aref_*lRef_*pDyn);
-
-            scalar Cl = sum(coeffs[0]);
-            scalar Cd = sum(coeffs[1]);
-            scalar Cm = sum(coeffs[2]);
-
-            scalar Clf = Cl/2.0 + Cm;
-            scalar Clr = Cl/2.0 - Cm;
-
-            file(0)
-                << obr_.time().value() << tab
-                << Cm << tab << Cd << tab << Cl << tab << Clf << tab << Clr
-                << endl;
-
-            if (log_)
             {
-                Info<< type() << " output:" << nl
+                for (label i = 1; i < coeffs[0].size(); i++)
-                    << "    Cm    = " << Cm << nl
-                    << "    Cd    = " << Cd << nl
-                    << "    Cl    = " << Cl << nl
-                    << "    Cl(f) = " << Clf << nl
-                    << "    Cl(r) = " << Clr << endl;
-            }
-
-            if (nBin_ > 1)
-            {
-                if (binCumulative_)
                 {
-                    for (label i = 1; i < coeffs[0].size(); i++)
+                    coeffs[0][i] += coeffs[0][i-1];
-                    {
+                    coeffs[1][i] += coeffs[1][i-1];
-                        coeffs[0][i] += coeffs[0][i-1];
+                    coeffs[2][i] += coeffs[2][i-1];
-                        coeffs[1][i] += coeffs[1][i-1];
-                        coeffs[2][i] += coeffs[2][i-1];
-                    }
                 }
-
-                file(1)<< obr_.time().value();
-
-                forAll(coeffs[0], i)
-                {
-                    file(1)
-                        << tab << coeffs[2][i]
-                        << tab << coeffs[1][i]
-                        << tab << coeffs[0][i];
-                }
-
-                file(1) << endl;
             }
 
-            if (log_)
+            file(1)<< obr_.time().value();
+
+            forAll(coeffs[0], i)
             {
-                Info<< endl;
+                file(1)
+                    << tab << coeffs[2][i]
+                    << tab << coeffs[1][i]
+                    << tab << coeffs[0][i];
             }
+
+            file(1) << endl;
+        }
+
+        if (log_)
+        {
+            Info<< endl;
         }
     }
 }

src/postProcessing/functionObjects/forces/forces/forces.C

@@ -130,6 +130,55 @@ void Foam::forces::writeFileHeader(const label i)
 }
 
 
+void Foam::forces::initialise()
+{
+    if (initialised_ || !active_)
+    {
+        return;
+    }
+
+    if (directForceDensity_)
+    {
+        if (!obr_.foundObject<volVectorField>(fDName_))
+        {
+            active_ = false;
+            WarningIn("void Foam::forces::initialise()")
+                << "Could not find " << fDName_ << " in database." << nl
+                << "    De-activating forces."
+                << endl;
+        }
+    }
+    else
+    {
+        if
+        (
+            !obr_.foundObject<volVectorField>(UName_)
+         || !obr_.foundObject<volScalarField>(pName_)
+         || (
+                rhoName_ != "rhoInf"
+             && !obr_.foundObject<volScalarField>(rhoName_)
+            )
+        )
+        {
+            active_ = false;
+
+            WarningIn("void Foam::forces::initialise()")
+                << "Could not find " << UName_ << ", " << pName_;
+
+            if (rhoName_ != "rhoInf")
+            {
+                Info<< " or " << rhoName_;
+            }
+
+            Info<< " in database." << nl
+                << "    De-activating forces." << endl;
+        }
+    }
+
+    initialised_ = true;
+}
+
+
 Foam::tmp<Foam::volSymmTensorField> Foam::forces::devRhoReff() const
 {
     typedef compressible::turbulenceModel cmpTurbModel;
@@ -474,7 +523,8 @@ Foam::forces::forces
     binMin_(GREAT),
     binPoints_(),
     binFormat_("undefined"),
-    binCumulative_(true)
+    binCumulative_(true),
+    initialised_(false)
 {
     // Check if the available mesh is an fvMesh otherise deactivate
     if (!isA<fvMesh>(obr_))
@@ -534,7 +584,8 @@ Foam::forces::forces
     binMin_(GREAT),
     binPoints_(),
     binFormat_("undefined"),
-    binCumulative_(true)
+    binCumulative_(true),
+    initialised_(false)
 {
     forAll(force_, i)
     {
@@ -556,6 +607,8 @@ void Foam::forces::read(const dictionary& dict)
 {
     if (active_)
     {
+        initialised_ = false;
+
         log_ = dict.lookupOrDefault<Switch>("log", false);
         directForceDensity_ = dict.lookupOrDefault("directForceDensity", false);
 
@@ -568,19 +621,6 @@ void Foam::forces::read(const dictionary& dict)
         {
             // Optional entry for fDName
             fDName_ = dict.lookupOrDefault<word>("fDName", "fD");
-
-            // Check whether fDName exists, if not deactivate forces
-            if
-            (
-                !obr_.foundObject<volVectorField>(fDName_)
-            )
-            {
-                active_ = false;
-                WarningIn("void forces::read(const dictionary&)")
-                    << "Could not find " << fDName_ << " in database." << nl
-                    << "    De-activating forces."
-                    << endl;
-            }
         }
         else
         {
@@ -589,32 +629,6 @@ void Foam::forces::read(const dictionary& dict)
             UName_ = dict.lookupOrDefault<word>("UName", "U");
             rhoName_ = dict.lookupOrDefault<word>("rhoName", "rho");
 
-            // Check whether UName, pName and rhoName exists,
-            // if not deactivate forces
-            if
-            (
-                !obr_.foundObject<volVectorField>(UName_)
-             || !obr_.foundObject<volScalarField>(pName_)
-             || (
-                    rhoName_ != "rhoInf"
-                 && !obr_.foundObject<volScalarField>(rhoName_)
-                )
-            )
-            {
-                active_ = false;
-
-                WarningIn("void forces::read(const dictionary&)")
-                    << "Could not find " << UName_ << ", " << pName_;
-
-                if (rhoName_ != "rhoInf")
-                {
-                    Info<< " or " << rhoName_;
-                }
-
-                Info<< " in database." << nl
-                    << "    De-activating forces." << endl;
-            }
-
             // Reference density needed for incompressible calculations
             rhoRef_ = readScalar(dict.lookup("rhoInf"));
 
@@ -745,13 +759,13 @@ void Foam::forces::timeSet()
 
 void Foam::forces::write()
 {
+    calcForcesMoment();
+
     if (!active_)
     {
         return;
     }
 
-    calcForcesMoment();
-
     if (Pstream::master())
     {
         functionObjectFile::write();
@@ -770,6 +784,13 @@ void Foam::forces::write()
 
 void Foam::forces::calcForcesMoment()
 {
+    initialise();
+
+    if (!active_)
+    {
+        return;
+    }
+
     force_[0] = vector::zero;
     force_[1] = vector::zero;
     force_[2] = vector::zero;

src/postProcessing/functionObjects/forces/forces/forces.H

@@ -32,7 +32,8 @@ Description
     pressure and skin-friction forces over a given list of patches.
 
     Member function forces::write() calculates the forces/moments and
-    writes the forces/moments into the file \<timeDir\>/forces.dat
+    writes the forces/moments into the file \<timeDir\>/forces.dat and bin
+    data (if selected) to the file \<timeDir\>/forces_bin.dat
 
     Example of function object specification:
     \verbatim
@@ -45,7 +46,6 @@ Description
         patches     (walls);
         nBin        20;
         binDir      (1 0 0);
-        binFormat   gnuplot;
     }
     \endverbatim
 
@@ -57,7 +57,6 @@ Description
         patches      | patches included in the forces calculation | yes |
         nBin         | number of data bins     | no          |
         binDir       | direction along which bins are defined | no |
-        binFormat    | output format for bin data | no       |
         pName        | pressure field name     | no          | p
         UName        | velocity field name     | no          | U
         rhoName      | density field name (see below) | no   | rho
@@ -216,6 +215,10 @@ protected:
                 bool binCumulative_;
 
 
+            //- Initialised flag
+            bool initialised_;
+
+
     // Protected Member Functions
 
         //- Create file names for forces and bins
@@ -224,6 +227,9 @@ protected:
         //- Output file header information
         virtual void writeFileHeader(const label i);
 
+        //- Initialise the fields
+        void initialise();
+
         //- Return the effective viscous stress (laminar + turbulent).
         tmp<volSymmTensorField> devRhoReff() const;