Merge branch 'feature/OversetMesh' into CumulativeDevelopment-VukoVukcevic

This commit is contained in:
Vuko Vukcevic 2018-02-28 13:42:22 +01:00
commit 3040021a85
751 changed files with 789307 additions and 0 deletions

View file

@ -0,0 +1,3 @@
icoDyMOversetFoam.C
EXE = $(FOAM_APPBIN)/icoDyMOversetFoam

View file

@ -0,0 +1,14 @@
EXE_INC = \
-I$(LIB_SRC)/finiteVolume/lnInclude \
-I$(LIB_SRC)/meshTools/lnInclude \
-I$(LIB_SRC)/dynamicMesh/dynamicFvMesh/lnInclude \
-I$(LIB_SRC)/dynamicMesh/dynamicMesh/lnInclude \
-I$(LIB_SRC)/oversetMesh/oversetMesh/lnInclude
EXE_LIBS = \
-lfiniteVolume \
-lmeshTools \
-ldynamicFvMesh \
-ldynamicMesh \
-loversetMesh \
-llduSolvers

View file

@ -0,0 +1,12 @@
// Momentum equation
fvVectorMatrix UEqn
(
fvm::ddt(U)
+ fvm::div(phi, U)
- fvm::laplacian(nu, U)
);
if (piso.momentumPredictor())
{
solve(UEqn == -fvc::grad(p));
}

View file

@ -0,0 +1,57 @@
{
wordList pcorrTypes(p.boundaryField().types());
for (label i = 0; i < p.boundaryField().size(); i++)
{
if (p.boundaryField()[i].fixesValue())
{
pcorrTypes[i] = fixedValueFvPatchScalarField::typeName;
}
}
volScalarField pcorr
(
IOobject
(
"pcorr",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::NO_WRITE
),
mesh,
dimensionedScalar("pcorr", p.dimensions(), 0),
pcorrTypes
);
// Initialise flux with interpolated velocity
phi = faceOversetMask*(fvc::interpolate(U) & mesh.Sf());
// Adjust overset fluxes
oversetAdjustPhi(phi, U); // Fringe flux adjustment
globalOversetAdjustPhi(phi, U, pcorr); // Global flux adjustment
while (piso.correctNonOrthogonal())
{
fvScalarMatrix pcorrEqn
(
fvm::laplacian(rAU, pcorr) == fvc::div(phi)
);
// Adjust non-orthogonal fluxes if necessary
om.correctNonOrthoFluxes(pcorrEqn, U);
pcorrEqn.setReference(pRefCell, pRefValue);
pcorrEqn.solve();
if (piso.finalNonOrthogonalIter())
{
phi -= pcorrEqn.flux();
}
// Fluxes are corrected to absolute velocity and further corrected
// later. HJ, 6/Feb/2009
}
# include "oversetContinuityErrs.H"
}

View file

@ -0,0 +1,11 @@
#include "createTimeControls.H"
bool correctPhi
(
piso.dict().lookupOrDefault("correctPhi", false)
);
bool checkMeshCourantNo
(
piso.dict().lookupOrDefault("checkMeshCourantNo", false)
);

View file

@ -0,0 +1,79 @@
Info<< "Reading transportProperties\n" << endl;
IOdictionary transportProperties
(
IOobject
(
"transportProperties",
runTime.constant(),
mesh,
IOobject::MUST_READ,
IOobject::NO_WRITE
)
);
dimensionedScalar nu
(
transportProperties.lookup("nu")
);
Info<< "Reading field p\n" << endl;
volScalarField p
(
IOobject
(
"p",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);
Info<< "Reading field U\n" << endl;
volVectorField U
(
IOobject
(
"U",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);
# include "createPhi.H"
// Set reference cell taking into account whether implicit or explicit
// overset aglorithm is used for field p.
const oversetMesh& om = oversetMesh::New(mesh);
label pRefCell = 0;
scalar pRefValue = 0.0;
om.setRefCell(p, piso.dict(), pRefCell, pRefValue);
mesh.schemesDict().setFluxRequired(p.name());
if (correctPhi)
{
mesh.schemesDict().setFluxRequired("pcorr");
}
Info<< "Reading field rAU if present\n" << endl;
volScalarField rAU
(
IOobject
(
"rAU",
runTime.timeName(),
mesh,
IOobject::READ_IF_PRESENT,
IOobject::AUTO_WRITE
),
mesh,
dimensionedScalar("dt", dimTime, 1.0),
zeroGradientFvPatchScalarField::typeName
);
# include "writeOversetMasks.H"

View file

@ -0,0 +1,167 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | foam-extend: Open Source CFD
\\ / O peration | Version: 4.0
\\ / A nd | Web: http://www.foam-extend.org
\\/ M anipulation | For copyright notice see file Copyright
-------------------------------------------------------------------------------
License
This file is part of foam-extend.
foam-extend 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.
foam-extend 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 foam-extend. If not, see <http://www.gnu.org/licenses/>.
Application
icoDyMOversetFoam
Description
Transient solver for incompressible, laminar flow of Newtonian fluids
with dynamic mesh and overset mesh support.
Author
Hrvoje Jasak, Wikki Ltd. All rights reserved
\*---------------------------------------------------------------------------*/
#include "fvCFD.H"
#include "dynamicFvMesh.H"
#include "oversetMesh.H"
#include "oversetFvPatchFields.H"
#include "oversetAdjustPhi.H"
#include "globalOversetAdjustPhi.H"
#include "pisoControl.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
int main(int argc, char *argv[])
{
# include "setRootCase.H"
# include "createTime.H"
# include "createDynamicFvMesh.H"
pisoControl piso(mesh);
# include "initContinuityErrs.H"
# include "initTotalVolume.H"
# include "createControls.H"
# include "createFields.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Info<< "\nStarting time loop\n" << endl;
while (runTime.run())
{
# include "readControls.H"
# include "checkTotalVolume.H"
// Make the fluxes absolute
fvc::makeAbsolute(phi, U);
runTime++;
Info<< "Time = " << runTime.timeName() << nl << endl;
bool meshChanged = mesh.update();
reduce(meshChanged, orOp<bool>());
# include "createOversetMasks.H"
if (correctPhi && (mesh.moving() || meshChanged))
{
# include "correctPhi.H"
}
// Make the fluxes relative to the mesh motion
fvc::makeRelative(phi, U);
# include "oversetCourantNo.H"
# include "setDeltaT.H"
if (mesh.moving() && checkMeshCourantNo)
{
# include "meshCourantNo.H"
}
# include "UEqn.H"
// --- PISO loop
while (piso.correct())
{
p.boundaryField().updateCoeffs();
rAU = 1.0/UEqn.A();
oversetFvPatchScalarField::oversetInterpolate(rAU); // Overset update
U = rAU*UEqn.H();
oversetFvPatchVectorField::oversetInterpolate(U); // Overset update
phi = fvc::interpolate(U) & mesh.Sf();
// Adjust overset fluxes
oversetAdjustPhi(phi, U); // Fringe flux adjustment
globalOversetAdjustPhi(phi, U, p); // Global flux adjustment
while (piso.correctNonOrthogonal())
{
fvScalarMatrix pEqn
(
fvm::laplacian(rAU, p) == fvc::div(phi)
);
// Adjust non-orthogonal fringe fluxes if necessary
om.correctNonOrthoFluxes(pEqn, U);
pEqn.setReference(pRefCell, pRefValue);
pEqn.solve
(
mesh.solutionDict().solver(p.select(piso.finalInnerIter()))
);
if (piso.finalNonOrthogonalIter())
{
phi -= pEqn.flux();
}
// Perform overset interpolation (after flux reconstruction)
oversetFvPatchScalarField::oversetInterpolate(p);
}
# include "oversetContinuityErrs.H"
U -= rAU*fvc::grad(p);
U.correctBoundaryConditions();
oversetFvPatchVectorField::oversetInterpolate(U); // Overset update
// Note: if implicit fringe is switched on, we are doing the
// interpolation twice (once in correctBoundaryConditions and once
// in oversetInterpolate). Reorganize. VV, 4/Oct/2016.
// Make the fluxes relative to the mesh motion
fvc::makeRelative(phi, U);
}
runTime.write();
Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
<< " ClockTime = " << runTime.elapsedClockTime() << " s"
<< nl << endl;
}
Info<< "End\n" << endl;
return(0);
}
// ************************************************************************* //

View file

@ -0,0 +1,5 @@
#include "readTimeControls.H"
correctPhi = piso.dict().lookupOrDefault("correctPhi", false);
checkMeshCourantNo = piso.dict().lookupOrDefault("checkMeshCourantNo", false);

View file

@ -0,0 +1,3 @@
icoOversetFoam.C
EXE = $(FOAM_APPBIN)/icoOversetFoam

View file

@ -0,0 +1,10 @@
EXE_INC = \
-I$(LIB_SRC)/finiteVolume/lnInclude \
-I$(LIB_SRC)/triSurface/lnInclude \
-I$(LIB_SRC)/meshTools/lnInclude \
-I$(LIB_SRC)/oversetMesh/oversetMesh/lnInclude
EXE_LIBS = \
-lfiniteVolume \
-loversetMesh \
-llduSolvers

View file

@ -0,0 +1,74 @@
Info<< "Reading transportProperties\n" << endl;
IOdictionary transportProperties
(
IOobject
(
"transportProperties",
runTime.constant(),
mesh,
IOobject::MUST_READ,
IOobject::NO_WRITE
)
);
dimensionedScalar nu
(
transportProperties.lookup("nu")
);
Info<< "Reading field p\n" << endl;
volScalarField p
(
IOobject
(
"p",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);
Info<< "Reading field U\n" << endl;
volVectorField U
(
IOobject
(
"U",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);
# include "createPhi.H"
// Set reference cell taking into account whether implicit or explicit
// overset aglorithm is used for field p.
label pRefCell = 0;
scalar pRefValue = 0.0;
om.setRefCell(p, piso.dict(), pRefCell, pRefValue);
mesh.schemesDict().setFluxRequired(p.name());
Info<< "Reading field rAU if present\n" << endl;
volScalarField rAU
(
IOobject
(
"rAU",
runTime.timeName(),
mesh,
IOobject::READ_IF_PRESENT,
IOobject::AUTO_WRITE
),
mesh,
dimensionedScalar("dt", dimTime, 1.0),
zeroGradientFvPatchScalarField::typeName
);
# include "writeOversetMasks.H"

View file

@ -0,0 +1,150 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | foam-extend: Open Source CFD
\\ / O peration | Version: 4.0
\\ / A nd | Web: http://www.foam-extend.org
\\/ M anipulation | For copyright notice see file Copyright
-------------------------------------------------------------------------------
License
This file is part of foam-extend.
foam-extend 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.
foam-extend 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 foam-extend. If not, see <http://www.gnu.org/licenses/>.
Application
icoOversetFoam
Description
Transient solver for incompressible, laminar flow of Newtonian fluids
with overset mesh support.
Author
Hrvoje Jasak, Wikki Ltd. All rights reserved
\*---------------------------------------------------------------------------*/
#include "fvCFD.H"
#include "pisoControl.H"
#include "oversetMesh.H"
#include "oversetFvPatchFields.H"
#include "oversetAdjustPhi.H"
#include "globalOversetAdjustPhi.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
int main(int argc, char *argv[])
{
# include "setRootCase.H"
# include "createTime.H"
# include "createMesh.H"
pisoControl piso(mesh);
# include "createOversetMasks.H"
# include "createFields.H"
# include "initContinuityErrs.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Info<< "\nStarting time loop\n" << endl;
while (runTime.loop())
{
Info<< "Time = " << runTime.timeName() << nl << endl;
# include "oversetCourantNo.H"
// Momentum equation
fvVectorMatrix UEqn
(
fvm::ddt(U)
+ fvm::div(phi, U)
- fvm::laplacian(nu, U)
);
if (piso.momentumPredictor())
{
solve(UEqn == -fvc::grad(p));
}
// --- PISO loop
while (piso.correct())
{
p.boundaryField().updateCoeffs();
rAU = 1.0/UEqn.A();
oversetFvPatchScalarField::oversetInterpolate(rAU); // Overset update
U = rAU*UEqn.H();
oversetFvPatchVectorField::oversetInterpolate(U); // Overset update
phi = fvc::interpolate(U) & mesh.Sf();
// Adjust overset fluxes
oversetAdjustPhi(phi, U); // Fringe flux adjustment
globalOversetAdjustPhi(phi, U, p); // Global flux adjustment
while (piso.correctNonOrthogonal())
{
fvScalarMatrix pEqn
(
fvm::laplacian(rAU, p) == fvc::div(phi)
);
// Adjust non-orthogonal fringe fluxes if necessary
om.correctNonOrthoFluxes(pEqn, U);
pEqn.setReference(pRefCell, pRefValue);
pEqn.solve
(
mesh.solutionDict().solver(p.select(piso.finalInnerIter()))
);
if (piso.finalNonOrthogonalIter())
{
phi -= pEqn.flux();
}
// Perform overset interpolation (after flux reconstruction)
oversetFvPatchScalarField::oversetInterpolate(p);
}
# include "oversetContinuityErrs.H"
U -= rAU*fvc::grad(p);
U.correctBoundaryConditions();
oversetFvPatchVectorField::oversetInterpolate(U); // Overset update
// Note: if implicit fringe is switched on, we are doing the
// interpolation twice (once in correctBoundaryConditions and once
// in oversetInterpolate). Reorganize. VV, 4/Oct/2016.
// Make the fluxes relative to the mesh motion
fvc::makeRelative(phi, U);
}
runTime.write();
Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
<< " ClockTime = " << runTime.elapsedClockTime() << " s"
<< nl << endl;
}
Info<< "End\n" << endl;
return(0);
}
// ************************************************************************* //

View file

@ -0,0 +1,3 @@
interOversetFoam.C
EXE = $(FOAM_APPBIN)/interOversetFoam

View file

@ -0,0 +1,21 @@
EXE_INC = \
-I$(LIB_SRC)/turbulenceModels \
-I$(LIB_SRC)/turbulenceModels/incompressible/RAS/RASModel \
-I$(LIB_SRC)/transportModels \
-I$(LIB_SRC)/transportModels/incompressible/lnInclude \
-I$(LIB_SRC)/transportModels/interfaceProperties/lnInclude \
-I$(LIB_SRC)/turbulenceModels/incompressible/turbulenceModel \
-I$(LIB_SRC)/finiteVolume/lnInclude \
-I$(LIB_SRC)/triSurface/lnInclude \
-I$(LIB_SRC)/meshTools/lnInclude \
-I$(LIB_SRC)/oversetMesh/oversetMesh/lnInclude
EXE_LIBS = \
-linterfaceProperties \
-lincompressibleTransportModels \
-lincompressibleTurbulenceModel \
-lincompressibleRASModels \
-lincompressibleLESModels \
-lfiniteVolume \
-loversetMesh \
-llduSolvers

View file

@ -0,0 +1,33 @@
surfaceScalarField muEff
(
"muEff",
twoPhaseProperties.muf()
+ fvc::interpolate(rho*turbulence->nut())
);
fvVectorMatrix UEqn
(
fvm::ddt(rho, U)
+ fvm::div(rhoPhi, U)
- fvm::laplacian(muEff, U)
- (fvc::grad(U) & fvc::grad(muEff))
);
UEqn.relax();
if (pimple.momentumPredictor())
{
solve
(
UEqn
==
fvc::reconstruct
(
(
fvc::interpolate(interface.sigmaK())*fvc::snGrad(alpha1)
- ghf*fvc::snGrad(rho)
- fvc::snGrad(pd)
) * mesh.magSf()
)
);
}

View file

@ -0,0 +1,60 @@
{
word alphaScheme("div(phi,alpha)");
word alpharScheme("div(phirb,alpha)");
// New formulation of phir: Co condition
// Note faceOversetMask. HJ, 16/Jun/2015
surfaceScalarField phir
(
"phir",
faceOversetMask*interface.cAlpha()*interface.nHatf()*
min
(
scalar(0.5)/ // This is ref compression Co number for cAlpha = 1
(
mesh.surfaceInterpolation::deltaCoeffs()*
mesh.time().deltaT()
),
mag(phi/mesh.magSf())
+ dimensionedScalar("small", dimVelocity, 1e-3)
)
);
fvScalarMatrix alpha1Eqn
(
fvm::ddt(alpha1)
+ fvm::div(phi, alpha1, alphaScheme)
+ fvm::div
(
-fvc::flux(-phir, scalar(1) - alpha1, alpharScheme),
alpha1,
alpharScheme
)
);
alpha1Eqn.relax();
alpha1Eqn.solve();
// Update mass flux
rhoPhi = alpha1Eqn.flux()*(rho1 - rho2) + phi*rho2;
Info<< "Liquid phase volume fraction = "
<< alpha1.weightedAverage(mesh.V()).value()
<< " Min(alpha1) = " << min(cellOversetMask*alpha1).value()
<< " Max(alpha1) = " << max(cellOversetMask*alpha1).value()
<< endl;
// Limit alpha to handle overset cells
alpha1.max(0);
alpha1.min(1);
// Update of interface and density
interface.correct();
twoPhaseProperties.correct();
// Calculate density using limited alpha1
rho = twoPhaseProperties.rho();
// Parallel update in snGrad. HJ, 8/Dec/2012
rho.correctBoundaryConditions();
}

View file

@ -0,0 +1,60 @@
{
// Set flux required
mesh.schemesDict().setFluxRequired("pcorr");
wordList pcorrTypes(p.boundaryField().types());
for (label i = 0; i < p.boundaryField().size(); i++)
{
if (p.boundaryField()[i].fixesValue())
{
pcorrTypes[i] = fixedValueFvPatchScalarField::typeName;
}
}
volScalarField pcorr
(
IOobject
(
"pcorr",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::NO_WRITE
),
mesh,
dimensionedScalar("pcorr", p.dimensions(), 0),
pcorrTypes
);
// Initialise flux with interpolated velocity
phi = faceOversetMask*(fvc::interpolate(U) & mesh.Sf());
// Adjust overset fluxes
oversetAdjustPhi(phi, U); // Fringe flux adjustment
globalOversetAdjustPhi(phi, U, pcorr); // Global flux adjustment
while (pimple.correctNonOrthogonal())
{
fvScalarMatrix pcorrEqn
(
fvm::laplacian(rAU, pcorr) == fvc::div(phi)
);
// Adjust non-orthogonal fluxes if necessary
om.correctNonOrthoFluxes(pcorrEqn, U);
pcorrEqn.setReference(pRefCell, pRefValue);
pcorrEqn.solve();
if (pimple.finalNonOrthogonalIter())
{
phi -= pcorrEqn.flux();
}
// Fluxes are corrected to absolute velocity and further corrected
// later. HJ, 6/Feb/2009
}
# include "oversetContinuityErrs.H"
}

View file

@ -0,0 +1,137 @@
Info<< "Reading field pd\n" << endl;
volScalarField pd
(
IOobject
(
"pd",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);
Info<< "Reading field alpha1\n" << endl;
volScalarField alpha1
(
IOobject
(
"alpha1",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);
Info<< "Reading field U\n" << endl;
volVectorField U
(
IOobject
(
"U",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);
# include "createPhi.H"
Info<< "Reading transportProperties\n" << endl;
twoPhaseMixture twoPhaseProperties(U, phi, "alpha1");
const dimensionedScalar& rho1 = twoPhaseProperties.rho1();
const dimensionedScalar& rho2 = twoPhaseProperties.rho2();
// Need to store rho for ddt(rho, U)
volScalarField rho
(
IOobject
(
"rho",
runTime.timeName(),
mesh,
IOobject::READ_IF_PRESENT
),
alpha1*rho1 + (scalar(1) - alpha1)*rho2,
alpha1.boundaryField().types()
);
rho.oldTime();
// Mass flux
// Initialisation does not matter because rhoPhi is reset after the
// alpha1 solution before it is used in the U equation.
surfaceScalarField rhoPhi
(
IOobject
(
"rho*phi",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::NO_WRITE
),
rho1*phi
);
Info<< "Calculating field g.h\n" << endl;
volScalarField gh("gh", g & mesh.C());
surfaceScalarField ghf("gh", g & mesh.Cf());
volScalarField p
(
IOobject
(
"p",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
pd + rho*gh
);
// Set reference cell taking into account whether implicit or explicit
// overset aglorithm is used for field p.
label pRefCell = 0;
scalar pRefValue = 0.0;
om.setRefCell(p, pimple.dict(), pRefCell, pRefValue);
mesh.schemesDict().setFluxRequired(pd.name());
mesh.schemesDict().setFluxRequired(alpha1.name());
// Construct interface from alpha1 distribution
interfaceProperties interface(alpha1, U, twoPhaseProperties);
// Construct incompressible turbulence model
autoPtr<incompressible::turbulenceModel> turbulence
(
incompressible::turbulenceModel::New(U, phi, twoPhaseProperties)
);
Info<< "Reading field rAU if present\n" << endl;
volScalarField rAU
(
IOobject
(
"rAU",
runTime.timeName(),
mesh,
IOobject::READ_IF_PRESENT,
IOobject::AUTO_WRITE
),
mesh,
runTime.deltaT()/rho1,
zeroGradientFvPatchScalarField::typeName
);
# include "writeOversetMasks.H"

View file

@ -0,0 +1,125 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | foam-extend: Open Source CFD
\\ / O peration | Version: 4.0
\\ / A nd | Web: http://www.foam-extend.org
\\/ M anipulation | For copyright notice see file Copyright
-------------------------------------------------------------------------------
License
This file is part of foam-extend.
foam-extend 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.
foam-extend 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 foam-extend. If not, see <http://www.gnu.org/licenses/>.
Application
interOversetFoam
Description
Solver for 2 incompressible, isothermal immiscible fluids using a VOF
(volume of fluid) phase-fraction based interface capturing approach,
with overset mesh support
The momentum and other fluid properties are of the "mixture" and a single
momentum equation is solved.
Turbulence modelling is generic, i.e. laminar, RAS or LES may be selected.
For a two-fluid approach see twoPhaseEulerFoam.
Author
Hrvoje Jasak, Wikki Ltd. All rights reserved
\*---------------------------------------------------------------------------*/
#include "fvCFD.H"
#include "interfaceProperties.H"
#include "twoPhaseMixture.H"
#include "turbulenceModel.H"
#include "pimpleControl.H"
#include "oversetMesh.H"
#include "oversetFvPatchFields.H"
#include "oversetAdjustPhi.H"
#include "globalOversetAdjustPhi.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
int main(int argc, char *argv[])
{
# include "setRootCase.H"
# include "createTime.H"
# include "createMesh.H"
pimpleControl pimple(mesh);
# include "readGravitationalAcceleration.H"
# include "initContinuityErrs.H"
# include "createOversetMasks.H"
# include "createFields.H"
# include "createTimeControls.H"
# include "correctPhi.H"
# include "CourantNo.H"
# include "setInitialDeltaT.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Info<< "\nStarting time loop\n" << endl;
while (runTime.run())
{
# include "readTimeControls.H"
# include "oversetCourantNo.H"
# include "setDeltaT.H"
runTime++;
Info<< "Time = " << runTime.timeName() << nl << endl;
// Pressure-velocity corrector
while (pimple.loop())
{
twoPhaseProperties.correct();
# include "alphaEqn.H"
# include "UEqn.H"
// --- PISO loop
while (pimple.correct())
{
# include "pdEqn.H"
}
# include "oversetContinuityErrs.H"
# include "limitU.H"
p = pd + cellOversetMask*rho*gh;
turbulence->correct();
}
runTime.write();
Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
<< " ClockTime = " << runTime.elapsedClockTime() << " s"
<< nl << endl;
}
Info<< "End\n" << endl;
return 0;
}
// ************************************************************************* //

View file

@ -0,0 +1,29 @@
{
scalar limitMagU = readScalar(pimple.dict().lookup("limitMagU"));
volScalarField magU(mag(U));
scalar maxMagU = max(magU).value();
Info<< "mag(U): max: " << maxMagU
<< " avg: " << magU.weightedAverage(mesh.V()).value();
if (maxMagU > limitMagU)
{
U.internalField() *=
neg(magU.internalField() - limitMagU)
+ pos(magU.internalField() - limitMagU)*
limitMagU/(magU.internalField() + SMALL);
U.correctBoundaryConditions();
oversetFvPatchVectorField::oversetInterpolate(U); // Overset update
// Note: if implicit fringe is switched on, we are doing the
// interpolation twice (once in correctBoundaryConditions and once in
// oversetInterpolate). Reorganize. VV, 4/Oct/2016.
Info << " ...limiting" << endl;
}
else
{
Info << endl;
}
}

View file

@ -0,0 +1,65 @@
{
pd.boundaryField().updateCoeffs();
rAU = 1.0/UEqn.A();
oversetFvPatchScalarField::oversetInterpolate(rAU); // Overset update
surfaceScalarField rAUf = fvc::interpolate(rAU);
U = rAU*UEqn.H();
oversetFvPatchVectorField::oversetInterpolate(U); // Overset update
surfaceScalarField phiU
(
"phiU",
fvc::interpolate(U) & mesh.Sf()
);
phi = phiU
+ (
fvc::interpolate(interface.sigmaK())*fvc::snGrad(alpha1)
- ghf*fvc::snGrad(rho)
)*rAUf*mesh.magSf();
// Adjust overset fluxes
oversetAdjustPhi(phi, U); // Fringe flux adjustment
globalOversetAdjustPhi(phi, U, p); // Global flux adjustment
while (pimple.correctNonOrthogonal())
{
fvScalarMatrix pdEqn
(
fvm::laplacian(rAUf, pd) == fvc::div(phi)
);
// Adjust non-orthogonal fringe fluxes if necessary
om.correctNonOrthoFluxes(pdEqn, U);
pdEqn.setReference(pRefCell, pRefValue);
pdEqn.solve
(
mesh.solutionDict().solver(pd.select(pimple.finalInnerIter()))
);
if (pimple.finalNonOrthogonalIter())
{
phi -= pdEqn.flux();
}
// Perform overset interpolation (after flux reconstruction)
oversetFvPatchScalarField::oversetInterpolate(p);
}
// Explicitly relax pressure for momentum corrector
if (!pimple.finalIter())
{
p.relax();
}
U += rAU*fvc::reconstruct((phi - phiU)/rAUf);
U.correctBoundaryConditions();
oversetFvPatchVectorField::oversetInterpolate(U); // Overset update
// Note: if implicit fringe is switched on, we are doing the interpolation
// twice (once in correctBoundaryConditions and once in oversetInterpolate)
// Reorganize. VV, 4/Oct/2016.
}

View file

@ -0,0 +1,3 @@
laplacianOversetFoam.C
EXE = $(FOAM_APPBIN)/laplacianOversetFoam

View file

@ -0,0 +1,8 @@
EXE_INC = \
-I$(LIB_SRC)/finiteVolume/lnInclude \
-I$(LIB_SRC)/oversetMesh/oversetMesh/lnInclude
EXE_LIBS = \
-lfiniteVolume \
-loversetMesh \
-llduSolvers

View file

@ -0,0 +1,37 @@
Info<< "Reading field T\n" << endl;
volScalarField T
(
IOobject
(
"T",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);
Info<< "Reading transportProperties\n" << endl;
IOdictionary transportProperties
(
IOobject
(
"transportProperties",
runTime.constant(),
mesh,
IOobject::MUST_READ,
IOobject::NO_WRITE
)
);
Info<< "Reading diffusivity DT\n" << endl;
dimensionedScalar DT
(
transportProperties.lookup("DT")
);

View file

@ -0,0 +1,99 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | foam-extend: Open Source CFD
\\ / O peration | Version: 4.0
\\ / A nd | Web: http://www.foam-extend.org
\\/ M anipulation | For copyright notice see file Copyright
-------------------------------------------------------------------------------
License
This file is part of foam-extend.
foam-extend 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.
foam-extend 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 foam-extend. If not, see <http://www.gnu.org/licenses/>.
Application
laplaciantOversetFoam
Description
Solves a simple Laplace equation, e.g. for thermal diffusion in a solid
with support for overset meshes.
Experimental: no overset mesh changes required
Author
Hrvoje Jasak, Wikki Ltd. All rights reserved
\*---------------------------------------------------------------------------*/
#include "fvCFD.H"
#include "simpleControl.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
int main(int argc, char *argv[])
{
# include "setRootCase.H"
# include "createTime.H"
# include "createMesh.H"
simpleControl simple(mesh);
# include "createFields.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Info<< "\nCalculating temperature distribution\n" << endl;
while (simple.loop())
{
Info<< "Time = " << runTime.timeName() << nl << endl;
T.correctBoundaryConditions();
Pout<< "T: " << T.internalField() << endl;
return 0;
fvScalarMatrix TEqn
(
fvm::laplacian(DT, T)
);
Pout<< "TEqn: " << TEqn << endl;
// TEqn.solve();
volScalarField residual
(
"residual",
T
);
residual.internalField() = TEqn.residual();
// residual.boundaryField() == 0;
residual.write();
Info<< "residual " << gSumMag(residual.internalField()) << endl;
# include "write.H"
Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
<< " ClockTime = " << runTime.elapsedClockTime() << " s"
<< nl << endl;
}
Info<< "End\n" << endl;
return 0;
}
// ************************************************************************* //

View file

@ -0,0 +1,46 @@
if (runTime.outputTime())
{
volVectorField gradT = fvc::grad(T);
volScalarField gradTx
(
IOobject
(
"gradTx",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
gradT.component(vector::X)
);
volScalarField gradTy
(
IOobject
(
"gradTy",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
gradT.component(vector::Y)
);
volScalarField gradTz
(
IOobject
(
"gradTz",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
gradT.component(vector::Z)
);
runTime.write();
}

View file

@ -0,0 +1,3 @@
pimpleDyMOversetFoam.C
EXE = $(FOAM_APPBIN)/pimpleDyMOversetFoam

View file

@ -0,0 +1,25 @@
EXE_INC = \
-I$(LIB_SRC)/dynamicMesh/dynamicFvMesh/lnInclude \
-I$(LIB_SRC)/dynamicMesh/dynamicMesh/lnInclude \
-I$(LIB_SRC)/meshTools/lnInclude \
-I$(LIB_SRC)/turbulenceModels/incompressible/turbulenceModel \
-I$(LIB_SRC)/transportModels \
-I$(LIB_SRC)/transportModels/incompressible/singlePhaseTransportModel \
-I$(LIB_SRC)/finiteVolume/lnInclude \
-I$(LIB_SRC)/oversetMesh/oversetMesh/lnInclude \
-I$(LIB_SRC)/oversetMesh/oversetDynamicFvMesh/lnInclude
EXE_LIBS = \
-ldynamicFvMesh \
-ltopoChangerFvMesh \
-ldynamicMesh \
-lmeshTools \
-lincompressibleTransportModels \
-lincompressibleTurbulenceModel \
-lincompressibleRASModels \
-lincompressibleLESModels \
-lfiniteVolume \
-llduSolvers \
-L$(MESQUITE_LIB_DIR) -lmesquite \
-loversetMesh \
-loversetDynamicFvMesh

View file

@ -0,0 +1,19 @@
fvVectorMatrix UEqn
(
fvm::ddt(U)
+ fvm::div(phi, U)
+ turbulence->divDevReff()
);
UEqn.relax
(
mesh.solutionDict().equationRelaxationFactor
(
U.select(pimple.finalIter())
)
);
if (pimple.momentumPredictor())
{
solve(UEqn == -fvc::grad(p));
}

View file

@ -0,0 +1,57 @@
{
wordList pcorrTypes(p.boundaryField().types());
for (label i = 0; i < p.boundaryField().size(); i++)
{
if (p.boundaryField()[i].fixesValue())
{
pcorrTypes[i] = fixedValueFvPatchScalarField::typeName;
}
}
volScalarField pcorr
(
IOobject
(
"pcorr",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::NO_WRITE
),
mesh,
dimensionedScalar("pcorr", p.dimensions(), 0),
pcorrTypes
);
// Initialise flux with interpolated velocity
phi = faceOversetMask*(fvc::interpolate(U) & mesh.Sf());
// Adjust overset fluxes
oversetAdjustPhi(phi, U); // Fringe flux adjustment
globalOversetAdjustPhi(phi, U, pcorr); // Global flux adjustment
while (pimple.correctNonOrthogonal())
{
fvScalarMatrix pcorrEqn
(
fvm::laplacian(rAU, pcorr) == fvc::div(phi)
);
// Adjust non-orthogonal fluxes if necessary
om.correctNonOrthoFluxes(pcorrEqn, U);
pcorrEqn.setReference(pRefCell, pRefValue);
pcorrEqn.solve();
if (pimple.finalNonOrthogonalIter())
{
phi -= pcorrEqn.flux();
}
// Fluxes are corrected to absolute velocity and further corrected
// later. HJ, 6/Feb/2009
}
# include "oversetContinuityErrs.H"
}

View file

@ -0,0 +1,11 @@
#include "createTimeControls.H"
bool correctPhi
(
pimple.dict().lookupOrDefault("correctPhi", false)
);
bool checkMeshCourantNo
(
pimple.dict().lookupOrDefault("checkMeshCourantNo", false)
);

View file

@ -0,0 +1,68 @@
Info<< "Reading field p\n" << endl;
volScalarField p
(
IOobject
(
"p",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);
Info<< "Reading field U\n" << endl;
volVectorField U
(
IOobject
(
"U",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);
# include "createPhi.H"
// Set reference cell taking into account whether implicit or explicit
// overset aglorithm is used for field p.
const oversetMesh& om = oversetMesh::New(mesh);
label pRefCell = 0;
scalar pRefValue = 0.0;
om.setRefCell(p, pimple.dict(), pRefCell, pRefValue);
mesh.schemesDict().setFluxRequired(p.name());
if (correctPhi)
{
mesh.schemesDict().setFluxRequired("pcorr");
}
singlePhaseTransportModel laminarTransport(U, phi);
autoPtr<incompressible::turbulenceModel> turbulence
(
incompressible::turbulenceModel::New(U, phi, laminarTransport)
);
Info<< "Reading field rAU if present\n" << endl;
volScalarField rAU
(
IOobject
(
"rAU",
runTime.timeName(),
mesh,
IOobject::READ_IF_PRESENT,
IOobject::AUTO_WRITE
),
mesh,
runTime.deltaT(),
zeroGradientFvPatchScalarField::typeName
);
# include "writeOversetMasks.H"

View file

@ -0,0 +1,76 @@
{
p.boundaryField().updateCoeffs();
rAU = 1.0/UEqn.A();
oversetFvPatchScalarField::oversetInterpolate(rAU); // Overset update
surfaceScalarField rAUf = fvc::interpolate(rAU);
U = rAU*UEqn.H();
oversetFvPatchVectorField::oversetInterpolate(U); // Overset update
surfaceScalarField phiU
(
"phiU",
fvc::interpolate(U) & mesh.Sf()
);
phi = phiU;
// Adjust overset fluxes
oversetAdjustPhi(phi, U); // Fringe flux adjustment
globalOversetAdjustPhi(phi, U, p); // Global flux adjustment
while (pimple.correctNonOrthogonal())
{
fvScalarMatrix pEqn
(
fvm::laplacian(rAUf, p) == fvc::div(phi)
);
// Adjust non-orthogonal fringe fluxes if necessary
om.correctNonOrthoFluxes(pEqn, U);
pEqn.setReference(pRefCell, pRefValue);
if
(
pimple.corrPISO() == pimple.nCorrPISO()
&& pimple.finalNonOrthogonalIter()
)
{
pEqn.solve(mesh.solutionDict().solver("pFinal"));
}
else
{
pEqn.solve();
}
if (pimple.finalNonOrthogonalIter())
{
phi -= pEqn.flux();
}
// Perform overset interpolation (after flux reconstruction)
oversetFvPatchScalarField::oversetInterpolate(p);
}
# include "oversetContinuityErrs.H"
// Explicitly relax pressure for momentum corrector
if (!pimple.finalIter())
{
p.relax();
}
# include "movingMeshContinuityErrs.H"
U += rAU*fvc::reconstruct((phi - phiU)/rAUf);
U.correctBoundaryConditions();
oversetFvPatchVectorField::oversetInterpolate(U); // Overset update
// Note: if implicit fringe is switched on, we are doing the interpolation
// twice (once in correctBoundaryConditions and once in oversetInterpolate)
// Reorganize. VV, 4/Oct/2016.
// Make the fluxes relative to the mesh motion
fvc::makeRelative(phi, U);
}

View file

@ -0,0 +1,129 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | foam-extend: Open Source CFD
\\ / O peration | Version: 4.0
\\ / A nd | Web: http://www.foam-extend.org
\\/ M anipulation | For copyright notice see file Copyright
-------------------------------------------------------------------------------
License
This file is part of foam-extend.
foam-extend 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.
foam-extend 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 foam-extend. If not, see <http://www.gnu.org/licenses/>.
Application
pimpleDyMOversetFoam.C
Description
Transient solver for incompressible, flow of Newtonian fluids
with with dynamic mesh and overset mesh using the PIMPLE
(merged PISO-SIMPLE) algorithm.
Turbulence modelling is generic, i.e. laminar, RAS or LES may be selected.
\*---------------------------------------------------------------------------*/
#include "fvCFD.H"
#include "singlePhaseTransportModel.H"
#include "turbulenceModel.H"
#include "dynamicFvMesh.H"
#include "pimpleControl.H"
#include "oversetMesh.H"
#include "oversetFvPatchFields.H"
#include "oversetAdjustPhi.H"
#include "globalOversetAdjustPhi.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
int main(int argc, char *argv[])
{
# include "setRootCase.H"
# include "createTime.H"
# include "createDynamicFvMesh.H"
pimpleControl pimple(mesh);
# include "initContinuityErrs.H"
# include "createControls.H"
# include "createFields.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Info<< "\nStarting time loop\n" << endl;
while (runTime.run())
{
# include "readControls.H"
// Make the fluxes absolute
fvc::makeAbsolute(phi, U);
runTime++;
Info<< "Time = " << runTime.timeName() << nl << endl;
bool meshChanged = mesh.update();
reduce(meshChanged, orOp<bool>());
# include "createOversetMasks.H"
# include "volContinuity.H"
if (checkMeshCourantNo)
{
# include "meshCourantNo.H"
}
// Mesh motion update
if (correctPhi && meshChanged)
{
// Fluxes will be corrected to absolute velocity
// HJ, 6/Feb/2009
# include "correctPhi.H"
}
// Make the fluxes relative to the mesh motion
fvc::makeRelative(phi, U);
# include "oversetCourantNo.H"
# include "setDeltaT.H"
// --- PIMPLE loop
while (pimple.loop())
{
# include "UEqn.H"
// --- PISO loop
while (pimple.correct())
{
# include "pEqn.H"
}
turbulence->correct();
}
runTime.write();
Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
<< " ClockTime = " << runTime.elapsedClockTime() << " s"
<< nl << endl;
}
Info<< "End\n" << endl;
return 0;
}
// ************************************************************************* //

View file

@ -0,0 +1,5 @@
#include "readTimeControls.H"
correctPhi = pimple.dict().lookupOrDefault("correctPhi", false);
checkMeshCourantNo = pimple.dict().lookupOrDefault("checkMeshCourantNo", false);

View file

@ -0,0 +1,3 @@
potentialDyMOversetFoam.C
EXE = $(FOAM_APPBIN)/potentialDyMOversetFoam

View file

@ -0,0 +1,18 @@
EXE_INC = \
-I$(LIB_SRC)/dynamicMesh/dynamicFvMesh/lnInclude \
-I$(LIB_SRC)/dynamicMesh/dynamicMesh/lnInclude \
-I$(LIB_SRC)/meshTools/lnInclude \
-I$(LIB_SRC)/finiteVolume/lnInclude \
-I$(LIB_SRC)/oversetMesh/oversetMesh/lnInclude \
-I$(LIB_SRC)/oversetMesh/oversetDynamicFvMesh/lnInclude
EXE_LIBS = \
-ldynamicFvMesh \
-ldynamicMesh \
-lengine \
-lmeshTools \
-lfiniteVolume \
-loversetMesh \
-loversetDynamicFvMesh \
-llduSolvers \
-L$(MESQUITE_LIB_DIR) -lmesquite

View file

@ -0,0 +1,66 @@
Info<< "Reading field p\n" << endl;
volScalarField p
(
IOobject
(
"p",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);
Info<< "Reading field U\n" << endl;
volVectorField U
(
IOobject
(
"U",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);
surfaceScalarField phi
(
IOobject
(
"phi",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
fvc::interpolate(U) & mesh.Sf()
);
// Write overset div phi for visualisation purposes
volScalarField oversetDivPhi
(
IOobject
(
"oversetDivPhi",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
fvc::div(phi)
);
// Set reference cell taking into account whether implicit or explicit
// overset algorithm is used for field p.
const oversetMesh& om = oversetMesh::New(mesh);
label pRefCell = 0;
scalar pRefValue = 0.0;
om.setRefCell(p, piso.dict(), pRefCell, pRefValue);
mesh.schemesDict().setFluxRequired(p.name());
// Write overset masks
# include "writeOversetMasks.H"

View file

@ -0,0 +1,190 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | foam-extend: Open Source CFD
\\ / O peration | Version: 4.0
\\ / A nd | Web: http://www.foam-extend.org
\\/ M anipulation | For copyright notice see file Copyright
-------------------------------------------------------------------------------
License
This file is part of foam-extend.
foam-extend 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.
foam-extend 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 foam-extend. If not, see <http://www.gnu.org/licenses/>.
Application
potentialDyMOversetFoam
Description
Transient solver for potential flow with dynamic overset mesh.
Author
Hrvoje Jasak, Wikki Ltd. All rights reserved.
Vuko Vukcevic, FMENA Zagreb. All rights reserved.
\*---------------------------------------------------------------------------*/
#include "fvCFD.H"
#include "dynamicFvMesh.H"
#include "pisoControl.H"
#include "oversetMesh.H"
#include "oversetFvPatchFields.H"
#include "oversetAdjustPhi.H"
#include "globalOversetAdjustPhi.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
int main(int argc, char *argv[])
{
argList::validOptions.insert("reconstructU", "");
# include "setRootCase.H"
# include "createTime.H"
# include "createDynamicFvMesh.H"
pisoControl piso(mesh);
# include "createFields.H"
# include "initTotalVolume.H"
# include "initContinuityErrs.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Info<< "\nStarting time loop\n" << endl;
while (runTime.loop())
{
# include "checkTotalVolume.H"
Info<< "Time = " << runTime.timeName() << nl << endl;
bool meshChanged = mesh.update();
reduce(meshChanged, orOp<bool>());
# include "createOversetMasks.H"
// Update moving wall velocity boundary condition and calculate the flux
U.correctBoundaryConditions();
// Forced overset update: make sure the overset interpolation is
// performed regardless of whether the coupledFringe is specified.
oversetFvPatchVectorField::oversetInterpolate(U);
phi == (linearInterpolate(U) & mesh.Sf());
// Resetting pressure field
p.internalField() = 0;
# include "volContinuity.H"
# include "meshCourantNo.H"
// Solve potential flow equations
// Adjust fluxes
oversetAdjustPhi(phi, U); // Fringe flux adjustment
globalOversetAdjustPhi(phi, U, p); // Global flux adjustment
while (piso.correctNonOrthogonal())
{
p.storePrevIter();
Info<< "Initial flux contour continuity error = "
<< mag(sum(phi.boundaryField()))
<< endl;
fvScalarMatrix pEqn
(
fvm::laplacian
(
dimensionedScalar
(
"1",
dimTime/p.dimensions()*dimensionSet(0, 2, -2, 0, 0),
1
),
p
)
==
fvc::div(phi)
);
// Adjust non-orthogonal fringe fluxes if necessary
om.correctNonOrthoFluxes(pEqn, U);
pEqn.setReference(pRefCell, pRefValue);
pEqn.solve();
if (piso.finalNonOrthogonalIter())
{
phi -= pEqn.flux();
# include "oversetContinuityErrs.H"
}
else
{
p.relax();
}
// Perform overset interpolation (after flux reconstruction)
oversetFvPatchScalarField::oversetInterpolate(p);
}
// Update div phi field for visualisation purposes
oversetDivPhi = cellOversetMask*fvc::div(phi);
if (args.optionFound("reconstructU"))
{
U = fvc::reconstruct(phi);
U.correctBoundaryConditions();
}
Info<< "Interpolated U error = "
<< (
sqrt
(
sum
(
sqr
(
faceOversetMask*
(
(fvc::interpolate(U) & mesh.Sf())
- phi
)
)
)
)/sum(mesh.magSf())
).value()
<< endl;
// Calculate velocity magnitude
{
volScalarField magU = mag(U);
Info<< "mag(U): max: " << gMax(magU.internalField())
<< " min: " << gMin(magU.internalField()) << endl;
}
runTime.write();
Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
<< " ClockTime = " << runTime.elapsedClockTime() << " s"
<< nl << endl;
}
Info<< "End\n" << endl;
return(0);
}
// ************************************************************************* //

View file

@ -0,0 +1,3 @@
potentialOversetFoam.C
EXE = $(FOAM_APPBIN)/potentialOversetFoam

View file

@ -0,0 +1,10 @@
EXE_INC = \
-I$(LIB_SRC)/finiteVolume/lnInclude \
-I$(LIB_SRC)/triSurface/lnInclude \
-I$(LIB_SRC)/meshTools/lnInclude \
-I$(LIB_SRC)/oversetMesh/oversetMesh/lnInclude
EXE_LIBS = \
-lfiniteVolume \
-loversetMesh \
-llduSolvers

View file

@ -0,0 +1,56 @@
Info<< "Reading field p\n" << endl;
volScalarField p
(
IOobject
(
"p",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::NO_WRITE
),
mesh
);
// Do not reset p
Info<< "Reading field U\n" << endl;
volVectorField U
(
IOobject
(
"U",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);
if (args.optionFound("resetU"))
{
U.internalField() = vector::zero;
}
surfaceScalarField phi
(
IOobject
(
"phi",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
fvc::interpolate(U) & mesh.Sf()
);
// Set reference cell taking into account whether implicit or explicit
// overset algorithm is used for field p.
label pRefCell = 0;
scalar pRefValue = 0.0;
om.setRefCell(p, simple.dict(), pRefCell, pRefValue);
mesh.schemesDict().setFluxRequired(p.name());

View file

@ -0,0 +1,262 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | foam-extend: Open Source CFD
\\ / O peration | Version: 4.0
\\ / A nd | Web: http://www.foam-extend.org
\\/ M anipulation | For copyright notice see file Copyright
-------------------------------------------------------------------------------
License
This file is part of foam-extend.
foam-extend 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.
foam-extend 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 foam-extend. If not, see <http://www.gnu.org/licenses/>.
Application
potentialOversetFoam
Description
Potential flow solver with overset mesh support.
Author
Hrvoje Jasak, Wikki Ltd. All rights reserved
\*---------------------------------------------------------------------------*/
#include "fvCFD.H"
#include "simpleControl.H"
#include "oversetMesh.H"
#include "oversetFvPatchFields.H"
#include "oversetAdjustPhi.H"
#include "globalOversetAdjustPhi.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
int main(int argc, char *argv[])
{
argList::validOptions.insert("resetU", "");
argList::validOptions.insert("writep", "");
# include "setRootCase.H"
# include "createTime.H"
# include "createMesh.H"
simpleControl simple(mesh);
# include "createOversetMasks.H"
# include "createFields.H"
# include "writeOversetMasks.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Info<< nl << "Calculating potential flow" << endl;
// Do correctors over the complete set
while (simple.correctNonOrthogonal())
{
phi = (linearInterpolate(U) & mesh.Sf());
// Adjust overset fluxes
oversetAdjustPhi(phi, U); // Fringe flux adjustment
globalOversetAdjustPhi(phi, U, p); // Global flux adjustment
Info<< "Initial flux contour continuity error = "
<< mag(sum(phi.boundaryField()))
<< endl;
p.storePrevIter();
fvScalarMatrix pEqn
(
fvm::laplacian
(
dimensionedScalar
(
"1",
dimTime/p.dimensions()*dimensionSet(0, 2, -2, 0, 0),
1
),
p
)
==
fvc::div(phi)
);
// Adjust non-orthogonal fringe fluxes if necessary
om.correctNonOrthoFluxes(pEqn, U);
pEqn.setReference(pRefCell, pRefValue);
pEqn.solve();
// Correct the flux
phi -= pEqn.flux();
volScalarField divFlux
(
"divFlux",
cellOversetMask*fvc::div(phi)
);
divFlux.write();
// Perform overset interpolation (after flux reconstruction)
oversetFvPatchScalarField::oversetInterpolate(p);
if (!simple.finalNonOrthogonalIter())
{
p.relax();
}
Info<< "p min " << gMin(p.internalField())
<< " max " << gMax(p.internalField())
<< " masked min "
<< gMin(cellOversetMask.internalField()*p.internalField())
<< " max "
<< gMax(cellOversetMask.internalField()*p.internalField())
<< endl;
Info<< "continuity error = "
<< mag
(
cellOversetMask*fvc::div(phi)
)().weightedAverage(mesh.V()).value()
<< endl;
Info<< "Contour continuity error = "
<< mag(sum(phi.boundaryField()))
<< endl;
U = fvc::reconstruct(phi);
U.correctBoundaryConditions();
oversetFvPatchVectorField::oversetInterpolate(U); // Overset update
// Note: if implicit fringe is switched on, we are doing the
// interpolation twice (once in correctBoundaryConditions and once
// in oversetInterpolate). Reorganize. VV, 4/Oct/2016.
Info<< "Interpolated U error = "
<< (
sqrt
(
sum
(
sqr
(
faceOversetMask*
(
(fvc::interpolate(U) & mesh.Sf())
- phi
)
)
)
)/sum(mesh.magSf())
).value()
<< endl;
}
// Calculate velocity magnitude
{
volScalarField magU = cellOversetMask*mag(U);
Info << "Masked mag(U): max: " << gMax(magU.internalField())
<< " min: " << gMin(magU.internalField()) << endl;
}
// Force the write
U.write();
p.write();
phi.write();
cellOversetMask.write();
if (args.optionFound("writep"))
{
// Find reference patch
label refPatch = -1;
scalar maxMagU = 0;
// Go through all velocity patches and find the one that fixes
// velocity to the largest value
forAll (U.boundaryField(), patchI)
{
const fvPatchVectorField& Upatch = U.boundaryField()[patchI];
if (Upatch.fixesValue())
{
// Calculate mean velocity
scalar u = sum(mag(Upatch));
label patchSize = Upatch.size();
reduce(u, sumOp<scalar>());
reduce(patchSize, sumOp<label>());
if (patchSize > 0)
{
scalar curMag = u/patchSize;
if (curMag > maxMagU)
{
refPatch = patchI;
maxMagU = curMag;
}
}
}
}
if (refPatch > -1)
{
// Calculate reference pressure
const fvPatchVectorField& Upatch = U.boundaryField()[refPatch];
const fvPatchScalarField& pPatch = p.boundaryField()[refPatch];
scalar patchE = sum(mag(pPatch + 0.5*magSqr(Upatch)));
label patchSize = Upatch.size();
reduce(patchE, sumOp<scalar>());
reduce(patchSize, sumOp<label>());
scalar e = patchE/patchSize;
Info<< "Using reference patch " << refPatch
<< " with mag(U) = " << maxMagU
<< " p + 0.5*U^2 = " << e << endl;
p.internalField() = e - 0.5*magSqr(U.internalField());
p.correctBoundaryConditions();
oversetFvPatchScalarField::oversetInterpolate(p); // Overset update
// Note: if implicit fringe is switched on, we are doing the
// interpolation twice (once in correctBoundaryConditions and once
// in oversetInterpolate). Reorganize. VV, 4/Oct/2016.
}
else
{
Info<< "No reference patch found. Writing potential function"
<< endl;
}
p.write();
}
Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
<< " ClockTime = " << runTime.elapsedClockTime() << " s"
<< nl << endl;
Info<< "End\n" << endl;
return 0;
}
// ************************************************************************* //

View file

@ -0,0 +1,3 @@
scalarTransportOversetFoam.C
EXE = $(FOAM_APPBIN)/scalarTransportOversetFoam

View file

@ -0,0 +1,8 @@
EXE_INC = \
-I$(LIB_SRC)/finiteVolume/lnInclude \
-I$(LIB_SRC)/oversetMesh/oversetMesh/lnInclude
EXE_LIBS = \
-lfiniteVolume \
-loversetMesh \
-llduSolvers

View file

@ -0,0 +1,75 @@
Info<< "Reading field T\n" << endl;
volScalarField T
(
IOobject
(
"T",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);
Info<< "Reading field V\n" << endl;
volVectorField V
(
IOobject
(
"V",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);
Info<< "Reading field U\n" << endl;
volVectorField U
(
IOobject
(
"U",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);
Info<< "Reading transportProperties\n" << endl;
IOdictionary transportProperties
(
IOobject
(
"transportProperties",
runTime.constant(),
mesh,
IOobject::MUST_READ,
IOobject::NO_WRITE
)
);
Info<< "Reading diffusivity DT\n" << endl;
dimensionedScalar DT
(
transportProperties.lookup("DT")
);
dimensionedScalar DV
(
transportProperties.lookup("DV")
);
# include "createPhi.H"

View file

@ -0,0 +1,120 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | foam-extend: Open Source CFD
\\ / O peration | Version: 4.0
\\ / A nd | Web: http://www.foam-extend.org
\\/ M anipulation | For copyright notice see file Copyright
-------------------------------------------------------------------------------
License
This file is part of foam-extend.
foam-extend 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.
foam-extend 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 foam-extend. If not, see <http://www.gnu.org/licenses/>.
Application
scalarTransportOversetFoam
Description
Solves a transport equation for a passive scalar with support for
overset meshes.
Experimental: no overset mesh changes required
Author
Hrvoje Jasak, Wikki Ltd. All rights reserved
\*---------------------------------------------------------------------------*/
#include "fvCFD.H"
#include "simpleControl.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
int main(int argc, char *argv[])
{
# include "setRootCase.H"
# include "createTime.H"
# include "createMesh.H"
simpleControl simple(mesh);
# include "createFields.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Info<< "\nCalculating scalar transport\n" << endl;
# include "CourantNo.H"
while (simple.loop())
{
Info<< "Time = " << runTime.timeName() << nl << endl;
fvScalarMatrix TEqn
(
fvm::ddt(T)
+ fvm::div(phi, T)
- fvm::laplacian(DT, T)
);
TEqn.solve();
// volScalarField Tresidual
// (
// "Tresidual",
// T
// );
// Tresidual.writeOpt() = IOobject::AUTO_WRITE;
// Tresidual.internalField() = TEqn.residual();
// Tresidual.boundaryField() == 0;
// Info<< "residual " << gSumMag(Tresidual.internalField()) << endl;
fvVectorMatrix VEqn
(
fvm::ddt(V)
+ fvm::div(phi, V)
- fvm::laplacian(DV, V)
);
VEqn.solve();
// return 0;
// VEqn.solve();
// volVectorField Vresidual
// (
// "Vresidual",
// V
// );
// Vresidual.writeOpt() = IOobject::AUTO_WRITE;
// Vresidual.internalField() = VEqn.residual();
// Vresidual.boundaryField() == vector::zero;
// Info<< "residual " << gSumMag(Vresidual.internalField()) << endl;
runTime.write();
}
Info<< "End\n" << endl;
return 0;
}
// ************************************************************************* //

View file

@ -0,0 +1,3 @@
simpleMRFOversetFoam.C
EXE = $(FOAM_APPBIN)/simpleMRFOversetFoam

View file

@ -0,0 +1,16 @@
EXE_INC = \
-I$(LIB_SRC)/turbulenceModels \
-I$(LIB_SRC)/turbulenceModels/incompressible/RAS/RASModel \
-I$(LIB_SRC)/transportModels \
-I$(LIB_SRC)/transportModels/incompressible/singlePhaseTransportModel \
-I$(LIB_SRC)/finiteVolume/lnInclude \
-I$(LIB_SRC)/triSurface/lnInclude \
-I$(LIB_SRC)/meshTools/lnInclude \
-I$(LIB_SRC)/oversetMesh/oversetMesh/lnInclude
EXE_LIBS = \
-lincompressibleRASModels \
-lincompressibleTransportModels \
-lfiniteVolume \
-loversetMesh \
-llduSolvers

View file

@ -0,0 +1,14 @@
// Solve the Momentum equation
tmp<fvVectorMatrix> UEqn
(
fvm::div(phi, U)
+ turbulence->divDevReff()
);
mrfZones.addCoriolis(UEqn());
UEqn().relax();
solve
(
UEqn() == -fvc::grad(p)
);

View file

@ -0,0 +1,83 @@
Info << "Reading field p\n" << endl;
volScalarField p
(
IOobject
(
"p",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);
Info<< "Reading field U\n" << endl;
volVectorField U
(
IOobject
(
"U",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);
# include "createPhi.H"
// Set reference cell taking into account whether implicit or explicit
// overset aglorithm is used for field p.
label pRefCell = 0;
scalar pRefValue = 0.0;
om.setRefCell(p, simple.dict(), pRefCell, pRefValue);
mesh.schemesDict().setFluxRequired(p.name());
singlePhaseTransportModel laminarTransport(U, phi);
autoPtr<incompressible::RASModel> turbulence
(
incompressible::RASModel::New(U, phi, laminarTransport)
);
// Create MRF zones
MRFZones mrfZones(mesh);
mrfZones.correctBoundaryVelocity(U);
// Create Urel as a permanent field to make it available for on-the-fly
// post-processing operations
volVectorField Urel
(
IOobject
(
"Urel",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
U
);
mrfZones.relativeVelocity(Urel);
Info<< "Reading field rAU if present\n" << endl;
volScalarField rAU
(
IOobject
(
"rAU",
runTime.timeName(),
mesh,
IOobject::READ_IF_PRESENT,
IOobject::AUTO_WRITE
),
mesh,
dimensionedScalar("dt", dimTime, 1.0),
zeroGradientFvPatchScalarField::typeName
);
# include "writeOversetMasks.H"

View file

@ -0,0 +1,59 @@
{
p.boundaryField().updateCoeffs();
rAU = 1.0/UEqn().A();
oversetFvPatchScalarField::oversetInterpolate(rAU); // Overset update
// Update boundary velocity for consistency with the flux
// This is needed for a ramped MRF
mrfZones.correctBoundaryVelocity(U);
U = rAU*UEqn().H();
oversetFvPatchVectorField::oversetInterpolate(U); // Overset update
UEqn.clear();
phi = (fvc::interpolate(U) & mesh.Sf());
mrfZones.relativeFlux(phi);
// Adjust overset fluxes
oversetAdjustPhi(phi, U); // Fringe flux adjustment
globalOversetAdjustPhi(phi, U, p); // Global flux adjustment
// Non-orthogonal pressure corrector loop
while (simple.correctNonOrthogonal())
{
fvScalarMatrix pEqn
(
fvm::laplacian(rAU, p) == fvc::div(phi)
);
// Adjust non-orthogonal fringe fluxes if necessary
om.correctNonOrthoFluxes(pEqn, U);
pEqn.setReference(pRefCell, pRefValue);
// Retain the residual from the first iteration
pEqn.solve();
if (simple.finalNonOrthogonalIter())
{
phi -= pEqn.flux();
}
// Perform overset interpolation (after flux reconstruction)
oversetFvPatchScalarField::oversetInterpolate(p);
}
# include "oversetContinuityErrs.H"
// Explicitly relax pressure for momentum corrector
p.relax();
// Momentum corrector
U -= rAU*fvc::grad(p);
U.correctBoundaryConditions();
oversetFvPatchVectorField::oversetInterpolate(U); // Overset update
// Note: if implicit fringe is switched on, we are doing the interpolation
// twice (once in correctBoundaryConditions and once in oversetInterpolate)
// Reorganize. VV, 4/Oct/2016.
}

View file

@ -0,0 +1,94 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | foam-extend: Open Source CFD
\\ / O peration | Version: 4.0
\\ / A nd | Web: http://www.foam-extend.org
\\/ M anipulation | For copyright notice see file Copyright
-------------------------------------------------------------------------------
License
This file is part of foam-extend.
foam-extend 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.
foam-extend 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 foam-extend. If not, see <http://www.gnu.org/licenses/>.
Application
simpleMRFOversetFoam
Description
Steady-state solver for incompressible, turbulent flow
with overset mesh support and MRF regions.
Author
Vuko Vukcevic, Wikki Ltd. All rights reserved.
\*---------------------------------------------------------------------------*/
#include "fvCFD.H"
#include "singlePhaseTransportModel.H"
#include "RASModel.H"
#include "simpleControl.H"
#include "MRFZones.H"
#include "oversetMesh.H"
#include "oversetFvPatchFields.H"
#include "oversetAdjustPhi.H"
#include "globalOversetAdjustPhi.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
int main(int argc, char *argv[])
{
# include "setRootCase.H"
# include "createTime.H"
# include "createMesh.H"
simpleControl simple(mesh);
# include "createOversetMasks.H"
# include "createFields.H"
# include "initContinuityErrs.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Info<< "\nStarting time loop\n" << endl;
while (simple.loop())
{
Info<< "Time = " << runTime.timeName() << nl << endl;
// Pressure-velocity SIMPLE corrector
{
# include "UEqn.H"
# include "pEqn.H"
}
// Calculate relative velocity
Urel == U;
mrfZones.relativeVelocity(Urel);
turbulence->correct();
runTime.write();
Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
<< " ClockTime = " << runTime.elapsedClockTime() << " s"
<< nl << endl;
}
Info<< "End\n" << endl;
return 0;
}
// ************************************************************************* //

View file

@ -0,0 +1,3 @@
simpleOversetFoam.C
EXE = $(FOAM_APPBIN)/simpleOversetFoam

View file

@ -0,0 +1,16 @@
EXE_INC = \
-I$(LIB_SRC)/turbulenceModels \
-I$(LIB_SRC)/turbulenceModels/incompressible/RAS/RASModel \
-I$(LIB_SRC)/transportModels \
-I$(LIB_SRC)/transportModels/incompressible/singlePhaseTransportModel \
-I$(LIB_SRC)/finiteVolume/lnInclude \
-I$(LIB_SRC)/triSurface/lnInclude \
-I$(LIB_SRC)/meshTools/lnInclude \
-I$(LIB_SRC)/oversetMesh/oversetMesh/lnInclude
EXE_LIBS = \
-lincompressibleRASModels \
-lincompressibleTransportModels \
-lfiniteVolume \
-loversetMesh \
-llduSolvers

View file

@ -0,0 +1,13 @@
// Solve the Momentum equation
tmp<fvVectorMatrix> UEqn
(
fvm::div(phi, U)
+ turbulence->divDevReff()
);
UEqn().relax();
solve
(
UEqn() == -fvc::grad(p)
);

View file

@ -0,0 +1,63 @@
Info << "Reading field p\n" << endl;
volScalarField p
(
IOobject
(
"p",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);
Info<< "Reading field U\n" << endl;
volVectorField U
(
IOobject
(
"U",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);
# include "createPhi.H"
// Set reference cell taking into account whether implicit or explicit
// overset aglorithm is used for field p.
label pRefCell = 0;
scalar pRefValue = 0.0;
om.setRefCell(p, simple.dict(), pRefCell, pRefValue);
mesh.schemesDict().setFluxRequired(p.name());
singlePhaseTransportModel laminarTransport(U, phi);
autoPtr<incompressible::RASModel> turbulence
(
incompressible::RASModel::New(U, phi, laminarTransport)
);
Info<< "Reading field rAU if present\n" << endl;
volScalarField rAU
(
IOobject
(
"rAU",
runTime.timeName(),
mesh,
IOobject::READ_IF_PRESENT,
IOobject::AUTO_WRITE
),
mesh,
dimensionedScalar("dt", dimTime, 1.0),
zeroGradientFvPatchScalarField::typeName
);
# include "writeOversetMasks.H"

View file

@ -0,0 +1,54 @@
{
p.boundaryField().updateCoeffs();
rAU = 1.0/UEqn().A();
oversetFvPatchScalarField::oversetInterpolate(rAU); // Overset update
U = rAU*UEqn().H();
oversetFvPatchVectorField::oversetInterpolate(U); // Overset update
UEqn.clear();
phi = (fvc::interpolate(U) & mesh.Sf());
// Adjust overset fluxes
oversetAdjustPhi(phi, U); // Fringe flux adjustment
globalOversetAdjustPhi(phi, U, p); // Global flux adjustment
// Non-orthogonal pressure corrector loop
while (simple.correctNonOrthogonal())
{
fvScalarMatrix pEqn
(
fvm::laplacian(rAU, p) == fvc::div(phi)
);
// Adjust non-orthogonal fringe fluxes if necessary
om.correctNonOrthoFluxes(pEqn, U);
pEqn.setReference(pRefCell, pRefValue);
// Retain the residual from the first iteration
pEqn.solve();
if (simple.finalNonOrthogonalIter())
{
phi -= pEqn.flux();
}
// Perform overset interpolation (after flux reconstruction)
oversetFvPatchScalarField::oversetInterpolate(p);
}
# include "oversetContinuityErrs.H"
// Explicitly relax pressure for momentum corrector
p.relax();
// Momentum corrector
U -= rAU*fvc::grad(p);
U.correctBoundaryConditions();
oversetFvPatchVectorField::oversetInterpolate(U); // Overset update
// Note: if implicit fringe is switched on, we are doing the interpolation
// twice (once in correctBoundaryConditions and once in oversetInterpolate)
// Reorganize. VV, 4/Oct/2016.
}

View file

@ -0,0 +1,89 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | foam-extend: Open Source CFD
\\ / O peration | Version: 4.0
\\ / A nd | Web: http://www.foam-extend.org
\\/ M anipulation | For copyright notice see file Copyright
-------------------------------------------------------------------------------
License
This file is part of foam-extend.
foam-extend 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.
foam-extend 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 foam-extend. If not, see <http://www.gnu.org/licenses/>.
Application
simpleOversetFoam
Description
Steady-state solver for incompressible, turbulent flow
with overset mesh support.
Author
Hrvoje Jasak, Wikki Ltd. All rights reserved
\*---------------------------------------------------------------------------*/
#include "fvCFD.H"
#include "singlePhaseTransportModel.H"
#include "RASModel.H"
#include "simpleControl.H"
#include "oversetMesh.H"
#include "oversetFvPatchFields.H"
#include "oversetAdjustPhi.H"
#include "globalOversetAdjustPhi.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
int main(int argc, char *argv[])
{
# include "setRootCase.H"
# include "createTime.H"
# include "createMesh.H"
simpleControl simple(mesh);
# include "createOversetMasks.H"
# include "createFields.H"
# include "initContinuityErrs.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Info<< "\nStarting time loop\n" << endl;
while (simple.loop())
{
Info<< "Time = " << runTime.timeName() << nl << endl;
// Pressure-velocity SIMPLE corrector
{
# include "UEqn.H"
# include "pEqn.H"
}
turbulence->correct();
runTime.write();
Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
<< " ClockTime = " << runTime.elapsedClockTime() << " s"
<< nl << endl;
}
Info<< "End\n" << endl;
return 0;
}
// ************************************************************************* //

View file

@ -0,0 +1,3 @@
calcOverset.C
EXE = $(FOAM_APPBIN)/calcOverset

View file

@ -0,0 +1,16 @@
EXE_INC = \
-I$(LIB_SRC)/turbulenceModels \
-I$(LIB_SRC)/turbulenceModels/incompressible/RAS/RASModel \
-I$(LIB_SRC)/transportModels \
-I$(LIB_SRC)/transportModels/incompressible/singlePhaseTransportModel \
-I$(LIB_SRC)/finiteVolume/lnInclude \
-I$(LIB_SRC)/triSurface/lnInclude \
-I$(LIB_SRC)/meshTools/lnInclude \
-I$(LIB_SRC)/oversetMesh/oversetMesh/lnInclude
EXE_LIBS = \
-lincompressibleRASModels \
-lincompressibleTransportModels \
-lfiniteVolume \
-loversetMesh \
-llduSolvers

View file

@ -0,0 +1,87 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright held by original author
\\/ 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 2 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, write to the Free Software Foundation,
Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
Application
calcOverset
Description
Calculates overset addressing and everything necessary for overset
interpolation. Used for parallel scaling tests of overset assembly.
Author
Vuko Vukcevic, FMENA Zagreb. All rights reserved.
\*---------------------------------------------------------------------------*/
#include "fvCFD.H"
#include "oversetMesh.H"
#include "oversetFvPatchFields.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
int main(int argc, char *argv[])
{
# include "setRootCase.H"
# include "createTime.H"
# include "createMesh.H"
volVectorField cellCentres
(
IOobject
(
"cellCentres",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
mesh.C(),
"zeroGradient"
);
volVectorField origCellCentres("origCellCentres", cellCentres);
origCellCentres.write();
// Create CPU time object
cpuTime oversetTime;
// Make sure everything necessary for overset interpolation has been
// calculated by calling a single overset interpolation
oversetFvPatchVectorField::oversetInterpolate(cellCentres);
Info<< "Elapsed time for overset assembly and single interpolation: "
<< oversetTime.cpuTimeIncrement() << endl;
// Write interpolated cell centres
cellCentres.write();
# include "writeOversetMasks.H"
Info<< "End\n" << endl;
return 0;
}
// ************************************************************************* //

View file

@ -0,0 +1,3 @@
oversetRegionIndex.C
EXE = $(FOAM_APPBIN)/oversetRegionIndex

View file

@ -0,0 +1,7 @@
EXE_INC = \
-I$(LIB_SRC)/finiteVolume/lnInclude \
-I$(LIB_SRC)/meshTools/lnInclude
EXE_LIBS = \
-lfiniteVolume \
-lmeshTools

View file

@ -0,0 +1,82 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright held by original author
\\/ 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 2 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, write to the Free Software Foundation,
Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
Application
oversetRegionID
Description
Visualise region split in overset mesh
Author
Hrvoje Jasak, Wikki Ltd. All rights reserved.
\*---------------------------------------------------------------------------*/
#include "fvCFD.H"
#include "regionSplit.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
int main(int argc, char *argv[])
{
# include "setRootCase.H"
# include "createTime.H"
# include "createMesh.H"
volScalarField regionIndex
(
IOobject
(
"regionIndex",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::NO_WRITE
),
mesh,
dimensionedScalar("zero", dimless, 0)
);
scalarField& rIn = regionIndex.internalField();
// Create region split. Region index depends only on mesh ordering
regionSplit rs(mesh);
forAll (rIn, cellI)
{
rIn[cellI] = rs[cellI];
}
Info<< "Number of regions: " << rs.nRegions() << endl;
Info << "Write overset region split ... ";
regionIndex.write();
Info << "done" << endl;
return 0;
}
// ************************************************************************* //

View file

@ -81,6 +81,9 @@ wmake libso immersedBoundary/immersedBoundaryTurbulence
wmake libso immersedBoundary/immersedBoundaryForce
wmake libso immersedBoundary/immersedBoundaryDynamicMesh
wmake libso overset/oversetMesh
wmake libso overset/oversetDynamicFvMesh
( cd cudaSolvers ; ./Allwmake )
# ----------------------------------------------------------------- end-of-file

View file

@ -0,0 +1,4 @@
movingOversetRegion/movingOversetRegion.C
oversetSolidBodyMotionFvMesh/oversetSolidBodyMotionFvMesh.C
LIB = $(FOAM_LIBBIN)/liboversetDynamicFvMesh

View file

@ -0,0 +1,19 @@
EXE_INC = \
-I$(LIB_SRC)/finiteVolume/lnInclude \
-I$(LIB_SRC)/meshTools/lnInclude \
-I$(LIB_SRC)/surfMesh/lnInclude \
-I$(LIB_SRC)/sampling/lnInclude \
-I$(LIB_SRC)/dynamicMesh/dynamicMesh/lnInclude \
-I$(LIB_SRC)/dynamicMesh/dynamicFvMesh/lnInclude \
-I$(LIB_SRC)/dynamicMesh/meshMotion/solidBodyMotion/lnInclude \
-I$(LIB_SRC)/oversetMesh/oversetMesh/lnInclude
EXE_LIBS = \
-lfiniteVolume \
-lmeshTools \
-lsurfMesh \
-lsampling \
-ldynamicMesh \
-ldynamicFvMesh \
-lsolidBodyMotion \
-loversetMesh

View file

@ -0,0 +1,153 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | foam-extend: Open Source CFD
\\ / O peration |
\\ / A nd | For copyright notice see file Copyright
\\/ M anipulation |
-------------------------------------------------------------------------------
License
This file is part of foam-extend.
foam-extend 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.
foam-extend 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 foam-extend. If not, see <http://www.gnu.org/licenses/>.
\*---------------------------------------------------------------------------*/
#include "movingOversetRegion.H"
#include "transformField.H"
// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
void Foam::movingOversetRegion::calcMotionMask() const
{
if (motionMaskPtr_)
{
FatalErrorIn("void movingOversetRegion::calcMotionMask() const")
<< "Motion mask for movingOversetRegion " << name()
<< " already calculated"
<< abort(FatalError);
}
motionMaskPtr_ = new scalarField(mesh().allPoints().size(), 0);
scalarField& mask = *motionMaskPtr_;
const labelListList& cp = mesh().cellPoints();
forAll (movingZoneNames_, mzI)
{
// Find moving cell zone
const label zoneID =
mesh().cellZones().findZoneID(movingZoneNames_[mzI]);
if (zoneID < 0)
{
// Zone not found. Consider if this is valid in parallel?
FatalErrorIn("void movingOversetRegion::calcMotionMask() const")
<< "Cannot find moving zone " << movingZoneNames_[mzI]
<< " for movingOversetRegion " << name()
<< abort(FatalError);
}
// Get zone cells and mark vertices
const cellZone& mz = mesh().cellZones()[zoneID];
forAll (mz, cellI)
{
const labelList& curCp = cp[mz[cellI]];
forAll (curCp, pointI)
{
mask[curCp[pointI]] = 1;
}
}
}
// if (debug)
{
// Count moving points
label nMovingPoints = 0;
forAll (mask, i)
{
if (mask[i] > 0)
{
nMovingPoints++;
}
}
Info<< "movingOversetRegion " << name() << ": "
<< returnReduce(nMovingPoints, sumOp<scalar>()) << " moving points"
<< endl;
}
}
void Foam::movingOversetRegion::clearOut()
{
deleteDemandDrivenData(motionMaskPtr_);
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
Foam::movingOversetRegion::movingOversetRegion
(
const word& name,
const fvMesh& mesh,
const dictionary& dict
)
:
name_(name),
mesh_(mesh),
sbmfPtr_(solidBodyMotionFunction::New(dict, mesh_.time())),
movingZoneNames_(dict.lookup("movingZones")),
motionMaskPtr_(NULL)
{}
// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
Foam::movingOversetRegion::~movingOversetRegion()
{
clearOut();
}
// * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
const Foam::scalarField& Foam::movingOversetRegion::motionMask() const
{
if (!motionMaskPtr_)
{
calcMotionMask();
}
return *motionMaskPtr_;
}
const Foam::tmp<Foam::pointField>
Foam::movingOversetRegion::motionIncrement
(
const pointField& undisplacedPoints
) const
{
return motionMask()*
(
transform(sbmfPtr_->transformation(), undisplacedPoints)
- undisplacedPoints
);
}
// ************************************************************************* //

View file

@ -0,0 +1,137 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | foam-extend: Open Source CFD
\\ / O peration |
\\ / A nd | For copyright notice see file Copyright
\\/ M anipulation |
-------------------------------------------------------------------------------
License
This file is part of foam-extend.
foam-extend 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.
foam-extend 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 foam-extend. If not, see <http://www.gnu.org/licenses/>.
Class
Foam::movingOversetRegion
Description
Moving region in overset motion hierarchy. Motion is prescribed using a
solid body motion function and moving cells are chosen as a set of
cell zones.
SourceFiles
movingOversetRegion.C
\*---------------------------------------------------------------------------*/
#ifndef movingOversetRegion_H
#define movingOversetRegion_H
#include "fvMesh.H"
#include "solidBodyMotionFunction.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
/*---------------------------------------------------------------------------*\
Class movingOversetRegion Declaration
\*---------------------------------------------------------------------------*/
class movingOversetRegion
{
// Private data
//- Name
const word name_;
//- Mesh reference
const fvMesh& mesh_;
//- Overset region motion control function
autoPtr<solidBodyMotionFunction> sbmfPtr_;
//- Moving cell zone names
wordList movingZoneNames_;
//- Region motion mask
mutable scalarField* motionMaskPtr_;
// Private Member Functions
//- Disallow default bitwise copy construct
movingOversetRegion(const movingOversetRegion&);
//- Disallow default bitwise assignment
void operator=(const movingOversetRegion&);
//- Calculate motion mask
void calcMotionMask() const;
// Clear storage
void clearOut();
public:
// Constructors
//- Construct from dictionary
movingOversetRegion
(
const word& name,
const fvMesh& mesh,
const dictionary& dict
);
//- Destructor
~movingOversetRegion();
// Member Functions
//- Return name
const word& name() const
{
return name_;
}
//- Return mesh
const fvMesh& mesh() const
{
return mesh_;
}
//- Return motion mask
const scalarField& motionMask() const;
//- Return motion increment
const tmp<pointField> motionIncrement
(
const pointField& undisplacedPoints
) const;
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //

View file

@ -0,0 +1,136 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright held by original author
\\/ 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 2 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, write to the Free Software Foundation,
Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
\*---------------------------------------------------------------------------*/
#include "oversetSolidBodyMotionFvMesh.H"
#include "addToRunTimeSelectionTable.H"
#include "volFields.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
namespace Foam
{
defineTypeNameAndDebug(oversetSolidBodyMotionFvMesh, 0);
addToRunTimeSelectionTable
(
dynamicFvMesh,
oversetSolidBodyMotionFvMesh,
IOobject
);
}
// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
Foam::oversetSolidBodyMotionFvMesh::oversetSolidBodyMotionFvMesh
(
const IOobject& io
)
:
dynamicFvMesh(io),
dynamicMeshCoeffs_
(
IOdictionary
(
IOobject
(
"dynamicMeshDict",
time().constant(),
*this,
IOobject::MUST_READ,
IOobject::NO_WRITE
)
).subDict(typeName + "Coeffs")
),
motionRegions_(),
undisplacedPoints_
(
IOobject
(
"points",
time().constant(),
meshSubDir,
*this,
IOobject::MUST_READ,
IOobject::NO_WRITE
)
),
curTimeIndex_(-1)
{
// Read motion function for all regions
PtrList<entry> motionDicts(dynamicMeshCoeffs_.lookup("motionFunctions"));
motionRegions_.setSize(motionDicts.size());
forAll (motionDicts, mI)
{
motionRegions_.set
(
mI,
new movingOversetRegion
(
motionDicts[mI].keyword(),
*this,
motionDicts[mI].dict()
)
);
}
}
// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
Foam::oversetSolidBodyMotionFvMesh::~oversetSolidBodyMotionFvMesh()
{}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
bool Foam::oversetSolidBodyMotionFvMesh::update()
{
// Move the mesh only once per time step as this is prescribed motion
if (curTimeIndex_ < this->time().timeIndex())
{
pointField newPoints = undisplacedPoints_;
forAll (motionRegions_, regionI)
{
newPoints +=
motionRegions_[regionI].motionIncrement(undisplacedPoints_);
}
fvMesh::movePoints(newPoints);
curTimeIndex_ = this->time().timeIndex();
}
return false;
}
// ************************************************************************* //

View file

@ -0,0 +1,115 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright held by original author
\\/ 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 2 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, write to the Free Software Foundation,
Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
Class
Foam::oversetSolidBodyMotionFvMesh
Description
Solid-body motion of the overset mesh specified by a run-time selectable
motion function.
SourceFiles
oversetSolidBodyMotionFvMesh.C
\*---------------------------------------------------------------------------*/
#ifndef oversetSolidBodyMotionFvMesh_H
#define oversetSolidBodyMotionFvMesh_H
#include "dynamicFvMesh.H"
#include "dictionary.H"
#include "movingOversetRegion.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
/*---------------------------------------------------------------------------*\
Class oversetSolidBodyMotionFvMesh Declaration
\*---------------------------------------------------------------------------*/
class oversetSolidBodyMotionFvMesh
:
public dynamicFvMesh
{
// Private data
//- Dictionary of motion control parameters
dictionary dynamicMeshCoeffs_;
//- Overset region motion control function
PtrList<movingOversetRegion> motionRegions_;
//- Reference points which are transformed
pointIOField undisplacedPoints_;
// Helper variables
//- Current time index
label curTimeIndex_;
// Private Member Functions
//- Disallow default bitwise copy construct
oversetSolidBodyMotionFvMesh(const oversetSolidBodyMotionFvMesh&);
//- Disallow default bitwise assignment
void operator=(const oversetSolidBodyMotionFvMesh&);
public:
//- Runtime type information
TypeName("oversetSolidBodyMotionFvMesh");
// Constructors
//- Construct from IOobject
explicit oversetSolidBodyMotionFvMesh(const IOobject& io);
// Destructor
virtual ~oversetSolidBodyMotionFvMesh();
// Member Functions
//- Update the mesh for both mesh motion and topology change
virtual bool update();
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //

View file

@ -0,0 +1,43 @@
oversetRegion/oversetRegion.C
oversetFringe/oversetFringe/oversetFringe.C
oversetFringe/oversetFringe/newOversetFringe.C
oversetFringe/manualFringe/manualFringe.C
oversetFringe/faceCellsFringe/faceCellsFringe.C
oversetFringe/overlapFringe/overlapFringe.C
oversetFringe/compositeFringe/compositeFringe.C
oversetFringe/overlapFringe/donorSuitability/donorSuitability/donorSuitability.C
oversetFringe/overlapFringe/donorSuitability/donorSuitability/newDonorSuitability.C
oversetFringe/overlapFringe/donorSuitability/noSuitability/noSuitability.C
oversetFringe/overlapFringe/donorSuitability/patchDistance/patchDistance.C
oversetFringe/overlapFringe/donorSuitability/cellVolumes/cellVolumes.C
oversetFringe/overlapFringe/donorSuitability/faceArea/faceArea.C
oversetFringe/overlapFringe/donorSuitability/cellBoundingBoxDiagonal/cellBoundingBoxDiagonal.C
oversetInterpolation/oversetInterpolation/oversetInterpolation.C
oversetInterpolation/oversetInterpolation/newOversetInterpolation.C
oversetInterpolation/injectionInterpolation/injectionInterpolation.C
oversetInterpolation/inverseDistanceInterpolation/inverseDistanceInterpolation.C
oversetInterpolation/averageValueInterpolation/averageValueInterpolation.C
oversetMesh/oversetMesh.C
oversetMesh/oversetMeshAddressing.C
oversetLduInterface/oversetLduInterface.C
oversetLduInterfaceField/oversetLduInterfaceField.C
oversetPolyPatch/oversetPolyPatch.C
oversetPointPatch/oversetPointPatch.C
oversetFvPatch/oversetFvPatch.C
oversetFvPatchField/oversetFvPatchFields.C
oversetFvPatchField/emptyOversetFvPatchField/emptyOversetFvPatchFields.C
oversetFvsPatchField/oversetFvsPatchFields.C
oversetAdjustPhi/oversetAdjustPhi.C
oversetAdjustPhi/globalOversetAdjustPhi.C
oversetAdjustPhi/regionWiseOversetAdjustPhi.C
LIB = $(FOAM_LIBBIN)/liboversetMesh

View file

@ -0,0 +1,13 @@
EXE_INC = \
-I$(LIB_SRC)/finiteVolume/lnInclude \
-I$(LIB_SRC)/meshTools/lnInclude \
-I$(LIB_SRC)/surfMesh/lnInclude \
-I$(LIB_SRC)/sampling/lnInclude \
-I$(LIB_SRC)/dynamicMesh/dynamicMesh/lnInclude
EXE_LIBS = \
-lfiniteVolume \
-lmeshTools \
-lsurfMesh \
-lsampling \
-ldynamicMesh

View file

@ -0,0 +1,423 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | foam-extend: Open Source CFD
\\ / O peration |
\\ / A nd | For copyright notice see file Copyright
\\/ M anipulation |
-------------------------------------------------------------------------------
License
This file is part of foam-extend.
foam-extend 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.
foam-extend 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 foam-extend. If not, see <http://www.gnu.org/licenses/>.
Class
Foam::donorAcceptor
Description
Class holds donor and acceptor data. Used for searching and communications
SourceFiles
Author
Hrvoje Jasak, Wikki Ltd. All rights reserved.
Contributor
Vuko Vukcevic, FMENA Zagreb. All rights reserved.
\*---------------------------------------------------------------------------*/
#ifndef donorAcceptor_H
#define donorAcceptor_H
#include "label.H"
#include "point.H"
#include "DynamicList.H"
#include "boolList.H"
#include "vectorField.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
// Forward declaration of classes
class Istream;
class Ostream;
// Forward declaration of friend functions and operators
class donorAcceptor;
inline bool operator==(const donorAcceptor&, const donorAcceptor&);
inline bool operator!=(const donorAcceptor&, const donorAcceptor&);
Istream& operator>>(Istream&, donorAcceptor&);
Ostream& operator<<(Ostream&, const donorAcceptor&);
/*---------------------------------------------------------------------------*\
Class donorAcceptor Declaration
\*---------------------------------------------------------------------------*/
class donorAcceptor
{
public:
// Public data types
typedef DynamicList<label, 10> DynamicLabelList;
typedef DynamicList<point, 10> DynamicPointList;
private:
// Private data
// Acceptor side
//- Acceptor cell number
label acceptorCell_;
//- Acceptor processor number
label acceptorProcNo_;
//- Acceptor location
point acceptorPoint_;
// Donor side
//- Acceptor cell number
label donorCell_;
//- Donor processor number
label donorProcNo_;
//- Donor location
point donorPoint_;
//- Extended donor cell numbers
DynamicLabelList extendedDonorCells_;
//- Extended donor points, used to construct weights in different
// overset interpolation schemes
DynamicPointList extendedDonorPoints_;
// Note that the processor number is the same for all extended
// donors
public:
// Constructors
//- Construct null
inline donorAcceptor()
{}
//- Construct from acceptor data
inline donorAcceptor
(
const label acceptorCell,
const label acceptorProcNo,
const point& acceptorPoint
)
:
acceptorCell_(acceptorCell),
acceptorProcNo_(acceptorProcNo),
acceptorPoint_(acceptorPoint),
donorCell_(-1),
donorProcNo_(-1),
donorPoint_(vector::zero),
extendedDonorCells_(),
extendedDonorPoints_()
{}
//- Construct from Istream
inline donorAcceptor(Istream& is)
:
acceptorCell_(readLabel(is)),
acceptorProcNo_(readLabel(is)),
acceptorPoint_(is),
donorCell_(readLabel(is)),
donorProcNo_(readLabel(is)),
donorPoint_(is),
extendedDonorCells_(is),
extendedDonorPoints_(is)
{}
//- Copy constructor - default
//- Destructor - default
// Member Functions
// Acceptor side
//- Return acceptor cell number
label acceptorCell() const
{
return acceptorCell_;
}
//- Return access to acceptor cell number
label& acceptorCell()
{
return acceptorCell_;
}
//- Return acceptor processor number
label acceptorProcNo() const
{
return acceptorProcNo_;
}
//- Return access to acceptor processor number
label& acceptorProcNo()
{
return acceptorProcNo_;
}
//- Return acceptor location
const point& acceptorPoint() const
{
return acceptorPoint_;
}
// Donor side
//- Has a donor been found?
bool donorFound() const
{
return donorCell_ > -1;
}
//- Return donor cell number
label donorCell() const
{
return donorCell_;
}
//- Return access to donor cell number
label& donorCell()
{
return donorCell_;
}
//- Return donor processor number
label donorProcNo() const
{
return donorProcNo_;
}
//- Return access to donor processor number
label& donorProcNo()
{
return donorProcNo_;
}
//- Return donor point
const point& donorPoint() const
{
if (donorFound())
{
return donorPoint_;
}
else
{
return vector::max;
}
}
//- Return extended donor cell numbers
const DynamicLabelList& extendedDonorCells() const
{
return extendedDonorCells_;
}
//- Return extended donor cell centres
const DynamicPointList& extendedDonorPoints() const
{
return extendedDonorPoints_;
}
//- Return donor to acceptor distance
scalar distance() const
{
if (donorFound())
{
return mag(acceptorPoint_ - donorPoint_);
}
else
{
return GREAT;
}
}
// Edit
//- Set hit: donor found
void setDonor
(
const label& donorCell,
const label& donorProcNo,
const point& donorPoint
)
{
donorCell_ = donorCell;
donorProcNo_ = donorProcNo;
donorPoint_ = donorPoint;
}
//- Set extended donors by going through neighbours of currently set
// "best" donor, taking into account whether the cell is eligible
// donor (not acceptor or hole).
// Note 1: we do not check whether "best" donor is actually set for
// better performance.
// Note 2: setting extended donors is not necessary for injection
// interpolation, consider optimising for this case.
void setExtendedDonors
(
const boolList& eligibleDonors,
const labelListList& cCells,
const vectorField& cellCentres
)
{
// Get neighbours of this donor
const labelList& nbrs = cCells[donorCell_];
// Loop through neighbours and mark eligible ones
forAll (nbrs, nbrI)
{
const label& nbrCellI = nbrs[nbrI];
if (eligibleDonors[nbrCellI])
{
extendedDonorCells_.append(nbrCellI);
extendedDonorPoints_.append(cellCentres[nbrCellI]);
}
}
// Note: no need to shrink DynamicLists for parallel efficiency
// as IOstream operators first cast into List, thus sending only
// necessary data (not the capacity of the DynamicList).
}
//- Copy extended donors from another candidate donor
void setExtendedDonors(const donorAcceptor& rd)
{
extendedDonorCells_ = rd.extendedDonorCells_;
extendedDonorPoints_ = rd.extendedDonorPoints_;
}
//- Reject donor: used in overlap minimisation algorithms via Donor
// Suitability Function
void rejectDonor()
{
donorCell_ = -1;
donorProcNo_ = -1;
donorPoint_ = vector::zero;
}
// Member Operators
void operator=(const donorAcceptor& rd)
{
acceptorCell_ = rd.acceptorCell_;
acceptorProcNo_ = rd.acceptorProcNo_;
acceptorPoint_ = rd.acceptorPoint_;
donorCell_ = rd.donorCell_;
donorProcNo_ = rd.donorProcNo_;
donorPoint_ = rd.donorPoint_;
extendedDonorCells_ = rd.extendedDonorCells_;
extendedDonorPoints_ = rd.extendedDonorPoints_;
}
// Friend Operators
friend bool operator==
(
const donorAcceptor& a,
const donorAcceptor& b
)
{
return
a.acceptorCell_ == b.acceptorCell_
&& a.acceptorProcNo_ == b.acceptorProcNo_
&& a.acceptorPoint_ == b.acceptorPoint_
&& a.donorCell_ == b.donorCell_
&& a.donorProcNo_ == b.donorProcNo_;
// Note: do not check whether extended neighbours are the same, we
// assume they will be if donor data is the same
}
friend bool operator!=
(
const donorAcceptor& a,
const donorAcceptor& b
)
{
return !(a == b);
}
// IOstream Operators
//- Istream operator
friend Istream& operator>>(Istream& is, donorAcceptor& rd)
{
is >> rd.acceptorCell_
>> rd.acceptorProcNo_
>> rd.acceptorPoint_
>> rd.donorCell_
>> rd.donorProcNo_
>> rd.donorPoint_
>> rd.extendedDonorCells_
>> rd.extendedDonorPoints_;
return is;
}
//- Ostream operator
friend Ostream& operator<<(Ostream& os, const donorAcceptor& rd)
{
os << rd.acceptorCell_ << token::SPACE
<< rd.acceptorProcNo_ << token::SPACE
<< rd.acceptorPoint_ << token::SPACE
<< rd.donorCell_ << token::SPACE
<< rd.donorProcNo_ << token::SPACE
<< rd.donorPoint_ << token::SPACE
<< rd.extendedDonorCells_ << token::SPACE
<< rd.extendedDonorPoints_;
return os;
}
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //

View file

@ -0,0 +1,52 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | foam-extend: Open Source CFD
\\ / O peration | Version: 3.2
\\ / A nd | Web: http://www.foam-extend.org
\\/ M anipulation | For copyright notice see file Copyright
-------------------------------------------------------------------------------
License
This file is part of foam-extend.
foam-extend 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.
foam-extend 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 foam-extend. If not, see <http://www.gnu.org/licenses/>.
Typedef
Foam::donorAcceptorList
Description
container classes
\*---------------------------------------------------------------------------*/
#ifndef donorAcceptorList_H
#define donorAcceptorList_H
#include "donorAcceptor.H"
#include "List.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
typedef List<donorAcceptor> donorAcceptorList;
typedef List<donorAcceptorList> donorAcceptorListList;
typedef DynamicList<donorAcceptor> donorAcceptorDynamicList;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //

View file

@ -0,0 +1,546 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | foam-extend: Open Source CFD
\\ / O peration | Version: 4.0
\\ / A nd | Web: http://www.foam-extend.org
\\/ M anipulation | For copyright notice see file Copyright
-------------------------------------------------------------------------------
License
This file is part of foam-extend.
foam-extend 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.
foam-extend 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 foam-extend. If not, see <http://www.gnu.org/licenses/>.
\*---------------------------------------------------------------------------*/
#include "oversetEulerDdtScheme.H"
#include "oversetMesh.H"
#include "surfaceInterpolate.H"
#include "fvcDiv.H"
#include "fvMatrices.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace fv
{
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
template<class Type>
tmp<GeometricField<Type, fvPatchField, volMesh> >
oversetEulerDdtScheme<Type>::fvcDdt
(
const dimensioned<Type>& dt
)
{
// Get mesh
const fvMesh& mesh = mesh:
const dimensionedScalar rDeltaT = 1.0/mesh.time().deltaT();
// Create the ddt field, initialised to zero for dimensioned<Type>
tmp<GeometricField<Type, fvPatchField, volMesh> > tdtdt
(
new GeometricField<Type, fvPatchField, volMesh>
(
IOobject
(
"ddt(" + dt.name() + ')',
mesh.time().timeName(),
mesh
),
mesh,
dimensioned<Type>
(
"0",
dt.dimensions()/dimTime,
pTraits<Type>::zero
)
)
);
// Take into account moving mesh if the mesh is moving
if (mesh.moving())
{
tdtdt().internalField() =
rDeltaT.value()*dt.value()*(1.0 - mesh.V0()/mesh.V());
}
// Correct for overset "mesh motion" due to change in acceptor/live cells
// Get overset mesh
const oversetMesh& om = oversetMesh::New(mesh);
// Get cells that have been converted from acceptors to live
const labelList& accToLive = om.acceptorsToLive();
// Get the ddt field
scalarField dtdt = tdtdt().internalField();
// Correct for local volume gain (acceptors to live)
forAll (accToLive, i)
{
// Get cell index
const label cellI = accToLive[i];
// Calculate the contribution due to abrupt volume gain
dtdt[cellI] += rDeltaT.value()*dt.value();
}
// Not sure what to do with the cells that have been live and now they are
// acceptors. It doesn't really make sense to add their contribution to
// acceptors because the equations for acceptors are replaced with
// interpolation. Maybe we need to take into account this volume change
// (volume loss) in neighbouring live cells across fringe faces?
// VV, 31/Jul/2017.
}
template<class Type>
tmp<GeometricField<Type, fvPatchField, volMesh> >
oversetEulerDdtScheme<Type>::fvcDdt
(
const GeometricField<Type, fvPatchField, volMesh>& vf
)
{
dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
IOobject ddtIOobject
(
"ddt("+vf.name()+')',
mesh().time().timeName(),
mesh()
);
if (mesh().moving())
{
return tmp<GeometricField<Type, fvPatchField, volMesh> >
(
new GeometricField<Type, fvPatchField, volMesh>
(
ddtIOobject,
mesh(),
rDeltaT.dimensions()*vf.dimensions(),
rDeltaT.value()*
(
vf.internalField()
- vf.oldTime().internalField()*mesh().V0()/mesh().V()
),
rDeltaT.value()*
(
vf.boundaryField() - vf.oldTime().boundaryField()
)
)
);
}
else
{
return tmp<GeometricField<Type, fvPatchField, volMesh> >
(
new GeometricField<Type, fvPatchField, volMesh>
(
ddtIOobject,
rDeltaT*(vf - vf.oldTime())
)
);
}
}
template<class Type>
tmp<GeometricField<Type, fvPatchField, volMesh> >
oversetEulerDdtScheme<Type>::fvcDdt
(
const dimensionedScalar& rho,
const GeometricField<Type, fvPatchField, volMesh>& vf
)
{
dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
IOobject ddtIOobject
(
"ddt("+rho.name()+','+vf.name()+')',
mesh().time().timeName(),
mesh()
);
if (mesh().moving())
{
return tmp<GeometricField<Type, fvPatchField, volMesh> >
(
new GeometricField<Type, fvPatchField, volMesh>
(
ddtIOobject,
mesh(),
rDeltaT.dimensions()*rho.dimensions()*vf.dimensions(),
rDeltaT.value()*rho.value()*
(
vf.internalField()
- vf.oldTime().internalField()*mesh().V0()/mesh().V()
),
rDeltaT.value()*rho.value()*
(
vf.boundaryField() - vf.oldTime().boundaryField()
)
)
);
}
else
{
return tmp<GeometricField<Type, fvPatchField, volMesh> >
(
new GeometricField<Type, fvPatchField, volMesh>
(
ddtIOobject,
rDeltaT*rho*(vf - vf.oldTime())
)
);
}
}
template<class Type>
tmp<GeometricField<Type, fvPatchField, volMesh> >
oversetEulerDdtScheme<Type>::fvcDdt
(
const volScalarField& rho,
const GeometricField<Type, fvPatchField, volMesh>& vf
)
{
dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
IOobject ddtIOobject
(
"ddt("+rho.name()+','+vf.name()+')',
mesh().time().timeName(),
mesh()
);
if (mesh().moving())
{
return tmp<GeometricField<Type, fvPatchField, volMesh> >
(
new GeometricField<Type, fvPatchField, volMesh>
(
ddtIOobject,
mesh(),
rDeltaT.dimensions()*rho.dimensions()*vf.dimensions(),
rDeltaT.value()*
(
rho.internalField()*vf.internalField()
- rho.oldTime().internalField()
*vf.oldTime().internalField()*mesh().V0()/mesh().V()
),
rDeltaT.value()*
(
rho.boundaryField()*vf.boundaryField()
- rho.oldTime().boundaryField()
*vf.oldTime().boundaryField()
)
)
);
}
else
{
return tmp<GeometricField<Type, fvPatchField, volMesh> >
(
new GeometricField<Type, fvPatchField, volMesh>
(
ddtIOobject,
rDeltaT*(rho*vf - rho.oldTime()*vf.oldTime())
)
);
}
}
template<class Type>
tmp<fvMatrix<Type> >
oversetEulerDdtScheme<Type>::fvmDdt
(
const GeometricField<Type, fvPatchField, volMesh>& vf
)
{
tmp<fvMatrix<Type> > tfvm
(
new fvMatrix<Type>
(
vf,
vf.dimensions()*dimVol/dimTime
)
);
fvMatrix<Type>& fvm = tfvm();
scalar rDeltaT = 1.0/mesh().time().deltaT().value();
fvm.diag() = rDeltaT*mesh().V();
if (mesh().moving())
{
fvm.source() = rDeltaT*vf.oldTime().internalField()*mesh().V0();
}
else
{
fvm.source() = rDeltaT*vf.oldTime().internalField()*mesh().V();
}
return tfvm;
}
template<class Type>
tmp<fvMatrix<Type> >
oversetEulerDdtScheme<Type>::fvmDdt
(
const dimensionedScalar& rho,
const GeometricField<Type, fvPatchField, volMesh>& vf
)
{
tmp<fvMatrix<Type> > tfvm
(
new fvMatrix<Type>
(
vf,
rho.dimensions()*vf.dimensions()*dimVol/dimTime
)
);
fvMatrix<Type>& fvm = tfvm();
scalar rDeltaT = 1.0/mesh().time().deltaT().value();
fvm.diag() = rDeltaT*rho.value()*mesh().V();
if (mesh().moving())
{
fvm.source() = rDeltaT
*rho.value()*vf.oldTime().internalField()*mesh().V0();
}
else
{
fvm.source() = rDeltaT
*rho.value()*vf.oldTime().internalField()*mesh().V();
}
return tfvm;
}
template<class Type>
tmp<fvMatrix<Type> >
oversetEulerDdtScheme<Type>::fvmDdt
(
const volScalarField& rho,
const GeometricField<Type, fvPatchField, volMesh>& vf
)
{
tmp<fvMatrix<Type> > tfvm
(
new fvMatrix<Type>
(
vf,
rho.dimensions()*vf.dimensions()*dimVol/dimTime
)
);
fvMatrix<Type>& fvm = tfvm();
scalar rDeltaT = 1.0/mesh().time().deltaT().value();
fvm.diag() = rDeltaT*rho.internalField()*mesh().V();
if (mesh().moving())
{
fvm.source() = rDeltaT
*rho.oldTime().internalField()
*vf.oldTime().internalField()*mesh().V0();
}
else
{
fvm.source() = rDeltaT
*rho.oldTime().internalField()
*vf.oldTime().internalField()*mesh().V();
}
return tfvm;
}
template<class Type>
tmp<typename oversetEulerDdtScheme<Type>::fluxFieldType>
oversetEulerDdtScheme<Type>::fvcDdtPhiCorr
(
const volScalarField& rA,
const GeometricField<Type, fvPatchField, volMesh>& U,
const fluxFieldType& phiAbs
)
{
dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
IOobject ddtIOobject
(
"ddtPhiCorr(" + rA.name() + ',' + U.name() + ',' + phiAbs.name() + ')',
mesh().time().timeName(),
mesh()
);
tmp<fluxFieldType> phiCorr =
phiAbs.oldTime() - (fvc::interpolate(U.oldTime()) & mesh().Sf());
return tmp<fluxFieldType>
(
new fluxFieldType
(
ddtIOobject,
this->fvcDdtPhiCoeff(U.oldTime(), phiAbs.oldTime(), phiCorr())
*fvc::interpolate(rDeltaT*rA)*phiCorr
)
);
}
template<class Type>
tmp<typename oversetEulerDdtScheme<Type>::fluxFieldType>
oversetEulerDdtScheme<Type>::fvcDdtPhiCorr
(
const volScalarField& rA,
const volScalarField& rho,
const GeometricField<Type, fvPatchField, volMesh>& U,
const fluxFieldType& phiAbs
)
{
dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
IOobject ddtIOobject
(
"ddtPhiCorr("
+ rA.name() + ','
+ rho.name() + ','
+ U.name() + ','
+ phiAbs.name() + ')',
mesh().time().timeName(),
mesh()
);
if
(
U.dimensions() == dimVelocity
&& phiAbs.dimensions() == dimVelocity*dimArea
)
{
return tmp<fluxFieldType>
(
new fluxFieldType
(
ddtIOobject,
rDeltaT
*this->fvcDdtPhiCoeff(U.oldTime(), phiAbs.oldTime())
*(
fvc::interpolate(rA*rho.oldTime())*phiAbs.oldTime()
- (fvc::interpolate(rA*rho.oldTime()*U.oldTime())
& mesh().Sf())
)
)
);
}
else if
(
U.dimensions() == dimVelocity
&& phiAbs.dimensions() == rho.dimensions()*dimVelocity*dimArea
)
{
return tmp<fluxFieldType>
(
new fluxFieldType
(
ddtIOobject,
rDeltaT
*this->fvcDdtPhiCoeff
(
U.oldTime(),
phiAbs.oldTime()/fvc::interpolate(rho.oldTime())
)
*(
fvc::interpolate(rA*rho.oldTime())
*phiAbs.oldTime()/fvc::interpolate(rho.oldTime())
- (
fvc::interpolate
(
rA*rho.oldTime()*U.oldTime()
) & mesh().Sf()
)
)
)
);
}
else if
(
U.dimensions() == rho.dimensions()*dimVelocity
&& phiAbs.dimensions() == rho.dimensions()*dimVelocity*dimArea
)
{
return tmp<fluxFieldType>
(
new fluxFieldType
(
ddtIOobject,
rDeltaT
*this->fvcDdtPhiCoeff
(
rho.oldTime(),
U.oldTime(),
phiAbs.oldTime()
)
*(
fvc::interpolate(rA)*phiAbs.oldTime()
- (fvc::interpolate(rA*U.oldTime()) & mesh().Sf())
)
)
);
}
else
{
FatalErrorIn
(
"oversetEulerDdtScheme<Type>::fvcDdtPhiCorr"
) << "dimensions of phiAbs are not correct"
<< abort(FatalError);
return fluxFieldType::null();
}
}
template<class Type>
tmp<surfaceScalarField> oversetEulerDdtScheme<Type>::meshPhi
(
const GeometricField<Type, fvPatchField, volMesh>&
)
{
return mesh().phi();
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace fv
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam
// ************************************************************************* //

View file

@ -0,0 +1,209 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | foam-extend: Open Source CFD
\\ / O peration | Version: 4.0
\\ / A nd | Web: http://www.foam-extend.org
\\/ M anipulation | For copyright notice see file Copyright
-------------------------------------------------------------------------------
License
This file is part of foam-extend.
foam-extend 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.
foam-extend 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 foam-extend. If not, see <http://www.gnu.org/licenses/>.
Class
Foam::fv::oversetEulerDdtScheme
Description
Basic first-order oversetEuler implicit/explicit ddt using only the current
and previous time-step values.
The scheme ensures that the Space Conservation Law is satisfied when using
overset grids.
The whole hierarchy should be reorganized in a better way, e.g. this class
could be derived from Euler ddt scheme, which requires changes in the
finite volume library. Another thing we should consider is using the pimpl
idiom.
SourceFiles
oversetEulerDdtScheme.C
Author
Vuko Vukcevic, Wikki Ltd. All rights reserved.
\*---------------------------------------------------------------------------*/
#ifndef oversetEulerDdtScheme_H
#define oversetEulerDdtScheme_H
#include "ddtScheme.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace fv
{
/*---------------------------------------------------------------------------*\
Class oversetEulerDdtScheme Declaration
\*---------------------------------------------------------------------------*/
template<class Type>
class oversetEulerDdtScheme
:
public ddtScheme<Type>
{
// Private Member Functions
//- Disallow default bitwise copy construct
oversetEulerDdtScheme(const oversetEulerDdtScheme&);
//- Disallow default bitwise assignment
void operator=(const oversetEulerDdtScheme&);
public:
//- Runtime type information
TypeName("oversetEuler");
// Constructors
//- Construct from mesh
oversetEulerDdtScheme(const fvMesh& mesh)
:
ddtScheme<Type>(mesh)
{}
//- Construct from mesh and Istream
oversetEulerDdtScheme(const fvMesh& mesh, Istream& is)
:
ddtScheme<Type>(mesh, is)
{}
// Member Functions
//- Return mesh reference
const fvMesh& mesh() const
{
return fv::ddtScheme<Type>::mesh();
}
tmp<GeometricField<Type, fvPatchField, volMesh> > fvcDdt
(
const dimensioned<Type>&
);
tmp<GeometricField<Type, fvPatchField, volMesh> > fvcDdt
(
const GeometricField<Type, fvPatchField, volMesh>&
);
tmp<GeometricField<Type, fvPatchField, volMesh> > fvcDdt
(
const dimensionedScalar&,
const GeometricField<Type, fvPatchField, volMesh>&
);
tmp<GeometricField<Type, fvPatchField, volMesh> > fvcDdt
(
const volScalarField&,
const GeometricField<Type, fvPatchField, volMesh>&
);
tmp<fvMatrix<Type> > fvmDdt
(
const GeometricField<Type, fvPatchField, volMesh>&
);
tmp<fvMatrix<Type> > fvmDdt
(
const dimensionedScalar&,
const GeometricField<Type, fvPatchField, volMesh>&
);
tmp<fvMatrix<Type> > fvmDdt
(
const volScalarField&,
const GeometricField<Type, fvPatchField, volMesh>&
);
typedef typename ddtScheme<Type>::fluxFieldType fluxFieldType;
tmp<fluxFieldType> fvcDdtPhiCorr
(
const volScalarField& rA,
const GeometricField<Type, fvPatchField, volMesh>& U,
const fluxFieldType& phi
);
tmp<fluxFieldType> fvcDdtPhiCorr
(
const volScalarField& rA,
const volScalarField& rho,
const GeometricField<Type, fvPatchField, volMesh>& U,
const fluxFieldType& phi
);
tmp<surfaceScalarField> meshPhi
(
const GeometricField<Type, fvPatchField, volMesh>&
);
};
template<>
tmp<surfaceScalarField> oversetEulerDdtScheme<scalar>::fvcDdtPhiCorr
(
const volScalarField& rA,
const volScalarField& U,
const surfaceScalarField& phi
);
template<>
tmp<surfaceScalarField> oversetEulerDdtScheme<scalar>::fvcDdtPhiCorr
(
const volScalarField& rA,
const volScalarField& rho,
const volScalarField& U,
const surfaceScalarField& phi
);
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace fv
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#ifdef NoRepository
# include "oversetEulerDdtScheme.C"
#endif
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //

View file

@ -0,0 +1,39 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | foam-extend: Open Source CFD
\\ / O peration | Version: 4.0
\\ / A nd | Web: http://www.foam-extend.org
\\/ M anipulation | For copyright notice see file Copyright
-------------------------------------------------------------------------------
License
This file is part of foam-extend.
foam-extend 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.
foam-extend 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 foam-extend. If not, see <http://www.gnu.org/licenses/>.
\*---------------------------------------------------------------------------*/
#include "oversetEulerDdtScheme.H"
#include "fvMesh.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
namespace fv
{
makeFvDdtScheme(oversetEulerDdtScheme)
}
}
// ************************************************************************* //

View file

@ -0,0 +1,13 @@
//
// createOversetMask.H
// ~~~~~~~~~~~~
const oversetMesh& om = oversetMesh::New(mesh);
const volScalarField& cellOversetMask = om.gamma();
const surfaceScalarField& faceOversetMask = om.sGamma();
if (runTime.outputTime())
{
# include "writeOversetMasks.H"
}

View file

@ -0,0 +1,69 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2011 Wikki 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, write to the Free Software Foundation,
Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
Global
oversetAlphaCourantNo
Description
Calculates and outputs the mean and maximum Courant numbers for interface
cells with overset mesh correction.
\*---------------------------------------------------------------------------*/
scalar maxAlphaCo
(
readScalar(runTime.controlDict().lookup("maxAlphaCo"))
);
scalar alphaCoNum = 0.0;
scalar meanAlphaCoNum = 0.0;
surfaceScalarField alpha1f =
fvc::interpolate(min(max(alpha1, scalar(0)), scalar(1)));
const dimensionedScalar alphaOffset("alphaOffset", dimless, dAlpha);
if (mesh.nInternalFaces())
{
surfaceScalarField magAlphaPhi
(
pos(alpha1f - alphaOffset)*
pos(scalar(1) - alphaOffset - alpha1f)*
mag(faceOversetMask*phi)
);
surfaceScalarField SfUfbyDelta =
mesh.surfaceInterpolation::deltaCoeffs()*magAlphaPhi;
const scalar deltaT = runTime.deltaT().value();
alphaCoNum = max(SfUfbyDelta/mesh.magSf()).value()*deltaT;
meanCoNum = (sum(SfUfbyDelta)/sum(mesh.magSf())).value()*deltaT;
}
Info<< "Interface Courant Number mean: " << meanAlphaCoNum
<< " max: " << alphaCoNum << endl;
// ************************************************************************* //

View file

@ -0,0 +1,53 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright held by original author
\\/ 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 2 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, write to the Free Software Foundation,
Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
Global
oversetContinuityErrs
Description
Calculates and prints the continuity errors with
overset mesh correction.
\*---------------------------------------------------------------------------*/
{
volScalarField contErr = cellOversetMask*fvc::div(phi);
sumLocalContErr = runTime.deltaT().value()*
mag(contErr)().weightedAverage(mesh.V()).value();
globalContErr = runTime.deltaT().value()*
contErr.weightedAverage(mesh.V()).value();
cumulativeContErr += globalContErr;
Info<< "Overset time step continuity errors : "
<< "max local = " << gMax(contErr)
<< ", sum local = " << sumLocalContErr
<< ", global = " << globalContErr
<< ", cumulative = " << cumulativeContErr
<< endl;
}
// ************************************************************************* //

View file

@ -0,0 +1,59 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright held by original author
\\/ 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 2 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, write to the Free Software Foundation,
Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
Global
CourantNo
Description
Calculates and outputs the mean and maximum Courant Numbers with overset
mesh correction.
\*---------------------------------------------------------------------------*/
scalar CoNum = 0.0;
scalar meanCoNum = 0.0;
scalar velMag = 0.0;
if (mesh.nInternalFaces())
{
surfaceScalarField magPhi = mag(faceOversetMask*phi);
surfaceScalarField SfUfbyDelta =
mesh.surfaceInterpolation::deltaCoeffs()*magPhi;
CoNum = max(SfUfbyDelta/mesh.magSf())
.value()*runTime.deltaT().value();
meanCoNum = (sum(SfUfbyDelta)/sum(mesh.magSf()))
.value()*runTime.deltaT().value();
velMag = max(magPhi/mesh.magSf()).value();
}
Info<< "Courant Number mean: " << meanCoNum
<< " max: " << CoNum
<< " velocity magnitude: " << velMag
<< endl;
// ************************************************************************* //

View file

@ -0,0 +1,46 @@
{
const oversetMesh& om = oversetMesh::New(mesh);
om.gamma().write();
om.gammaExt().write();
om.sGamma().write();
om.oversetTypes().write();
// Create region ID field
volScalarField regionIndex
(
IOobject
(
"regionIndex",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::NO_WRITE
),
mesh,
dimless
);
scalarField& regionIndexIn = regionIndex.internalField();
forAll (om.regions(), regionI)
{
const cellZone& cz = om.regions()[regionI].zone();
forAll (cz, cellI)
{
regionIndexIn[cz[cellI]] = regionI;
}
}
// Update boundary values
forAll (regionIndex.boundaryField(), patchI)
{
if (!regionIndex.boundaryField()[patchI].coupled())
{
regionIndex.boundaryField()[patchI] =
regionIndex.boundaryField()[patchI].patchInternalField();
}
}
// regionIndex.write();
}

View file

@ -0,0 +1,222 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | foam-extend: Open Source CFD
\\ / O peration | Version: 3.2
\\ / A nd | Web: http://www.foam-extend.org
\\/ M anipulation | For copyright notice see file Copyright
-------------------------------------------------------------------------------
License
This file is part of foam-extend.
foam-extend 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.
foam-extend 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 foam-extend. If not, see <http://www.gnu.org/licenses/>.
\*---------------------------------------------------------------------------*/
#include "globalOversetAdjustPhi.H"
#include "volFields.H"
#include "surfaceFields.H"
#include "processorFvsPatchFields.H"
#include "inletOutletFvPatchFields.H"
#include "fvc.H"
#include "tolerancesSwitch.H"
#include "oversetMesh.H"
#include "emptyOversetFvPatchFields.H"
// * * * * * * * * * * * * * * * Global Functions * * * * * * * * * * * * * //
bool Foam::globalOversetAdjustPhi
(
surfaceScalarField& phi,
const volVectorField& U,
volScalarField& p
)
{
if (p.needReference())
{
p.boundaryField().updateCoeffs();
scalar massIn = 0.0;
scalar fixedMassOut = 0.0;
scalar adjustableMassOut = 0.0;
// If the mesh is moving, adjustment needs to be calculated on
// relative fluxes. HJ, 13/Feb/2009
if (phi.mesh().moving())
{
fvc::makeRelative(phi, U);
}
forAll (phi.boundaryField(), patchi)
{
const fvPatchVectorField& Up = U.boundaryField()[patchi];
const fvsPatchScalarField& phip = phi.boundaryField()[patchi];
// Bug fix: All coupled patches should also be unaffected
// HJ, 12/Feb/2010
if
(
!Up.coupled()
&& !isA<emptyOversetFvPatchVectorField>(Up)
)
{
if
(
Up.fixesValue()
&& !isA<inletOutletFvPatchVectorField>(Up)
)
{
forAll (phip, i)
{
if (phip[i] < 0.0)
{
massIn -= phip[i];
}
else
{
fixedMassOut += phip[i];
}
}
}
else
{
forAll(phip, i)
{
if (phip[i] < 0.0)
{
massIn -= phip[i];
}
else
{
adjustableMassOut += phip[i];
}
}
}
}
}
// HR 16.03.10: Bug fix for moving meshes. Changing domain volume needs
// to be taken into account.
if (phi.mesh().moving())
{
dimensionedScalar Vdiff =
sum(phi.mesh().V()) - sum(phi.mesh().V0());
fixedMassOut += Vdiff.value()/phi.time().deltaT().value();
}
reduce(massIn, sumOp<scalar>());
reduce(fixedMassOut, sumOp<scalar>());
reduce(adjustableMassOut, sumOp<scalar>());
scalar massCorr = 1.0;
static const debug::tolerancesSwitch closedDomainTol
(
"closedDomainTol",
1e-10
);
if (mag(adjustableMassOut) > SMALL)
{
massCorr = (massIn - fixedMassOut)/adjustableMassOut;
}
else if
(
mag(fixedMassOut - massIn)
> closedDomainTol()*Foam::max(1.0, mag(massIn))
)
{
phi.write();
U.write();
p.write();
// Cannot adjust
FatalErrorIn
(
"globalOversetAdjustPhi\n"
"(\n"
" surfaceScalarField& phi,\n"
" const volVectorField& U,\n"
" const volScalarField& p\n"
")"
) << "Continuity error cannot be removed by adjusting the"
" outflow.\nPlease check the velocity boundary conditions"
" and/or run potentialFoam to initialise the outflow." << nl
<< "Specified mass inflow : " << massIn << nl
<< "Specified mass outflow : " << fixedMassOut << nl
<< "Difference : "
<< mag(fixedMassOut - massIn) << nl
<< "Adjustable mass outflow : " << adjustableMassOut << nl
<< exit(FatalError);
}
if (oversetMesh::debug)
{
Info<< "bool Foam::globalOversetAdjustPhi(...) massIn: " << massIn
<< " fixedMassOut: " << fixedMassOut
<< " adjustableMassOut: " << adjustableMassOut
<< " mass corr: " << massCorr
<< endl;
}
forAll (phi.boundaryField(), patchi)
{
const fvPatchVectorField& Up = U.boundaryField()[patchi];
fvsPatchScalarField& phip = phi.boundaryField()[patchi];
// Bug fix: All coupled patches should also be unaffected
// HJ, 12/Feb/2010
if
(
!Up.coupled()
&& !isA<emptyOversetFvPatchVectorField>(Up)
)
{
if
(
!Up.fixesValue()
|| isA<inletOutletFvPatchVectorField>(Up)
)
{
forAll (phip, i)
{
if (phip[i] > 0.0)
{
phip[i] *= massCorr;
}
}
}
}
}
// If the mesh is moving, adjustment needs to be calculated on
// relative fluxes. Now reverting to absolute fluxes. HJ, 13/Feb/2009
if (phi.mesh().moving())
{
fvc::makeAbsolute(phi, U);
}
return mag(massIn) < SMALL
&& mag(fixedMassOut) < SMALL
&& mag(adjustableMassOut) < SMALL;
}
else
{
return false;
}
}
// ************************************************************************* //

View file

@ -0,0 +1,65 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | foam-extend: Open Source CFD
\\ / O peration | Version: 3.2
\\ / A nd | Web: http://www.foam-extend.org
\\/ M anipulation | For copyright notice see file Copyright
-------------------------------------------------------------------------------
License
This file is part of foam-extend.
foam-extend 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.
foam-extend 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 foam-extend. If not, see <http://www.gnu.org/licenses/>.
InNamespace
Foam
Description
Adjust the global balance of fluxes for cases which do not have a pressure
boundary. Global flux calculation/adjustment is done on all patches except
emptyOversetFvPatch.
The function returns true if the domain is closed.
SourceFiles
globalOversetAdjustPhi.C
\*---------------------------------------------------------------------------*/
#ifndef globalOversetAdjustPhi_H
#define globalOversetAdjustPhi_H
#include "volFieldsFwd.H"
#include "surfaceFieldsFwd.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
//- Adjust the balance of fluxes to obey continuity.
// For cases which do not have a pressure boundary.
// Return true if the domain is closed.
bool globalOversetAdjustPhi
(
surfaceScalarField& phi,
const volVectorField& U,
volScalarField& p
);
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //

View file

@ -0,0 +1,304 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright held by original author
\\/ 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 2 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, write to the Free Software Foundation,
Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
\*---------------------------------------------------------------------------*/
#include "oversetAdjustPhi.H"
#include "volFields.H"
#include "surfaceFields.H"
#include "oversetMesh.H"
#include "fvc.H"
#include "processorPolyPatch.H"
// * * * * * * * * * * * * * * * Global Functions * * * * * * * * * * * * * //
void Foam::oversetAdjustPhi
(
surfaceScalarField& phi,
const volVectorField& U
)
{
const fvMesh& mesh = phi.mesh();
// If the mesh is moving, adjustment needs to be calculated on
// relative fluxes. HJ, 13/Feb/2009
if (mesh.moving())
{
fvc::makeRelative(phi, U);
}
// Get overset mesh
const oversetMesh& om = oversetMesh::New(mesh);
// Get addressing to fringe faces
const labelList& fringeFaces = om.fringeFaces();
const boolList& fringeFaceFlips = om.fringeFaceFlips();
// Get internal owner-neighbour addressing
const unallocLabelList& owner = mesh.owner();
const unallocLabelList& neighbour = mesh.neighbour();
// Get region split to identify separate mesh components
const labelList& regionID = om.regionID();
// Sum up incoming and outgoing flux
scalar fringeIn = 0;
scalar fringeOut = 0;
// Adjust fluxes per region
scalarField& phiIn = phi.internalField();
// Note: fringe faces contain both internal and processor boundary
// faces. Make sure processor faces are not counted double.
// HJ, 1/May/2015
forAll (fringeFaces, ffI)
{
const label curFace = fringeFaces[ffI];
const bool curFlip = fringeFaceFlips[ffI];
if (mesh.isInternalFace(curFace))
{
// Internal face
// Debug check
if (regionID[owner[curFace]] != regionID[neighbour[curFace]])
{
FatalErrorIn
(
"void Foam::oversetAdjustPhi\n"
"(\n"
" surfaceScalarField& phi,\n"
" volVectorField& U\n"
")"
) << "Region index different for owner and neighbour "
<< "of face " << curFace
<< abort(FatalError);
}
// Get reference to the current face flux
const scalar& curPhi = phiIn[curFace];
if (curFlip)
{
if (curPhi > 0.0)
{
// Flux going out of the fringe.
fringeOut += curPhi;
}
else
{
// Flux coming into the fringe. Note reverse sign.
fringeIn -= curPhi;
}
}
else
{
if (curPhi > 0.0)
{
// Flux coming into the fringe.
fringeIn += curPhi;
}
else
{
// Flux going out of the fringe. Note reverse sign.
fringeOut -= curPhi;
}
}
}
else
{
// Processor boundary fringe face
// Find patch and face
const label patchI = mesh.boundaryMesh().whichPatch(curFace);
const label faceI = mesh.boundaryMesh()[patchI].whichFace(curFace);
if (patchI < 0)
{
FatalErrorIn
(
"void Foam::oversetAdjustPhi\n"
"(\n"
" surfaceScalarField& phi,\n"
" volVectorField& U\n"
")"
) << "Cannot find patch for fringe face" << curFace
<< abort(FatalError);
}
// Only account for processor face from the owner processor side
if (isA<processorPolyPatch>(mesh.boundaryMesh()[patchI]))
{
const processorPolyPatch& procPatch =
refCast<const processorPolyPatch>
(
mesh.boundaryMesh()[patchI]
);
if (procPatch.owner())
{
// Processor patch, master side
// Get reference to the current face flux
const scalar& curPhi = phi.boundaryField()[patchI][faceI];
if (curFlip)
{
if (curPhi > 0.0)
{
// Flux going out of the fringe.
fringeOut += curPhi;
}
else
{
// Flux coming into the fringe. Note reverse sign.
fringeIn -= curPhi;
}
}
else
{
if (curPhi > 0.0)
{
// Flux coming into the fringe.
fringeIn += curPhi;
}
else
{
// Flux going out of the fringe. Note reverse sign.
fringeOut -= curPhi;
}
}
}
}
else
{
FatalErrorIn
(
"void Foam::oversetAdjustPhi\n"
"(\n"
" surfaceScalarField& phi,\n"
" volVectorField& U\n"
")"
) << "Patch for fringe face" << curFace
<< " is not of processor type"
<< abort(FatalError);
}
}
}
// Do a global reduce of fluxes
reduce(fringeIn, sumOp<scalar>());
reduce(fringeOut, sumOp<scalar>());
// Calculate region flux correction
const scalar fluxScale = fringeIn/(fringeOut + SMALL);
if (oversetMesh::debug)
{
Info<< "Fringe balance for " << phi.name() << ": in = " << fringeIn
<< ", out = " << fringeOut
<< ", balance = " << fringeOut - fringeIn
<< ", fluxScale = " << fluxScale
<< endl;
}
// Go through all fringe faces on each region and balance the fluxes
forAll (fringeFaces, ffI)
{
const label curFace = fringeFaces[ffI];
const bool curFlip = fringeFaceFlips[ffI];
if (mesh.isInternalFace(curFace))
{
// Internal face
// Get reference to the flux for scaling
scalar& curPhi = phiIn[curFace];
// Scale outgoing flux to match the incoming one
if (curFlip && (curPhi > 0))
{
curPhi *= fluxScale;
}
else if (!curFlip && (curPhi < 0))
{
curPhi *= fluxScale;
}
}
else
{
// Processor boundary fringe face
// Find patch and face
const label patchI = mesh.boundaryMesh().whichPatch(curFace);
const label faceI = mesh.boundaryMesh()[patchI].whichFace(curFace);
// We need to scale both owner and neighbour fluxes because they
// represent the same flux
const processorPolyPatch& procPatch =
refCast<const processorPolyPatch>
(
mesh.boundaryMesh()[patchI]
);
// Get reference to the flux for scaling
scalar& curPhi = phi.boundaryField()[patchI][faceI];
if (procPatch.owner()) // Owner side
{
// Scale outgoing flux to match the incoming one
if (curFlip && (curPhi > 0))
{
curPhi *= fluxScale;
}
else if (!curFlip && (curPhi < 0))
{
curPhi *= fluxScale;
}
}
else // Neighbouring processor side
{
// Scale outgoing flux to match the incoming one
// Note: change in curPhi sign and curFlip
if (!curFlip && (curPhi < 0))
{
curPhi *= fluxScale;
}
else if (curFlip && (curPhi > 0))
{
curPhi *= fluxScale;
}
}
}
}
// If the mesh is moving, adjustment needs to be calculated on
// relative fluxes. Now reverting to absolute fluxes. HJ, 13/Feb/2009
if (mesh.moving())
{
fvc::makeAbsolute(phi, U);
}
}
// ************************************************************************* //

View file

@ -0,0 +1,60 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright held by original author
\\/ 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 2 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, write to the Free Software Foundation,
Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
InNamespace
Foam
Description
Adjust overset mesh fringe fluxes to obey continuity
SourceFiles
oversetAdjustPhi.C
\*---------------------------------------------------------------------------*/
#ifndef oversetAdjustPhi_H
#define oversetAdjustPhi_H
#include "volFieldsFwd.H"
#include "surfaceFieldsFwd.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
//- Adjust the fluxes on overset mesh fringe to obey continuity.
void oversetAdjustPhi
(
surfaceScalarField& phi,
const volVectorField& U
);
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //

View file

@ -0,0 +1,405 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright held by original author
\\/ 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 2 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, write to the Free Software Foundation,
Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
\*---------------------------------------------------------------------------*/
#include "regionWiseOversetAdjustPhi.H"
#include "volFields.H"
#include "surfaceFields.H"
#include "oversetMesh.H"
#include "fvc.H"
#include "processorPolyPatch.H"
#include "emptyOversetFvPatchFields.H"
// * * * * * * * * * * * * * * * Global Functions * * * * * * * * * * * * * //
void Foam::regionWiseOversetAdjustPhi
(
surfaceScalarField& phi,
const volVectorField& U
)
{
const fvMesh& mesh = phi.mesh();
// If the mesh is moving, adjustment needs to be calculated on
// relative fluxes. HJ, 13/Feb/2009
if (mesh.moving())
{
fvc::makeRelative(phi, U);
}
// Get overset mesh
const oversetMesh& om = oversetMesh::New(mesh);
// Get addressing to fringe faces
const labelList& fringeFaces = om.fringeFaces();
const boolList& fringeFaceFlips = om.fringeFaceFlips();
// Get internal owner-neighbour addressing
const unallocLabelList& owner = mesh.owner();
const unallocLabelList& neighbour = mesh.neighbour();
// Get region split to identify separate mesh components
const labelList& regionID = om.regionID();
// Incoming and outgoing region fluxes
scalarField regionIn(om.regions().size(), 0);
scalarField regionOut(om.regions().size(), 0);
// Incoming and outgoing region fluxes through the fringe
scalarField regionFringeIn(om.regions().size(), 0);
scalarField regionFringeOut(om.regions().size(), 0);
// Get flux internal field (for fringe faces later on)
scalarField& phiIn = phi.internalField();
// STAGE 1
// Loop through boundary fields and collect the incoming/outgoing fluxes
forAll (phi.boundaryField(), patchI)
{
// Get necessary references
const fvPatchVectorField& Up = U.boundaryField()[patchI];
const fvsPatchScalarField& phip = phi.boundaryField()[patchI];
const unallocLabelList& fc = Up.patch().faceCells();
// All coupled and emptyOverset patches should not be taken into account
if
(
!Up.coupled()
&& !isA<emptyOversetFvPatchVectorField>(Up)
)
{
// Note: no need to distinguish between fixed value fluxes or
// adjustable fluxes as only the fringe fluxes are going to be
// scaled
forAll (phip, i)
{
// Get current region index
const label curRegion = regionID[fc[i]];
if (phip[i] < 0.0)
{
regionIn[curRegion] -= phip[i];
}
else
{
regionOut[curRegion] += phip[i];
}
}
}
}
// STAGE 2
// Loop through fringe faces and collect incoming/outgoing fluxes
// Note: fringe faces contain both internal and processor boundary
// faces. Make sure processor faces are not counted twice.
forAll (fringeFaces, ffI)
{
// Get face information
const label curFace = fringeFaces[ffI];
const bool curFlip = fringeFaceFlips[ffI];
if (mesh.isInternalFace(curFace))
{
// Internal face
// Get region index
const label curRegion = regionID[owner[curFace]];
// Check whether owner and neighbour belong to the same region
if (curRegion != regionID[neighbour[curFace]])
{
FatalErrorIn
(
"void Foam::regionWiseOversetAdjustPhi\n"
"(\n"
" surfaceScalarField& phi,\n"
" volVectorField& U\n"
")"
) << "Region index different for owner and neighbour "
<< "of face " << curFace
<< abort(FatalError);
}
// Get current face flux
const scalar& curPhi = phiIn[curFace];
if (curFlip)
{
if (curPhi > 0.0)
{
// Flux going into the region (out of the fringe).
// Note that positive sign is kept.
regionFringeIn[curRegion] += curPhi;
}
else
{
// Flux coming out of the region (into the fringe).
// Note reverted sign.
regionFringeOut[curRegion] -= curPhi;
}
}
else
{
if (curPhi > 0.0)
{
// Flux going out of the region (into the fringe).
// Note that positive sign is kept.
regionFringeOut[curRegion] += curPhi;
}
else
{
// Flux going into the region (out of the fringe).
// Note reverted sign.
regionFringeIn[curRegion] -= curPhi;
}
}
}
else
{
// Processor boundary fringe face
// Find patch and face
const label patchI = mesh.boundaryMesh().whichPatch(curFace);
const label faceI = mesh.boundaryMesh()[patchI].whichFace(curFace);
// Sanity check
if (patchI < 0)
{
FatalErrorIn
(
"void Foam::regionWiseOversetAdjustPhi\n"
"(\n"
" surfaceScalarField& phi,\n"
" volVectorField& U\n"
")"
) << "Cannot find patch for fringe face" << curFace
<< abort(FatalError);
}
// Only account for processor face from the owner processor side
if (isA<processorPolyPatch>(mesh.boundaryMesh()[patchI]))
{
const processorPolyPatch& procPatch =
refCast<const processorPolyPatch>
(
mesh.boundaryMesh()[patchI]
);
if (procPatch.owner())
{
// Processor patch, master side
// Get region index
const label curRegion =
regionID[mesh.boundary()[patchI].faceCells()[faceI]];
// Get reference to the current face flux
const scalar& curPhi = phi.boundaryField()[patchI][faceI];
if (curFlip)
{
if (curPhi > 0.0)
{
// Flux going into the region (out of the fringe).
// Note that positive sign is kept.
regionFringeIn[curRegion] += curPhi;
}
else
{
// Flux coming out of the region (into the fringe).
// Note reverted sign.
regionFringeOut[curRegion] -= curPhi;
}
}
else
{
if (curPhi > 0.0)
{
// Flux going out of the region (into the fringe).
// Note that positive sign is kept.
regionFringeOut[curRegion] += curPhi;
}
else
{
// Flux going into the region (out of the fringe).
// Note reverted sign.
regionFringeIn[curRegion] -= curPhi;
}
}
}
}
else
{
FatalErrorIn
(
"void Foam::regionWiseOversetAdjustPhi\n"
"(\n"
" surfaceScalarField& phi,\n"
" volVectorField& U\n"
")"
) << "Patch for fringe face" << curFace
<< " is not of processor type"
<< abort(FatalError);
}
}
}
// Do a global reduce of fluxes
reduce(regionIn, sumOp<scalarField>());
reduce(regionOut, sumOp<scalarField>());
reduce(regionFringeIn, sumOp<scalarField>());
reduce(regionFringeOut, sumOp<scalarField>());
// This method will inevitably fail when regionFringeOut is zero, which can
// happen if we do not have fully enclosed mesh inside a background mesh
// (i.e. when we have a mesh which consists of left and right regions
// connected through a planar fringe layer)
if (min(regionFringeOut) < SMALL)
{
WarningIn
(
"void Foam::regionWiseOversetAdjustPhi\n"
"(\n"
" surfaceScalarField& phi,\n"
" volVectorField& U\n"
")"
) << "Outflow through fringe in a certain region is zero. "
<< "Skipping flux adjustment for: " << phi.name()
<< endl;
return;
}
// Calculate region flux correction. Note: only fluxes going out of the
// region through the fringe will be scaled.
scalarField regionFringeFluxScale =
(regionIn - regionOut + regionFringeIn)/regionFringeOut;
if (oversetMesh::debug)
{
const scalarField totalRegionIn = regionIn + regionFringeIn;
const scalarField totalRegionOut = regionOut + regionFringeOut;
Info<< "Region fringe balance for " << phi.name()
<< ": in = " << totalRegionIn
<< ", out = " << totalRegionOut
<< ", fringeIn = " << regionFringeIn
<< ", fringeOut = " << regionFringeOut
<< ", balance = " << totalRegionOut - totalRegionIn
<< ", fluxScale = " << regionFringeFluxScale
<< endl;
}
// Go through all fringe faces on each region and balance the fluxes
forAll (fringeFaces, ffI)
{
const label curFace = fringeFaces[ffI];
const bool curFlip = fringeFaceFlips[ffI];
if (mesh.isInternalFace(curFace))
{
// Internal face
// Get region index
const label curRegion = regionID[owner[curFace]];
// Get reference to the flux for scaling
scalar& curPhi = phiIn[curFace];
// Scale fringe flux going out of the region
if (curFlip && (curPhi < 0))
{
curPhi *= regionFringeFluxScale[curRegion];
}
else if (!curFlip && (curPhi > 0))
{
curPhi *= regionFringeFluxScale[curRegion];
}
}
else
{
// Processor boundary fringe face
// Find patch and face
const label patchI = mesh.boundaryMesh().whichPatch(curFace);
const label faceI = mesh.boundaryMesh()[patchI].whichFace(curFace);
// We need to scale both owner and neighbour fluxes because they
// represent the same flux
const processorPolyPatch& procPatch =
refCast<const processorPolyPatch>
(
mesh.boundaryMesh()[patchI]
);
// Get reference to the flux for scaling
scalar& curPhi = phi.boundaryField()[patchI][faceI];
if (procPatch.owner()) // Owner side
{
// Get region index
const label curRegion =
regionID[mesh.boundary()[patchI].faceCells()[faceI]];
// Scale fringe flux going out of the region
if (curFlip && (curPhi < 0))
{
curPhi *= regionFringeFluxScale[curRegion];
}
else if (!curFlip && (curPhi > 0))
{
curPhi *= regionFringeFluxScale[curRegion];
}
}
else // Neighbouring processor side
{
// Get region index
const label curRegion =
regionID[mesh.boundary()[patchI].faceCells()[faceI]];
// Scale fringe flux going out of the region
// Note: change in curPhi sign and curFlip
if (!curFlip && (curPhi > 0))
{
curPhi *= regionFringeFluxScale[curRegion];
}
else if (curFlip && (curPhi < 0))
{
curPhi *= regionFringeFluxScale[curRegion];
}
}
}
}
// If the mesh is moving, adjustment needs to be calculated on
// relative fluxes. Now reverting to absolute fluxes. HJ, 13/Feb/2009
if (mesh.moving())
{
fvc::makeAbsolute(phi, U);
}
}
// ************************************************************************* //

View file

@ -0,0 +1,65 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright held by original author
\\/ 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 2 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, write to the Free Software Foundation,
Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
InNamespace
Foam
Description
Adjust overset mesh fluxes to obey continuity region wise. All patches
except coupled and emptyOverset are visited and the ingoing/outgoing fluxes
are calculated correspondinly. The fringe fluxes are then scaled to obey
continuity in all regions.
SourceFiles
regionWiseOversetAdjustPhi.C
Author
Vuko Vukcevic, FMENA Zagreb. All rights reserved
\*---------------------------------------------------------------------------*/
#ifndef regionWiseOversetAdjustPhi_H
#define regionWiseOversetAdjustPhi_H
#include "volFieldsFwd.H"
#include "surfaceFieldsFwd.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
void regionWiseOversetAdjustPhi
(
surfaceScalarField& phi,
const volVectorField& U
);
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //

View file

@ -0,0 +1,252 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright held by original author
\\/ 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 2 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, write to the Free Software Foundation,
Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
\*---------------------------------------------------------------------------*/
#include "compositeFringe.H"
#include "oversetRegion.H"
#include "addToRunTimeSelectionTable.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
namespace Foam
{
defineTypeNameAndDebug(compositeFringe, 0);
addToRunTimeSelectionTable(oversetFringe, compositeFringe, dictionary);
}
// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
void Foam::compositeFringe::readBaseFringes(const dictionary& dict)
{
if (!baseFringes_.empty())
{
baseFringes_.clear();
}
PtrList<entry> baseFringeEntries(dict.lookup("baseFringes"));
baseFringes_.setSize(baseFringeEntries.size());
forAll (baseFringes_, bfI)
{
baseFringes_.set
(
bfI,
oversetFringe::New
(
this->mesh(),
this->region(),
baseFringeEntries[bfI].dict()
)
);
}
}
void Foam::compositeFringe::calcAddressing() const
{
if (fringeHolesPtr_ || acceptorsPtr_)
{
FatalErrorIn
(
"void Foam::compositeFringe::calcAddressing() const"
) << "Fringe addressing already calculated"
<< abort(FatalError);
}
// Get reference to region cell zone
const cellZone& rcz = region().zone();
// Make a hash set to collect acceptor points
labelHashSet acceptorSet;
// Make a hash set to collect hole points
labelHashSet fringeHoleSet;
// Go through all base fringes and record holes and acceptors
forAll (baseFringes_, bfI)
{
const labelList& ch = baseFringes_[bfI].fringeHoles();
forAll (ch, chI)
{
// Check if cell is in region zone
if (rcz.whichCell(ch[chI]) > -1)
{
// Found acceptor
acceptorSet.insert(ch[chI]);
}
}
const labelList& ca = baseFringes_[bfI].candidateAcceptors();
forAll (ca, caI)
{
// Check if cell is in region zone
if (rcz.whichCell(ca[caI]) > -1)
{
// Found acceptor
acceptorSet.insert(ca[caI]);
}
}
}
// Collect fringeHoles
fringeHolesPtr_ = new labelList(fringeHoleSet.sortedToc());
// Collect acceptors
acceptorsPtr_ = new labelList(acceptorSet.sortedToc());
}
void Foam::compositeFringe::clearAddressing() const
{
deleteDemandDrivenData(fringeHolesPtr_);
deleteDemandDrivenData(acceptorsPtr_);
deleteDemandDrivenData(finalDonorAcceptorsPtr_);
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
// Construct from dictionary
Foam::compositeFringe::compositeFringe
(
const fvMesh& mesh,
const oversetRegion& region,
const dictionary& dict
)
:
oversetFringe(mesh, region, dict),
fringeHolesPtr_(NULL),
acceptorsPtr_(NULL),
finalDonorAcceptorsPtr_(NULL)
{
// Read fringes
readBaseFringes(dict);
}
// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
Foam::compositeFringe::~compositeFringe()
{
clearAddressing();
}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
bool Foam::compositeFringe::updateIteration
(
donorAcceptorList& donorAcceptorRegionData
) const
{
// If the donorAcceptor list has been allocated, something went wrong with
// the iteration procedure (not-updated flag): this function has been called
// more than once, which should not happen for compositeFringe
if (finalDonorAcceptorsPtr_)
{
FatalErrorIn("compositeFringe::updateIteration(donorAcceptorList&")
<< "finalDonorAcceptorPtr_ already allocated. Something went "
<< "wrong with the iteration procedure (flag was not updated)."
<< nl << "This should not happen for compositeFringe."
<< abort(FatalError);
}
// Allocate the list by reusing the argument list
finalDonorAcceptorsPtr_ = new donorAcceptorList
(
donorAcceptorRegionData,
true
);
// Set the flag to true and return
updateSuitableOverlapFlag(true);
return foundSuitableOverlap();
}
const Foam::labelList& Foam::compositeFringe::fringeHoles() const
{
if (!fringeHolesPtr_)
{
calcAddressing();
}
return *fringeHolesPtr_;
}
const Foam::labelList& Foam::compositeFringe::candidateAcceptors() const
{
if (!acceptorsPtr_)
{
calcAddressing();
}
return *acceptorsPtr_;
}
Foam::donorAcceptorList& Foam::compositeFringe::finalDonorAcceptors() const
{
if (!finalDonorAcceptorsPtr_)
{
FatalErrorIn("compositeFringe::finalDonorAcceptors()")
<< "finalDonorAcceptorPtr_ not allocated. Make sure you have "
<< "called compositeFringe::updateIteration() before asking for "
<< "final set of donor/acceptor pairs."
<< abort(FatalError);
}
if (!foundSuitableOverlap())
{
FatalErrorIn("compositeFringe::finalDonorAcceptors()")
<< "Attemted to access finalDonorAcceptors but suitable overlap "
<< "has not been found. This is not allowed. "
<< abort(FatalError);
}
return *finalDonorAcceptorsPtr_;
}
void Foam::compositeFringe::update() const
{
Info<< "void compositeFringe::update() const" << endl;
forAll (baseFringes_, bfI)
{
// Update current base fringe
baseFringes_[bfI].update();
}
// Update flag for composite fringe
updateSuitableOverlapFlag(false);
}
// ************************************************************************* //

View file

@ -0,0 +1,150 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright held by original author
\\/ 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 2 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, write to the Free Software Foundation,
Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
Class
compositeFringe
Description
Composite fringe selection algorithm, consisting of multiple
fringes on a single region
Author
Hrvoje Jasak, Wikki Ltd. All rights reserved.
Contributor
Vuko Vukcevic, Wikki Ltd. All rights reserved.
SourceFiles
compositeFringe.C
\*---------------------------------------------------------------------------*/
#ifndef compositeFringe_H
#define compositeFringe_H
#include "oversetFringe.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
/*---------------------------------------------------------------------------*\
Class compositeFringe Declaration
\*---------------------------------------------------------------------------*/
class compositeFringe
:
public oversetFringe
{
// Private data
//- List of primitive selection algorithms
PtrList<oversetFringe> baseFringes_;
//- Fringe hole cells
mutable labelList* fringeHolesPtr_;
//- Acceptor cells
mutable labelList* acceptorsPtr_;
//- Final donor/acceptor pairs for this region (fringe)
mutable donorAcceptorList* finalDonorAcceptorsPtr_;
// Private Member Functions
//- Disallow default bitwise copy construct
compositeFringe(const compositeFringe&);
//- Disallow default bitwise assignment
void operator=(const compositeFringe&);
//- Read base fringes
void readBaseFringes(const dictionary& dict);
// Calculate hole and acceptor addressing
void calcAddressing() const;
// Clear addressing
void clearAddressing() const;
public:
//- Runtime type information
TypeName("composite");
// Constructors
//- Construct from dictionary
compositeFringe
(
const fvMesh& mesh,
const oversetRegion& region,
const dictionary& dict
);
// Destructor
virtual ~compositeFringe();
// Member Functions
//- Update iteration. Note: invalidates parameter
virtual bool updateIteration
(
donorAcceptorList& donorAcceptorRegionData
) const;
//- Return list of deactivated (hole) cells
// Fringe hole cells are collected in addition to geometric hole
// cells, which fall outside of all donor regions
virtual const labelList& fringeHoles() const;
//- Return list of acceptor cells
virtual const labelList& candidateAcceptors() const;
//- Return list of final donor acceptor pairs. Note: caller may
// invalidate finalDonorAcceptorsPtr_ for optimisation purposes
virtual donorAcceptorList& finalDonorAcceptors() const;
//- Update the fringe
virtual void update() const;
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //

View file

@ -0,0 +1,217 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright held by original author
\\/ 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 2 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, write to the Free Software Foundation,
Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
\*---------------------------------------------------------------------------*/
#include "donorBasedLayeredFringe.H"
#include "oversetRegion.H"
#include "polyPatchID.H"
#include "addToRunTimeSelectionTable.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
namespace Foam
{
defineTypeNameAndDebug(donorBasedLayeredFringe, 0);
addToRunTimeSelectionTable
(
oversetFringe,
donorBasedLayeredFringe,
dictionary
);
}
// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
void Foam::donorBasedLayeredFringe::calcAddressing() const
{
if (fringeHolesPtr_ || acceptorsPtr_)
{
FatalErrorIn
(
"void Foam::donorBasedLayeredFringe::calcAddressing() const"
) << "Fringe addressing already calculated"
<< abort(FatalError);
}
// We need to wait for the fringe from donor region to finish the
// assembly. Simply create empty acceptor and fringe hole lists
acceptorsPtr_ = new labelList();
fringeHolesPtr_ = new labelList();
}
void Foam::donorBasedLayeredFringe::clearAddressing() const
{
deleteDemandDrivenData(fringeHolesPtr_);
deleteDemandDrivenData(acceptorsPtr_);
deleteDemandDrivenData(finalDonorAcceptorsPtr_);
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
// Construct from dictionary
Foam::donorBasedLayeredFringe::donorBasedLayeredFringe
(
const fvMesh& mesh,
const oversetRegion& region,
const dictionary& dict
)
:
oversetFringe(mesh, region, dict),
fringeHolesPtr_(NULL),
acceptorsPtr_(NULL),
finalDonorAcceptorsPtr_(NULL),
nLayers_(readLabel(dict.lookup("nLayers")))
{
// Sanity check for number of layers
if (nLayers_ < 1)
{
FatalErrorIn
(
"donorBasedLayeredFringe::donorBasedLayeredFringe\n"
"(\n"
" const fvMesh& mesh,\n"
" const oversetRegion& region,\n"
" const dictionary& dict,\n"
")\n"
) << "Invalid number of layers specified in donor based layered"
<< "fringe dictionary."
<< nl
<< "Please specify value between greater than 0."
<< endl;
}
}
// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
Foam::donorBasedLayeredFringe::~donorBasedLayeredFringe()
{
clearAddressing();
}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
bool Foam::donorBasedLayeredFringe::updateIteration
(
donorAcceptorList& donorAcceptorRegionData
) const
{
if (finalDonorAcceptorsPtr_)
{
FatalErrorIn("donorBasedLayeredFringe::updateIteration(donorAcceptorList&")
<< "finalDonorAcceptorPtr_ already allocated. Something went "
<< "wrong with the iteration procedure (flag was not updated)."
<< abort(FatalError);
}
// Get overset mesh
const oversetMesh& om = region().overset();
// Get list of donor regions of this overset region
const labelList& donorRegions = region().donorRegions();
// Check whether all donor regions have finished with the donor search
// Allocate the list by reusing the argument list
finalDonorAcceptorsPtr_ = new donorAcceptorList
(
donorAcceptorRegionData,
true
);
// Set the flag to true and return
updateSuitableOverlapFlag(true);
return foundSuitableOverlap();
}
const Foam::labelList& Foam::donorBasedLayeredFringe::fringeHoles() const
{
if (!fringeHolesPtr_)
{
calcAddressing();
}
return *fringeHolesPtr_;
}
const Foam::labelList& Foam::donorBasedLayeredFringe::candidateAcceptors() const
{
if (!acceptorsPtr_)
{
calcAddressing();
}
return *acceptorsPtr_;
}
Foam::donorAcceptorList& Foam::donorBasedLayeredFringe::finalDonorAcceptors() const
{
if (!finalDonorAcceptorsPtr_)
{
FatalErrorIn("donorBasedLayeredFringe::finalDonorAcceptors()")
<< "finalDonorAcceptorPtr_ not allocated. Make sure you have "
<< "called donorBasedLayeredFringe::updateIteration() before asking for "
<< "final set of donor/acceptor pairs."
<< abort(FatalError);
}
if (!foundSuitableOverlap())
{
FatalErrorIn("donorBasedLayeredFringe::finalDonorAcceptors()")
<< "Attemted to access finalDonorAcceptors but suitable overlap "
<< "has not been found. This is not allowed. "
<< abort(FatalError);
}
return *finalDonorAcceptorsPtr_;
}
void Foam::donorBasedLayeredFringe::update() const
{
if (updateFringe_)
{
Info<< "donorBasedLayeredFringe::update() const" << endl;
// Clear out
clearAddressing();
}
// Set flag to false and clear final donor/acceptors only
deleteDemandDrivenData(finalDonorAcceptorsPtr_);
updateSuitableOverlapFlag(false);
}
// ************************************************************************* //

View file

@ -0,0 +1,148 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright held by original author
\\/ 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 2 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, write to the Free Software Foundation,
Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
Class
donorBasedLayeredFringeFringe
Description
Fringe algorithm based on layered approach, where the acceptors of the
region are defined as neighbours of the donors toward the interior of the
region.
This is achieved by waiting for the acceptors of the donor region to be
found and then defining acceptors.
Author
Vuko Vukcevic, Wikki Ltd. All rights reserved.
SourceFiles
donorBasedLayeredFringeFringe.C
\*---------------------------------------------------------------------------*/
#ifndef donorBasedLayeredFringeFringe_H
#define donorBasedLayeredFringeFringe_H
#include "oversetFringe.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
/*---------------------------------------------------------------------------*\
Class donorBasedLayeredFringeFringe Declaration
\*---------------------------------------------------------------------------*/
class donorBasedLayeredFringeFringe
:
public oversetFringe
{
// Private data
//- Fringe hole cells
mutable labelList* fringeHolesPtr_;
//- Acceptor cells
mutable labelList* acceptorsPtr_;
//- Final donor/acceptor pairs for this region (fringe)
mutable donorAcceptorList* finalDonorAcceptorsPtr_;
//- Number of layers from donors to acceptors of this region
label nLayers_;
// Private Member Functions
//- Disallow default bitwise copy construct
donorBasedLayeredFringeFringe(const donorBasedLayeredFringeFringe&);
//- Disallow default bitwise assignment
void operator=(const donorBasedLayeredFringeFringe&);
// Calculate hole and acceptor addressing
void calcAddressing() const;
// Clear addressing
void clearAddressing() const;
public:
//- Runtime type information
TypeName("donorBasedLayeredFringe");
// Constructors
//- Construct from dictionary
donorBasedLayeredFringeFringe
(
const fvMesh& mesh,
const oversetRegion& region,
const dictionary& dict
);
// Destructor
virtual ~donorBasedLayeredFringeFringe();
// Member Functions
//- Update iteration. Note: invalidates parameter
virtual bool updateIteration
(
donorAcceptorList& donorAcceptorRegionData
) const;
//- Return list of deactivated (hole) cells
// Fringe hole cells are collected in addition to geometric hole
// cells, which fall outside of all donor regions
virtual const labelList& fringeHoles() const;
//- Return list of acceptor cells
virtual const labelList& candidateAcceptors() const;
//- Return list of final donor acceptor pairs. Note: caller may
// invalidate finalDonorAcceptorsPtr_ for optimisation purposes
virtual donorAcceptorList& finalDonorAcceptors() const;
//- Update the fringe
virtual void update() const;
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //

View file

@ -0,0 +1,234 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright held by original author
\\/ 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 2 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, write to the Free Software Foundation,
Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
\*---------------------------------------------------------------------------*/
#include "faceCellsFringe.H"
#include "oversetRegion.H"
#include "polyPatchID.H"
#include "addToRunTimeSelectionTable.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
namespace Foam
{
defineTypeNameAndDebug(faceCellsFringe, 0);
addToRunTimeSelectionTable(oversetFringe, faceCellsFringe, dictionary);
}
// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
void Foam::faceCellsFringe::calcAddressing() const
{
if (acceptorsPtr_)
{
FatalErrorIn
(
"void Foam::faceCellsFringe::calcAddressing() const"
) << "Addressing already calculated"
<< abort(FatalError);
}
// Get reference to region cell zone
const cellZone& rcz = region().zone();
// Make a hash set to collect acceptor points
// Note: 2 patches connecting at the corner may create a duplicate,
// which is filtered on insertion
labelHashSet acceptorSet;
// Find patches and mark cells
forAll (patchNames_, nameI)
{
const polyPatchID curFringePatch
(
patchNames_[nameI],
mesh().boundaryMesh()
);
if (!curFringePatch.active())
{
FatalErrorIn
(
"void faceCellsFringe::calcAddressing() const"
) << "Fringe patch " << patchNames_[nameI]
<< " cannot be found"
<< abort(FatalError);
}
const unallocLabelList& curFaceCells =
mesh().boundaryMesh()[curFringePatch.index()].faceCells();
forAll (curFaceCells, fcI)
{
// Check if cell is in region zone
if (rcz.whichCell(curFaceCells[fcI]) > -1)
{
// Found acceptor
acceptorSet.insert(curFaceCells[fcI]);
}
}
}
// Collect acceptors
acceptorsPtr_ = new labelList(acceptorSet.sortedToc());
// Holes are empty for this fringe
fringeHolesPtr_ = new labelList();
}
void Foam::faceCellsFringe::clearAddressing() const
{
deleteDemandDrivenData(fringeHolesPtr_);
deleteDemandDrivenData(acceptorsPtr_);
deleteDemandDrivenData(finalDonorAcceptorsPtr_);
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
// Construct from dictionary
Foam::faceCellsFringe::faceCellsFringe
(
const fvMesh& mesh,
const oversetRegion& region,
const dictionary& dict
)
:
oversetFringe(mesh, region, dict),
patchNames_(dict.lookup("patches")),
fringeHolesPtr_(NULL),
acceptorsPtr_(NULL),
finalDonorAcceptorsPtr_(NULL),
updateFringe_
(
dict.lookupOrDefault<Switch>("updateAcceptors", false)
)
{}
// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
Foam::faceCellsFringe::~faceCellsFringe()
{
clearAddressing();
}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
bool Foam::faceCellsFringe::updateIteration
(
donorAcceptorList& donorAcceptorRegionData
) const
{
// If the donorAcceptor list has been allocated, something went wrong with
// the iteration procedure (not-updated flag): this function has been called
// more than once, which should not happen for faceCellsFringe
if (finalDonorAcceptorsPtr_)
{
FatalErrorIn("faceCellsFringe::updateIteration(donorAcceptorList&")
<< "finalDonorAcceptorPtr_ already allocated. Something went "
<< "wrong with the iteration procedure (flag was not updated)."
<< nl << "This should not happen for faceCellsFringe."
<< abort(FatalError);
}
// Allocate the list by reusing the argument list
finalDonorAcceptorsPtr_ = new donorAcceptorList
(
donorAcceptorRegionData,
true
);
// Set the flag to true and return
updateSuitableOverlapFlag(true);
return foundSuitableOverlap();
}
const Foam::labelList& Foam::faceCellsFringe::fringeHoles() const
{
if (!fringeHolesPtr_)
{
calcAddressing();
}
return *fringeHolesPtr_;
}
const Foam::labelList& Foam::faceCellsFringe::candidateAcceptors() const
{
if (!acceptorsPtr_)
{
calcAddressing();
}
return *acceptorsPtr_;
}
Foam::donorAcceptorList& Foam::faceCellsFringe::finalDonorAcceptors() const
{
if (!finalDonorAcceptorsPtr_)
{
FatalErrorIn("faceCellsFringe::finalDonorAcceptors()")
<< "finalDonorAcceptorPtr_ not allocated. Make sure you have "
<< "called faceCellsFringe::updateIteration() before asking for "
<< "final set of donor/acceptor pairs."
<< abort(FatalError);
}
if (!foundSuitableOverlap())
{
FatalErrorIn("faceCellsFringe::finalDonorAcceptors()")
<< "Attemted to access finalDonorAcceptors but suitable overlap "
<< "has not been found. This is not allowed. "
<< abort(FatalError);
}
return *finalDonorAcceptorsPtr_;
}
void Foam::faceCellsFringe::update() const
{
if (updateFringe_)
{
Info<< "faceCellsFringe::update() const" << endl;
// Clear out
clearAddressing();
}
// Set flag to false and clear final donor/acceptors only
deleteDemandDrivenData(finalDonorAcceptorsPtr_);
updateSuitableOverlapFlag(false);
}
// ************************************************************************* //

View file

@ -0,0 +1,150 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright held by original author
\\/ 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 2 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, write to the Free Software Foundation,
Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
Class
faceCellsFringe
Description
Face-cells fringe algorithm
Author
Hrvoje Jasak, Wikki Ltd. All rights reserved.
Contributor
Vuko Vukcevic, Wikki Ltd. All rights reserved.
SourceFiles
faceCellsFringe.C
\*---------------------------------------------------------------------------*/
#ifndef faceCellsFringe_H
#define faceCellsFringe_H
#include "oversetFringe.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
/*---------------------------------------------------------------------------*\
Class faceCellsFringe Declaration
\*---------------------------------------------------------------------------*/
class faceCellsFringe
:
public oversetFringe
{
// Private data
//- Patch names
wordList patchNames_;
//- Fringe hole cells
mutable labelList* fringeHolesPtr_;
//- Acceptor cells
mutable labelList* acceptorsPtr_;
//- Final donor/acceptor pairs for this region (fringe)
mutable donorAcceptorList* finalDonorAcceptorsPtr_;
//- Switch whether to update acceptors/donors for possible changes in
// patch cells (switched off by default)
const Switch updateFringe_;
// Private Member Functions
//- Disallow default bitwise copy construct
faceCellsFringe(const faceCellsFringe&);
//- Disallow default bitwise assignment
void operator=(const faceCellsFringe&);
// Calculate hole and acceptor addressing
void calcAddressing() const;
// Clear addressing
void clearAddressing() const;
public:
//- Runtime type information
TypeName("faceCells");
// Constructors
//- Construct from dictionary
faceCellsFringe
(
const fvMesh& mesh,
const oversetRegion& region,
const dictionary& dict
);
// Destructor
virtual ~faceCellsFringe();
// Member Functions
//- Update iteration. Note: invalidates parameter
virtual bool updateIteration
(
donorAcceptorList& donorAcceptorRegionData
) const;
//- Return list of deactivated (hole) cells
// Fringe hole cells are collected in addition to geometric hole
// cells, which fall outside of all donor regions
virtual const labelList& fringeHoles() const;
//- Return list of acceptor cells
virtual const labelList& candidateAcceptors() const;
//- Return list of final donor acceptor pairs. Note: caller may
// invalidate finalDonorAcceptorsPtr_ for optimisation purposes
virtual donorAcceptorList& finalDonorAcceptors() const;
//- Update the fringe
virtual void update() const;
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //

View file

@ -0,0 +1,248 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright held by original author
\\/ 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 2 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, write to the Free Software Foundation,
Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
\*---------------------------------------------------------------------------*/
#include "manualFringe.H"
#include "oversetRegion.H"
#include "foamTime.H"
#include "cellSet.H"
#include "addToRunTimeSelectionTable.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
namespace Foam
{
defineTypeNameAndDebug(manualFringe, 0);
addToRunTimeSelectionTable(oversetFringe, manualFringe, dictionary);
}
// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
void Foam::manualFringe::calcAddressing() const
{
if (fringeHolesPtr_ || acceptorsPtr_)
{
FatalErrorIn
(
"void Foam::manualFringe::calcAddressing() const"
) << "Fringe addressing already calculated"
<< abort(FatalError);
}
fringeHolesPtr_ = new labelList
(
cellSet
(
mesh(),
holesSetName_,
IOobject::MUST_READ,
IOobject::NO_WRITE
).toc()
);
acceptorsPtr_ = new labelList
(
cellSet
(
mesh(),
acceptorsSetName_,
IOobject::MUST_READ,
IOobject::NO_WRITE
).toc()
);
// Debug
// Get reference to region cell zone
const cellZone& rcz = region().zone();
// Check holes
const labelList& h = *fringeHolesPtr_;
forAll (h, holeI)
{
if (rcz.whichCell(h[holeI]) < 0)
{
FatalErrorIn
(
"void Foam::manualFringe::calcAddressing() const"
) << "Invalid hole cell for region " << region().name()
<< ": cell " << h[holeI] << " does not belong to this region"
<< abort(FatalError);
}
}
// Check acceptors
const labelList& a = *acceptorsPtr_;
forAll (a, accI)
{
if (rcz.whichCell(a[accI]) < 0)
{
FatalErrorIn
(
"void Foam::manualFringe::calcAddressing() const"
) << "Invalid acceptor cell for region " << region().name()
<< ": cell " << a[accI] << " does not belong to this region"
<< abort(FatalError);
}
}
}
void Foam::manualFringe::clearAddressing() const
{
deleteDemandDrivenData(fringeHolesPtr_);
deleteDemandDrivenData(acceptorsPtr_);
deleteDemandDrivenData(finalDonorAcceptorsPtr_);
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
// Construct from dictionary
Foam::manualFringe::manualFringe
(
const fvMesh& mesh,
const oversetRegion& region,
const dictionary& dict
)
:
oversetFringe(mesh, region, dict),
holesSetName_(dict.lookup("holes")),
acceptorsSetName_(dict.lookup("acceptors")),
fringeHolesPtr_(NULL),
acceptorsPtr_(NULL),
finalDonorAcceptorsPtr_(NULL),
updateFringe_
(
dict.lookupOrDefault<Switch>("updateAcceptors", false)
)
{}
// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
Foam::manualFringe::~manualFringe()
{
clearAddressing();
}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
bool Foam::manualFringe::updateIteration
(
donorAcceptorList& donorAcceptorRegionData
) const
{
// If the donorAcceptor list has been allocated, something went wrong with
// the iteration procedure (not-updated flag): this function has been called
// more than once, which should not happen for manualFringe
if (finalDonorAcceptorsPtr_)
{
FatalErrorIn("manualFringe::updateIteration(donorAcceptorList&)")
<< "finalDonorAcceptorPtr_ already allocated. Something went "
<< "wrong with the iteration procedure (flag was not updated)."
<< nl << "This should not happen for manualFringe."
<< abort(FatalError);
}
// Allocate the list by reusing the argument list
finalDonorAcceptorsPtr_ = new donorAcceptorList
(
donorAcceptorRegionData,
true
);
// Set the flag to true and return
updateSuitableOverlapFlag(true);
return foundSuitableOverlap();
}
const Foam::labelList& Foam::manualFringe::fringeHoles() const
{
if (!fringeHolesPtr_)
{
calcAddressing();
}
return *fringeHolesPtr_;
}
const Foam::labelList& Foam::manualFringe::candidateAcceptors() const
{
if (!acceptorsPtr_)
{
calcAddressing();
}
return *acceptorsPtr_;
}
Foam::donorAcceptorList& Foam::manualFringe::finalDonorAcceptors() const
{
if (!finalDonorAcceptorsPtr_)
{
FatalErrorIn("manualFringe::finalDonorAcceptors()")
<< "finalDonorAcceptorPtr_ not allocated. Make sure you have "
<< "called manualFringe::updateIteration() before asking for "
<< "final set of donor/acceptor pairs."
<< abort(FatalError);
}
if (!foundSuitableOverlap())
{
FatalErrorIn("manualFringe::finalDonorAcceptors()")
<< "Attemted to access finalDonorAcceptors but suitable overlap "
<< "has not been found. This is not allowed. "
<< abort(FatalError);
}
return *finalDonorAcceptorsPtr_;
}
void Foam::manualFringe::update() const
{
if (updateFringe_)
{
Info<< "manualFringe::update() const" << endl;
// Clear out
clearAddressing();
}
// Set flag to false and clear final donor/acceptors only
deleteDemandDrivenData(finalDonorAcceptorsPtr_);
updateSuitableOverlapFlag(false);
}
// ************************************************************************* //

View file

@ -0,0 +1,153 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright held by original author
\\/ 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 2 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, write to the Free Software Foundation,
Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
Class
manualFringe
Description
Manual fringe selection algorithm
Author
Hrvoje Jasak, Wikki Ltd. All rights reserved.
Contributor
Vuko Vukcevic, Wikki Ltd. All rights reserved.
SourceFiles
manualFringe.C
\*---------------------------------------------------------------------------*/
#ifndef manualFringe_H
#define manualFringe_H
#include "oversetFringe.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
/*---------------------------------------------------------------------------*\
Class manualFringe Declaration
\*---------------------------------------------------------------------------*/
class manualFringe
:
public oversetFringe
{
// Private data
//- Cell set of deactivated (hole) cell
word holesSetName_;
//- Cell set of acceptor cell
word acceptorsSetName_;
//- Hole cells
mutable labelList* fringeHolesPtr_;
//- Acceptor cells
mutable labelList* acceptorsPtr_;
//- Final donor/acceptor pairs for this region (fringe)
mutable donorAcceptorList* finalDonorAcceptorsPtr_;
//- Switch whether to update acceptors/donors for possible changes in
// cell sets (switched off by default)
const Switch updateFringe_;
// Private Member Functions
//- Disallow default bitwise copy construct
manualFringe(const manualFringe&);
//- Disallow default bitwise assignment
void operator=(const manualFringe&);
// Calculate hole and acceptor addressing
void calcAddressing() const;
// Clear addressing
void clearAddressing() const;
public:
//- Runtime type information
TypeName("manual");
// Constructors
//- Construct from dictionary
manualFringe
(
const fvMesh& mesh,
const oversetRegion& region,
const dictionary& dict
);
// Destructor
virtual ~manualFringe();
// Member Functions
//- Update iteration. Note: invalidates parameter
virtual bool updateIteration
(
donorAcceptorList& donorAcceptorRegionData
) const;
//- Return list of deactivated (hole) cells
// Fringe hole cells are collected in addition to geometric hole
// cells, which fall outside of all donor regions
virtual const labelList& fringeHoles() const;
//- Return list of acceptor cells
virtual const labelList& candidateAcceptors() const;
//- Return list of final donor acceptor pairs. Note: caller may
// invalidate finalDonorAcceptorsPtr_ for optimisation purposes
virtual donorAcceptorList& finalDonorAcceptors() const;
//- Update the fringe
virtual void update() const;
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //

View file

@ -0,0 +1,84 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright held by original author
\\/ 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 2 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, write to the Free Software Foundation,
Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
\*---------------------------------------------------------------------------*/
#include "cellBoundingBoxDiagonal.H"
#include "overlapFringe.H"
#include "oversetRegion.H"
#include "addToRunTimeSelectionTable.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
namespace donorSuitability
{
defineTypeNameAndDebug(cellBoundingBoxDiagonal, 0);
addToRunTimeSelectionTable
(
donorSuitability,
cellBoundingBoxDiagonal,
dictionary
);
}
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
Foam::donorSuitability::cellBoundingBoxDiagonal::cellBoundingBoxDiagonal
(
const overlapFringe& overlapFringeAlgorithm,
const dictionary& dict
)
:
donorSuitability(overlapFringeAlgorithm, dict)
{
// Get reference to fvMesh
const fvMesh& mesh = overlapFringeAlgorithm.mesh();
// Get necessary mesh data
const cellList& cells = mesh.cells();
const faceList& faces = mesh.faces();
const pointField& points = mesh.points();
// Create local donor suitability function
scalarField localDsf(mesh.nCells(), 0);
// Loop through cells and calculate the bounding box diagonal of each cell
forAll (cells, cellI)
{
const boundBox bb(cells[cellI].points(faces, points), false);
localDsf[cellI] = bb.mag();
}
// Combine donor suitability function data across processors for parallel
// run
this->combineDonorSuitabilityFunction(localDsf);
}
// ************************************************************************* //

View file

@ -0,0 +1,97 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright held by original author
\\/ 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 2 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, write to the Free Software Foundation,
Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
Class
Foam::donorSuitability::cellBoundingBoxDiagonal
Description
Cell bounding box diagonal suitability. Donor suitability function
corresponds to the diagonal of a bounding box encompassing a cell.
Minimises the overlap by shifting it towards cells with similar bounding box
diagonal (example: a tet cell with vertices on a hex cell faces (opposite
faces in a hex, have the same bounding box). This might be suitable for
unstructured tet meshes.
SourceFiles
cellBoundingBoxDiagonal.C
Author
Vuko Vukcevic, FMENA Zagreb. All rights reserved.
\*---------------------------------------------------------------------------*/
#ifndef cellBoundingBoxDiagonal_H
#define cellBoundingBoxDiagonal_H
#include "donorSuitability.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
namespace donorSuitability
{
/*---------------------------------------------------------------------------*\
Class cellBoundingBoxDiagonal Declaration
\*---------------------------------------------------------------------------*/
class cellBoundingBoxDiagonal
:
public donorSuitability
{
public:
//- Runtime type information
TypeName("cellBoundingBoxDiagonal");
// Constructors
//- Construct from components
cellBoundingBoxDiagonal
(
const overlapFringe& overlapFringeAlgorithm,
const dictionary& dict
);
//- Destructor
virtual ~cellBoundingBoxDiagonal()
{}
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace cellBoundingBoxDiagonal
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //

Some files were not shown because too many files have changed in this diff Show more