Commit 81cea099 authored by Henry Weller's avatar Henry Weller Committed by Andrew Heather
Browse files

ENH: pimpleDyMFoam: Improved efficiency and consistency when running on a static mesh

Now pimpleDyMFoam is exactly equivalent to pimpleFoam when running on a
staticFvMesh.  Also when the constant/dynamicMeshDict is not present a
staticFvMesh is automatically constructed so that the pimpleDyMFoam solver can
run any pimpleFoam case without change.

pimpleDyMFoam: Store Uf as an autoPtr for better error handling

pimpleFoam: Set initial deltaT from the Courant number

for improved stability on start-up and compatibility with pimpleDyMFoam

ENH: pimpleFoam: Merged dynamic mesh functionality of pimpleDyMFoam into pimpleFoam

and replaced pimpleDyMFoam with a script which reports this change.

The pimpleDyMFoam tutorials have been moved into the pimpleFoam directory.

This change is the first of a set of developments to merge dynamic mesh
functionality into the standard solvers to improve consistency, usability,
flexibility and maintainability of these solvers.

Henry G. Weller
CFD Direct Ltd.

tutorials/incomp...
parent be9d388a
......@@ -4,14 +4,19 @@ EXE_INC = \
-I$(LIB_SRC)/transportModels \
-I$(LIB_SRC)/transportModels/incompressible/singlePhaseTransportModel \
-I$(LIB_SRC)/finiteVolume/lnInclude \
-I$(LIB_SRC)/sampling/lnInclude \
-I$(LIB_SRC)/dynamicFvMesh/lnInclude \
-I$(LIB_SRC)/dynamicMesh/lnInclude \
-I$(LIB_SRC)/meshTools/lnInclude \
-I$(LIB_SRC)/sampling/lnInclude
EXE_LIBS = \
-lturbulenceModels \
-lincompressibleTurbulenceModels \
-lincompressibleTransportModels \
-lfiniteVolume \
-lmeshTools \
-lfvOptions \
-lsampling
-lsampling \
-ldynamicFvMesh \
-ltopoChangerFvMesh \
-ldynamicMesh \
-lmeshTools
......@@ -4,12 +4,17 @@ surfaceScalarField phiHbyA
(
"phiHbyA",
fvc::flux(HbyA)
+ fvc::interpolate(rAU)*fvc::ddtCorr(U, phi)
+ fvc::interpolate(rAU)*fvc::ddtCorr(U, phi, Uf)
);
MRF.makeRelative(phiHbyA);
adjustPhi(phiHbyA, U, p);
if (p.needReference())
{
fvc::makeRelative(phiHbyA, U);
adjustPhi(phiHbyA, U, p);
fvc::makeAbsolute(phiHbyA, U);
}
tmp<volScalarField> rAtU(rAU);
......@@ -32,7 +37,6 @@ constrainPressure(p, U, phiHbyA, rAtU(), MRF);
// Non-orthogonal pressure corrector loop
while (pimple.correctNonOrthogonal())
{
// Pressure corrector
fvScalarMatrix pEqn
(
fvm::laplacian(rAtU(), p) == fvc::div(phiHbyA)
......@@ -53,6 +57,12 @@ while (pimple.correctNonOrthogonal())
// Explicitly relax pressure for momentum corrector
p.relax();
U = HbyA - rAtU()*fvc::grad(p);
U = HbyA - rAtU*fvc::grad(p);
U.correctBoundaryConditions();
fvOptions.correct(U);
// Correct Uf if the mesh is moving
fvc::correctUf(Uf, U, phi);
// Make the fluxes relative to the mesh motion
fvc::makeRelative(phi, U);
EXE_INC = \
-I. \
-I.. \
-I$(LIB_SRC)/TurbulenceModels/turbulenceModels/lnInclude \
-I$(LIB_SRC)/TurbulenceModels/incompressible/lnInclude \
-I$(LIB_SRC)/transportModels \
-I$(LIB_SRC)/transportModels/incompressible/singlePhaseTransportModel \
-I$(LIB_SRC)/finiteVolume/lnInclude \
-I$(LIB_SRC)/sampling/lnInclude \
-I$(LIB_SRC)/dynamicFvMesh/lnInclude \
-I$(LIB_SRC)/dynamicMesh/lnInclude \
-I$(LIB_SRC)/meshTools/lnInclude \
EXE_LIBS = \
-lturbulenceModels \
-lincompressibleTurbulenceModels \
-lincompressibleTransportModels \
-lfiniteVolume \
-lfvOptions \
-lsampling \
-ldynamicFvMesh \
-ltopoChangerFvMesh \
-ldynamicMesh \
-lmeshTools
volScalarField rAU(1.0/UEqn.A());
volVectorField HbyA(constrainHbyA(rAU*UEqn.H(), U, p));
surfaceScalarField phiHbyA
(
"phiHbyA",
fvc::flux(HbyA)
+ fvc::interpolate(rAU)*fvc::ddtCorr(U, Uf)
);
MRF.makeRelative(phiHbyA);
if (p.needReference())
{
fvc::makeRelative(phiHbyA, U);
adjustPhi(phiHbyA, U, p);
fvc::makeAbsolute(phiHbyA, U);
}
tmp<volScalarField> rAtU(rAU);
if (pimple.consistent())
{
rAtU = 1.0/max(1.0/rAU - UEqn.H1(), 0.1/rAU);
phiHbyA +=
fvc::interpolate(rAtU() - rAU)*fvc::snGrad(p)*mesh.magSf();
HbyA -= (rAU - rAtU())*fvc::grad(p);
}
if (pimple.nCorrPISO() <= 1)
{
tUEqn.clear();
}
// Update the pressure BCs to ensure flux consistency
constrainPressure(p, U, phiHbyA, rAtU(), MRF);
// Non-orthogonal pressure corrector loop
while (pimple.correctNonOrthogonal())
{
fvScalarMatrix pEqn
(
fvm::laplacian(rAtU(), p) == fvc::div(phiHbyA)
);
pEqn.setReference(pRefCell, pRefValue);
pEqn.solve(mesh.solver(p.select(pimple.finalInnerIter())));
if (pimple.finalNonOrthogonalIter())
{
phi = phiHbyA - pEqn.flux();
}
}
#include "continuityErrs.H"
// Explicitly relax pressure for momentum corrector
p.relax();
U = HbyA - rAtU*fvc::grad(p);
U.correctBoundaryConditions();
fvOptions.correct(U);
{
Uf = fvc::interpolate(U);
surfaceVectorField n(mesh.Sf()/mesh.magSf());
Uf += n*(phi/mesh.magSf() - (n & Uf));
}
// Make the fluxes relative to the mesh motion
fvc::makeRelative(phi, U);
......@@ -22,14 +22,14 @@ License
along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
Application
pimpleFoam
pimpleDyMFoam.C
Group
grpIncompressibleSolvers
Description
Large time-step transient solver for incompressible, turbulent flow, using
the PIMPLE (merged PISO-SIMPLE) algorithm.
Transient solver for incompressible, turbulent flow of Newtonian fluids
on a moving mesh.
\heading Solver details
The solver uses the PIMPLE (merged PISO-SIMPLE) algorithm to solve the
......@@ -68,9 +68,11 @@ Description
\*---------------------------------------------------------------------------*/
#include "fvCFD.H"
#include "dynamicFvMesh.H"
#include "singlePhaseTransportModel.H"
#include "turbulentTransportModel.H"
#include "pimpleControl.H"
#include "CorrectPhi.H"
#include "fvOptions.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
......@@ -81,11 +83,13 @@ int main(int argc, char *argv[])
#include "setRootCase.H"
#include "createTime.H"
#include "createMesh.H"
#include "createControl.H"
#include "createTimeControls.H"
#include "createFields.H"
#include "createDynamicFvMesh.H"
#include "initContinuityErrs.H"
#include "createControls.H"
#include "createFields.H"
#include "createUfIfPresent.H"
#include "CourantNo.H"
#include "setInitialDeltaT.H"
turbulence->validate();
......@@ -95,7 +99,7 @@ int main(int argc, char *argv[])
while (runTime.run())
{
#include "readTimeControls.H"
#include "readControls.H"
#include "CourantNo.H"
#include "setDeltaT.H"
......@@ -103,6 +107,31 @@ int main(int argc, char *argv[])
Info<< "Time = " << runTime.timeName() << nl << endl;
mesh.update();
#include "updateUf.H"
if (mesh.changing())
{
MRF.update();
if (correctPhi)
{
// Calculate absolute flux from the mapped surface velocity
phi = mesh.Sf() & Uf();
#include "correctPhi.H"
// Make the flux relative to the mesh motion
fvc::makeRelative(phi, U);
}
if (checkMeshCourantNo)
{
#include "meshCourantNo.H"
}
}
// --- Pressure-velocity PIMPLE corrector loop
while (pimple.loop())
{
......
#!/bin/sh
#------------------------------------------------------------------------------
# ========= |
# \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
# \\ / O peration |
# \\ / A nd | Copyright (C) 2017 OpenFOAM Foundation
# \\/ 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/>.
#
# Script
# supersededByPostProcess
#
# Description
# Script to suggest using the new "postProcess" utility.
#
#------------------------------------------------------------------------------
DyMsolver=${0##*/}
solver=$(echo $DyMsolver | sed 's/DyM//')
echo "The dynamic mesh functionality in $DyMsolver has been merged into $solver"
echo "and the $DyMsolver tutorials moved into the $solver tutorials directory"
#------------------------------------------------------------------------------
#!/bin/sh
#------------------------------------------------------------------------------
# ========= |
# \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
# \\ / O peration |
# \\ / A nd | Copyright (C) 2017 OpenFOAM Foundation
# \\/ 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/>.
#
# Script
# supersededByPostProcess
#
# Description
# Script to suggest using the new "postProcess" utility.
#
#------------------------------------------------------------------------------
DyMsolver=${0##*/}
solver=$(echo $DyMsolver | sed 's/DyM//')
echo "The dynamic mesh functionality in $DyMsolver has been merged into $solver"
echo "and the $DyMsolver tutorials moved into the $solver tutorials directory"
#------------------------------------------------------------------------------
......@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2011 OpenFOAM Foundation
\\ / A nd | Copyright (C) 2011-2017 OpenFOAM Foundation
\\/ M anipulation | Copyright (C) 2018 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
......@@ -37,8 +37,6 @@ SourceFiles
#define dynamicFvMesh_H
#include "fvMesh.H"
#include "autoPtr.H"
#include "runTimeSelectionTables.H"
#define MOVING_MESH
......@@ -117,7 +115,9 @@ public:
// Selectors
//- Select null constructed
//- Select, construct and return the dynamicFvMesh
// If the constant/dynamicMeshDict does not exist
// a staticFvMesh is returned
static autoPtr<dynamicFvMesh> New(const IOobject& io);
......@@ -131,13 +131,16 @@ public:
virtual bool update() = 0;
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#include "staticFvMesh.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //
......@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2011-2015 OpenFOAM Foundation
\\ / A nd | Copyright (C) 2011-2017 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
......@@ -24,8 +24,6 @@ License
\*---------------------------------------------------------------------------*/
#include "dynamicFvMesh.H"
#include "Time.H"
#include "dlLibraryTable.H"
// * * * * * * * * * * * * * * * * Selectors * * * * * * * * * * * * * * * * //
......@@ -36,51 +34,56 @@ Foam::autoPtr<Foam::dynamicFvMesh> Foam::dynamicFvMesh::New(const IOobject& io)
// - defaultRegion (region0) gets loaded from constant, other ones
// get loaded from constant/<regionname>. Normally we'd use
// polyMesh::dbDir() but we haven't got a polyMesh yet ...
IOdictionary dict
IOobject dictHeader
(
IOobject
(
"dynamicMeshDict",
io.time().constant(),
(io.name() == polyMesh::defaultRegion ? "" : io.name()),
io.db(),
IOobject::MUST_READ_IF_MODIFIED,
IOobject::NO_WRITE,
false
)
);
const word dynamicFvMeshTypeName(dict.lookup("dynamicFvMesh"));
Info<< "Selecting dynamicFvMesh " << dynamicFvMeshTypeName << endl;
const_cast<Time&>(io.time()).libs().open
(
dict,
"dynamicFvMeshLibs",
IOobjectConstructorTablePtr_
"dynamicMeshDict",
io.time().constant(),
(io.name() == polyMesh::defaultRegion ? "" : io.name()),
io.db(),
IOobject::MUST_READ_IF_MODIFIED,
IOobject::NO_WRITE,
false
);
if (!IOobjectConstructorTablePtr_)
if (dictHeader.typeHeaderOk<IOdictionary>(true))
{
FatalErrorInFunction
<< "dynamicFvMesh table is empty"
<< exit(FatalError);
}
IOdictionary dict(dictHeader);
auto cstrIter = IOobjectConstructorTablePtr_->cfind(dynamicFvMeshTypeName);
const word dynamicFvMeshTypeName(dict.lookup("dynamicFvMesh"));
if (!cstrIter.found())
{
FatalErrorInFunction
<< "Unknown dynamicFvMesh type "
<< dynamicFvMeshTypeName << nl << nl
<< "Valid dynamicFvMesh types :" << endl
<< IOobjectConstructorTablePtr_->sortedToc()
<< exit(FatalError);
Info<< "Selecting dynamicFvMesh " << dynamicFvMeshTypeName << endl;
const_cast<Time&>(io.time()).libs().open
(
dict,
"dynamicFvMeshLibs",
IOobjectConstructorTablePtr_
);
if (!IOobjectConstructorTablePtr_)
{
FatalErrorInFunction
<< "dynamicFvMesh table is empty"
<< exit(FatalError);
}
IOobjectConstructorTable::iterator cstrIter =
IOobjectConstructorTablePtr_->find(dynamicFvMeshTypeName);
if (cstrIter == IOobjectConstructorTablePtr_->end())
{
FatalErrorInFunction
<< "Unknown dynamicFvMesh type "
<< dynamicFvMeshTypeName << nl << nl
<< "Valid dynamicFvMesh types are :" << endl
<< IOobjectConstructorTablePtr_->sortedToc()
<< exit(FatalError);
}
return autoPtr<dynamicFvMesh>(cstrIter()(io));
}
return autoPtr<dynamicFvMesh>(cstrIter()(io));
return autoPtr<dynamicFvMesh>(new staticFvMesh(io));
}
......
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2017 OpenFOAM Foundation
\\/ 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/>.
Global
createUf
Description
Creates and initialises the velocity field Uf if present.
\*---------------------------------------------------------------------------*/
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
autoPtr<surfaceVectorField> Uf;
IOobject UfHeader
(
"Uf",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
);
if (UfHeader.typeHeaderOk<surfaceVectorField>(true))
{
Info<< "Reading face velocity Uf\n" << endl;
Uf = new surfaceVectorField(UfHeader, mesh);
}
// ************************************************************************* //
......@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2011-2017 OpenFOAM Foundation
\\ / A nd | Copyright (C) 2017 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
......@@ -21,110 +21,35 @@ License
You should have received a copy of the GNU General Public License
along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
Application
pimpleDyMFoam.C
Group
grpIncompressibleSolvers grpMovingMeshSolvers
Global
updateUf
Description
Transient solver for incompressible, turbulent flow of Newtonian fluids
on a moving mesh.
Turbulence modelling is generic, i.e. laminar, RAS or LES may be selected.
Constructs the face velocity field Uf if not already constructed.
\*---------------------------------------------------------------------------*/
#include "fvCFD.H"
#include "dynamicFvMesh.H"
#include "singlePhaseTransportModel.H"
#include "turbulentTransportModel.H"
#include "pimpleControl.H"
#include "CorrectPhi.H"
#include "fvOptions.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
int main(int argc, char *argv[])
if (mesh.changing())
{
#include "postProcess.H"
#include "setRootCase.H"
#include "createTime.H"
#include "createDynamicFvMesh.H"
#include "initContinuityErrs.H"
#include "createControls.H"
#include "createFields.H"
#include "createUf.H"
#include "CourantNo.H"
#include "setInitialDeltaT.H"