Commit b98e0609 authored by Andrew Heather's avatar Andrew Heather
Browse files

Merge branch 'feature-ami' into 'develop'

AMI code enhancements

See merge request !367
parents b61dd6fd 6e716c66
volScalarField rAU(1.0/UEqn.A());
volVectorField HbyA(constrainHbyA(rAU*UEqn.H(), U, p));
surfaceScalarField phiHbyA
(
"phiHbyA",
fvc::flux(HbyA)
+ MRF.zeroFilter(fvc::interpolate(rAU)*fvc::ddtCorr(U, phi, Uf))
);
surfaceScalarField phiHbyA("phiHbyA", fvc::flux(HbyA));
if (pimple.ddtCorr())
{
phiHbyA += MRF.zeroFilter(fvc::interpolate(rAU)*fvc::ddtCorr(U, phi, Uf));
}
MRF.makeRelative(phiHbyA);
......
......@@ -252,11 +252,7 @@ int main(int argc, char *argv[])
meshToMesh::interpolationMethod method =
meshToMesh::interpolationMethodNames_[mapMethod];
patchMapMethod =
AMIPatchToPatchInterpolation::interpolationMethodNames_
[
meshToMesh::interpolationMethodAMI(method)
];
patchMapMethod = meshToMesh::interpolationMethodAMI(method);
}
word procMapMethod =
......
......@@ -45,6 +45,7 @@ void Foam::polyMesh::updateMesh(const mapPolyMesh& mpm)
<< "Updating addressing and (optional) pointMesh/pointFields" << endl;
// Update boundaryMesh (note that patches themselves already ok)
// boundary_.updateMesh(mpm);
boundary_.updateMesh();
// Update zones
......
......@@ -63,6 +63,7 @@ void Foam::polyPatch::movePoints(PstreamBuffers&, const pointField& p)
primitivePatch::movePoints(p);
}
void Foam::polyPatch::updateMesh(PstreamBuffers&)
{
primitivePatch::clearGeom();
......
......@@ -55,6 +55,7 @@ namespace Foam
// Forward declarations
class polyBoundaryMesh;
class polyPatch;
class polyTopoChange;
class PstreamBuffers;
Ostream& operator<<(Ostream&, const polyPatch&);
......@@ -419,6 +420,19 @@ public:
labelList& rotation
) const;
//- For dynamic mesh cases - return true if this patch will change the
//- topology
virtual bool changeTopology() const
{
return false;
}
//- Collect topology changes in a polyTopoChange object
virtual bool setTopology(polyTopoChange&)
{
return false;
}
// Member operators
......
......@@ -47,6 +47,7 @@ License
#include "gravityMeshObject.H"
#include "turbulentTransportModel.H"
#include "demandDrivenData.H"
#include "unitConversion.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
......
......@@ -10,4 +10,9 @@ dynamicMotionSolverListFvMesh/dynamicMotionSolverListFvMesh.C
simplifiedDynamicFvMesh/simplifiedDynamicFvMeshes.C
simplifiedDynamicFvMesh/simplifiedDynamicFvMesh.C
dynamicMotionSolverFvMeshAMI/dynamicMotionSolverFvMeshAMI.C
LIB = $(FOAM_LIBBIN)/libdynamicFvMesh
......@@ -92,6 +92,9 @@ const Foam::motionSolver& Foam::dynamicMotionSolverFvMesh::motion() const
bool Foam::dynamicMotionSolverFvMesh::update()
{
// Scan through AMI patches and update
fvMesh::movePoints(motionPtr_->newPoints());
volVectorField* Uptr = getObjectPtr<volVectorField>("U");
......
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | www.openfoam.com
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2019-2020 OpenCFD Ltd.
-------------------------------------------------------------------------------
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/>.
\*---------------------------------------------------------------------------*/
#include "dynamicMotionSolverFvMeshAMI.H"
#include "addToRunTimeSelectionTable.H"
#include "motionSolver.H"
#include "volFields.H"
#include "surfaceFields.H"
#include "cyclicAMIPolyPatch.H"
#include "polyTopoChange.H"
#include "MeshObject.H"
#include "lduMesh.H"
#include "surfaceInterpolate.H"
#include "processorFvPatch.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
namespace Foam
{
defineTypeNameAndDebug(dynamicMotionSolverFvMeshAMI, 0);
addToRunTimeSelectionTable
(
dynamicFvMesh,
dynamicMotionSolverFvMeshAMI,
IOobject
);
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
Foam::dynamicMotionSolverFvMeshAMI::dynamicMotionSolverFvMeshAMI
(
const IOobject& io
)
:
dynamicFvMesh(io),
motionPtr_(motionSolver::New(*this))
{}
Foam::dynamicMotionSolverFvMeshAMI::dynamicMotionSolverFvMeshAMI
(
const IOobject& io,
pointField&& points,
faceList&& faces,
labelList&& allOwner,
labelList&& allNeighbour,
const bool syncPar
)
:
dynamicFvMesh
(
io,
std::move(points),
std::move(faces),
std::move(allOwner),
std::move(allNeighbour),
syncPar
),
motionPtr_(motionSolver::New(*this))
{}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
const Foam::motionSolver& Foam::dynamicMotionSolverFvMeshAMI::motion() const
{
return *motionPtr_;
}
bool Foam::dynamicMotionSolverFvMeshAMI::update()
{
// Mesh not moved/changed yet
moving(false);
topoChanging(false);
if (debug)
{
for (const fvPatch& fvp : boundary())
{
if (!isA<processorFvPatch>(fvp))
{
Info<< "1 --- patch:" << fvp.patch().name()
<< " area:" << gSum(fvp.magSf()) << endl;
}
}
}
pointField newPoints(motionPtr_->curPoints());
polyBoundaryMesh& pbm = const_cast<polyBoundaryMesh&>(boundaryMesh());
// Scan all patches and see if we want to apply a mesh topology update
bool changeRequired = false;
for (polyPatch& pp : pbm)
{
DebugInfo
<< "pre-topology change: patch " << pp.name()
<< " size:" << returnReduce(pp.size(), sumOp<label>())
<< " mag(faceAreas):" << gSum(mag(pp.faceAreas())) << endl;
//changeRequired = pp.changeTopology(newPoints) || changeRequired;
changeRequired = pp.changeTopology() || changeRequired;
}
reduce(changeRequired, orOp<bool>());
if (changeRequired)
{
polyTopoChange polyTopo(*this);
// Set new point positions in polyTopo object
polyTopo.movePoints(newPoints);
// Accumulate the patch-based mesh changes on the current mesh
// Note:
// - updates the AMIs using the new points
// - creates a topo change object that removes old added faces and
// adds the new faces
for (polyPatch& pp : pbm)
{
pp.setTopology(polyTopo);
}
// Update geometry
// Note
// - changeMesh leads to polyMesh::resetPrimitives which will also
// trigger polyBoundaryMesh::updateMesh (init and update) and
// ::calcGeometry (with topoChanging = false)
// - BUT: mesh still corresponds to original (non-extended mesh) so
// we want to bypass these calls...
// - after changes topoChanging = true
autoPtr<mapPolyMesh> map =
polyTopo.changeMesh
(
*this,
true // We will be calling movePoints after this update
);
// Apply topology change - update fv geometry and map fields
// - polyMesh::updateMesh
// - fires initUpdateMesh and updateMesh in AMI BCs - called before
// mapFields
// - AMI addressing must be up-to-date - used by, e.g. FaceCellWave
// - will trigger (again) polyBoundaryMesh::updateMesh (init and update)
updateMesh(map());
// Move points and update derived properties
// Note:
// - resets face areas based on raw point locations!
// - polyBoundaryMesh::updateMesh (init and update)
// Note:
// - processorPolyPatches will trigger calculation of faceCentres
// (and therefore cell volumes), so need to update faceAreas in
// initMovePoints since proc patches will be evaluated later than
// AMI patches
if (map().hasMotionPoints())
{
movePoints(map().preMotionPoints());
}
}
else
{
fvMesh::movePoints(newPoints);
}
volVectorField* Uptr = getObjectPtr<volVectorField>("U");
if (Uptr)
{
Uptr->correctBoundaryConditions();
surfaceVectorField* UfPtr = getObjectPtr<surfaceVectorField>("Uf");
if (UfPtr)
{
*UfPtr = fvc::interpolate(*Uptr);
}
}
if (debug)
{
for (const fvPatch& fvp : boundary())
{
if (!isA<processorFvPatch>(fvp))
{
Info<< "2 --- patch:" << fvp.patch().name()
<< " area:" << gSum(fvp.magSf()) << endl;
}
}
}
return true;
}
// ************************************************************************* //
......@@ -5,7 +5,8 @@
\\ / A nd | www.openfoam.com
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2011-2012 OpenFOAM Foundation
Copyright (C) 2011-2017 OpenFOAM Foundation
Copyright (C) 2020 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
......@@ -23,16 +24,98 @@ License
You should have received a copy of the GNU General Public License
along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
Class
Foam::dynamicMotionSolverFvMeshAMI
Description
The dynamicMotionSolverFvMeshAMI
SourceFiles
dynamicMotionSolverFvMeshAMI.C
\*---------------------------------------------------------------------------*/
#include "AMIInterpolation.H"
#ifndef dynamicMotionSolverFvMeshAMI_H
#define dynamicMotionSolverFvMeshAMI_H
#include "dynamicFvMesh.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
defineTypeNameAndDebug(AMIInterpolationName, 0);
}
class motionSolver;
/*---------------------------------------------------------------------------*\
Class dynamicMotionSolverFvMeshAMI Declaration
\*---------------------------------------------------------------------------*/
class dynamicMotionSolverFvMeshAMI
:
public dynamicFvMesh
{
// Private data
autoPtr<motionSolver> motionPtr_;
// Private Member Functions
//- No copy construct
dynamicMotionSolverFvMeshAMI
(
const dynamicMotionSolverFvMeshAMI&
) = delete;
//- No copy assignment
void operator=(const dynamicMotionSolverFvMeshAMI&) = delete;
public:
//- Runtime type information
TypeName("dynamicMotionSolverFvMeshAMI");
// Constructors
//- Construct from IOobject
dynamicMotionSolverFvMeshAMI(const IOobject& io);
//- Construct from components without boundary.
// Boundary is added using addFvPatches() member function
dynamicMotionSolverFvMeshAMI
(
const IOobject& io,
pointField&& points,
faceList&& faces,
labelList&& allOwner,
labelList&& allNeighbour,
const bool syncPar = true
);
//- Destructor
virtual ~dynamicMotionSolverFvMeshAMI() = default;
// Member Functions
//- Return the motionSolver
const motionSolver& motion() const;
//- Update the mesh for both mesh motion and topology change
virtual bool update();
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //
EXE_INC = \
-I$(LIB_SRC)/fileFormats/lnInclude \
-I$(LIB_SRC)/surfMesh/lnInclude \
-I$(LIB_SRC)/meshTools/lnInclude
-I$(LIB_SRC)/meshTools/lnInclude \
-I$(LIB_SRC)/dynamicMesh/lnInclude
LIB_LIBS = \
-lOpenFOAM \
......
......@@ -52,6 +52,7 @@ bool Foam::pimpleControl::read()
pimpleDict.getOrDefault("turbOnFinalIterOnly", true);
finalOnLastPimpleIterOnly_ =
pimpleDict.getOrDefault("finalOnLastPimpleIterOnly", false);
ddtCorr_ = pimpleDict.getOrDefault("ddtCorr", true);
return true;
}
......@@ -155,8 +156,9 @@ Foam::pimpleControl::pimpleControl
corrPISO_(0),
SIMPLErho_(false),
turbOnFinalIterOnly_(true),
converged_(false),
finalOnLastPimpleIterOnly_(false)
finalOnLastPimpleIterOnly_(false),
ddtCorr_(true),
converged_(false)
{
read();
......
......@@ -85,19 +85,22 @@ protected:
label corrPISO_;
//- Flag to indicate whether to update density in SIMPLE
// rather than PISO mode
//- rather than PISO mode
bool SIMPLErho_;
//- Flag to indicate whether to only solve turbulence on final iter
bool turbOnFinalIterOnly_;
//- Converged flag
bool converged_;
//- Flag to indicate wheter the final solver is used only on the
// final pimple iter
//- final pimple iter
bool finalOnLastPimpleIterOnly_;
//- Flag to indicate that ddtCorr should be applied; default = yes
bool ddtCorr_;
//- Converged flag
bool converged_;
// Protected Member Functions
......@@ -181,6 +184,9 @@ public:
//- Return true to solve for turbulence
inline bool turbCorr();
//- Return true to apply ddtCorr
inline bool ddtCorr() const;
};
......
......@@ -143,4 +143,10 @@ inline bool Foam::pimpleControl::turbCorr()
}
inline bool Foam::pimpleControl::ddtCorr() const
{
return ddtCorr_;
}
// ************************************************************************* //
......@@ -226,10 +226,9 @@ void Foam::cyclicACMIFvPatchField<Type>::updateInterfaceMatrix
{
const cyclicACMIPolyPatch& cpp = cyclicACMIPatch_.cyclicACMIPatch();
// note: only applying coupled contribution
// Note: only applying coupled contribution
const labelUList& nbrFaceCellsCoupled =
cpp.neighbPatch().faceCells();
const labelUList& nbrFaceCellsCoupled = cpp.neighbPatch().faceCells();
solveScalarField pnf(psiInternal, nbrFaceCellsCoupled);
......@@ -254,7 +253,7 @@ void Foam::cyclicACMIFvPatchField<Type>::updateInterfaceMatrix
{
const cyclicACMIPolyPatch& cpp = cyclicACMIPatch_.cyclicACMIPatch();
// note: only applying coupled contribution
// Note: only applying coupled contribution
const labelUList& nbrFaceCellsCoupled = cpp.neighbPatch().faceCells();
......@@ -277,7 +276,7 @@ void Foam::cyclicACMIFvPatchField<Type>::manipulateMatrix
{
const scalarField& mask = cyclicACMIPatch_.cyclicACMIPatch().mask();
// nothing to be done by the AMI, but re-direct to non-overlap patch
// Nothing to be done by the AMI, but re-direct to non-overlap patch
// with non-overlap patch weights
const fvPatchField<Type>& npf = nonOverlapPatchField();
......
......@@ -173,7 +173,7 @@ public:
//- Return true if this patch field fixes a value
// Needed to check if a level has to be specified while solving
// Poissons equations
// Poisson equations
virtual bool fixesValue() const
{
const scalarField& mask =
......
......@@ -726,6 +726,8 @@ void Foam::fvMesh::mapFields(const mapPolyMesh& meshMap)
Foam::tmp<Foam::scalarField> Foam::fvMesh::movePoints(const pointField& p)
{