exaFOAM/foam-extend4.1-coherent-io
Archived
1
1
Fork 0
Code Issues Pull requests Projects Releases Packages Wiki Activity

Removed old immersed boundary solvers

Browse source
This commit is contained in:
Hrvoje Jasak 2017-12-30 09:24:50 +00:00
parent 2c9b34ac61
commit 13efece5db
34 changed files with 0 additions and 1682 deletions
Show stats Download patch file Download diff file Expand all files Collapse all files

49
applications/solvers/immersedBoundary/SOLVER_COOKBOOK
View file

@ -1,49 +0,0 @@
Make/options: add
-I$(LIB_SRC)/triSurface/lnInclude \
-I$(LIB_SRC)/meshTools/lnInclude \
-I../immersedBoundary/lnInclude
-ltriSurface \
-lmeshTools \
-L$(FOAM_USER_LIBBIN) -limmersedBoundary \
createFields.H: add
# include "createIbMasks.H"
solver:
1) add immersed boundary headers
#include "immersedBoundaryFvPatch.H"
#include "immersedBoundaryAdjustPhi.H"
2) on calculation of face fluxes (or their components), mask with faceIbMask
3) before to adjustPhi, add:
// Adjust immersed boundary fluxes
immersedBoundaryAdjustPhi(phi, U);
4) on explicit updates, add correctBoundaryConditions();
eg.
U = fvc::reconstruct(phi);
U.correctBoundaryConditions();
5) On reports of continuity error, add masking:
Info<< "IB-masked continuity error = "
<< mag(cellIbMask*fvc::div(phi))().weightedAverage(mesh.V()).value()
<< endl;
or use immersedBoundaryContinuityErrs.H
5) Chenge Courant number check to be IB-sensitive, using
immersedBoundaryCourantNo.H

3
applications/solvers/immersedBoundary/icoDyMIbFoam/Make/files
View file

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

15
applications/solvers/immersedBoundary/icoDyMIbFoam/Make/options
View file

@ -1,15 +0,0 @@
EXE_INC = \
-I$(LIB_SRC)/finiteVolume/lnInclude \
-I$(LIB_SRC)/meshTools/lnInclude \
-I$(LIB_SRC)/immersedBoundary/immersedBoundary/lnInclude \
-I$(LIB_SRC)/dynamicMesh/dynamicFvMesh/lnInclude \
-I$(LIB_SRC)/dynamicMesh/dynamicMesh/lnInclude \
-I./lnInclude
EXE_LIBS = \
-lfiniteVolume \
-lmeshTools \
-limmersedBoundary \
-ldynamicFvMesh \
-ldynamicMesh \
-llduSolvers

56
applications/solvers/immersedBoundary/icoDyMIbFoam/createFields.H
View file

@ -1,56 +0,0 @@
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"
label pRefCell = 0;
scalar pRefValue = 0.0;
setRefCell(p, pimple.dict(), pRefCell, pRefValue);
mesh.schemesDict().setFluxRequired(p.name());

152
applications/solvers/immersedBoundary/icoDyMIbFoam/icoDyMIbFoam.C
View file

@ -1,152 +0,0 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / 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 immersed boundary mesh support.
Author
Hrvoje Jasak, Wikki Ltd. All rights reserved
\*---------------------------------------------------------------------------*/
#include "fvCFD.H"
#include "dynamicFvMesh.H"
#include "immersedBoundaryFvPatch.H"
#include "immersedBoundaryAdjustPhi.H"
#include "pimpleControl.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
int main(int argc, char *argv[])
{
# include "setRootCase.H"
# include "createTime.H"
# include "createDynamicFvMesh.H"
pimpleControl pimple(mesh);
# include "createFields.H"
# include "initContinuityErrs.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Info<< "\nStarting time loop\n" << endl;
while (runTime.loop())
{
Info<< "Time = " << runTime.timeName() << nl << endl;
// Make the fluxes absolute
fvc::makeAbsolute(phi, U);
bool meshChanged = mesh.update();
reduce(meshChanged, orOp<bool>());
Info<< "Mesh update" << meshChanged << endl;
# include "createIbMasks.H"
// Make the fluxes relative to the mesh motion
fvc::makeRelative(phi, U);
# include "CourantNo.H"
// Pressure-velocity corrector
while (pimple.loop())
{
fvVectorMatrix UEqn
(
fvm::ddt(U)
+ fvm::div(phi, U)
- fvm::laplacian(nu, U)
);
if (pimple.momentumPredictor())
{
solve(UEqn == -fvc::grad(p));
}
// --- PISO loop
while (pimple.correct())
{
volScalarField rUA = 1.0/UEqn.A();
U = rUA*UEqn.H();
// Immersed boundary update
U.correctBoundaryConditions();
phi = faceIbMask*(fvc::interpolate(U) & mesh.Sf());
// Adjust immersed boundary fluxes
immersedBoundaryAdjustPhi(phi, U);
adjustPhi(phi, U, p);
// Non-orthogonal pressure corrector loop
while (pimple.correctNonOrthogonal())
{
fvScalarMatrix pEqn
(
fvm::laplacian(rUA, p) == fvc::div(phi)
);
pEqn.setReference(pRefCell, pRefValue);
pEqn.solve
(
mesh.solutionDict().solver
(
p.select(pimple.finalInnerIter())
)
);
if (pimple.finalNonOrthogonalIter())
{
phi -= pEqn.flux();
}
}
# include "immersedBoundaryContinuityErrs.H"
// Make the fluxes relative to the mesh motion
fvc::makeRelative(phi, U);
U -= rUA*fvc::grad(p);
U.correctBoundaryConditions();
}
}
runTime.write();
Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
<< " ClockTime = " << runTime.elapsedClockTime() << " s"
<< nl << endl;
}
Info<< "End\n" << endl;
return(0);
}
// ************************************************************************* //

3
applications/solvers/immersedBoundary/icoIbFoam/Make/files
View file

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

13
applications/solvers/immersedBoundary/icoIbFoam/Make/options
View file

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

56
applications/solvers/immersedBoundary/icoIbFoam/createFields.H
View file

@ -1,56 +0,0 @@
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"
label pRefCell = 0;
scalar pRefValue = 0.0;
setRefCell(p, pimple.dict(), pRefCell, pRefValue);
mesh.schemesDict().setFluxRequired(p.name());

139
applications/solvers/immersedBoundary/icoIbFoam/icoIbFoam.C
View file

@ -1,139 +0,0 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / 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
icoIbFoam
Description
Transient solver for incompressible, laminar flow of Newtonian fluids
with immersed boundary support.
Author
Hrvoje Jasak, Wikki Ltd. All rights reserved
\*---------------------------------------------------------------------------*/
#include "fvCFD.H"
#include "immersedBoundaryFvPatch.H"
#include "immersedBoundaryAdjustPhi.H"
#include "pimpleControl.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
int main(int argc, char *argv[])
{
# include "setRootCase.H"
# include "createTime.H"
# include "createMesh.H"
pimpleControl pimple(mesh);
# include "createIbMasks.H"
# include "createFields.H"
# include "initContinuityErrs.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Info<< "\nStarting time loop\n" << endl;
while (runTime.loop())
{
Info<< "Time = " << runTime.timeName() << nl << endl;
# include "CourantNo.H"
// Pressure-velocity corrector
while (pimple.loop())
{
fvVectorMatrix UEqn
(
fvm::ddt(U)
+ fvm::div(phi, U)
- fvm::laplacian(nu, U)
);
if (pimple.momentumPredictor())
{
solve(UEqn == -fvc::grad(p));
}
// --- PISO loop
while (pimple.correct())
{
volScalarField rUA = 1.0/UEqn.A();
U = rUA*UEqn.H();
// Immersed boundary update
U.correctBoundaryConditions();
phi = faceIbMask*(fvc::interpolate(U) & mesh.Sf());
// Adjust immersed boundary fluxes
immersedBoundaryAdjustPhi(phi, U);
adjustPhi(phi, U, p);
// Non-orthogonal pressure corrector loop
while (pimple.correctNonOrthogonal())
{
fvScalarMatrix pEqn
(
fvm::laplacian(rUA, p) == fvc::div(phi)
);
pEqn.setReference(pRefCell, pRefValue);
pEqn.solve
(
mesh.solutionDict().solver
(
p.select(pimple.finalInnerIter())
)
);
if (pimple.finalNonOrthogonalIter())
{
phi -= pEqn.flux();
}
}
# include "immersedBoundaryContinuityErrs.H"
U -= rUA*fvc::grad(p);
U.correctBoundaryConditions();
}
}
runTime.write();
Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
<< " ClockTime = " << runTime.elapsedClockTime() << " s"
<< nl << endl;
}
Info<< "End\n" << endl;
return 0;
}
// ************************************************************************* //

3
applications/solvers/immersedBoundary/interIbFoam/Make/files
View file

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

22
applications/solvers/immersedBoundary/interIbFoam/Make/options
View file

@ -1,22 +0,0 @@
EXE_INC = \
-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)/meshTools/lnInclude \
-I$(LIB_SRC)/immersedBoundary/immersedBoundary/lnInclude
EXE_LIBS = \
-linterfaceProperties \
-lincompressibleTransportModels \
-lincompressibleTurbulenceModel \
-lincompressibleRASModels \
-lincompressibleLESModels \
-lfiniteVolume \
-lmeshTools \
-lsurfMesh \
-lsampling \
-ldynamicMesh \
-limmersedBoundary \
-llduSolvers

27
applications/solvers/immersedBoundary/interIbFoam/UEqn.H
View file

@ -1,27 +0,0 @@
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
==
- gh*fvc::grad(rho)
- fvc::grad(pd)
);
}

56
applications/solvers/immersedBoundary/interIbFoam/alphaEqn.H
View file

@ -1,56 +0,0 @@
{
word alphaScheme("div(phi,alpha)");
word alpharScheme("div(phirb,alpha)");
// New formulation of phir: Co condition
surfaceScalarField phir
(
"phir",
faceIbMask*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();
Info<< "Liquid phase volume fraction = "
<< alpha1.weightedAverage(mesh.V()).value()
<< " Min(alpha1) = " << min(cellIbMask*alpha1).value()
<< " Max(alpha1) = " << max(cellIbMask*alpha1).value()
<< endl;
rhoPhi = faceIbMask*(alpha1Eqn.flux()*(rho1 - rho2) + phi*rho2);
// Limit alpha to handle IB cells
// alpha1.max(0);
// alpha1.min(1);
// Update of interface and density
interface.correct();
twoPhaseProperties.correct();
// Calculate density using limited alpha1
rho = twoPhaseProperties.rho();
rho.correctBoundaryConditions();
}

58
applications/solvers/immersedBoundary/interIbFoam/correctPhi.H
View file

@ -1,58 +0,0 @@
{
# include "continuityErrs.H"
wordList pcorrTypes
(
pd.boundaryField().size(),
zeroGradientFvPatchScalarField::typeName
);
for (label i=0; i<pd.boundaryField().size(); i++)
{
if (pd.boundaryField()[i].fixesValue())
{
pcorrTypes[i] = fixedValueFvPatchScalarField::typeName;
}
}
volScalarField pcorr
(
IOobject
(
"pcorr",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::NO_WRITE
),
mesh,
dimensionedScalar("pcorr", pd.dimensions(), 0.0),
pcorrTypes
);
dimensionedScalar rUAf("(1|A(U))", dimTime/rho.dimensions(), 1.0);
// Adjust immersed boundary fluxes
immersedBoundaryAdjustPhi(phi, U);
adjustPhi(phi, U, pcorr);
mesh.schemesDict().setFluxRequired(pcorr.name());
while(pimple.correctNonOrthogonal())
{
fvScalarMatrix pcorrEqn
(
fvm::laplacian(rUAf, pcorr) == fvc::div(phi)
);
pcorrEqn.setReference(pdRefCell, pdRefValue);
pcorrEqn.solve();
if (pimple.finalNonOrthogonalIter())
{
phi -= pcorrEqn.flux();
}
}
# include "immersedBoundaryContinuityErrs.H"
}

118
applications/solvers/immersedBoundary/interIbFoam/createFields.H
View file

@ -1,118 +0,0 @@
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
);
label pdRefCell = 0;
scalar pdRefValue = 0.0;
setRefCell(pd, pimple.dict(), pdRefCell, pdRefValue);
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)
);

124
applications/solvers/immersedBoundary/interIbFoam/interIbFoam.C
View file

@ -1,124 +0,0 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / 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
interIbFoam
Description
Solver for 2 incompressible, isothermal immiscible fluids using a VOF
(volume of fluid) phase-fraction based interface capturing approach,
with immersed boundary 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 "immersedBoundaryFvPatch.H"
#include "immersedBoundaryAdjustPhi.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 "createFields.H"
# include "createIbMasks.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 "immersedBoundaryCourantNo.H"
# include "setDeltaT.H"
runTime++;
Info<< "Time = " << runTime.timeName() << nl << endl;
// Pressure-velocity corrector
while (pimple.loop())
{
twoPhaseProperties.correct();
# include "UEqn.H"
// --- PISO loop
while (pimple.correct())
{
# include "pEqn.H"
}
# include "immersedBoundaryContinuityErrs.H"
# include "limitU.H"
p = pd + rho*gh;
p.correctBoundaryConditions();
# include "alphaEqn.H"
turbulence->correct();
}
runTime.write();
Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
<< " ClockTime = " << runTime.elapsedClockTime() << " s"
<< nl << endl;
}
Info<< "End\n" << endl;
return 0;
}
// ************************************************************************* //

24
applications/solvers/immersedBoundary/interIbFoam/limitU.H
View file

@ -1,24 +0,0 @@
{
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();
Info << " ...limiting" << endl;
}
else
{
Info << endl;
}
}

53
applications/solvers/immersedBoundary/interIbFoam/pEqn.H
View file

@ -1,53 +0,0 @@
{
volScalarField rUA = 1.0/UEqn.A();
surfaceScalarField rUAf = fvc::interpolate(rUA);
U = rUA*UEqn.H();
// Immersed boundary update
U.correctBoundaryConditions();
# include "limitU.H"
surfaceScalarField phiU
(
"phiU",
faceIbMask*(fvc::interpolate(U) & mesh.Sf())
);
phi = phiU
+ faceIbMask*
(
fvc::interpolate(interface.sigmaK())*fvc::snGrad(alpha1)
- ghf*fvc::snGrad(rho)
)*rUAf*mesh.magSf();
// Adjust immersed boundary fluxes
immersedBoundaryAdjustPhi(phi, U);
adjustPhi(phi, U, pd);
while (pimple.correctNonOrthogonal())
{
fvScalarMatrix pdEqn
(
fvm::laplacian(rUAf, pd, "laplacian(rAU,pd)") == fvc::div(phi)
);
pdEqn.setReference(pdRefCell, pdRefValue);
pdEqn.solve
(
mesh.solutionDict().solver(pd.select(pimple.finalInnerIter()))
);
if (pimple.finalNonOrthogonalIter())
{
phi -= pdEqn.flux();
}
}
// Explicitly relax pressure
pd.relax();
U += rUA*fvc::reconstruct((phi - phiU)/rUAf);
U.correctBoundaryConditions();
}

3
applications/solvers/immersedBoundary/porousSimpleIbFoam/Make/files
View file

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

20
applications/solvers/immersedBoundary/porousSimpleIbFoam/Make/options
View file

@ -1,20 +0,0 @@
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)/meshTools/lnInclude \
-I$(LIB_SRC)/immersedBoundary/immersedBoundary/lnInclude
EXE_LIBS = \
-lincompressibleRASModels \
-lincompressibleTransportModels \
-lfiniteVolume \
-lmeshTools \
-lsurfMesh \
-lsampling \
-ldynamicMesh \
-limmersedBoundary \
-limmersedBoundaryTurbulence \
-llduSolvers

12
applications/solvers/immersedBoundary/porousSimpleIbFoam/UEqn.H
View file

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

44
applications/solvers/immersedBoundary/porousSimpleIbFoam/createFields.H
View file

@ -1,44 +0,0 @@
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"
label pRefCell = 0;
scalar pRefValue = 0.0;
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)
);
porousZones pZones(mesh);

41
applications/solvers/immersedBoundary/porousSimpleIbFoam/pEqn.H
View file

@ -1,41 +0,0 @@
{
volScalarField AU = UEqn().A();
U = UEqn().H()/AU;
// Immersed boundary update
U.correctBoundaryConditions();
UEqn.clear();
phi = faceIbMask*(fvc::interpolate(U) & mesh.Sf());
// Adjust immersed boundary fluxes
immersedBoundaryAdjustPhi(phi, U);
adjustPhi(phi, U, p);
// Non-orthogonal pressure corrector loop
while (simple.correctNonOrthogonal())
{
fvScalarMatrix pEqn
(
fvm::laplacian(1.0/AU, p) == fvc::div(phi)
);
pEqn.setReference(pRefCell, pRefValue);
pEqn.solve();
if (simple.finalNonOrthogonalIter())
{
phi -= pEqn.flux();
}
}
# include "immersedBoundaryContinuityErrs.H"
// Explicitly relax pressure for momentum corrector
p.relax();
// Momentum corrector
U -= fvc::grad(p)/AU;
U.correctBoundaryConditions();
}

88
applications/solvers/immersedBoundary/porousSimpleIbFoam/porousSimpleIbFoam.C
View file

@ -1,88 +0,0 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / 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
simpleIbFoam
Description
Steady-state solver for incompressible, turbulent flow
with porous zones and immersed boundary support.
Author
Hrvoje Jasak, Wikki Ltd. All rights reserved
\*---------------------------------------------------------------------------*/
#include "fvCFD.H"
#include "singlePhaseTransportModel.H"
#include "RASModel.H"
#include "porousZones.H"
#include "simpleControl.H"
#include "immersedBoundaryFvPatch.H"
#include "immersedBoundaryAdjustPhi.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
int main(int argc, char *argv[])
{
# include "setRootCase.H"
# include "createTime.H"
# include "createMesh.H"
simpleControl simple(mesh);
# include "createIbMasks.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;
}
// ************************************************************************* //

3
applications/solvers/immersedBoundary/potentialIbFoam/Make/files
View file

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

16
applications/solvers/immersedBoundary/potentialIbFoam/Make/options
View file

@ -1,16 +0,0 @@
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)/immersedBoundary/immersedBoundary/lnInclude
EXE_LIBS = \
-lfiniteVolume \
-lmeshTools \
-lsurfMesh \
-lsampling \
-ldynamicMesh \
-ldynamicFvMesh \
-limmersedBoundary \
-llduSolvers

50
applications/solvers/immersedBoundary/potentialIbFoam/createFields.H
View file

@ -1,50 +0,0 @@
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
);
// Do not reset U
surfaceScalarField phi
(
IOobject
(
"phi",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
fvc::interpolate(U) & mesh.Sf()
);
label pRefCell = 0;
scalar pRefValue = 0.0;
setRefCell(p, simple.dict(), pRefCell, pRefValue);
mesh.schemesDict().setFluxRequired(p.name());

231
applications/solvers/immersedBoundary/potentialIbFoam/potentialIbFoam.C
View file

@ -1,231 +0,0 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / 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
potentialIbFoam
Description
Potential flow solver with immersed boundary support.
\*---------------------------------------------------------------------------*/
#include "fvCFD.H"
#include "immersedBoundaryFvPatch.H"
#include "immersedBoundaryAdjustPhi.H"
#include "simpleControl.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
int main(int argc, char *argv[])
{
argList::validOptions.insert("writep", "");
# include "setRootCase.H"
# include "createTime.H"
# include "createMesh.H"
simpleControl simple(mesh);
# include "createIbMasks.H"
# include "createFields.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Info<< nl << "Calculating potential flow" << endl;
// Do correctors over the complete set
while (simple.correctNonOrthogonal())
{
phi = faceIbMask*(linearInterpolate(U) & mesh.Sf());
// Adjust immersed boundary fluxes
immersedBoundaryAdjustPhi(phi, U);
// Adjust fluxes
adjustPhi(phi, U, p);
p.storePrevIter();
fvScalarMatrix pEqn
(
fvm::laplacian
(
dimensionedScalar
(
"1",
dimTime/p.dimensions()*dimensionSet(0, 2, -2, 0, 0),
1
),
p
)
==
fvc::div(phi)
);
pEqn.setReference(pRefCell, pRefValue);
pEqn.solve();
// Correct the flux
phi -= pEqn.flux();
if (!simple.finalNonOrthogonalIter())
{
p.relax();
}
Info<< "p min " << gMin(p.internalField())
<< " max " << gMax(p.internalField())
<< " masked min "
<< gMin(cellIbMask.internalField()*p.internalField())
<< " max "
<< gMax(cellIbMask.internalField()*p.internalField())
<< endl;
Info<< "continuity error = "
<< mag
(
fvc::div(faceIbMask*phi)
)().weightedAverage(mesh.V()).value()
<< endl;
Info<< "Contour continuity error = "
<< mag(sum(phi.boundaryField()))
<< endl;
U = fvc::reconstruct(phi);
U.correctBoundaryConditions();
Info<< "Interpolated U error = "
<< (
sqrt
(
sum
(
sqr
(
faceIbMask*
(
fvc::interpolate(U) & mesh.Sf()
)
- phi
)
)
)/sum(mesh.magSf())
).value()
<< endl;
}
// Calculate velocity magnitude
{
volScalarField magU = cellIbMask*mag(U);
Info << "IB-masked mag(U): max: " << gMax(magU.internalField())
<< " min: " << gMin(magU.internalField()) << endl;
}
// Force the write
U.write();
phi.write();
cellIbMask.write();
cellIbMaskExt.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();
}
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;
}
// ************************************************************************* //

3
applications/solvers/immersedBoundary/simpleIbFoam/Make/files
View file

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

20
applications/solvers/immersedBoundary/simpleIbFoam/Make/options
View file

@ -1,20 +0,0 @@
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)/meshTools/lnInclude \
-I$(LIB_SRC)/immersedBoundary/immersedBoundary/lnInclude
EXE_LIBS = \
-lincompressibleRASModels \
-lincompressibleTransportModels \
-lfiniteVolume \
-lmeshTools \
-lsurfMesh \
-lsampling \
-ldynamicMesh \
-limmersedBoundary \
-limmersedBoundaryTurbulence \
-llduSolvers

10
applications/solvers/immersedBoundary/simpleIbFoam/UEqn.H
View file

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

42
applications/solvers/immersedBoundary/simpleIbFoam/createFields.H
View file

@ -1,42 +0,0 @@
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"
label pRefCell = 0;
scalar pRefValue = 0.0;
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)
);

41
applications/solvers/immersedBoundary/simpleIbFoam/pEqn.H
View file

@ -1,41 +0,0 @@
{
volScalarField AU = UEqn().A();
U = UEqn().H()/AU;
// Immersed boundary update
U.correctBoundaryConditions();
UEqn.clear();
phi = faceIbMask*(fvc::interpolate(U) & mesh.Sf());
// Adjust immersed boundary fluxes
immersedBoundaryAdjustPhi(phi, U);
adjustPhi(phi, U, p);
// Non-orthogonal pressure corrector loop
while (simple.correctNonOrthogonal())
{
fvScalarMatrix pEqn
(
fvm::laplacian(1.0/AU, p) == fvc::div(phi)
);
pEqn.setReference(pRefCell, pRefValue);
pEqn.solve();
if (simple.finalNonOrthogonalIter())
{
phi -= pEqn.flux();
}
}
# include "immersedBoundaryContinuityErrs.H"
// Explicitly relax pressure for momentum corrector
p.relax();
// Momentum corrector
U -= fvc::grad(p)/AU;
U.correctBoundaryConditions();
}

87
applications/solvers/immersedBoundary/simpleIbFoam/simpleIbFoam.C
View file

@ -1,87 +0,0 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / 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
simpleIbFoam
Description
Steady-state solver for incompressible, turbulent flow
with immersed boundary support.
Author
Hrvoje Jasak, Wikki Ltd. All rights reserved
\*---------------------------------------------------------------------------*/
#include "fvCFD.H"
#include "singlePhaseTransportModel.H"
#include "RASModel.H"
#include "simpleControl.H"
#include "immersedBoundaryFvPatch.H"
#include "immersedBoundaryAdjustPhi.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
int main(int argc, char *argv[])
{
# include "setRootCase.H"
# include "createTime.H"
# include "createMesh.H"
simpleControl simple(mesh);
# include "createIbMasks.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;
}
// ************************************************************************* //