ICODE-ADC/CFDEMcoupling-PUBLIC
4
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

release on 2015-02-25_18-52-32

Browse Source
This commit is contained in:
Christoph Goniva, DCS Computing GmbH
2015-02-25 18:52:32 +01:00
parent f231cd4928
commit 027a03fbcb
208 changed files with 5511 additions and 1458 deletions
Download Patch File Download Diff File Expand all files Collapse all files

26
README
View File

@ -25,16 +25,17 @@ License
Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
Description
This code is designed to realize coupled CFD-DEM simulations using LIGGGHTS
and OpenFOAM. Note: this code is not part of OpenFOAM (see DISCLAIMER).
This code provides models and solvers to realize coupled CFD-DEM simulations
using LIGGGHTS and OpenFOAM.
Note: this code is not part of OpenFOAM (see DISCLAIMER).
\*---------------------------------------------------------------------------*/
CFDEM coupling provides an open source parallel coupled CFD-DEM framework
combining the strengths of LIGGGHTS DEM code and the Open Source
CFD package OpenFOAM(R)(*). The CFDEMcoupling toolbox allows to expand
CFDEM(R) coupling provides an open source parallel coupled CFD-DEM framework
combining the strengths of LIGGGHTS(R) DEM code and the Open Source
CFD package OpenFOAM(R)(*). The CFDEM(R)coupling toolbox allows to expand
standard CFD solvers of OpenFOAM(R)(*) to include a coupling to the DEM
code LIGGGHTS. In this toolbox the particle representation within the
code LIGGGHTS(R). In this toolbox the particle representation within the
CFD solver is organized by "cloud" classes. Key functionalities are organised
in sub-models (e.g. force models, data exchange models, etc.) which can easily
be selected and combined by dictionary settings.
@ -54,7 +55,7 @@ The file structure:
- "src" directory including the source files of the coupling toolbox and models
- "applications" directory including the solver files for coupled CFD-DEM simulations
- "doc" directory including the documentation of CFDEMcoupling
- "doc" directory including the documentation of CFDEM(R)coupling
- "tutorials" directory including basic tutorial cases showing the functionality
@ -64,18 +65,17 @@ Details on installation are given on the "www.cfdem.com"
The functionality of this CFD-DEM framwork is described via "tutorial cases" showing
how to use different solvers and models.
CFDEMcoupling stands for Computational Fluid Dynamics (CFD) -
CFDEM(R)coupling stands for Computational Fluid Dynamics (CFD) -
Discrete Element Method (DEM) coupling.
CFDEMcoupling is an open-source code, distributed freely under the terms of the
CFDEM(R)coupling is an open-source code, distributed freely under the terms of the
GNU Public License (GPL).
Core development of CFDEMcoupling is done by
Core development of CFDEM(R)coupling is done by
Christoph Goniva and Christoph Kloss, both at DCS Computing GmbH, 2012
\*---------------------------------------------------------------------------*/
(*) "OpenFOAM(R)"_of is a registered trade mark of the ESI Group.
This offering is not affiliated, approved or endorsed by ESI Group,
the producer of the OpenFOAM® software and owner of the OpenFOAM® trade mark.
(*) "OpenFOAM(R)"_of is a registered trade mark of OpenCFD Limited, a wholly owned subsidiary of the ESI Group.
This offering is not approved or endorsed by OpenCFD Limited, the producer of the OpenFOAM software and owner of the OPENFOAM® and OpenCFD® trade marks.
\*---------------------------------------------------------------------------*/

7
applications/.gitignore vendored Normal file
View File

@ -0,0 +1,7 @@
*.o
*.d
*.a
*.dep
log_*
log.*
*~

8
applications/solvers/cfdemSolverIB/Make/options
View File

@ -1,3 +1,5 @@
include $(CFDEM_ADD_LIBS_DIR)/additionalLibs
EXE_INC = \
-I$(LIB_SRC)/turbulenceModels/incompressible/turbulenceModel \
-I$(LIB_SRC)/transportModels \
@ -11,7 +13,6 @@ EXE_INC = \
-I$(LIB_SRC)/dynamicMesh/dynamicFvMesh/lnInclude \
-I$(LIB_SRC)/dynamicMesh/dynamicMesh/lnInclude \
-I$(LIB_SRC)/fvOptions/lnInclude
EXE_LIBS = \
-L$(CFDEM_LIB_DIR)\
-lincompressibleRASModels \
@ -21,4 +22,7 @@ EXE_LIBS = \
-ldynamicFvMesh \
-ldynamicMesh \
-lfvOptions \
-l$(CFDEM_LIB_NAME)
-l$(CFDEM_LIB_NAME) \
$(CFDEM_ADD_LIB_PATHS) \
$(CFDEM_ADD_LIBS)

6
applications/solvers/cfdemSolverIB/cfdemSolverIB.C
View File

@ -72,7 +72,7 @@ int main(int argc, char *argv[])
#include "initContinuityErrs.H"
#if defined(version22)
#include "createFvOptions.H"
#include "createFvOptions.H"
#endif
// create cfdemCloud
@ -128,10 +128,12 @@ int main(int argc, char *argv[])
for (int corr=0; corr<nCorr; corr++)
{
volScalarField rUA = 1.0/UEqn.A();
surfaceScalarField rUAf(fvc::interpolate(rUA));
U = rUA*UEqn.H();
#ifdef version23
phi = (fvc::interpolate(U) & mesh.Sf()); // there is a new version in 23x
phi = (fvc::interpolate(U) & mesh.Sf())
+ rUAf*fvc::ddtCorr(U, phi);
#else
phi = (fvc::interpolate(U) & mesh.Sf())
+ fvc::ddtPhiCorr(rUA, U, phi);

6
applications/solvers/cfdemSolverPiso/Make/options
View File

@ -1,3 +1,5 @@
include $(CFDEM_ADD_LIBS_DIR)/additionalLibs
EXE_INC = \
-I$(LIB_SRC)/turbulenceModels/incompressible/turbulenceModel \
-I$(LIB_SRC)/transportModels \
@ -12,4 +14,6 @@ EXE_LIBS = \
-lincompressibleLESModels \
-lincompressibleTransportModels \
-lfiniteVolume \
-l$(CFDEM_LIB_NAME)
-l$(CFDEM_LIB_NAME) \
$(CFDEM_ADD_LIB_PATHS) \
$(CFDEM_ADD_LIBS)

52
applications/solvers/cfdemSolverPiso/cfdemSolverPiso.C
View File

@ -42,6 +42,7 @@ Description
#include "implicitCouple.H"
#include "clockModel.H"
#include "smoothingModel.H"
#include "forceModel.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@ -59,11 +60,10 @@ int main(int argc, char *argv[])
#include "checkModelType.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Info<< "\nStarting time loop\n" << endl;
while (runTime.loop())
{
Info<< "\nStarting time loop\n" << endl;
particleCloud.clockM().start(1,"Global");
particleCloud.clockM().start(1,"Global");
Info<< "Time = " << runTime.timeName() << nl << endl;
@ -72,13 +72,24 @@ int main(int argc, char *argv[])
// do particle stuff
particleCloud.clockM().start(2,"Coupling");
particleCloud.evolve(voidfraction,Us,U);
bool hasEvolved = particleCloud.evolve(voidfraction,Us,U);
if(hasEvolved)
{
particleCloud.smoothingM().smoothen(particleCloud.forceM(0).impParticleForces());
}
Info << "update Ksl.internalField()" << endl;
Ksl = particleCloud.momCoupleM(0).impMomSource();
particleCloud.smoothingM().smoothen(Ksl);
Ksl.correctBoundaryConditions();
//Force Checks
vector fTotal(0,0,0);
vector fImpTotal = sum(mesh.V()*Ksl.internalField()*(Us.internalField()-U.internalField()));
reduce(fImpTotal, sumOp<vector>());
Info << "TotalForceExp: " << fTotal << endl;
Info << "TotalForceImp: " << fImpTotal << endl;
#include "solverDebugInfo.H"
particleCloud.clockM().stop("Coupling");
@ -91,23 +102,19 @@ int main(int argc, char *argv[])
// Momentum predictor
fvVectorMatrix UEqn
(
fvm::ddt(voidfraction,U) + fvm::Sp(fvc::ddt(voidfraction),U)
+ fvm::div(phi,U) + fvm::Sp(fvc::div(phi),U)
fvm::ddt(voidfraction,U) - fvm::Sp(fvc::ddt(voidfraction),U)
+ fvm::div(phi,U) - fvm::Sp(fvc::div(phi),U)
// + turbulence->divDevReff(U)
+ particleCloud.divVoidfractionTau(U, voidfraction)
==
- fvm::Sp(Ksl/rho,U)
);
if (modelType=="B")
UEqn == - fvc::grad(p) + Ksl/rho*Us;
else
UEqn == - voidfraction*fvc::grad(p) + Ksl/rho*Us;
UEqn.relax();
if (momentumPredictor)
solve(UEqn);
if (momentumPredictor && (modelType=="B" || modelType=="Bfull"))
solve(UEqn == - fvc::grad(p) + Ksl/rho*Us);
else if (momentumPredictor)
solve(UEqn == - voidfraction*fvc::grad(p) + Ksl/rho*Us);
// --- PISO loop
@ -120,11 +127,17 @@ int main(int argc, char *argv[])
surfaceScalarField rUAf("(1|A(U))", fvc::interpolate(rUA));
volScalarField rUAvoidfraction("(voidfraction2|A(U))",rUA*voidfraction);
surfaceScalarField rUAfvoidfraction("(voidfraction2|A(U)F)", fvc::interpolate(rUAvoidfraction));
U = rUA*UEqn.H();
phi = (fvc::interpolate(U*voidfraction) & mesh.Sf() );
//+ fvc::ddtPhiCorr(rUAvoidfraction, U, phi);
#ifdef version23
phi = ( fvc::interpolate(U*voidfraction) & mesh.Sf() )
+ rUAfvoidfraction*fvc::ddtCorr(U, phi);
#else
phi = ( fvc::interpolate(U*voidfraction) & mesh.Sf() )
+ fvc::ddtPhiCorr(rUAvoidfraction, U, phi);
#endif
surfaceScalarField phiS(fvc::interpolate(Us*voidfraction) & mesh.Sf());
surfaceScalarField phiGes = phi + rUAf*(fvc::interpolate(Ksl/rho) * phiS);
@ -157,13 +170,14 @@ int main(int argc, char *argv[])
if (nonOrth == nNonOrthCorr)
{
phiGes -= pEqn.flux();
phi = phiGes;
}
} // end non-orthogonal corrector loop
#include "continuityErrorPhiPU.H"
if (modelType=="B")
if (modelType=="B" || modelType=="Bfull")
U -= rUA*fvc::grad(p) - Ksl/rho*Us*rUA;
else
U -= voidfraction*rUA*fvc::grad(p) - Ksl/rho*Us*rUA;
@ -191,7 +205,7 @@ int main(int argc, char *argv[])
}
Info<< "End\n" << endl;
return 0;
}

6
applications/solvers/cfdemSolverPisoScalar/Make/options
View File

@ -1,3 +1,5 @@
include $(CFDEM_ADD_LIBS_DIR)/additionalLibs
EXE_INC = \
-I$(LIB_SRC)/turbulenceModels/incompressible/turbulenceModel \
-I$(LIB_SRC)/transportModels \
@ -12,4 +14,6 @@ EXE_LIBS = \
-lincompressibleLESModels \
-lincompressibleTransportModels \
-lfiniteVolume \
-l$(CFDEM_LIB_NAME)
-l$(CFDEM_LIB_NAME) \
$(CFDEM_ADD_LIB_PATHS) \
$(CFDEM_ADD_LIBS)

211
applications/solvers/cfdemSolverPisoScalar/cfdemSolverPisoScalar.C
View File

@ -30,7 +30,7 @@ Application
Description
Transient solver for incompressible flow.
Turbulence modelling is generic, i.e. laminar, RAS or LES may be selected.
The code is an evolution of the solver pisoFoam in OpenFOAM(R) 1.6,
The code is an evolution of the solver pisoFoam in OpenFOAM(R) 1.6,
where additional functionality for CFD-DEM coupling is added.
\*---------------------------------------------------------------------------*/
@ -40,30 +40,26 @@ Description
#include "cfdemCloud.H"
#include "implicitCouple.H"
#include "forceModel.H"
#include "smoothingModel.H"
#include "forceModel.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
int main(int argc, char *argv[])
{
#include "setRootCase.H"
#include "createTime.H"
#include "createMesh.H"
#include "createFields.H"
#include "initContinuityErrs.H"
// create cfdemCloud
#include "readGravitationalAcceleration.H"
cfdemCloud particleCloud(mesh);
#include "checkModelType.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Info<< "\nStarting time loop\n" << endl;
while (runTime.loop())
{
Info<< "Time = " << runTime.timeName() << nl << endl;
@ -72,11 +68,15 @@ int main(int argc, char *argv[])
#include "CourantNo.H"
// do particle stuff
Info << "- evolve()" << endl;
particleCloud.evolve(voidfraction,Us,U);
bool hasEvolved = particleCloud.evolve(voidfraction,Us,U);
Ksl.internalField() = particleCloud.momCoupleM(0).impMomSource();
particleCloud.smoothingM().smoothen(Ksl);
if(hasEvolved)
{
particleCloud.smoothingM().smoothen(particleCloud.forceM(0).impParticleForces());
}
Info << "update Ksl.internalField()" << endl;
Ksl = particleCloud.momCoupleM(0).impMomSource();
Ksl.correctBoundaryConditions();
@ -87,107 +87,116 @@ int main(int argc, char *argv[])
Tsource.correctBoundaryConditions();
// solve scalar transport equation
phi = fvc::interpolate(U*voidfraction) & mesh.Sf();
solve
fvScalarMatrix TEqn
(
fvm::ddt(voidfraction,T)
+ fvm::div(phi, T)
- fvm::laplacian(DT*voidfraction, T)
==
Tsource
fvm::ddt(voidfraction,T) - fvm::Sp(fvc::ddt(voidfraction),T)
+ fvm::div(phi, T) - fvm::Sp(fvc::div(phi),T)
- fvm::laplacian(DT*voidfraction, T)
==
Tsource
);
TEqn.relax();
TEqn.solve();
// Pressure-velocity PISO corrector
if(particleCloud.solveFlow())
{
// Momentum predictor
fvVectorMatrix UEqn
(
fvm::ddt(voidfraction,U)
+ fvm::div(phi, U)
+ turbulence->divDevReff(U)
==
- fvm::Sp(Ksl/rho,U)
);
UEqn.relax();
if (momentumPredictor)
// Pressure-velocity PISO corrector
{
//solve UEqn
if (modelType=="B")
// Momentum predictor
fvVectorMatrix UEqn
(
fvm::ddt(voidfraction,U) - fvm::Sp(fvc::ddt(voidfraction),U)
+ fvm::div(phi,U) - fvm::Sp(fvc::div(phi),U)
// + turbulence->divDevReff(U)
+ particleCloud.divVoidfractionTau(U, voidfraction)
==
- fvm::Sp(Ksl/rho,U)
);
UEqn.relax();
if (momentumPredictor && (modelType=="B" || modelType=="Bfull"))
solve(UEqn == - fvc::grad(p) + Ksl/rho*Us);
else
else if (momentumPredictor)
solve(UEqn == - voidfraction*fvc::grad(p) + Ksl/rho*Us);
// --- PISO loop
//for (int corr=0; corr<nCorr; corr++)
int nCorrSoph = nCorr + 5 * pow((1-particleCloud.dataExchangeM().timeStepFraction()),1);
for (int corr=0; corr<nCorrSoph; corr++)
{
volScalarField rUA = 1.0/UEqn.A();
surfaceScalarField rUAf("(1|A(U))", fvc::interpolate(rUA));
volScalarField rUAvoidfraction("(voidfraction2|A(U))",rUA*voidfraction);
surfaceScalarField rUAfvoidfraction("(voidfraction2|A(U)F)", fvc::interpolate(rUAvoidfraction));
U = rUA*UEqn.H();
#ifdef version23
phi = ( fvc::interpolate(U*voidfraction) & mesh.Sf() )
+ rUAfvoidfraction*fvc::ddtCorr(U, phi);
#else
phi = ( fvc::interpolate(U*voidfraction) & mesh.Sf() )
+ fvc::ddtPhiCorr(rUAvoidfraction, U, phi);
#endif
surfaceScalarField phiS(fvc::interpolate(Us*voidfraction) & mesh.Sf());
surfaceScalarField phiGes = phi + rUAf*(fvc::interpolate(Ksl/rho) * phiS);
if (modelType=="A")
rUAvoidfraction = volScalarField("(voidfraction2|A(U))",rUA*voidfraction*voidfraction);
// Non-orthogonal pressure corrector loop
for (int nonOrth=0; nonOrth<=nNonOrthCorr; nonOrth++)
{
// Pressure corrector
fvScalarMatrix pEqn
(
fvm::laplacian(rUAvoidfraction, p) == fvc::div(phiGes) + particleCloud.ddtVoidfraction()
);
pEqn.setReference(pRefCell, pRefValue);
if
(
corr == nCorr-1
&& nonOrth == nNonOrthCorr
)
{
pEqn.solve(mesh.solver("pFinal"));
}
else
{
pEqn.solve();
}
if (nonOrth == nNonOrthCorr)
{
phiGes -= pEqn.flux();
phi = phiGes;
}
} // end non-orthogonal corrector loop
#include "continuityErrorPhiPU.H"
if (modelType=="B" || modelType=="Bfull")
U -= rUA*fvc::grad(p) - Ksl/rho*Us*rUA;
else
U -= voidfraction*rUA*fvc::grad(p) - Ksl/rho*Us*rUA;
U.correctBoundaryConditions();
} // end piso loop
}
// --- PISO loop
//for (int corr=0; corr<nCorr; corr++)
int nCorrSoph = nCorr + 5 * pow((1-particleCloud.dataExchangeM().timeStepFraction()),1);
for (int corr=0; corr<nCorrSoph; corr++)
{
volScalarField rUA = 1.0/UEqn.A();
surfaceScalarField rUAf("(1|A(U))", fvc::interpolate(rUA));
volScalarField rUAvoidfraction("(voidfraction2|A(U))",rUA*voidfraction);
U = rUA*UEqn.H();
#ifdef version23
phi = ( fvc::interpolate(U*voidfraction) & mesh.Sf() );
#else
phi = ( fvc::interpolate(U*voidfraction) & mesh.Sf() )
+ fvc::ddtPhiCorr(rUAvoidfraction, U, phi);
#endif
surfaceScalarField phiS(fvc::interpolate(Us*voidfraction) & mesh.Sf());
surfaceScalarField phiGes = phi + rUAf*(fvc::interpolate(Ksl/rho) * phiS);
if (modelType=="A")
rUAvoidfraction = volScalarField("(voidfraction2|A(U))",rUA*voidfraction*voidfraction);
// Non-orthogonal pressure corrector loop
for (int nonOrth=0; nonOrth<=nNonOrthCorr; nonOrth++)
{
// Pressure corrector
fvScalarMatrix pEqn
(
fvm::laplacian(rUAvoidfraction, p) == fvc::div(phiGes) + particleCloud.ddtVoidfraction()
);
pEqn.setReference(pRefCell, pRefValue);
if
(
corr == nCorr-1
&& nonOrth == nNonOrthCorr
)
{
pEqn.solve(mesh.solver("pFinal"));
}
else
{
pEqn.solve();
}
if (nonOrth == nNonOrthCorr)
{
phiGes -= pEqn.flux();
}
} // end non-orthogonal corrector loop
#include "continuityErrorPhiPU.H"
if (modelType=="B")
U -= rUA*fvc::grad(p) - Ksl/rho*Us*rUA;
else
U -= voidfraction*rUA*fvc::grad(p) - Ksl/rho*Us*rUA;
U.correctBoundaryConditions();
} // end piso loop
turbulence->correct();
}// end solveFlow
else
{
Info << "skipping flow solution." << endl;
}
turbulence->correct();
runTime.write();
Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"

4
applications/utilities/cfdemPostproc/Make/options
View File

@ -1,3 +1,5 @@
include $(CFDEM_ADD_LIBS_DIR)/additionalLibs
EXE_INC = \
-I$(LIB_SRC)/turbulenceModels/incompressible/turbulenceModel \
-I$(LIB_SRC)/transportModels \
@ -14,4 +16,6 @@ EXE_LIBS = \
-lincompressibleTransportModels \
-lfiniteVolume \
-l$(CFDEM_LIB_NAME) \
$(CFDEM_ADD_LIB_PATHS) \
$(CFDEM_ADD_LIBS)

5
applications/utilities/cfdemPostproc/cfdemPostproc.C
View File

@ -71,6 +71,7 @@ int main(int argc, char *argv[])
double **voidfractions_;
double **particleWeights_;
double **particleVolumes_;
double **particleV_;
double **cellIDs_;
particleCloud.dataExchangeM().allocateArray(positions_,0.,3);
@ -80,6 +81,7 @@ int main(int argc, char *argv[])
particleCloud.dataExchangeM().allocateArray(voidfractions_,0.,1);
particleCloud.dataExchangeM().allocateArray(particleWeights_,0.,1);
particleCloud.dataExchangeM().allocateArray(particleVolumes_,0.,1);
particleCloud.dataExchangeM().allocateArray(particleV_,0.,1);
particleCloud.get_cellIDs(cellIDs_); // get ref to cellIDs
//particleCloud.dataExchangeM().allocateArray(cellIDs_,0.,1);
@ -105,7 +107,7 @@ int main(int argc, char *argv[])
particleCloud.locateM().findCell(NULL,positions_,cellIDs_,particleCloud.numberOfParticles());
particleCloud.setPos(positions_);
particleCloud.voidFractionM().setvoidFraction(NULL,voidfractions_,particleWeights_,particleVolumes_);
particleCloud.voidFractionM().setvoidFraction(NULL,voidfractions_,particleWeights_,particleVolumes_,particleV_);
voidfraction.internalField() = particleCloud.voidFractionM().voidFractionInterp();
voidfraction.correctBoundaryConditions();
@ -135,6 +137,7 @@ int main(int argc, char *argv[])
particleCloud.dataExchangeM().destroy(voidfractions_,1);
particleCloud.dataExchangeM().destroy(particleWeights_,1);
particleCloud.dataExchangeM().destroy(particleVolumes_,1);
particleCloud.dataExchangeM().destroy(particleV_,1);
//particleCloud.dataExchangeM().destroy(cellIDs_); // destroyed in cloud
Info<< "End\n" << endl;

4
applications/utilities/cfdemPostproc/createFields.H
View File

@ -26,7 +26,7 @@
IOobject::NO_WRITE
),
mesh,
dimensionedVector("0", dimensionSet(0, 1, -1, 0, 0), vector::zero)
dimensionedVector("0", dimensionSet(0, 1, -1, 0, 0), Foam::vector::zero)
);
//========================
@ -62,7 +62,7 @@
IOobject::AUTO_WRITE
),
mesh,
dimensionedVector("0", dimensionSet(0, 1, -1, 0, 0), vector::zero)
dimensionedVector("0", dimensionSet(0, 1, -1, 0, 0), Foam::vector::zero)
);
//========================

7
doc/.gitignore vendored Normal file
View File

@ -0,0 +1,7 @@
*.o
*.d
*.a
*.dep
log_*
log.*
*~

10
doc/CFDEMcoupling_Manual.html
View File

@ -164,7 +164,7 @@ In order to get the latest code version, please use the git repository at http:/
</P>
<PRE>modelType
</PRE>
<P>"modelType" refers to the formulation of the equations to be solved. Choose "A" or "B", according to Zhou et al. (2010): "Discrete particle simulation of particle-fluid flow: model formulations and their applicability", JFM. "A" requires the use of the force models gradPForce and viscForce, whereas "B" requires the force model "Archimedes".
<P>"modelType" refers to the formulation of the equations to be solved. Choose "A", "B" or "Bfull", according to Zhou et al. (2010): "Discrete particle simulation of particle-fluid flow: model formulations and their applicability", JFM. "A" requires the use of the force models gradPForce and viscForce, whereas "B" requires the force model "Archimedes". "Bfull" refers to model type I.
</P>
<PRE>couplingInterval
</PRE>
@ -217,9 +217,11 @@ listing below of styles within certain commands.
<TR ALIGN="center"><TD ><A HREF = "forceModel_DiFeliceDrag.html">forceModel_DiFeliceDrag</A></TD><TD ><A HREF = "forceModel_GidaspowDrag.html">forceModel_GidaspowDrag</A></TD></TR>
<TR ALIGN="center"><TD ><A HREF = "forceModel_KochHillDrag.html">forceModel_KochHillDrag</A></TD><TD ><A HREF = "forceModel_LaEuScalarTemp.html">forceModel_LaEuScalarTemp</A></TD></TR>
<TR ALIGN="center"><TD ><A HREF = "forceModel_MeiLift.html">forceModel_MeiLift</A></TD><TD ><A HREF = "forceModel_SchillerNaumannDrag.html">forceModel_SchillerNaumannDrag</A></TD></TR>
<TR ALIGN="center"><TD ><A HREF = "forceModel_ShirgaonkarIB.html">forceModel_ShirgaonkarIB</A></TD><TD ><A HREF = "forceModel_gradPForce.html">forceModel_gradPForce</A></TD></TR>
<TR ALIGN="center"><TD ><A HREF = "forceModel_noDrag.html">forceModel_noDrag</A></TD><TD ><A HREF = "forceModel_particleCellVolume.html">forceModel_particleCellVolume</A></TD></TR>
<TR ALIGN="center"><TD ><A HREF = "forceModel_virtualMassForce.html">forceModel_virtualMassForce</A></TD><TD ><A HREF = "forceModel_viscForce.html">forceModel_viscForce</A></TD></TR>
<TR ALIGN="center"><TD ><A HREF = "forceModel_ShirgaonkarIB.html">forceModel_ShirgaonkarIB</A></TD><TD ><A HREF = "forceModel_fieldStore.html">forceModel_fieldStore</A></TD></TR>
<TR ALIGN="center"><TD ><A HREF = "forceModel_gradPForce.html">forceModel_gradPForce</A></TD><TD ><A HREF = "forceModel_noDrag.html">forceModel_noDrag</A></TD></TR>
<TR ALIGN="center"><TD ><A HREF = "forceModel_particleCellVolume.html">forceModel_particleCellVolume</A></TD><TD ><A HREF = "forceModel_virtualMassForce.html">forceModel_virtualMassForce</A></TD></TR>
<TR ALIGN="center"><TD ><A HREF = "forceModel_viscForce.html">forceModel_viscForce</A></TD><TD ><A HREF = "forceSubModel.html">forceSubModel</A></TD></TR>
<TR ALIGN="center"><TD ><A HREF = "forceSubModel_ImEx.html">forceSubModel_ImEx</A></TD><TD ><A HREF = "forceSubModel_ImExCorr.html">forceSubModel_ImExCorr</A></TD></TR>
<TR ALIGN="center"><TD ><A HREF = "liggghtsCommandModel.html">liggghtsCommandModel</A></TD><TD ><A HREF = "liggghtsCommandModel_execute.html">liggghtsCommandModel_execute</A></TD></TR>
<TR ALIGN="center"><TD ><A HREF = "liggghtsCommandModel_readLiggghtsData.html">liggghtsCommandModel_readLiggghtsData</A></TD><TD ><A HREF = "liggghtsCommandModel_runLiggghts.html">liggghtsCommandModel_runLiggghts</A></TD></TR>
<TR ALIGN="center"><TD ><A HREF = "liggghtsCommandModel_writeLiggghts.html">liggghtsCommandModel_writeLiggghts</A></TD><TD ><A HREF = "locateModel.html">locateModel</A></TD></TR>

BIN
doc/CFDEMcoupling_Manual.pdf
View File

Binary file not shown.

10
doc/CFDEMcoupling_Manual.txt
View File

@ -133,7 +133,7 @@ Reasonable example settings for the "couplingProperties" dictionary are given in
modelType :pre
"modelType" refers to the formulation of the equations to be solved. Choose "A" or "B", according to Zhou et al. (2010): "Discrete particle simulation of particle-fluid flow: model formulations and their applicability", JFM. "A" requires the use of the force models gradPForce and viscForce, whereas "B" requires the force model "Archimedes".
"modelType" refers to the formulation of the equations to be solved. Choose "A", "B" or "Bfull", according to Zhou et al. (2010): "Discrete particle simulation of particle-fluid flow: model formulations and their applicability", JFM. "A" requires the use of the force models gradPForce and viscForce, whereas "B" requires the force model "Archimedes". "Bfull" refers to model type I.
couplingInterval :pre
@ -221,6 +221,10 @@ Reasonable example settings for the "liggghtsCommands" dictionary are given in t
@ -261,11 +265,15 @@ listing below of styles within certain commands.
"forceModel_MeiLift"_forceModel_MeiLift.html,
"forceModel_SchillerNaumannDrag"_forceModel_SchillerNaumannDrag.html,
"forceModel_ShirgaonkarIB"_forceModel_ShirgaonkarIB.html,
"forceModel_fieldStore"_forceModel_fieldStore.html,
"forceModel_gradPForce"_forceModel_gradPForce.html,
"forceModel_noDrag"_forceModel_noDrag.html,
"forceModel_particleCellVolume"_forceModel_particleCellVolume.html,
"forceModel_virtualMassForce"_forceModel_virtualMassForce.html,
"forceModel_viscForce"_forceModel_viscForce.html,
"forceSubModel"_forceSubModel.html,
"forceSubModel_ImEx"_forceSubModel_ImEx.html,
"forceSubModel_ImExCorr"_forceSubModel_ImExCorr.html,
"liggghtsCommandModel"_liggghtsCommandModel.html,
"liggghtsCommandModel_execute"_liggghtsCommandModel_execute.html,
"liggghtsCommandModel_readLiggghtsData"_liggghtsCommandModel_readLiggghtsData.html,

4
doc/cfdemSolverPisoScalar.html
View File

@ -11,12 +11,14 @@
</H3>
<P><B>Description:</B>
</P>
<P>"cfdemSolverPisoScalar" is a coupled CFD-DEM solver using CFDEMcoupling, an open source parallel coupled CFD-DEM framework. Based on pisoFoam(R)(*), a finite volume based solver for turbulent Navier-Stokes equations applying PISO algorithm, "cfdemSolverPisoScalar" has additional functionality for a coupling to the DEM code "LIGGGHTS" as well as a scalar transport equation. The volume averaged Navier-Stokes Equations are solved accounting for momentum exchange and volume displacement of discrete particles whose trajectories are calculated in the DEM code LIGGGHTS. The scalar transport equation is coupled to scalar properties of the particle phase, thus convective heat transfer in a fluid granular system can be modeled with "cfdemSolverPisoScalar".
<P>"cfdemSolverPisoScalar" is a coupled CFD-DEM solver using CFDEMcoupling, an open source parallel coupled CFD-DEM framework. Based on pisoFoam(R)(*), a finite volume based solver for turbulent Navier-Stokes equations applying PISO algorithm, "cfdemSolverPisoScalar" has additional functionality for a coupling to the DEM code "LIGGGHTS" as well as a scalar transport equation. The volume averaged Navier-Stokes Equations are solved accounting for momentum exchange and volume displacement of discrete particles, whose trajectories are calculated in the DEM code LIGGGHTS. The scalar transport equation is coupled to scalar properties of the particle phase, thus convective heat transfer in a fluid granular system can be modeled with "cfdemSolverPisoScalar".
</P>
<P>see:
</P>
<P>GONIVA, C., KLOSS, C., HAGER,A. and PIRKER, S. (2010): "An Open Source CFD-DEM Perspective", Proc. of OpenFOAM Workshop, Göteborg, June 22.-24.
</P>
<P>The heat transfer equation is implemented according to Nield & Bejan (2013), Convection in Porous Media, DOI 10.1007/978-1-4614-5541-7_2, Springer
</P>
<HR>
<P>(*) This offering is not approved or endorsed by OpenCFD Limited, the producer of the OpenFOAM software and owner of the OPENFOAM® and OpenCFD® trade marks. OPENFOAM® is a registered trade mark of OpenCFD Limited, a wholly owned subsidiary of the ESI Group.

5
doc/cfdemSolverPisoScalar.txt
View File

@ -9,12 +9,15 @@ cfdemSolverPisoScalar command :h3
[Description:]
"cfdemSolverPisoScalar" is a coupled CFD-DEM solver using CFDEMcoupling, an open source parallel coupled CFD-DEM framework. Based on pisoFoam(R)(*), a finite volume based solver for turbulent Navier-Stokes equations applying PISO algorithm, "cfdemSolverPisoScalar" has additional functionality for a coupling to the DEM code "LIGGGHTS" as well as a scalar transport equation. The volume averaged Navier-Stokes Equations are solved accounting for momentum exchange and volume displacement of discrete particles whose trajectories are calculated in the DEM code LIGGGHTS. The scalar transport equation is coupled to scalar properties of the particle phase, thus convective heat transfer in a fluid granular system can be modeled with "cfdemSolverPisoScalar".
"cfdemSolverPisoScalar" is a coupled CFD-DEM solver using CFDEMcoupling, an open source parallel coupled CFD-DEM framework. Based on pisoFoam(R)(*), a finite volume based solver for turbulent Navier-Stokes equations applying PISO algorithm, "cfdemSolverPisoScalar" has additional functionality for a coupling to the DEM code "LIGGGHTS" as well as a scalar transport equation. The volume averaged Navier-Stokes Equations are solved accounting for momentum exchange and volume displacement of discrete particles, whose trajectories are calculated in the DEM code LIGGGHTS. The scalar transport equation is coupled to scalar properties of the particle phase, thus convective heat transfer in a fluid granular system can be modeled with "cfdemSolverPisoScalar".
see:
GONIVA, C., KLOSS, C., HAGER,A. and PIRKER, S. (2010): "An Open Source CFD-DEM Perspective", Proc. of OpenFOAM Workshop, Göteborg, June 22.-24.
The heat transfer equation is implemented according to Nield & Bejan (2013), Convection in Porous Media, DOI 10.1007/978-1-4614-5541-7_2, Springer
:line
(*) This offering is not approved or endorsed by OpenCFD Limited, the producer of the OpenFOAM software and owner of the OPENFOAM® and OpenCFD® trade marks. OPENFOAM® is a registered trade mark of OpenCFD Limited, a wholly owned subsidiary of the ESI Group.

2
doc/clockModel_standardClock.html
View File

@ -23,7 +23,7 @@
</P>
<P>The "standardClock" model is a basic clockModel model which measures the run time between every ".start(int arrayPos,string name)" and ".stop(string name)" statement placed in the code. If a ".start(name)" is called more than once (e.g. in a loop) the accumulated times are calculated. After the simulation has finished, the data is stored in $caseDir/CFD/clockData/$startTime/*.txt .
Since the measurements are stored in an array, it is necessary to put a variable <I>arrayPos</I> (type integer) at the start command. Those do not need to be in ascending order and positions may be omitted. The standard size of this array is 30 and can be changed at the initialization of the standardClock class. If <I>arrayPos</I> is out of bounds, the array size will be doubled. The stop command does not need <I>arrayPos</I>, since the class remembers the positions. The string name is intended for easier evaluation afterwards an may be omitted like ".start(int arrayPos)" and ".stop()". The command ".stop(string name)" is a safety feature, because if the name is not equal to the started name, output will be produced for information.
After the case ran you may use the matPlot.py script located in $CFDEM_UT_DIR/vizClock/ to produce a graphical output of your measurements. The usage is like 'python < matPlot.py' and you have to be in the directory of the desired time step, where there is a file called "timeEvalFull.txt", which contains averaged and maximum data with respect to the number of processes.
After the case ran you may use the matPlot.py script located in $CFDEM_UT_DIR/vizClock/ to produce a graphical output of your measurements. The usage is like 'python < matPlot.py' and you have to be in the directory of the desired time step, where there is a file called "timeEvalFull.txt", which contains averaged and maximum data with respect to the number of processes. There is an alias called "vizClock" to run this python routine for visualizing the data.
</P>
<P><B>Restrictions:</B> none.
</P>

2
doc/clockModel_standardClock.txt
View File

@ -21,7 +21,7 @@ clockModel standardClock; :pre
The "standardClock" model is a basic clockModel model which measures the run time between every ".start(int arrayPos,string name)" and ".stop(string name)" statement placed in the code. If a ".start(name)" is called more than once (e.g. in a loop) the accumulated times are calculated. After the simulation has finished, the data is stored in $caseDir/CFD/clockData/$startTime/*.txt .
Since the measurements are stored in an array, it is necessary to put a variable {arrayPos} (type integer) at the start command. Those do not need to be in ascending order and positions may be omitted. The standard size of this array is 30 and can be changed at the initialization of the standardClock class. If {arrayPos} is out of bounds, the array size will be doubled. The stop command does not need {arrayPos}, since the class remembers the positions. The string name is intended for easier evaluation afterwards an may be omitted like ".start(int arrayPos)" and ".stop()". The command ".stop(string name)" is a safety feature, because if the name is not equal to the started name, output will be produced for information.
After the case ran you may use the matPlot.py script located in $CFDEM_UT_DIR/vizClock/ to produce a graphical output of your measurements. The usage is like 'python < matPlot.py' and you have to be in the directory of the desired time step, where there is a file called "timeEvalFull.txt", which contains averaged and maximum data with respect to the number of processes.
After the case ran you may use the matPlot.py script located in $CFDEM_UT_DIR/vizClock/ to produce a graphical output of your measurements. The usage is like 'python < matPlot.py' and you have to be in the directory of the desired time step, where there is a file called "timeEvalFull.txt", which contains averaged and maximum data with respect to the number of processes. There is an alias called "vizClock" to run this python routine for visualizing the data.
[Restrictions:] none.

2
doc/forceModel.html
View File

@ -33,7 +33,7 @@
</P>
<P><B>Description:</B>
</P>
<P>The force model performs the calculation of forces (e.g. fluid-particle interaction forces) acting on each DEM particle. All force models selected are executed sequentially and the forces on the particles are superposed.
<P>The force model performs the calculation of forces (e.g. fluid-particle interaction forces) acting on each DEM particle. All force models selected are executed sequentially and the forces on the particles are superposed. If the fluid density field is needed, by default a field named "rho" will be used. Via the forceSubModel an alternative field can be chosen.
</P>
<P><B>Restrictions:</B>
</P>

2
doc/forceModel.txt
View File

@ -31,7 +31,7 @@ Note: This examples list might not be complete - please look for other models (f
[Description:]
The force model performs the calculation of forces (e.g. fluid-particle interaction forces) acting on each DEM particle. All force models selected are executed sequentially and the forces on the particles are superposed.
The force model performs the calculation of forces (e.g. fluid-particle interaction forces) acting on each DEM particle. All force models selected are executed sequentially and the forces on the particles are superposed. If the fluid density field is needed, by default a field named "rho" will be used. Via the forceSubModel an alternative field can be chosen.
[Restrictions:]

6
doc/forceModel_Archimedes.html
View File

@ -19,13 +19,10 @@
);
ArchimedesProps
{
densityFieldName "density";
gravityFieldName "gravity";
};
</PRE>
<UL><LI><I>density</I> = name of the finite volume density field
<LI><I>gravity</I> = name of the finite volume gravity field
<UL><LI><I>gravity</I> = name of the finite volume gravity field
</UL>
@ -37,7 +34,6 @@ ArchimedesProps
);
ArchimedesProps
{
densityFieldName "rho";
gravityFieldName "g";
}
</PRE>

5
doc/forceModel_Archimedes.txt
View File

@ -17,12 +17,10 @@ forceModels
);
ArchimedesProps
\{
densityFieldName "density";
gravityFieldName "gravity";
\}; :pre
{density} = name of the finite volume density field :ulb,l
{gravity} = name of the finite volume gravity field :l
{gravity} = name of the finite volume gravity field :ulb,l
:ule
[Examples:]
@ -33,7 +31,6 @@ forceModels
);
ArchimedesProps
\{
densityFieldName "rho";
gravityFieldName "g";
\} :pre

6
doc/forceModel_ArchimedesIB.html
View File

@ -19,14 +19,11 @@
);
ArchimedesIBProps
{
densityFieldName "density";
gravityFieldName "gravity";
voidfractionFieldName "voidfraction";
};
</PRE>
<UL><LI><I>density</I> = name of the finite volume density field
<LI><I>gravity</I> = name of the finite volume gravity field
<UL><LI><I>gravity</I> = name of the finite volume gravity field
<LI><I>voidfraction</I> = name of the finite volume voidfraction field
@ -40,7 +37,6 @@ ArchimedesIBProps
);
ArchimedesIBProps
{
densityFieldName "rho";
gravityFieldName "g";
voidfractionFieldName "voidfractionNext";
}

5
doc/forceModel_ArchimedesIB.txt
View File

@ -17,13 +17,11 @@ forceModels
);
ArchimedesIBProps
\{
densityFieldName "density";
gravityFieldName "gravity";
voidfractionFieldName "voidfraction";
\}; :pre
{density} = name of the finite volume density field :ulb,l
{gravity} = name of the finite volume gravity field :l
{gravity} = name of the finite volume gravity field :ulb,l
{voidfraction} = name of the finite volume voidfraction field :l
:ule
@ -35,7 +33,6 @@ forceModels
);
ArchimedesIBProps
\{
densityFieldName "rho";
gravityFieldName "g";
voidfractionFieldName "voidfractionNext";
\} :pre

10
doc/forceModel_DiFeliceDrag.html
View File

@ -20,15 +20,12 @@
DiFeliceDragProps
{
velFieldName "U";
densityFieldName "density";
interpolation;
interpolation switch1;
};
</PRE>
<UL><LI><I>U</I> = name of the finite volume fluid velocity field
<LI><I>density</I> = name of the finite volume gravity field
<LI><I>interpolation</I> = flag to use interpolated voidfraction and velocity values (normally off)
<LI><I>switch1</I> = flag to use interpolated voidfraction and velocity values (normally off)
</UL>
@ -41,8 +38,7 @@ DiFeliceDragProps
DiFeliceDragProps
{
velFieldName "U";
densityFieldName "rho";
interpolation;
interpolation true;
}
</PRE>
<P><B>Description:</B>

9
doc/forceModel_DiFeliceDrag.txt
View File

@ -18,13 +18,11 @@ forceModels
DiFeliceDragProps
\{
velFieldName "U";
densityFieldName "density";
interpolation;
interpolation switch1;
\}; :pre
{U} = name of the finite volume fluid velocity field :ulb,l
{density} = name of the finite volume gravity field :l
{interpolation} = flag to use interpolated voidfraction and velocity values (normally off) :l
{switch1} = flag to use interpolated voidfraction and velocity values (normally off) :l
:ule
[Examples:]
@ -36,8 +34,7 @@ forceModels
DiFeliceDragProps
\{
velFieldName "U";
densityFieldName "rho";
interpolation;
interpolation true;
\} :pre
[Description:]

16
doc/forceModel_GidaspowDrag.html
View File

@ -20,24 +20,24 @@
GidaspowDragProps
{
velFieldName "U";
densityFieldName "density";
voidfractionFieldName "voidfraction";
granVelFieldName "Us";
phi "scalar";
interpolation;
implDEM;
interpolation switch1;
implForceDEM switch2;
};
</PRE>
<UL><LI><I>U</I> = name of the finite volume fluid velocity field
<LI><I>density</I> = name of the finite volume gravity field
<LI><I>voidfraction</I> = name of the finite volume voidfraction field
<LI><I>Us</I> = name of the finite volume cell averaged particle velocity field
<LI><I>phi</I> = drag correction factor (in doubt 1)
<LI><I>interpolation</I> = (optional, normally off) flag to use interpolated voidfraction and fluid velocity values
<LI><I>switch1</I> = (optional, normally off) flag to use interpolated voidfraction and fluid velocity values
<I>implDEM</I> = (optional, normally off) flag to use implicit formulation of drag on DEM side:l
<I>switch2</I> = (optional, normally off) flag to use implicit formulation of drag on DEM side:l
</UL>
<P><B>Examples:</B>
@ -49,8 +49,8 @@ GidaspowDragProps
GidaspowDragProps
{
velFieldName "U";
densityFieldName "rho";
voidfractionFieldName "voidfraction";
granVelFieldName "Us";
}
</PRE>
<P><B>Description:</B>

14
doc/forceModel_GidaspowDrag.txt
View File

@ -18,19 +18,19 @@ forceModels
GidaspowDragProps
\{
velFieldName "U";
densityFieldName "density";
voidfractionFieldName "voidfraction";
granVelFieldName "Us";
phi "scalar";
interpolation;
implDEM;
interpolation switch1;
implForceDEM switch2;
\}; :pre
{U} = name of the finite volume fluid velocity field :ulb,l
{density} = name of the finite volume gravity field :l
{voidfraction} = name of the finite volume voidfraction field :l
{Us} = name of the finite volume cell averaged particle velocity field :l
{phi} = drag correction factor (in doubt 1) :l
{interpolation} = (optional, normally off) flag to use interpolated voidfraction and fluid velocity values :l
{implDEM} = (optional, normally off) flag to use implicit formulation of drag on DEM side:l
{switch1} = (optional, normally off) flag to use interpolated voidfraction and fluid velocity values :l
{switch2} = (optional, normally off) flag to use implicit formulation of drag on DEM side:l
:ule
[Examples:]
@ -42,8 +42,8 @@ forceModels
GidaspowDragProps
\{
velFieldName "U";
densityFieldName "rho";
voidfractionFieldName "voidfraction";
granVelFieldName "Us";
\} :pre
[Description:]

12
doc/forceModel_KochHillDrag.html
View File

@ -20,21 +20,18 @@
KochHillDragProps
{
velFieldName "U";
densityFieldName "density";
voidfractionFieldName "voidfraction";
interpolation;
implDEM;
interpolation "bool1";
implForceDEM "bool2";
};
</PRE>
<UL><LI><I>U</I> = name of the finite volume fluid velocity field
<LI><I>density</I> = name of the finite volume gravity field
<LI><I>voidfraction</I> = name of the finite volume voidfraction field
<LI><I>interpolation</I> = (optional, normally off) flag to use interpolated voidfraction and fluid velocity values
<LI><I>bool1</I> = (optional, normally off) flag to use interpolated voidfraction and fluid velocity values
<I>implDEM</I> = (optional, normally off) flag to use implicit formulation of drag on DEM side:l
<I>bool2</I> = (optional, normally off) flag to use implicit formulation of drag on DEM side:l
</UL>
<P><B>Examples:</B>
@ -46,7 +43,6 @@ KochHillDragProps
KochHillDragProps
{
velFieldName "U";
densityFieldName "rho";
voidfractionFieldName "voidfraction";
}
</PRE>

11
doc/forceModel_KochHillDrag.txt
View File

@ -18,17 +18,15 @@ forceModels
KochHillDragProps
\{
velFieldName "U";
densityFieldName "density";
voidfractionFieldName "voidfraction";
interpolation;
implDEM;
interpolation "bool1";
implForceDEM "bool2";
\}; :pre
{U} = name of the finite volume fluid velocity field :ulb,l
{density} = name of the finite volume gravity field :l
{voidfraction} = name of the finite volume voidfraction field :l
{interpolation} = (optional, normally off) flag to use interpolated voidfraction and fluid velocity values :l
{implDEM} = (optional, normally off) flag to use implicit formulation of drag on DEM side:l
{bool1} = (optional, normally off) flag to use interpolated voidfraction and fluid velocity values :l
{bool2} = (optional, normally off) flag to use implicit formulation of drag on DEM side:l
:ule
[Examples:]
@ -40,7 +38,6 @@ forceModels
KochHillDragProps
\{
velFieldName "U";
densityFieldName "rho";
voidfractionFieldName "voidfraction";
\} :pre

18
doc/forceModel_LaEuScalarTemp.html
View File

@ -21,23 +21,19 @@ LaEuScalarTempProps
{
velFieldName "U";
tempFieldName "T";
tempSourceFieldName "Tsource";
voidfractionFieldName "voidfraction";
partTempName "Temp";
partHeatFluxName "convectiveHeatFlux";
lambda value;
Cp value1;
densityFieldName "density";
interpolation;
verbose;
interpolation "switch1";
verbose "switch2";
};
</PRE>
<UL><LI><I>U</I> = name of the finite volume fluid velocity field
<LI><I>T</I> = name of the finite volume scalar temperature field
<LI><I>Tsource</I> = name of the finite volume scalar temperature source field
<LI><I>voidfraction</I> = name of the finite volume voidfraction field
<LI><I>Temp</I> = name of the DEM data representing the particles temperature
@ -48,11 +44,9 @@ LaEuScalarTempProps
<LI><I>value1</I> = fluid specific heat capacity [W*s/(kg*K)]
<LI><I>density</I> = name of the finite volume fluid density field
<LI><I>switch1</I> = (optional, normally off) flag to use interpolated voidfraction and fluid velocity values
<LI><I>interpolation</I> = (optional, normally off) flag to use interpolated voidfraction and fluid velocity values
<LI><I>verbose</I> = (normally off) for verbose run
<LI><I>switch2</I> = (normally off) for verbose run
</UL>
@ -66,13 +60,11 @@ LaEuScalarTempProps
{
velFieldName "U";
tempFieldName "T";
tempSourceFieldName "Tsource";
voidfractionFieldName "voidfraction";
partTempName "Temp";
partHeatFluxName "convectiveHeatFlux";
lambda 0.0256;
Cp 1007;
densityFieldName "rho";
}
</PRE>
<P><B>Description:</B>
@ -81,7 +73,7 @@ LaEuScalarTempProps
</P>
<P><B>Restrictions:</B>
</P>
<P>Goes only with cfdemSolverScalar.
<P>Goes only with cfdemSolverScalar. The force model has to be the second (!!!) model in the forces list.
</P>
<P><B>Related commands:</B>
</P>

16
doc/forceModel_LaEuScalarTemp.txt
View File

@ -19,28 +19,24 @@ LaEuScalarTempProps
\{
velFieldName "U";
tempFieldName "T";
tempSourceFieldName "Tsource";
voidfractionFieldName "voidfraction";
partTempName "Temp";
partHeatFluxName "convectiveHeatFlux";
lambda value;
Cp value1;
densityFieldName "density";
interpolation;
verbose;
interpolation "switch1";
verbose "switch2";
\}; :pre
{U} = name of the finite volume fluid velocity field :ulb,l
{T} = name of the finite volume scalar temperature field :l
{Tsource} = name of the finite volume scalar temperature source field :l
{voidfraction} = name of the finite volume voidfraction field :l
{Temp} = name of the DEM data representing the particles temperature :l
{convectiveHeatFlux} = name of the DEM data representing the particle-fluid convective heat flux :l
{value} = fluid thermal conductivity \[W/(m*K)\] :l
{value1} = fluid specific heat capacity \[W*s/(kg*K)\] :l
{density} = name of the finite volume fluid density field :l
{interpolation} = (optional, normally off) flag to use interpolated voidfraction and fluid velocity values :l
{verbose} = (normally off) for verbose run :l
{switch1} = (optional, normally off) flag to use interpolated voidfraction and fluid velocity values :l
{switch2} = (normally off) for verbose run :l
:ule
[Examples:]
@ -53,13 +49,11 @@ LaEuScalarTempProps
\{
velFieldName "U";
tempFieldName "T";
tempSourceFieldName "Tsource";
voidfractionFieldName "voidfraction";
partTempName "Temp";
partHeatFluxName "convectiveHeatFlux";
lambda 0.0256;
Cp 1007;
densityFieldName "rho";
\} :pre
[Description:]
@ -68,7 +62,7 @@ This "forceModel" does not influence the particles or the fluid flow! Using the
[Restrictions:]
Goes only with cfdemSolverScalar.
Goes only with cfdemSolverScalar. The force model has to be the second (!!!) model in the forces list.
[Related commands:]

16
doc/forceModel_MeiLift.html
View File

@ -20,21 +20,18 @@
MeiLiftProps
{
velFieldName "U";
densityFieldName "density";
useSecondOrderTerms;
interpolation;
verbose;
interpolation "switch1";
verbose "switch2";
};
</PRE>
<UL><LI><I>U</I> = name of the finite volume fluid velocity field
<LI><I>density</I> = name of the finite volume fluid density field
<LI><I>useSecondOrderTerms</I> = switch to activate second order terms in the lift force model
<LI><I>interpolation</I> = switch to activate tri-linear interpolation of the flow quantities at the particle position
<LI><I>switch1</I> = switch to activate tri-linear interpolation of the flow quantities at the particle position
<LI><I>verbose</I> = switch to activate the report of per-particle quantities to the screen
<LI><I>switch2</I> = switch to activate the report of per-particle quantities to the screen
</UL>
@ -47,10 +44,9 @@ MeiLiftProps
MeiLiftProps
{
velFieldName "U";
densityFieldName "rho";
useSecondOrderTerms;
interpolation;
verbose;
interpolation true;
verbose true;
}
</PRE>
<P><B>Description:</B>

15
doc/forceModel_MeiLift.txt
View File

@ -18,17 +18,15 @@ forceModels
MeiLiftProps
\{
velFieldName "U";
densityFieldName "density";
useSecondOrderTerms;
interpolation;
verbose;
interpolation "switch1";
verbose "switch2";
\}; :pre
{U} = name of the finite volume fluid velocity field :ulb,l
{density} = name of the finite volume fluid density field :l
{useSecondOrderTerms} = switch to activate second order terms in the lift force model :l
{interpolation} = switch to activate tri-linear interpolation of the flow quantities at the particle position :l
{verbose} = switch to activate the report of per-particle quantities to the screen :l
{switch1} = switch to activate tri-linear interpolation of the flow quantities at the particle position :l
{switch2} = switch to activate the report of per-particle quantities to the screen :l
:ule
[Examples:]
@ -40,10 +38,9 @@ forceModels
MeiLiftProps
\{
velFieldName "U";
densityFieldName "rho";
useSecondOrderTerms;
interpolation;
verbose;
interpolation true;
verbose true;
\} :pre
[Description:]

4
doc/forceModel_SchillerNaumannDrag.html
View File

@ -20,13 +20,10 @@
SchillerNaumannDragProps
{
velFieldName "U";
densityFieldName "density";
};
</PRE>
<UL><LI><I>U</I> = name of the finite volume fluid velocity field
<LI><I>density</I> = name of the finite volume gravity field
</UL>
<P><B>Examples:</B>
@ -38,7 +35,6 @@ SchillerNaumannDragProps
SchillerNaumannDragProps
{
velFieldName "U";
densityFieldName "rho";
}
</PRE>
<P><B>Description:</B>

3
doc/forceModel_SchillerNaumannDrag.txt
View File

@ -18,11 +18,9 @@ forceModels
SchillerNaumannDragProps
\{
velFieldName "U";
densityFieldName "density";
\}; :pre
{U} = name of the finite volume fluid velocity field :ulb,l
{density} = name of the finite volume gravity field :l
:ule
[Examples:]
@ -34,7 +32,6 @@ forceModels
SchillerNaumannDragProps
\{
velFieldName "U";
densityFieldName "rho";
\} :pre
[Description:]

4
doc/forceModel_ShirgaonkarIB.html
View File

@ -20,14 +20,11 @@
ShirgaonkarIBProps
{
velFieldName "U";
densityFieldName "density";
pressureFieldName "pressure";
};
</PRE>
<UL><LI><I>U</I> = name of the finite volume fluid velocity field
<LI><I>density</I> = name of the finite volume density field
<LI><I>pressure</I> = name of the finite volume pressure field
@ -41,7 +38,6 @@ ShirgaonkarIBProps
ShirgaonkarIBProps
{
velFieldName "U";
densityFieldName "rho";
pressureFieldName "p";
}
</PRE>

3
doc/forceModel_ShirgaonkarIB.txt
View File

@ -18,12 +18,10 @@ forceModels
ShirgaonkarIBProps
\{
velFieldName "U";
densityFieldName "density";
pressureFieldName "pressure";
\}; :pre
{U} = name of the finite volume fluid velocity field :ulb,l
{density} = name of the finite volume density field :l
{pressure} = name of the finite volume pressure field :l
:ule
@ -36,7 +34,6 @@ forceModels
ShirgaonkarIBProps
\{
velFieldName "U";
densityFieldName "rho";
pressureFieldName "p";
\} :pre

68
doc/forceModel_fieldStore.html Normal file
View File

@ -0,0 +1,68 @@
<HTML>
<CENTER><A HREF = "http://www.cfdem.com">CFDEMproject WWW Site</A> - <A HREF = "CFDEMcoupling_Manual.html#comm">CFDEM Commands</A>
</CENTER>
<HR>
<H3>forceModel_fieldStore command
</H3>
<P><B>Syntax:</B>
</P>
<P>Defined in couplingProperties dictionary.
</P>
<PRE>forceModels
(
fieldStore
);
fieldStoreProps
{
scalarFieldNames
(
"scalarField"
);
vectorFieldNames
(
"vectorField"
);
};
</PRE>
<UL><LI><I>scalarField</I> = names of the finite volume scalar fields to be stored
<LI><I>vectorField</I> = names of the finite volume vector fields to be stored
</UL>
<P><B>Examples:</B>
</P>
<PRE>forceModels
(
fieldStore
);
fieldStoreProps
{
scalarFieldNames
(
"voidfraction"
);
vectorFieldNames
(
"U"
);
}
</PRE>
<P><B>Description:</B>
</P>
<P>This "forceModel" does not influence the particles or the flow - it is a tool to store a scalar/vector field! This is especially useful if you use a boundary condition which cannot interpreted correctly in your postporcessor (e.g. paraview).
</P>
<P><B>Restrictions:</B>
</P>
<P>none.
</P>
<P><B>Related commands:</B>
</P>
<P><A HREF = "forceModel.html">forceModel</A>
</P>
</HTML>

63
doc/forceModel_fieldStore.txt Normal file
View File

@ -0,0 +1,63 @@
"CFDEMproject WWW Site"_lws - "CFDEM Commands"_lc :c
:link(lws,http://www.cfdem.com)
:link(lc,CFDEMcoupling_Manual.html#comm)
:line
forceModel_fieldStore command :h3
[Syntax:]
Defined in couplingProperties dictionary.
forceModels
(
fieldStore
);
fieldStoreProps
\{
scalarFieldNames
(
"scalarField"
);
vectorFieldNames
(
"vectorField"
);
\}; :pre
{scalarField} = names of the finite volume scalar fields to be stored :ulb,l
{vectorField} = names of the finite volume vector fields to be stored :l
:ule
[Examples:]
forceModels
(
fieldStore
);
fieldStoreProps
\{
scalarFieldNames
(
"voidfraction"
);
vectorFieldNames
(
"U"
);
\} :pre
[Description:]
This "forceModel" does not influence the particles or the flow - it is a tool to store a scalar/vector field! This is especially useful if you use a boundary condition which cannot interpreted correctly in your postporcessor (e.g. paraview).
[Restrictions:]
none.
[Related commands:]
"forceModel"_forceModel.html

10
doc/forceModel_gradPForce.html
View File

@ -20,18 +20,15 @@
gradPForceProps
{
pFieldName "pressure";
densityFieldName "density";
velocityFieldName "U";
interpolation;
interpolation switch1;
};
</PRE>
<UL><LI><I>pressure</I> = name of the finite volume fluid pressure field
<LI><I>density</I> = name of the finite volume gravity field
<LI><I>U</I> = name of the finite volume fluid velocity field
<LI><I>interpolation</I> = flag to use interpolated pressure values (normally off)
<LI><I>switch1</I> = flag to use interpolated pressure values (normally off)
</UL>
@ -44,9 +41,8 @@ gradPForceProps
gradPForceProps
{
pFieldName "p";
densityFieldName "rho";
velocityFieldName "U";
interpolation;
interpolation true;
}
</PRE>
<P><B>Description:</B>

9
doc/forceModel_gradPForce.txt
View File

@ -18,15 +18,13 @@ forceModels
gradPForceProps
\{
pFieldName "pressure";
densityFieldName "density";
velocityFieldName "U";
interpolation;
interpolation switch1;
\}; :pre
{pressure} = name of the finite volume fluid pressure field :ulb,l
{density} = name of the finite volume gravity field :l
{U} = name of the finite volume fluid velocity field :l
{interpolation} = flag to use interpolated pressure values (normally off) :l
{switch1} = flag to use interpolated pressure values (normally off) :l
:ule
[Examples:]
@ -38,9 +36,8 @@ forceModels
gradPForceProps
\{
pFieldName "p";
densityFieldName "rho";
velocityFieldName "U";
interpolation;
interpolation true;
\} :pre
[Description:]

4
doc/forceModel_virtualMassForce.html
View File

@ -20,13 +20,10 @@
virtualMassForceProps
{
velFieldName "U";
densityFieldName "density";
};
</PRE>
<UL><LI><I>U</I> = name of the finite volume fluid velocity field
<LI><I>density</I> = name of the finite volume fluid density field
</UL>
<P><B>Examples:</B>
@ -38,7 +35,6 @@ virtualMassForceProps
virtualMassForceProps
{
velFieldName "U";
densityFieldName "rho";
}
</PRE>
<P><B>Description:</B>

3
doc/forceModel_virtualMassForce.txt
View File

@ -18,11 +18,9 @@ forceModels
virtualMassForceProps
\{
velFieldName "U";
densityFieldName "density";
\}; :pre
{U} = name of the finite volume fluid velocity field :ulb,l
{density} = name of the finite volume fluid density field :l
:ule
[Examples:]
@ -34,7 +32,6 @@ forceModels
virtualMassForceProps
\{
velFieldName "U";
densityFieldName "rho";
\} :pre
[Description:]

8
doc/forceModel_viscForce.html
View File

@ -20,15 +20,12 @@
viscForceProps
{
velocityFieldName "U";
densityFieldName "density";
interpolation;
interpolation "switch";
};
</PRE>
<UL><LI><I>U</I> = name of the finite volume fluid velocity field
<LI><I>density</I> = name of the finite volume gravity field
<LI><I>interpolation</I> = flag to use interpolated stress values (normally off)
<LI><I>switch</I> = flag to use interpolated stress values (normally off)
</UL>
@ -41,7 +38,6 @@ viscForceProps
viscForceProps
{
velocityFieldName "U";
densityFieldName "density";
}
</PRE>
<P><B>Description:</B>

7
doc/forceModel_viscForce.txt
View File

@ -18,13 +18,11 @@ forceModels
viscForceProps
\{
velocityFieldName "U";
densityFieldName "density";
interpolation;
interpolation "switch";
\}; :pre
{U} = name of the finite volume fluid velocity field :ulb,l
{density} = name of the finite volume gravity field :l
{interpolation} = flag to use interpolated stress values (normally off) :l
{switch} = flag to use interpolated stress values (normally off) :l
:ule
[Examples:]
@ -36,7 +34,6 @@ forceModels
viscForceProps
\{
velocityFieldName "U";
densityFieldName "density";
\} :pre
[Description:]

49
doc/forceSubModel.html Normal file
View File

@ -0,0 +1,49 @@
<HTML>
<CENTER><A HREF = "http://www.cfdem.com">CFDEMproject WWW Site</A> - <A HREF = "CFDEMcoupling_Manual.html#comm">CFDEM Commands</A>
</CENTER>
<HR>
<H3>forceSubModel command
</H3>
<P><B>Syntax:</B>
</P>
<P>Defined in couplingProperties sub-dictionary of the force model in use. If no force sub-model is applied ImEx is used as default. If the keyword "forceSubModels" is provided, a choice of sub model is demanded.
</P>
<PRE>forceSubModels
(
model_x
model_y
);
</PRE>
<UL><LI>model = name of force sub-model to be applied
</UL>
<P><B>Examples:</B>
</P>
<PRE>forceSubModels
(
ImEx
);
</PRE>
<P>Note: This examples list might not be complete - please look for other models (forceSubModel_XY) in this documentation.
</P>
<P><B>Description:</B>
</P>
<P>The force sub model is designed to hold the settings a force model can have. For now it handles the treatExplicit, treatDEM and implDEM option.
</P>
<P><B>Restrictions:</B>
</P>
<P>None.
</P>
<P><B>Related commands:</B>
</P>
<P><A HREF = "forceSubModel_ImEx.html">ImEx</A>
</P>
<P>Note: This examples list may be incomplete - please look for other models (forceSubModel_XY) in this documentation.
</P>
<P><B>Default:</B> none.
</P>
</HTML>

45
doc/forceSubModel.txt Normal file
View File

@ -0,0 +1,45 @@
"CFDEMproject WWW Site"_lws - "CFDEM Commands"_lc :c
:link(lws,http://www.cfdem.com)
:link(lc,CFDEMcoupling_Manual.html#comm)
:line
forceSubModel command :h3
[Syntax:]
Defined in couplingProperties sub-dictionary of the force model in use. If no force sub-model is applied ImEx is used as default. If the keyword "forceSubModels" is provided, a choice of sub model is demanded.
forceSubModels
(
model_x
model_y
); :pre
model = name of force sub-model to be applied :ul
[Examples:]
forceSubModels
(
ImEx
); :pre
Note: This examples list might not be complete - please look for other models (forceSubModel_XY) in this documentation.
[Description:]
The force sub model is designed to hold the settings a force model can have. For now it handles the treatExplicit, treatDEM and implDEM option.
[Restrictions:]
None.
[Related commands:]
"ImEx"_forceSubModel_ImEx.html
Note: This examples list may be incomplete - please look for other models (forceSubModel_XY) in this documentation.
[Default:] none.

45
doc/forceSubModel_ImEx.html Normal file
View File

@ -0,0 +1,45 @@
<HTML>
<CENTER><A HREF = "http://www.cfdem.com">CFDEMproject WWW Site</A> - <A HREF = "CFDEMcoupling_Manual.html#comm">CFDEM Commands</A>
</CENTER>
<HR>
<H3>forceSubModel_ImEx command
</H3>
<P><B>Syntax:</B>
</P>
<P>Defined in couplingProperties sub-dictionary of the force model in use.
</P>
<P>forceSubModels
(
ImEx;
);
</P>
<P>treatExplicit true; // optional for some force models.
treatDEM true; // optional for some force models.
implDEM true; // optional for some force models.
</P>
<P><B>Examples:</B>
</P>
<P>forceSubModels
(
ImEx;
);
treatExplicit true; // optional for some force models.
</P>
<P><B>Description:</B>
</P>
<P> If no force sub-model is applied ImEx is used as default. If the keyword "forceSubModels" is provided, a choice of sub model is demanded. Depending on the force model different keywords are read and can therefrore be set (see the log file). If the keyword is provided, its value is used.
</P>
<P><B>Restrictions:</B>
</P>
<P>none.
</P>
<P><B>Related commands:</B>
</P>
<P><A HREF = "forceSubModel.html">forceSubModel</A>
</P>
</HTML>

42
doc/forceSubModel_ImEx.txt Normal file
View File

@ -0,0 +1,42 @@
"CFDEMproject WWW Site"_lws - "CFDEM Commands"_lc :c
:link(lws,http://www.cfdem.com)
:link(lc,CFDEMcoupling_Manual.html#comm)
:line
forceSubModel_ImEx command :h3
[Syntax:]
Defined in couplingProperties sub-dictionary of the force model in use.
forceSubModels
(
ImEx;
);
treatExplicit true; // optional for some force models.
treatDEM true; // optional for some force models.
implDEM true; // optional for some force models.
[Examples:]
forceSubModels
(
ImEx;
);
treatExplicit true; // optional for some force models.
[Description:]
If no force sub-model is applied ImEx is used as default. If the keyword "forceSubModels" is provided, a choice of sub model is demanded. Depending on the force model different keywords are read and can therefrore be set (see the log file). If the keyword is provided, its value is used.
[Restrictions:]
none.
[Related commands:]
"forceSubModel"_forceSubModel.html

46
doc/forceSubModel_ImExCorr.html Normal file
View File

@ -0,0 +1,46 @@
<HTML>
<CENTER><A HREF = "http://www.cfdem.com">CFDEMproject WWW Site</A> - <A HREF = "CFDEMcoupling_Manual.html#comm">CFDEM Commands</A>
</CENTER>
<HR>
<H3>forceSubModel_ImExCorr command
</H3>
<P><B>Syntax:</B>
</P>
<P>Defined in couplingProperties sub-dictionary of the force model in use.
</P>
<P>forceSubModels
(
ImExCorr;
);
</P>
<P>treatExplicit true; // optional for some force models.
treatDEM true; // optional for some force models.
implDEM true; // optional for some force models.
explicitCorr true; // optional for some force models.
</P>
<P><B>Examples:</B>
</P>
<P>forceSubModels
(
ImExCorr;
);
treatExplicit true; // optional for some force models.
</P>
<P><B>Description:</B>
</P>
<P> Same as ImEx, but it additionally reads "explicitCorr" to correct the error steming from interpolation of Ufluid and averaging of Uparticles.
</P>
<P><B>Restrictions:</B>
</P>
<P>none.
</P>
<P><B>Related commands:</B>
</P>
<P><A HREF = "forceSubModel.html">forceSubModel</A>
</P>
</HTML>

43
doc/forceSubModel_ImExCorr.txt Normal file
View File

@ -0,0 +1,43 @@
"CFDEMproject WWW Site"_lws - "CFDEM Commands"_lc :c
:link(lws,http://www.cfdem.com)
:link(lc,CFDEMcoupling_Manual.html#comm)
:line
forceSubModel_ImExCorr command :h3
[Syntax:]
Defined in couplingProperties sub-dictionary of the force model in use.
forceSubModels
(
ImExCorr;
);
treatExplicit true; // optional for some force models.
treatDEM true; // optional for some force models.
implDEM true; // optional for some force models.
explicitCorr true; // optional for some force models.
[Examples:]
forceSubModels
(
ImExCorr;
);
treatExplicit true; // optional for some force models.
[Description:]
Same as ImEx, but it additionally reads "explicitCorr" to correct the error steming from interpolation of Ufluid and averaging of Uparticles.
[Restrictions:]
none.
[Related commands:]
"forceSubModel"_forceSubModel.html

BIN
doc/githubAccess_public.pdf
View File

Binary file not shown.

10
doc/liggghtsCommandModel_writeLiggghts.html
View File

@ -21,18 +21,20 @@
writeLiggghtsProps
{
writeLast switch1;
path "path";
writeName "name";
overwrite switch2;
verbose;
}
</PRE>
<UL><LI><I>switch1</I> = switch (choose on/off) to select if only last step is stored or every write step.
<UL><LI><I>switch1</I> = switch (choose on/off) to select if only last step is stored or every write step (default on).
<LI><I>name</I> = name of the restart file to be written in /$caseDir/DEM/ default default "liggghts.restartCFDEM"
<LI><I>path</I> = optionally an alternative path for saving the restart file can be defined (default "../DEM").
<LI><I>switch2</I> = switch (choose on/off) to select if only one restart file $name or many files $name_$timeStamp are written
<LI><I>name</I> = name of the restart file to be written in /$caseDir/DEM/ default (default "liggghts.restartCFDEM")
<LI><I>verbose</I> = (normally off) for verbose run
<LI><I>switch2</I> = switch (choose on/off) to select if only one restart file $name or many files $name_$timeStamp are written (default off):l
<I>verbose</I> = (default off) for verbose run
</UL>

10
doc/liggghtsCommandModel_writeLiggghts.txt
View File

@ -19,15 +19,17 @@ liggghtsCommandModels
writeLiggghtsProps
\{
writeLast switch1;
path "path";
writeName "name";
overwrite switch2;
verbose;
\} :pre
{switch1} = switch (choose on/off) to select if only last step is stored or every write step. :ulb,l
{name} = name of the restart file to be written in /$caseDir/DEM/ default default "liggghts.restartCFDEM" :l
{switch2} = switch (choose on/off) to select if only one restart file $name or many files $name_$timeStamp are written :l
{verbose} = (normally off) for verbose run :l
{switch1} = switch (choose on/off) to select if only last step is stored or every write step (default on). :ulb,l
{path} = optionally an alternative path for saving the restart file can be defined (default "../DEM"). :l
{name} = name of the restart file to be written in /$caseDir/DEM/ default (default "liggghts.restartCFDEM") :l
{switch2} = switch (choose on/off) to select if only one restart file $name or many files $name_$timeStamp are written (default off):l
{verbose} = (default off) for verbose run :l
:ule
[Examples:]

3
doc/momCoupleModel.html
View File

@ -33,7 +33,8 @@
</P>
<P>Note that the variable "imExSplitFactor" can be set in the couplingProperties in order to treat implicitly defined forces (in the implementation of the force model) as explicit ones. "imExSplitFactor 1.0;" is set by default, meaning that all implicit forces will be considered implicitly, whereas "imExSplitFactor 0.0;" would mean that implicitly defined forces will be treated in an explicit fashion.
</P>
<P><B>Description:</B>
<P>Note that the switch "treatVoidCellsAsExplicitForce true;" can be set in the couplingProperties in order to change the treatment of cells which are void of particles. This is only relevant if (i) smoothing is used, and (ii) implicit force coupling is performed. By default, the particle veloctiy field (Us) will be smoothed to obtain a meaningful reference quantity for the implicit force coupling. In case "treatVoidCellsAsExplicitForce true;" is set, however, Us will not be smoothed and implicit forces (after the smoothing has been performed) in cells void of particles be treated as explicit ones. This avoids the problem of defining Us in cells that are void of particles, but for which an implicit coupling force is obtained in the smoothing process.
<B>Description:</B>
</P>
<P>The momCoupleModel is the base class for momentum exchange between DEM and CFD simulation.
</P>

1
doc/momCoupleModel.txt
View File

@ -31,6 +31,7 @@ Forces can be coupled in an implicit way to the fluid solver (i.e., when solving
Note that the variable "imExSplitFactor" can be set in the couplingProperties in order to treat implicitly defined forces (in the implementation of the force model) as explicit ones. "imExSplitFactor 1.0;" is set by default, meaning that all implicit forces will be considered implicitly, whereas "imExSplitFactor 0.0;" would mean that implicitly defined forces will be treated in an explicit fashion.
Note that the switch "treatVoidCellsAsExplicitForce true;" can be set in the couplingProperties in order to change the treatment of cells which are void of particles. This is only relevant if (i) smoothing is used, and (ii) implicit force coupling is performed. By default, the particle veloctiy field (Us) will be smoothed to obtain a meaningful reference quantity for the implicit force coupling. In case "treatVoidCellsAsExplicitForce true;" is set, however, Us will not be smoothed and implicit forces (after the smoothing has been performed) in cells void of particles be treated as explicit ones. This avoids the problem of defining Us in cells that are void of particles, but for which an implicit coupling force is obtained in the smoothing process.
[Description:]
The momCoupleModel is the base class for momentum exchange between DEM and CFD simulation.

7
src/.gitignore vendored Normal file
View File

@ -0,0 +1,7 @@
*.o
*.d
*.a
*.dep
log_*
log.*
*~

38
src/lagrangian/cfdemParticle/cfdTools/checkModelType.H
View File

@ -4,8 +4,7 @@
word modelType = particleCloud.modelType();
//Warning << "model type not being checked" << endl;
if (modelType=="B"){
if (modelType=="Bfull"){
Info << "solving volume averaged Navier Stokes equations of type B\n"<< endl;
// check if Archimedes is used
@ -18,6 +17,41 @@
if(!found)
FatalError <<"Archimedes model not found!\n" << abort(FatalError);
// check if gradPForce is used
found=false;
forAll(particleCloud.forceModels(),i)
{
if(particleCloud.forceModels()[i]=="gradPForce")
found=true;
}
if(!found)
FatalError <<"gradPForce model not found!\n" << abort(FatalError);
// check if viscForce is used
found=false;
forAll(particleCloud.forceModels(),i)
{
if(particleCloud.forceModels()[i]=="viscForce")
found=true;
}
if(!found)
FatalError <<"viscForce model not found!\n" << abort(FatalError);
}else if(modelType=="B"){
Info << "solving volume averaged Navier Stokes equations of type B\n"<< endl;
// check if Archimedes is used
bool found=false;
forAll(particleCloud.forceModels(),i)
{
if(particleCloud.forceModels()[i]=="Archimedes")
found=true;
}
if(!found)
FatalError <<"Archimedes model not found!\n" << abort(FatalError);
// check if gradP and viscForce are used
found=false;
forAll(particleCloud.forceModels(),i)

53
src/lagrangian/cfdemParticle/cfdTools/compressibleContinuityErrsPU.H Normal file
View File

@ -0,0 +1,53 @@
/*---------------------------------------------------------------------------*\
CFDEMcoupling - Open Source CFD-DEM coupling
CFDEMcoupling is part of the CFDEMproject
www.cfdem.com
Christoph Goniva, christoph.goniva@cfdem.com
Copyright (C) 1991-2009 OpenCFD Ltd.
Copyright (C) 2009-2012 JKU, Linz
Copyright (C) 2012- DCS Computing GmbH,Linz
-------------------------------------------------------------------------------
License
This file is part of CFDEMcoupling.
CFDEMcoupling 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.
CFDEMcoupling 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 CFDEMcoupling. If not, see <http://www.gnu.org/licenses/>.
Global
continuityErrs
Description
Calculates and prints the continuity errors.
The code is an evolution of compressibleContinuityErrs.H in OpenFOAM(R) 2.3.x,
where additional functionality for CFD-DEM coupling is added.
\*---------------------------------------------------------------------------*/
{
dimensionedScalar totalMass = fvc::domainIntegrate(rho*voidfraction);
scalar sumLocalContErr =
(fvc::domainIntegrate(mag(rho - thermo.rho())*voidfraction)/totalMass).value();
scalar globalContErr =
(fvc::domainIntegrate((rho - thermo.rho())*voidfraction)/totalMass).value();
cumulativeContErr += globalContErr;
Info<< "time step continuity errors : sum local = " << sumLocalContErr
<< ", global = " << globalContErr
<< ", cumulative = " << cumulativeContErr
<< endl;
}
// ************************************************************************* //

80
src/lagrangian/cfdemParticle/cfdTools/global.C Normal file
View File

@ -0,0 +1,80 @@
/*---------------------------------------------------------------------------*\
CFDEMcoupling - Open Source CFD-DEM coupling
CFDEMcoupling is part of the CFDEMproject
www.cfdem.com
Christoph Goniva, christoph.goniva@cfdem.com
Copyright 2009-2012 JKU Linz
Copyright 2012- DCS Computing GmbH, Linz
-------------------------------------------------------------------------------
License
This file is part of CFDEMcoupling.
CFDEMcoupling 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.
CFDEMcoupling 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 CFDEMcoupling; if not, write to the Free Software Foundation,
Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
Description
This code is designed to realize coupled CFD-DEM simulations using LIGGGHTS
and OpenFOAM(R). Note: this code is not part of OpenFOAM(R) (see DISCLAIMER).
\*---------------------------------------------------------------------------*/
#include "error.H"
#include "global.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
defineTypeNameAndDebug(global, 0);
defineRunTimeSelectionTable(global, dictionary);
// * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
void global::info()
{
Info << "\nYou are currently using:" << endl;
Info << "OF version: " << FOAMversion << endl;
Info << "OF build: " << FOAMbuild << endl;
Info << "CFDEM build: " << CFDEMversion << "\n" << endl;
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
// Construct from components
global::global
(
const dictionary& dict,
cfdemCloud& sm
)
:
dict_(dict),
particleCloud_(sm),
CFDEMversion(GITVERSION)
{}
// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
global::~global()
{}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam
// ************************************************************************* //

128
src/lagrangian/cfdemParticle/cfdTools/global.H Normal file
View File

@ -0,0 +1,128 @@
/*---------------------------------------------------------------------------*\
CFDEMcoupling - Open Source CFD-DEM coupling
CFDEMcoupling is part of the CFDEMproject
www.cfdem.com
Christoph Goniva, christoph.goniva@cfdem.com
Copyright 2009-2012 JKU Linz
Copyright 2012- DCS Computing GmbH, Linz
-------------------------------------------------------------------------------
License
This file is part of CFDEMcoupling.
CFDEMcoupling 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.
CFDEMcoupling 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 CFDEMcoupling; if not, write to the Free Software Foundation,
Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
Description
This code is designed to realize coupled CFD-DEM simulations using LIGGGHTS
and OpenFOAM(R). Note: this code is not part of OpenFOAM(R) (see DISCLAIMER).
Class
global
SourceFiles
global.Cver
\*---------------------------------------------------------------------------*/
#ifndef global_H
#define global_H
#include "fvCFD.H"
#include "cfdemCloud.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
/*---------------------------------------------------------------------------*\
Class global Declaration
\*---------------------------------------------------------------------------*/
class global
{
protected:
// Protected data
const dictionary& dict_;
cfdemCloud& particleCloud_;
const char* const CFDEMversion;
// Protected member functions
public:
//- Runtime type information
TypeName("global");
// Declare runtime constructor selection table
declareRunTimeSelectionTable
(
autoPtr,
global,
dictionary,
(
const dictionary& dict,
cfdemCloud& sm
),
(dict,sm)
);
// Constructors
//- Construct from components
global
(
const dictionary& dict,
cfdemCloud& sm
);
// Destructor
virtual ~global();
// Selector
static autoPtr<global> New
(
const dictionary& dict,
cfdemCloud& sm
);
// Member Function
void info();
// Access
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //

265
src/lagrangian/cfdemParticle/cfdTools/mathExtra.H Normal file
View File

@ -0,0 +1,265 @@
/* ----------------------------------------------------------------------
CFDEMcoupling - Open Source CFD-DEM coupling
CFDEMcoupling is part of the CFDEMproject
www.cfdem.com
Christoph Goniva, christoph.goniva@cfdem.com
Copyright 2009-2012 JKU Linz
Copyright 2012- DCS Computing GmbH, Linz
-------------------------------------------------------------------------------
License
This file is part of CFDEMcoupling.
CFDEMcoupling 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.
CFDEMcoupling 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 CFDEMcoupling; if not, write to the Free Software Foundation,
Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
Description
This code is designed to realize coupled CFD-DEM simulations using LIGGGHTS
and OpenFOAM(R). Note: this code is not part of OpenFOAM(R) (see DISCLAIMER).
Copyright of this contribution:
Copyright 2014- TU Graz, IPPT
------------------------------------------------------------------------- */
#ifndef CFDEM_MATH_EXTRA_H
#define CFDEM_MATH_EXTRA_H
#include <vector>
//#include "math.h"
#include "stdio.h"
#include "string.h"
#include "error.h"
#include "ctype.h"
#define TOLERANCE_ORTHO 1e-10
namespace MathExtra
{
// inline void outerProduct(double *vec1, double *vec2, double **m);
// inline double spheroidGeometry(int index, double bi, double ai);
//--------------------------------------------------------------------
// Outer Product of two vectors
inline void outerProduct(double *vec1, double *vec2, double **m)
{
int i, j;
//debug output
for( i = 0; i < 3; ++i )
// printf("OUTER PRODUCT: Input: vec1 element %d = %g", i, vec1[i]);
for( i = 0; i < 3; ++i )
// printf("OUTER PRODUCT: Input: vec2 element %d=%g", i, vec2[i]);
//calculation
for( i = 0; i < 3; ++i )
for( j = 0; j < 3; ++j )
{
m[i][j] = vec1[i] * vec2[j];
printf("OUTER PRODUCT: Result: m[%d][%d]=%g", i, j, m[i][j]);
}
}
//--------------------------------------------------------------------
// Compute the major, minor axis and eccentricity parameters of a prolate spheroid
inline bool spheroidGeometry(double radius, double aspectRatio, //inputs
double& ai, double& bi, //outputs
double& ei, double& Le //outputs
) //
{
//INPUT
// radius ...volume-equivalent radius of the spheroid
// aspectRatio ...major/minor aspect ratio
//OUTPUT
// ai ...
// bi ...
// ei ...
// Le ...
if(radius<=0.0) //avoid troubles in case radius is 0 or negative
return false;
ai = radius * std::pow(aspectRatio*aspectRatio,0.33333333333333333333333);
bi = ai / aspectRatio;
ei = std::sqrt(
1.0
- 1.0 / (aspectRatio*aspectRatio)
);
Le = std::log(
(1.0+ei)
/(1.0-ei)
);
return true;
}
//--------------------------------------------------------------------
// Compute the major, minor axis and eccentricity parameters of a prolate spheroid
inline double Pi()
{
return 3.1415926535897932384626433832795;
}
//--------------------------------------------------------------------
// Compute the major, minor axis and eccentricity parameters of a prolate spheroid
inline bool spheroidGeometry2(double radius, double aspectRatio, //inputs
double& ai, double& bi, //outputs
double& XAe, double& YAe, //outputs
double& XCe, double& YCe, //outputs
double& YHe
) //
{
//INPUT
// radius ...volume-equivalent radius of the spheroid
// aspectRatio ...major/minor aspect ratio
//OUTPUT
// XAe ...Eccentricity dependet parameter
// YAe ...Eccentricity dependet parameter
// XCe ...Eccentricity dependet parameter
// XCe ...Eccentricity dependet parameter
// YCe ...Eccentricity dependet parameter
// YHe ...Eccentricity dependet parameter
double ei(0.0), Le(0.0);
bool result =
spheroidGeometry(radius, aspectRatio, //inputs
ai, bi, //outputs
ei, Le //outputs
);
if(!result)
return false;
XAe= 2.6666666666666666666666667
*ei*ei*ei
/(-2.0*ei+(1.0+ei*ei)*Le);
YAe= 5.333333333333333333333333333
*ei*ei*ei
/(2.0*ei+(3*ei*ei-1.0)*Le);
XCe= 1.333333333333333333333333333
*ei*ei*ei
*(1.0-ei*ei)
/(2.0*ei-(1.0-ei*ei)*Le);
YCe= 1.3333333333333333333333
*ei*ei*ei
*(2.0-ei*ei)
/(-2.0*ei+(1.0+ei*ei)*Le);
YHe= 1.3333333333333333333333
*ei*ei*ei*ei*ei
/(-2.0*ei+(1.0+ei*ei)*Le);
return true;
}
//--------------------------------------------------------------------
// zeroize a 3x3x3 tensor
inline void zeroize333(double tensor[3][3][3] )
{
for(int iX=0; iX<3; iX++)
for(int iY=0; iY<3; iY++)
for(int iZ=0; iZ<3; iZ++)
tensor[iX][iY][iZ] = 0.0;
}
//--------------------------------------------------------------------
// zeroize a 3x3 tensor
inline void zeroize33(double tensor[3][3] )
{
for(int iX=0; iX<3; iX++)
for(int iY=0; iY<3; iY++)
tensor[iX][iY] = 0.0;
}
//--------------------------------------------------------------------
// multiply a 3x3x3 tensor with a scalar
inline void multiply333(double scalar, double tensor[3][3][3] )
{
for(int iX=0; iX<3; iX++)
for(int iY=0; iY<3; iY++)
for(int iZ=0; iZ<3; iZ++)
tensor[iX][iY][iZ] *= scalar;
}
//--------------------------------------------------------------------
// Compute a dot and dyadic product of with a vector
inline void permutationTensor(double tensor[3][3][3] )
{
zeroize333(tensor);
tensor[0][1][2] = 1.0;
tensor[1][2][0] = 1.0;
tensor[2][0][1] = 1.0;
tensor[0][2][1] =-1.0;
tensor[2][1][0] =-1.0;
tensor[1][0][2] =-1.0;
}
//--------------------------------------------------------------------
// Compute a dot product of the permutation tensor and
// then a dyadic product of with a vector
inline bool permutationDotDyadic(
double vector[3],
double tensor[3][3][3]
)
{
//Generate permutation tensor
double permutationT[3][3][3];
permutationTensor(permutationT);
//Step 1: compute dot prodcut of permutation tensor
double permutationDotProd[3][3];
zeroize33(permutationDotProd);
for(int iX=0; iX<3; iX++)
for(int iY=0; iY<3; iY++)
for(int iZ=0; iZ<3; iZ++)
permutationDotProd[iX][iY] += permutationT[iX][iY][iZ]
* vector[iZ];
for(int iX=0; iX<3; iX++)
for(int iY=0; iY<3; iY++)
for(int iZ=0; iZ<3; iZ++)
tensor[iX][iY][iZ] = permutationDotProd[iX][iY]
* vector[iZ];
return true;
}
//--------------------------------------------------------------------
// Compute a dot and dyadic product of with a vector
inline bool doubleDotTensor333Tensor33(double tensor333[3][3][3],
double tensor33[3][3],
double result[3]
)
{
result[0]=0.0;result[1]=0.0;result[2]=0.0;
for(int iX=0; iX<3; iX++)
for(int iY=0; iY<3; iY++)
for(int iZ=0; iZ<3; iZ++)
result[iX] += tensor333[iX][iY][iZ] * tensor33[iY][iZ];
return true;
}
}; //end of namespace
#endif

84
src/lagrangian/cfdemParticle/cfdTools/newGlobal.C Normal file
View File

@ -0,0 +1,84 @@
/*---------------------------------------------------------------------------*\
CFDEMcoupling - Open Source CFD-DEM coupling
CFDEMcoupling is part of the CFDEMproject
www.cfdem.com
Christoph Goniva, christoph.goniva@cfdem.com
Copyright 2009-2012 JKU Linz
Copyright 2012- DCS Computing GmbH, Linz
-------------------------------------------------------------------------------
License
This file is part of CFDEMcoupling.
CFDEMcoupling 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.
CFDEMcoupling 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 CFDEMcoupling; if not, write to the Free Software Foundation,
Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
Description
This code is designed to realize coupled CFD-DEM simulations using LIGGGHTS
and OpenFOAM(R). Note: this code is not part of OpenFOAM(R) (see DISCLAIMER).
\*---------------------------------------------------------------------------*/
#include "error.H"
#include "global.H"
#include "dilute.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
autoPtr<global> global::New
(
const dictionary& dict,
cfdemCloud& sm
)
{
word globalType
(
dict.lookup("global")
);
Info<< "Selecting global "
<< globalType << endl;
dictionaryConstructorTable::iterator cstrIter =
dictionaryConstructorTablePtr_->find(globalType);
if (cstrIter == dictionaryConstructorTablePtr_->end())
{
FatalError
<< "global::New(const dictionary&, const spray&) : "
<< endl
<< " unknown globalType type "
<< globalType
<< ", constructor not in hash table" << endl << endl
<< " Valid global types are :"
<< endl;
Info<< dictionaryConstructorTablePtr_->toc()
<< abort(FatalError);
}
return autoPtr<global>(cstrIter()(dict,sm));
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam
// ************************************************************************* //

47
src/lagrangian/cfdemParticle/cfdTools/versionInfo.H
View File

@ -1,9 +1,44 @@
word CFDEMversion="cfdem-2.7.1";
word compatibleLIGGGHTSversion="3.0.2";
/*---------------------------------------------------------------------------*\
CFDEMcoupling - Open Source CFD-DEM coupling
CFDEMcoupling is part of the CFDEMproject
www.cfdem.com
Christoph Goniva, christoph.goniva@cfdem.com
Copyright 2009-2012 JKU Linz
Copyright 2012- DCS Computing GmbH, Linz
-------------------------------------------------------------------------------
License
This file is part of CFDEMcoupling.
CFDEMcoupling 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.
CFDEMcoupling 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 CFDEMcoupling; if not, write to the Free Software Foundation,
Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
Description
This code is designed to realize coupled CFD-DEM simulations using LIGGGHTS
and OpenFOAM(R). Note: this code is not part of OpenFOAM(R) (see DISCLAIMER).
\*---------------------------------------------------------------------------*/
#ifndef versionInfo_H
#define versionInfo_H
word CFDEMversion="cfdem-2.9.0";
word compatibleLIGGGHTSversion="3.1.0";
word OFversion="2.3.x-commit-4d6f4a3115ff76ec4154c580eb041bc95ba4ec09";
Info << "\nCFDEMcoupling version: " << CFDEMversion << "\n" << endl;
Info << "\n, compatible to LIGGGHTS version: " << compatibleLIGGGHTSversion << "\n" << endl;
Info << "\n, compatible to OF version: " << OFversion << "\n" << endl;
Info << "\nCFDEMcoupling version: " << CFDEMversion << endl;
Info << ", compatible to LIGGGHTS version: " << compatibleLIGGGHTSversion << endl;
Info << ", compatible to OF version and build: " << OFversion << endl;
#endif

77
src/lagrangian/cfdemParticle/cfdemCloud/cfdemCloud.C
View File

@ -31,6 +31,7 @@ Description
#include "fileName.H"
#include "cfdemCloud.H"
#include "global.H"
#include "forceModel.H"
#include "locateModel.H"
#include "momCoupleModel.H"
@ -80,6 +81,7 @@ Foam::cfdemCloud::cfdemCloud
positions_(NULL),
velocities_(NULL),
fluidVel_(NULL),
fAcc_(NULL),
impForces_(NULL),
expForces_(NULL),
DEMForces_(NULL),
@ -89,7 +91,9 @@ Foam::cfdemCloud::cfdemCloud
cellIDs_(NULL),
particleWeights_(NULL),
particleVolumes_(NULL),
particleV_(NULL),
numberOfParticles_(0),
d32_(-1),
numberOfParticlesChanged_(false),
arraysReallocated_(false),
forceModels_(couplingProperties_.lookup("forceModels")),
@ -99,7 +103,9 @@ Foam::cfdemCloud::cfdemCloud
cg_(1.),
cgOK_(true),
impDEMdrag_(false),
impDEMdragAcc_(false),
imExSplitFactor_(1.0),
treatVoidCellsAsExplicitForce_(false),
useDDTvoidfraction_(false),
ddtVoidfraction_
(
@ -117,7 +123,7 @@ Foam::cfdemCloud::cfdemCloud
turbulence_
(
#if defined(version21) || defined(version16ext)
#ifdef comp
#ifdef compre
mesh.lookupObject<compressible::turbulenceModel>
#else
mesh.lookupObject<incompressible::turbulenceModel>
@ -213,12 +219,16 @@ Foam::cfdemCloud::cfdemCloud
)
{
#include "versionInfo.H"
global buildInfo(couplingProperties_,*this);
buildInfo.info();
Info << "If BC are important, please provide volScalarFields -imp/expParticleForces-" << endl;
if (couplingProperties_.found("solveFlow"))
solveFlow_=Switch(couplingProperties_.lookup("solveFlow"));
if (couplingProperties_.found("imExSplitFactor"))
imExSplitFactor_ = readScalar(couplingProperties_.lookup("imExSplitFactor"));
if (couplingProperties_.found("treatVoidCellsAsExplicitForce"))
treatVoidCellsAsExplicitForce_ = readBool(couplingProperties_.lookup("treatVoidCellsAsExplicitForce"));
if (couplingProperties_.found("verbose")) verbose_=true;
if (couplingProperties_.found("ignore")) ignore_=true;
if (turbulenceModelType_=="LESProperties")
@ -276,6 +286,7 @@ Foam::cfdemCloud::~cfdemCloud()
dataExchangeM().destroy(positions_,3);
dataExchangeM().destroy(velocities_,3);
dataExchangeM().destroy(fluidVel_,3);
dataExchangeM().destroy(fAcc_,3);
dataExchangeM().destroy(impForces_,3);
dataExchangeM().destroy(expForces_,3);
dataExchangeM().destroy(DEMForces_,3);
@ -285,6 +296,7 @@ Foam::cfdemCloud::~cfdemCloud()
dataExchangeM().destroy(cellIDs_,1);
dataExchangeM().destroy(particleWeights_,1);
dataExchangeM().destroy(particleVolumes_,1);
dataExchangeM().destroy(particleV_,1);
}
// * * * * * * * * * * * * * * * private Member Functions * * * * * * * * * * * * * //
void Foam::cfdemCloud::getDEMdata()
@ -292,6 +304,9 @@ void Foam::cfdemCloud::getDEMdata()
dataExchangeM().getData("radius","scalar-atom",radii_);
dataExchangeM().getData("x","vector-atom",positions_);
dataExchangeM().getData("v","vector-atom",velocities_);
if(impDEMdragAcc_)
dataExchangeM().getData("dragAcc","vector-atom",fAcc_); // array is used twice - might be necessary to clean it first
}
void Foam::cfdemCloud::giveDEMdata()
@ -302,6 +317,7 @@ void Foam::cfdemCloud::giveDEMdata()
if(impDEMdrag_)
{
if(verbose_) Info << "sending Ksl and uf" << endl;
dataExchangeM().giveData("Ksl","scalar-atom",Cds_);
dataExchangeM().giveData("uf","vector-atom",fluidVel_);
}
@ -441,6 +457,16 @@ label Foam::cfdemCloud::liggghtsCommandModelIndex(word name)
return index;
}
std::vector< std::vector<double*> >* Foam::cfdemCloud::getVprobe()
{
return probeModel_->getVprobe();
}
std::vector< std::vector<double> >* Foam::cfdemCloud::getSprobe()
{
return probeModel_->getSprobe();
}
// * * * * * * * * * * * * * * * WRITE * * * * * * * * * * * * * //
// * * * write cfdemCloud internal data * * * //
@ -458,9 +484,10 @@ bool Foam::cfdemCloud::evolve
if(!ignore())
{
if (dataExchangeM().couple())
if (dataExchangeM().doCoupleNow())
{
Info << "\n Coupling..." << endl;
dataExchangeM().couple(0);
doCouple=true;
// reset vol Fields
@ -496,7 +523,7 @@ bool Foam::cfdemCloud::evolve
// set void fraction field
clockM().start(19,"setvoidFraction");
if(verbose_) Info << "- setvoidFraction()" << endl;
voidFractionM().setvoidFraction(NULL,voidfractions_,particleWeights_,particleVolumes_);
voidFractionM().setvoidFraction(NULL,voidfractions_,particleWeights_,particleVolumes_,particleV_);
if(verbose_) Info << "setvoidFraction done." << endl;
clockM().stop("setvoidFraction");
@ -508,7 +535,12 @@ bool Foam::cfdemCloud::evolve
//Smoothen "next" fields
smoothingM().dSmoothing();
smoothingM().smoothen(voidFractionM().voidFractionNext());
smoothingM().smoothenReferenceField(averagingM().UsNext());
//only smoothen if we use implicit force coupling in cells void of particles
//because we need unsmoothened Us field to detect cells for explicit
//force coupling
if(!treatVoidCellsAsExplicitForce())
smoothingM().smoothenReferenceField(averagingM().UsNext());
clockM().stop("setVectorAverage");
}
@ -560,6 +592,8 @@ bool Foam::cfdemCloud::evolve
clockM().start(23,"giveDEMdata");
giveDEMdata();
clockM().stop("giveDEMdata");
dataExchangeM().couple(1);
}//end dataExchangeM().couple()
@ -584,15 +618,40 @@ bool Foam::cfdemCloud::reAllocArrays() const
dataExchangeM().allocateArray(positions_,0.,3);
dataExchangeM().allocateArray(velocities_,0.,3);
dataExchangeM().allocateArray(fluidVel_,0.,3);
dataExchangeM().allocateArray(fAcc_,0.,3);
dataExchangeM().allocateArray(impForces_,0.,3);
dataExchangeM().allocateArray(expForces_,0.,3);
dataExchangeM().allocateArray(DEMForces_,0.,3);
dataExchangeM().allocateArray(Cds_,0.,1);
dataExchangeM().allocateArray(radii_,0.,1);
dataExchangeM().allocateArray(voidfractions_,1.,voidFractionM().maxCellsPerParticle());
dataExchangeM().allocateArray(cellIDs_,0.,voidFractionM().maxCellsPerParticle());
dataExchangeM().allocateArray(cellIDs_,-1.,voidFractionM().maxCellsPerParticle());
dataExchangeM().allocateArray(particleWeights_,0.,voidFractionM().maxCellsPerParticle());
dataExchangeM().allocateArray(particleVolumes_,0.,voidFractionM().maxCellsPerParticle());
dataExchangeM().allocateArray(particleV_,0.,1);
arraysReallocated_ = true;
return true;
}
return false;
}
bool Foam::cfdemCloud::reAllocArrays(int nP, bool forceRealloc) const
{
if( (numberOfParticlesChanged_ && !arraysReallocated_) || forceRealloc)
{
// get arrays of new length
dataExchangeM().allocateArray(positions_,0.,3,nP);
dataExchangeM().allocateArray(velocities_,0.,3,nP);
dataExchangeM().allocateArray(fluidVel_,0.,3,nP);
dataExchangeM().allocateArray(impForces_,0.,3,nP);
dataExchangeM().allocateArray(expForces_,0.,3,nP);
dataExchangeM().allocateArray(DEMForces_,0.,3,nP);
dataExchangeM().allocateArray(Cds_,0.,1,nP);
dataExchangeM().allocateArray(radii_,0.,1,nP);
dataExchangeM().allocateArray(voidfractions_,1.,voidFractionM().maxCellsPerParticle(),nP);
dataExchangeM().allocateArray(cellIDs_,0.,voidFractionM().maxCellsPerParticle(),nP);
dataExchangeM().allocateArray(particleWeights_,0.,voidFractionM().maxCellsPerParticle(),nP);
dataExchangeM().allocateArray(particleVolumes_,0.,voidFractionM().maxCellsPerParticle(),nP);
arraysReallocated_ = true;
return true;
}
@ -604,7 +663,7 @@ tmp<fvVectorMatrix> cfdemCloud::divVoidfractionTau(volVectorField& U,volScalarFi
return
(
- fvm::laplacian(voidfractionNuEff(voidfraction), U)
- fvc::div(voidfractionNuEff(voidfraction)*dev(fvc::grad(U)().T()))
- fvc::div(voidfractionNuEff(voidfraction)*dev2(fvc::grad(U)().T()))
);
}
@ -637,11 +696,11 @@ void cfdemCloud::calcDdtVoidfraction(volScalarField& voidfraction) const
tmp<volScalarField> cfdemCloud::voidfractionNuEff(volScalarField& voidfraction) const
{
if (modelType_=="B")
if (modelType_=="B" || modelType_=="Bfull")
{
return tmp<volScalarField>
(
#ifdef comp
#ifdef compre
new volScalarField("viscousTerm", (turbulence_.mut() + turbulence_.mu()))
#else
new volScalarField("viscousTerm", (turbulence_.nut() + turbulence_.nu()))
@ -652,7 +711,7 @@ tmp<volScalarField> cfdemCloud::voidfractionNuEff(volScalarField& voidfraction)
{
return tmp<volScalarField>
(
#ifdef comp
#ifdef compre
new volScalarField("viscousTerm", voidfraction*(turbulence_.mut() + turbulence_.mu()))
#else
new volScalarField("viscousTerm", voidfraction*(turbulence_.nut() + turbulence_.nu()))

41
src/lagrangian/cfdemParticle/cfdemCloud/cfdemCloud.H
View File

@ -44,6 +44,7 @@ SourceFiles
// choose version
#include "OFversion.H"
#include <vector>
#include "fvCFD.H"
#include "IFstream.H"
@ -102,6 +103,8 @@ protected:
mutable double **fluidVel_;
mutable double **fAcc_;
mutable double **impForces_;
mutable double **expForces_;
@ -120,8 +123,12 @@ protected:
mutable double **particleVolumes_;
mutable double **particleV_;
int numberOfParticles_;
scalar d32_;
bool numberOfParticlesChanged_;
mutable bool arraysReallocated_;
@ -140,14 +147,18 @@ protected:
bool impDEMdrag_;
bool impDEMdragAcc_;
mutable scalar imExSplitFactor_;
mutable bool treatVoidCellsAsExplicitForce_; //will treat the coupling force in cells with no Us data explicitly
bool useDDTvoidfraction_;
mutable volScalarField ddtVoidfraction_;
#if defined(version21) || defined(version16ext)
#ifdef comp
#ifdef compre
const compressible::turbulenceModel& turbulence_;
#else
const incompressible::turbulenceModel& turbulence_;
@ -200,6 +211,7 @@ public:
friend class dataExchangeModel;
friend class voidFractionModel;
friend class forceModel;
friend class forceSubModel;
// Constructors
@ -243,8 +255,12 @@ public:
inline const bool& impDEMdrag() const;
inline const bool& impDEMdragAcc() const;
inline const scalar& imExSplitFactor() const;
inline const bool& treatVoidCellsAsExplicitForce() const;
inline const bool& ignore() const;
inline const fvMesh& mesh() const;
@ -261,6 +277,8 @@ public:
inline double ** fluidVels() const;
inline double ** fAccs() const;
inline double ** impForces() const;
inline double ** expForces() const;
@ -289,6 +307,7 @@ public:
virtual inline double d(int);
inline scalar d32(bool recalc=true);
virtual inline double dMin() {return -1;}
virtual inline double dMax() {return -1;}
virtual inline int minType() {return -1;}
@ -298,9 +317,19 @@ public:
virtual inline double ** particleDensity() const {return NULL;};
virtual inline int ** particleTypes() const {return NULL;};
virtual label particleType(label index) const {return -1;};
//access to the particle's rotation and torque data
virtual inline double ** DEMTorques() const {return NULL;};
virtual inline double ** omegaArray() const {return NULL;};
virtual vector omega(int) const {return Foam::vector(0,0,0);};
//access to the particles' orientation information
virtual inline double ** exArray() const {return NULL;};
virtual vector ex(int) const {return Foam::vector(0,0,0);};
//Detector if SRF module is enable or not
virtual inline bool SRFOn(){return false;}
inline int numberOfParticles() const;
inline bool numberOfParticlesChanged() const;
@ -336,7 +365,7 @@ public:
inline autoPtr<liggghtsCommandModel>* liggghtsCommand() const;
#if defined(version21) || defined(version16ext)
#ifdef comp
#ifdef compre
inline const compressible::turbulenceModel& turbulence() const;
#else
inline const incompressible::turbulenceModel& turbulence() const;
@ -352,6 +381,9 @@ public:
virtual bool reAllocArrays() const;
virtual bool reAllocArrays(int nP, bool forceRealloc) const; //force number of particles during reallocation
// IO
void writeScalarFieldToTerminal(double**&);
@ -369,6 +401,11 @@ public:
tmp<volScalarField> voidfractionNuEff(volScalarField&) const;
void resetArray(double**&,int,int,double resetVal=0.);
std::vector< std::vector<double*> >* getVprobe();
std::vector< std::vector<double> >* getSprobe();
};

38
src/lagrangian/cfdemParticle/cfdemCloud/cfdemCloudI.H
View File

@ -55,11 +55,21 @@ inline const bool& cfdemCloud::impDEMdrag() const
return impDEMdrag_;
};
inline const bool& cfdemCloud::impDEMdragAcc() const
{
return impDEMdragAcc_;
};
inline const scalar& cfdemCloud::imExSplitFactor() const
{
return imExSplitFactor_;
};
inline const bool& cfdemCloud::treatVoidCellsAsExplicitForce() const
{
return treatVoidCellsAsExplicitForce_;
}
inline const scalar& cfdemCloud::cg() const
{
return cg_;
@ -105,6 +115,11 @@ inline double ** cfdemCloud::fluidVels() const
return fluidVel_;
}
inline double ** cfdemCloud::fAccs() const
{
return fAcc_;
}
inline double ** cfdemCloud::impForces() const
{
return impForces_;
@ -169,7 +184,7 @@ inline label Foam::cfdemCloud::body(int index)
inline double cfdemCloud::particleVolume(int index)
{
return particleVolumes_[index][0];
return particleV_[index][0];
}
inline scalar cfdemCloud::radius(int index)
@ -182,6 +197,25 @@ inline double cfdemCloud::d(int index)
return 2*radii_[index][0];
}
inline double cfdemCloud::d32(bool recalc)
{
if(d32_<0 || recalc)
{
scalar Ntot(0);
scalar Dtot(0);
scalar r(0);
for(int index = 0;index < numberOfParticles(); ++index)
{
r=radii_[index][0];
Ntot+=2*r*r*r;
Dtot+=r*r;
}
d32_=Ntot/Dtot;
}
return d32_;
}
inline int cfdemCloud::numberOfParticles() const
{
return numberOfParticles_;
@ -269,7 +303,7 @@ inline autoPtr<liggghtsCommandModel>* cfdemCloud::liggghtsCommand() const
}
#if defined(version21) || defined(version16ext)
#ifdef comp
#ifdef compre
inline const compressible::turbulenceModel& cfdemCloud::turbulence() const
#else
inline const incompressible::turbulenceModel& cfdemCloud::turbulence() const

8
src/lagrangian/cfdemParticle/derived/cfdemCloudIB/cfdemCloudIB.C
View File

@ -101,9 +101,10 @@ bool Foam::cfdemCloudIB::evolve()
arraysReallocated_=false;
bool doCouple=false;
if (dataExchangeM().couple())
if (dataExchangeM().doCoupleNow())
{
Info << "\n timeStepFraction() = " << dataExchangeM().timeStepFraction() << endl;
dataExchangeM().couple(0);
doCouple=true;
// Info << "skipLagrangeToEulerMapping_: " << skipLagrangeToEulerMapping_
@ -121,7 +122,7 @@ bool Foam::cfdemCloudIB::evolve()
// set void fraction field
if(verbose_) Info << "- setvoidFraction()" << endl;
voidFractionM().setvoidFraction(NULL,voidfractions_,particleWeights_,particleVolumes_);
voidFractionM().setvoidFraction(NULL,voidfractions_,particleWeights_,particleVolumes_,particleV_);
if(verbose_) Info << "setvoidFraction done." << endl;
}
@ -140,6 +141,8 @@ bool Foam::cfdemCloudIB::evolve()
// write DEM data
if(verbose_) Info << " -giveDEMdata()" << endl;
giveDEMdata();
dataExchangeM().couple(1);
haveEvolvedOnce_=true;
}
@ -190,6 +193,7 @@ void Foam::cfdemCloudIB::calcVelocityCorrection
}
//}
}
U.correctBoundaryConditions();
// make field divergence free - set reference value in case it is needed
fvScalarMatrix phiIBEqn

8
src/lagrangian/cfdemParticle/etc/additionalLibs Normal file
View File

@ -0,0 +1,8 @@
# paths for additional libraries
CFDEM_ADD_LIB_PATHS = \
# additional libraries to be linked to solvers
CFDEM_ADD_LIBS = \
# additional static libraries to be linked to lagrangian library
CFDEM_ADD_STATICLIBS = \

38
src/lagrangian/cfdemParticle/etc/bashrc
View File

@ -37,12 +37,32 @@
#- export environment variables (adapt to your paths)
#------------------------------------------------------------------------------
#check if default lammps lib path should be used
if [[ $CFDEM_LAMMPS_LIB_DIR == "" ]]; then
export CFDEM_LAMMPS_LIB_DIR=$CFDEM_LIGGGHTS_SRC_DIR/../lib/
else
echo "using CFDEM_LAMMPS_LIB_DIR=$CFDEM_LAMMPS_LIB_DIR defined by user."
fi
#- LIGGGHTS lib name
export CFDEM_LIGGGHTS_LIB_NAME=lmp_$CFDEM_LIGGGHTS_MAKEFILE_NAME
#- CFDEM lib name
export CFDEM_LIB_NAME=lagrangianCFDEM-$CFDEM_VERSION-$WM_PROJECT_VERSION
#- CFDEM compressible lib name
export CFDEM_LIB_COMP_NAME=lagrangianCFDEMcomp-$CFDEM_VERSION-$WM_PROJECT_VERSION
#check if additional libraries should be compiled together with solvers
if [[ $CFDEM_ADD_LIBS_DIR == "" ]]; then
export CFDEM_ADD_LIBS_DIR=$CFDEM_SRC_DIR/lagrangian/cfdemParticle/etc
else
echo "using CFDEM_ADD_LIBS_DIR=$CFDEM_ADD_LIBS_DIR defined by user."
fi
#-----------------------------------------------------
# additional libraries
#- LMP Many2Many lib path and makefile
export CFDEM_Many2ManyLIB_PATH=$CFDEM_SRC_DIR/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library
export CFDEM_Many2ManyLIB_MAKEFILENAME=$CFDEM_LIGGGHTS_MAKEFILE_NAME
@ -52,9 +72,17 @@ export CFDEM_M2MLIB_PATH=$CFDEM_SRC_DIR/lagrangian/cfdemParticle/subModels/dataE
export CFDEM_M2MLIB_MAKEFILENAME=$CFDEM_LIGGGHTS_MAKEFILE_NAME
#- LMP POEMS lib path and makefile
export CFDEM_POEMSLIB_PATH=$CFDEM_LIGGGHTS_SRC_DIR/../lib/poems
export CFDEM_POEMSLIB_PATH=$CFDEM_LAMMPS_LIB_DIR/poems
export CFDEM_POEMSLIB_MAKEFILENAME=g++
#- LMP ASPHERE lib path and makefile
export CFDEM_ASPHERELIB_PATH=$CFDEM_LAMMPS_LIB_DIR/poems
export CFDEM_ASPHERELIB_MAKEFILENAME=g++
#-C3PO library
export C3PO_SRC_DIR=$CFDEM_SRC_DIR/c3po
#-----------------------------------------------------
#- path to test harness
export CFDEM_TEST_HARNESS_PATH=$CFDEM_PROJECT_USER_DIR/log/logFilesCFDEM-$CFDEM_VERSION-$WM_PROJECT_VERSION
@ -145,6 +173,9 @@ alias cfdemCompCFDEMuti='bash $CFDEM_SRC_DIR/lagrangian/cfdemParticle/etc/compil
#- shortcut to test basic tutorials
alias cfdemTestTUT='bash $CFDEM_SRC_DIR/lagrangian/cfdemParticle/etc/testTutorials.sh'
#- shortcut to visualize the clock model data
alias vizClock='python $CFDEM_UT_DIR/vizClock/matPlot.py'
#- recursive touch of current directory
alias touchRec='find ./* -exec touch {} \;'
@ -160,6 +191,11 @@ export -f cfdemLiggghtsPar
cfdemGrep() { grep -rl "$1" ./* | xargs gedit; }
export -f cfdemGrep
#- shortcut lo list files in a directory
#cfdemListFiles() { find $1 | sed s:""$1"":: > listOfFiles.txt; } #leave out the dir iteslf in list
cfdemListFiles() { find $1 > listOfFiles.txt; } #keep the dir in list
export -f cfdemListFiles
# check if the run directory exists
if [ -d "$CFDEM_PROJECT_USER_DIR" ]; then
:

26
src/lagrangian/cfdemParticle/etc/cfdemSystemTest.sh
View File

@ -11,10 +11,15 @@ source $CFDEM_SRC_DIR/lagrangian/cfdemParticle/etc/functions.sh
#- show gcc settings
checkGPP="true"
#- sys check for add on
checkAddOn="true"
#- system settings
echo "*******************"
echo "system settings:"
echo "*******************"
printHeader
echo "*********************************"
echo "CFDEM(R)coupling system settings:"
echo "*********************************"
echo "CFDEM_VERSION=$CFDEM_VERSION"
echo "couple to OF_VERSION=$WM_PROJECT_VERSION"
@ -29,8 +34,10 @@ checkDirComment "$CFDEM_SOLVER_DIR" '$CFDEM_SOLVER_DIR' "yes"
checkDirComment "$CFDEM_TUT_DIR" '$CFDEM_TUT_DIR' "yes"
checkDirComment "$CFDEM_LIGGGHTS_SRC_DIR" '$CFDEM_LIGGGHTS_SRC_DIR' "yes"
checkDirComment "$CFDEM_LPP_DIR" '$CFDEM_LPP_DIR' "yes"
checkDirComment "$CFDEM_ADD_LIBS_DIR" '$CFDEM_ADD_LIBS_DIR' "yes"
checkDirComment "$CFDEM_PIZZA_DIR" '$CFDEM_PIZZA_DIR' "no"
checkDirComment "$CFDEM_TEST_HARNESS_PATH" '$CFDEM_TEST_HARNESS_PATH' "no"
checkDirComment "$C3PO_SRC_DIR" '$C3PO_SRC_DIR' "no"
echo ""
echo "library names"
@ -61,3 +68,16 @@ if [ $checkGPP == "true" ]
mpirun --version
fi
echo "**********************"
echo "additional packages..."
if [ $checkAddOn == "true" ]
then
packageName=c3po
filePath=$CFDEM_SRC_DIR/$packageName
if [ $(checkDir $filePath) == "true" ]; then
source $filePath/etc/$packageName"SystemTest.sh"
else
echo "$packageName does not exist."
fi
fi

19
src/lagrangian/cfdemParticle/etc/compileCFDEMcoupling_sol.sh
View File

@ -51,15 +51,16 @@ else
logpath="$(dirname "$(readlink -f ${BASH_SOURCE[0]})")/$logDir"
##number of solvers compiled at a time
nsteps=$WM_NCOMPPROCS
echo "do compilation on $nsteps procs"
nchunk=`echo $njobs/$nsteps+1 | bc`
if [[ $WM_NCOMPPROCS == "" ]]; then
echo "do compilation in serial"
if [[ $WM_NCOMPPROCS == "" ]] || [ $WM_NCOMPPROCS -eq 1 ]; then
nsteps=1
nchunk=1
else
echo "do compilation on $nsteps procs in $nchunk chunks"
let nchunk=$njobs+1 # +1, to wait for the last compilation too
echo "do compilation in serial"
else
nsteps=$WM_NCOMPPROCS
nchunk=`echo $njobs/$nsteps+1 | bc`
echo "do compilation on $nsteps procs in $nchunk chunks"
let nchunk++ # +1, to wait for the last compilation too
fi
counter=0
@ -108,6 +109,8 @@ else
let counter++
fi
done
sleep 1 # wait a second until compilation starts
done
echo "compilation done."

84
src/lagrangian/cfdemParticle/etc/compileCFDEMcoupling_src.sh
View File

@ -20,48 +20,48 @@ mkdir -p $logDir
#================================================================================#
# compile src
#================================================================================#
whitelist="$CFDEM_SRC_DIR/lagrangian/cfdemParticle/etc/library-list.txt"
echo ""
echo "Please provide the libraries to be compiled in the $CWD/$whitelist file."
whitelist="$CFDEM_SRC_DIR/lagrangian/cfdemParticle/etc/library-list.txt"
echo ""
echo "Please provide the libraries to be compiled in the $CWD/$whitelist file."
if [ ! -f "$CWD/$whitelist" ];then
echo "$whitelist does not exist in $CWD. Nothing will be done."
NLINES=0
COUNT=0
else
NLINES=`wc -l < $CWD/$whitelist`
COUNT=0
fi
if [ ! -f "$CWD/$whitelist" ];then
echo "$whitelist does not exist in $CWD. Nothing will be done."
NLINES=0
COUNT=0
else
NLINES=`wc -l < $CWD/$whitelist`
COUNT=0
fi
while [ $COUNT -lt $NLINES ]
do
let COUNT++
LINE=`head -n $COUNT $CWD/$whitelist | tail -1`
# white lines
if [[ "$LINE" == "" ]]; then
echo "compile $LINE"
continue
# comments
elif [[ "$LINE" == \#* ]]; then
continue
# paths
elif [[ "$LINE" == */dir ]]; then
echo "will change path..."
LINE=$(echo "${LINE%????}")
path="$CFDEM_SRC_DIR/$LINE"
cd $path
#continue
fi
#--------------------------------------------------------------------------------#
#- define variables
logpath="$(dirname "$(readlink -f ${BASH_SOURCE[0]})")/$logDir"
logfileName="log_compileCFDEMcoupling_"$(basename $LINE)""
casePath="$path"
headerText="$logfileName""-$NOW"
#--------------------------------------------------------------------------------#
# remove old log file
rm "$logpath/$logfileName"*
compileLib $logpath $logfileName $casePath $headerText
done
while [ $COUNT -lt $NLINES ]
do
let COUNT++
LINE=`head -n $COUNT $CWD/$whitelist | tail -1`
# white lines
if [[ "$LINE" == "" ]]; then
echo "compile $LINE"
continue
# comments
elif [[ "$LINE" == \#* ]]; then
continue
# paths
elif [[ "$LINE" == */dir ]]; then
echo "will change path..."
LINE=$(echo "${LINE%????}")
path="$CFDEM_SRC_DIR/$LINE"
cd $path
#continue
fi
#--------------------------------------------------------------------------------#
#- define variables
logpath="$(dirname "$(readlink -f ${BASH_SOURCE[0]})")/$logDir"
logfileName="log_compileCFDEMcoupling_"$(basename $LINE)""
casePath="$path"
headerText="$logfileName""-$NOW"
#--------------------------------------------------------------------------------#
# remove old log file
rm "$logpath/$logfileName"*
compileLib $logpath $logfileName $casePath $headerText
done

117
src/lagrangian/cfdemParticle/etc/compileCFDEMcoupling_uti.sh
View File

@ -1,36 +1,115 @@
#!/bin/bash
#===================================================================#
# compile routine for CFDEMcoupling solvers, part of CFDEMproject
# compile routine for CFDEMcoupling utilities, part of CFDEMproject
# Christoph Goniva - May. 2012, DCS Computing GmbH
#===================================================================#
whitelist="utilities-list.txt"
#- include functions
source $CFDEM_SRC_DIR/lagrangian/cfdemParticle/etc/functions.sh
NOW="$(date +"%Y-%m-%d-%H:%M")"
logDir="log"
cd $CFDEM_SRC_DIR/lagrangian/cfdemParticle/etc
mkdir -p $logDir
#================================================================================#
# compile utilities
#================================================================================#
CWD="$(dirname "$(readlink -f ${BASH_SOURCE[0]})")"
NOW="$(date +"%Y-%m-%d-%H:%M")"
echo ""
echo "This routine will compile the utilities specified in utilities-list.txt"
echo ""
#echo "Are the variables CFDEM_UT_DIR=$CFDEM_UT_DIR"
#echo "and CFDEM_SRC_DIR=$CFDEM_SRC_DIR/lagrangian/cfdemParticle correct? (y/n)"
#read YN
#if [ "$YN" != "y" ];then
# echo "Aborted by user."
# exit 1
#fi
echo ""
echo "Please provide the utilities to be compiled in the $CWD/$whitelist file."
echo "structure:"
echo "path to provide the path relative to CFDEM_UT_DIR"
echo ""
echo "example:"
echo "cfdemPostproc/dir"
echo ""
if [ ! -f "$CWD/$whitelist" ];then
echo "$whitelist does not exist in $CWD"
else
njobs=`wc -l < $CWD/$whitelist`
echo ""
echo "running compilation in pseudo-parallel mode of $njobs utilities"
for utName in "cfdemPostproc"
do
#--------------------------------------------------------------------------------#
#- define variables
logpath="$(dirname "$(readlink -f ${BASH_SOURCE[0]})")/$logDir"
logfileName="log_compileCFDEMcoupling""_$utName"
casePath="$CFDEM_UT_DIR/$utName"
headerText="$logfileName""_$utName""-$NOW"
#--------------------------------------------------------------------------------#
compileSolver $logpath $logfileName $casePath $headerText
done
echo "Note: the list of utilities compiled might be incomplete."
echo "please check $CFDEM_UT_DIR for more utilities available"
##number of utilities compiled at a time
if [[ $WM_NCOMPPROCS == "" ]] || [ $WM_NCOMPPROCS -eq 1 ]; then
nsteps=1
let nchunk=$njobs+1 # +1, to wait for the last compilation too
echo "do compilation in serial"
else
nsteps=$WM_NCOMPPROCS
nchunk=`echo $njobs/$nsteps+1 | bc`
echo "do compilation on $nsteps procs in $nchunk chunks"
let nchunk++ # +1, to wait for the last compilation too
fi
counter=0
for i in `seq $nchunk`
do
#wait until prev. compilation is finished
echo "waiting..."
until [ `ps -C make | wc -l` -eq 1 ];
do
sleep 2
done
for j in `seq $nsteps`
do
let solNr=($i-1)*$nsteps+$j
LINE=`head -n $solNr $CWD/$whitelist | tail -1`
# white lines
if [[ "$LINE" == "" ]]; then
continue
# comments
elif [[ "$LINE" == \#* ]]; then
continue
# paths
elif [[ "$LINE" == */dir ]]; then
#echo "change path"
LINE=$(echo "${LINE%????}")
path="$CFDEM_UT_DIR/$LINE"
#cd $path
let solNr++
fi
if [[ "$counter" -lt "$njobs" ]]; then
#--------------------------------------------------------------------------------#
#- define variables
#logpath="$(dirname "$(readlink -f ${BASH_SOURCE[0]})")/$logDir"
logfileName="log_compileCFDEMcoupling""_$LINE"
casePath="$CFDEM_UT_DIR/$LINE"
headerText="$logfileName""_$LINE""-$NOW"
parallel="true"
#--------------------------------------------------------------------------------#
#echo "compiling $LINE"
compileSolver $logpath $logfileName $casePath $headerText $parallel
let counter++
fi
done
sleep 1 # wait a second until compilation starts
done
echo "compilation done."
fi

4
src/lagrangian/cfdemParticle/etc/compileLIGGGHTS_lib.sh
View File

@ -33,6 +33,8 @@ mkdir -p $logDir
COUNT=0
fi
logpath="$(dirname "$(readlink -f ${BASH_SOURCE[0]})")/$logDir"
while [ $COUNT -lt $NLINES ]
do
let COUNT++
@ -57,7 +59,7 @@ mkdir -p $logDir
#--------------------------------------------------------------------------------#
#- define variables
logpath="$(dirname "$(readlink -f ${BASH_SOURCE[0]})")/$logDir"
#logpath="$(dirname "$(readlink -f ${BASH_SOURCE[0]})")/$logDir"
logfileName="log_compile$LINE""lib"
headerText="$logfileName""-$NOW"
libVarMakefileName="CFDEM_$LINE""LIB_MAKEFILENAME"

23
src/lagrangian/cfdemParticle/etc/cshrc
View File

@ -37,12 +37,32 @@
#- export environment variables (adapt to your paths)
#------------------------------------------------------------------------------
#check if default lammps lib path should be used
if ( ! ($?CFDEM_LAMMPS_LIB_DIR) ) then
setenv CFDEM_LAMMPS_LIB_DIR $CFDEM_LIGGGHTS_SRC_DIR"/../lib/"
else
echo "using CFDEM_LAMMPS_LIB_DIR=$CFDEM_LAMMPS_LIB_DIR defined by user."
endif
#- LIGGGHTS lib name
setenv CFDEM_LIGGGHTS_LIB_NAME lmp_$CFDEM_LIGGGHTS_MAKEFILE_NAME
#- CFDEM lib name
setenv CFDEM_LIB_NAME lagrangianCFDEM-$CFDEM_VERSION-$WM_PROJECT_VERSION
#- CFDEM compressible lib name
setenv CFDEM_LIB_COMP_NAME lagrangianCFDEMcomp-$CFDEM_VERSION-$WM_PROJECT_VERSION
#check if additional libraries should be compiled together with solvers
if ( ! ($?CFDEM_ADD_LIBS_DIR) ) then
setenv CFDEM_ADD_LIBS_DIR $CFDEM_SRC_DIR/lagrangian/cfdemParticle/etc
else
echo "using CFDEM_ADD_LIBS_DIR=$CFDEM_ADD_LIBS_DIR defined by user."
endif
#-----------------------------------------------------
# additional libraries
#- LMP Many2Many lib path and makefile
setenv CFDEM_Many2ManyLIB_PATH $CFDEM_SRC_DIR/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library
setenv CFDEM_Many2ManyLIB_MAKEFILENAME $CFDEM_LIGGGHTS_MAKEFILE_NAME
@ -143,6 +163,9 @@ alias cfdemCompCFDEMuti 'bash $CFDEM_SRC_DIR/lagrangian/cfdemParticle/etc/compil
#- shortcut to test basic tutorials
alias cfdemTestTUT 'bash $CFDEM_SRC_DIR/lagrangian/cfdemParticle/etc/testTutorials.sh'
#- shortcut to visualize the clock model data
alias vizClock 'python $CFDEM_UT_DIR/vizClock/matPlot.py'
#- shortcut to run liggghts in serial
alias cfdemLiggghts '$CFDEM_LIGGGHTS_SRC_DIR/lmp_$CFDEM_LIGGGHTS_MAKEFILE_NAME'

82
src/lagrangian/cfdemParticle/etc/functions.sh
View File

@ -71,6 +71,21 @@ compileLib()
#- wclean and wmake
#if [ $doClean != "noClean" ]; then
# check library to compile is compressible
str=$casePath
i=$((${#str}-4))
ending=${str:$i:4}
if [[ $ending == "Comp" ]]; then
echo "Compiling a compressible library - so doing an rmdepall of incomp library first."
echo "Please make sure to have the compressible libraries first in the library-list.txt!"
cd $CFDEM_SRC_DIR/lagrangian/cfdemParticle
echo "changing to $PWD"
rmdepall 2>&1 | tee -a $logpath/$logfileName
cd $casePath
echo "changing to $PWD"
else
echo "Compiling a incompressible library."
fi
rmdepall 2>&1 | tee -a $logpath/$logfileName
wclean 2>&1 | tee -a $logpath/$logfileName
#fi
@ -139,6 +154,7 @@ compileLIGGGHTS()
logpath="$1"
logfileName="$2"
headerText="$3"
clean="$4"
#--------------------------------------------------------------------------------#
#- clean up old log file
@ -157,9 +173,14 @@ compileLIGGGHTS()
echo 2>&1 | tee -a $logpath/$logfileName
#- wclean and wmake
rm $CFDEM_LIGGGHTS_SRC_DIR/"lmp_"$CFDEM_LIGGGHTS_MAKEFILE_NAME
rm $CFDEM_LIGGGHTS_SRC_DIR/"lib"$CFDEM_LIGGGHTS_LIB_NAME".a"
make clean-all 2>&1 | tee -a $logpath/$logfileName
if [[ $clean == "false" ]]; then
echo "not cleaning LIGGGHTS"
else
rm $CFDEM_LIGGGHTS_SRC_DIR/"lmp_"$CFDEM_LIGGGHTS_MAKEFILE_NAME
rm $CFDEM_LIGGGHTS_SRC_DIR/"lib"$CFDEM_LIGGGHTS_LIB_NAME".a"
make clean-all 2>&1 | tee -a $logpath/$logfileName
echo "cleaning LIGGGHTS"
fi
if [[ $WM_NCOMPPROCS == "" ]]; then
echo "compiling LIGGGHTS on one CPU"
make $CFDEM_LIGGGHTS_MAKEFILE_NAME 2>&1 | tee -a $logpath/$logfileName
@ -189,7 +210,13 @@ compileLMPlib()
rm $logpath/$logfileName
#- change path
cd $libraryPath
if [ -d "$libraryPath" ]; then
cd $libraryPath
else
echo ""
echo "lib path $libraryPath does not exist - check settings in .../etc/bashrc."
read
fi
#- header
echo 2>&1 | tee -a $logpath/$logfileName
@ -550,6 +577,7 @@ parCFDDEMrun()
machineFileName="$7"
debugMode="$8"
reconstuctCase="$9"
cleanCase="$10"
#--------------------------------------------------------------------------------#
if [ $debugMode == "on" ]; then
@ -903,6 +931,52 @@ checkDirComment()
fi
}
#========================================#
#- function to check if a variable exits
checkEnv()
{
#--------------------------------------------------------------------------------#
#- define variables
var="$1"
#--------------------------------------------------------------------------------#
if [[ $var == "" ]]; then
echo "false"
else
echo "true"
fi
}
#========================================#
#- function to check if a variable exits
checkEnvComment()
{
#--------------------------------------------------------------------------------#
#- define variables
var="$1"
varName="$2"
critical="$3"
#--------------------------------------------------------------------------------#
if [ $(checkEnv $var) == "true" ]; then
echo "valid:yes critical:$critical - $varName = $var"
else
echo "valid:NO critical:$critical - $varName = $var variable not set!"
fi
}
#========================================#
#- function to print a header to terminal
printHeader()
{
echo ""
echo "*********************************************"
echo "* C F D E M (R) c o u p l i n g *"
echo "* *"
echo "* by DCS Computing GmbH *"
echo "* www.dcs-computing.com *"
echo "*********************************************"
echo ""
}
#========================================#
#- track memory usage
trackMem()

9
src/lagrangian/cfdemParticle/subModels/IOModel/IOModel/IOModel.C
View File

@ -50,6 +50,15 @@ int IOModel::dumpDEMdata() const
return -1;
}
bool IOModel::dumpNow() const
{
//bool dmp(false);
//if (time_.value()+SMALL > time_.endTime().value()-time_.deltaT().value() || time_.outputTime())
// dmp=true;
return time_.outputTime();
}
fileName IOModel::createTimeDir(fileName path) const
{
fileName timeDirPath(path/time_.timeName());

2
src/lagrangian/cfdemParticle/subModels/IOModel/IOModel/IOModel.H
View File

@ -113,6 +113,8 @@ public:
virtual int dumpDEMdata() const;
bool dumpNow() const;
fileName createTimeDir(fileName) const;
fileName createLagrangianDir(fileName) const;

2
src/lagrangian/cfdemParticle/subModels/IOModel/basicIO/basicIO.C
View File

@ -84,7 +84,7 @@ basicIO::~basicIO()
int basicIO::dumpDEMdata() const
{
if (time_.outputTime())
if (dumpNow())
{
// make time directory
if (parOutput_) lagPath_=buildFilePath(dirName_);

2
src/lagrangian/cfdemParticle/subModels/IOModel/sophIO/sophIO.C
View File

@ -79,7 +79,7 @@ int sophIO::dumpDEMdata() const
{
int npProcs(-1);
if (time_.outputTime())
if (dumpNow())
{
npProcs=basicIO::dumpDEMdata();

2
src/lagrangian/cfdemParticle/subModels/IOModel/trackIO/trackIO.C
View File

@ -77,7 +77,7 @@ int trackIO::dumpDEMdata() const
{
int npProcs(-1);
if (time_.outputTime())
if (dumpNow())
{
npProcs = sophIO::dumpDEMdata();

3
src/lagrangian/cfdemParticle/subModels/averagingModel/dense/dense.C
View File

@ -150,6 +150,9 @@ void dense::setVectorAverage
{
for(int i=0;i<3;i++) valueVec[i] = value[index][i];
weightP = weight[index][subCell];
if(weightP<SMALL) Warning << "error in dense::setVectorAverage" << endl;
// not yet implemented:
//if(weightWithParticleDensity) weightP *= particleCloud_.particleDensity()[index][0];

6
src/lagrangian/cfdemParticle/subModels/clockModel/clockModel/clockModel.C
View File

@ -181,7 +181,7 @@ void Foam::clockModel::evalFile() const
{
std::ofstream outFile;
std::string fileName(path_/"timeEval.txt");
outFile.open(fileName.data(),ios_base::app);
outFile.open(fileName.data(),ios_base::trunc);
outFile << "Time Evaluation"<<nl;
outFile << eval();
outFile.close();
@ -201,7 +201,7 @@ void Foam::clockModel::evalPar() const
strs << myrank;
fileName.append(strs.str());
fileName.append(".txt");
outFile.open(fileName.data(),ios_base::app);
outFile.open(fileName.data(),ios_base::trunc);
outFile << "Time Evaluation for Processor Nr." << myrank <<nl;
outFile << eval();
outFile.close();
@ -280,7 +280,7 @@ void Foam::clockModel::evalPar() const
if(myrank == 0)
{
std::string fileName(path_/"timeEvalFull.txt");
outFile.open(fileName.data(),ios_base::app);
outFile.open(fileName.data(),ios_base::trunc);
outFile << msg;
outFile.close();
}

14
src/lagrangian/cfdemParticle/subModels/clockModel/clockModel/clockModel.H
View File

@ -120,17 +120,17 @@ public:
// Member Functions
virtual void start(int pos) const; //start measurement with custom string identifier
virtual void start(int pos,std::string identifier) const; //start measurement with custom string identifier
virtual void stop() const; //stop last started measurement
virtual void stop(std::string identifier) const; //stop last started measurement with check if identifier is equal
virtual void start(int pos) const; //start measurement with custom string identifier
virtual void start(int pos,std::string identifier) const; //start measurement with custom string identifier
virtual void stop() const; //stop last started measurement
virtual void stop(std::string identifier) const; //stop last started measurement with check if identifier is equal
virtual std::string eval() const;
virtual void evalFile() const;
virtual void evalPar() const;
void initElems();
std::vector<int> calcShift() const; //detects empty indices in vector, when times are evaluated
void Hist() const; //calc Histogram
virtual void normHist() const; //calc normalized Histogram
std::vector<int> calcShift() const; //detects empty indices in vector, when times are evaluated
void Hist() const; //calc Histogram
virtual void normHist() const; //calc normalized Histogram
void plotHist(double,std::string,int,int) const; //plot histogramm to terminal
void getRAMUsage() const;

17
src/lagrangian/cfdemParticle/subModels/dataExchangeModel/dataExchangeModel/dataExchangeModel.C
View File

@ -137,7 +137,7 @@ void Foam::dataExchangeModel::allocateArray
{
int len=0;
if(strcmp(length,"nparticles")==0) len = particleCloud_.numberOfParticles();
else if (strcmp(length,"nbodies")==0) len = nClumpTypes_;
else if (strcmp(length,"nbodies")==0) len = particleCloud_.numberOfClumps();
else FatalError<<"call allocateArray with length, nparticles or nbodies!\n" << abort(FatalError);
allocateArray(array,initVal,width,len);
}
@ -203,7 +203,7 @@ void Foam::dataExchangeModel::destroy(double* array) const
//====
bool Foam::dataExchangeModel::couple() const
bool Foam::dataExchangeModel::couple(int i) const
{
bool coupleNow = false;
if (doCoupleNow())
@ -226,7 +226,7 @@ scalar Foam::dataExchangeModel::timeStepFraction() const
}
int Foam::dataExchangeModel::getNumberOfParticles() const
{
Warning << "ask for nr of clumps - which is not supported for this dataExchange model" << endl;
Warning << "ask for nr of particles - which is not supported for this dataExchange model" << endl;
return -1;
}
@ -235,7 +235,17 @@ int Foam::dataExchangeModel::getNumberOfClumps() const
Warning << "ask for nr of clumps - which is not supported for this dataExchange model" << endl;
return -1;
}
int Foam::dataExchangeModel::getNumberOfTypes() const
{
Warning << "ask for nr of types - which is not supported for this dataExchange model" << endl;
return -1;
}
double* Foam::dataExchangeModel::getTypeVol() const
{
Warning << "ask for type volume - which is not supported for this dataExchange model" << endl;
return NULL;
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
// Construct from components
@ -248,7 +258,6 @@ dataExchangeModel::dataExchangeModel
dict_(dict),
particleCloud_(sm),
maxNumberOfParticles_(0),
nClumpTypes_(1),
couplingStep_(0),
DEMts_(-1.),
couplingInterval_(readScalar(dict_.lookup("couplingInterval")))

11
src/lagrangian/cfdemParticle/subModels/dataExchangeModel/dataExchangeModel/dataExchangeModel.H
View File

@ -62,8 +62,6 @@ protected:
int maxNumberOfParticles_;
int nClumpTypes_;
mutable int couplingStep_;
scalar DEMts_;
@ -71,7 +69,6 @@ protected:
int couplingInterval_;
// Protected member functions
void setNumberOfParticles(int) const;
public:
@ -118,9 +115,9 @@ public:
// Member Function
inline const int& maxNumberOfParticles() const {return maxNumberOfParticles_;};
void setNumberOfParticles(int) const;
inline int nClumpTypes() const {return nClumpTypes_;};
inline const int& maxNumberOfParticles() const {return maxNumberOfParticles_;};
template <typename T>
void getData
@ -179,7 +176,7 @@ public:
virtual void destroy(double*) const;
//====
virtual bool couple() const;
virtual bool couple(int) const;
virtual scalar timeStepFraction() const;
@ -240,6 +237,8 @@ public:
virtual int getNumberOfParticles() const;
virtual int getNumberOfClumps() const;
virtual int getNumberOfTypes() const;
virtual double* getTypeVol() const;
inline void setPositions(label n,double* pos) const
{

39
src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMPI/twoWayMPI.C
View File

@ -112,8 +112,9 @@ twoWayMPI::twoWayMPI
// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
twoWayMPI::~twoWayMPI()
{}
{
if (liggghts == 1) delete lmp;
}
// * * * * * * * * * * * * * * * private Member Functions * * * * * * * * * * * * * //
char* twoWayMPI::wordToChar(word& inWord) const
@ -124,6 +125,7 @@ char* twoWayMPI::wordToChar(word& inWord) const
// * * * * * * * * * * * * * * * public Member Functions * * * * * * * * * * * * * //
void twoWayMPI::getData
(
word name,
@ -251,10 +253,10 @@ void Foam::twoWayMPI::destroy(double* array) const
}
//============
bool Foam::twoWayMPI::couple() const
bool Foam::twoWayMPI::couple(int i) const
{
bool coupleNow = false;
if (doCoupleNow())
if (i==0)
{
couplingStep_++;
coupleNow = true;
@ -306,9 +308,9 @@ bool Foam::twoWayMPI::couple() const
DEMstepsToInterrupt[ind] -= DEMstepsToInterrupt[ind-1];
}
Info << "Foam::twoWayMPI::couple(): interruptTimes=" << interruptTimes << endl;
Info << "Foam::twoWayMPI::couple(): DEMstepsToInterrupt=" << DEMstepsToInterrupt << endl;
Info << "Foam::twoWayMPI::couple(): lcModel=" << lcModel << endl;
Info << "Foam::twoWayMPI::couple(i): interruptTimes=" << interruptTimes << endl;
Info << "Foam::twoWayMPI::couple(i): DEMstepsToInterrupt=" << DEMstepsToInterrupt << endl;
Info << "Foam::twoWayMPI::couple(i): lcModel=" << lcModel << endl;
}
if(particleCloud_.liggghtsCommand()[i]().type()=="runLiggghts")
@ -391,7 +393,6 @@ bool Foam::twoWayMPI::couple() const
// give nr of particles to cloud
double newNpart = liggghts_get_maxtag(lmp);
setNumberOfParticles(newNpart);
// re-allocate arrays of cloud
@ -414,9 +415,29 @@ int Foam::twoWayMPI::getNumberOfClumps() const
return liggghts_get_maxtag_ms(lmp);
#endif
Warning << "liggghts_get_maxtag_ms(lmp) is commented here!" << endl;
Warning << "liggghts_get_maxtag_ms(lmp) is not available here!" << endl;
return -1;
}
int Foam::twoWayMPI::getNumberOfTypes() const
{
#ifdef multisphere
return liggghts_get_ntypes_ms(lmp);
#endif
Warning << "liggghts_get_maxtag_ms(lmp) is not available here!" << endl;
return -1;
}
double* Foam::twoWayMPI::getTypeVol() const
{
#ifdef multisphere
return liggghts_get_vclump_ms(lmp);
#endif
Warning << "liggghts_get_vclump_ms(lmp) is not available here!" << endl;
return NULL;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam

10
src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMPI/twoWayMPI.H
View File

@ -91,12 +91,13 @@ private:
MPI_Comm comm_liggghts;
LAMMPS_NS::LAMMPS *lmp;
// private member functions
char* wordToChar(word&) const;
protected:
LAMMPS_NS::LAMMPS *lmp;
public:
//- Runtime type information
@ -118,6 +119,7 @@ public:
// Member Functions
void getData
(
word name,
@ -160,10 +162,12 @@ public:
void destroy(int*) const;
//==============
bool couple() const;
bool couple(int) const;
int getNumberOfParticles() const;
int getNumberOfClumps() const;
int getNumberOfTypes() const;
double* getTypeVol() const;
word myType() const{return typeName; };

48
src/lagrangian/cfdemParticle/subModels/forceModel/Archimedes/Archimedes.C
View File

@ -63,8 +63,6 @@ Archimedes::Archimedes
forceModel(dict,sm),
propsDict_(dict.subDict(typeName + "Props")),
twoDimensional_(false),
densityFieldName_(propsDict_.lookup("densityFieldName")),
rho_(sm.mesh().lookupObject<volScalarField> (densityFieldName_)),
gravityFieldName_(propsDict_.lookup("gravityFieldName")),
#if defined(version21) || defined(version16ext)
g_(sm.mesh().lookupObject<uniformDimensionedVectorField> (gravityFieldName_))
@ -86,14 +84,25 @@ Archimedes::Archimedes
Info << "2-dimensional simulation - make sure DEM side is 2D" << endl;
}
if (propsDict_.found("treatExplicit")) treatExplicit_=true;
if (modelType_=="A"){
treatDEM_=true;
Info << "accounting for Archimedes only on DEM side!" << endl;
// init force sub model
setForceSubModels(propsDict_);
// define switches which can be read from dict
forceSubM(0).setSwitchesList(0,true); // activate treatExplicit switch
forceSubM(0).setSwitchesList(1,true); // activate treatForceDEM switch
forceSubM(0).readSwitches();
if (modelType_=="A" || modelType_=="Bfull"){
if(!forceSubM(0).switches()[1]) // treatDEM != true
{
Warning << "Usually model type A and Bfull need Archimedes only on DEM side only (treatForceDEM=true)! are you sure about your settings?" << endl;
}
}
if (modelType_=="B"){
treatDEM_=false;
Info << "accounting for Archimedes on DEM and CFD side!" << endl;
if(forceSubM(0).switches()[1]) // treatDEM = true
{
Warning << "Usually model type B needs Archimedes on CFD and DEM side (treatForceDEM=false)! are you sure about your settings?" << endl;
}
}
particleCloud_.checkCG(true);
@ -127,12 +136,21 @@ void Archimedes::setForce() const
if(twoDimensional_)
{
force = -g_.value()*rho_[cellI]*pow(dp,2)/4*M_PI;
force = -g_.value()*forceSubM(0).rhoField()[cellI]*pow(dp,2)/4*M_PI;
Warning << "Archimedes::setForce() : this functionality is not tested!" << endl;
}else{
force = -g_.value()*rho_[cellI]*particleCloud_.particleVolume(index);
force = -g_.value()*forceSubM(0).rhoField()[cellI]*particleCloud_.particleVolume(index);
}
//if(index >=0 && index <100)
//{
// Pout << "cellI = " << cellI << endl;
// Pout << "index = " << index << endl;
// Pout << "forceSubM(0).rhoField()[cellI] = " << forceSubM(0).rhoField()[cellI] << endl;
// Pout << "particleCloud_.particleVolume(index) = " << particleCloud_.particleVolume(index) << endl;
// Pout << "force = " << force << endl;
//}
//Set value fields and write the probe
if(probeIt_)
{
@ -143,14 +161,8 @@ void Archimedes::setForce() const
}
}
if(!treatDEM_)
{
if(treatExplicit_)
for(int j=0;j<3;j++) expForces()[index][j] += force[j];
else
for(int j=0;j<3;j++) impForces()[index][j] += force[j];
}
for(int j=0;j<3;j++) DEMForces()[index][j] += force[j];
// write particle based data to global array
forceSubM(0).partToArray(index,force,vector::zero);
//}
}
}

4
src/lagrangian/cfdemParticle/subModels/forceModel/Archimedes/Archimedes.H
View File

@ -63,10 +63,6 @@ private:
bool twoDimensional_;
word densityFieldName_;
const volScalarField& rho_; // ref to fluid density
word gravityFieldName_;
#ifdef version21

28
src/lagrangian/cfdemParticle/subModels/forceModel/ArchimedesIB/ArchimedesIB.C
View File

@ -63,8 +63,6 @@ ArchimedesIB::ArchimedesIB
forceModel(dict,sm),
propsDict_(dict.subDict(typeName + "Props")),
twoDimensional_(false),
densityFieldName_(propsDict_.lookup("densityFieldName")),
rho_(sm.mesh().lookupObject<volScalarField> (densityFieldName_)),
voidfractionFieldName_(propsDict_.lookup("voidfractionFieldName")), //mod by alice
voidfractions_(sm.mesh().lookupObject<volScalarField> (voidfractionFieldName_)),//mod by alice
gravityFieldName_(propsDict_.lookup("gravityFieldName")),
@ -85,9 +83,18 @@ ArchimedesIB::ArchimedesIB
Info << "2-dimensional simulation - make sure DEM side is 2D" << endl;
}
if (propsDict_.found("treatExplicit")) treatExplicit_=true;
treatDEM_=true;
// init force sub model
setForceSubModels(propsDict_);
// define switches which can be read from dict
forceSubM(0).setSwitchesList(0,true); // activate treatExplicit switch
// read those switches defined above, if provided in dict
forceSubM(0).readSwitches();
forceSubM(0).setSwitches(1,true); // treatDEM = true
Info << "accounting for Archimedes only on DEM side!" << endl;
particleCloud_.checkCG(true);
}
@ -116,8 +123,8 @@ void ArchimedesIB::setForce() const
label cellI = particleCloud_.cellIDs()[index][subCell];
if (cellI > -1) // particle Found
{
//force += -g_.value()*rho_[cellI]*rho_.mesh().V()[cellI]*(1-particleCloud_.voidfractions()[index][subCell]);//mod by alice
force += -g_.value()*rho_[cellI]*rho_.mesh().V()[cellI]*(1-voidfractions_[cellI]);//mod by alice
//force += -g_.value()*forceSubM(0).rhoField()[cellI]*forceSubM(0).rhoField().mesh().V()[cellI]*(1-particleCloud_.voidfractions()[index][subCell]);//mod by alice
force += -g_.value()*forceSubM(0).rhoField()[cellI]*particleCloud_.mesh().V()[cellI]*(1-voidfractions_[cellI]);//mod by alice
}
}
@ -131,12 +138,9 @@ void ArchimedesIB::setForce() const
// set force on particle
if(twoDimensional_) Warning<<"ArchimedesIB model doesn't work for 2D right now!!\n"<< endl;
if(!treatDEM_)
{
if(treatExplicit_) for(int j=0;j<3;j++) expForces()[index][j] += force[j];
else for(int j=0;j<3;j++) impForces()[index][j] += force[j];
}
for(int j=0;j<3;j++) DEMForces()[index][j] += force[j];
// write particle based data to global array
forceSubM(0).partToArray(index,force,vector::zero);
//}
}
}

4
src/lagrangian/cfdemParticle/subModels/forceModel/ArchimedesIB/ArchimedesIB.H
View File

@ -65,10 +65,6 @@ private:
bool twoDimensional_;
word densityFieldName_;
const volScalarField& rho_; // ref to fluid density
word voidfractionFieldName_;
const volScalarField& voidfractions_;

92
src/lagrangian/cfdemParticle/subModels/forceModel/DiFeliceDrag/DiFeliceDrag.C
View File

@ -64,15 +64,10 @@ DiFeliceDrag::DiFeliceDrag
:
forceModel(dict,sm),
propsDict_(dict.subDict(typeName + "Props")),
verbose_(false),
velFieldName_(propsDict_.lookup("velFieldName")),
U_(sm.mesh().lookupObject<volVectorField> (velFieldName_)),
densityFieldName_(propsDict_.lookup("densityFieldName")),
rho_(sm.mesh().lookupObject<volScalarField> (densityFieldName_)),
voidfractionFieldName_(propsDict_.lookup("voidfractionFieldName")),
voidfraction_(sm.mesh().lookupObject<volScalarField> (voidfractionFieldName_)),
interpolation_(false),
splitImplicitExplicit_(false),
UsFieldName_(propsDict_.lookup("granVelFieldName")),
UsField_(sm.mesh().lookupObject<volVectorField> (UsFieldName_)),
scaleDia_(1.),
@ -87,25 +82,24 @@ DiFeliceDrag::DiFeliceDrag
particleCloud_.probeM().scalarFields_.append("voidfraction"); //other are debug
particleCloud_.probeM().writeHeader();
if (propsDict_.found("verbose")) verbose_=true;
if (propsDict_.found("treatExplicit")) treatExplicit_=true;
if (propsDict_.found("interpolation"))
{
Info << "using interpolated value of U." << endl;
interpolation_=true;
}
if (propsDict_.found("splitImplicitExplicit"))
{
Info << "will split implicit / explicit force contributions." << endl;
splitImplicitExplicit_ = true;
if(!interpolation_)
Info << "WARNING: will only consider fluctuating particle velocity in implicit / explicit force split!" << endl;
}
particleCloud_.checkCG(true);
if (propsDict_.found("scale"))
scaleDia_=scalar(readScalar(propsDict_.lookup("scale")));
if (propsDict_.found("scaleDrag"))
scaleDrag_=scalar(readScalar(propsDict_.lookup("scaleDrag")));
// init force sub model
setForceSubModels(propsDict_);
// define switches which can be read from dict
forceSubM(0).setSwitchesList(0,true); // activate treatExplicit switch
forceSubM(0).setSwitchesList(2,true); // activate implDEM switch
forceSubM(0).setSwitchesList(3,true); // activate search for verbose switch
forceSubM(0).setSwitchesList(4,true); // activate search for interpolate switch
forceSubM(0).setSwitchesList(8,true); // activate scalarViscosity switch
// read those switches defined above, if provided in dict
forceSubM(0).readSwitches();
}
@ -125,18 +119,16 @@ void DiFeliceDrag::setForce() const
scaleDia_=particleCloud_.cg();
Info << "DiFeliceDrag using scale from liggghts cg = " << scaleDia_ << endl;
}
// get viscosity field
#ifdef comp
const volScalarField nufField = particleCloud_.turbulence().mu() / rho_;
#else
const volScalarField& nufField = particleCloud_.turbulence().nu();
#endif
const volScalarField& nufField = forceSubM(0).nuField();
const volScalarField& rhoField = forceSubM(0).rhoField();
vector position(0,0,0);
scalar voidfraction(1);
vector Ufluid(0,0,0);
vector drag(0,0,0);
vector dragExplicit(0,0,0);
scalar dragCoefficient(0);
label cellI=0;
vector Us(0,0,0);
vector Ur(0,0,0);
@ -147,11 +139,6 @@ void DiFeliceDrag::setForce() const
scalar Rep(0);
scalar Cd(0);
vector UfluidFluct(0,0,0);
vector UsFluct(0,0,0);
vector dragExplicit(0,0,0);
scalar dragCoefficient(0);
interpolationCellPoint<scalar> voidfractionInterpolator_(voidfraction_);
interpolationCellPoint<vector> UInterpolator_(U_);
@ -159,15 +146,15 @@ void DiFeliceDrag::setForce() const
for(int index = 0;index < particleCloud_.numberOfParticles(); index++)
{
//if(mask[index][0])
//{
cellI = particleCloud_.cellIDs()[index][0];
drag = vector(0,0,0);
dragExplicit = vector(0,0,0);
dragCoefficient=0;
Ufluid =vector(0,0,0);
if (cellI > -1) // particle Found
{
if(interpolation_)
if(forceSubM(0).interpolation())
{
position = particleCloud_.position(index);
voidfraction = voidfractionInterpolator_.interpolate(position,cellI);
@ -182,11 +169,10 @@ void DiFeliceDrag::setForce() const
Ur = Ufluid-Us;
ds = 2*particleCloud_.radius(index);
nuf = nufField[cellI];
rho = rho_[cellI];
rho = rhoField[cellI];
magUr = mag(Ur);
Rep = 0;
Cd = 0;
dragCoefficient = 0;
if (magUr > 0)
{
@ -212,16 +198,10 @@ void DiFeliceDrag::setForce() const
drag = dragCoefficient*Ur; //total drag force!
//Split forces
if(splitImplicitExplicit_)
{
UfluidFluct = Ufluid - U_[cellI];
UsFluct = Us - UsField_[cellI];
dragExplicit = dragCoefficient*(UfluidFluct - UsFluct); //explicit part of force
}
forceSubM(0).explicitCorr(drag,dragExplicit,dragCoefficient,Ufluid,U_[cellI],Us,UsField_[cellI],forceSubM(0).verbose(),index);
}
if(verbose_ && index >-1 && index <102)
if(forceSubM(0).verbose() && index >-1 && index <102)
{
Pout << "index = " << index << endl;
Pout << "Us = " << Us << endl;
@ -233,12 +213,6 @@ void DiFeliceDrag::setForce() const
Pout << "Rep = " << Rep << endl;
Pout << "Cd = " << Cd << endl;
Pout << "drag (total) = " << drag << endl;
if(splitImplicitExplicit_)
{
Pout << "UfluidFluct = " << UfluidFluct << endl;
Pout << "UsFluct = " << UsFluct << endl;
Pout << "dragExplicit = " << dragExplicit << endl;
}
}
//Set value fields and write the probe
@ -253,20 +227,10 @@ void DiFeliceDrag::setForce() const
particleCloud_.probeM().writeProbe(index, sValues, vValues);
}
}
// set force on particle
if(treatExplicit_) for(int j=0;j<3;j++) expForces()[index][j] += drag[j];
else //implicit treatment, taking explicit force contribution into account
{
for(int j=0;j<3;j++)
{
impForces()[index][j] += drag[j] - dragExplicit[j]; //only consider implicit part!
expForces()[index][j] += dragExplicit[j];
}
}
for(int j=0;j<3;j++) DEMForces()[index][j] += drag[j];
// write particle based data to global array
forceSubM(0).partToArray(index,drag,dragExplicit,Ufluid,dragCoefficient);
}
//}
}

10
src/lagrangian/cfdemParticle/subModels/forceModel/DiFeliceDrag/DiFeliceDrag.H
View File

@ -60,24 +60,14 @@ class DiFeliceDrag
private:
dictionary propsDict_;
bool verbose_;
word velFieldName_;
const volVectorField& U_;
word densityFieldName_;
const volScalarField& rho_;
word voidfractionFieldName_;
const volScalarField& voidfraction_;
bool interpolation_; // use interpolated U field values
bool splitImplicitExplicit_; // use splitting of implicit and explict force contribution
word UsFieldName_;
const volVectorField& UsField_; // the average particle velocity field (for implicit/expliti force split)

106
src/lagrangian/cfdemParticle/subModels/forceModel/GidaspowDrag/GidaspowDrag.C
View File

@ -63,17 +63,16 @@ GidaspowDrag::GidaspowDrag
:
forceModel(dict,sm),
propsDict_(dict.subDict(typeName + "Props")),
verbose_(false),
velFieldName_(propsDict_.lookup("velFieldName")),
U_(sm.mesh().lookupObject<volVectorField> (velFieldName_)),
densityFieldName_(propsDict_.lookup("densityFieldName")),
rho_(sm.mesh().lookupObject<volScalarField> (densityFieldName_)),
voidfractionFieldName_(propsDict_.lookup("voidfractionFieldName")),
voidfraction_(sm.mesh().lookupObject<volScalarField> (voidfractionFieldName_)),
phi_(readScalar(propsDict_.lookup("phi"))),
interpolation_(false),
UsFieldName_(propsDict_.lookup("granVelFieldName")),
UsField_(sm.mesh().lookupObject<volVectorField> (UsFieldName_)),
scaleDia_(1.),
scaleDrag_(1.)
scaleDrag_(1.),
switchingVoidfraction_(0.8)
{
//Append the field names to be probed
particleCloud_.probeM().initialize(typeName, "gidaspowDrag.logDat");
@ -84,23 +83,24 @@ GidaspowDrag::GidaspowDrag
particleCloud_.probeM().scalarFields_.append("voidfraction");
particleCloud_.probeM().writeHeader();
if (propsDict_.found("verbose")) verbose_=true;
if (propsDict_.found("treatExplicit")) treatExplicit_=true;
if (propsDict_.found("interpolation")) interpolation_=true;
if (propsDict_.found("implDEM"))
{
treatExplicit_=false;
implDEM_=true;
setImpDEMdrag();
Info << "Using implicit DEM drag formulation." << endl;
}
// init force sub model
setForceSubModels(propsDict_);
// define switches which can be read from dict
forceSubM(0).setSwitchesList(0,true); // activate treatExplicit switch
forceSubM(0).setSwitchesList(2,true); // activate implDEM switch
forceSubM(0).setSwitchesList(3,true); // activate search for verbose switch
forceSubM(0).setSwitchesList(4,true); // activate search for interpolate switch
forceSubM(0).setSwitchesList(8,true); // activate scalarViscosity switch
forceSubM(0).readSwitches();
particleCloud_.checkCG(true);
if (propsDict_.found("scale"))
scaleDia_=scalar(readScalar(propsDict_.lookup("scale")));
if (propsDict_.found("scaleDrag"))
scaleDrag_=scalar(readScalar(propsDict_.lookup("scaleDrag")));
Info << "Gidaspow - interpolation switch: " << interpolation_ << endl;
if (propsDict_.found("switchingVoidfraction"))
switchingVoidfraction_ = readScalar(propsDict_.lookup("switchingVoidfraction"));
}
@ -121,12 +121,8 @@ void GidaspowDrag::setForce() const
Info << "Gidaspow using scale from liggghts cg = " << scaleDia_ << endl;
}
// get viscosity field
#ifdef comp
const volScalarField nufField = particleCloud_.turbulence().mu() / rho_;
#else
const volScalarField& nufField = particleCloud_.turbulence().nu();
#endif
const volScalarField& nufField = forceSubM(0).nuField();
const volScalarField& rhoField = forceSubM(0).rhoField();
vector position(0,0,0);
scalar voidfraction(1);
@ -145,9 +141,11 @@ void GidaspowDrag::setForce() const
scalar localPhiP(0);
scalar CdMagUrLag(0); //Cd of the very particle
scalar KslLag(0); //momentum exchange of the very particle (per unit volume)
scalar betaP(0); //momentum exchange of the very particle
vector dragExplicit(0,0,0);
scalar dragCoefficient(0);
interpolationCellPoint<scalar> voidfractionInterpolator_(voidfraction_);
interpolationCellPoint<vector> UInterpolator_(U_);
@ -159,15 +157,17 @@ void GidaspowDrag::setForce() const
//{
cellI = particleCloud_.cellIDs()[index][0];
drag = vector(0,0,0);
dragExplicit = vector(0,0,0);
betaP = 0;
Vs = 0;
Ufluid =vector(0,0,0);
voidfraction=0;
dragCoefficient = 0;
if (cellI > -1) // particle Found
{
if(interpolation_)
if( forceSubM(0).interpolation() )
{
position = particleCloud_.position(index);
voidfraction = voidfractionInterpolator_.interpolate(position,cellI);
@ -186,46 +186,47 @@ void GidaspowDrag::setForce() const
Us = particleCloud_.velocity(index);
Ur = Ufluid-Us;
magUr = mag(Ur);
ds = 2*particleCloud_.radius(index)*phi_;
rho = rho_[cellI];
ds = 2*particleCloud_.radius(index);
rho = rhoField[cellI];
nuf = nufField[cellI];
Rep=0.0;
localPhiP = 1.0f-voidfraction+SMALL;
Vs = ds*ds*ds*M_PI/6;
//Compute specific drag coefficient (i.e., Force per unit slip velocity and per m³ SUSPENSION)
if(voidfraction > 0.8) //dilute
// calc particle's drag coefficient (i.e., Force per unit slip velocity and per m³ PARTICLE)
if(voidfraction > switchingVoidfraction_) //dilute
{
Rep=ds/scaleDia_*voidfraction*magUr/nuf;
CdMagUrLag = (24.0*nuf/(ds/scaleDia_*voidfraction)) //1/magUr missing here, but compensated in expression for KslLag!
*(scalar(1)+0.15*Foam::pow(Rep, 0.687));
CdMagUrLag = (24.0*nuf/(ds/scaleDia_*voidfraction)) //1/magUr missing here, but compensated in expression for betaP!
*(scalar(1.0)+0.15*Foam::pow(Rep, 0.687));
KslLag = 0.75*(
rho*localPhiP*voidfraction*CdMagUrLag
betaP = 0.75*( //this is betaP = beta / localPhiP!
rho*voidfraction*CdMagUrLag
/
(ds/scaleDia_*Foam::pow(voidfraction,2.65))
);
}
else //dense
{
KslLag = (150*Foam::pow(localPhiP,2)*nuf*rho)/
(voidfraction*ds/scaleDia_*ds/scaleDia_+SMALL)
betaP = (150 * localPhiP*nuf*rho) //this is betaP = beta / localPhiP!
/ (voidfraction*ds/scaleDia_*phi_*ds/scaleDia_*phi_)
+
(1.75*(localPhiP) * magUr * rho)/
((ds/scaleDia_));
(1.75 * magUr * rho)
/((ds/scaleDia_*phi_));
}
// calc particle's drag coefficient (i.e., Force per unit slip velocity and per m³ PARTICLE)
betaP = KslLag / localPhiP;
// calc particle's drag
drag = Vs * betaP * Ur * scaleDrag_;
dragCoefficient = Vs*betaP*scaleDrag_;
if (modelType_=="B")
drag /= voidfraction;
dragCoefficient /= voidfraction;
if(verbose_ && index >=0 && index <2)
drag = dragCoefficient * Ur;
// explicitCorr
forceSubM(0).explicitCorr(drag,dragExplicit,dragCoefficient,Ufluid,U_[cellI],Us,UsField_[cellI],forceSubM(0).verbose());
if(forceSubM(0).verbose() && index >=0 && index <2)
{
Pout << "cellI = " << cellI << endl;
Pout << "index = " << index << endl;
@ -255,23 +256,8 @@ void GidaspowDrag::setForce() const
}
}
// set force on particle
if(treatExplicit_) for(int j=0;j<3;j++) expForces()[index][j] += drag[j];
else for(int j=0;j<3;j++) impForces()[index][j] += drag[j];
// set Cd
if(implDEM_)
{
for(int j=0;j<3;j++) fluidVel()[index][j]=Ufluid[j];
if (modelType_=="B" && cellI > -1)
Cds()[index][0] = Vs*betaP/voidfraction*scaleDrag_;
else
Cds()[index][0] = Vs*betaP*scaleDrag_;
}else{
for(int j=0;j<3;j++) DEMForces()[index][j] += drag[j];
}
// write particle based data to global array
forceSubM(0).partToArray(index,drag,dragExplicit,Ufluid,dragCoefficient);
//}// end if mask
}// end loop particles

14
src/lagrangian/cfdemParticle/subModels/forceModel/GidaspowDrag/GidaspowDrag.H
View File

@ -31,6 +31,8 @@ Description
Gidaspow drag law
- only valid for low-Reynolds number systems (Re_p<1000)
- including interpolation of the velocity to the exact position
- splits off explicit drag component due to fluctuation in fluid and particle
velocity (optional via forceSubModel "ImExCorr")
Class
GidaspowDrag
@ -62,28 +64,26 @@ class GidaspowDrag
private:
dictionary propsDict_;
bool verbose_;
word velFieldName_;
const volVectorField& U_;
word densityFieldName_;
const volScalarField& rho_;
word voidfractionFieldName_;
const volScalarField& voidfraction_;
const scalar phi_;
bool interpolation_; // use interpolated field values
word UsFieldName_;
const volVectorField& UsField_; // the average particle velocity field
mutable scalar scaleDia_;
mutable scalar scaleDrag_;
mutable scalar switchingVoidfraction_; //voidfraction above which dilute formulation will be used
public:
//- Runtime type information

94
src/lagrangian/cfdemParticle/subModels/forceModel/KochHillDrag/KochHillDrag.C
View File

@ -64,14 +64,12 @@ KochHillDrag::KochHillDrag
:
forceModel(dict,sm),
propsDict_(dict.subDict(typeName + "Props")),
verbose_(false),
velFieldName_(propsDict_.lookup("velFieldName")),
U_(sm.mesh().lookupObject<volVectorField> (velFieldName_)),
densityFieldName_(propsDict_.lookup("densityFieldName")),
rho_(sm.mesh().lookupObject<volScalarField> (densityFieldName_)),
voidfractionFieldName_(propsDict_.lookup("voidfractionFieldName")),
voidfraction_(sm.mesh().lookupObject<volScalarField> (voidfractionFieldName_)),
interpolation_(false),
UsFieldName_(propsDict_.lookupOrDefault("granVelFieldName",word("Us"))),
UsField_(sm.mesh().lookupObject<volVectorField> (UsFieldName_)),
scaleDia_(1.),
scaleDrag_(1.)
{
@ -84,16 +82,20 @@ KochHillDrag::KochHillDrag
particleCloud_.probeM().scalarFields_.append("voidfraction"); //other are debug
particleCloud_.probeM().writeHeader();
if (propsDict_.found("verbose")) verbose_=true;
if (propsDict_.found("treatExplicit")) treatExplicit_=true;
if (propsDict_.found("interpolation")) interpolation_=true;
if (propsDict_.found("implDEM"))
{
treatExplicit_=false;
implDEM_=true;
setImpDEMdrag();
Info << "Using implicit DEM drag formulation." << endl;
}
// init force sub model
setForceSubModels(propsDict_);
// define switches which can be read from dict
forceSubM(0).setSwitchesList(0,true); // activate treatExplicit switch
forceSubM(0).setSwitchesList(2,true); // activate implDEM switch
forceSubM(0).setSwitchesList(3,true); // activate search for verbose switch
forceSubM(0).setSwitchesList(4,true); // activate search for interpolate switch
forceSubM(0).setSwitchesList(7,true); // activate implForceDEMacc switch
forceSubM(0).setSwitchesList(8,true); // activate scalarViscosity switch
// read those switches defined above, if provided in dict
forceSubM(0).readSwitches();
particleCloud_.checkCG(true);
if (propsDict_.found("scale"))
@ -120,17 +122,15 @@ void KochHillDrag::setForce() const
Info << "KochHill using scale from liggghts cg = " << scaleDia_ << endl;
}
// get viscosity field
#ifdef comp
const volScalarField nufField = particleCloud_.turbulence().mu()/rho_;
#else
const volScalarField& nufField = particleCloud_.turbulence().nu();
#endif
const volScalarField& nufField = forceSubM(0).nuField();
const volScalarField& rhoField = forceSubM(0).rhoField();
vector position(0,0,0);
scalar voidfraction(1);
vector Ufluid(0,0,0);
vector drag(0,0,0);
vector dragExplicit(0,0,0);
scalar dragCoefficient(0);
label cellI=0;
vector Us(0,0,0);
@ -144,6 +144,8 @@ void KochHillDrag::setForce() const
scalar volumefraction(0);
scalar betaP(0);
int couplingInterval(particleCloud_.dataExchangeM().couplingInterval());
interpolationCellPoint<scalar> voidfractionInterpolator_(voidfraction_);
interpolationCellPoint<vector> UInterpolator_(U_);
@ -151,10 +153,10 @@ void KochHillDrag::setForce() const
for(int index = 0;index < particleCloud_.numberOfParticles(); index++)
{
//if(mask[index][0])
//{
cellI = particleCloud_.cellIDs()[index][0];
drag = vector(0,0,0);
dragExplicit = vector(0,0,0);
dragCoefficient=0;
betaP = 0;
Vs = 0;
Ufluid =vector(0,0,0);
@ -162,7 +164,7 @@ void KochHillDrag::setForce() const
if (cellI > -1) // particle Found
{
if(interpolation_)
if(forceSubM(0).interpolation())
{
position = particleCloud_.position(index);
voidfraction = voidfractionInterpolator_.interpolate(position,cellI);
@ -181,11 +183,11 @@ void KochHillDrag::setForce() const
Ur = Ufluid-Us;
ds = particleCloud_.d(index);
nuf = nufField[cellI];
rho = rho_[cellI];
rho = rhoField[cellI];
magUr = mag(Ur);
Rep = 0;
Vs = ds*ds*ds*M_PI/6;
volumefraction = 1-voidfraction+SMALL;
volumefraction = max(SMALL,min(1-SMALL,1-voidfraction));
if (magUr > 0)
{
@ -209,20 +211,32 @@ void KochHillDrag::setForce() const
// calc model coefficient F3
scalar F3 = 0.0673+0.212*volumefraction+0.0232/pow(voidfraction,5);
//Calculate F
//Calculate F (the factor 0.5 is introduced, since Koch and Hill, ARFM 33:619–47, use the radius
//to define Rep, and we use the particle diameter, see vanBuijtenen et al., CES 66:2368–2376.
scalar F = voidfraction * (F0 + 0.5*F3*Rep);
// calc drag model coefficient betaP
betaP = 18.*nuf*rho/(ds/scaleDia_*ds/scaleDia_)*voidfraction*F;
// calc particle's drag
drag = Vs*betaP*Ur*scaleDrag_;
dragCoefficient = Vs*betaP*scaleDrag_;
if (modelType_=="B")
drag /= voidfraction;
dragCoefficient /= voidfraction;
if(forceSubM(0).switches()[7]) // implForceDEMaccumulated=true
{
//get drag from the particle itself
for (int j=0 ; j<3 ; j++) drag[j] = particleCloud_.fAccs()[index][j]/couplingInterval;
}else
{
drag = dragCoefficient * Ur;
// explicitCorr
forceSubM(0).explicitCorr(drag,dragExplicit,dragCoefficient,Ufluid,U_[cellI],Us,UsField_[cellI],forceSubM(0).verbose());
}
}
if(verbose_ && index >=0 && index <2)
if(forceSubM(0).verbose() && index >=0 && index <2)
{
Pout << "cellI = " << cellI << endl;
Pout << "index = " << index << endl;
@ -250,25 +264,9 @@ void KochHillDrag::setForce() const
particleCloud_.probeM().writeProbe(index, sValues, vValues);
}
}
// set force on particle
if(treatExplicit_) for(int j=0;j<3;j++) expForces()[index][j] += drag[j];
else for(int j=0;j<3;j++) impForces()[index][j] += drag[j];
// set Cd
if(implDEM_)
{
for(int j=0;j<3;j++) fluidVel()[index][j]=Ufluid[j];
if (modelType_=="B" && cellI > -1)
Cds()[index][0] = Vs*betaP/voidfraction*scaleDrag_;
else
Cds()[index][0] = Vs*betaP*scaleDrag_;
}else{
for(int j=0;j<3;j++) DEMForces()[index][j] += drag[j];
}
//}
// write particle based data to global array
forceSubM(0).partToArray(index,drag,dragExplicit,Ufluid,dragCoefficient);
}
}

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