potentialFoam: Solve for velocity potential named Phi rather than using the pressure field for this purpose

The Phi field is read if available otherwise created automatically with boundary conditions obtained automatically from the pressure field if available (with optional name) otherwise inferred from the velocity field. Phi Laplacian scheme and solver specification are required. See tutorials for examples.
2025-12-08 06:57:46 +00:00 · 2015-02-14 11:03:37 +00:00
parent 2edacf0ceb
commit cdadf4865d
16 changed files with 129 additions and 98 deletions
--- a/applications/solvers/basic/potentialFoam/createFields.H
+++ b/applications/solvers/basic/potentialFoam/createFields.H
@ -1,21 +1,4 @@
-    Info<< "Reading field p\n" << endl;
-    volScalarField p
-    (
-        IOobject
-        (
-            "p",
-            runTime.timeName(),
-            mesh,
-            IOobject::MUST_READ,
-            IOobject::NO_WRITE
-        ),
-        mesh
-    );
-
-    p = dimensionedScalar("zero", p.dimensions(), 0.0);
-
-
-    Info<< "Reading field U\n" << endl;
+    Info<< "Reading velocity field U\n" << endl;
    volVectorField U
    (
        IOobject
@ -51,12 +34,65 @@
    }


-    label pRefCell = 0;
-    scalar pRefValue = 0.0;
+    // Default name for the pressure field
+    word pName("p");
+
+    // Update name of the pressure field from the command-line option
+    args.optionReadIfPresent("pName", pName);
+
+    // Infer the pressure BCs from the velocity BCs
+    wordList pBCTypes
+    (
+        U.boundaryField().size(),
+        fixedValueFvPatchScalarField::typeName
+    );
+
+    forAll(U.boundaryField(), patchi)
+    {
+        if (U.boundaryField()[patchi].fixesValue())
+        {
+            pBCTypes[patchi] = zeroGradientFvPatchScalarField::typeName;
+        }
+    }
+
+    Info<< "Constructing pressure field " << pName << nl << endl;
+    volScalarField p
+    (
+        IOobject
+        (
+            pName,
+            runTime.timeName(),
+            mesh,
+            IOobject::READ_IF_PRESENT,
+            IOobject::NO_WRITE
+        ),
+        mesh,
+        dimensionedScalar(pName, sqr(dimVelocity), 0),
+        pBCTypes
+    );
+
+    Info<< "Constructing velocity potential field Phi\n" << endl;
+    volScalarField Phi
+    (
+        IOobject
+        (
+            "Phi",
+            runTime.timeName(),
+            mesh,
+            IOobject::READ_IF_PRESENT,
+            IOobject::NO_WRITE
+        ),
+        mesh,
+        dimensionedScalar("Phi", dimLength*dimVelocity, 0),
+        p.boundaryField().types()
+    );
+
+    label PhiRefCell = 0;
+    scalar PhiRefValue = 0;
    setRefCell
    (
-        p,
-        potentialFlow,
-        pRefCell,
-        pRefValue
+        Phi,
+        potentialFlow.dict(),
+        PhiRefCell,
+        PhiRefValue
    );
--- a/applications/solvers/basic/potentialFoam/potentialFoam.C
+++ b/applications/solvers/basic/potentialFoam/potentialFoam.C
@ -2,7 +2,7 @@
  =========                 |
  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2015 OpenFOAM Foundation
     \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@ -25,29 +25,48 @@ Application
    potentialFoam

 Description
-    Simple potential flow solver which can be used to generate starting fields
-    for full Navier-Stokes codes.
+    Potential flow solver which solves for the velocity potential
+    from which the flux-field is obtained and velocity field by reconstructing
+    the flux.
+
+    This application is particularly useful to generate starting fields for
+    Navier-Stokes codes.

 \*---------------------------------------------------------------------------*/

 #include "fvCFD.H"
 #include "fvIOoptionList.H"
+#include "pisoControl.H"

 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

 int main(int argc, char *argv[])
 {
-    argList::addBoolOption("writep", "write the final pressure field");
+    argList::addOption
+    (
+        "pName",
+        "pName",
+        "Name of the pressure field"
+    );
+
    argList::addBoolOption
    (
        "initialiseUBCs",
-        "initialise U boundary conditions"
+        "Initialise U boundary conditions"
+    );
+
+    argList::addBoolOption
+    (
+        "writePhi",
+        "Write the velocity potential field"
    );

    #include "setRootCase.H"
    #include "createTime.H"
    #include "createMesh.H"
-    #include "readControls.H"
+
+    pisoControl potentialFlow(mesh, "potentialFlow");
+
    #include "createFields.H"
    #include "createFvOptions.H"

@ -63,55 +82,47 @@ int main(int argc, char *argv[])

    adjustPhi(phi, U, p);

-
-    for (int nonOrth=0; nonOrth<=nNonOrthCorr; nonOrth++)
+    // Non-orthogonal velocity potential corrector loop
+    while (potentialFlow.correctNonOrthogonal())
    {
-        fvScalarMatrix pEqn
+        fvScalarMatrix PhiEqn
        (
-            fvm::laplacian
-            (
-                dimensionedScalar
-                (
-                    "1",
-                    dimTime/p.dimensions()*dimensionSet(0, 2, -2, 0, 0),
-                    1
-                ),
-                p
-            )
+            fvm::laplacian(dimensionedScalar("1", dimless, 1), Phi)
         ==
            fvc::div(phi)
        );

-        pEqn.setReference(pRefCell, pRefValue);
-        pEqn.solve();
+        PhiEqn.setReference(PhiRefCell, PhiRefValue);
+        PhiEqn.solve();

-        if (nonOrth == nNonOrthCorr)
+        if (potentialFlow.finalNonOrthogonalIter())
        {
-            phi -= pEqn.flux();
+            phi -= PhiEqn.flux();
        }
    }

    fvOptions.makeAbsolute(phi);

-    Info<< "continuity error = "
+    Info<< "Continuity error = "
        << mag(fvc::div(phi))().weightedAverage(mesh.V()).value()
        << endl;

    U = fvc::reconstruct(phi);
    U.correctBoundaryConditions();

-    Info<< "Interpolated U error = "
+    Info<< "Interpolated velocity error = "
        << (sqrt(sum(sqr((fvc::interpolate(U) & mesh.Sf()) - phi)))
          /sum(mesh.magSf())).value()
        << endl;

-    // Force the write
+    // Write U and phi
    U.write();
    phi.write();

-    if (args.optionFound("writep"))
+    // Optionally write Phi
+    if (args.optionFound("writePhi"))
    {
-        p.write();
+        Phi.write();
    }

    runTime.functionObjects().end();
--- a/tutorials/basic/potentialFoam/cylinder/Allrun
+++ b/tutorials/basic/potentialFoam/cylinder/Allrun
@ -6,26 +6,6 @@ cd ${0%/*} || exit 1    # Run from this directory

 application=`getApplication`

-
-# This case uses the #codeStream which is disabled by default. Enable for
-# just this case.
-MAIN_CONTROL_DICT=`foamEtcFile controlDict`
-if [ -f "$MAIN_CONTROL_DICT" ]
-then
-    echo "Modifying ${MAIN_CONTROL_DICT} to enable allowSystemOperations"
-
-    # Clean up on termination and on Ctrl-C
-    trap 'mv ${MAIN_CONTROL_DICT}.$$ ${MAIN_CONTROL_DICT} 2>/dev/null; exit 0' \
-        EXIT TERM INT
-    cp ${MAIN_CONTROL_DICT} ${MAIN_CONTROL_DICT}.$$
-
-    echo "Enabling allowSystemOperations in ${MAIN_CONTROL_DICT}."
-
-    sed \
-    -e s/"\(allowSystemOperations[ \t]*\)\([0-9]\);"/"\1 1;"/g \
-    ${MAIN_CONTROL_DICT}.$$ > ${MAIN_CONTROL_DICT}
-fi
-
 cp -r 0.org 0 > /dev/null 2>&1
 runApplication blockMesh
 runApplication $application
--- a/tutorials/basic/potentialFoam/cylinder/system/fvSchemes
+++ b/tutorials/basic/potentialFoam/cylinder/system/fvSchemes
@ -48,7 +48,7 @@ snGradSchemes
 fluxRequired
 {
    default         no;
-    p               ;
+    Phi             ;
 }


--- a/tutorials/basic/potentialFoam/cylinder/system/fvSolution
+++ b/tutorials/basic/potentialFoam/cylinder/system/fvSolution
@ -17,7 +17,7 @@ FoamFile

 solvers
 {
-    p
+    Phi
    {
        solver          PCG;
        preconditioner  DIC;
--- a/tutorials/basic/potentialFoam/pitzDaily/system/fvSchemes
+++ b/tutorials/basic/potentialFoam/pitzDaily/system/fvSchemes
@ -48,7 +48,7 @@ snGradSchemes
 fluxRequired
 {
    default         no;
-    p               ;
+    Phi             ;
 }


--- a/tutorials/basic/potentialFoam/pitzDaily/system/fvSolution
+++ b/tutorials/basic/potentialFoam/pitzDaily/system/fvSolution
@ -17,7 +17,7 @@ FoamFile

 solvers
 {
-    p
+    Phi
    {
        solver          PCG;
        preconditioner  DIC;
--- a/tutorials/incompressible/pisoFoam/les/motorBike/motorBike/system/fvSchemes
+++ b/tutorials/incompressible/pisoFoam/les/motorBike/motorBike/system/fvSchemes
@ -55,7 +55,8 @@ snGradSchemes
 fluxRequired
 {
    default         no;
-    p               ;
+    p;
+    Phi;
 }

 wallDist
--- a/tutorials/incompressible/pisoFoam/les/motorBike/motorBike/system/fvSolution
+++ b/tutorials/incompressible/pisoFoam/les/motorBike/motorBike/system/fvSolution
@ -38,6 +38,11 @@ solvers
        relTol          0;
    };

+    Phi
+    {
+        $p;
+    }
+
    "(U|nuTilda)"
    {
        solver          smoothSolver;
--- a/tutorials/incompressible/simpleFoam/motorBike/system/fvSchemes
+++ b/tutorials/incompressible/simpleFoam/motorBike/system/fvSchemes
@ -53,6 +53,7 @@ fluxRequired
 {
    default         no;
    p;
+    Phi;
 }

 wallDist
--- a/tutorials/incompressible/simpleFoam/motorBike/system/fvSolution
+++ b/tutorials/incompressible/simpleFoam/motorBike/system/fvSolution
@ -30,6 +30,11 @@ solvers
        mergeLevels      1;
    }

+    Phi
+    {
+        $p;
+    }
+
    U
    {
        solver           smoothSolver;
--- a/tutorials/incompressible/simpleFoam/pipeCyclic/Allrun
+++ b/tutorials/incompressible/simpleFoam/pipeCyclic/Allrun
@ -7,26 +7,6 @@ cd ${0%/*} || exit 1    # Run from this directory
 # Get application directory
 application=`getApplication`

-# This case uses the #codeStream which is disabled by default. Enable for
-# just this case.
-MAIN_CONTROL_DICT=`foamEtcFile controlDict`
-if [ -f "$MAIN_CONTROL_DICT" ]
-then
-    echo "Modifying ${MAIN_CONTROL_DICT} to enable allowSystemOperations"
-
-    # Clean up on termination and on Ctrl-C
-    trap 'mv ${MAIN_CONTROL_DICT}.$$ ${MAIN_CONTROL_DICT} 2>/dev/null; exit 0' \
-        EXIT TERM INT
-    cp ${MAIN_CONTROL_DICT} ${MAIN_CONTROL_DICT}.$$
-
-    echo "Enabling allowSystemOperations in ${MAIN_CONTROL_DICT}."
-
-    sed \
-    -e s/"\(allowSystemOperations[ \t]*\)\([0-9]\);"/"\1 1;"/g \
-    ${MAIN_CONTROL_DICT}.$$ > ${MAIN_CONTROL_DICT}
-fi
-
-
 runApplication blockMesh
 runApplication topoSet
 runApplication refineHexMesh c0 -overwrite
--- a/tutorials/lagrangian/reactingParcelFoam/verticalChannel/system/fvSchemes
+++ b/tutorials/lagrangian/reactingParcelFoam/verticalChannel/system/fvSchemes
@ -58,6 +58,7 @@ fluxRequired
 {
    default         no;
    p;
+    Phi;
 }

 wallDist
--- a/tutorials/lagrangian/reactingParcelFoam/verticalChannel/system/fvSolution
+++ b/tutorials/lagrangian/reactingParcelFoam/verticalChannel/system/fvSolution
@ -68,6 +68,11 @@ solvers
        relTol          0;
    }

+    Phi
+    {
+        $p;
+    }
+
    "(Yi|O2|N2|H2O)"
    {
        solver          PBiCG;
@ -88,7 +93,7 @@ solvers

 potentialFlow
 {
-    // used for potentialFoam initialisation
+    // Used for potentialFoam initialisation
    nNonOrthogonalCorrectors 5;
 }

--- a/tutorials/lagrangian/simpleReactingParcelFoam/verticalChannel/system/fvSchemes
+++ b/tutorials/lagrangian/simpleReactingParcelFoam/verticalChannel/system/fvSchemes
@ -59,6 +59,7 @@ fluxRequired
 {
    default         no;
    p;
+    Phi;
 }

 wallDist
--- a/tutorials/lagrangian/simpleReactingParcelFoam/verticalChannel/system/fvSolution
+++ b/tutorials/lagrangian/simpleReactingParcelFoam/verticalChannel/system/fvSolution
@ -37,6 +37,11 @@ solvers
        maxIter          50;
    };

+    Phi
+    {
+        $p;
+    }
+
    "(U|Yi|h|k|omega)"
    {
        solver          smoothSolver;