From ac4fa1374587b1a31ea845604978125f38a79e8e Mon Sep 17 00:00:00 2001 From: Thomas Lichtenegger Date: Mon, 4 Sep 2017 13:38:52 +0200 Subject: [PATCH 1/4] Monitor mass changes. --- .../cfdemSolverRhoPimpleChem/cfdemSolverRhoPimpleChem.C | 8 +++----- .../solvers/cfdemSolverRhoPimpleChem/monitorMass.H | 7 +++++++ 2 files changed, 10 insertions(+), 5 deletions(-) create mode 100644 applications/solvers/cfdemSolverRhoPimpleChem/monitorMass.H diff --git a/applications/solvers/cfdemSolverRhoPimpleChem/cfdemSolverRhoPimpleChem.C b/applications/solvers/cfdemSolverRhoPimpleChem/cfdemSolverRhoPimpleChem.C index c769134e..a7482e81 100644 --- a/applications/solvers/cfdemSolverRhoPimpleChem/cfdemSolverRhoPimpleChem.C +++ b/applications/solvers/cfdemSolverRhoPimpleChem/cfdemSolverRhoPimpleChem.C @@ -137,17 +137,15 @@ int main(int argc, char *argv[]) #include "molConc.H" #include "pEqn.H" } - m=gSum(rhoeps*1.0*rhoeps.mesh().V()); - if(counter==0) m0=m; - counter++; - Info << "\ncurrent gas mass = " << m << "\n" << endl; - Info << "\ncurrent added gas mass = " << m-m0 << "\n" << endl; + if (pimple.turbCorr()) { turbulence->correct(); } } + #include "monitorMass.H" + particleCloud.clockM().start(31,"postFlow"); particleCloud.postFlow(); particleCloud.clockM().stop("postFlow"); diff --git a/applications/solvers/cfdemSolverRhoPimpleChem/monitorMass.H b/applications/solvers/cfdemSolverRhoPimpleChem/monitorMass.H new file mode 100644 index 00000000..5940864f --- /dev/null +++ b/applications/solvers/cfdemSolverRhoPimpleChem/monitorMass.H @@ -0,0 +1,7 @@ +{ + m=gSum(rhoeps*1.0*rhoeps.mesh().V()); + if(counter==0) m0=m; + counter++; + Info << "\ncurrent gas mass = " << m << "\n" << endl; + Info << "\ncurrent added gas mass = " << m-m0 << "\n" << endl; +} From 4bd7f6529e881ee4388ccdf6e7d08322a5f5b89a Mon Sep 17 00:00:00 2001 From: Thomas Lichtenegger Date: Mon, 4 Sep 2017 13:40:16 +0200 Subject: [PATCH 2/4] Fix original pressure equation, add alternative formulation. --- .../solvers/cfdemSolverRhoPimpleChem/pEqn.H | 22 +--- .../pEqn_alternative.H | 109 ++++++++++++++++++ 2 files changed, 114 insertions(+), 17 deletions(-) create mode 100644 applications/solvers/cfdemSolverRhoPimpleChem/pEqn_alternative.H diff --git a/applications/solvers/cfdemSolverRhoPimpleChem/pEqn.H b/applications/solvers/cfdemSolverRhoPimpleChem/pEqn.H index cd843446..5d73d7fe 100644 --- a/applications/solvers/cfdemSolverRhoPimpleChem/pEqn.H +++ b/applications/solvers/cfdemSolverRhoPimpleChem/pEqn.H @@ -3,12 +3,6 @@ rho = max(rho, rhoMin); rho = min(rho, rhoMax); rho.relax(); -rhoeps = rho * voidfraction; - -// Thermodynamic density needs to be updated by psi*d(p) after the -// pressure solution - done in 2 parts. Part 1: -thermo.rho() -= psi*p; - volScalarField rAU(1.0/UEqn.A()); surfaceScalarField rhorAUf("rhorAUf", fvc::interpolate(rhoeps*rAU)); if (modelType=="A") @@ -44,22 +38,20 @@ else // Update the pressure BCs to ensure flux consistency constrainPressure(p, rhoeps, U, phi, rhorAUf); - + volScalarField SmbyP(particleCloud.chemistryM(0).Sm() / p); + while (pimple.correctNonOrthogonal()) { // Pressure corrector fvScalarMatrix pEqn ( - // fvm::ddt(psi*voidfraction, p) - fvc::ddt(rhoeps) + psi*correction(fvm::ddt(voidfraction,p)) + fvm::ddt(voidfraction, psi, p) + fvc::div(phi) - fvm::laplacian(rhorAUf, p) == - // particleCloud.chemistryM(0).Sm() fvm::Sp(SmbyP, p) + fvOptions(psi, p, rho.name()) - ); pEqn.solve(mesh.solver(p.select(pimple.finalInnerIter()))); @@ -76,21 +68,17 @@ else // Explicitly relax pressure for momentum corrector p.relax(); -// Second part of thermodynamic density update -thermo.rho() += psi*p; // Recalculate density from the relaxed pressure rho = thermo.rho(); - rho = max(rho, rhoMin); rho = min(rho, rhoMax); rho.relax(); - -rhoeps = rho * voidfraction; - Info<< "rho max/min : " << max(rho).value() << " " << min(rho).value() << endl; +rhoeps = rho * voidfraction; + if (modelType=="A") { U = HbyA - rAU*(voidfraction*fvc::grad(p)-Ksl*Us); diff --git a/applications/solvers/cfdemSolverRhoPimpleChem/pEqn_alternative.H b/applications/solvers/cfdemSolverRhoPimpleChem/pEqn_alternative.H new file mode 100644 index 00000000..cd843446 --- /dev/null +++ b/applications/solvers/cfdemSolverRhoPimpleChem/pEqn_alternative.H @@ -0,0 +1,109 @@ +rho = thermo.rho(); +rho = max(rho, rhoMin); +rho = min(rho, rhoMax); +rho.relax(); + +rhoeps = rho * voidfraction; + +// Thermodynamic density needs to be updated by psi*d(p) after the +// pressure solution - done in 2 parts. Part 1: +thermo.rho() -= psi*p; + +volScalarField rAU(1.0/UEqn.A()); +surfaceScalarField rhorAUf("rhorAUf", fvc::interpolate(rhoeps*rAU)); +if (modelType=="A") +{ + rhorAUf *= fvc::interpolate(voidfraction); +} +volVectorField HbyA(constrainHbyA(rAU*UEqn.H(), U, p)); + +surfaceScalarField phiUs("phiUs", fvc::interpolate(rhoeps*rAU*Ksl*Us)& mesh.Sf()); + +if (pimple.nCorrPISO() <= 1) +{ + tUEqn.clear(); +} + +if (pimple.transonic()) +{ +// transonic version not implemented yet +} +else +{ + surfaceScalarField phiHbyA + ( + "phiHbyA", + ( + fvc::flux(rhoeps*HbyA) + // + rhorAUf*fvc::ddtCorr(rho, U, phi) + ) + ); + + // flux without pressure gradient contribution + phi = phiHbyA + phiUs; + + // Update the pressure BCs to ensure flux consistency + constrainPressure(p, rhoeps, U, phi, rhorAUf); + + volScalarField SmbyP(particleCloud.chemistryM(0).Sm() / p); + while (pimple.correctNonOrthogonal()) + { + // Pressure corrector + fvScalarMatrix pEqn + ( + // fvm::ddt(psi*voidfraction, p) + fvc::ddt(rhoeps) + psi*correction(fvm::ddt(voidfraction,p)) + + fvc::div(phi) + - fvm::laplacian(rhorAUf, p) + == + // particleCloud.chemistryM(0).Sm() + fvm::Sp(SmbyP, p) + + fvOptions(psi, p, rho.name()) + + ); + + pEqn.solve(mesh.solver(p.select(pimple.finalInnerIter()))); + + if (pimple.finalNonOrthogonalIter()) + { + phi += pEqn.flux(); + } + } +} + +#include "rhoEqn.H" +#include "compressibleContinuityErrsPU.H" + +// Explicitly relax pressure for momentum corrector +p.relax(); +// Second part of thermodynamic density update +thermo.rho() += psi*p; + +// Recalculate density from the relaxed pressure +rho = thermo.rho(); + +rho = max(rho, rhoMin); +rho = min(rho, rhoMax); +rho.relax(); + +rhoeps = rho * voidfraction; + +Info<< "rho max/min : " << max(rho).value() + << " " << min(rho).value() << endl; + +if (modelType=="A") +{ + U = HbyA - rAU*(voidfraction*fvc::grad(p)-Ksl*Us); +} +else +{ + U = HbyA - rAU*(fvc::grad(p)-Ksl*Us); +} +U.correctBoundaryConditions(); +fvOptions.correct(U); +K = 0.5*magSqr(U); + +if (thermo.dpdt()) +{ + dpdt = fvc::ddt(voidfraction,p); +} From f14ce1b7d207603daec39992f192563982c15281 Mon Sep 17 00:00:00 2001 From: Thomas Lichtenegger Date: Mon, 4 Sep 2017 13:40:55 +0200 Subject: [PATCH 3/4] Output mass source in separate line. --- .../cfdemParticle/subModels/chemistryModel/species/species.C | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/lagrangian/cfdemParticle/subModels/chemistryModel/species/species.C b/src/lagrangian/cfdemParticle/subModels/chemistryModel/species/species.C index fe2ce18b..e64bdad4 100644 --- a/src/lagrangian/cfdemParticle/subModels/chemistryModel/species/species.C +++ b/src/lagrangian/cfdemParticle/subModels/chemistryModel/species/species.C @@ -334,7 +334,7 @@ void species::execute() } massSourceCurr_ = gSum(changeOfGasMassField_*1.0*changeOfGasMassField_.mesh().V() * Nevery_ * timestep); massSourceTot_ += massSourceCurr_; - Info << "total conversion of mass:\tcurrent source = " << massSourceCurr_ << ", total source = " << massSourceTot_ << "\n" << endl; + Info << "total conversion of mass:\n\tcurrent source = " << massSourceCurr_ << "\n\ttotal source = " << massSourceTot_ << "\n" << endl; Info << "get data done" << endl; } From 4e7246381dbe931e2f972d5a79b30c24feadb3da Mon Sep 17 00:00:00 2001 From: Thomas Lichtenegger Date: Mon, 4 Sep 2017 13:42:53 +0200 Subject: [PATCH 4/4] Some changes for tutorial case with single shrinking particle. --- .../CFD/0/CO | 2 +- .../CFD/0/O2 | 2 +- .../Shrinking_particle_model_multispecies/CFD/0/T | 15 ++++++++++----- .../CFD/constant/couplingProperties | 2 +- .../CFD/system/controlDict | 5 +++-- .../CFD/system/fvOptions | 4 ++-- .../CFD/system/fvSolution | 10 +++++----- .../CFD/system/probesDict | 2 +- .../DEM/in.liggghts_run | 10 +++++----- .../parCFDDEMrun.sh | 2 +- 10 files changed, 30 insertions(+), 24 deletions(-) diff --git a/tutorials/cfdemSolverRhoPimpleChem/Shrinking_particle_model_multispecies/CFD/0/CO b/tutorials/cfdemSolverRhoPimpleChem/Shrinking_particle_model_multispecies/CFD/0/CO index 1d8cec8a..e655f71b 100644 --- a/tutorials/cfdemSolverRhoPimpleChem/Shrinking_particle_model_multispecies/CFD/0/CO +++ b/tutorials/cfdemSolverRhoPimpleChem/Shrinking_particle_model_multispecies/CFD/0/CO @@ -17,7 +17,7 @@ FoamFile dimensions [0 0 0 0 0 0 0]; -internalField uniform 0.3; +internalField uniform 0.0; boundaryField { diff --git a/tutorials/cfdemSolverRhoPimpleChem/Shrinking_particle_model_multispecies/CFD/0/O2 b/tutorials/cfdemSolverRhoPimpleChem/Shrinking_particle_model_multispecies/CFD/0/O2 index d1a9870f..6b0f56c8 100644 --- a/tutorials/cfdemSolverRhoPimpleChem/Shrinking_particle_model_multispecies/CFD/0/O2 +++ b/tutorials/cfdemSolverRhoPimpleChem/Shrinking_particle_model_multispecies/CFD/0/O2 @@ -17,7 +17,7 @@ FoamFile dimensions [0 0 0 0 0 0 0]; -internalField uniform 0.7; +internalField uniform 1.0; boundaryField { diff --git a/tutorials/cfdemSolverRhoPimpleChem/Shrinking_particle_model_multispecies/CFD/0/T b/tutorials/cfdemSolverRhoPimpleChem/Shrinking_particle_model_multispecies/CFD/0/T index dc46d191..03bc87ab 100644 --- a/tutorials/cfdemSolverRhoPimpleChem/Shrinking_particle_model_multispecies/CFD/0/T +++ b/tutorials/cfdemSolverRhoPimpleChem/Shrinking_particle_model_multispecies/CFD/0/T @@ -23,27 +23,32 @@ boundaryField { top { - type zeroGradient; + type fixedValue; + value uniform 1293.15; } bottom { - type zeroGradient; + type fixedValue; + value uniform 1293.15; } side-walls { - type zeroGradient; + type fixedValue; + value uniform 1293.15; } inlet { - type zeroGradient; + type fixedValue; + value uniform 1293.15; } outlet { - type zeroGradient; + type fixedValue; + value uniform 1293.15; } } diff --git a/tutorials/cfdemSolverRhoPimpleChem/Shrinking_particle_model_multispecies/CFD/constant/couplingProperties b/tutorials/cfdemSolverRhoPimpleChem/Shrinking_particle_model_multispecies/CFD/constant/couplingProperties index c286cb4e..6684af60 100644 --- a/tutorials/cfdemSolverRhoPimpleChem/Shrinking_particle_model_multispecies/CFD/constant/couplingProperties +++ b/tutorials/cfdemSolverRhoPimpleChem/Shrinking_particle_model_multispecies/CFD/constant/couplingProperties @@ -31,7 +31,7 @@ FoamFile modelType "A"; // A or B -couplingInterval 100;//1000; +couplingInterval 10;//1000; voidFractionModel divided;//centre;// diff --git a/tutorials/cfdemSolverRhoPimpleChem/Shrinking_particle_model_multispecies/CFD/system/controlDict b/tutorials/cfdemSolverRhoPimpleChem/Shrinking_particle_model_multispecies/CFD/system/controlDict index e0c9228f..c7d9d501 100644 --- a/tutorials/cfdemSolverRhoPimpleChem/Shrinking_particle_model_multispecies/CFD/system/controlDict +++ b/tutorials/cfdemSolverRhoPimpleChem/Shrinking_particle_model_multispecies/CFD/system/controlDict @@ -25,11 +25,11 @@ stopAt endTime; endTime 5.0; -deltaT 0.0001; +deltaT 0.00001; writeControl timeStep; -writeInterval 50; +writeInterval 1000; purgeWrite 0; @@ -104,6 +104,7 @@ functions ( rhoeps rho + molarConc ); } diff --git a/tutorials/cfdemSolverRhoPimpleChem/Shrinking_particle_model_multispecies/CFD/system/fvOptions b/tutorials/cfdemSolverRhoPimpleChem/Shrinking_particle_model_multispecies/CFD/system/fvOptions index ee7ab4cb..b47a13ee 100644 --- a/tutorials/cfdemSolverRhoPimpleChem/Shrinking_particle_model_multispecies/CFD/system/fvOptions +++ b/tutorials/cfdemSolverRhoPimpleChem/Shrinking_particle_model_multispecies/CFD/system/fvOptions @@ -22,7 +22,7 @@ source1 { active yes; selectionMode all; - Tmin 1293.1; - Tmax 1293.2; + Tmin 1293.15; + Tmax 2000; } } diff --git a/tutorials/cfdemSolverRhoPimpleChem/Shrinking_particle_model_multispecies/CFD/system/fvSolution b/tutorials/cfdemSolverRhoPimpleChem/Shrinking_particle_model_multispecies/CFD/system/fvSolution index 5ff3e6d2..2967bf16 100644 --- a/tutorials/cfdemSolverRhoPimpleChem/Shrinking_particle_model_multispecies/CFD/system/fvSolution +++ b/tutorials/cfdemSolverRhoPimpleChem/Shrinking_particle_model_multispecies/CFD/system/fvSolution @@ -58,14 +58,14 @@ solvers { solver smoothSolver; smoother symGaussSeidel; - tolerance 1e-05; + tolerance 1e-06; relTol 0.1; } "(U|h|R|k|epsilon)Final" { $U; - tolerance 1e-05; + tolerance 1e-07; relTol 0; } @@ -77,7 +77,7 @@ solvers "(Yi|CO2|O2)Final" { $Yi; - tolerance 1e-06; + tolerance 1e-08; relTol 0; } } @@ -98,13 +98,13 @@ relaxationFactors { fields { - ".*" 1; + p 0.3; } equations { ".*" 1; } } - +*/ // ************************************************************************* // diff --git a/tutorials/cfdemSolverRhoPimpleChem/Shrinking_particle_model_multispecies/CFD/system/probesDict b/tutorials/cfdemSolverRhoPimpleChem/Shrinking_particle_model_multispecies/CFD/system/probesDict index 5c180b82..5a78f453 100644 --- a/tutorials/cfdemSolverRhoPimpleChem/Shrinking_particle_model_multispecies/CFD/system/probesDict +++ b/tutorials/cfdemSolverRhoPimpleChem/Shrinking_particle_model_multispecies/CFD/system/probesDict @@ -21,7 +21,7 @@ fields rho p T - molC + molarConc O2 O CO2 diff --git a/tutorials/cfdemSolverRhoPimpleChem/Shrinking_particle_model_multispecies/DEM/in.liggghts_run b/tutorials/cfdemSolverRhoPimpleChem/Shrinking_particle_model_multispecies/DEM/in.liggghts_run index 110f436d..bd8992d3 100644 --- a/tutorials/cfdemSolverRhoPimpleChem/Shrinking_particle_model_multispecies/DEM/in.liggghts_run +++ b/tutorials/cfdemSolverRhoPimpleChem/Shrinking_particle_model_multispecies/DEM/in.liggghts_run @@ -45,21 +45,21 @@ fix zwalls2 all wall/gran model hertz tangential history primitive type 1 zp ############################################### # cfd coupling -fix cfd all couple/cfd couple_every 100 mpi +fix cfd all couple/cfd couple_every 10 mpi fix cfd2 all couple/cfd/force # this should invoke chemistry fix cfd3 all couple/cfd/chemistry n_species 5 species_names O2 CO2 N2 CO O # this should shrink the particle -#fix cfd4 all chem/shrink speciesA O2 molMassA 31.99 speciesC CO2 molMassC 44.01 molMassB 12.01 k 2.5e1 rmin 0.005 nevery 100 screen yes -#fix cfd5 all chem/shrink speciesA CO2 molMassA 44.01 speciesC CO molMassC 28.01 nuC 2 molMassB 12.01 k 2.5e1 rmin 0.005 nevery 100 +fix cfd4 all chem/shrink speciesA O2 molMassA 31.99 speciesC CO2 molMassC 44.01 molMassB 12.01 k 0.1 rmin 0.005 nevery 10 screen yes +fix cfd5 all chem/shrink speciesA CO2 molMassA 44.01 speciesC CO molMassC 28.01 nuC 2 molMassB 12.01 k 0.02 rmin 0.005 nevery 10 screen yes # values for k0 and T0 according to Shen et al., Fuel (2011) with a coke density of 800 kg / m^3 # C + O2 -> CO2 k0 = 4e3 m/s T0 = 10855 K # C + CO2 -> 2CO k0 = 6e6 m/s T0 = 29018 K -fix cfd4 all chem/shrink/Arrhenius speciesA O2 molMassA 31.99 speciesC CO2 molMassC 44.01 molMassB 12.01 rmin 0.005 nevery 100 screen yes k 4.0e3 T 10855 -fix cfd5 all chem/shrink/Arrhenius speciesA CO2 molMassA 44.01 speciesC CO molMassC 28.01 nuC 2 molMassB 12.01 rmin 0.005 nevery 100 screen yes k 6e6 T 29018 +#fix cfd4 all chem/shrink/Arrhenius speciesA O2 molMassA 31.99 speciesC CO2 molMassC 44.01 molMassB 12.01 rmin 0.005 nevery 10 screen yes k 4.0e3 T 10855 +#fix cfd5 all chem/shrink/Arrhenius speciesA CO2 molMassA 44.01 speciesC CO molMassC 28.01 nuC 2 molMassB 12.01 rmin 0.005 nevery 10 screen yes k 6e6 T 29018 # apply nve integration to all particles that are inserted as single particles diff --git a/tutorials/cfdemSolverRhoPimpleChem/Shrinking_particle_model_multispecies/parCFDDEMrun.sh b/tutorials/cfdemSolverRhoPimpleChem/Shrinking_particle_model_multispecies/parCFDDEMrun.sh index 6674ef7f..28b95f78 100755 --- a/tutorials/cfdemSolverRhoPimpleChem/Shrinking_particle_model_multispecies/parCFDDEMrun.sh +++ b/tutorials/cfdemSolverRhoPimpleChem/Shrinking_particle_model_multispecies/parCFDDEMrun.sh @@ -10,7 +10,7 @@ . ~/.bashrc #- include functions -source $CFDEM_SRC_DIR/lagrangian/cfdemParticle/etc/functions.sh +source $CFDEM_PROJECT_DIR/etc/functions.sh #--------------------------------------------------------------------------------# #- define variables