- 修改配置文件，使不输出警告；

- 测试案例默认不调用octave进行后处理，默认启用clean

.gitignore

@@ -1,6 +1,8 @@
 *.o
 *.dep
-log_*
-log.*
+log*
+.vscode
+platforms
 *~
+testCase*
 

README

@@ -1,107 +0,0 @@
-/*---------------------------------------------------------------------------*\
-    CFDEMcoupling - Open Source CFD-DEM coupling
-
-    CFD-DEM coupling toolbox, released by 
-    DCS Computing GmbH, Linz, Austria
-    www.dcs-computing.com, office@dcs-computing.com
-
-    CFDEMcoupling is part of the CFDEMproject
-    www.cfdem.com
-
-    CFDEMcoupling is free software; you can redistribute it and/or modify it
-    under the terms of the GNU General Public License as published by the
-    Free Software Foundation; either version 3 of the License, or (at your
-    option) any later version.
-
-    CFDEMcoupling is part of CFDEMproject: 
-    www.liggghts.com | www.cfdem.com
-
-    Core developer and main author:
-    Christoph Goniva, christoph.goniva@dcs-computing.com
-
-    LIGGGHTS® and CFDEM® are registered trade marks of DCS Computing GmbH, 
-    the producer of the LIGGGHTS® software and the CFDEM®coupling software
-    See http://www.cfdem.com/terms-trademark-policy for details.
-/*---------------------------------------------------------------------------*\
-Copyright
-
-    Copyright 2012-         DCS Computing GmbH, Linz
-    Copyright 2009-2015     JKU Linz
-
-    Some parts of CFDEMcoupling are based on OpenFOAM® and Copyright on these
-    parts is held by the OpenFOAM® Foundation (www.openFoam.org)
-    and potentially other parties.
-    Some parts of CFDEMcoupling are contributied by other parties, which are
-    holding the Copyright. This is listed in each file of the distribution.
-/*---------------------------------------------------------------------------*\
-License
-
-    This file is part of CFDEMcoupling.
-
-    CFDEMcoupling is free software; you can redistribute it and/or modify it
-    under the terms of the GNU General Public License as published by the
-    Free Software Foundation; either version 3 of the License, or (at your
-    option) any later version.
-
-    CFDEMcoupling is distributed in the hope that it will be useful, but WITHOUT
-    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
-    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
-    for more details.
-
-    You should have received a copy of the GNU General Public License
-    along with CFDEMcoupling; if not, write to the Free Software Foundation,
-    Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
-/*---------------------------------------------------------------------------*\
-Description
-    This code provides models and solvers to realize coupled CFD-DEM simulations
-    using LIGGGHTS and OpenFOAM® technology.
-    Note:   this code is not part of OpenFOAM®.
-            This offering is not approved or endorsed by OpenCFD Limited, 
-            the producer of the OpenFOAM® software and owner of the OPENFOAM®
-            and OpenCFD®  trade marks.
-\*---------------------------------------------------------------------------*/
-
-
-CFDEM® coupling provides an open source parallel coupled CFD-DEM framework 
-combining the strengths of LIGGGHTS® DEM code and the Open Source 
-CFD package OpenFOAM®(*). The CFDEM®coupling toolbox allows to expand 
-standard CFD solvers of OpenFOAM®(*) to include a coupling to the DEM 
-code LIGGGHTS®. In this toolbox the particle representation within the 
-CFD solver is organized by "cloud" classes. Key functionalities are organised 
-in sub-models (e.g. force models, data exchange models, etc.), which can easily 
-be selected and combined by dictionary settings.
-
-The coupled solvers run fully parallel on distributed-memory clusters. 
-
-Features are:
-
-- its modular approach allows users to easily implement new models
-- its MPI parallelization enables to use it for large scale problems
-- the "forum" (www.cfdem.com) on CFD-DEM gives the possibility 
-  to exchange with other users / developers
-- the use of GIT allows to easily update to the latest version
-- basic documentation is provided
-
-The CFDEM®coupling distribution includes the following files and directories:
-
-- "README" file (this file)
-- "COPYING" file the GNU General Public License (GPL)
-- "DISCLAIMER" file
-- "src" directory including the source files of the coupling toolbox and models
-- "applications" directory including the solver files for coupled CFD-DEM simulations
-- "doc" directory including the documentation of CFDEM®coupling
-- "tutorials" directory including basic tutorial cases showing the functionality
-
-Details on installation are given on the "www.cfdem.com" website.
-
-The functionality of this CFD-DEM framwork is described via "tutorial cases" showing 
-how to use different solvers and models.
-
-CFDEM®coupling stands for Computational Fluid Dynamics (CFD) -
-Discrete Element Method (DEM) coupling.
-
-\*---------------------------------------------------------------------------*/
-(*) This offering is not approved or endorsed by OpenCFD Limited,
-    the producer of the OpenFOAM software and 
-    owner of the OPENFOAM® and OpenCFD® trade marks.
-\*---------------------------------------------------------------------------*/

README.md

@@ -0,0 +1,88 @@
+<div align ="center">
+
+<p float="left">
+  <img src="doc/cfdem_logo.jpg" width="300" />
+  <img src="doc/cfdem_logo2.png" width="300" /> 
+</p>
+  
+</div>
+
+## The extended version of CFDEM®coupling is constantly developed and maintained by DCS Computing GmbH - you can find more information [here](https://www.aspherix-dem.com/software/cfdemcoupling-functionalities/)
+
+This repository contains older, publicly available versions of CFDEM®coupling (referred to as "CFDEM®coupling-PUBLIC") and will not be updated in the foreseeable future.
+
+CFDEM®coupling-PUBLIC is compatible with [LIGGGHTS-PUBLIC](https://github.com/CFDEMproject/LIGGGHTS-PUBLIC).
+
+The extended version of CFDEM®coupling is compatible with the DEM software [Aspherix®](https://www.aspherix-dem.com/), the "successor" of LIGGGHTS® developed and distributed by DCS Computing GmbH, Linz, Austria.
+
+</div>
+
+## CFDEM®coupling-PUBLIC
+
+CFDEM®coupling-PUBLIC is an Open Source coupled CFD-DEM framework combining the strengths of LIGGGHTS® DEM code and the Open Source CFD package OpenFOAM® [^1]. CFDEM®coupling stands for Computational Fluid Dynamics (CFD) - Discrete Element Method (DEM) coupling.
+
+LIGGGHTS® and CFDEM® are registered trade marks of DCS Computing GmbH, the producer of the LIGGGHTS® software and the CFDEM®coupling software; see http://www.cfdem.com/terms-trademark-policy for details.
+
+CFDEM®coupling-PUBLIC  is open-source, distributed under the terms of the GNU Public License, version 3 or later. CFDEM®coupling-PUBLIC is part of CFDEM®project: www.liggghts.com | www.cfdem.com. Core developer and main author: Christoph Goniva, christoph.goniva@dcs-computing.com
+
+[^1]: This offering is not approved or endorsed by OpenCFD Limited, the producer of the OpenFOAM software and owner of the OPENFOAM® and OpenCFD® trade marks.
+
+</div>
+
+## Features
+
+* Modular approach that allows users to easily implement new models
+* MPI parallelization that enables to use it for large scale problems
+* The "forum" (www.cfdem.com) on CFD-DEM gives the possibility to exchange with other users / developers
+* The use of GIT allows to easily update to the latest version
+* Basic documentation is provided
+
+In this toolbox the particle representation within the CFD solver is organized by "cloud" classes. Key functionalities are organised in sub-models (e.g. force models, data exchange models, etc.), which can easily be selected and combined by dictionary settings.
+
+</div>
+
+## Structure
+The CFDEM®coupling-PUBLIC distribution includes the following files and directories:
+
+* `README`: this file
+* `COPYING`: the GNU General Public License (GPL) file
+* `DISCLAIMER`
+* `src`: directory including the source files of the coupling toolbox and models
+* `applications`: directory including the solver files for coupled CFD-DEM simulations
+* `doc`: directory including the documentation of CFDEM®coupling
+* `tutorials`: directory including basic tutorial cases showing the functionality
+
+Details on installation are given on the www.cfdem.com website.
+
+</div>
+
+## License
+
+CFDEM®coupling-PUBLIC 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.
+
+CFDEM®coupling is distributed in the hope that it will be useful, but **without any warranty**; without even the implied warranty of **merchantabiity** or **fitness for a particular purpouse**; see the GNU General Public License for more details.
+
+You should have received a copy of the GNU General Public License along with CFDEM®coupling-PUBLIC; if not, write to the Free Software Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA.
+
+</div>
+
+## Copyrights
+* Copyright 2012-now: DCS Computing GmbH, Linz
+* Copyright 2009-2015: JKU Linz
+
+Some parts of CFDEM®coupling are based on OpenFOAM® and Copyright on these parts is held by the OpenFOAM® Foundation (www.openFoam.org) and potentially other parties. Some parts of CFDEMcoupling are contributied by other parties, which are holding the Copyright. This is listed in each file of the distribution.
+
+## How to cite CFDEM®coupling-PUBLIC
+If you are using PhasicFlow in your research or industrial work, cite the following [article](https://doi.org/10.1016/j.partic.2012.05.002):
+```
+@article{goniva2012influence,
+  title={Influence of rolling friction on single spout fluidized bed simulation},
+  author={Goniva, Christoph and Kloss, Christoph and Deen, Niels G and Kuipers, Johannes AM and Pirker, Stefan},
+  journal={Particuology},
+  volume={10},
+  number={5},
+  pages={582--591},
+  year={2012},
+  publisher={Elsevier}
+}
+```

applications/solvers/cfdemSolverPiso/cfdemSolverPiso.C

@@ -37,209 +37,150 @@ Description
 #include "fvCFD.H"
 #include "singlePhaseTransportModel.H"
 
-#include "OFversion.H"
-#if defined(version30)
-    #include "turbulentTransportModel.H"
-    #include "pisoControl.H"
-#else
-    #include "turbulenceModel.H"
-#endif
-#if defined(versionv1606plus) || defined(version40)
-    #include "fvOptions.H"
-#else
-    #include "fvIOoptionList.H"
-#endif
+#include "turbulentTransportModel.H"
+#include "pisoControl.H"
+
+#include "fvOptions.H"
 #include "fixedFluxPressureFvPatchScalarField.H"
 #include "cfdemCloud.H"
 
-#if defined(anisotropicRotation)
-    #include "cfdemCloudRotation.H"
-#endif
-#if defined(superquadrics_flag)
-    #include "cfdemCloudRotationSuperquadric.H"
-#endif
 #include "implicitCouple.H"
 #include "clockModel.H"
 #include "smoothingModel.H"
 #include "forceModel.H"
 
-
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
 int main(int argc, char *argv[])
 {
-    #include "setRootCase.H"
-    #include "createTime.H"
-    #include "createMesh.H"
-    #if defined(version30)
-        pisoControl piso(mesh);
-        #include "createTimeControls.H"
-    #endif
-    #include "createFields.H"
-    #include "createFvOptions.H"
-    #include "initContinuityErrs.H"
+    // int i = 0;
+    // while(1){
+    //     if(i!=0){
+    //         break;
+    //     }
+    //     sleep(1);
+    // }
+    #include "setRootCase.H" // 更新listOptions; 定义args; 更新listOutput
+    #include "createTime.H"  // 定义Foam::Time runTime(Time::controlDictname, args)
+    #include "createMesh.H"  // 定义Foam::FvMesh mesh
+    pisoControl piso(mesh);
+    // #include "createTimeControls.H" // 读取adjustTimeStep, maxCo, maxDeltaT
+    #include "createFields.H"       // 初始化物理场
+    #include "createFvOptions.H"    // 定义fvOptions
+    #include "initContinuityErrs.H" // 定义cumulativeContErr
 
     // create cfdemCloud
-    #include "readGravitationalAcceleration.H"
-    #include "checkImCoupleM.H"
-    #if defined(anisotropicRotation)
-        cfdemCloudRotation particleCloud(mesh);
-    #elif defined(superquadrics_flag)
-        cfdemCloudRotationSuperquadric particleCloud(mesh);
-    #else
-        cfdemCloud particleCloud(mesh);
-    #endif
-    #include "checkModelType.H"
+    // #include "readGravitationalAcceleration.H"
+    #include "createGravity.H"                  // 定于重力加速度g
+    #include "checkImCoupleM.H"                 // 
+    cfdemCloud particleCloud(mesh);             // 定义particleCloud
+    #include "checkModelType.H"                 // CFD-DEM模型类型, referring to: Zhou et al. 2010,JFM
 
     // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-    Info<< "\nStarting time loop\n" << endl;
+    Info << "\nStarting time loop\n" << endl;
     while (runTime.loop())
     {
-        Info<< "Time = " << runTime.timeName() << nl << endl;
+        Info << "Time = " << runTime.timeName() << nl << endl;
 
-        #if defined(version30)
-            #include "readTimeControls.H"
-            #include "CourantNo.H"
-            #include "setDeltaT.H"
-        #else
-            #include "readPISOControls.H"
-            #include "CourantNo.H"
-        #endif
+        // #include "readTimeControls.H"
+        #include "CourantNo.H"
+        // #include "setDeltaT.H"
 
         // do particle stuff
-        particleCloud.clockM().start(1,"Global");
-        particleCloud.clockM().start(2,"Coupling");
-        bool hasEvolved = particleCloud.evolve(voidfraction,Us,U);
+        particleCloud.clockM().start(1, "Global");
+        particleCloud.clockM().start(2, "Coupling");
+        bool hasEvolved = particleCloud.evolve(voidfraction, Us, U);
 
-        if(hasEvolved)
+        if (hasEvolved)
         {
             particleCloud.smoothingM().smoothenAbsolutField(particleCloud.forceM(0).impParticleForces());
         }
-    
+
         Ksl = particleCloud.momCoupleM(particleCloud.registryM().getProperty("implicitCouple_index")).impMomSource();
         Ksl.correctBoundaryConditions();
 
         surfaceScalarField voidfractionf = fvc::interpolate(voidfraction);
-        phi = voidfractionf*phiByVoidfraction;
+        phi = voidfractionf * phiByVoidfraction;
 
-        //Force Checks
+        // Force Checks
         #include "forceCheckIm.H"
 
         #include "solverDebugInfo.H"
         particleCloud.clockM().stop("Coupling");
 
-        particleCloud.clockM().start(26,"Flow");
+        particleCloud.clockM().start(26, "Flow");
 
-        if(particleCloud.solveFlow())
+        if (particleCloud.solveFlow())
         {
             // Pressure-velocity PISO corrector
             {
                 // Momentum predictor
-                fvVectorMatrix UEqn
-                (
-                    fvm::ddt(voidfraction,U) - fvm::Sp(fvc::ddt(voidfraction),U)
-                  + fvm::div(phi,U) - fvm::Sp(fvc::div(phi),U)
-//                + turbulence->divDevReff(U)
-                  + particleCloud.divVoidfractionTau(U, voidfraction)
-                  ==
-                  - fvm::Sp(Ksl/rho,U)
-                  + fvOptions(U)
-                );
+                fvVectorMatrix UEqn(
+                    fvm::ddt(voidfraction, U) 
+                    - fvm::Sp(fvc::ddt(voidfraction), U) 
+                    + fvm::div(phi, U) 
+                    - fvm::Sp(fvc::div(phi), U)
+                //  + turbulence->divDevReff(U)
+                    + particleCloud.divVoidfractionTau(U, voidfraction) 
+                    == -fvm::Sp(Ksl / rho, U) 
+                    + fvOptions(U));
 
                 UEqn.relax();
                 fvOptions.constrain(UEqn);
 
-                #if defined(version30)
-                    if (piso.momentumPredictor())
-                #else
-                    if (momentumPredictor)
-                #endif
+                if (piso.momentumPredictor())
                 {
-                    if (modelType=="B" || modelType=="Bfull")
-                        solve(UEqn == - fvc::grad(p) + Ksl/rho*Us);
+                    if (modelType == "B" || modelType == "Bfull")
+                        solve(UEqn == -fvc::grad(p) + Ksl / rho * Us);
                     else
-                        solve(UEqn == - voidfraction*fvc::grad(p) + Ksl/rho*Us);
+                        solve(UEqn == -voidfraction * fvc::grad(p) + Ksl / rho * Us);
 
                     fvOptions.correct(U);
                 }
-
                 // --- PISO loop
-                #if defined(version30)
-                    while (piso.correct())
-                #else
-                    for (int corr=0; corr<nCorr; corr++)
-                #endif
+                while (piso.correct())
                 {
-                    volScalarField rUA = 1.0/UEqn.A();
+                    volScalarField rUA = 1.0 / UEqn.A();
 
                     surfaceScalarField rUAf("(1|A(U))", fvc::interpolate(rUA));
-                    volScalarField rUAvoidfraction("(voidfraction2|A(U))",rUA*voidfraction);
+                    volScalarField rUAvoidfraction("(voidfraction2|A(U))", rUA * voidfraction);
                     surfaceScalarField rUAfvoidfraction("(voidfraction2|A(U)F)", fvc::interpolate(rUAvoidfraction));
 
-                    U = rUA*UEqn.H();
+                    U = rUA * UEqn.H();
 
-                    #ifdef version23
-                        phi = ( fvc::interpolate(U) & mesh.Sf() )
-                            + rUAfvoidfraction*fvc::ddtCorr(U, phiByVoidfraction);
-                    #else
-                        phi = ( fvc::interpolate(U) & mesh.Sf() )
-                            + fvc::ddtPhiCorr(rUAvoidfraction, U, phiByVoidfraction);
-                    #endif
+                    phi = (fvc::interpolate(U) & mesh.Sf()) + rUAfvoidfraction * fvc::ddtCorr(U, phiByVoidfraction);
                     surfaceScalarField phiS(fvc::interpolate(Us) & mesh.Sf());
-                    phi += rUAf*(fvc::interpolate(Ksl/rho) * phiS);
+                    phi += rUAf * (fvc::interpolate(Ksl / rho) * phiS);
 
-                    if (modelType=="A")
-                        rUAvoidfraction = volScalarField("(voidfraction2|A(U))",rUA*voidfraction*voidfraction);
+                    if (modelType == "A")
+                        rUAvoidfraction = volScalarField("(voidfraction2|A(U))", rUA * voidfraction * voidfraction);
 
                     // Update the fixedFluxPressure BCs to ensure flux consistency
                     #include "fixedFluxPressureHandling.H"
-                    
 
                     // Non-orthogonal pressure corrector loop
-                    #if defined(version30)
-                        while (piso.correctNonOrthogonal())
-                    #else
-                        for (int nonOrth=0; nonOrth<=nNonOrthCorr; nonOrth++)
-                    #endif
+                    while (piso.correctNonOrthogonal())
                     {
                         // Pressure corrector
-                        fvScalarMatrix pEqn
-                        (
-                            fvm::laplacian(rUAvoidfraction, p) == fvc::div(voidfractionf*phi) + particleCloud.ddtVoidfraction()
-                        );
+                        fvScalarMatrix pEqn(
+                            fvm::laplacian(rUAvoidfraction, p) == fvc::div(voidfractionf * phi) + particleCloud.ddtVoidfraction());
                         pEqn.setReference(pRefCell, pRefValue);
 
-                        #if defined(version30)
-                            pEqn.solve(mesh.solver(p.select(piso.finalInnerIter())));
-                            if (piso.finalNonOrthogonalIter())
-                            {
-                                phiByVoidfraction = phi - pEqn.flux()/voidfractionf;
-                            }
-                        #else
-                            if( corr == nCorr-1 && nonOrth == nNonOrthCorr )
-                                #if defined(versionExt32)
-                                    pEqn.solve(mesh.solutionDict().solver("pFinal"));
-                                #else
-                                    pEqn.solve(mesh.solver("pFinal"));
-                                #endif
-                            else
-                                pEqn.solve();
-
-                            if (nonOrth == nNonOrthCorr)
-                            {
-                                phiByVoidfraction = phi - pEqn.flux()/voidfractionf;
-                            }
-                        #endif
+                        pEqn.solve(mesh.solver(p.select(piso.finalInnerIter())));
+                        if (piso.finalNonOrthogonalIter())
+                        {
+                            phiByVoidfraction = phi - pEqn.flux() / voidfractionf;
+                        }
 
                     } // end non-orthogonal corrector loop
 
-                    phi = voidfractionf*phiByVoidfraction;
+                    phi = voidfractionf * phiByVoidfraction;
                     #include "continuityErrorPhiPU.H"
 
-                    if (modelType=="B" || modelType=="Bfull")
-                        U -= rUA*fvc::grad(p) - Ksl/rho*Us*rUA;
+                    if (modelType == "B" || modelType == "Bfull")
+                        U -= rUA * fvc::grad(p) - Ksl / rho * Us * rUA;
                     else
-                        U -= voidfraction*rUA*fvc::grad(p) - Ksl/rho*Us*rUA;
+                        U -= voidfraction * rUA * fvc::grad(p) - Ksl / rho * Us * rUA;
 
                     U.correctBoundaryConditions();
                     fvOptions.correct(U);
@@ -249,7 +190,7 @@ int main(int argc, char *argv[])
 
             laminarTransport.correct();
             turbulence->correct();
-        }// end solveFlow
+        } // end solveFlow
         else
         {
             Info << "skipping flow solution." << endl;
@@ -257,18 +198,18 @@ int main(int argc, char *argv[])
 
         runTime.write();
 
-        Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
-            << "  ClockTime = " << runTime.elapsedClockTime() << " s"
-            << nl << endl;
+        Info << "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
+             << "  ClockTime = " << runTime.elapsedClockTime() << " s"
+             << nl << endl;
 
         particleCloud.clockM().stop("Flow");
         particleCloud.clockM().stop("Global");
     }
 
-    Info<< "End\n" << endl;
+    Info << "End\n"
+         << endl;
 
     return 0;
 }
 
-
 // ************************************************************************* //

applications/solvers/cfdemSolverPiso/cfdemSolverPiso.C.bk

@@ -0,0 +1,274 @@
+/*---------------------------------------------------------------------------*\
+    CFDEMcoupling - Open Source CFD-DEM coupling
+
+    CFDEMcoupling is part of the CFDEMproject
+    www.cfdem.com
+                                Christoph Goniva, christoph.goniva@cfdem.com
+                                Copyright (C) 1991-2009 OpenCFD Ltd.
+                                Copyright (C) 2009-2012 JKU, Linz
+                                Copyright (C) 2012-     DCS Computing GmbH,Linz
+-------------------------------------------------------------------------------
+License
+    This file is part of CFDEMcoupling.
+
+    CFDEMcoupling is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    CFDEMcoupling is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with CFDEMcoupling.  If not, see <http://www.gnu.org/licenses/>.
+
+Application
+    cfdemSolverPiso
+
+Description
+    Transient solver for incompressible flow.
+    Turbulence modelling is generic, i.e. laminar, RAS or LES may be selected.
+    The code is an evolution of the solver pisoFoam in OpenFOAM(R) 1.6,
+    where additional functionality for CFD-DEM coupling is added.
+\*---------------------------------------------------------------------------*/
+
+#include "fvCFD.H"
+#include "singlePhaseTransportModel.H"
+
+#include "OFversion.H"
+#if defined(version30)
+    #include "turbulentTransportModel.H"
+    #include "pisoControl.H"
+#else
+    #include "turbulenceModel.H"
+#endif
+#if defined(versionv1606plus) || defined(version40)
+    #include "fvOptions.H"
+#else
+    #include "fvIOoptionList.H"
+#endif
+#include "fixedFluxPressureFvPatchScalarField.H"
+#include "cfdemCloud.H"
+
+#if defined(anisotropicRotation)
+    #include "cfdemCloudRotation.H"
+#endif
+#if defined(superquadrics_flag)
+    #include "cfdemCloudRotationSuperquadric.H"
+#endif
+#include "implicitCouple.H"
+#include "clockModel.H"
+#include "smoothingModel.H"
+#include "forceModel.H"
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+int main(int argc, char *argv[])
+{
+    #include "setRootCase.H"
+    #include "createTime.H"
+    #include "createMesh.H"
+    #if defined(version30)
+        pisoControl piso(mesh);
+        #include "createTimeControls.H"
+    #endif
+    #include "createFields.H"
+    #include "createFvOptions.H"
+    #include "initContinuityErrs.H"
+
+    // create cfdemCloud
+    #include "readGravitationalAcceleration.H"
+    #include "checkImCoupleM.H"
+    #if defined(anisotropicRotation)
+        cfdemCloudRotation particleCloud(mesh);
+    #elif defined(superquadrics_flag)
+        cfdemCloudRotationSuperquadric particleCloud(mesh);
+    #else
+        cfdemCloud particleCloud(mesh);
+    #endif
+    #include "checkModelType.H"
+
+    // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+    Info<< "\nStarting time loop\n" << endl;
+    while (runTime.loop())
+    {
+        Info<< "Time = " << runTime.timeName() << nl << endl;
+
+        #if defined(version30)
+            #include "readTimeControls.H"
+            #include "CourantNo.H"
+            #include "setDeltaT.H"
+        #else
+            #include "readPISOControls.H"
+            #include "CourantNo.H"
+        #endif
+
+        // do particle stuff
+        particleCloud.clockM().start(1,"Global");
+        particleCloud.clockM().start(2,"Coupling");
+        bool hasEvolved = particleCloud.evolve(voidfraction,Us,U);
+
+        if(hasEvolved)
+        {
+            particleCloud.smoothingM().smoothenAbsolutField(particleCloud.forceM(0).impParticleForces());
+        }
+    
+        Ksl = particleCloud.momCoupleM(particleCloud.registryM().getProperty("implicitCouple_index")).impMomSource();
+        Ksl.correctBoundaryConditions();
+
+        surfaceScalarField voidfractionf = fvc::interpolate(voidfraction);
+        phi = voidfractionf*phiByVoidfraction;
+
+        //Force Checks
+        #include "forceCheckIm.H"
+
+        #include "solverDebugInfo.H"
+        particleCloud.clockM().stop("Coupling");
+
+        particleCloud.clockM().start(26,"Flow");
+
+        if(particleCloud.solveFlow())
+        {
+            // Pressure-velocity PISO corrector
+            {
+                // Momentum predictor
+                fvVectorMatrix UEqn
+                (
+                    fvm::ddt(voidfraction,U) - fvm::Sp(fvc::ddt(voidfraction),U)
+                  + fvm::div(phi,U) - fvm::Sp(fvc::div(phi),U)
+//                + turbulence->divDevReff(U)
+                  + particleCloud.divVoidfractionTau(U, voidfraction)
+                  ==
+                  - fvm::Sp(Ksl/rho,U)
+                  + fvOptions(U)
+                );
+
+                UEqn.relax();
+                fvOptions.constrain(UEqn);
+
+                #if defined(version30)
+                    if (piso.momentumPredictor())
+                #else
+                    if (momentumPredictor)
+                #endif
+                {
+                    if (modelType=="B" || modelType=="Bfull")
+                        solve(UEqn == - fvc::grad(p) + Ksl/rho*Us);
+                    else
+                        solve(UEqn == - voidfraction*fvc::grad(p) + Ksl/rho*Us);
+
+                    fvOptions.correct(U);
+                }
+
+                // --- PISO loop
+                #if defined(version30)
+                    while (piso.correct())
+                #else
+                    for (int corr=0; corr<nCorr; corr++)
+                #endif
+                {
+                    volScalarField rUA = 1.0/UEqn.A();
+
+                    surfaceScalarField rUAf("(1|A(U))", fvc::interpolate(rUA));
+                    volScalarField rUAvoidfraction("(voidfraction2|A(U))",rUA*voidfraction);
+                    surfaceScalarField rUAfvoidfraction("(voidfraction2|A(U)F)", fvc::interpolate(rUAvoidfraction));
+
+                    U = rUA*UEqn.H();
+
+                    #ifdef version23
+                        phi = ( fvc::interpolate(U) & mesh.Sf() )
+                            + rUAfvoidfraction*fvc::ddtCorr(U, phiByVoidfraction);
+                    #else
+                        phi = ( fvc::interpolate(U) & mesh.Sf() )
+                            + fvc::ddtPhiCorr(rUAvoidfraction, U, phiByVoidfraction);
+                    #endif
+                    surfaceScalarField phiS(fvc::interpolate(Us) & mesh.Sf());
+                    phi += rUAf*(fvc::interpolate(Ksl/rho) * phiS);
+
+                    if (modelType=="A")
+                        rUAvoidfraction = volScalarField("(voidfraction2|A(U))",rUA*voidfraction*voidfraction);
+
+                    // Update the fixedFluxPressure BCs to ensure flux consistency
+                    #include "fixedFluxPressureHandling.H"
+                    
+
+                    // Non-orthogonal pressure corrector loop
+                    #if defined(version30)
+                        while (piso.correctNonOrthogonal())
+                    #else
+                        for (int nonOrth=0; nonOrth<=nNonOrthCorr; nonOrth++)
+                    #endif
+                    {
+                        // Pressure corrector
+                        fvScalarMatrix pEqn
+                        (
+                            fvm::laplacian(rUAvoidfraction, p) == fvc::div(voidfractionf*phi) + particleCloud.ddtVoidfraction()
+                        );
+                        pEqn.setReference(pRefCell, pRefValue);
+
+                        #if defined(version30)
+                            pEqn.solve(mesh.solver(p.select(piso.finalInnerIter())));
+                            if (piso.finalNonOrthogonalIter())
+                            {
+                                phiByVoidfraction = phi - pEqn.flux()/voidfractionf;
+                            }
+                        #else
+                            if( corr == nCorr-1 && nonOrth == nNonOrthCorr )
+                                #if defined(versionExt32)
+                                    pEqn.solve(mesh.solutionDict().solver("pFinal"));
+                                #else
+                                    pEqn.solve(mesh.solver("pFinal"));
+                                #endif
+                            else
+                                pEqn.solve();
+
+                            if (nonOrth == nNonOrthCorr)
+                            {
+                                phiByVoidfraction = phi - pEqn.flux()/voidfractionf;
+                            }
+                        #endif
+
+                    } // end non-orthogonal corrector loop
+
+                    phi = voidfractionf*phiByVoidfraction;
+                    #include "continuityErrorPhiPU.H"
+
+                    if (modelType=="B" || modelType=="Bfull")
+                        U -= rUA*fvc::grad(p) - Ksl/rho*Us*rUA;
+                    else
+                        U -= voidfraction*rUA*fvc::grad(p) - Ksl/rho*Us*rUA;
+
+                    U.correctBoundaryConditions();
+                    fvOptions.correct(U);
+
+                } // end piso loop
+            }
+
+            laminarTransport.correct();
+            turbulence->correct();
+        }// end solveFlow
+        else
+        {
+            Info << "skipping flow solution." << endl;
+        }
+
+        runTime.write();
+
+        Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
+            << "  ClockTime = " << runTime.elapsedClockTime() << " s"
+            << nl << endl;
+
+        particleCloud.clockM().stop("Flow");
+        particleCloud.clockM().stop("Global");
+    }
+
+    Info<< "End\n" << endl;
+
+    return 0;
+}
+
+
+// ************************************************************************* //

applications/solvers/cfdemSolverPiso/createGravity.H

@@ -0,0 +1,12 @@
+    Info<< "\nReading g" << endl;
+    uniformDimensionedVectorField g
+    (
+        IOobject
+        (
+            "g",
+            runTime.constant(),
+            mesh,
+            IOobject::MUST_READ,
+            IOobject::NO_WRITE
+        )
+    );

doc/_build/html/CFDEMcoupling_Manual.html

@@ -193,7 +193,7 @@ A general formulation with a solid phase:</p>
 <li>Calculate the voidfraction field from particle positions, particle volumes and cell volumes. (voidFractionModel)</li>
 <li>Calculate the forces on the CFD field and the particle forces. These may be assembled by various forceModels, e.g. buyoancy and a drag-law. (forceModel and forceSubModel)</li>
 <li>Calculate the implicit and if necessary explicit source terms from the force fields for the Navier-Stokes-Equations, the solver needs to solve. (momCoupleModel)</li>
-<li>Calculate the &#8220;Us&#8221; field from particle data, map the particle data to the CFD-field. (averagingModel)</li>
+<li>Calculate the &#8220;Us&#8221; field from particle data, map the particle data to the CFD field. (averagingModel)</li>
 </ul>
 <p>In the standard setting the momentum coupling model takes the calculated forces of forceModels and produces the implicit momentum source term Ksl with</p>
 <img alt="_images/Ksl.png" class="align-center" src="_images/Ksl.png" />

doc/_build/html/_sources/CFDEMcoupling_Manual.txt

@@ -100,7 +100,7 @@ All cfdemSolvers use a "cfdemCloud" object, which takes care of many things. The
 * Calculate the voidfraction field from particle positions, particle volumes and cell volumes. (voidFractionModel)
 * Calculate the forces on the CFD field and the particle forces. These may be assembled by various forceModels, e.g. buyoancy and a drag-law. (forceModel and forceSubModel)
 * Calculate the implicit and if necessary explicit source terms from the force fields for the Navier-Stokes-Equations, the solver needs to solve. (momCoupleModel)
-* Calculate the "Us" field from particle data, map the particle data to the CFD-field. (averagingModel)
+* Calculate the "Us" field from particle data, map the particle data to the CFD field. (averagingModel)
 
 In the standard setting the momentum coupling model takes the calculated forces of forceModels and produces the implicit momentum source term Ksl with
 

doc/release-notes-PUBLIC-3.8.0.txt

@@ -0,0 +1,65 @@
+
+===============================================================
+Release notes for PUBLIC, version 3.8.0 
+
+src/lagrangian/cfdemParticle/subModels/IOModel/trackIO/trackIO.C |  Josef Kerbl (DCS) reworked the trackIO IOModel to transfer the LIGGGHTS® particle IDs to the CFDEMcoupling output 
+src/lagrangian/cfdemParticle/subModels/forceModel/checkCouplingInterval/checkCouplingInterval.C |  Josef Kerbl (DCS) changed the default maximum values of CFL and particle CFL numbers to 1 
+applications/solvers/cfdemSolverIB/cfdemSolverIB.C |  Josef Kerbl (DCS) added compatibility for non-Newtonian viscosity models 
+applications/solvers/cfdemSolverPiso/cfdemSolverPiso.C |  Josef Kerbl (DCS) added compatibility for non-Newtonian viscosity models 
+applications/solvers/cfdemSolverPisoSTM/cfdemSolverPisoSTM.C |  Josef Kerbl (DCS) added compatibility for non-Newtonian viscosity models 
+applications/solvers/cfdemSolverPisoScalar/cfdemSolverPisoScalar.C |  Josef Kerbl (DCS) added compatibility for non-Newtonian viscosity models 
+tutorials/.gitignore |  Alice Hager (DCS) moved all blockMeshDict files from constant/polyMesh to system 
+doc/_build/html/liggghtsCommandModel_execute.html |  Alice Hager (DCS) updated the documentation of new features and keywords (handling of numbers and labels, defined keywords, noBlanks, ...) 
+doc/_build/html/liggghtsCommandModel_setDEMGravity.html |  Alice Hager (DCS) added keyword unfix to documentation 
+doc/_build/html/CFDEMcoupling_Manual.html |  Alice Hager (DCS) added information about the solveFlow keyword. 
+applications/solvers/cfdemSolverIB/cfdemSolverIB.C |  Alice Hager (DCS) added the solveFlow keyword 
+cfdemCloudIB.C |  Alice Hager (DCS) removed the DEMTorques_ from the getDEMdata function (bugfix) 
+src/lagrangian/cfdemParticle/cfdemCloud/cfdemCloud.C |  Josef Kerbl (DCS) removed the unused coupleForce_ switch 
+src/lagrangian/cfdemParticle/subModels/forceModel/forceModel/forceModel.C |  Josef Kerbl (DCS) removed the unused coupleForce_ switch 
+src/lagrangian/cfdemParticle/subModels/forceModel/forceModel/forceModel.H |  Josef Kerbl (DCS) removed the unused coupleForce_ switch 
+src/lagrangian/cfdemParticle/subModels/forceModel/noDrag/noDrag.C |  Josef Kerbl (DCS) removed the unused coupleForce_ switch 
+src/lagrangian/cfdemParticle/subModels/dataExchangeModel/dataExchangeModel/dataExchangeModel.H |  Josef Kerbl (DCS) elevated the warning of non-matching DEM-timestep, couplingInterval and CFD-timestep to error 
+checkCouplingInterval.C |  Josef Kerbl (DCS) added checkCG to forceModels missing this 
+fieldStore.C |  Josef Kerbl (DCS) added checkCG to forceModels missing this 
+fieldTimeAverage.C |  Josef Kerbl (DCS) added checkCG to forceModels missing this 
+noDrag.C |  Josef Kerbl (DCS) added checkCG to forceModels missing this 
+volWeightedAverage.C |  Josef Kerbl (DCS) added checkCG to forceModels missing this 
+Archimedes.C |  Josef Kerbl (DCS) added support for typeSpecificCG to some force models 
+checkCouplingInterval.C |  Josef Kerbl (DCS) added support for typeSpecificCG to some force models 
+cfdemCloud.C |  Josef Kerbl (DCS) changed the default value of cgWarnOnly to false 
+oneWayVTK.C |  Christoph Goniva (DCS) bugfix to enable correct reading of v field (thx to Moritz Höfert) 
+applications/utilities/cfdemPostproc/cfdemPostproc.C |  Christoph Goniva (DCS) enable calculating an writing Us field (thx to Moritz Höfert) 
+meanSupVelocityForce.C |  Christoph Goniva (DCS) enable using cellZones 
+forceSubModel.C |  Josef Kerbl (DCS) added check, if scaleDrag > 0 
+bashrc |  Josef Kerbl (DCS) added support to automatically compile LIGGGHTS in debug mode for the Makefile.auto, if WM_COMPILE_OPTION is set to Debug 
+cshrc |  Josef Kerbl (DCS) added support to automatically compile LIGGGHTS in debug mode for the Makefile.auto, if WM_COMPILE_OPTION is set to Debug 
+functions.sh |  Josef Kerbl (DCS) added support to automatically compile LIGGGHTS in debug mode for the Makefile.auto, if WM_COMPILE_OPTION is set to Debug 
+cfdemSolverPisoSTM.C |  Christoph Goniva (DCS) added fvOptions and fixedFluxPressureBC for 5.x 
+cfdemSolverPiso.C |  Alice Hager (DCS) ported 5.x 
+cfdemSolverPisoSTM.C |  Christoph Goniva (DCS) - fvOptions added 
+applications/utilities/cfdemPostproc/cfdemPostproc.C |  Alice Hager (DCS) ported 5.x 
+meanSupVelocityForce.C |  Alice Hager (DCS) ported 5.x & major revision 
+eulerianScalarField.C |  Alice Hager (DCS) ported 5.x 
+cfdemCloud.C |  Alice Hager (DCS) ported 5.x 
+cfdemCloud.C |  Christoph Goniva (DCS) minor revision 
+cfdemCloudIB.C |  Alexander Podlozhnyuk (DCS) adaptions for periodic bc 
+IOModel.C |  Josef Kerbl (DCS) - revision ofor 5.x 
+dataExchangeModel.C |  Christoph Goniva (DCS) improved calculation of time step fraction 
+oneWayVTK.C |  Josef Kerbl (DCS) bugfix 
+LaEuScalarTemp.C |  Alexander Podlozhnyuk (DCS) extension for compressible flow 
+MeiLift.H |  Stefan Radl (TUG) reference added 
+checkCouplingInterval.C |  Christoph Goniva (DCS) consider scaleDia 
+checkCouplingInterval.C |  Christoph Goniva (DCS) default maxCFL=1 
+checkCouplingInterval.C |  Christoph Goniva (DCS) default maxPCFL=1 
+forceSubModel.C |  Christoph Goniva (DCS) new option scaleDH 
+gradPForce.C |  Christoph Goniva (DCS) verbose option 
+noDrag.C |  Christoph Goniva (DCS) bugfix 
+particleCellVolume.C |  Josef Kerbl (DCS) added header to output 
+scalarGeneralExchange.C |  Stefan Radl (TUG) major revisions 
+viscForce.C |  Christoph Goniva (DCS) verbose option 
+engineSearch.C |  Josef Kerbl (DCS) bugfix for MS use 
+engineSearchIB.C |  Alexander Podlozhnyuk (DCS) Stefan Radl (TUG) improved periodic 
+standardSearch.C |  Josef Kerbl (DCS) bugfix for MS use 
+IBVoidFraction.C |  Alexander Podlozhnyuk (DCS) Stefan Radl (TUG) improved periodic 
+===============================================================
+

src/eulerian/fvOptionsCFDEM/sources/derived/meanSupVelocityForce/meanSupVelocityForce.C

@@ -103,13 +103,13 @@ void Foam::fv::meanSupVelocityForce::correct(volVectorField& U)
 {
     if (twoPhase_)
     {
-        word alphaName = coeffs_.lookup("alphaField");
+        word alphaName = coeffs_.get<word>("alphaField");
         alpha_ = mesh_.lookupObject<volScalarField>(alphaName);
     }
 
     if (coupled_)
     {
-        word voidfractionName = coeffs_.lookup("voidfractionField");
+        word voidfractionName = coeffs_.get<word>("voidfractionField");
         voidfraction_ = mesh_.lookupObject<volScalarField>(voidfractionName);
     }
     

src/eulerian/scalarTransportModelsCFDEM/scalarTransportModel/newScalarTransportModel.C

@@ -64,7 +64,7 @@ autoPtr<scalarTransportModel> scalarTransportModel::New
     Info<< "Selecting scalarTransportModel "
 		 << scalarTransportModelType << endl;
 
-    dictionaryConstructorTable::iterator cstrIter =
+    auto cstrIter =
         dictionaryConstructorTablePtr_->find(scalarTransportModelType);
 
     if (cstrIter == dictionaryConstructorTablePtr_->end())

src/lagrangian/cfdemParticle/cfdTools/fixedFluxPressureHandling.H

@@ -1,41 +1,41 @@
-#ifdef version50
+// #ifdef version50
     if (modelType=="A")
     {
         volScalarField rUsed = rUA*voidfraction;
         constrainPressure(p, U, phi, rUsed);
     }
     else constrainPressure(p, U, phi, rUA);
-#endif
-#ifndef versionExt32
-#ifndef version40
-    if (modelType=="A")
-    {
-        setSnGrad<fixedFluxPressureFvPatchScalarField>
-        (
-            #ifdef versionv1612plus
-                p.boundaryFieldRef(),
-            #else
-                p.boundaryField(),
-            #endif
-            (
-                phi.boundaryField()
-              - (mesh.Sf().boundaryField() & U.boundaryField())
-            )/(mesh.magSf().boundaryField()*rUAf.boundaryField()*voidfractionf.boundaryField())
-        );
-    }else
-    {
-        setSnGrad<fixedFluxPressureFvPatchScalarField>
-        (
-            #ifdef versionv1612plus
-                p.boundaryFieldRef(),
-            #else
-                p.boundaryField(),
-            #endif
-            (
-                phi.boundaryField()
-              - (mesh.Sf().boundaryField() & U.boundaryField())
-            )/(mesh.magSf().boundaryField()*rUAf.boundaryField())
-        );
-    }
-#endif
-#endif
+// #endif
+// #ifndef versionExt32
+// #ifndef version40
+//     if (modelType=="A")
+//     {
+//         setSnGrad<fixedFluxPressureFvPatchScalarField>
+//         (
+//             #ifdef versionv1612plus
+//                 p.boundaryFieldRef(),
+//             #else
+//                 p.boundaryField(),
+//             #endif
+//             (
+//                 phi.boundaryField()
+//               - (mesh.Sf().boundaryField() & U.boundaryField())
+//             )/(mesh.magSf().boundaryField()*rUAf.boundaryField()*voidfractionf.boundaryField())
+//         );
+//     }else
+//     {
+//         setSnGrad<fixedFluxPressureFvPatchScalarField>
+//         (
+//             #ifdef versionv1612plus
+//                 p.boundaryFieldRef(),
+//             #else
+//                 p.boundaryField(),
+//             #endif
+//             (
+//                 phi.boundaryField()
+//               - (mesh.Sf().boundaryField() & U.boundaryField())
+//             )/(mesh.magSf().boundaryField()*rUAf.boundaryField())
+//         );
+//     }
+// #endif
+// #endif

src/lagrangian/cfdemParticle/cfdTools/forceCheckEx.H

@@ -9,11 +9,11 @@ reduce(fExpSource, sumOp<dimensionedVector>()); //must reduce, since sum in func
 //calc. total explicit force
 dimensionedVector fExpTotal = -gSum( 
                                     mesh.V()
-                                #if defined(version40)
+                                // #if defined(version40)
                                     * f.primitiveField()
-                                #else
-                                    * f.internalField()
-                                #endif
+                                // #else
+                                //     * f.internalField()
+                                // #endif
                                     //
                                   ); //negative because we want the force on the fluid
 

src/lagrangian/cfdemParticle/cfdTools/global.H

@@ -41,6 +41,9 @@ SourceFiles
 
 #include "fvCFD.H"
 #include "cfdemCloud.H"
+#include "foamVersion.H"
+// #include "runTimeSelectionTables.H"
+// #include "addToRunTimeSelectionTable.H"
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
 namespace Foam

src/lagrangian/cfdemParticle/cfdTools/newGlobal.C

@@ -56,7 +56,7 @@ autoPtr<global> global::New
          << globalType << endl;
 
 
-    dictionaryConstructorTable::iterator cstrIter =
+    auto cstrIter =
         dictionaryConstructorTablePtr_->find(globalType);
 
     if (cstrIter == dictionaryConstructorTablePtr_->end())

src/lagrangian/cfdemParticle/cfdTools/versionInfo.H

@@ -33,9 +33,10 @@ Description
 #ifndef versionInfo_H
 #define versionInfo_H
 
-word CFDEMversion="cfdem-3.8.0";
+word CFDEMversion="cfdem-3.8.1";
 word compatibleLIGGGHTSversion="3.8.0";
-word OFversion="5.x-commit-538044ac05c4672b37c7df607dca1116fa88df88";
+// word OFversion="5.x-commit-538044ac05c4672b37c7df607dca1116fa88df88";
+word OFversion="v2406";
 
 Info << "\nCFDEMcoupling version: " << CFDEMversion <<  endl;
 Info << ", compatible to LIGGGHTS version: " << compatibleLIGGGHTSversion <<  endl;

src/lagrangian/cfdemParticle/cfdemCloud/cfdemCloud.C

@@ -46,109 +46,65 @@ Description
 #include "smoothingModel.H"
 #include "liggghtsCommandModel.H"
 #include "cyclicPolyPatch.H"
+#include "processorPolyPatch.H"
 
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
-Foam::cfdemCloud::cfdemCloud
-(
-    const fvMesh& mesh
-)
-:
-    mesh_(mesh),
-    couplingProperties_
-    (
-        IOobject
-        (
-            "couplingProperties",
-            mesh_.time().constant(),
-            mesh_,
-            IOobject::MUST_READ,
-            IOobject::NO_WRITE
-        )
-    ),
-    liggghtsCommandDict_
-    (
-        IOobject
-        (
-            "liggghtsCommands",
-            mesh_.time().constant(),
-            mesh_,
-            IOobject::MUST_READ,
-            IOobject::NO_WRITE
-        )
-    ),
-    allowAdjustTimeStep_(couplingProperties_.lookupOrDefault<Switch>("allowAdjustTimeStep", false)),
-    solveFlow_(true),
-    solveScalarTransport_(true),
-    verbose_(couplingProperties_.lookupOrDefault<Switch>("verbose", false)),
-    expCorrDeltaUError_(false),
-    debug_(false),
-    allowCFDsubTimestep_(true),
-    ignore_(false),
-    modelType_(couplingProperties_.lookup("modelType")),
-    positions_(NULL),
-    velocities_(NULL),
-    fluidVel_(NULL),
-    fAcc_(NULL),
-    impForces_(NULL),
-    expForces_(NULL),
-    DEMForces_(NULL),
-    Cds_(NULL),
-    radii_(NULL),
-    voidfractions_(NULL),
-    cellIDs_(NULL),
-    particleWeights_(NULL),
-    particleVolumes_(NULL),
-    particleV_(NULL),
-    dragPrev_(NULL),
-    numberOfParticles_(0),
-    d32_(-1),
-    numberOfParticlesChanged_(false),
-    arraysReallocated_(false),
-    forceModels_(couplingProperties_.lookup("forceModels")),
-    momCoupleModels_(couplingProperties_.lookup("momCoupleModels")),
-    liggghtsCommandModelList_(liggghtsCommandDict_.lookup("liggghtsCommandModels")),
-    turbulenceModelType_(couplingProperties_.lookup("turbulenceModelType")),
-    isLES_(false),
-    cg_(1.),
-    cgOK_(true),
-    impDEMdrag_(false),
-    impDEMdragAcc_(false),
-    imExSplitFactor_(1.0),
-    treatVoidCellsAsExplicitForce_(false),
-    useDDTvoidfraction_("off"),
-    ddtVoidfraction_
-    (   
-        IOobject
-        (
-            "ddtVoidfraction",
-            mesh.time().timeName(),
-            mesh,
-            IOobject::NO_READ,
-            IOobject::AUTO_WRITE
-        ),
-        mesh,
-        dimensionedScalar("zero", dimensionSet(0,0,-1,0,0), 0)  // 1/s
-    ),
-    checkPeriodicCells_(false),
-    wall_periodicityCheckRange_(vector(1,1,1)),
-    wall_periodicityCheckTolerance_(1e-07),
-    meshHasUpdated_(false),
-    turbulence_
-    (
-        #if defined(version24Dev)
-            mesh.lookupObject<turbulenceModel>
-        #elif defined(version21) || defined(version16ext)
-            #ifdef compre
-                mesh.lookupObject<compressible::turbulenceModel>
-            #else
-                mesh.lookupObject<incompressible::turbulenceModel>
-            #endif
-        #elif defined(version15)
-            mesh.lookupObject<incompressible::RASModel>
-        #endif
-        (
-            turbulenceModelType_
-        )
+Foam::cfdemCloud::cfdemCloud(const fvMesh &mesh)
+    : mesh_(mesh),
+      couplingProperties_(IOobject("couplingProperties",
+                                   mesh_.time().constant(), mesh_,
+                                   IOobject::MUST_READ, IOobject::NO_WRITE)),
+      liggghtsCommandDict_(IOobject("liggghtsCommands", mesh_.time().constant(),
+                                    mesh_, IOobject::MUST_READ,
+                                    IOobject::NO_WRITE)),
+      allowAdjustTimeStep_(couplingProperties_.lookupOrDefault<Switch>(
+          "allowAdjustTimeStep", false)),
+      solveFlow_(true), solveScalarTransport_(true),
+      verbose_(couplingProperties_.lookupOrDefault<Switch>("verbose", false)),
+      expCorrDeltaUError_(false), debug_(false), allowCFDsubTimestep_(true),
+      ignore_(false), writeTimePassed_(false),
+      modelType_(couplingProperties_.lookup("modelType")), positions_(NULL),
+      velocities_(NULL), fluidVel_(NULL), fAcc_(NULL), impForces_(NULL),
+      expForces_(NULL), DEMForces_(NULL), Cds_(NULL), radii_(NULL),
+      voidfractions_(NULL), cellIDs_(NULL), particleWeights_(NULL),
+      particleVolumes_(NULL), particleV_(NULL), dragPrev_(NULL),
+      numberOfParticles_(0), d32_(-1), numberOfParticlesChanged_(false),
+      arraysReallocated_(false),
+      forceModels_(couplingProperties_.lookup("forceModels")),
+      momCoupleModels_(couplingProperties_.lookup("momCoupleModels")),
+      liggghtsCommandModelList_(
+          liggghtsCommandDict_.lookup("liggghtsCommandModels")),
+      turbulenceModelType_(couplingProperties_.lookup("turbulenceModelType")),
+      isLES_(false), cg_(1.), cgOK_(true), impDEMdrag_(false),
+      impDEMdragAcc_(false), imExSplitFactor_(1.0),
+      treatVoidCellsAsExplicitForce_(false), useDDTvoidfraction_("off"),
+      ddtVoidfraction_(
+          IOobject("ddtVoidfraction", mesh.time().timeName(), mesh,
+                   IOobject::NO_READ, IOobject::AUTO_WRITE),
+          mesh,
+          dimensionedScalar("zero", dimensionSet(0, 0, -1, 0, 0), 0) // 1/s
+          ),
+      checkPeriodicCells_(false), wall_periodicityCheckRange_(vector(1, 1, 1)),
+      wall_periodicityCheckTolerance_(1e-07), meshHasUpdated_(false),
+      turbulence_(
+          // #if defined(version24Dev)
+          //     mesh.lookupObject<turbulenceModel>
+          // #elif defined(version21) || defined(version16ext)
+          //     #ifdef compre
+          //         mesh.lookupObject<compressible::turbulenceModel>
+          //     #else
+          //         mesh.lookupObject<incompressible::turbulenceModel>
+          //     #endif
+          // #elif defined(version15)
+          //     mesh.lookupObject<incompressible::RASModel>
+          // #endif
+          // (
+          //     turbulenceModelType_
+          // )
+          mesh.lookupObject<turbulenceModel>
+          (
+               turbulenceModelType_
+          )
     ),
     turbulenceMultiphase_
     (   
@@ -300,7 +256,7 @@ Foam::cfdemCloud::cfdemCloud
 
     if (couplingProperties_.found("treatVoidCellsAsExplicitForce"))
         treatVoidCellsAsExplicitForce_ = readBool(couplingProperties_.lookup("treatVoidCellsAsExplicitForce"));
-    if (couplingProperties_.found("ignore")) ignore_=true;
+    if (couplingProperties_.found("ignore")) ignore_=true; // 指定ignore就会被忽略
     if (turbulenceModelType_=="LESProperties")
     {
         isLES_ = true;
@@ -383,17 +339,17 @@ Foam::cfdemCloud::cfdemCloud
     forAll(patches, patchI)
     {
         const polyPatch& pp = patches[patchI];
-        #if defined(versionExt32)
-        if (isA<cyclicPolyPatch>(pp))
-            nPatchesCyclic++;
-        else if (!isA<processorPolyPatch>(pp))
-            nPatchesNonCyclic++;
-        #else
+        // #if defined(versionExt32)
+        // if (isA<cyclicPolyPatch>(pp))
+        //     nPatchesCyclic++;
+        // else if (!isA<processorPolyPatch>(pp))
+        //     nPatchesNonCyclic++;
+        // #else
         if (isA<cyclicPolyPatch>(pp) || isA<cyclicAMIPolyPatch>(pp))
             nPatchesCyclic++;
         else if (!isA<processorPolyPatch>(pp))
             nPatchesNonCyclic++;
-        #endif
+        // #endif
     }
     if(nPatchesNonCyclic==0)
         checkPeriodicCells_=true;
@@ -411,6 +367,7 @@ Foam::cfdemCloud::cfdemCloud
     //Check if user attempts to change fluid time step
     if( mesh_.time().controlDict().lookupOrDefault<Switch>("adjustTimeStep", false) && !allowAdjustTimeStep_ )
     {
+        // NOTE(zhyang): 不支持自适应时间步长
         FatalError << "cfdemCloud:: you want to adjustTimeStep in controlDict. This is not allowed in this version of CFDEM."
                    << abort(FatalError);
     }
@@ -705,10 +662,11 @@ bool Foam::cfdemCloud::evolve
 
     if(!ignore())
     {
+        if(!writeTimePassed_ && mesh_.time().outputTime()) writeTimePassed_=true;
         if (dataExchangeM().doCoupleNow())
         {
             Info << "\n Coupling..." << endl;
-            dataExchangeM().couple(0);
+            dataExchangeM().couple(0); // NOTE(zhyang)：耦合
             doCouple=true;
 
             // reset vol Fields

src/lagrangian/cfdemParticle/cfdemCloud/cfdemCloud.H

@@ -49,11 +49,11 @@ SourceFiles
 #include "fvCFD.H"
 #include "IFstream.H"
 
-#if defined(version21) || defined(version16ext)
-    #include "turbulenceModel.H"
-#elif defined(version15)
-    #include "RASModel.H"
-#endif
+// #if defined(version21) || defined(version16ext)
+#include "turbulenceModel.H"
+// #elif defined(version15)
+//     #include "RASModel.H"
+// #endif
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
@@ -106,6 +106,8 @@ protected:
 
     bool ignore_;
 
+    bool writeTimePassed_;
+
     const word modelType_;
 
     mutable double **positions_;
@@ -187,17 +189,17 @@ protected:
 
     bool meshHasUpdated_;
 
-    #if defined(version24Dev)
-        const turbulenceModel& turbulence_;
-    #elif defined(version21) || defined(version16ext)
-        #ifdef compre
-            const compressible::turbulenceModel& turbulence_;
-        #else
-            const incompressible::turbulenceModel& turbulence_;
-        #endif
-    #elif defined(version15)
-        const incompressible::RASModel& turbulence_;
-    #endif
+    // #if defined(version24Dev)
+    const turbulenceModel& turbulence_;
+    // #elif defined(version21) || defined(version16ext)
+    //     #ifdef compre
+    //         const compressible::turbulenceModel& turbulence_;
+    //     #else
+    //         const incompressible::turbulenceModel& turbulence_;
+    //     #endif
+    // #elif defined(version15)
+    //     const incompressible::RASModel& turbulence_;
+    // #endif
   
     //Multiphase Turbulence (e.g., slip-induced turbulence)
     mutable volScalarField turbulenceMultiphase_;
@@ -332,6 +334,8 @@ public:
 
        inline const bool& ignore() const;
 
+       inline const bool& writeTimePassed() const;
+
        inline const fvMesh& mesh() const;
 
        inline bool allowAdjustTimeStep() const;
@@ -513,17 +517,17 @@ public:
 
        inline autoPtr<liggghtsCommandModel>* liggghtsCommand() const;
 
-       #if defined(version24Dev)
+      //  #if defined(version24Dev)
             inline const turbulenceModel& turbulence() const;
-       #elif defined(version21) || defined(version16ext)
-            #ifdef compre
-                inline const compressible::turbulenceModel& turbulence() const;
-            #else
-                inline const incompressible::turbulenceModel& turbulence() const;
-            #endif
-       #elif defined(version15)
-           inline const incompressible::RASModel& turbulence() const;
-       #endif
+      //  #elif defined(version21) || defined(version16ext)
+      //       #ifdef compre
+      //           inline const compressible::turbulenceModel& turbulence() const;
+      //       #else
+      //           inline const incompressible::turbulenceModel& turbulence() const;
+      //       #endif
+      //  #elif defined(version15)
+      //      inline const incompressible::RASModel& turbulence() const;
+      //  #endif
 
         inline void makeSpecific(volScalarField&);
 

src/lagrangian/cfdemParticle/cfdemCloud/cfdemCloudI.H

@@ -101,6 +101,11 @@ inline const bool& cfdemCloud::ignore() const
     return ignore_;
 }
 
+inline const bool& cfdemCloud::writeTimePassed() const
+{
+    return writeTimePassed_;
+}
+
 inline const fvMesh& cfdemCloud::mesh() const
 {
     return mesh_;
@@ -345,17 +350,17 @@ inline autoPtr<liggghtsCommandModel>* cfdemCloud::liggghtsCommand() const
     return liggghtsCommand_;
 }
 
-#if defined(version24Dev)
-    inline const turbulenceModel& cfdemCloud::turbulence() const
-#elif defined(version21) || defined(version16ext)
-    #ifdef compre
-        inline const compressible::turbulenceModel& cfdemCloud::turbulence() const
-    #else
-        inline const incompressible::turbulenceModel& cfdemCloud::turbulence() const
-    #endif
-#elif defined(version15)
-    inline const incompressible::RASModel& cfdemCloud::turbulence() const
-#endif
+// #if defined(version24Dev)
+inline const turbulenceModel& cfdemCloud::turbulence() const
+// #elif defined(version21) || defined(version16ext)
+//     #ifdef compre
+//         inline const compressible::turbulenceModel& cfdemCloud::turbulence() const
+//     #else
+//         inline const incompressible::turbulenceModel& cfdemCloud::turbulence() const
+//     #endif
+// #elif defined(version15)
+//     inline const incompressible::RASModel& cfdemCloud::turbulence() const
+// #endif
 {
     return turbulence_;
 }

src/lagrangian/cfdemParticle/derived/cfdemCloudIB/cfdemCloudIB.C

@@ -84,7 +84,7 @@ cfdemCloudIB::cfdemCloudIB
     if(this->couplingProperties().found("tortuosity"))
     {
         calculateTortuosity_ = true;
-        flowDir_ = this->couplingProperties().subDict("tortuosity").lookup("flowDirection");
+        flowDir_ = this->couplingProperties().subDict("tortuosity").get<vector>("flowDirection");
         flowDir_ = flowDir_ / mag(flowDir_);
         Info << "Will calculate tortuosity in the mean flow direction ("<<flowDir_[0]<<" "<<flowDir_[1]<<" "<<flowDir_[2]<<")"<< endl;
     }
@@ -165,6 +165,7 @@ bool Foam::cfdemCloudIB::evolve
         skipLagrangeToEulerMapping_=true;
     }
 
+    if(!writeTimePassed_ && mesh_.time().outputTime()) writeTimePassed_=true;
     if (dataExchangeM().doCoupleNow())
     {
         Info << "\n timeStepFraction() = " << dataExchangeM().timeStepFraction() << endl;

src/lagrangian/cfdemParticle/etc/OFversion/OFversion.H

@@ -14,6 +14,8 @@
     #define versionv1612plus
 #elif(CFDEMWMPROJECTVERSION == 1706)
     #define versionv1706
+#elif(CFDEMWMPROJECTVERSION == 2406)
+    #define versionv2406
 #endif
 
 //define anisotropicRotation cloud models
@@ -48,6 +50,10 @@
     #define versionv1606plus
 #endif
 
+#if defined(versionv2406plus)
+    #define versionv2406plus
+#endif
+
 // features of 3.0 work also in v1606+
 #if defined(versionv1606plus)
     #define version30

src/lagrangian/cfdemParticle/etc/addLibs_universal/additionalLibs_v2406

@@ -0,0 +1,56 @@
+# Specify additional include and library paths, as well as libraries for the compilation
+#
+#CFDEM_ADD_INC = \
+#CFDEM_ADD_LIB_PATHS = \
+#CFDEM_ADD_LIBS = \
+
+# additional static libraries to be linked to lagrangian library
+CFDEM_ADD_STATICLIBS = \
+
+# include flags for compiling with SQ
+include $(CFDEM_ADD_LIBS_DIR)/additionalLibs_superquadric
+
+#################################################################
+## SETTINGS FOR 3.0.x                                          ##
+#################################################################
+#----------------------------------------------------------------
+# incompressible turbulence model settings
+#----------------------------------------------------------------
+# paths for incompressible turbulence models to use
+CFDEM_ADD_INCOMPTURBMOD_PATHS = \
+-I$(LIB_SRC)/TurbulenceModels/turbulenceModels/lnInclude \
+-I$(LIB_SRC)/TurbulenceModels/incompressible/lnInclude \
+-I$(LIB_SRC)/fvOptions/lnInclude \
+
+# libs for turbulence models to use
+CFDEM_ADD_INCOMPTURBMOD_LIBS = \
+-lturbulenceModels \
+-lincompressibleTurbulenceModels \
+-lfvOptions \
+
+CFDEM_TRI_SURF = \
+-lsurfMesh
+
+# 对应的是libthermophysicalProperties，当前没有被使用
+CFDEM_SPRAY_LIBS = \
+    -lliquidProperties \
+    -lliquidMixtureProperties \
+    -lsolidProperties \
+    -lsolidMixtureProperties \
+    -lthermophysicalFunctions 
+
+#----------------------------------------------------------------
+# compressible turbulence model settings
+#----------------------------------------------------------------
+# paths for compressible turbulence models to use
+CFDEM_ADD_COMPTURBMOD_PATHS = \
+-I$(LIB_SRC)/TurbulenceModels/turbulenceModels/lnInclude \
+-I$(LIB_SRC)/TurbulenceModels/compressible/lnInclude \
+-I$(LIB_SRC)/transportModels/compressible/lnInclude \
+-I$(LIB_SRC)/thermophysicalModels/radiation/lnInclude \
+
+# libs for turbulence models to use
+CFDEM_ADD_COMPTURBMOD_LIBS = \
+-lturbulenceModels \
+-lcompressibleTurbulenceModels \
+#################################################################

src/lagrangian/cfdemParticle/etc/bashrc

@@ -17,13 +17,13 @@
 ## STANDARD BLOCK FOR REGULAR INSTALLATIONS
 ##================================================#
 ##- source cfdem env vars
-#export CFDEM_VERSION=PUBLIC
-#export CFDEM_PROJECT_DIR=$HOME/CFDEM/CFDEMcoupling-$CFDEM_VERSION-$WM_PROJECT_VERSION
-#export CFDEM_PROJECT_USER_DIR=$HOME/CFDEM/$LOGNAME-$CFDEM_VERSION-$WM_PROJECT_VERSION
-#export CFDEM_bashrc=$CFDEM_PROJECT_DIR/src/lagrangian/cfdemParticle/etc/bashrc
-#export CFDEM_LIGGGHTS_SRC_DIR=$HOME/LIGGGHTS/LIGGGHTS-PUBLIC/src
-#export CFDEM_LIGGGHTS_MAKEFILE_NAME=auto
-#export CFDEM_LPP_DIR=$HOME/LIGGGHTS/mylpp/src
+export CFDEM_VERSION=PUBLIC
+export CFDEM_PROJECT_DIR=$HOME/CFDEM/CFDEMcoupling-$CFDEM_VERSION-$WM_PROJECT_VERSION
+export CFDEM_PROJECT_USER_DIR=$HOME/CFDEM/$LOGNAME-$CFDEM_VERSION-$WM_PROJECT_VERSION
+export CFDEM_bashrc=$CFDEM_PROJECT_DIR/src/lagrangian/cfdemParticle/etc/bashrc
+export CFDEM_LIGGGHTS_SRC_DIR=$HOME/LIGGGHTS/LIGGGHTS-PUBLIC/src
+export CFDEM_LIGGGHTS_MAKEFILE_NAME=auto
+export CFDEM_LPP_DIR=$HOME/LIGGGHTS/LPP/src
 #. $CFDEM_bashrc
 #================================================#
 #------------------------------------------------------------------------------
@@ -356,6 +356,9 @@ elif [[ $WM_PROJECT_VERSION == v1706 ]]; then
     cfdemecho "In doubt switch to OpenFOAM(R) 5.x."
     cfdemecho "************************************"
     sleep 1.5
+elif [[ $WM_PROJECT_VERSION == v2406 ]]; then 
+    export CFDEM_WM_PROJECT_VERSION=2406
+    sleep 1.5
 elif [[ $WM_PROJECT_VERSION == 2.4.* ]]; then
     export CFDEM_WM_PROJECT_VERSION=24
     cfdemecho "************************************"

src/lagrangian/cfdemParticle/etc/functions.sh

@@ -80,7 +80,7 @@ compileLib()
                 echo "Please make sure to have the incompressible libraries first in the library-list.txt!"
                 cd $CFDEM_SRC_DIR/lagrangian/cfdemParticle
                 echo "changing to $PWD"
-                if [[ $WM_PROJECT_VERSION == dev || $WM_PROJECT_VERSION == 3.0.* || $WM_PROJECT_VERSION == 4.*  || $WM_PROJECT_VERSION == 5.* ]]; then
+                if [[ $WM_PROJECT_VERSION == v2406 || $WM_PROJECT_VERSION == dev || $WM_PROJECT_VERSION == 3.0.* || $WM_PROJECT_VERSION == 4.*  || $WM_PROJECT_VERSION == 5.* ]]; then
                     wrmdep 2>&1 | tee -a $logpath/$logfileName
                 else
                     rmdepall 2>&1 | tee -a $logpath/$logfileName
@@ -90,7 +90,7 @@ compileLib()
             else
                 echo "Compiling a incompressible library."
         fi
-        if [[ $WM_PROJECT_VERSION == dev || $WM_PROJECT_VERSION == 3.0.* || $WM_PROJECT_VERSION == 4.*  || $WM_PROJECT_VERSION == 5.* ]]; then
+        if [[ $WM_PROJECT_VERSION == v2406 || $WM_PROJECT_VERSION == dev || $WM_PROJECT_VERSION == 3.0.* || $WM_PROJECT_VERSION == 4.*  || $WM_PROJECT_VERSION == 5.* ]]; then
             wrmdep 2>&1 | tee -a $logpath/$logfileName
         else
             rmdepall 2>&1 | tee -a $logpath/$logfileName

src/lagrangian/cfdemParticle/subModels/IOModel/IOModel/newIOModel.C

@@ -54,7 +54,7 @@ autoPtr<IOModel> IOModel::New
     Info<< "Selecting IOModel "
          << IOModelType << endl;
 
-    dictionaryConstructorTable::iterator cstrIter =
+    auto cstrIter =
         dictionaryConstructorTablePtr_->find(IOModelType);
 
     if (cstrIter == dictionaryConstructorTablePtr_->end())

src/lagrangian/cfdemParticle/subModels/averagingModel/averagingModel/newAveragingModel.C

@@ -56,7 +56,7 @@ autoPtr<averagingModel> averagingModel::New
          << averagingModelType << endl;
 
 
-    dictionaryConstructorTable::iterator cstrIter =
+    auto cstrIter =
         dictionaryConstructorTablePtr_->find(averagingModelType);
 
     if (cstrIter == dictionaryConstructorTablePtr_->end())

src/lagrangian/cfdemParticle/subModels/clockModel/clockModel/newClockModel.C

@@ -56,7 +56,7 @@ autoPtr<clockModel> clockModel::New
          << clockModelType << endl;
 
 
-    dictionaryConstructorTable::iterator cstrIter =
+    auto cstrIter =
         dictionaryConstructorTablePtr_->find(clockModelType);
 
     if (cstrIter == dictionaryConstructorTablePtr_->end())

src/lagrangian/cfdemParticle/subModels/dataExchangeModel/dataExchangeModel/dataExchangeModel.H

@@ -246,7 +246,7 @@ public:
             checkTSsize();
         }
 
-        inline bool doCoupleNow() const
+        inline bool doCoupleNow() const // NOTE(zhyang): 耦合步判断
         {
             // couple at first sub TS
             if ((particleCloud_.mesh().time().timeIndex()-timeIndexOffset_)*particleCloud_.mesh().time().deltaT().value()-SMALL

src/lagrangian/cfdemParticle/subModels/dataExchangeModel/dataExchangeModel/newDataExchangeModel.C

@@ -56,7 +56,7 @@ autoPtr<dataExchangeModel> dataExchangeModel::New
          << dataExchangeModelType << endl;
 
 
-    dictionaryConstructorTable::iterator cstrIter =
+    auto cstrIter =
         dictionaryConstructorTablePtr_->find(dataExchangeModelType);
 
     if (cstrIter == dictionaryConstructorTablePtr_->end())

src/lagrangian/cfdemParticle/subModels/dataExchangeModel/oneWayVTK/oneWayVTK.C

@@ -111,8 +111,8 @@ void oneWayVTK::getData
         // set file pointer
         string HH=string(H);
         const char * paricleFilePath=HH.c_str();
-        ifstream* inputPtr;
-        inputPtr = new ifstream(paricleFilePath);
+        std::ifstream* inputPtr;
+        inputPtr = new std::ifstream(paricleFilePath);
         if(!*inputPtr) FatalError << "File not found!, " << H << "\n" << abort(FatalError);
 
         if (name == "radius")
@@ -156,8 +156,8 @@ void oneWayVTK::getData
         // set file pointer
         string HH=string(H);
         const char * paricleFilePath=HH.c_str();
-        ifstream* inputPtr;
-        inputPtr = new ifstream(paricleFilePath);
+        std::ifstream* inputPtr;
+        inputPtr = new std::ifstream(paricleFilePath);
         if(!*inputPtr) FatalError << "File not found!, " << H << "\n" << abort(FatalError);
 
         // read position data from VTK file

src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMPI/twoWayMPI.C

@@ -228,14 +228,14 @@ bool Foam::twoWayMPI::couple(int i) const
 
             // check if liggghtsCommandModels with exaxt timing are being run
             bool exactTiming(false);
-            int runComNr = -10;
-            DynamicList<scalar> interruptTimes(0);
-            DynamicList<int> DEMstepsToInterrupt(0);
-            DynamicList<int> lcModel(0);
+            int runComNr = -10; // runCommand index
+            DynamicList<scalar> interruptTimes(0);   // 时间列表
+            DynamicList<int> DEMstepsToInterrupt(0); // 时间步数列表
+            DynamicList<int> lcModel(0); // NOTE(zhyang): liggghtsCommandModel
 
             forAll(particleCloud_.liggghtsCommandModelList(),i)
             {
-                // Check if exact timing is needed
+                // Check if exact timing is needed NOTE(zhyang): 确定执行的时间列表
                 // get time for execution
                 // store time for execution in list
                 if(particleCloud_.liggghtsCommand()[i]().exactTiming())
@@ -290,7 +290,7 @@ bool Foam::twoWayMPI::couple(int i) const
                     particleCloud_.liggghtsCommand()[runComNr]().set(DEMstepsToInterrupt[j]);
                     const char* command = particleCloud_.liggghtsCommand()[runComNr]().command(0);
                     Info << "Executing run command: '"<< command <<"'"<< endl;
-                    lmp->input->one(command);
+                    lmp->input->one(command); // NOTE(zhyang): 调用liggghts
     
                     // run liggghts command with exact timing
                     command = particleCloud_.liggghtsCommand()[lcModel[j]]().command(0);
@@ -338,7 +338,7 @@ bool Foam::twoWayMPI::couple(int i) const
                         {
                             const char* command = particleCloud_.liggghtsCommand()[i]().command(j);
                             Info << "Executing command: '"<< command <<"'"<< endl;
-                            lmp->input->one(command);
+                            lmp->input->one(command); // NOTE(zhyang): run ... 例如run 100, 根据couplingInterval
                         }
                     }
                 }
@@ -355,7 +355,7 @@ bool Foam::twoWayMPI::couple(int i) const
 
         // re-allocate arrays of cloud
         particleCloud_.clockM().start(4,"LIGGGHTS_reallocArrays");
-        particleCloud_.reAllocArrays();
+        particleCloud_.reAllocArrays(); // NOTE(zhyang): 根据颗粒数重新分配粒子相关容器（数组）的大小
         particleCloud_.clockM().stop("LIGGGHTS_reallocArrays");
     }
 

src/lagrangian/cfdemParticle/subModels/forceModel/Archimedes/Archimedes.C

@@ -64,11 +64,7 @@ Archimedes::Archimedes
     propsDict_(dict.subDict(typeName + "Props")),
     twoDimensional_(false),
     gravityFieldName_(propsDict_.lookup("gravityFieldName")),
-    #if defined(version21) || defined(version16ext)
-        g_(sm.mesh().lookupObject<uniformDimensionedVectorField> (gravityFieldName_))
-    #elif defined(version15)
-        g_(dimensionedVector(sm.mesh().lookupObject<IOdictionary>("environmentalProperties").lookup(gravityFieldName_)).value())
-    #endif
+    g_(sm.mesh().lookupObject<uniformDimensionedVectorField> (gravityFieldName_))
 {
 
     //Append the field names to be probed

src/lagrangian/cfdemParticle/subModels/forceModel/Archimedes/Archimedes.H

@@ -65,11 +65,7 @@ private:
 
     word gravityFieldName_;
 
-    #ifdef version21
-        const uniformDimensionedVectorField& g_;    // ref to gravity
-    #elif defined(version16ext) || defined(version15)
-        const dimensionedVector& g_;    // ref to gravity
-    #endif
+    const uniformDimensionedVectorField& g_;    // ref to gravity
 
 
 public:

src/lagrangian/cfdemParticle/subModels/forceModel/ArchimedesIB/ArchimedesIB.C

@@ -66,11 +66,7 @@ ArchimedesIB::ArchimedesIB
     voidfractionFieldName_(propsDict_.lookup("voidfractionFieldName")), //mod by alice
     voidfractions_(sm.mesh().lookupObject<volScalarField> (voidfractionFieldName_)),//mod by alice
     gravityFieldName_(propsDict_.lookup("gravityFieldName")),
-    #if defined(version21) || defined(version16ext)
-        g_(sm.mesh().lookupObject<uniformDimensionedVectorField> (gravityFieldName_))
-    #elif  defined(version15)
-        g_(dimensionedVector(sm.mesh().lookupObject<IOdictionary>("environmentalProperties").lookup(gravityFieldName_)).value())
-    #endif
+    g_(sm.mesh().lookupObject<uniformDimensionedVectorField> (gravityFieldName_))
 {
     //Append the field names to be probed
     particleCloud_.probeM().initialize(typeName, typeName+".logDat");

src/lagrangian/cfdemParticle/subModels/forceModel/ArchimedesIB/ArchimedesIB.H

@@ -71,11 +71,7 @@ private:
 
     word gravityFieldName_;
 
-    #ifdef version21
-        const uniformDimensionedVectorField& g_;    // ref to gravity
-    #elif defined(version16ext) || defined(version15)
-        const dimensionedVector& g_;    // ref to gravity
-    #endif
+    const uniformDimensionedVectorField& g_;    // ref to gravity
 
 
 public:

src/lagrangian/cfdemParticle/subModels/forceModel/GidaspowDrag/GidaspowDrag.C

@@ -119,13 +119,14 @@ void GidaspowDrag::setForce() const
 {
 
     const volScalarField& rhoField = forceSubM(0).rhoField();
-    #if defined(version24Dev)
-       // there seems to have been a change in the return value of 
-       // particleCloud_.turbulence().nu() used by forceSubM(0).nuField();
-       const volScalarField& nufField = particleCloud_.turbulence().nu();
-    #else
+    // NOTE(zhyang)：需要check
+    // #if defined(version24Dev)
+    //    // there seems to have been a change in the return value of 
+    //    // particleCloud_.turbulence().nu() used by forceSubM(0).nuField();
+    //    const volScalarField& nufField = particleCloud_.turbulence().nu();
+    // #else
         const volScalarField& nufField = forceSubM(0).nuField();
-    #endif
+    // #endif
 
 
     //update force submodels to prepare for loop

src/lagrangian/cfdemParticle/subModels/forceModel/forceModel/forceModel.C

@@ -95,26 +95,26 @@ forceModel::forceModel
     implicitRotation_(false),
     forceSubModels_(0),
     forceSubModel_(new autoPtr<forceSubModel>[nrForceSubModels()]),
-    voidfractionInterpolator_(NULL),
-    UInterpolator_(NULL),
-    vorticityInterpolator_(NULL),
-    gradPInterpolator_(NULL),
-    gradUInterpolator_(NULL),
-    gradVoidfractionInterpolator_(NULL),
-    Up1Interpolator_(NULL),
-    Up2Interpolator_(NULL),
-    dSauterInterpolator_(NULL),
-    phiP1Interpolator_(NULL),
-    phiP2Interpolator_(NULL),
-    alphaInterpolator_(NULL),
-    gradAlphaInterpolator_(NULL),
-    TInterpolator_(NULL),
-    UsInterpolator_(NULL),
-    fluidScalarFieldInterpolator_(NULL),
-    gradPsolidInterpolator_(NULL),
-    shearRateInterpolator_(NULL),
-    DDtUInterpolator_(NULL),
-    divTauInterpolator_(NULL)
+    voidfractionInterpolator_(nullptr),
+    UInterpolator_(nullptr),
+    vorticityInterpolator_(nullptr),
+    gradPInterpolator_(nullptr),
+    gradUInterpolator_(nullptr),
+    gradVoidfractionInterpolator_(nullptr),
+    Up1Interpolator_(nullptr),
+    Up2Interpolator_(nullptr),
+    dSauterInterpolator_(nullptr),
+    phiP1Interpolator_(nullptr),
+    phiP2Interpolator_(nullptr),
+    alphaInterpolator_(nullptr),
+    gradAlphaInterpolator_(nullptr),
+    TInterpolator_(nullptr),
+    UsInterpolator_(nullptr),
+    fluidScalarFieldInterpolator_(nullptr),
+    gradPsolidInterpolator_(nullptr),
+    shearRateInterpolator_(nullptr),
+    DDtUInterpolator_(nullptr),
+    divTauInterpolator_(nullptr)
 {}
 
 

src/lagrangian/cfdemParticle/subModels/forceModel/forceModel/newForceModel.C

@@ -51,7 +51,7 @@ autoPtr<forceModel> forceModel::New
     Info<< "Selecting forceModel "
          << forceType << endl;
 
-    dictionaryConstructorTable::iterator cstrIter =
+    auto cstrIter =
         dictionaryConstructorTablePtr_->find(forceType);
 
     if (cstrIter == dictionaryConstructorTablePtr_->end())

src/lagrangian/cfdemParticle/subModels/forceModel/forceSubModels/forceSubModel/newForceSubModel.C

@@ -51,7 +51,7 @@ autoPtr<forceSubModel> forceSubModel::New
     Info<< " Selecting forceSubModel "
          << forceType << endl;
 
-    dictionaryConstructorTable::iterator cstrIter =
+    auto cstrIter =
         dictionaryConstructorTablePtr_->find(forceType);
 
     if (cstrIter == dictionaryConstructorTablePtr_->end())

src/lagrangian/cfdemParticle/subModels/forceModel/scalarGeneralExchange/scalarGeneralExchange.C

@@ -615,7 +615,7 @@ void scalarGeneralExchange::setupModel() const
     if(partHeatFluxName_!="none")
     {
         validPartFlux_=true;
-        particleCloud_.registerNamesFieldsUserCFDEMToExt(propsDict_.lookup("partHeatFluxName"), partHeatFluxPositionInRegister_);
+        particleCloud_.registerNamesFieldsUserCFDEMToExt(propsDict_.get<word>("partHeatFluxName"), partHeatFluxPositionInRegister_);
     }
     bool validspeciesFlux = 
         particleCloud_.checkAndregisterNamesFieldsUserCFDEMToExt(partSpeciesFluxNames_,

src/lagrangian/cfdemParticle/subModels/liggghtsCommandModel/liggghtsCommandModel/liggghtsCommandModel.C

@@ -176,7 +176,7 @@ bool liggghtsCommandModel::runThisCommand(int couplingStep)
     bool runIt=false;
     if(
        (!runEveryWriteStep_ && firstCouplingStep_  <= couplingStep && lastCouplingStep_  >= couplingStep)  ||
-       (runEveryWriteStep_  && particleCloud_.mesh().time().outputTime())
+       (runEveryWriteStep_  && particleCloud_.writeTimePassed())
       )
     {
         if(couplingStep >= nextRun_)

src/lagrangian/cfdemParticle/subModels/liggghtsCommandModel/liggghtsCommandModel/newLiggghtsCommandModel.C

@@ -52,7 +52,7 @@ autoPtr<liggghtsCommandModel> liggghtsCommandModel::New
     Info<< "Selecting liggghtsCommandModel "
          << liggghtsCommandType << " ,provide dicts, numbered from 0 to n" << endl;
 
-    dictionaryConstructorTable::iterator cstrIter =
+    auto cstrIter =
         dictionaryConstructorTablePtr_->find(liggghtsCommandType);
 
     if (cstrIter == dictionaryConstructorTablePtr_->end())

src/lagrangian/cfdemParticle/subModels/liggghtsCommandModel/readLiggghtsData/readLiggghtsData.C

@@ -75,8 +75,8 @@ readLiggghtsData::readLiggghtsData
     myName_=word(typeName + "Props" + h);
     propsDict_=dictionary(dict.subDict(myName_));
 
-    if (propsDict_.found("exactTiming"))
-        exactTiming_=true;
+    if (propsDict_.found("exactTiming")) // NOTE(zhyang)：只要找到这个关键字，结果就设置true; 不设置就是false
+        exactTiming_ = true;
     Info << "exactTiming==" << exactTiming_ << endl;
 
     if (propsDict_.found("verbose")) verbose_=true;

src/lagrangian/cfdemParticle/subModels/liggghtsCommandModel/setDEMGravity/setDEMGravity.C

@@ -69,11 +69,7 @@ setDEMGravity::setDEMGravity
     filePathList_(100),
     scalarList_(40),
     labelList_(1),
-    #if defined(version21) || defined(version16ext)
-        g_(sm.mesh().lookupObject<uniformDimensionedVectorField> ("g")),
-    #elif defined(version15)
-        g_(dimensionedVector(sm.mesh().lookupObject<IOdictionary>("environmentalProperties").lookup("g")).value()),
-    #endif
+    g_(sm.mesh().lookupObject<uniformDimensionedVectorField> ("g")),
     unfix_(false)
 {
     // define dictionary

src/lagrangian/cfdemParticle/subModels/liggghtsCommandModel/setDEMGravity/setDEMGravity.H

@@ -70,11 +70,7 @@ private:
 
     labelList labelList_;
 
-    #ifdef version21
-        const uniformDimensionedVectorField& g_;    // ref to gravity
-    #elif defined(version16ext) || defined(version15)
-        const dimensionedVector& g_;    // ref to gravity
-    #endif
+    const uniformDimensionedVectorField& g_;    // ref to gravity
 
     Switch unfix_;
 

src/lagrangian/cfdemParticle/subModels/locateModel/engineSearch/engineSearch.C

@@ -63,13 +63,7 @@ engineSearch::engineSearch
     locateModel(dict,sm),
     propsDict_(dict.subDict(name == "" ? typeName + "Props" : name + "Props")),
     treeSearch_(propsDict_.lookupOrDefault<Switch>("treeSearch", true)),
-    #if defined(version30)
-        searchEngine_(particleCloud_.mesh(),polyMesh::FACE_PLANES)
-    #elif defined(version21)
-        searchEngine_(particleCloud_.mesh(),polyMesh::FACEPLANES) // FACEPLANES or FACECENTRETETS; FACEDIAGTETS not stable
-    #elif defined(version16ext)
-        searchEngine_(particleCloud_.mesh(),false) //(particleCloud_.mesh(),faceDecomp_)
-    #endif
+    searchEngine_(particleCloud_.mesh(),polyMesh::FACE_PLANES)
 {}
 
 

src/lagrangian/cfdemParticle/subModels/locateModel/locateModel/newLocateModel.C

@@ -56,7 +56,7 @@ autoPtr<locateModel> locateModel::New
          << locateModelType << endl;
 
 
-    dictionaryConstructorTable::iterator cstrIter =
+    auto cstrIter =
         dictionaryConstructorTablePtr_->find(locateModelType);
 
     if (cstrIter == dictionaryConstructorTablePtr_->end())

src/lagrangian/cfdemParticle/subModels/locateModel/standardSearch/standardSearch.C

@@ -90,13 +90,7 @@ label standardSearch::findCell
             for(int i=0;i<3;i++) position[i] = positions[index][i];
 
             // find cell
-            #if defined(version30)
-                cellIDs[index][0] = particleCloud_.mesh().findCell(position, polyMesh::FACE_PLANES);
-            #elif defined(version21)
-                cellIDs[index][0] = particleCloud_.mesh().findCell(position, polyMesh::FACEPLANES);
-            #elif defined(version16ext)
-                cellIDs[index][0] = particleCloud_.mesh().findCell(position);
-            #endif
+            cellIDs[index][0] = particleCloud_.mesh().findCell(position, polyMesh::FACE_PLANES);
         }
         else cellIDs[index][0]=-1;
     }
@@ -111,13 +105,7 @@ label standardSearch::findSingleCell
 ) const
 {
     // find cell
-    #if defined(version30)
-        return particleCloud_.mesh().findCell(position, polyMesh::FACE_PLANES);
-    #elif defined(version21)
-        return particleCloud_.mesh().findCell(position, polyMesh::FACEPLANES);
-    #elif defined(version16ext)
-        return particleCloud_.mesh().findCell(position);
-    #endif
+    return particleCloud_.mesh().findCell(position, polyMesh::FACE_PLANES);
 }
 
 

src/lagrangian/cfdemParticle/subModels/meshMotionModel/meshMotionModel/newMeshMotionModel.C

@@ -56,7 +56,7 @@ autoPtr<meshMotionModel> meshMotionModel::New
          << meshMotionModelType << endl;
 
 
-    dictionaryConstructorTable::iterator cstrIter =
+    auto cstrIter =
         dictionaryConstructorTablePtr_->find(meshMotionModelType);
 
     if (cstrIter == dictionaryConstructorTablePtr_->end())

src/lagrangian/cfdemParticle/subModels/momCoupleModel/momCoupleModel/newMomCoupleModel.C

@@ -52,7 +52,7 @@ autoPtr<momCoupleModel> momCoupleModel::New
          << momCoupleType << endl;
 
 
-    dictionaryConstructorTable::iterator cstrIter =
+    auto cstrIter =
         dictionaryConstructorTablePtr_->find(momCoupleType);
 
     if (cstrIter == dictionaryConstructorTablePtr_->end())

src/lagrangian/cfdemParticle/subModels/probeModel/probeModel/newProbeModel.C

@@ -58,7 +58,7 @@ autoPtr<probeModel> probeModel::New
          << probeModelType << endl;
 
 
-    dictionaryConstructorTable::iterator cstrIter =
+    auto cstrIter =
         dictionaryConstructorTablePtr_->find(probeModelType);
 
     if (cstrIter == dictionaryConstructorTablePtr_->end())

src/lagrangian/cfdemParticle/subModels/registryModel/registryModel/newRegistryModel.C

@@ -56,7 +56,7 @@ autoPtr<registryModel> registryModel::New
          << registryModelType << endl;*/
 
 
-    dictionaryConstructorTable::iterator cstrIter =
+    auto cstrIter =
         dictionaryConstructorTablePtr_->find(registryModelType);
 
     if (cstrIter == dictionaryConstructorTablePtr_->end())

src/lagrangian/cfdemParticle/subModels/smoothingModel/constDiffSmoothing/constDiffSmoothing.C

@@ -199,7 +199,7 @@ void Foam::constDiffSmoothing::smoothenReferenceField(volVectorField& fieldSrc)
             IOobject::NO_WRITE
         ),
         particleCloud_.mesh(),
-        0.0
+        dimensionedScalar(dimless, Zero)
     );
 
 

src/lagrangian/cfdemParticle/subModels/smoothingModel/smoothingModel/newSmoothingModel.C

@@ -56,7 +56,7 @@ autoPtr<smoothingModel> smoothingModel::New
          << smoothingModelType << endl;
 
 
-    dictionaryConstructorTable::iterator cstrIter =
+    auto cstrIter =
         dictionaryConstructorTablePtr_->find(smoothingModelType);
 
     if (cstrIter == dictionaryConstructorTablePtr_->end())

src/lagrangian/cfdemParticle/subModels/voidFractionModel/voidFractionModel/newVoidFractionModel.C

@@ -56,7 +56,7 @@ autoPtr<voidFractionModel> voidFractionModel::New
          << voidFractionModelType << endl;
 
 
-    dictionaryConstructorTable::iterator cstrIter =
+    auto cstrIter =
         dictionaryConstructorTablePtr_->find(voidFractionModelType);
 
     if (cstrIter == dictionaryConstructorTablePtr_->end())

tutorials/cfdemSolverPiso/ErgunTestMPI/parCFDDEMrun.sh

@@ -23,7 +23,7 @@ nrProcs="4"
 machineFileName="none"   # yourMachinefileName | none
 debugMode="off"          # on | off| strict | profile
 testHarnessPath="$CFDEM_TEST_HARNESS_PATH"
-runOctave="true"
+runOctave="false"
 postproc="false"
 #--------------------------------------------------------------------------------#
 

tutorials/cfdemSolverPiso/ErgunTestMPI_cgs/parCFDDEMrun.sh

@@ -23,7 +23,7 @@ nrProcs="4"
 machineFileName="none"   # yourMachinefileName | none
 debugMode="off"          # on | off| strict
 testHarnessPath="$CFDEM_TEST_HARNESS_PATH"
-runOctave="true"
+runOctave="false"
 postproc="false"
 #--------------------------------------------------------------------------------#
 

tutorials/cfdemSolverPiso/periodicChannel/CFD/constant/couplingProperties

@@ -48,6 +48,8 @@ clockModel off;//standardClock;//
 
 smoothingModel off;// localPSizeDiffSmoothing;// constDiffSmoothing; //
 
+verbose on;
+
 forceModels
 (
     KochHillDrag

tutorials/cfdemSolverPiso/periodicChannel/CFD/system/controlDict

@@ -57,30 +57,27 @@ functions
     probes
     {
         type        probes;
-        // Where to load it from
-        functionObjectLibs ( "libsampling.so" );  
-        // Name of the directory for probe data
-        name        probes;
+        libs        ( "libsampling.so" );  // Where to load it from
+        name        probes; // Name of the directory for probe data
         probeLocations
         (
             (0.00003 0.050001 0.050001)
             (0.39999 0.050001 0.050001)
         );
 
-        // Fields to be probed
-        fields ( p U voidfraction );
+        fields          ( p U voidfraction ); // Fields to be probed
 
         // Write at same frequency as fields
-        outputControl   timeStep;//outputTime;
-        outputInterval  1;
+        writeControl    timeStep;
+        writeInterval   1;
     }
 
     volFlow_inlet
     {
         type            surfaceFieldValue;
-        functionObjectLibs ("libfieldFunctionObjects.so");
-        outputControl   timeStep;
-        outputInterval  1;
+        libs            ("libfieldFunctionObjects.so");
+        writeControl    timeStep;
+        writeInterval   1;
         log             false;
         writeFields     false;
         valueOutput     false;
@@ -96,9 +93,9 @@ functions
     volFlow_outlet
     {
         type            surfaceFieldValue;
-        functionObjectLibs ("libfieldFunctionObjects.so");
-        outputControl   timeStep;
-        outputInterval  1;
+        libs            ("libfieldFunctionObjects.so");
+        writeControl    timeStep;
+        writeInterval   1;
         log             false;
         writeFields     false;
         valueOutput     false;
@@ -114,9 +111,9 @@ functions
     volFlow_wall
     {
         type            surfaceFieldValue;
-        functionObjectLibs ("libfieldFunctionObjects.so");
-        outputControl   timeStep;
-        outputInterval  1;
+        libs            ("libfieldFunctionObjects.so");
+        writeControl    timeStep;
+        writeInterval   1;
         log             false;
         writeFields     false;
         valueOutput     false;
@@ -129,4 +126,3 @@ functions
         );
     }
 }
-// ************************************************************************* //

tutorials/cfdemSolverPiso/periodicChannel/parCFDDEMrun.sh

@@ -23,7 +23,7 @@ nrProcs="4"
 machineFileName="none"   # yourMachinefileName | none
 debugMode="off"          # on | off| strict
 testHarnessPath="$CFDEM_TEST_HARNESS_PATH"
-runOctave="true"
+runOctave="false"
 postproc="false"
 cleanCase="true"
 #--------------------------------------------------------------------------------#

tutorials/cfdemSolverPiso/settlingTestMPI/CFD/constant/couplingProperties

@@ -49,6 +49,8 @@ smoothingModel  constDiffSmoothing; //off;//
 
 //useDDTvoidfraction;
 
+verbose on;
+
 forceModels
 (
     DiFeliceDrag

tutorials/cfdemSolverPiso/settlingTestMPI/parCFDDEMrun.sh

@@ -23,7 +23,7 @@ nrProcs="2"
 machineFileName="none"   # yourMachinefileName | none
 debugMode="off"          # on | off| strict
 testHarnessPath="$CFDEM_TEST_HARNESS_PATH"
-runOctave="true"
+runOctave="false"
 cleanUp="true"
 postproc="false"
 #--------------------------------------------------------------------------------#

tutorials/cfdemSolverPiso/voidfractionTest/parCFDDEMrun.sh

@@ -24,7 +24,7 @@ machineFileName="none"   # yourMachinefileName | none
 debugMode="off"          # on | off| strict
 reconstuctCase="false"
 testHarnessPath="$CFDEM_TEST_HARNESS_PATH"
-runOctave="true"
+runOctave="false"
 cleanUp="true"
 postproc="false"
 #--------------------------------------------------------------------------------#