ENH: potentialFoam tutorial: use 'coded' functionObject to generate analytical solution

applications/solvers/basic/potentialFoam/potentialFoam.C

@@ -49,6 +49,10 @@ int main(int argc, char *argv[])
 
     Info<< nl << "Calculating potential flow" << endl;
 
+    // Since solver contains no time loop it would never execute
+    // function objects so do it ourselves.
+    runTime.functionObjects().start();
+
     adjustPhi(phi, U, p);
 
     for (int nonOrth=0; nonOrth<=nNonOrthCorr; nonOrth++)
@@ -99,6 +103,9 @@ int main(int argc, char *argv[])
         p.write();
     }
 
+    runTime.functionObjects().end();
+
+
     Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
         << "  ClockTime = " << runTime.elapsedClockTime() << " s"
         << nl << endl;

tutorials/basic/potentialFoam/cylinder/analyticalCylinder/Make/files

-analyticalCylinder.C
-
-EXE = $(FOAM_USER_APPBIN)/analyticalCylinder

tutorials/basic/potentialFoam/cylinder/analyticalCylinder/Make/options

-EXE_INC = \
-    -I$(LIB_SRC)/finiteVolume/lnInclude
-
-EXE_LIBS = \
-    -lfiniteVolume

tutorials/basic/potentialFoam/cylinder/analyticalCylinder/analyticalCylinder.C

-/*---------------------------------------------------------------------------*\
-  =========                 |
-  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
-   \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2004-2010 OpenCFD Ltd.
-     \\/     M anipulation  |
--------------------------------------------------------------------------------
-License
-    This file is part of OpenFOAM.
-
-    OpenFOAM is free software: you can redistribute it and/or modify it
-    under the terms of the GNU General Public License as published by
-    the Free Software Foundation, either version 3 of the License, or
-    (at your option) any later version.
-
-    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
-    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
-    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
-    for more details.
-
-    You should have received a copy of the GNU General Public License
-    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
-
-Application
-    analyticalCylinder
-
-Description
-    Generates an analytical solution for potential flow around a cylinder.
-    Can be compared with the solution from the potentialFlow/cylinder example.
-
-\*---------------------------------------------------------------------------*/
-
-#include "fvCFD.H"
-
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-int main(int argc, char *argv[])
-{
-
-#   include "setRootCase.H"
-
-#   include "createTime.H"
-#   include "createMesh.H"
-#   include "createFields.H"
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-    Info<< "\nEvaluating analytical solution" << endl;
-
-    volVectorField centres = UA.mesh().C();
-    volScalarField magCentres = mag(centres);
-    volScalarField theta = acos((centres & vector(1,0,0))/magCentres);
-
-    volVectorField cs2theta =
-        cos(2*theta)*vector(1,0,0)
-      + sin(2*theta)*vector(0,1,0);
-
-    UA = uInfX*(dimensionedVector(vector(1,0,0))
-      - pow((radius/magCentres),2)*cs2theta);
-
-    // Force writing of UA (since time has not changed)
-    UA.write();
-
-    Info<< "end" << endl;
-
-    return 0;
-}
-
-// ************************************************************************* //

tutorials/basic/potentialFoam/cylinder/analyticalCylinder/createFields.H

-Info<< "Reading field U\n" << endl;
-volVectorField U
-(
-    IOobject
-    (
-        "U",
-        runTime.timeName(),
-        mesh,
-        IOobject::MUST_READ,
-        IOobject::NO_WRITE
-    ),
-    mesh
-);
-
-Info<< "Reading inlet velocity  uInfX\n" << endl;
-
-dimensionedScalar uInfX
-(
-    "uInfx",
-    dimensionSet(0, 1, -1, 0, 0),
-    U.boundaryField()[3][0].x()
-);
-Info << "U at inlet = " << uInfX.value() << " m/s" << endl;
-
-dimensionedScalar radius
-( 
-    "radius",
-    dimensionSet(0, 1, 0, 0, 0),
-    mag(U.mesh().boundary()[4].Cf()[0])
-);
-    
-Info << "Cylinder radius = " << radius.value() << " m" << endl;
-
-volVectorField UA
-(
-    IOobject
-    (
-        "UA",
-        runTime.timeName(),
-        mesh,
-        IOobject::NO_READ,
-        IOobject::AUTO_WRITE
-    ),
-    U
-);

tutorials/basic/potentialFoam/cylinder/system/controlDict

@@ -45,5 +45,74 @@ timePrecision   6;
 
 runTimeModifiable true;
 
+functions
+{
+    difference
+    {
+        functionObjectLibs ("libutilityFunctionObjects.so");
+        type coded;
+        redirectType error;
+        code
+        #{
+            // Lookup U
+            Info<< "Looking up field U\n" << endl;
+            const volVectorField& U = mesh().lookupObject<volVectorField>("U");
+
+            Info<< "Reading inlet velocity  uInfX\n" << endl;
+
+            dimensionedScalar uInfX
+            (
+                "uInfx",
+                dimensionSet(0, 1, -1, 0, 0),
+                U.boundaryField()[3][0].x()
+            );
+            Info << "U at inlet = " << uInfX.value() << " m/s" << endl;
+
+            dimensionedScalar radius
+            (
+                "radius",
+                dimensionSet(0, 1, 0, 0, 0),
+                mag(U.mesh().boundary()[4].Cf()[0])
+            );
+
+            Info << "Cylinder radius = " << radius.value() << " m" << endl;
+
+            volVectorField UA
+            (
+                IOobject
+                (
+                    "UA",
+                    mesh().time().timeName(),
+                    U.mesh(),
+                    IOobject::NO_READ,
+                    IOobject::AUTO_WRITE
+                ),
+                U
+            );
+
+            Info<< "\nEvaluating analytical solution" << endl;
+
+            volVectorField centres = UA.mesh().C();
+            volScalarField magCentres = mag(centres);
+            volScalarField theta = acos((centres & vector(1,0,0))/magCentres);
+
+            volVectorField cs2theta =
+                cos(2*theta)*vector(1,0,0)
+              + sin(2*theta)*vector(0,1,0);
+
+            UA = uInfX*(dimensionedVector(vector(1,0,0))
+              - pow((radius/magCentres),2)*cs2theta);
+
+            // Force writing of UA (since time has not changed)
+            UA.write();
+
+            volScalarField error("error", mag(U-UA)/mag(UA));
+            Info<<"Writing relative error in U to " << error.objectPath()
+                << endl;
+            error.write();
+        #};
+    }
+}
+
 
 // ************************************************************************* //