diff --git a/.circleci/config.yml b/.circleci/config.yml index 8ad53ebb..6a8c6c01 100644 --- a/.circleci/config.yml +++ b/.circleci/config.yml @@ -21,18 +21,18 @@ jobs: - run: name: Make project and user dir - command: mkdir -p /root/CFDEM/CFDEMcoupling && mkdir -p /root/CFDEM/-4.1 + command: mkdir -p /root/CFDEM/CFDEMcoupling && mkdir -p /root/CFDEM/-6 - checkout: path: /root/CFDEM/CFDEMcoupling - run: name: Add OpenFOAM package repository - command: sudo apt-get install -y software-properties-common wget apt-transport-https && sudo add-apt-repository http://dl.openfoam.org/ubuntu && sudo sh -c "wget -O - http://dl.openfoam.org/gpg.key | apt-key add -" + command: sudo apt-get install -y software-properties-common wget apt-transport-https && sudo sh -c "wget -O - http://dl.openfoam.org/gpg.key | apt-key add -" && sudo add-apt-repository http://dl.openfoam.org/ubuntu - run: - name: Install OpenFOAM 4.1 - command: sudo apt-get update && sudo apt-get -y install openfoam4 + name: Install OpenFOAM 6 + command: sudo apt-get update && sudo apt-get -y install openfoam6 - run: name: Clone LIGGGHTS repository @@ -42,7 +42,7 @@ jobs: name: Build LIGGGHTS command: > shopt -s expand_aliases && - source /opt/openfoam4/etc/bashrc && + source /opt/openfoam6/etc/bashrc && source /root/CFDEM/CFDEMcoupling/etc/bashrc && bash /root/CFDEM/CFDEMcoupling/etc/compileLIGGGHTS.sh no_output_timeout: 30m @@ -51,7 +51,7 @@ jobs: name: Build CFDEMcoupling library command: > shopt -s expand_aliases && - source /opt/openfoam4/etc/bashrc && + source /opt/openfoam6/etc/bashrc && source /root/CFDEM/CFDEMcoupling/etc/bashrc && bash /root/CFDEM/CFDEMcoupling/etc/compileCFDEMcoupling_src.sh @@ -59,7 +59,7 @@ jobs: name: Build CFDEMcoupling solvers command: > shopt -s expand_aliases && - source /opt/openfoam4/etc/bashrc && + source /opt/openfoam6/etc/bashrc && source /root/CFDEM/CFDEMcoupling/etc/bashrc && bash /root/CFDEM/CFDEMcoupling/etc/compileCFDEMcoupling_sol.sh @@ -67,6 +67,6 @@ jobs: name: Build CFDEMcoupling utilities command: > shopt -s expand_aliases && - source /opt/openfoam4/etc/bashrc && + source /opt/openfoam6/etc/bashrc && source /root/CFDEM/CFDEMcoupling/etc/bashrc && bash /root/CFDEM/CFDEMcoupling/etc/compileCFDEMcoupling_uti.sh diff --git a/applications/solvers/cfdemSolverMultiphase/Make/options b/applications/solvers/cfdemSolverMultiphase/Make/options index 0cfccaff..6244f98f 100644 --- a/applications/solvers/cfdemSolverMultiphase/Make/options +++ b/applications/solvers/cfdemSolverMultiphase/Make/options @@ -1,6 +1,10 @@ +FOAM_VERSION_MAJOR := $(word 1,$(subst ., ,$(WM_PROJECT_VERSION))) +PFLAGS+= -DOPENFOAM_VERSION_MAJOR=$(FOAM_VERSION_MAJOR) + include $(CFDEM_ADD_LIBS_DIR)/additionalLibs EXE_INC = \ + $(PFLAGS) \ -I$(CFDEM_OFVERSION_DIR) \ -ImultiphaseMixture/lnInclude \ -I$(LIB_SRC)/transportModels \ diff --git a/applications/solvers/cfdemSolverMultiphase/cfdemSolverMultiphase.C b/applications/solvers/cfdemSolverMultiphase/cfdemSolverMultiphase.C index 42e514da..312b9168 100644 --- a/applications/solvers/cfdemSolverMultiphase/cfdemSolverMultiphase.C +++ b/applications/solvers/cfdemSolverMultiphase/cfdemSolverMultiphase.C @@ -46,6 +46,11 @@ Description int main(int argc, char *argv[]) { + +#if OPENFOAM_VERSION_MAJOR >= 6 + FatalError << "cfdemSolverMultiphase requires OpenFOAM 4.x or 5.x to work properly" << exit(FatalError); +#endif + #include "postProcess.H" #include "setRootCase.H" #include "createTime.H" diff --git a/applications/solvers/cfdemSolverMultiphase/multiphaseMixture/Make/options b/applications/solvers/cfdemSolverMultiphase/multiphaseMixture/Make/options index f8ffa1cf..15704fe2 100644 --- a/applications/solvers/cfdemSolverMultiphase/multiphaseMixture/Make/options +++ b/applications/solvers/cfdemSolverMultiphase/multiphaseMixture/Make/options @@ -1,4 +1,8 @@ +FOAM_VERSION_MAJOR := $(word 1,$(subst ., ,$(WM_PROJECT_VERSION))) +PFLAGS+= -DOPENFOAM_VERSION_MAJOR=$(FOAM_VERSION_MAJOR) + EXE_INC = \ + $(PFLAGS) \ -IalphaContactAngle \ -I$(LIB_SRC)/transportModels \ -I$(LIB_SRC)/transportModels/incompressible/lnInclude \ diff --git a/applications/solvers/cfdemSolverMultiphase/multiphaseMixture/multiphaseMixture.C b/applications/solvers/cfdemSolverMultiphase/multiphaseMixture/multiphaseMixture.C index a1ea1537..e152ae4a 100644 --- a/applications/solvers/cfdemSolverMultiphase/multiphaseMixture/multiphaseMixture.C +++ b/applications/solvers/cfdemSolverMultiphase/multiphaseMixture/multiphaseMixture.C @@ -679,8 +679,13 @@ void Foam::multiphaseMixture::solveAlphas alphaPhiCorr, zeroField(), zeroField(), +#if OPENFOAM_VERSION_MAJOR < 6 1, 0, +#else + oneField(), + zeroField(), +#endif true ); diff --git a/applications/solvers/cfdemSolverRhoPimple/Make/options b/applications/solvers/cfdemSolverRhoPimple/Make/options index 0377ece5..cc67835d 100644 --- a/applications/solvers/cfdemSolverRhoPimple/Make/options +++ b/applications/solvers/cfdemSolverRhoPimple/Make/options @@ -1,5 +1,7 @@ include $(CFDEM_ADD_LIBS_DIR)/additionalLibs +FOAM_VERSION_MAJOR := $(word 1,$(subst ., ,$(WM_PROJECT_VERSION))) +PFLAGS+= -DOPENFOAM_VERSION_MAJOR=$(FOAM_VERSION_MAJOR) PFLAGS+= -Dcompre EXE_INC = \ diff --git a/applications/solvers/cfdemSolverRhoPimple/cfdemSolverRhoPimple.C b/applications/solvers/cfdemSolverRhoPimple/cfdemSolverRhoPimple.C index eb4d2c47..f7d3d217 100644 --- a/applications/solvers/cfdemSolverRhoPimple/cfdemSolverRhoPimple.C +++ b/applications/solvers/cfdemSolverRhoPimple/cfdemSolverRhoPimple.C @@ -92,7 +92,7 @@ int main(int argc, char *argv[]) particleCloud.clockM().start(2,"Coupling"); bool hasEvolved = particleCloud.evolve(voidfraction,Us,U); - if(hasEvolved) + if(hasEvolved && smoothenForces) { particleCloud.smoothingM().smoothen(particleCloud.forceM(0).impParticleForces()); } @@ -113,7 +113,11 @@ int main(int argc, char *argv[]) particleCloud.clockM().start(26,"Flow"); +#if OPENFOAM_VERSION_MAJOR < 6 if (pimple.nCorrPIMPLE() <= 1) +#else + if (pimple.nCorrPimple() <= 1) +#endif { #include "rhoEqn.H" } diff --git a/applications/solvers/cfdemSolverRhoPimple/createFields.H b/applications/solvers/cfdemSolverRhoPimple/createFields.H index 35364cbd..23ddbab4 100644 --- a/applications/solvers/cfdemSolverRhoPimple/createFields.H +++ b/applications/solvers/cfdemSolverRhoPimple/createFields.H @@ -177,6 +177,17 @@ Info<< "Reading thermophysical properties\n" << endl; ) ); + bool smoothenForces + ( + pimple.dict().lookupOrDefault + ( + "smoothenForces", + false + ) + ); + if (smoothenForces) Info << "Smoothening implicit particle forces.\n" << endl; + else Info << "Not smoothening implicit particle forces.\n" << endl; + Info<< "Creating turbulence model\n" << endl; autoPtr turbulence ( diff --git a/applications/solvers/cfdemSolverRhoPimpleChem/Make/options b/applications/solvers/cfdemSolverRhoPimpleChem/Make/options index aa2a2273..8f25c877 100644 --- a/applications/solvers/cfdemSolverRhoPimpleChem/Make/options +++ b/applications/solvers/cfdemSolverRhoPimpleChem/Make/options @@ -1,5 +1,7 @@ include $(CFDEM_ADD_LIBS_DIR)/additionalLibs +FOAM_VERSION_MAJOR := $(word 1,$(subst ., ,$(WM_PROJECT_VERSION))) +PFLAGS+= -DOPENFOAM_VERSION_MAJOR=$(FOAM_VERSION_MAJOR) PFLAGS+= -Dcompre EXE_INC = \ @@ -19,12 +21,8 @@ EXE_INC = \ -I$(LIB_SRC)/thermophysicalModels/specie/lnInclude \ -I$(LIB_SRC)/transportModels/compressible/lnInclude \ -I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \ - -I$(LIB_SRC)/thermophysicalModels/properties/liquidProperties/lnInclude \ - -I$(LIB_SRC)/thermophysicalModels/properties/liquidMixtureProperties/lnInclude \ - -I$(LIB_SRC)/thermophysicalModels/thermophysicalFunctions/lnInclude \ -I$(LIB_SRC)/thermophysicalModels/reactionThermo/lnInclude \ -I$(LIB_SRC)/thermophysicalModels/chemistryModel/lnInclude \ - -I$(LIB_SRC)/thermophysicalModels/radiationModels/lnInclude \ -I$(LIB_SRC)/regionModels/regionModel/lnInclude \ -I$(LIB_SRC)/regionModels/surfaceFilmModels/lnInclude \ -I$(LIB_SRC)/ODE/lnInclude \ @@ -48,9 +46,6 @@ EXE_LIBS = \ -l$(CFDEM_LIB_COMP_NAME) \ $(CFDEM_ADD_LIB_PATHS) \ $(CFDEM_ADD_LIBS) \ - -lliquidProperties \ - -lliquidMixtureProperties \ - -lthermophysicalFunctions \ -lreactionThermophysicalModels \ -lchemistryModel \ -lradiationModels \ diff --git a/applications/solvers/cfdemSolverRhoPimpleChem/YEqn.H b/applications/solvers/cfdemSolverRhoPimpleChem/YEqn.H index 3c82ddc8..e4d58344 100644 --- a/applications/solvers/cfdemSolverRhoPimpleChem/YEqn.H +++ b/applications/solvers/cfdemSolverRhoPimpleChem/YEqn.H @@ -13,7 +13,11 @@ tmp > mvConvection { combustion->correct(); +#if OPENFOAM_VERSION_MAJOR < 5 dQ = combustion->dQ(); +#else + Qdot = combustion->Qdot(); +#endif label inertIndex = -1; volScalarField Yt(0.0*Y[0]); @@ -72,4 +76,5 @@ tmp > mvConvection Y[inertIndex].max(0.0); } } + particleCloud.clockM().stop("Y"); diff --git a/applications/solvers/cfdemSolverRhoPimpleChem/cfdemSolverRhoPimpleChem.C b/applications/solvers/cfdemSolverRhoPimpleChem/cfdemSolverRhoPimpleChem.C index a634ef7a..c4da39ff 100644 --- a/applications/solvers/cfdemSolverRhoPimpleChem/cfdemSolverRhoPimpleChem.C +++ b/applications/solvers/cfdemSolverRhoPimpleChem/cfdemSolverRhoPimpleChem.C @@ -30,7 +30,12 @@ Description #include "fvCFD.H" #include "turbulentFluidThermoModel.H" +#if OPENFOAM_VERSION_MAJOR < 6 #include "rhoCombustionModel.H" +#else +#include "rhoReactionThermo.H" +#include "CombustionModel.H" +#endif #include "bound.H" #include "pimpleControl.H" #include "fvOptions.H" @@ -115,7 +120,11 @@ int main(int argc, char *argv[]) particleCloud.clockM().start(26,"Flow"); +#if OPENFOAM_VERSION_MAJOR < 6 if (pimple.nCorrPIMPLE() <= 1) +#else + if (pimple.nCorrPimple() <= 1) +#endif { #include "rhoEqn.H" } diff --git a/applications/solvers/cfdemSolverRhoPimpleChem/createFields.H b/applications/solvers/cfdemSolverRhoPimpleChem/createFields.H index 85624d06..4cda15ee 100644 --- a/applications/solvers/cfdemSolverRhoPimpleChem/createFields.H +++ b/applications/solvers/cfdemSolverRhoPimpleChem/createFields.H @@ -1,20 +1,28 @@ // thermodynamics, chemistry +#if OPENFOAM_VERSION_MAJOR < 6 Info<< "Creating combustion model\n" << endl; - autoPtr combustion ( combustionModels::rhoCombustionModel::New(mesh) ); - rhoReactionThermo& thermo = combustion->thermo(); +#else + Info<< "Reading thermophysical properties\n" << endl; + autoPtr pThermo(rhoReactionThermo::New(mesh)); + rhoReactionThermo& thermo = pThermo(); +#endif thermo.validate(args.executable(), "h", "e"); basicSpecieMixture& composition = thermo.composition(); PtrList& Y = composition.Y(); // read molecular weight +#if OPENFOAM_VERSION_MAJOR < 6 volScalarField W(composition.W()); +#else + volScalarField W(thermo.W()); +#endif bool propagateInertSpecie = true; @@ -198,6 +206,14 @@ ) ); +#if OPENFOAM_VERSION_MAJOR >= 6 + Info<< "Creating combustion model\n" << endl; + autoPtr> combustion + ( + CombustionModel::New(thermo, turbulence()) + ); +#endif + Info<< "Creating field dpdt\n" << endl; volScalarField dpdt ( @@ -214,6 +230,7 @@ Info<< "Creating field kinetic energy K\n" << endl; volScalarField K("K", 0.5*magSqr(U)); +#if OPENFOAM_VERSION_MAJOR < 5 volScalarField dQ ( IOobject @@ -227,6 +244,21 @@ mesh, dimensionedScalar("dQ", dimEnergy/dimTime, 0.0) ); +#else + volScalarField Qdot + ( + IOobject + ( + "Qdot", + runTime.timeName(), + mesh, + IOobject::NO_READ, + IOobject::AUTO_WRITE + ), + mesh, + dimensionedScalar("Qdot", dimEnergy/dimVolume/dimTime, 0.0) + ); +#endif Info<< "\nReading momentum exchange field Ksl\n" << endl; volScalarField Ksl diff --git a/applications/solvers/rStatAnalysis/Make/options b/applications/solvers/rStatAnalysis/Make/options index f91d7947..0d54afbc 100644 --- a/applications/solvers/rStatAnalysis/Make/options +++ b/applications/solvers/rStatAnalysis/Make/options @@ -23,4 +23,4 @@ EXE_LIBS = \ -lmeshTools \ -l$(CFDEM_LIB_NAME) \ $(CFDEM_ADD_LIB_PATHS) \ -$(CFDEM_ADD_LIBS) + $(CFDEM_ADD_LIBS) diff --git a/applications/solvers/rStatAnalysis/rStatAnalysis.C b/applications/solvers/rStatAnalysis/rStatAnalysis.C index 38356182..72523bfd 100644 --- a/applications/solvers/rStatAnalysis/rStatAnalysis.C +++ b/applications/solvers/rStatAnalysis/rStatAnalysis.C @@ -29,11 +29,6 @@ Description \*---------------------------------------------------------------------------*/ -// #include "fvCFD.H" -// #include "singlePhaseTransportModel.H" -// #include "turbulentTransportModel.H" -// #include "fvOptions.H" - #include "recBase.H" #include "recStatAnalysis.H" diff --git a/applications/solvers/rcfdemSolverHeattransfer/Make/files b/applications/solvers/rcfdemSolverHeattransfer/Make/files deleted file mode 100644 index cfbe9475..00000000 --- a/applications/solvers/rcfdemSolverHeattransfer/Make/files +++ /dev/null @@ -1,3 +0,0 @@ -rcfdemSolverHeattransfer.C - -EXE=$(CFDEM_APP_DIR)/rcfdemSolverHeattransfer diff --git a/applications/solvers/rcfdemSolverHeattransfer/Make/options b/applications/solvers/rcfdemSolverHeattransfer/Make/options deleted file mode 100644 index 9443c2f8..00000000 --- a/applications/solvers/rcfdemSolverHeattransfer/Make/options +++ /dev/null @@ -1,27 +0,0 @@ -include $(CFDEM_ADD_LIBS_DIR)/additionalLibs - -EXE_INC = \ - -I$(CFDEM_OFVERSION_DIR) \ - -I$(LIB_SRC)/finiteVolume/lnInclude \ - -I$(LIB_SRC)/meshTools/lnInclude \ - -I$(LIB_SRC)/TurbulenceModels/turbulenceModels/lnInclude \ - -I$(LIB_SRC)/TurbulenceModels/incompressible/lnInclude \ - -I$(LIB_SRC)/transportModels \ - -I$(LIB_SRC)/transportModels/incompressible/singlePhaseTransportModel \ - -I$(CFDEM_SRC_DIR)/lagrangian/cfdemParticle/lnInclude \ - -I$(CFDEM_SRC_DIR)/lagrangian/cfdemParticle/cfdTools \ - -I$(CFDEM_SRC_DIR)/recurrence/lnInclude \ - -I$(CFDEM_SRC_DIR)/lagrangian/cfdemParticle/derived/cfdemCloudRec \ - -EXE_LIBS = \ - -L$(CFDEM_LIB_DIR)\ - -lrecurrence \ - -lturbulenceModels \ - -lincompressibleTurbulenceModels \ - -lincompressibleTransportModels \ - -lfiniteVolume \ - -lmeshTools \ - -lfvOptions \ - -l$(CFDEM_LIB_NAME) \ - $(CFDEM_ADD_LIB_PATHS) \ -$(CFDEM_ADD_LIBS) diff --git a/applications/solvers/rcfdemSolverHeattransfer/TEqImp.H b/applications/solvers/rcfdemSolverHeattransfer/TEqImp.H deleted file mode 100644 index 104c3d87..00000000 --- a/applications/solvers/rcfdemSolverHeattransfer/TEqImp.H +++ /dev/null @@ -1,40 +0,0 @@ - volScalarField rhoeps = rhoRec*voidfractionRec; - - particleCloud.energyContributions(Qsource); - - particleCloud.energyCoefficients(QCoeff); - - K = 0.5*magSqr(URec); - - addSource = fvc::ddt(rhoeps, K) + fvc::div(phiRec, K) + - fvc::div - ( - fvc::absolute(phiRec/fvc::interpolate(rhoRec), voidfractionRec*URec), pRec - ); - - - fvScalarMatrix TEqn = - ( - fvm::ddt(rhoeps, T) - + fvm::div(phiRec, T) - + addSource/Cv - - fvm::laplacian(voidfractionRec*thCond/Cv, T) - - Qsource/Cv - - fvm::Sp(QCoeff/Cv, T) - == - fvOptions(rhoeps, T) // no fvOptions support yet - ); - - //TEqn.relax(relaxCoeff); - - fvOptions.constrain(TEqn); // no fvOptions support yet - - TEqn.solve(); - - particleCloud.clockM().start(31,"postFlow"); - counter++; - - if((counter - couplingSubStep) % dtDEM2dtCFD == 0) - particleCloud.postFlow(); - - particleCloud.clockM().stop("postFlow"); diff --git a/applications/solvers/rcfdemSolverHeattransfer/createFields.H b/applications/solvers/rcfdemSolverHeattransfer/createFields.H deleted file mode 100644 index cade4a1a..00000000 --- a/applications/solvers/rcfdemSolverHeattransfer/createFields.H +++ /dev/null @@ -1,206 +0,0 @@ - // dummy fields - Info << "\nCreating dummy pressure field\n" << endl; - volScalarField pRec - ( - IOobject - ( - "pRec", - runTime.timeName(), - mesh, - IOobject::MUST_READ, - IOobject::AUTO_WRITE - ), - mesh - ); - - // recurrence fields - Info << "\nCreating recurrence fields.\n" << endl; - - volScalarField rhoRec - ( - IOobject - ( - "rhoRec", - runTime.timeName(), - mesh, - IOobject::MUST_READ, - IOobject::AUTO_WRITE - ), - mesh//, - //dimensionedScalar("rhoRec", dimensionSet(1, -3, 0, 0, 0), 1.0) - ); - - volVectorField URec - ( - IOobject - ( - "URec", - runTime.timeName(), - mesh, - IOobject::MUST_READ, - IOobject::AUTO_WRITE - ), - mesh - ); - - volScalarField voidfractionRec - ( - IOobject - ( - "voidfractionRec", - runTime.timeName(), - mesh, - IOobject::MUST_READ, - IOobject::AUTO_WRITE - ), - mesh - ); - - volVectorField UsRec - ( - IOobject - ( - "UsRec", - runTime.timeName(), - mesh, - IOobject::MUST_READ, - IOobject::AUTO_WRITE - ), - mesh - ); - - // heat transfer fields - Info << "\nCreating heat transfer fields.\n" << endl; - - volScalarField Qsource - ( - IOobject - ( - "Qsource", - runTime.timeName(), - mesh, - IOobject::NO_READ, - IOobject::AUTO_WRITE - ), - mesh, - dimensionedScalar("zero", dimensionSet(1,-1,-3,0,0,0,0), 0.0) - ); - - volScalarField QCoeff - ( - IOobject - ( - "QCoeff", - runTime.timeName(), - mesh, - IOobject::NO_READ, - IOobject::AUTO_WRITE - ), - mesh, - dimensionedScalar("zero", dimensionSet(1,-1,-3,-1,0,0,0), 0.0) - ); - - volScalarField thCond - ( - IOobject - ( - "thCond", - runTime.timeName(), - mesh, - IOobject::READ_IF_PRESENT, - IOobject::AUTO_WRITE - ), - mesh, - dimensionedScalar("zero", dimensionSet(1,1,-3,-1,0,0,0), 0.0), - "zeroGradient" - ); - - volScalarField T - ( - IOobject - ( - "T", - runTime.timeName(), - mesh, - IOobject::MUST_READ, - IOobject::AUTO_WRITE - ), - mesh - ); - - // calculated fields - Info << "\nCreating fields subject to calculation\n" << endl; - volScalarField voidfraction - ( - IOobject - ( - "voidfraction", - runTime.timeName(), - mesh, - IOobject::NO_READ, - IOobject::AUTO_WRITE - ), - voidfractionRec - ); - - volVectorField Us - ( - IOobject - ( - "Us", - runTime.timeName(), - mesh, - IOobject::NO_READ, - IOobject::AUTO_WRITE - ), - UsRec - ); - - // write fields for t=t_start - voidfraction.write(); - Us.write(); -//=============================== - - - Info << "Calculating face flux field phiRec\n" << endl; - surfaceScalarField phiRec - ( - IOobject - ( - "phiRec", - runTime.timeName(), - mesh, - IOobject::READ_IF_PRESENT, - IOobject::AUTO_WRITE - ), - linearInterpolate(URec*voidfractionRec*rhoRec) & mesh.Sf() - ); - phiRec.write(); - - - singlePhaseTransportModel laminarTransport(URec, phiRec); - - autoPtr turbulence - ( - incompressible::turbulenceModel::New(URec, phiRec, laminarTransport) - ); - - const IOdictionary& transportProps = mesh.lookupObject("transportProperties"); - dimensionedScalar Cv(transportProps.lookup("Cv")); - - volScalarField addSource - ( - IOobject - ( - "addSource", - runTime.timeName(), - mesh, - IOobject::NO_READ, - IOobject::AUTO_WRITE - ), - mesh, - dimensionedScalar("zero", dimensionSet(1,-1,-3,0,0,0,0), 0.0) - ); - - Info << "Creating field kinetic energy K\n" << endl; - volScalarField K("K", 0.5*magSqr(URec)); diff --git a/applications/solvers/rcfdemSolverHeattransfer/rcfdemSolverHeattransfer.C b/applications/solvers/rcfdemSolverHeattransfer/rcfdemSolverHeattransfer.C deleted file mode 100644 index 7fa36bff..00000000 --- a/applications/solvers/rcfdemSolverHeattransfer/rcfdemSolverHeattransfer.C +++ /dev/null @@ -1,128 +0,0 @@ -/*---------------------------------------------------------------------------*\ - CFDEMcoupling academic - Open Source CFD-DEM coupling - - Contributing authors: - Thomas Lichtenegger - Copyright (C) 2015- Johannes Kepler University, Linz -------------------------------------------------------------------------------- -License - This file is part of CFDEMcoupling academic. - - CFDEMcoupling academic 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 academic 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 academic. If not, see . - -Application - rcfdemSolverHeattransfer - -Description - Solves heat transfer between fluid and particles based on rCFD - -\*---------------------------------------------------------------------------*/ - -#include "fvCFD.H" -#include "singlePhaseTransportModel.H" -#include "turbulentTransportModel.H" -#include "fvOptions.H" - -#include "cfdemCloudRec.H" -#include "recBase.H" -#include "recModel.H" - -#include "cfdemCloudEnergy.H" -#include "clockModel.H" -#include "thermCondModel.H" -#include "energyModel.H" - - - -// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // - -int main(int argc, char *argv[]) -{ - #include "postProcess.H" - #include "setRootCase.H" - #include "createTime.H" - #include "createMesh.H" - #include "createControl.H" - #include "createFields.H" - #include "createFvOptions.H" - - cfdemCloudRec particleCloud(mesh); - recBase recurrenceBase(mesh); - - // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // - - Info << "\nCalculating particle trajectories based on recurrence statistics\n" << endl; - - label recTimeIndex = 0; - scalar recTimeStep = recurrenceBase.recM().recTimeStep(); - scalar startTime = runTime.startTime().value(); - - // control coupling behavior in case of substepping - // assumes constant timestep size - label counter = 0; - label couplingSubStep = recurrenceBase.couplingSubStep();//5;//3; - double dtProp = particleCloud.dataExchangeM().couplingTime() / runTime.deltaTValue(); - label dtDEM2dtCFD = int(dtProp + 0.5); - Info << "deltaT_DEM / deltaT_CFD = " << dtDEM2dtCFD << endl; - if (dtDEM2dtCFD > 1) - Info << "coupling at substep " << couplingSubStep << endl; - - - while (runTime.run()) - { - runTime++; - - // do stuff (every lagrangian time step) - particleCloud.clockM().start(1,"Global"); - - Info << "Time = " << runTime.timeName() << nl << endl; - - particleCloud.clockM().start(2,"Coupling"); - - particleCloud.evolve(voidfraction,Us,URec); - - particleCloud.clockM().stop("Coupling"); - - particleCloud.clockM().start(26,"Flow"); - #include "TEqImp.H" - particleCloud.clockM().stop("Flow"); - - particleCloud.clockM().start(32,"ReadFields"); - if ( runTime.timeOutputValue() - startTime - (recTimeIndex+1)*recTimeStep + 1.0e-5 > 0.0 ) - { - recurrenceBase.updateRecFields(); - #include "updateFields.H" - recTimeIndex++; - } - particleCloud.clockM().stop("ReadFields"); - - particleCloud.clockM().start(33,"Output"); - runTime.write(); - particleCloud.clockM().stop("Output"); - - particleCloud.clockM().stop("Global"); - - Info << "ExecutionTime = " << runTime.elapsedCpuTime() << " s" - << " ClockTime = " << runTime.elapsedClockTime() << " s" - << nl << endl; - - } - - Info << "End\n" << endl; - - return 0; -} - - -// ************************************************************************* // diff --git a/applications/solvers/rcfdemSolverHeattransfer/updateFields.H b/applications/solvers/rcfdemSolverHeattransfer/updateFields.H deleted file mode 100644 index 5e9a8ab1..00000000 --- a/applications/solvers/rcfdemSolverHeattransfer/updateFields.H +++ /dev/null @@ -1,14 +0,0 @@ -recurrenceBase.recM().exportVolScalarField("voidfraction",voidfractionRec); -recurrenceBase.recM().exportVolScalarField("rho",rhoRec); -recurrenceBase.recM().exportVolScalarField("p",pRec); -recurrenceBase.recM().exportVolVectorField("Us",UsRec); -recurrenceBase.recM().exportVolVectorField("U",URec); -recurrenceBase.recM().exportSurfaceScalarField("phi",phiRec); - -{ - volScalarField& NuField(const_cast(mesh.lookupObject ("NuField"))); - recurrenceBase.recM().exportVolScalarField("NuField",NuField); -} - - -#include "updateRho.H" \ No newline at end of file diff --git a/applications/solvers/rcfdemSolverHeattransfer/updateRho.H b/applications/solvers/rcfdemSolverHeattransfer/updateRho.H deleted file mode 100644 index 159a8654..00000000 --- a/applications/solvers/rcfdemSolverHeattransfer/updateRho.H +++ /dev/null @@ -1,32 +0,0 @@ -// work-around for transient properties -// needs to be specified for each case - -// case 1 - -forAll(rhoRec,cellI) -{ - if (mesh.C()[cellI].z() < 0.00228) - rhoRec[cellI] = 1.18+(1.085-1.18)*Foam::exp(-0.065*runTime.timeOutputValue()); - else if (mesh.C()[cellI].z() < 0.00456) - rhoRec[cellI] = 1.18+(1.01-1.18)*Foam::exp(-0.05*runTime.timeOutputValue()); - else if (mesh.C()[cellI].z() < 0.00684) - rhoRec[cellI] = 1.18+(0.98-1.18)*Foam::exp(-0.0425*runTime.timeOutputValue()); - else - rhoRec[cellI] = 1.18+(0.955-1.18)*Foam::exp(-0.0425*runTime.timeOutputValue()); -} - - -// case 2 -/* -forAll(rhoRec,cellI) -{ - if (mesh.C()[cellI].z() < 0.00228) - rhoRec[cellI] = 1.18+(1.115-1.18)*Foam::exp(-0.065*runTime.timeOutputValue()); - else if (mesh.C()[cellI].z() < 0.00456) - rhoRec[cellI] = 1.18+(1.04-1.18)*Foam::exp(-0.05*runTime.timeOutputValue()); - else if (mesh.C()[cellI].z() < 0.00684) - rhoRec[cellI] = 1.18+(1.005-1.18)*Foam::exp(-0.0425*runTime.timeOutputValue()); - else - rhoRec[cellI] = 1.18+(0.96-1.18)*Foam::exp(-0.0425*runTime.timeOutputValue()); -} -*/ \ No newline at end of file diff --git a/applications/solvers/rcfdemSolverRhoSteadyPimple/Make/options b/applications/solvers/rcfdemSolverRhoSteadyPimple/Make/options index b2b37fe7..a97830ba 100644 --- a/applications/solvers/rcfdemSolverRhoSteadyPimple/Make/options +++ b/applications/solvers/rcfdemSolverRhoSteadyPimple/Make/options @@ -1,5 +1,7 @@ include $(CFDEM_ADD_LIBS_DIR)/additionalLibs +FOAM_VERSION_MAJOR := $(word 1,$(subst ., ,$(WM_PROJECT_VERSION))) +PFLAGS+= -DOPENFOAM_VERSION_MAJOR=$(FOAM_VERSION_MAJOR) PFLAGS+= -Dcompre EXE_INC = \ diff --git a/applications/solvers/rcfdemSolverRhoSteadyPimple/rcfdemSolverRhoSteadyPimple.C b/applications/solvers/rcfdemSolverRhoSteadyPimple/rcfdemSolverRhoSteadyPimple.C index 6d2fe6e1..5eaeed0b 100644 --- a/applications/solvers/rcfdemSolverRhoSteadyPimple/rcfdemSolverRhoSteadyPimple.C +++ b/applications/solvers/rcfdemSolverRhoSteadyPimple/rcfdemSolverRhoSteadyPimple.C @@ -69,15 +69,15 @@ int main(int argc, char *argv[]) #include "createFvOptions.H" // create cfdemCloud - // #include "readGravitationalAcceleration.H" + //#include "readGravitationalAcceleration.H" cfdemCloudRec particleCloud(mesh); #include "checkModelType.H" recBase recurrenceBase(mesh); #include "updateFields.H" turbulence->validate(); - // #include "compressibleCourantNo.H" - // #include "setInitialDeltaT.H" + //#include "compressibleCourantNo.H" + //#include "setInitialDeltaT.H" // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // @@ -108,8 +108,8 @@ int main(int argc, char *argv[]) particleCloud.clockM().start(2,"Coupling"); bool hasEvolved = particleCloud.evolve(voidfraction,Us,U); -//voidfraction = voidfractionRec; -//Us = UsRec; + //voidfraction = voidfractionRec; + //Us = UsRec; if(hasEvolved) { @@ -137,7 +137,11 @@ int main(int argc, char *argv[]) while (pimple.loop()) { // if needed, perform drag update here +#if OPENFOAM_VERSION_MAJOR < 6 if (pimple.nCorrPIMPLE() <= 1) +#else + if (pimple.nCorrPimple() <= 1) +#endif { #include "rhoEqn.H" } diff --git a/doc/CFDEMcoupling_install.txt b/doc/CFDEMcoupling_install.txt index b4634e3e..8c599862 100644 --- a/doc/CFDEMcoupling_install.txt +++ b/doc/CFDEMcoupling_install.txt @@ -20,7 +20,7 @@ In the CFDEMcoupling repository take a look at the file src/lagrangian/cfdemParticle/cfdTools/versionInfo.H :pre to find out the latest tested version of LIGGGHTS and OpenFOAM that work with -CFDEMcoupling. As of this writing the version of OpenFOAM to be used is 4.x. +CFDEMcoupling. As of this writing the version of OpenFOAM to be used is 6. You can then basically follow the instructions at "openfoam.org"_https://openfoam.org/download/source/, cloning OpenFOAM from the @@ -29,11 +29,11 @@ git repository. cd $HOME mkdir OpenFOAM cd OpenFOAM -git clone https://github.com/OpenFOAM/OpenFOAM-4.x.git :pre +git clone https://github.com/OpenFOAM/OpenFOAM-6.git :pre Clone the corresponding third party packages to the OpenFOAM folder. -git clone https://github.com/OpenFOAM/ThirdParty-4.x.git :pre +git clone https://github.com/OpenFOAM/ThirdParty-6.git :pre Switch to root user with sudo @@ -42,9 +42,8 @@ sudo su - :pre Install dependent packages required for OpenFOAM on Ubuntu by executing the following commands: -apt-get install build-essential flex bison cmake zlib1g-dev libboost-system-dev libboost-thread-dev libopenmpi-dev openmpi-bin gnuplot libreadline-dev libncurses-dev libxt-dev -apt-get install qt4-dev-tools libqt4-dev libqt4-opengl-dev freeglut3-dev libqtwebkit-dev -apt-get install libcgal-dev :pre +apt-get install build-essential flex bison git-core cmake zlib1g-dev libboost-system-dev libboost-thread-dev libopenmpi-dev openmpi-bin gnuplot libreadline-dev libncurses-dev libxt-dev +apt-get install libqt5x11extras5-dev libxt-dev qt5-default qttools5-dev curl :pre 2.1.2 Setup the environment :h5 @@ -56,7 +55,7 @@ gedit ~/.bashrc :pre and add the following lines: -source $HOME/OpenFOAM/OpenFOAM-4.x/etc/bashrc +source $HOME/OpenFOAM/OpenFOAM-6/etc/bashrc export WM_NCOMPPROCS=4 :pre Save the file and reload it: @@ -67,7 +66,7 @@ source ~/.bashrc :pre [Additional check] -Open ~/OpenFOAM/OpenFOAM-4.x/etc/bashrc and make sure that {WM_MPLIB} is set +Open ~/OpenFOAM/OpenFOAM-6/etc/bashrc and make sure that {WM_MPLIB} is set correctly: export WM_MPLIB=SYSTEMOPENMPI :pre @@ -89,7 +88,7 @@ cd $WM_THIRD_PARTY_DIR Paraview is a third-party software provided for graphical post-processing in OpenFOAM. Its compilation is automated using a script called makeParaView in the -ThirdParty-4.x directory. +ThirdParty-6 directory. Before installing Paraview, check the version of cmake that is installed on the system. This can be done by typing diff --git a/doc/CFDEMcoupling_models.txt b/doc/CFDEMcoupling_models.txt index 800e9eb4..2ce1c291 100644 --- a/doc/CFDEMcoupling_models.txt +++ b/doc/CFDEMcoupling_models.txt @@ -27,7 +27,8 @@ The "averagingModel"_averagingModel.html keyword entry defines the model used to map the Lagrangian data to Eulerian values. "dense"_averagingModel_dense.html, -"dilute"_averagingModel_dilute.html :tb(c=2,ea=c) +"dilute"_averagingModel_dilute.html, +off :tb(c=2,ea=c) 6.3 Chemistry models :h4 @@ -69,8 +70,9 @@ that performs the data exchange between the DEM code and the CFD code. The {energyModels} keyword specifies a list of energy models used for e.g. compressible, reacting flows. +heatTransferGranConduction, heatTransferGunn, -heatTransferGunnPartField, +heatTransferRanzMarshall, reactionHeat :tb(c=2,ea=c) @@ -92,15 +94,23 @@ Fines, "MeiLift"_forceModel_MeiLift.html, "SchillerNaumannDrag"_forceModel_SchillerNaumannDrag.html, "ShirgaonkarIB"_forceModel_ShirgaonkarIB.html, +dSauter, +deactivateForce, +directedDiffusiveRelaxation, +evaluateFluctuations, "fieldStore"_forceModel_fieldStore.html, "fieldTimeAverage"_forceModel_fieldTimeAverage.html, +freeStreaming, "gradPForce"_forceModel_gradPForce.html, "gradPForceSmooth"_forceModel_gradPForceSmooth.html, granKineticEnergy, "interface"_forceModel_interface.html, +isotropicFluctuations, "noDrag"_forceModel_noDrag.html, "particleCellVolume"_forceModel_particleCellVolume.html, +particleDeformation, "pdCorrelation"_forceModel_pdCorrelation.html, +potentialRelaxation, "surfaceTensionForce"_forceModel_surfaceTensionForce.html, "virtualMassForce"_forceModel_virtualMassForce.html, "viscForce"_forceModel_viscForce.html, @@ -153,6 +163,7 @@ used to manipulate the CFD mesh according to the DEM mesh motion. The "momCoupleModels"_momCoupleModel.html keyword specifies a list of models used for momentum exchange between DEM and CFD simulation +deactivateCouple, "explicitCouple"_momCoupleModel_explicitCouple.html, "implicitCouple"_momCoupleModel_implicitCouple.html, "off"_momCoupleModel_noCouple.html :tb(c=2,ea=c) @@ -214,5 +225,7 @@ accounting for the volume of the particles in the CFD domain. "IB"_voidFractionModel_IBVoidFraction.html, "bigParticle"_voidFractionModel_bigParticleVoidFraction.html, "centre"_voidFractionModel_centreVoidFraction.html, -"divided"_voidFractionModel_dividedVoidFraction.html :tb(c=2,ea=c) +"divided"_voidFractionModel_dividedVoidFraction.html, +off, +trilinear :tb(c=2,ea=c) diff --git a/doc/CFDEMcoupling_solvers.txt b/doc/CFDEMcoupling_solvers.txt index acbca05a..7fad7b71 100644 --- a/doc/CFDEMcoupling_solvers.txt +++ b/doc/CFDEMcoupling_solvers.txt @@ -15,5 +15,11 @@ This section lists all CFDEMcoupling solvers alphabetically. "cfdemSolverPisoScalar"_cfdemSolverPisoScalar.html, "cfdemSolverRhoPimple"_cfdemSolverRhoPimple.html, "cfdemSolverRhoPimpleChem"_cfdemSolverRhoPimpleChem.html, -"cfdemSolverRhoSimple"_cfdemSolverRhoSimple.html :tb(c=2,ea=c) +rStatAnalysis, +rcfdemSolverBase, +rcfdemSolverCoupledHeattransfer, +rcfdemSolverRhoSteadyPimple, +recSolverTurbTransport, +rtfmSolverSpecies, +testTwoFluidRecurrenceTurbulence :tb(c=2,ea=c) diff --git a/doc/cfdemSolverMultiphase.txt b/doc/cfdemSolverMultiphase.txt index 6f7a08a3..10be3bbf 100644 --- a/doc/cfdemSolverMultiphase.txt +++ b/doc/cfdemSolverMultiphase.txt @@ -29,6 +29,8 @@ END_RST --> For more details, see "Vångö et al. (2018)"_#Vångö2018. +IMPORTANT NOTE: This solver requires OpenFOAM 4.x or 5.x to work properly. + :line :link(Vångö2018) diff --git a/doc/cfdemSolverRhoSimple.txt b/doc/cfdemSolverRhoSimple.txt deleted file mode 100644 index d0ac3e00..00000000 --- a/doc/cfdemSolverRhoSimple.txt +++ /dev/null @@ -1,59 +0,0 @@ - - - - - -"CFDEMproject Website"_lws - "Main Page"_main :c - -:link(lws,http://www.cfdem.com) -:link(main,CFDEMcoupling_Manual.html) - -:line - -cfdemSolverRhoSimple command :h3 - -[Description:] - - -"cfdemSolverRhoSimple" is a coupled CFD-DEM solver using the CFDEMcoupling -framework. Based on the OpenFOAM®(*) solver rhoSimpleFoam, this is a -steady-state solver for turbulent flow of compressible fluids coupled with the -DEM code LIGGGHTS for solid particles. - - - - - -:line - - -NOTE: -(*) This offering is not approved or endorsed by OpenCFD Limited, producer and -distributor of the OpenFOAM software via www.openfoam.com, and owner of the -OPENFOAM® and OpenCFD® trade marks. -OPENFOAM® is a registered trade mark of OpenCFD Limited, producer and -distributor of the OpenFOAM software via www.openfoam.com. - - - diff --git a/doc/forceModel_particleCellVolume.txt b/doc/forceModel_particleCellVolume.txt index 08162efb..e8065699 100644 --- a/doc/forceModel_particleCellVolume.txt +++ b/doc/forceModel_particleCellVolume.txt @@ -22,11 +22,14 @@ particleCellVolumeProps lowerThreshold scalar2; startTime scalar3; verbose; + writeToFile switch1; \} :pre {scalar1} = only cells with a field value (magnitude) lower than this upper threshold are considered :l {scalar2} = only cells with a field value (magnitude) greater than this lower threshold are considered :l {scalar3} = (optional, default 0) start time of volume calculation and output :l +{verbose} = (optional, default false) keyword only (mostly used for debugging) :l +{switch1} = (optional, default false) switch for file output :l :ule [Examples:] diff --git a/doc/forceModel_volWeightedAverage.txt b/doc/forceModel_volWeightedAverage.txt index a5773317..0570bc64 100644 --- a/doc/forceModel_volWeightedAverage.txt +++ b/doc/forceModel_volWeightedAverage.txt @@ -30,6 +30,7 @@ volWeightedAverageProps upperThreshold scalar1; lowerThreshold scalar2; verbose; + writeToFile switch1; \} :pre {time} = (optional, default 0.) time to start the averaging :ulb,l @@ -38,6 +39,7 @@ volWeightedAverageProps {scalar1} = only cells with a field value (magnitude) lower than this upper threshold are considered :l {scalar2} = only cells with a field value (magnitude) greater than this lower threshold are considered :l {verbose} = (optional, default false) keyword only (mostly used for debugging) :l +{switch1} = (optional, default false) switch for file output :l :ule [Examples:] diff --git a/etc/bashrc b/etc/bashrc index 97057592..9c425502 100755 --- a/etc/bashrc +++ b/etc/bashrc @@ -17,7 +17,7 @@ #------------------------------------------------------------------------------ export CFDEM_PROJECT=CFDEM -export CFDEM_VERSION=19.09 +export CFDEM_VERSION=20.05 ################################################################################ # USER EDITABLE PART: Changes made here may be lost with the next upgrade @@ -25,7 +25,7 @@ export CFDEM_VERSION=19.09 # Please set to the appropriate path if the default is not correct. # # activate compatible OpenFOAM version -. $HOME/OpenFOAM/OpenFOAM-4.x/etc/bashrc +. $HOME/OpenFOAM/OpenFOAM-6/etc/bashrc # # CFDEMcoupling export CFDEM_INST_DIR=$HOME/$CFDEM_PROJECT @@ -199,6 +199,9 @@ alias cfdemCleanCFDEM='bash $CFDEM_PROJECT_DIR/etc/cleanCFDEMcoupling.sh' #- shortcut to compile LIGGGHTS + sublibraries alias cfdemCompLIG='bash $CFDEM_PROJECT_DIR/etc/compileLIGGGHTS.sh' +#- shortcut to compile LIGGGHTS sublibraries +alias cfdemCompLIGlib='bash $CFDEM_PROJECT_DIR/etc/compileLIGGGHTS_lib.sh' + #- shortcut to compile CFDEMcoupling +LIGGGHTS alias cfdemCompCFDEMall='bash $CFDEM_PROJECT_DIR/etc/compileCFDEMcoupling_all.sh' diff --git a/etc/cshrc b/etc/cshrc index feb0022c..e6d73425 100755 --- a/etc/cshrc +++ b/etc/cshrc @@ -15,7 +15,7 @@ #------------------------------------------------------------------------------ setenv CFDEM_PROJECT CFDEM -setenv CFDEM_VERSION 19.09 +setenv CFDEM_VERSION 20.05 ################################################################################ # USER EDITABLE PART: Changes made here may be lost with the next upgrade @@ -23,7 +23,7 @@ setenv CFDEM_VERSION 19.09 # Please set to the appropriate path if the default is not correct. # # activate compatible OpenFOAM version -. $HOME/OpenFOAM/OpenFOAM-4.x/etc/cshrc +. $HOME/OpenFOAM/OpenFOAM-6/etc/cshrc # # CFDEMcoupling setenv CFDEM_INST_DIR $HOME/$CFDEM_PROJECT @@ -232,6 +232,9 @@ alias cfdemCleanCFDEM 'bash $CFDEM_PROJECT_DIR/etc/cleanCFDEMcoupling.sh' #- shortcut to compile LIGGGHTS + sublibraries alias cfdemCompLIG 'bash $CFDEM_PROJECT_DIR/etc/compileLIGGGHTS.sh' +#- shortcut to compile LIGGGHTS sublibraries +alias cfdemCompLIGlib 'bash $CFDEM_PROJECT_DIR/etc/compileLIGGGHTS_lib.sh' + #- shortcut to compile CFDEMcoupling +LIGGGHTS alias cfdemCompCFDEMall 'bash $CFDEM_PROJECT_DIR/etc/compileCFDEMcoupling_all.sh' diff --git a/etc/solver-list.txt b/etc/solver-list.txt index 5eeca1d2..4c794769 100644 --- a/etc/solver-list.txt +++ b/etc/solver-list.txt @@ -1,5 +1,4 @@ rcfdemSolverRhoSteadyPimple/dir -rcfdemSolverHeattransfer/dir rcfdemSolverCoupledHeattransfer/dir rStatAnalysis/dir rcfdemSolverBase/dir @@ -7,7 +6,6 @@ rtfmSolverSpecies/dir cfdemSolverPisoMS/dir cfdemSolverPiso/dir cfdemSolverRhoPimple/dir -cfdemSolverRhoSimple/dir cfdemSolverIB/dir cfdemSolverPisoScalar/dir cfdemSolverRhoPimpleChem/dir diff --git a/src/lagrangian/cfdemParticle/Make/options b/src/lagrangian/cfdemParticle/Make/options index 8a84bb35..f8a3be9b 100644 --- a/src/lagrangian/cfdemParticle/Make/options +++ b/src/lagrangian/cfdemParticle/Make/options @@ -3,6 +3,8 @@ sinclude $(RULES)/mplib$(WM_MPLIB) GIT_VERSION := $(shell git describe --dirty --always --tags) PFLAGS+= -DGITVERSION=\"$(GIT_VERSION)\" +FOAM_VERSION_MAJOR := $(word 1,$(subst ., ,$(WM_PROJECT_VERSION))) +PFLAGS+= -DOPENFOAM_VERSION_MAJOR=$(FOAM_VERSION_MAJOR) include $(CFDEM_ADD_LIBS_DIR)/additionalLibs diff --git a/src/lagrangian/cfdemParticle/cfdTools/versionInfo.H b/src/lagrangian/cfdemParticle/cfdTools/versionInfo.H index 7fb67161..64c92d4e 100755 --- a/src/lagrangian/cfdemParticle/cfdTools/versionInfo.H +++ b/src/lagrangian/cfdemParticle/cfdTools/versionInfo.H @@ -34,9 +34,9 @@ Description #ifndef versionInfo_H #define versionInfo_H -word CFDEMversion="PFM 19.09"; -word compatibleLIGGGHTSversion="PFM 19.09"; -word OFversion="4.x"; +word CFDEMversion="PFM 20.05"; +word compatibleLIGGGHTSversion="PFM 20.05"; +word OFversion="6"; Info << "\nCFDEMcoupling version: " << CFDEMversion << endl; Info << "compatible to LIGGGHTS version: " << compatibleLIGGGHTSversion << endl; diff --git a/src/lagrangian/cfdemParticle/cfdemCloud/cfdemCloud.C b/src/lagrangian/cfdemParticle/cfdemCloud/cfdemCloud.C index 754f2b23..df719b60 100644 --- a/src/lagrangian/cfdemParticle/cfdemCloud/cfdemCloud.C +++ b/src/lagrangian/cfdemParticle/cfdemCloud/cfdemCloud.C @@ -84,10 +84,12 @@ cfdemCloud::cfdemCloud ignore_(couplingProperties_.found("ignore")), allowCFDsubTimestep_(true), limitDEMForces_(couplingProperties_.found("limitDEMForces")), + phaseInForces_(couplingProperties_.found("phaseInForcesTime")), getParticleDensities_(couplingProperties_.lookupOrDefault("getParticleDensities",false)), getParticleEffVolFactors_(couplingProperties_.lookupOrDefault("getParticleEffVolFactors",false)), getParticleTypes_(couplingProperties_.lookupOrDefault("getParticleTypes",false)), maxDEMForce_(0.), + phaseInForcesTime_(couplingProperties_.lookupOrDefault("phaseInForcesTime",0.0)), modelType_(couplingProperties_.lookup("modelType")), positions_(NULL), velocities_(NULL), @@ -368,6 +370,12 @@ cfdemCloud::cfdemCloud if (verbose_) Info << "nPatchesNonCyclic=" << nPatchesNonCyclic << ", nPatchesCyclic=" << nPatchesCyclic << endl; Warning << "Periodic handing is disabled because the domain is not fully periodic!\n" << endl; } + + // check if phasing-in time existing and is meaningful + if (phaseInForces_ && phaseInForcesTime_ < SMALL) + { + FatalError << "phasing-in time too small" << endl; + } } // * * * * * * * * * * * * * * * * Destructors * * * * * * * * * * * * * * // @@ -468,6 +476,19 @@ void cfdemCloud::setForces() } Info << "largest particle-fluid interaction on particle: " << maxF << endl; } + + if (phaseInForces_) + { + scalar tfrac = (mesh_.time().timeOutputValue()-mesh_.time().startTime().value())/phaseInForcesTime_; + if (tfrac <= 1.0) + { + for (int index = 0;index < numberOfParticles(); ++index) + { + for(int i=0;i<3;i++) DEMForces_[index][i] *= tfrac; + Cds_[index][0] *= tfrac; + } + } + } } void cfdemCloud::setParticleForceField() diff --git a/src/lagrangian/cfdemParticle/cfdemCloud/cfdemCloud.H b/src/lagrangian/cfdemParticle/cfdemCloud/cfdemCloud.H index 94a1e779..0b0d02ee 100644 --- a/src/lagrangian/cfdemParticle/cfdemCloud/cfdemCloud.H +++ b/src/lagrangian/cfdemParticle/cfdemCloud/cfdemCloud.H @@ -96,6 +96,8 @@ protected: const bool limitDEMForces_; + const bool phaseInForces_; + const bool getParticleDensities_; const bool getParticleEffVolFactors_; @@ -104,6 +106,8 @@ protected: scalar maxDEMForce_; + scalar phaseInForcesTime_; + const word modelType_; double **positions_; // particle positions diff --git a/src/lagrangian/cfdemParticle/subModels/IOModel/IOModel/IOModel.C b/src/lagrangian/cfdemParticle/subModels/IOModel/IOModel/IOModel.C index 3a773c16..ce63280d 100644 --- a/src/lagrangian/cfdemParticle/subModels/IOModel/IOModel/IOModel.C +++ b/src/lagrangian/cfdemParticle/subModels/IOModel/IOModel/IOModel.C @@ -30,6 +30,7 @@ Description \*---------------------------------------------------------------------------*/ #include "error.H" +#include "particle.H" #include "IOModel.H" // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // @@ -99,13 +100,22 @@ void IOModel::streamDataToPath(const fileName& path, const double* const* array, OFstream fileStream(path/name); fileStream - << "FoamFile\n{\n" - << " version " << fileStream.version() << ";\n" - << " format " << fileStream.format() << ";\n" - << " class " << className << ";\n" - << " location " << 0 << ";\n" - << " object " << name << ";\n" - << "}" << nl; + << "/*--------------------------------*- C++ -*----------------------------------*\\" << nl + << " ========= |" << nl + << " \\\\ / F ield | OpenFOAM: The Open Source CFD Toolbox" << nl + << " \\\\ / O peration | Website: https://openfoam.org" << nl + << " \\\\ / A nd | Version: " << FOAMversion << nl + << " \\\\/ M anipulation |" << nl + << "\\*---------------------------------------------------------------------------*/" << nl + << "FoamFile" << nl + << "{" << nl + << " version " << fileStream.version() << ";" << nl + << " format " << fileStream.format() << ";" << nl + << " class " << className << ";" << nl + << " location " << 0 << ";" << nl + << " object " << name << ";" << nl + << "}" << nl + << "// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //" << nl << nl; fileStream << nPProc << "\n"; @@ -129,7 +139,19 @@ void IOModel::streamDataToPath(const fileName& path, const double* const* array, } else if (type == "position") { - fileStream << "( "<< array[index][0] << " " << array[index][1] << " " << array[index][2] << " ) " << cellIDs[index][0] << " \n"; +#if OPENFOAM_VERSION_MAJOR < 5 + fileStream << "( " << array[index][0] << " " << array[index][1] << " " << array[index][2] << " ) " << cellIDs[index][0] << nl; +#else + particle part + ( + particleCloud_.mesh(), + vector(array[index][0],array[index][1],array[index][2]), + cellIDs[index][0] + ); + + part.writePosition(fileStream); + fileStream << nl; +#endif } else if (type == "vector") { diff --git a/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMPI/twoWayMPI.C b/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMPI/twoWayMPI.C index 6f7d0677..45af8a5c 100644 --- a/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMPI/twoWayMPI.C +++ b/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMPI/twoWayMPI.C @@ -68,6 +68,18 @@ twoWayMPI::twoWayMPI : dataExchangeModel(dict,sm), propsDict_(dict.subDict(typeName + "Props")), + DEMVariableNames_(propsDict_.lookupOrDefault("DEMVariables",wordList(0))), + DEMVariables_ + ( + IOobject + ( + "DEMVariables", + sm.mesh().time().timeName(), + sm.mesh(), + IOobject::NO_READ, + IOobject::NO_WRITE + ) + ), lmp(NULL) { Info << "Starting up LIGGGHTS for first time execution" << endl; @@ -92,10 +104,38 @@ twoWayMPI::twoWayMPI twoWayMPI::~twoWayMPI() { + if (propsDict_.found("liggghtsEndOfRunPath")) + { + const fileName liggghtsEndOfRunPath(propsDict_.lookup("liggghtsEndOfRunPath")); + lmp->input->file(liggghtsEndOfRunPath.c_str()); + } delete lmp; } + +// * * * * * * * * * * * * * * * private Member Functions * * * * * * * * * * * * * // + +void twoWayMPI::updateDEMVariables() +{ + scalar variablevalue = 0.0; + word variablename(""); + DEMVariables_.clear(); + for (label i=0; i #include "dataExchangeModel.H" +#include "scalarIOList.H" //=================================// //LAMMPS/LIGGGHTS @@ -71,7 +72,14 @@ private: dictionary propsDict_; MPI_Comm comm_liggghts; + wordList DEMVariableNames_; + + scalarIOList DEMVariables_; + // private member functions + double DEMVariableValue(word); + + void updateDEMVariables(); protected: LAMMPS_NS::LAMMPS *lmp; @@ -146,7 +154,6 @@ public: int getNumberOfTypes() const; double* getTypeVol() const; - scalar getCG() const { return lmp->force->cg(); } }; diff --git a/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/twoWayMany2Many.C b/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/twoWayMany2Many.C index 79d50ed5..da2df875 100644 --- a/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/twoWayMany2Many.C +++ b/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/twoWayMany2Many.C @@ -163,6 +163,12 @@ twoWayMany2Many::~twoWayMany2Many() delete lmp2foam_vec_; delete foam2lmp_vec_; delete foam2lmp_; + + if (propsDict_.found("liggghtsEndOfRunPath")) + { + const fileName liggghtsEndOfRunPath(propsDict_.lookup("liggghtsEndOfRunPath")); + lmp->input->file(liggghtsEndOfRunPath.c_str()); + } delete lmp; } diff --git a/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayOne2One/twoWayOne2One.C b/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayOne2One/twoWayOne2One.C index 0aa2ab91..b8c314f6 100644 --- a/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayOne2One/twoWayOne2One.C +++ b/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayOne2One/twoWayOne2One.C @@ -242,6 +242,11 @@ twoWayOne2One::~twoWayOne2One() destroy(prev_cell_ids_); destroy(dbl_cell_ids_); + if (propsDict_.found("liggghtsEndOfRunPath")) + { + const fileName liggghtsEndOfRunPath(propsDict_.lookup("liggghtsEndOfRunPath")); + lmp->input->file(liggghtsEndOfRunPath.c_str()); + } delete lmp; } @@ -961,6 +966,10 @@ void twoWayOne2One::extractCollected(T*& src, T**& dst, int width) const } } +int twoWayOne2One::getNumberOfParticles() const +{ + return particleCloud_.numberOfParticles(); +} // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // diff --git a/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayOne2One/twoWayOne2One.H b/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayOne2One/twoWayOne2One.H index 9f436b50..a6ee34da 100644 --- a/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayOne2One/twoWayOne2One.H +++ b/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayOne2One/twoWayOne2One.H @@ -208,6 +208,8 @@ public: template void extractCollected(T*&, T**&, int width=1) const; + int getNumberOfParticles() const; + scalar getCG() const { return lmp->force->cg(); } }; diff --git a/src/lagrangian/cfdemParticle/subModels/energyModel/heatTransferGunn/heatTransferGunn.C b/src/lagrangian/cfdemParticle/subModels/energyModel/heatTransferGunn/heatTransferGunn.C index 5570b610..ece81a57 100644 --- a/src/lagrangian/cfdemParticle/subModels/energyModel/heatTransferGunn/heatTransferGunn.C +++ b/src/lagrangian/cfdemParticle/subModels/energyModel/heatTransferGunn/heatTransferGunn.C @@ -51,6 +51,8 @@ heatTransferGunn::heatTransferGunn implicit_(propsDict_.lookupOrDefault("implicit",true)), calcTotalHeatFlux_(propsDict_.lookupOrDefault("calcTotalHeatFlux",true)), initPartTemp_(propsDict_.lookupOrDefault("initPartTemp",false)), + Tmin_(propsDict_.lookupOrDefault("Tmin",0.0)), + Tmax_(propsDict_.lookupOrDefault("Tmax",1e6)), totalHeatFlux_(0.0), NusseltScalingFactor_(1.0), QPartFluidName_(propsDict_.lookupOrDefault("QPartFluidName","QPartFluid")), @@ -148,7 +150,8 @@ heatTransferGunn::heatTransferGunn partRe_(NULL), partNu_(NULL), scaleDia_(1.), - typeCG_(propsDict_.lookupOrDefault("coarseGrainingFactors",scalarList(1,1.0))) + typeCG_(propsDict_.lookupOrDefault("coarseGrainingFactors",scalarList(1,1.0))), + maxTypeCG_(typeCG_.size()) { allocateMyArrays(); @@ -207,7 +210,10 @@ heatTransferGunn::heatTransferGunn scaleDia_=scalar(readScalar(propsDict_.lookup("scale"))); typeCG_[0] = scaleDia_; } - else if (typeCG_.size()>1) multiTypes_ = true; + else if (typeCG_.size()>1) + { + multiTypes_ = true; + } if (initPartTemp_ && !partTempField_.headerOk()) { @@ -338,6 +344,10 @@ void heatTransferGunn::calcEnergyContribution() if (multiTypes_) { partType = particleCloud_.particleType(index); + if (partType > maxTypeCG_) + { + FatalError<< "Too few coarse-graining factors provided." << abort(FatalError); + } cg = typeCG_[partType - 1]; scaleDia3 = cg*cg*cg; } @@ -596,12 +606,22 @@ void heatTransferGunn::partTempField() void heatTransferGunn::initPartTemp() { label cellI = 0; + scalar T = 0.0; for(int index = 0;index < particleCloud_.numberOfParticles(); ++index) { cellI = particleCloud_.cellIDs()[index][0]; if(cellI >= 0) { - partTemp_[index][0] = partTempField_[cellI]; + T = partTempField_[cellI]; + if (T < Tmin_) + { + T = Tmin_; + } + else if (T > Tmax_) + { + T = Tmax_; + } + partTemp_[index][0] = T; } } } diff --git a/src/lagrangian/cfdemParticle/subModels/energyModel/heatTransferGunn/heatTransferGunn.H b/src/lagrangian/cfdemParticle/subModels/energyModel/heatTransferGunn/heatTransferGunn.H index f9aa8e80..06b9a5ac 100644 --- a/src/lagrangian/cfdemParticle/subModels/energyModel/heatTransferGunn/heatTransferGunn.H +++ b/src/lagrangian/cfdemParticle/subModels/energyModel/heatTransferGunn/heatTransferGunn.H @@ -60,6 +60,10 @@ protected: bool initPartTemp_; + scalar Tmin_; + + scalar Tmax_; + scalar totalHeatFlux_; scalar NusseltScalingFactor_; @@ -122,6 +126,8 @@ protected: scalarList typeCG_; + const label maxTypeCG_; + void allocateMyArrays() const; void partTempField(); diff --git a/src/lagrangian/cfdemParticle/subModels/energyModel/heatTransferRanzMarshall/heatTransferRanzMarshall.C b/src/lagrangian/cfdemParticle/subModels/energyModel/heatTransferRanzMarshall/heatTransferRanzMarshall.C index 80938e60..fdfabc13 100644 --- a/src/lagrangian/cfdemParticle/subModels/energyModel/heatTransferRanzMarshall/heatTransferRanzMarshall.C +++ b/src/lagrangian/cfdemParticle/subModels/energyModel/heatTransferRanzMarshall/heatTransferRanzMarshall.C @@ -49,6 +49,9 @@ heatTransferRanzMarshall::heatTransferRanzMarshall verbose_(propsDict_.lookupOrDefault("verbose",false)), implicit_(propsDict_.lookupOrDefault("implicit",true)), calcTotalHeatFlux_(propsDict_.lookupOrDefault("calcTotalHeatFlux",true)), + initPartTemp_(propsDict_.lookupOrDefault("initPartTemp",false)), + Tmin_(propsDict_.lookupOrDefault("Tmin",0.0)), + Tmax_(propsDict_.lookupOrDefault("Tmax",1e6)), totalHeatFlux_(0.0), NusseltScalingFactor_(1.0), QPartFluidName_(propsDict_.lookupOrDefault("QPartFluidName","QPartFluid")), @@ -146,7 +149,8 @@ heatTransferRanzMarshall::heatTransferRanzMarshall partRe_(NULL), partNu_(NULL), scaleDia_(1.), - typeCG_(propsDict_.lookupOrDefault("coarseGrainingFactors",scalarList(1,1.0))) + typeCG_(propsDict_.lookupOrDefault("coarseGrainingFactors",scalarList(1,1.0))), + maxTypeCG_(typeCG_.size()) { allocateMyArrays(); @@ -205,7 +209,15 @@ heatTransferRanzMarshall::heatTransferRanzMarshall scaleDia_=scalar(readScalar(propsDict_.lookup("scale"))); typeCG_[0] = scaleDia_; } - else if (typeCG_.size()>1) multiTypes_ = true; + else if (typeCG_.size()>1) + { + multiTypes_ = true; + } + + if (initPartTemp_ && !partTempField_.headerOk()) + { + FatalError <<"Trying to initialize particle temperatures, but no field found.\n" << abort(FatalError); + } } @@ -252,8 +264,18 @@ void heatTransferRanzMarshall::calcEnergyContribution() // reset Scalar field QPartFluid_.primitiveFieldRef() = 0.0; - // get DEM data - particleCloud_.dataExchangeM().getData(partTempName_,"scalar-atom",partTemp_); + if (initPartTemp_) + { + // if particle temperatures are to be initialized from field, do a one-time push to DEM + initPartTemp(); + particleCloud_.dataExchangeM().giveData("Temp","scalar-atom", partTemp_); + initPartTemp_ = false; + } + else + { + // get DEM data + particleCloud_.dataExchangeM().getData(partTempName_,"scalar-atom",partTemp_); + } #ifdef compre const volScalarField mufField = particleCloud_.turbulence().mu(); @@ -322,6 +344,10 @@ void heatTransferRanzMarshall::calcEnergyContribution() if (multiTypes_) { partType = particleCloud_.particleType(index); + if (partType > maxTypeCG_) + { + FatalError<< "Too few coarse-graining factors provided." << abort(FatalError); + } cg = typeCG_[partType - 1]; scaleDia3 = cg*cg*cg; } @@ -564,6 +590,29 @@ void heatTransferRanzMarshall::partTempField() Info << "heatTransferRanzMarshall: average part. temp = " << partTempAve_.value() << endl; } + +void heatTransferRanzMarshall::initPartTemp() +{ + label cellI = 0; + scalar T = 0.0; + for(int index = 0;index < particleCloud_.numberOfParticles(); ++index) + { + cellI = particleCloud_.cellIDs()[index][0]; + if(cellI >= 0) + { + T = partTempField_[cellI]; + if (T < Tmin_) + { + T = Tmin_; + } + else if (T > Tmax_) + { + T = Tmax_; + } + partTemp_[index][0] = T; + } + } +} // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // } // End namespace Foam diff --git a/src/lagrangian/cfdemParticle/subModels/energyModel/heatTransferRanzMarshall/heatTransferRanzMarshall.H b/src/lagrangian/cfdemParticle/subModels/energyModel/heatTransferRanzMarshall/heatTransferRanzMarshall.H index f627dc12..17baedb8 100644 --- a/src/lagrangian/cfdemParticle/subModels/energyModel/heatTransferRanzMarshall/heatTransferRanzMarshall.H +++ b/src/lagrangian/cfdemParticle/subModels/energyModel/heatTransferRanzMarshall/heatTransferRanzMarshall.H @@ -57,6 +57,12 @@ protected: bool calcTotalHeatFlux_; + bool initPartTemp_; + + scalar Tmin_; + + scalar Tmax_; + scalar totalHeatFlux_; scalar NusseltScalingFactor_; @@ -119,6 +125,8 @@ protected: scalarList typeCG_; + const label maxTypeCG_; + void allocateMyArrays() const; void partTempField(); @@ -129,6 +137,8 @@ protected: virtual void heatFlux(label, scalar, scalar, scalar, scalar cg3 = 1.0); + virtual void initPartTemp(); + public: //- Runtime type information diff --git a/src/lagrangian/cfdemParticle/subModels/forceModel/BeetstraDrag/BeetstraDrag.C b/src/lagrangian/cfdemParticle/subModels/forceModel/BeetstraDrag/BeetstraDrag.C index 18381aa3..e4eef410 100644 --- a/src/lagrangian/cfdemParticle/subModels/forceModel/BeetstraDrag/BeetstraDrag.C +++ b/src/lagrangian/cfdemParticle/subModels/forceModel/BeetstraDrag/BeetstraDrag.C @@ -58,8 +58,9 @@ BeetstraDrag::BeetstraDrag minVoidfraction_(propsDict_.lookupOrDefault("minVoidfraction",0.1)), UsFieldName_(propsDict_.lookup("granVelFieldName")), UsField_(sm.mesh().lookupObject (UsFieldName_)), - scaleDia_(1.), + maxTypeCG_(1), typeCG_(propsDict_.lookupOrDefault("coarseGrainingFactors",scalarList(1,1.0))), + scaleDia_(1.), scaleDrag_(1.), rhoP_(0.), rho_(0.), @@ -102,7 +103,11 @@ BeetstraDrag::BeetstraDrag scaleDrag_=scalar(readScalar(propsDict_.lookup("scaleDrag"))); } - if (typeCG_.size()>1) multiTypes_ = true; + if (typeCG_.size()>1) + { + multiTypes_ = true; + maxTypeCG_ = typeCG_.size(); + } if (propsDict_.found("useFilteredDragModel")) { @@ -213,12 +218,14 @@ void BeetstraDrag::setForce() const voidfraction = voidfraction_[cellI]; Ufluid = U_[cellI]; } - // in case a fines phase is present, void fraction needs to be adapted - adaptVoidfraction(voidfraction, cellI); if (multiTypes_) { partType = particleCloud_.particleType(index); + if (partType > maxTypeCG_) + { + FatalError<< "Too few coarse-graining factors provided." << abort(FatalError); + } cg = typeCG_[partType - 1]; scaleDia3 = cg*cg*cg; } diff --git a/src/lagrangian/cfdemParticle/subModels/forceModel/BeetstraDrag/BeetstraDrag.H b/src/lagrangian/cfdemParticle/subModels/forceModel/BeetstraDrag/BeetstraDrag.H index cfeeea39..d66b81f3 100644 --- a/src/lagrangian/cfdemParticle/subModels/forceModel/BeetstraDrag/BeetstraDrag.H +++ b/src/lagrangian/cfdemParticle/subModels/forceModel/BeetstraDrag/BeetstraDrag.H @@ -64,10 +64,12 @@ protected: const volVectorField& UsField_; - mutable scalar scaleDia_; + mutable label maxTypeCG_; mutable scalarList typeCG_; + mutable scalar scaleDia_; + mutable scalar scaleDrag_; mutable scalar rhoP_; @@ -88,8 +90,6 @@ protected: bool usePC_; - virtual void adaptVoidfraction(double&, label) const {} - virtual scalar effDiameter(double d, label cellI, label index) const {return d;} virtual scalar meanSauterDiameter(double d, label cellI) const {return d;} diff --git a/src/lagrangian/cfdemParticle/subModels/forceModel/BeetstraDragPoly/BeetstraDragPoly.C b/src/lagrangian/cfdemParticle/subModels/forceModel/BeetstraDragPoly/BeetstraDragPoly.C index 04e12ed9..ae6e05e7 100644 --- a/src/lagrangian/cfdemParticle/subModels/forceModel/BeetstraDragPoly/BeetstraDragPoly.C +++ b/src/lagrangian/cfdemParticle/subModels/forceModel/BeetstraDragPoly/BeetstraDragPoly.C @@ -49,23 +49,26 @@ BeetstraDragPoly::BeetstraDragPoly : BeetstraDrag(dict,sm), fines_(propsDict_.lookupOrDefault("fines",false)), - dFine_(1.0) + dFine_(1.0), + alphaP_(NULL), + alphaSt_(NULL), + dSauter_(NULL) { // if fines are present, take mixture dSauter, otherwise normal dSauter if (fines_) { dFine_ = readScalar(propsDict_.lookup("dFine")); volScalarField& alphaP(const_cast(sm.mesh().lookupObject ("alphaP"))); - alphaP_.set(&alphaP); + alphaP_ = &alphaP; volScalarField& alphaSt(const_cast(sm.mesh().lookupObject ("alphaSt"))); - alphaSt_.set(&alphaSt); + alphaSt_ = &alphaSt; volScalarField& dSauter(const_cast(sm.mesh().lookupObject ("dSauterMix"))); - dSauter_.set(&dSauter); + dSauter_ = &dSauter; } else { volScalarField& dSauter(const_cast(sm.mesh().lookupObject ("dSauter"))); - dSauter_.set(&dSauter); + dSauter_ = &dSauter; } } @@ -77,22 +80,15 @@ BeetstraDragPoly::~BeetstraDragPoly() // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // - -void BeetstraDragPoly::adaptVoidfraction(double& voidfraction, label cellI) const -{ - if (fines_) voidfraction -= alphaSt_()[cellI]; - if (voidfraction < minVoidfraction_) voidfraction = minVoidfraction_; -} - scalar BeetstraDragPoly::effDiameter(double d, label cellI, label index) const { - scalar dS = dSauter_()[cellI]; + scalar dS = (*dSauter_)[cellI]; scalar effD = d*d / dS + 0.064*d*d*d*d / (dS*dS*dS); if (fines_) { scalar fineCorr = dFine_*dFine_ / dS + 0.064*dFine_*dFine_*dFine_*dFine_ / (dS*dS*dS); - fineCorr *= d*d*d / (dFine_*dFine_*dFine_) * alphaSt_()[cellI] / alphaP_()[cellI]; + fineCorr *= d*d*d / (dFine_*dFine_*dFine_) * (*alphaSt_)[cellI] / (*alphaP_)[cellI]; effD += fineCorr; } @@ -107,7 +103,7 @@ scalar BeetstraDragPoly::effDiameter(double d, label cellI, label index) const scalar BeetstraDragPoly::meanSauterDiameter(double d, label cellI) const { - return dSauter_()[cellI]; + return (*dSauter_)[cellI]; } // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // diff --git a/src/lagrangian/cfdemParticle/subModels/forceModel/BeetstraDragPoly/BeetstraDragPoly.H b/src/lagrangian/cfdemParticle/subModels/forceModel/BeetstraDragPoly/BeetstraDragPoly.H index 88965010..8c91b7cf 100644 --- a/src/lagrangian/cfdemParticle/subModels/forceModel/BeetstraDragPoly/BeetstraDragPoly.H +++ b/src/lagrangian/cfdemParticle/subModels/forceModel/BeetstraDragPoly/BeetstraDragPoly.H @@ -46,13 +46,11 @@ protected: scalar dFine_; - autoPtr alphaP_; + volScalarField* alphaP_; - autoPtr alphaSt_; + volScalarField* alphaSt_; - autoPtr dSauter_; - - void adaptVoidfraction(double&, label) const; + volScalarField* dSauter_; scalar effDiameter(double, label, label) const; diff --git a/src/lagrangian/cfdemParticle/subModels/forceModel/Fines/Fines.C b/src/lagrangian/cfdemParticle/subModels/forceModel/Fines/Fines.C index 5c4e6efa..ce779474 100644 --- a/src/lagrangian/cfdemParticle/subModels/forceModel/Fines/Fines.C +++ b/src/lagrangian/cfdemParticle/subModels/forceModel/Fines/Fines.C @@ -50,6 +50,11 @@ Fines::Fines forceModel(dict,sm), finesFields_(dict,sm) { + // safety check: fines should be first force model in list because they correct the voidfraction field + if (particleCloud_.forceModels()[0] != "dSauter" || particleCloud_.forceModels()[1] != "Fines") + { + FatalError <<"Force models dSauter and Fines need to be first in list.\n" << abort(FatalError); + } } diff --git a/src/lagrangian/cfdemParticle/subModels/forceModel/Fines/FinesFields.C b/src/lagrangian/cfdemParticle/subModels/forceModel/Fines/FinesFields.C index 1147bbe6..d17aabb5 100644 --- a/src/lagrangian/cfdemParticle/subModels/forceModel/Fines/FinesFields.C +++ b/src/lagrangian/cfdemParticle/subModels/forceModel/Fines/FinesFields.C @@ -45,11 +45,20 @@ FinesFields::FinesFields propsDict_(dict.subDict("FinesFieldsProps")), smoothing_(propsDict_.lookupOrDefault("smoothing",false)), verbose_(propsDict_.lookupOrDefault("verbose",false)), + clogKin_(propsDict_.lookupOrDefault("kineticClogging",false)), + clogStick_(propsDict_.lookupOrDefault("stickyClogging",false)), + movingBed_(propsDict_.lookupOrDefault("movingBed",true)), + useDepositionLength_(false), + fluxFieldName_(propsDict_.lookupOrDefault("fluxFieldName","phi")), + phi_(sm.mesh().lookupObject (fluxFieldName_)), velFieldName_(propsDict_.lookupOrDefault("velFieldName","U")), U_(sm.mesh().lookupObject (velFieldName_)), voidfractionFieldName_(propsDict_.lookupOrDefault("voidfractionFieldName","voidfraction")), - // this is probably really bad - voidfraction_(const_cast(sm.mesh().lookupObject (voidfractionFieldName_))), +#if OPENFOAM_VERSION_MAJOR < 5 + voidfraction_(const_cast(sm.mesh().lookupObject(voidfractionFieldName_))), +#else + voidfraction_(sm.mesh().lookupObjectRef (voidfractionFieldName_)), +#endif UsFieldName_(propsDict_.lookupOrDefault("granVelFieldName","Us")), UsField_(sm.mesh().lookupObject (UsFieldName_)), pFieldName_(propsDict_.lookupOrDefault("pFieldName","p")), @@ -200,7 +209,6 @@ FinesFields::FinesFields ), sm.mesh(), dimensionedScalar("zero", dimensionSet(1,0,-1,0,0), 0) - //dimensionedVector("zero", dimensionSet(1,-2,-1,0,0), vector::zero) ), uDyn_ ( IOobject @@ -219,54 +227,134 @@ FinesFields::FinesFields rhoFine_("rhoFine",dimensionSet(1,-3,0,0,0),0.0), g_("g",dimensionSet(0,1,-2,0,0),vector(0,0,-9.81)), alphaDynMax_(0.1), - alphaMax_(readScalar(propsDict_.lookup ("alphaMax"))), - critVoidfraction_(readScalar(propsDict_.lookup ("critVoidfraction"))), - depRate_(readScalar(propsDict_.lookup ("depRate"))), + alphaMax_(propsDict_.lookupOrDefault("alphaMax",0.95)), + alphaMinClog_(propsDict_.lookupOrDefault("alphaMinClog",0.1)), + critVoidfraction_(propsDict_.lookupOrDefault("critVoidfraction", 0.05)), + deltaT_(voidfraction_.mesh().time().deltaTValue()), + depositionLength_(0.0), exponent_(-1.33), - nCrit_(readScalar(propsDict_.lookup ("nCrit"))), - poresizeWidth_(readScalar(propsDict_.lookup ("poresizeWidth"))), + nCrit_(0.0), + poresizeWidth_(0.0), prefactor_(10.5), - ratioHydraulicPore_(1.5) + ratioHydraulicPore_(1.5), + tauDeposition_(0.0), + tauRelease_(readScalar(propsDict_.lookup ("tauRelease"))), + uBindHigh_(propsDict_.lookupOrDefault("uBindHigh",0.0)), + uBindLow_(propsDict_.lookupOrDefault("uBindLow",0.0)), + uMin_(0.001) { Sds_.write(); if (propsDict_.found("prefactor")) - prefactor_=readScalar(propsDict_.lookup ("prefactor")); - if (propsDict_.found("exponent")) - exponent_=readScalar(propsDict_.lookup ("exponent")); - if (propsDict_.found("dFine")) - dFine_.value()=readScalar(propsDict_.lookup ("dFine")); - else - FatalError <<"Please specify dFine.\n" << abort(FatalError); - if (propsDict_.found("diffCoeff")) - diffCoeff_.value()=readScalar(propsDict_.lookup ("diffCoeff")); - if (propsDict_.found("rhoFine")) - rhoFine_.value()=readScalar(propsDict_.lookup ("rhoFine")); - else - FatalError <<"Please specify rhoFine.\n" << abort(FatalError); - if (propsDict_.found("nuAve")) - nuAve_.value()=readScalar(propsDict_.lookup ("nuAve")); - if (propsDict_.found("alphaDynMax")) - alphaDynMax_=readScalar(propsDict_.lookup ("alphaDynMax")); - - if(verbose_) { - alphaG_.writeOpt() = IOobject::AUTO_WRITE; - alphaG_.write(); - alphaP_.writeOpt() = IOobject::AUTO_WRITE; - alphaP_.write(); - deltaAlpha_.writeOpt() = IOobject::AUTO_WRITE; - deltaAlpha_.write(); - dHydMix_.writeOpt() = IOobject::AUTO_WRITE; - dHydMix_.write(); - DragCoeff_.writeOpt() = IOobject::AUTO_WRITE; - DragCoeff_.write(); - dSauterMix_.writeOpt() = IOobject::AUTO_WRITE; - dSauterMix_.write(); - FanningCoeff_.writeOpt() = IOobject::AUTO_WRITE; - FanningCoeff_.write(); - Froude_.writeOpt() = IOobject::AUTO_WRITE; - Froude_.write(); + prefactor_=readScalar(propsDict_.lookup ("prefactor")); + } + if (propsDict_.found("exponent")) + { + exponent_=readScalar(propsDict_.lookup ("exponent")); + } + if (propsDict_.found("dFine")) + { + dFine_.value()=readScalar(propsDict_.lookup ("dFine")); + } + else + { + FatalError <<"Please specify dFine.\n" << abort(FatalError); + } + if (propsDict_.found("diffCoeff")) + { + diffCoeff_.value()=readScalar(propsDict_.lookup ("diffCoeff")); + } + if (propsDict_.found("rhoFine")) + { + rhoFine_.value()=readScalar(propsDict_.lookup ("rhoFine")); + } + else + { + FatalError <<"Please specify rhoFine.\n" << abort(FatalError); + } + if (propsDict_.found("nuAve")) + { + nuAve_.value()=readScalar(propsDict_.lookup ("nuAve")); + } + if (propsDict_.found("alphaDynMax")) + { + alphaDynMax_=readScalar(propsDict_.lookup ("alphaDynMax")); + } + + if (clogKin_) + { + nCrit_ = readScalar(propsDict_.lookup ("nCrit")); + poresizeWidth_ = readScalar(propsDict_.lookup ("poresizeWidth")); + } + + if (clogStick_) + { + if (uBindHigh_ - uBindLow_ < SMALL) + { + FatalError <<"No reasonable values for critical binding velocities.\n" << abort(FatalError); + } + + uReconstructed_.set + ( + new volVectorField + ( + IOobject + ( + "uReconstructed", + sm.mesh().time().timeName(), + sm.mesh(), + IOobject::NO_READ, + IOobject::NO_WRITE + ), + U_//sm.mesh() + ) + ); + + } + + if (propsDict_.found("depositionLength") && propsDict_.found("tauDeposition")) + { + FatalError <<"You cannot specify both a deposition length and time.\n" << abort(FatalError); + } + else if (propsDict_.found("depositionLength")) + { + depositionLength_ = readScalar(propsDict_.lookup ("depositionLength")); + useDepositionLength_ = true; + } + else if (propsDict_.found("tauDeposition")) + { + tauDeposition_ = readScalar(propsDict_.lookup ("tauDeposition")); + useDepositionLength_ = false; + } + else + { + FatalError <<"Specify either a deposition length or time.\n" << abort(FatalError); + } + + if (verbose_) + { + alphaG_.writeOpt() = IOobject::AUTO_WRITE; + alphaG_.write(); + alphaP_.writeOpt() = IOobject::AUTO_WRITE; + alphaP_.write(); + deltaAlpha_.writeOpt() = IOobject::AUTO_WRITE; + deltaAlpha_.write(); + dHydMix_.writeOpt() = IOobject::AUTO_WRITE; + dHydMix_.write(); + DragCoeff_.writeOpt() = IOobject::AUTO_WRITE; + DragCoeff_.write(); + dSauterMix_.writeOpt() = IOobject::AUTO_WRITE; + dSauterMix_.write(); + FanningCoeff_.writeOpt() = IOobject::AUTO_WRITE; + FanningCoeff_.write(); + Froude_.writeOpt() = IOobject::AUTO_WRITE; + Froude_.write(); + if (clogStick_) + { + uReconstructed_().writeOpt() = IOobject::AUTO_WRITE; + uReconstructed_().write(); + } } } @@ -281,76 +369,100 @@ FinesFields::~FinesFields() void FinesFields::update() { - if(verbose_) Info << "FinesFields: Updating alphaP.\n" << endl; + if (verbose_) Info << "FinesFields: Updating alphaP.\n" << endl; updateAlphaP(); - if(verbose_) Info << "FinesFields: Updating alphaG.\n" << endl; + if (verbose_) Info << "FinesFields: Updating alphaG.\n" << endl; updateAlphaG(); - if(verbose_) Info << "FinesFields: Calculating source terms.\n" << endl; + if (clogStick_) + { + if (verbose_) Info << "FinesFields: Updating uReconstructed.\n" << endl; + updateUReconstructed(); + } + if (verbose_) Info << "FinesFields: Calculating source terms.\n" << endl; calcSource(); - if(verbose_) Info << "FinesFields: Updating dSauter.\n" << endl; + if (verbose_) Info << "FinesFields: Updating dSauter.\n" << endl; updateDSauter(); - if(verbose_) Info << "FinesFields: Updating dHydMix.\n" << endl; + if (verbose_) Info << "FinesFields: Updating dHydMix.\n" << endl; updateDHydMix(); - if(verbose_) Info << "FinesFields: Updating Froude.\n" << endl; + if (verbose_) Info << "FinesFields: Updating Froude.\n" << endl; updateFroude(); - if(verbose_) Info << "FinesFields: Updating FanningCoeff.\n" << endl; + if (verbose_) Info << "FinesFields: Updating FanningCoeff.\n" << endl; updateFanningCoeff(); - if(verbose_) Info << "FinesFields: Updating DragCoeff.\n" << endl; + if (verbose_) Info << "FinesFields: Updating DragCoeff.\n" << endl; updateDragCoeff(); - if(verbose_) Info << "FinesFields: Updating uDyn.\n" << endl; + if (verbose_) Info << "FinesFields: Updating uDyn.\n" << endl; updateUDyn(); - if(verbose_) Info << "FinesFields: Integrating alphas.\n" << endl; + if (verbose_) Info << "FinesFields: Integrating alphas.\n" << endl; integrateFields(); - if(verbose_) Info << "FinesFields: Update finished.\n" << endl; + if (verbose_) Info << "FinesFields: Update finished.\n" << endl; } void FinesFields::calcSource() { + if (!clogKin_ & !clogStick_) return; Sds_.primitiveFieldRef() = 0; deltaAlpha_.primitiveFieldRef() = 0.0; - scalar f(0.0); - scalar critpore(0.0); - scalar dmean(0.0); - scalar d1(0.0); - scalar d2(0.0); + scalar fKin = 0.0; + scalar fStick = 0.0; + scalar critpore = 0.0; + scalar dmean = 0.0; + scalar d1 = 0.0; + scalar d2 = 0.0; + scalar magU = 0.0; + scalar tauDeposition = 0.0; forAll(Sds_,cellI) { - // calculate everything in units auf dSauter - critpore = nCrit_*dFine_.value()/dSauter_[cellI]; - // pore size from hydraulic radius - dmean = 2 * (1 - alphaP_[cellI]) / ( (1 + poresizeWidth_*poresizeWidth_/3) * 3 * alphaP_[cellI] ); - // Sweeney and Martin, Acta Materialia 51 (2003): ratio of hydraulic to pore throat radius - dmean /= ratioHydraulicPore_; - d1 = dmean * (1 - poresizeWidth_); - d2 = dmean * (1 + poresizeWidth_); + fKin = 0.0; + fStick = 0.0; + if (clogKin_ && alphaP_[cellI] > alphaMinClog_) // no cloggig in dilute regions + { + // calculate everything in units auf dSauter + critpore = nCrit_*dFine_.value()/dSauter_[cellI]; + // pore size from hydraulic radius + dmean = 2 * (1 - alphaP_[cellI]) / ( (1 + poresizeWidth_*poresizeWidth_/3) * 3 * alphaP_[cellI] ); + // Sweeney and Martin, Acta Materialia 51 (2003): ratio of hydraulic to pore throat radius + dmean /= ratioHydraulicPore_; + d1 = dmean * (1 - poresizeWidth_); + d2 = dmean * (1 + poresizeWidth_); - f = (critpore*critpore*critpore - d1 * d1 * d1) / (d2 * d2 * d2 - d1 * d1 * d1); - if (f < 0) - { - f = 0.0; + fKin = (critpore*critpore*critpore - d1 * d1 * d1) / (d2 * d2 * d2 - d1 * d1 * d1); + if (fKin < 0) fKin = 0.0; + else if (fKin > 1.0) fKin = 1.0; } - else if (f > 1.0) + + if (clogStick_ && alphaP_[cellI] > alphaMinClog_) // no cloggig in dilute regions { - f = 1.0; + magU = mag(uReconstructed_()[cellI]); // use U reconstructed from phi to suppress oscillations at interfaces + // fStick = 1.0 / ( 1.0 + magU/uBind_) * alphaP_[cellI] / 0.65; + if (magU < uBindLow_) fStick = 1.0; + else if (magU > uBindHigh_) fStick = 0.0; + else fStick = 1.0 - (magU - uBindLow_) / (uBindHigh_ - uBindLow_); + // if (fStick > 1.0) fStick = 1.0; + fStick *= alphaP_[cellI] / 0.65; } // at this point, voidfraction is still calculated from the true particle sizes - deltaAlpha_[cellI] = f * (alphaMax_ - alphaP_[cellI]) - alphaSt_[cellI]; - // too much volume occupied: release it (50% per time step) + deltaAlpha_[cellI] = max(fKin,fStick) * (alphaMax_ - alphaP_[cellI]) - alphaSt_[cellI]; + // too much volume occupied: release it if (deltaAlpha_[cellI] < 0.0) { - Sds_[cellI] = 0.5*deltaAlpha_[cellI]; - } - // volume too occupy available: deposit at most 80% of dyn hold up - else if (depRate_ * deltaAlpha_[cellI] > 0.8 * alphaDyn_[cellI]) - { - Sds_[cellI] = 0.8 * alphaDyn_[cellI]; + Sds_[cellI] = deltaAlpha_[cellI] / tauRelease_; } + // volume to occupy available else { - Sds_[cellI] = depRate_ * deltaAlpha_[cellI]; + if (useDepositionLength_) + { + tauDeposition = depositionLength_ / max(mag(U_[cellI]),uMin_); + } + else + { + tauDeposition = tauDeposition_; + } + if (tauDeposition < 4.0*deltaT_) tauDeposition = 4.0*deltaT_; // do not deposit all dyn hold-up within less than 3 time steps + Sds_[cellI] = alphaDyn_[cellI] / tauDeposition; } } } @@ -364,11 +476,12 @@ void FinesFields::integrateFields() fvScalarMatrix alphaStEqn ( - fvm::ddt(alphaSt_) - + fvm::div(phiSt,alphaSt_) + fvm::ddt(alphaSt_) == Sds_ ); + if (movingBed_) alphaStEqn += fvm::div(phiSt,alphaSt_); + fvScalarMatrix alphaDynEqn ( fvm::ddt(alphaDyn_) @@ -377,10 +490,11 @@ void FinesFields::integrateFields() == -Sds_ ); + alphaStEqn.solve(); alphaDynEqn.solve(); - if(smoothing_) + if (smoothing_) particleCloud_.smoothingM().smoothen(alphaDyn_); // limit hold-ups, should be done more elegantly @@ -422,13 +536,14 @@ void FinesFields::integrateFields() } -void FinesFields::updateAlphaG() +void FinesFields::updateAlphaG() // called after updateAlphaP() - correct voidfraction by removing space occupied by fines { alphaG_ = max(voidfraction_ - alphaSt_ - alphaDyn_, critVoidfraction_); + voidfraction_ = alphaG_; } -void FinesFields::updateAlphaP() +void FinesFields::updateAlphaP() // called first in the update cycle - voidfraction_ is current with DEM data { alphaP_ = 1.0 - voidfraction_ + SMALL; } @@ -439,7 +554,7 @@ void FinesFields::updateDHydMix() forAll(dHydMix_,cellI) { scalar aPSt = alphaP_[cellI] + alphaSt_[cellI]; - if(aPSt < SMALL || aPSt > 1 - SMALL) + if (aPSt < SMALL || aPSt > 1 - SMALL) dHydMix_[cellI] = SMALL; else dHydMix_[cellI] = 2*(1 - aPSt) / (3*aPSt ) * dSauterMix_[cellI]; @@ -459,11 +574,11 @@ void FinesFields::updateDragCoeff() forAll(DragCoeff_,cellI) { Ref1 = Ref[cellI]; - if(Ref1 <= SMALL) + if (Ref1 <= SMALL) Cd = 24.0 / SMALL; - else if(Ref1 <= 1.0) + else if (Ref1 <= 1.0) Cd = 24.0 / Ref1; - else if(Ref1 <= 1000) + else if (Ref1 <= 1000) Cd = 24 * (1.0 + 0.15 * Foam::pow(Ref1,0.687) ) / Ref1; else Cd = 0.44; @@ -475,11 +590,11 @@ void FinesFields::updateDragCoeff() forAll(DragCoeff_.boundaryField()[patchI], faceI) { Ref1 = Ref.boundaryField()[patchI][faceI]; - if(Ref1 <= SMALL) + if (Ref1 <= SMALL) Cd = 24.0 / SMALL; - else if(Ref1 <= 1.0) + else if (Ref1 <= 1.0) Cd = 24.0 / Ref1; - else if(Ref1 <= 1000) + else if (Ref1 <= 1000) Cd = 24 * (1.0 + 0.15 * Foam::pow(Ref1,0.687) ) / Ref1; else Cd = 0.44; @@ -496,9 +611,9 @@ void FinesFields::updateDSauter() { scalar aP = alphaP_[cellI]; scalar aSt = alphaSt_[cellI]; - if(aSt < SMALL) + if (aSt < SMALL) dSauterMix_[cellI] = dSauter_[cellI]; - else if(aP < SMALL) + else if (aP < SMALL) dSauterMix_[cellI] = dFine_.value(); else dSauterMix_[cellI] = (aP + aSt) / (aP / dSauter_[cellI] + aSt / dFine_.value() ); @@ -541,7 +656,7 @@ void FinesFields::updateUDyn() { scalar mU(mag(U_[cellI])); scalar muDyn(mag(uDyn_[cellI])); - if(muDyn > mU && muDyn > SMALL) + if (muDyn > mU && muDyn > SMALL) { uDyn_[cellI] *= mU / muDyn; } @@ -552,13 +667,24 @@ void FinesFields::updateUDyn() { scalar mU(mag(U_.boundaryField()[patchI][faceI])); scalar muDyn(mag(uDyn_.boundaryField()[patchI][faceI])); - if(muDyn > mU && muDyn > SMALL) + if (muDyn > mU && muDyn > SMALL) { uDyn_.boundaryFieldRef()[patchI][faceI] *= mU / muDyn; } } } +void FinesFields::updateUReconstructed() +{ + if (phi_.dimensions() == dimensionSet(1, 0, -1, 0, 0)) // compressible + { + uReconstructed_() = fvc::reconstruct(phi_) / (rhoG_ * voidfraction_); + } + else + { + uReconstructed_() = fvc::reconstruct(phi_) / voidfraction_; + } +} // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // } // End namespace Foam diff --git a/src/lagrangian/cfdemParticle/subModels/forceModel/Fines/FinesFields.H b/src/lagrangian/cfdemParticle/subModels/forceModel/Fines/FinesFields.H index 3098ff2f..ae673387 100644 --- a/src/lagrangian/cfdemParticle/subModels/forceModel/Fines/FinesFields.H +++ b/src/lagrangian/cfdemParticle/subModels/forceModel/Fines/FinesFields.H @@ -46,6 +46,18 @@ private: bool verbose_; + bool clogKin_; + + bool clogStick_; + + bool movingBed_; + + bool useDepositionLength_; + + word fluxFieldName_; + + const surfaceScalarField& phi_; + word velFieldName_; const volVectorField& U_; @@ -90,11 +102,12 @@ private: volScalarField Sds_; - //volVectorField massFluxDyn_; surfaceScalarField massFluxDyn_; volVectorField uDyn_; + autoPtr uReconstructed_; // velocity field can show oscillations at porosity steps + dimensionedScalar dFine_; dimensionedScalar diffCoeff_; @@ -109,9 +122,13 @@ private: scalar alphaMax_; + scalar alphaMinClog_; + scalar critVoidfraction_; - scalar depRate_; + scalar deltaT_; + + scalar depositionLength_; scalar exponent_; @@ -123,6 +140,16 @@ private: scalar ratioHydraulicPore_; + scalar tauDeposition_; + + scalar tauRelease_; + + scalar uBindHigh_; + + scalar uBindLow_; + + scalar uMin_; + void calcSource(); void integrateFields(); @@ -143,6 +170,8 @@ private: void updateUDyn(); + void updateUReconstructed(); + public: //- Runtime type information diff --git a/src/lagrangian/cfdemParticle/subModels/forceModel/KochHillRWDrag/KochHillRWDrag.C b/src/lagrangian/cfdemParticle/subModels/forceModel/KochHillRWDrag/KochHillRWDrag.C index 513f2cf7..496556d2 100755 --- a/src/lagrangian/cfdemParticle/subModels/forceModel/KochHillRWDrag/KochHillRWDrag.C +++ b/src/lagrangian/cfdemParticle/subModels/forceModel/KochHillRWDrag/KochHillRWDrag.C @@ -291,7 +291,11 @@ void KochHillRWDrag::setForce() const // modify current fluid velocity for (int dim=0; dim<3; dim++) { +#if OPENFOAM_VERSION_MAJOR < 6 partUfluct_[index][dim] = RanGen_.GaussNormal()*sqrt(2.*k/3.); +#else + partUfluct_[index][dim] = RanGen_.scalarNormal()*sqrt(2.*k/3.); +#endif //Pout << "RW-TEST: Ufluid[" << dim << "] = " << Ufluid[dim] << " Ufluct = " << partUfluct_[index][dim] << " k = " << k << endl; // TEST-Output Ufluid[dim] = Ufluid[dim] + partUfluct_[index][dim]; } diff --git a/src/lagrangian/cfdemParticle/subModels/forceModel/SchillerNaumannDrag/SchillerNaumannDrag.C b/src/lagrangian/cfdemParticle/subModels/forceModel/SchillerNaumannDrag/SchillerNaumannDrag.C index afe0e546..c274e462 100644 --- a/src/lagrangian/cfdemParticle/subModels/forceModel/SchillerNaumannDrag/SchillerNaumannDrag.C +++ b/src/lagrangian/cfdemParticle/subModels/forceModel/SchillerNaumannDrag/SchillerNaumannDrag.C @@ -66,7 +66,9 @@ SchillerNaumannDrag::SchillerNaumannDrag propsDict_(dict.subDict(typeName + "Props")), verbose_(propsDict_.found("verbose")), velFieldName_(propsDict_.lookup("velFieldName")), - U_(sm.mesh().lookupObject (velFieldName_)) + U_(sm.mesh().lookupObject (velFieldName_)), + scaleDia_(1.), + scaleDrag_(1.) { //Append the field names to be probed particleCloud_.probeM().initialize(typeName, typeName+".logDat"); @@ -85,7 +87,12 @@ SchillerNaumannDrag::SchillerNaumannDrag // read those switches defined above, if provided in dict forceSubM(0).readSwitches(); - particleCloud_.checkCG(false); + particleCloud_.checkCG(true); + + if (propsDict_.found("scale")) + scaleDia_ = scalar(readScalar(propsDict_.lookup("scale"))); + if (propsDict_.found("scaleDrag")) + scaleDrag_ = scalar(readScalar(propsDict_.lookup("scaleDrag"))); } @@ -99,6 +106,16 @@ SchillerNaumannDrag::~SchillerNaumannDrag() void SchillerNaumannDrag::setForce() const { + if (scaleDia_ > 1.0) + { + Info << "SchillerNaumann using scale = " << scaleDia_ << endl; + } + else if (particleCloud_.cg() > 1.0) + { + scaleDia_ = particleCloud_.cg(); + Info << "SchillerNaumann using scale from liggghts cg = " << scaleDia_ << endl; + } + #include "setupProbeModel.H" const volScalarField& nufField = forceSubM(0).nuField(); @@ -117,6 +134,7 @@ void SchillerNaumannDrag::setForce() const vector Us = particleCloud_.velocity(index); vector Ur = U_[cellI]-Us; scalar ds = 2*particleCloud_.radius(index); + scalar ds_scaled = ds/scaleDia_; scalar nuf = nufField[cellI]; scalar rho = rhoField[cellI]; scalar voidfraction = particleCloud_.voidfraction(index); @@ -127,30 +145,33 @@ void SchillerNaumannDrag::setForce() const if (magUr > 0) { // calc particle Re Nr - Rep = ds*magUr/nuf; + Rep = ds_scaled*magUr/nuf; // calc fluid drag Coeff Cd = max(0.44,24.0/Rep*(1.0+0.15*pow(Rep,0.687))); // calc particle's drag - drag = 0.125*Cd*rho*M_PI*ds*ds*magUr*Ur; + drag = 0.125*Cd*rho*M_PI*ds*ds*scaleDia_*magUr*Ur*scaleDrag_; if (modelType_=="B") drag /= voidfraction; } - if(verbose_ && index >100 && index <102) + if(verbose_ && index >=100 && index <102) { - Info << "index = " << index << endl; - Info << "Us = " << Us << endl; - Info << "Ur = " << Ur << endl; - Info << "ds = " << ds << endl; - Info << "rho = " << rho << endl; - Info << "nuf = " << nuf << endl; - Info << "voidfraction = " << voidfraction << endl; - Info << "Rep = " << Rep << endl; - Info << "Cd = " << Cd << endl; - Info << "drag = " << drag << endl; + Pout << "cellI = " << cellI << endl; + Pout << "index = " << index << endl; + Pout << "Ufluid = " << U_[cellI] << endl; + Pout << "Us = " << Us << endl; + Pout << "Ur = " << Ur << endl; + Pout << "ds = " << ds << endl; + Pout << "ds/scale = " << ds_scaled << endl; + Pout << "rho = " << rho << endl; + Pout << "nuf = " << nuf << endl; + Pout << "voidfraction = " << voidfraction << endl; + Pout << "Rep = " << Rep << endl; + Pout << "Cd = " << Cd << endl; + Pout << "drag = " << drag << endl; } //Set value fields and write the probe diff --git a/src/lagrangian/cfdemParticle/subModels/forceModel/SchillerNaumannDrag/SchillerNaumannDrag.H b/src/lagrangian/cfdemParticle/subModels/forceModel/SchillerNaumannDrag/SchillerNaumannDrag.H index 2fa743d5..a0f79aa6 100644 --- a/src/lagrangian/cfdemParticle/subModels/forceModel/SchillerNaumannDrag/SchillerNaumannDrag.H +++ b/src/lagrangian/cfdemParticle/subModels/forceModel/SchillerNaumannDrag/SchillerNaumannDrag.H @@ -66,6 +66,10 @@ private: const volVectorField& U_; + mutable scalar scaleDia_; + + mutable scalar scaleDrag_; + public: //- Runtime type information diff --git a/src/lagrangian/cfdemParticle/subModels/forceModel/dSauter/dSauter.C b/src/lagrangian/cfdemParticle/subModels/forceModel/dSauter/dSauter.C index 9f1a555b..88fed5d4 100644 --- a/src/lagrangian/cfdemParticle/subModels/forceModel/dSauter/dSauter.C +++ b/src/lagrangian/cfdemParticle/subModels/forceModel/dSauter/dSauter.C @@ -56,6 +56,7 @@ dSauter::dSauter multiTypes_(false), d2_(NULL), d3_(NULL), + maxTypeCG_(1), typeCG_(propsDict_.lookupOrDefault("coarseGrainingFactors",scalarList(1,1.0))), d2Field_ ( IOobject @@ -95,7 +96,11 @@ dSauter::dSauter "zeroGradient" ) { - if (typeCG_.size()>1) multiTypes_ = true; + if (typeCG_.size()>1) + { + multiTypes_ = true; + maxTypeCG_ = typeCG_.size(); + } allocateMyArrays(); dSauter_.write(); @@ -153,6 +158,10 @@ void dSauter::setForce() const if (multiTypes_) { partType = particleCloud_.particleType(index); + if (partType > maxTypeCG_) + { + FatalError<< "Too few coarse-graining factors provided." << abort(FatalError); + } cg = typeCG_[partType - 1]; } ds = particleCloud_.d(index); diff --git a/src/lagrangian/cfdemParticle/subModels/forceModel/dSauter/dSauter.H b/src/lagrangian/cfdemParticle/subModels/forceModel/dSauter/dSauter.H index f1e1ca6e..eb875ca5 100644 --- a/src/lagrangian/cfdemParticle/subModels/forceModel/dSauter/dSauter.H +++ b/src/lagrangian/cfdemParticle/subModels/forceModel/dSauter/dSauter.H @@ -49,6 +49,8 @@ private: mutable double **d3_; + label maxTypeCG_; + scalarList typeCG_; mutable volScalarField d2Field_; diff --git a/src/lagrangian/cfdemParticle/subModels/forceModel/forceModel/forceModel.C b/src/lagrangian/cfdemParticle/subModels/forceModel/forceModel/forceModel.C index 3bca651c..8c3154ae 100644 --- a/src/lagrangian/cfdemParticle/subModels/forceModel/forceModel/forceModel.C +++ b/src/lagrangian/cfdemParticle/subModels/forceModel/forceModel/forceModel.C @@ -69,7 +69,8 @@ forceModel::forceModel IOobject::AUTO_WRITE ), sm.mesh(), - dimensionedVector("zero", dimensionSet(1,1,-2,0,0), vector::zero) // N + dimensionedVector("zero", dimensionSet(1,1,-2,0,0), vector::zero), + "zeroGradient" ), expParticleForces_ ( IOobject @@ -81,7 +82,8 @@ forceModel::forceModel IOobject::AUTO_WRITE ), sm.mesh(), - dimensionedVector("zero", dimensionSet(1,1,-2,0,0), vector::zero) // N + dimensionedVector("zero", dimensionSet(1,1,-2,0,0), vector::zero), + "zeroGradient" ), coupleForce_(true), modelType_(sm.modelType()), diff --git a/src/lagrangian/cfdemParticle/subModels/forceModel/isotropicFluctuations/isotropicFluctuations.C b/src/lagrangian/cfdemParticle/subModels/forceModel/isotropicFluctuations/isotropicFluctuations.C index 6c9c5579..9b0be1df 100644 --- a/src/lagrangian/cfdemParticle/subModels/forceModel/isotropicFluctuations/isotropicFluctuations.C +++ b/src/lagrangian/cfdemParticle/subModels/forceModel/isotropicFluctuations/isotropicFluctuations.C @@ -1,6 +1,6 @@ /*---------------------------------------------------------------------------*\ CFDEMcoupling academic - Open Source CFD-DEM coupling - + Contributing authors: Thomas Lichtenegger Copyright (C) 2015- Johannes Kepler University, Linz @@ -82,8 +82,9 @@ isotropicFluctuations::isotropicFluctuations sm.mesh(), dimensionedScalar("D0", dimensionSet(0,0,0,0,0,0,0), D0_) ), + maxDisplacement_(propsDict_.lookupOrDefault("maxDisplacement", -1.0)), dtDEM_(particleCloud_.dataExchangeM().DEMts()), - ranGen_(osRandomInteger()) + ranGen_(clock::getTime()+pid()) { if(ignoreCellsName_ != "none") { @@ -122,10 +123,13 @@ void isotropicFluctuations::setForce() const vector flucU(0,0,0); label cellI=0; scalar relVolfractionExcess(0.0); - + interpolationCellPoint voidfractionInterpolator_(voidfraction_); interpolationCellPoint voidfractionRecInterpolator_(voidfractionRec_); - + + scalar maxVel = maxDisplacement_ / dtDEM_; + scalar magFlucU = 0.0; + for(int index = 0;index < particleCloud_.numberOfParticles(); ++index) { cellI = particleCloud_.cellIDs()[index][0]; @@ -151,14 +155,22 @@ void isotropicFluctuations::setForce() const voidfractionRec = voidfractionRec_[cellI]; } // write particle based data to global array - + deltaVoidfrac=voidfractionRec-voidfraction; relVolfractionExcess=deltaVoidfrac/(1-voidfraction+SMALL); if(deltaVoidfrac>0) { D = D0Field_[cellI]; - flucU=unitFlucDir()*fluctuationMag(relVolfractionExcess,D); - } + magFlucU = fluctuationMag(relVolfractionExcess,D); + if (maxVel > 0.0 && magFlucU > maxVel) + { + flucU=unitFlucDir()*maxVel; + } + else + { + flucU=unitFlucDir()*fluctuationMag(relVolfractionExcess,D); + } + } // write particle based data to global array for(int j=0;j<3;j++) @@ -168,11 +180,11 @@ void isotropicFluctuations::setForce() const } } } - + if (measureDiff_) { dimensionedScalar diff( fvc::domainIntegrate( sqr( voidfraction_ - voidfractionRec_ ) ) ); - scalar t = particleCloud_.mesh().time().timeOutputValue(); + scalar t = particleCloud_.mesh().time().timeOutputValue(); recErrorFile_ << t << "\t" << diff.value() << endl; } } diff --git a/src/lagrangian/cfdemParticle/subModels/forceModel/isotropicFluctuations/isotropicFluctuations.H b/src/lagrangian/cfdemParticle/subModels/forceModel/isotropicFluctuations/isotropicFluctuations.H index 47c74e58..360b128f 100644 --- a/src/lagrangian/cfdemParticle/subModels/forceModel/isotropicFluctuations/isotropicFluctuations.H +++ b/src/lagrangian/cfdemParticle/subModels/forceModel/isotropicFluctuations/isotropicFluctuations.H @@ -1,6 +1,6 @@ /*---------------------------------------------------------------------------*\ CFDEMcoupling academic - Open Source CFD-DEM coupling - + Contributing authors: Thomas Lichtenegger Copyright (C) 2015- Johannes Kepler University, Linz @@ -46,24 +46,24 @@ class isotropicFluctuations { protected: dictionary propsDict_; - + bool interpolate_; - + bool measureDiff_; - + mutable OFstream recErrorFile_; // ignore fluctuations in region word ignoreCellsName_; autoPtr ignoreCells_; - + bool existIgnoreCells_; word voidfractionFieldName_; - + const volScalarField& voidfraction_; - + word voidfractionRecFieldName_; const volScalarField& voidfractionRec_; @@ -76,6 +76,8 @@ protected: volScalarField D0Field_; + scalar maxDisplacement_; + const scalar dtDEM_; virtual vector unitFlucDir() const; diff --git a/src/lagrangian/cfdemParticle/subModels/forceModel/particleCellVolume/particleCellVolume.C b/src/lagrangian/cfdemParticle/subModels/forceModel/particleCellVolume/particleCellVolume.C index e09af175..c8aca4b0 100644 --- a/src/lagrangian/cfdemParticle/subModels/forceModel/particleCellVolume/particleCellVolume.C +++ b/src/lagrangian/cfdemParticle/subModels/forceModel/particleCellVolume/particleCellVolume.C @@ -94,8 +94,17 @@ particleCellVolume::particleCellVolume ), upperThreshold_(readScalar(propsDict_.lookup("upperThreshold"))), lowerThreshold_(readScalar(propsDict_.lookup("lowerThreshold"))), - verbose_(propsDict_.found("verbose")) + verbose_(propsDict_.found("verbose")), + writeToFile_(propsDict_.lookupOrDefault("writeToFile",false)), + filePtr_() { + // create the path and output file + if(writeToFile_) + { + fileName path(particleCloud_.IOM().createTimeDir("postProcessing/particleCellVolume")); + filePtr_.set(new OFstream(path/"particleCellVolume.txt")); + filePtr_() << "# time | total particle volume in cells | total volume of cells with particles | average volume fraction | min(voidfraction) | max(voidfraction)" << endl; + } } @@ -120,6 +129,8 @@ void particleCellVolume::setForce() const scalar fieldValue=-1; scalar cellVol=-1; + scalar minFieldVal=1e18; + scalar maxFieldVal=-1e18; forAll(field,cellI) { @@ -129,6 +140,8 @@ void particleCellVolume::setForce() const cellVol = mesh_.V()[cellI]; scalarField_[cellI] = (1-fieldValue) * cellVol; scalarField2_[cellI] = cellVol; + minFieldVal = min(minFieldVal, fieldValue); + maxFieldVal = max(maxFieldVal, fieldValue); } else { @@ -138,6 +151,8 @@ void particleCellVolume::setForce() const } scalarField_.ref() = gSum(scalarField_); scalarField2_.ref() = gSum(scalarField2_); + reduce(minFieldVal, minOp()); + reduce(maxFieldVal, maxOp()); if(verbose_) { @@ -147,8 +162,20 @@ void particleCellVolume::setForce() const << ", and > " << lowerThreshold_ << ",\n the total volume of cells holding particles = " << scalarField2_[0] << ",\n this results in an average volume fraction of:" << scalarField_[0]/(scalarField2_[0]+SMALL) + << ",\n the min occurring " << scalarFieldName_ << " is:" << minFieldVal + << ",\n the max occurring " << scalarFieldName_ << " is:" << maxFieldVal << endl; } + + if(writeToFile_) + { + filePtr_() << mesh_.time().value() << " " + << scalarField_[0] << " " + << scalarField2_[0] << " " + << scalarField_[0]/(scalarField2_[0]+SMALL) << " " + << minFieldVal << " " + << maxFieldVal << endl; + } }// end if time >= startTime_ } diff --git a/src/lagrangian/cfdemParticle/subModels/forceModel/particleCellVolume/particleCellVolume.H b/src/lagrangian/cfdemParticle/subModels/forceModel/particleCellVolume/particleCellVolume.H index b70f5329..f0bf2b94 100644 --- a/src/lagrangian/cfdemParticle/subModels/forceModel/particleCellVolume/particleCellVolume.H +++ b/src/lagrangian/cfdemParticle/subModels/forceModel/particleCellVolume/particleCellVolume.H @@ -74,6 +74,10 @@ private: const Switch verbose_; + const Switch writeToFile_; + + mutable autoPtr filePtr_; + public: //- Runtime type information diff --git a/src/lagrangian/cfdemParticle/subModels/forceModel/particleDeformation/particleDeformation.C b/src/lagrangian/cfdemParticle/subModels/forceModel/particleDeformation/particleDeformation.C index ac5f9446..ed1444c7 100644 --- a/src/lagrangian/cfdemParticle/subModels/forceModel/particleDeformation/particleDeformation.C +++ b/src/lagrangian/cfdemParticle/subModels/forceModel/particleDeformation/particleDeformation.C @@ -20,6 +20,7 @@ License #include "particleDeformation.H" #include "addToRunTimeSelectionTable.H" #include "dataExchangeModel.H" +#include "OFstream.H" // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // @@ -52,9 +53,13 @@ particleDeformation::particleDeformation initialExec_(true), refFieldName_(propsDict_.lookup("refFieldName")), refField_(), - partType_(propsDict_.lookupOrDefault