diff --git a/README b/README
deleted file mode 100755
index fd2ca756..00000000
--- a/README
+++ /dev/null
@@ -1,81 +0,0 @@
-/*---------------------------------------------------------------------------*\
-    CFDEMcoupling - Open Source CFD-DEM coupling
-
-    CFDEMcoupling is part of the CFDEMproject
-    www.cfdem.com
-                                Christoph Goniva, christoph.goniva@cfdem.com
-                                Copyright 2009-2012 JKU Linz
-                                Copyright 2012-     DCS Computing GmbH, Linz
--------------------------------------------------------------------------------
-License
-    This file is part of CFDEMcoupling.
-
-    CFDEMcoupling is free software; you can redistribute it and/or modify it
-    under the terms of the GNU General Public License as published by the
-    Free Software Foundation; either version 3 of the License, or (at your
-    option) any later version.
-
-    CFDEMcoupling is distributed in the hope that it will be useful, but WITHOUT
-    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
-    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
-    for more details.
-
-    You should have received a copy of the GNU General Public License
-    along with CFDEMcoupling; if not, write to the Free Software Foundation,
-    Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
-
-Description
-    This code is designed to realize coupled CFD-DEM simulations using LIGGGHTS
-    and OpenFOAM. Note: this code is not part of OpenFOAM (see DISCLAIMER).
-\*---------------------------------------------------------------------------*/
-
-
-CFDEM coupling provides an open source parallel coupled CFD-DEM framework 
-combining the strengths of LIGGGHTS DEM code and the Open Source 
-CFD package OpenFOAM(R)(*). The CFDEMcoupling toolbox allows to expand 
-standard CFD solvers of OpenFOAM(R)(*) 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"_lws 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 file structure:
-
-- "src" directory including the source files of the coupling toolbox and models
-- "applications" directory including the solver files for coupled CFD-DEM simulations
-- "doc" directory including the documentation of CFDEMcoupling
-- "tutorials" directory including basic tutorial cases showing the functionality
-
-
-
-Details on installation are given on the "www.cfdem.com"
-
-The functionality of this CFD-DEM framwork is described via "tutorial cases" showing 
-how to use different solvers and models.
-
-CFDEMcoupling stands for Computational Fluid Dynamics (CFD) -
-Discrete Element Method (DEM) coupling.
-
-CFDEMcoupling is an open-source code, distributed freely under the terms of the 
-GNU Public License (GPL).
-
-Core development of CFDEMcoupling is done by 
-Christoph Goniva and Christoph Kloss, both at DCS Computing GmbH, 2012
-
-
-\*---------------------------------------------------------------------------*/
-(*) "OpenFOAM(R)"_of is a registered trade mark of the ESI Group.
-This offering is not affiliated, approved or endorsed by ESI Group, 
-the producer of the OpenFOAM® software and owner of the OpenFOAM® trade mark.
-\*---------------------------------------------------------------------------*/
diff --git a/doc/CFDEMcoupling_Manual.html b/doc/CFDEMcoupling_Manual.html
index 0fc3dcb6..1f776e41 100644
--- a/doc/CFDEMcoupling_Manual.html
+++ b/doc/CFDEMcoupling_Manual.html
@@ -212,7 +212,7 @@ listing below of styles within certain commands.
 <TR ALIGN="center"><TD ><A HREF = "clockModel_noClock.html">clockModel_noClock</A></TD><TD ><A HREF = "clockModel_standardClock.html">clockModel_standardClock</A></TD></TR>
 <TR ALIGN="center"><TD ><A HREF = "dataExchangeModel.html">dataExchangeModel</A></TD><TD ><A HREF = "dataExchangeModel_noDataExchange.html">dataExchangeModel_noDataExchange</A></TD></TR>
 <TR ALIGN="center"><TD ><A HREF = "dataExchangeModel_oneWayVTK.html">dataExchangeModel_oneWayVTK</A></TD><TD ><A HREF = "dataExchangeModel_twoWayFiles.html">dataExchangeModel_twoWayFiles</A></TD></TR>
-<TR ALIGN="center"><TD ><A HREF = "dataExchangeModel_twoWayM2M.html">dataExchangeModel_twoWayM2M</A></TD><TD ><A HREF = "dataExchangeModel_twoWayMPI.html">dataExchangeModel_twoWayMPI</A></TD></TR>
+<TR ALIGN="center"><TD ><A HREF = "dataExchangeModel_twoWayMPI.html">dataExchangeModel_twoWayMPI</A></TD><TD ><A HREF = "dataExchangeModel_twoWayMany2Many.html">dataExchangeModel_twoWayMany2Many</A></TD></TR>
 <TR ALIGN="center"><TD ><A HREF = "forceModel.html">forceModel</A></TD><TD ><A HREF = "forceModel_Archimedes.html">forceModel_Archimedes</A></TD></TR>
 <TR ALIGN="center"><TD ><A HREF = "forceModel_ArchimedesIB.html">forceModel_ArchimedesIB</A></TD><TD ><A HREF = "forceModel_DiFeliceDrag.html">forceModel_DiFeliceDrag</A></TD></TR>
 <TR ALIGN="center"><TD ><A HREF = "forceModel_GidaspowDrag.html">forceModel_GidaspowDrag</A></TD><TD ><A HREF = "forceModel_KochHillDrag.html">forceModel_KochHillDrag</A></TD></TR>
@@ -225,17 +225,16 @@ listing below of styles within certain commands.
 <TR ALIGN="center"><TD ><A HREF = "liggghtsCommandModel_runLiggghts.html">liggghtsCommandModel_runLiggghts</A></TD><TD ><A HREF = "liggghtsCommandModel_writeLiggghts.html">liggghtsCommandModel_writeLiggghts</A></TD></TR>
 <TR ALIGN="center"><TD ><A HREF = "locateModel.html">locateModel</A></TD><TD ><A HREF = "locateModel_engineSearch.html">locateModel_engineSearch</A></TD></TR>
 <TR ALIGN="center"><TD ><A HREF = "locateModel_engineSearchIB.html">locateModel_engineSearchIB</A></TD><TD ><A HREF = "locateModel_standardSearch.html">locateModel_standardSearch</A></TD></TR>
-<TR ALIGN="center"><TD ><A HREF = "locateModel_turboEngineSearchM2M.html">locateModel_turboEngineM2MSearch</A></TD><TD ><A HREF = "locateModel_turboEngineSearch.html">locateModel_turboEngineSearch</A></TD></TR>
-<TR ALIGN="center"><TD ><A HREF = "meshMotionModel.html">meshMotionModel</A></TD><TD ><A HREF = "meshMotionModel_noMeshMotion.html">meshMotionModel_noMeshMotion</A></TD></TR>
-<TR ALIGN="center"><TD ><A HREF = "momCoupleModel.html">momCoupleModel</A></TD><TD ><A HREF = "momCoupleModel_explicitCouple.html">momCoupleModel_explicitCouple</A></TD></TR>
-<TR ALIGN="center"><TD ><A HREF = "momCoupleModel_implicitCouple.html">momCoupleModel_implicitCouple</A></TD><TD ><A HREF = "momCoupleModel_noCouple.html">momCoupleModel_noCouple</A></TD></TR>
-<TR ALIGN="center"><TD ><A HREF = "probeModel.html">probeModel</A></TD><TD ><A HREF = "probeModel_noProbe.html">probeModel_noProbe</A></TD></TR>
-<TR ALIGN="center"><TD ><A HREF = "regionModel.html">regionModel</A></TD><TD ><A HREF = "regionModel_allRegion.html">regionModel_allRegion</A></TD></TR>
-<TR ALIGN="center"><TD ><A HREF = "smoothingModel.html">smoothingModel</A></TD><TD ><A HREF = "smoothingModel_constDiffSmoothing.html">smoothingModel_constDiffSmoothing</A></TD></TR>
-<TR ALIGN="center"><TD ><A HREF = "smoothingModel_noSmoothing.html">smoothingModel_noSmoothing</A></TD><TD ><A HREF = "voidFractionModel.html">voidfractionModel</A></TD></TR>
-<TR ALIGN="center"><TD ><A HREF = "voidFractionModel_GaussVoidFraction.html">voidfractionModel_GaussVoidFraction</A></TD><TD ><A HREF = "voidFractionModel_IBVoidFraction.html">voidfractionModel_IBVoidFraction</A></TD></TR>
-<TR ALIGN="center"><TD ><A HREF = "voidFractionModel_bigParticleVoidFraction.html">voidfractionModel_bigParticleVoidFraction</A></TD><TD ><A HREF = "voidFractionModel_centreVoidFraction.html">voidfractionModel_centreVoidFraction</A></TD></TR>
-<TR ALIGN="center"><TD ><A HREF = "voidFractionModel_dividedVoidFraction.html">voidfractionModel_dividedVoidFraction</A> 
+<TR ALIGN="center"><TD ><A HREF = "locateModel_turboEngineSearch.html">locateModel_turboEngineSearch</A></TD><TD ><A HREF = "meshMotionModel.html">meshMotionModel</A></TD></TR>
+<TR ALIGN="center"><TD ><A HREF = "meshMotionModel_noMeshMotion.html">meshMotionModel_noMeshMotion</A></TD><TD ><A HREF = "momCoupleModel.html">momCoupleModel</A></TD></TR>
+<TR ALIGN="center"><TD ><A HREF = "momCoupleModel_explicitCouple.html">momCoupleModel_explicitCouple</A></TD><TD ><A HREF = "momCoupleModel_implicitCouple.html">momCoupleModel_implicitCouple</A></TD></TR>
+<TR ALIGN="center"><TD ><A HREF = "momCoupleModel_noCouple.html">momCoupleModel_noCouple</A></TD><TD ><A HREF = "probeModel.html">probeModel</A></TD></TR>
+<TR ALIGN="center"><TD ><A HREF = "probeModel_noProbe.html">probeModel_noProbe</A></TD><TD ><A HREF = "regionModel.html">regionModel</A></TD></TR>
+<TR ALIGN="center"><TD ><A HREF = "regionModel_allRegion.html">regionModel_allRegion</A></TD><TD ><A HREF = "smoothingModel.html">smoothingModel</A></TD></TR>
+<TR ALIGN="center"><TD ><A HREF = "smoothingModel_constDiffSmoothing.html">smoothingModel_constDiffSmoothing</A></TD><TD ><A HREF = "smoothingModel_noSmoothing.html">smoothingModel_noSmoothing</A></TD></TR>
+<TR ALIGN="center"><TD ><A HREF = "voidFractionModel.html">voidfractionModel</A></TD><TD ><A HREF = "voidFractionModel_GaussVoidFraction.html">voidfractionModel_GaussVoidFraction</A></TD></TR>
+<TR ALIGN="center"><TD ><A HREF = "voidFractionModel_IBVoidFraction.html">voidfractionModel_IBVoidFraction</A></TD><TD ><A HREF = "voidFractionModel_bigParticleVoidFraction.html">voidfractionModel_bigParticleVoidFraction</A></TD></TR>
+<TR ALIGN="center"><TD ><A HREF = "voidFractionModel_centreVoidFraction.html">voidfractionModel_centreVoidFraction</A></TD><TD ><A HREF = "voidFractionModel_dividedVoidFraction.html">voidfractionModel_dividedVoidFraction</A> 
 </TD></TR></TABLE></DIV>
 
 </HTML>
diff --git a/doc/CFDEMcoupling_Manual.pdf b/doc/CFDEMcoupling_Manual.pdf
index 9eeb4519..0bf492b9 100644
Binary files a/doc/CFDEMcoupling_Manual.pdf and b/doc/CFDEMcoupling_Manual.pdf differ
diff --git a/doc/CFDEMcoupling_Manual.txt b/doc/CFDEMcoupling_Manual.txt
index 69f8f64e..612fbedb 100644
--- a/doc/CFDEMcoupling_Manual.txt
+++ b/doc/CFDEMcoupling_Manual.txt
@@ -250,8 +250,8 @@ listing below of styles within certain commands.
 "dataExchangeModel_noDataExchange"_dataExchangeModel_noDataExchange.html,
 "dataExchangeModel_oneWayVTK"_dataExchangeModel_oneWayVTK.html,
 "dataExchangeModel_twoWayFiles"_dataExchangeModel_twoWayFiles.html,
-"dataExchangeModel_twoWayM2M"_dataExchangeModel_twoWayM2M.html,
 "dataExchangeModel_twoWayMPI"_dataExchangeModel_twoWayMPI.html,
+"dataExchangeModel_twoWayMany2Many"_dataExchangeModel_twoWayMany2Many.html,
 "forceModel"_forceModel.html,
 "forceModel_Archimedes"_forceModel_Archimedes.html,
 "forceModel_ArchimedesIB"_forceModel_ArchimedesIB.html,
@@ -276,7 +276,6 @@ listing below of styles within certain commands.
 "locateModel_engineSearch"_locateModel_engineSearch.html,
 "locateModel_engineSearchIB"_locateModel_engineSearchIB.html,
 "locateModel_standardSearch"_locateModel_standardSearch.html,
-"locateModel_turboEngineM2MSearch"_locateModel_turboEngineSearchM2M.html,
 "locateModel_turboEngineSearch"_locateModel_turboEngineSearch.html,
 "meshMotionModel"_meshMotionModel.html,
 "meshMotionModel_noMeshMotion"_meshMotionModel_noMeshMotion.html,
diff --git a/doc/dataExchangeModel_twoWayM2M.txt b/doc/dataExchangeModel_twoWayM2M.txt
deleted file mode 100644
index 94527637..00000000
--- a/doc/dataExchangeModel_twoWayM2M.txt
+++ /dev/null
@@ -1,42 +0,0 @@
-"CFDEMproject WWW Site"_lws - "CFDEM Commands"_lc :c
-
-:link(lws,http://www.cfdem.com)
-:link(lc,CFDEMcoupling_Manual.html#comm)
-
-:line
-
-dataExchangeModel_twoWayM2M command :h3
-
-[Syntax:]
-
-Defined in couplingProperties dictionary.
-
-dataExchangeModel twoWayM2M;
-twoWayM2MProps
-\{
-    liggghtsPath "path";
-\}; :pre
-
-{path} = path to the DEM simulation input file :ulb,l
-:ule
-
-[Examples:]
-
-dataExchangeModel twoWayM2M;
-twoWayM2MProps
-\{
-    liggghtsPath "../DEM/in.liggghts_init";
-\} :pre
-
-[Description:]
-
-The data exchange model performs the data exchange between the DEM code and the CFD code. The twoWayM2M model is a model that can exchange particle properties from DEM to CFD and from CFD to DEM. Data is exchanged via MPI technique using the many to many mapping scheme. The DEM run is executed by the coupling model, via a liggghtsCommandModel object.
-
-[Restrictions:]
-
-Should be used in combination with the turboEngineSearchM2M locate model to achieve best performance!
-
-[Related commands:]
-
-"dataExchangeModel"_dataExchangeModel.html
-
diff --git a/doc/dataExchangeModel_twoWayM2M.html b/doc/dataExchangeModel_twoWayMany2Many.html
similarity index 54%
rename from doc/dataExchangeModel_twoWayM2M.html
rename to doc/dataExchangeModel_twoWayMany2Many.html
index a80fd095..db5b6592 100644
--- a/doc/dataExchangeModel_twoWayM2M.html
+++ b/doc/dataExchangeModel_twoWayMany2Many.html
@@ -7,14 +7,14 @@
 
 <HR>
 
-<H3>dataExchangeModel_twoWayM2M command 
+<H3>dataExchangeModel_twoWayMany2Many command 
 </H3>
 <P><B>Syntax:</B>
 </P>
 <P>Defined in couplingProperties dictionary.
 </P>
-<PRE>dataExchangeModel twoWayM2M;
-twoWayM2MProps
+<PRE>dataExchangeModel twoWayMany2Many;
+twoWayMany2ManyProps
 {
     liggghtsPath "path";
 }; 
@@ -25,19 +25,19 @@ twoWayM2MProps
 </UL>
 <P><B>Examples:</B>
 </P>
-<PRE>dataExchangeModel twoWayM2M;
-twoWayM2MProps
+<PRE>dataExchangeModel twoWayMany2Many;
+twoWayMany2ManyProps
 {
     liggghtsPath "../DEM/in.liggghts_init";
 } 
 </PRE>
 <P><B>Description:</B>
 </P>
-<P>The data exchange model performs the data exchange between the DEM code and the CFD code. The twoWayM2M model is a model that can exchange particle properties from DEM to CFD and from CFD to DEM. Data is exchanged via MPI technique using the many to many mapping scheme. The DEM run is executed by the coupling model, via a liggghtsCommandModel object.
+<P>The data exchange model performs the data exchange between the DEM code and the CFD code. The twoWayMany2Many model is a model that can exchange particle properties from DEM to CFD and from CFD to DEM. Data is exchanged via MPI technique using the many to many mapping scheme. The DEM run is executed by the coupling model, via a liggghtsCommandModel object.
 </P>
 <P><B>Restrictions:</B>
 </P>
-<P>Should be used in combination with the turboEngineSearchM2M locate model to achieve best performance!
+<P>Must be used in combination with the engineSearchMany2Many locate model!
 </P>
 <P><B>Related commands:</B>
 </P>
diff --git a/doc/dataExchangeModel_twoWayMany2Many.txt b/doc/dataExchangeModel_twoWayMany2Many.txt
new file mode 100644
index 00000000..31f87e67
--- /dev/null
+++ b/doc/dataExchangeModel_twoWayMany2Many.txt
@@ -0,0 +1,42 @@
+"CFDEMproject WWW Site"_lws - "CFDEM Commands"_lc :c
+
+:link(lws,http://www.cfdem.com)
+:link(lc,CFDEMcoupling_Manual.html#comm)
+
+:line
+
+dataExchangeModel_twoWayMany2Many command :h3
+
+[Syntax:]
+
+Defined in couplingProperties dictionary.
+
+dataExchangeModel twoWayMany2Many;
+twoWayMany2ManyProps
+\{
+    liggghtsPath "path";
+\}; :pre
+
+{path} = path to the DEM simulation input file :ulb,l
+:ule
+
+[Examples:]
+
+dataExchangeModel twoWayMany2Many;
+twoWayMany2ManyProps
+\{
+    liggghtsPath "../DEM/in.liggghts_init";
+\} :pre
+
+[Description:]
+
+The data exchange model performs the data exchange between the DEM code and the CFD code. The twoWayMany2Many model is a model that can exchange particle properties from DEM to CFD and from CFD to DEM. Data is exchanged via MPI technique using the many to many mapping scheme. The DEM run is executed by the coupling model, via a liggghtsCommandModel object.
+
+[Restrictions:]
+
+Must be used in combination with the engineSearchMany2Many locate model!
+
+[Related commands:]
+
+"dataExchangeModel"_dataExchangeModel.html
+
diff --git a/doc/githubAccess_public.pdf b/doc/githubAccess_public.pdf
index a12c932e..8daebd98 100644
Binary files a/doc/githubAccess_public.pdf and b/doc/githubAccess_public.pdf differ
diff --git a/doc/locateModel_turboEngineSearchM2M.html b/doc/locateModel_turboEngineSearchM2M.html
deleted file mode 100644
index 5082aaa7..00000000
--- a/doc/locateModel_turboEngineSearchM2M.html
+++ /dev/null
@@ -1,58 +0,0 @@
-<HTML>
-<CENTER><A HREF = "http://www.cfdem.com">CFDEMproject WWW Site</A> - <A HREF = "CFDEMcoupling_Manual.html#comm">CFDEM Commands</A> 
-</CENTER>
-
-
-
-
-<HR>
-
-<H3>locateModel_turboEngineM2MSearch command 
-</H3>
-<P><B>Syntax:</B>
-</P>
-<P>Defined in couplingProperties dictionary.
-</P>
-<PRE>locateModel turboEngineM2M;
-turboEngineM2MProps
-{
-    turboEngineProps
-    {
-        treeSearch switch1;
-    }
-} 
-</PRE>
-<UL><LI><I>switch1</I> = switch to use tree search algorithm 
-
-
-</UL>
-<P><B>Examples:</B>
-</P>
-<PRE>locateModel turboEngineM2M;
-turboEngineM2MProps
-{
-    turboEngineProps
-    {
-        treeSearch true;
-    }
-} 
-</PRE>
-<P><B>Description:</B>
-</P>
-<P>The locateModel "turboEngineM2M" locates the CFD cell and cellID corresponding to a given position. The algorithm is improved compared to engine search to show better parallel performance.
-</P>
-<P>The turboEngineM2MSearch locate Model can be used with different settings to use different algorithms:
-</P>
-<UL><LI>faceDecomp false; treeSearch false;     will execute some geometric (linear) search using the last known cellID (recommended) 
-
-<LI>faceDecomp false; treeSearch true;      will use a recursive tree structure to find the cell. 
-
-
-</UL>
-<P><B>Restrictions:</B> This model can only be used with many to many data exchange model!
-</P>
-<P><B>Related commands:</B>
-</P>
-<P><A HREF = "locateModel.html">locateModel</A>
-</P>
-</HTML>
diff --git a/doc/locateModel_turboEngineSearchM2M.txt b/doc/locateModel_turboEngineSearchM2M.txt
deleted file mode 100644
index c126f9f3..00000000
--- a/doc/locateModel_turboEngineSearchM2M.txt
+++ /dev/null
@@ -1,52 +0,0 @@
-"CFDEMproject WWW Site"_lws - "CFDEM Commands"_lc :c
-
-:link(lws,http://www.cfdem.com)
-:link(lc,CFDEMcoupling_Manual.html#comm)
-
-:line
-
-locateModel_turboEngineM2MSearch command :h3
-
-[Syntax:]
-
-Defined in couplingProperties dictionary.
-
-locateModel turboEngineM2M;
-turboEngineM2MProps
-\{
-    turboEngineProps
-    \{
-        treeSearch switch1;
-    \}
-\} :pre
-
-{switch1} = switch to use tree search algorithm :ulb,l
-:ule
-
-[Examples:]
-
-locateModel turboEngineM2M;
-turboEngineM2MProps
-\{
-    turboEngineProps
-    \{
-        treeSearch true;
-    \}
-\} :pre
-
-[Description:]
-
-The locateModel "turboEngineM2M" locates the CFD cell and cellID corresponding to a given position. The algorithm is improved compared to engine search to show better parallel performance.
-
-The turboEngineM2MSearch locate Model can be used with different settings to use different algorithms:
-
-faceDecomp false; treeSearch false;     will execute some geometric (linear) search using the last known cellID (recommended) :ulb,l
-faceDecomp false; treeSearch true;      will use a recursive tree structure to find the cell. :l
-:ule
-
-[Restrictions:] This model can only be used with many to many data exchange model!
-
-[Related commands:]
-
-"locateModel"_locateModel.html
-
diff --git a/src/lagrangian/cfdemParticle/Make/files b/src/lagrangian/cfdemParticle/Make/files
index e8798d13..3ce050bb 100644
--- a/src/lagrangian/cfdemParticle/Make/files
+++ b/src/lagrangian/cfdemParticle/Make/files
@@ -3,10 +3,10 @@ forceModels = subModels/forceModel
 forceModelsMS = subModels/forceModelMS
 IOModels = subModels/IOModel
 voidFractionModels = subModels/voidFractionModel
+voidFractionModelsMS = subModels/voidFractionModelMS
 locateModels = subModels/locateModel
 meshMotionModels = subModels/meshMotionModel
 momCoupleModels = subModels/momCoupleModel
-regionModels = subModels/regionModel
 dataExchangeModels = subModels/dataExchangeModel
 averagingModels = subModels/averagingModel
 clockModels = subModels/clockModel
@@ -15,6 +15,7 @@ smoothingModels = subModels/smoothingModel
 probeModels = subModels/probeModel
 
 $(cfdemCloud)/cfdemCloud.C
+derived/cfdemCloudBiDisperse/cfdemCloudBiDisperse.C
 derived/cfdemCloudIB/cfdemCloudIB.C
 derived/cfdemCloudMS/cfdemCloudMS.C
 
@@ -23,36 +24,59 @@ $(forceModels)/forceModel/newForceModel.C
 $(forceModels)/noDrag/noDrag.C
 $(forceModels)/checkCouplingInterval/checkCouplingInterval.C
 $(forceModels)/DiFeliceDrag/DiFeliceDrag.C
+$(forceModels)/DiFeliceDragNLift/DiFeliceDragNLift.C
 $(forceModels)/GidaspowDrag/GidaspowDrag.C
 $(forceModels)/SchillerNaumannDrag/SchillerNaumannDrag.C
 $(forceModels)/Archimedes/Archimedes.C
 $(forceModels)/ArchimedesIB/ArchimedesIB.C
 $(forceModels)/interface/interface.C
 $(forceModels)/ShirgaonkarIB/ShirgaonkarIB.C
+$(forceModels)/interfaceParticleProbe/interfaceParticleProbe.C
+$(forceModels)/fieldStore/fieldStore.C
+$(forceModels)/fieldTimeAverage/fieldTimeAverage.C
+$(forceModels)/fieldBound/fieldBound.C
+$(forceModels)/volWeightedAverage/volWeightedAverage.C
+$(forceModels)/totalMomentumExchange/totalMomentumExchange.C
 $(forceModels)/KochHillDrag/KochHillDrag.C
 $(forceModels)/KochHillRWDrag/KochHillRWDrag.C
+$(forceModels)/BeetstraDrag/multiphaseFlowBasic/multiphaseFlowBasic.C
+$(forceModels)/BeetstraDrag/BeetstraDrag.C
+$(forceModels)/LaEuScalarLiquid/LaEuScalarLiquid.C
 $(forceModels)/LaEuScalarTemp/LaEuScalarTemp.C
+$(forceModels)/LaEuScalarDust/LaEuScalarDust.C
 $(forceModels)/virtualMassForce/virtualMassForce.C
 $(forceModels)/gradPForce/gradPForce.C
+$(forceModels)/gradULiftForce/gradULiftForce.C
+$(forceModels)/HollowayDrag/HollowayDrag.C
 $(forceModels)/viscForce/viscForce.C
 $(forceModels)/MeiLift/MeiLift.C
+$(forceModels)/melting/melting.C
+$(forceModels)/KochHillDragNLift/KochHillDragNLift.C
+/*$(forceModels)/stokesSpheroidDrag/stokesSpheroidDrag.C*/
+$(forceModels)/solidsPressureForce/solidsPressureForce.C
+$(forceModels)/periodicPressure/periodicPressure.C
+$(forceModels)/periodicPressureControl/periodicPressureControl.C
+$(forceModels)/averageSlipVel/averageSlipVel.C
 $(forceModels)/particleCellVolume/particleCellVolume.C
-$(forceModels)/fieldTimeAverage/fieldTimeAverage.C
-$(forceModels)/volWeightedAverage/volWeightedAverage.C
 
 $(forceModelsMS)/forceModelMS/forceModelMS.C
 $(forceModelsMS)/forceModelMS/newForceModelMS.C
 $(forceModelsMS)/DiFeliceDragMS/DiFeliceDragMS.C
+$(forceModelsMS)/GidaspowDragMS/GidaspowDragMS.C
+$(forceModelsMS)/noDragMS/noDragMS.C
 
 $(probeModels)/probeModel/probeModel.C
 $(probeModels)/probeModel/newProbeModel.C
 $(probeModels)/noProbe/noProbe.C
 $(probeModels)/particleProbe/particleProbe.C
+/*$(probeModels)/interfaceParticleProbe/interfaceParticleProbe.C*/
 
 $(IOModels)/IOModel/IOModel.C
 $(IOModels)/IOModel/newIOModel.C
 $(IOModels)/noIO/noIO.C
 $(IOModels)/basicIO/basicIO.C
+$(IOModels)/tempIO/tempIO.C
+$(IOModels)/colorIO/colorIO.C
 $(IOModels)/trackIO/trackIO.C
 $(IOModels)/sophIO/sophIO.C
 
@@ -60,34 +84,39 @@ $(voidFractionModels)/voidFractionModel/voidFractionModel.C
 $(voidFractionModels)/voidFractionModel/newVoidFractionModel.C
 $(voidFractionModels)/centreVoidFraction/centreVoidFraction.C
 $(voidFractionModels)/dividedVoidFraction/dividedVoidFraction.C
+$(voidFractionModels)/dividedVoidFractionBiDi/dividedVoidFractionBiDi.C
 $(voidFractionModels)/dividedVoidFractionMS/dividedVoidFractionMS.C
 $(voidFractionModels)/bigParticleVoidFraction/bigParticleVoidFraction.C
 $(voidFractionModels)/GaussVoidFraction/GaussVoidFraction.C
 $(voidFractionModels)/IBVoidFraction/IBVoidFraction.C
+$(voidFractionModels)/weightedNeigbhorsVoidFraction/weightedNeigbhorsVoidFraction.C
+
+$(voidFractionModelsMS)/voidFractionModelMS/voidFractionModelMS.C
+$(voidFractionModelsMS)/voidFractionModelMS/newVoidFractionModelMS.C
 
 $(locateModels)/locateModel/locateModel.C
 $(locateModels)/locateModel/newLocateModel.C
 $(locateModels)/standardSearch/standardSearch.C
 $(locateModels)/engineSearch/engineSearch.C
+$(locateModels)/engineSearchMany2Many/engineSearchMany2Many.C
 $(locateModels)/turboEngineSearch/turboEngineSearch.C
 $(locateModels)/turboEngineSearchM2M/turboEngineSearchM2M.C
 $(locateModels)/engineSearchIB/engineSearchIB.C
-
+$(locateModels)/hyperEngineSearch/hyperEngineSearch.C
+$(locateModels)/ijkSearch/ijkSearch.C
 
 $(meshMotionModels)/meshMotionModel/meshMotionModel.C
 $(meshMotionModels)/meshMotionModel/newMeshMotionModel.C
 $(meshMotionModels)/noMeshMotion/noMeshMotion.C
+$(meshMotionModels)/DEMdrivenMeshMotion/DEMdrivenMeshMotion.C
 
 $(momCoupleModels)/momCoupleModel/momCoupleModel.C
 $(momCoupleModels)/momCoupleModel/newMomCoupleModel.C
 $(momCoupleModels)/explicitCouple/explicitCouple.C
+$(momCoupleModels)/explicitCoupleSource/explicitCoupleSource.C
 $(momCoupleModels)/implicitCouple/implicitCouple.C
 $(momCoupleModels)/noCouple/noCouple.C
 
-$(regionModels)/regionModel/regionModel.C
-$(regionModels)/regionModel/newRegionModel.C
-$(regionModels)/allRegion/allRegion.C
-
 $(dataExchangeModels)/dataExchangeModel/dataExchangeModel.C
 $(dataExchangeModels)/dataExchangeModel/newDataExchangeModel.C
 $(dataExchangeModels)/oneWayVTK/oneWayVTK.C
@@ -95,11 +124,13 @@ $(dataExchangeModels)/twoWayFiles/twoWayFiles.C
 $(dataExchangeModels)/noDataExchange/noDataExchange.C
 $(dataExchangeModels)/twoWayMPI/twoWayMPI.C
 $(dataExchangeModels)/twoWayM2M/twoWayM2M.C
+$(dataExchangeModels)/twoWayMany2Many/twoWayMany2Many.C
 
 $(averagingModels)/averagingModel/averagingModel.C
 $(averagingModels)/averagingModel/newAveragingModel.C
 $(averagingModels)/dilute/dilute.C
 $(averagingModels)/dense/dense.C
+$(averagingModels)/denseBiDi/denseBiDi.C
 
 $(clockModels)/clockModel/clockModel.C
 $(clockModels)/clockModel/newClockModel.C
@@ -108,6 +139,7 @@ $(clockModels)/noClock/noClock.C
 
 $(liggghtsCommandModels)/liggghtsCommandModel/liggghtsCommandModel.C
 $(liggghtsCommandModels)/liggghtsCommandModel/newLiggghtsCommandModel.C
+$(liggghtsCommandModels)/colorParticles/colorParticles.C
 $(liggghtsCommandModels)/execute/execute.C
 $(liggghtsCommandModels)/runLiggghts/runLiggghts.C
 $(liggghtsCommandModels)/writeLiggghts/writeLiggghts.C
@@ -117,5 +149,6 @@ $(smoothingModels)/smoothingModel/smoothingModel.C
 $(smoothingModels)/smoothingModel/newSmoothingModel.C
 $(smoothingModels)/noSmoothing/noSmoothing.C
 $(smoothingModels)/constDiffSmoothing/constDiffSmoothing.C
+$(smoothingModels)/localPSizeDiffSmoothing/localPSizeDiffSmoothing.C
 
 LIB = $(CFDEM_LIB_DIR)/lib$(CFDEM_LIB_NAME)
diff --git a/src/lagrangian/cfdemParticle/Make/options b/src/lagrangian/cfdemParticle/Make/options
index e3ef2b01..4d01d3d8 100644
--- a/src/lagrangian/cfdemParticle/Make/options
+++ b/src/lagrangian/cfdemParticle/Make/options
@@ -14,11 +14,12 @@ EXE_INC = \
     -I$(LIB_SRC)/OpenFOAM/containers/HashTables/labelHashSet \
     -I$(CFDEM_LIGGGHTS_SRC_DIR) \
     -I$(CFDEM_M2MLIB_PATH) \
+    -I$(CFDEM_Many2ManyLIB_PATH) \
     -I$(CFDEM_SRC_DIR)/cfdTools \
 
 LIB_LIBS = \
      $(PLIBS) \
-    -L$(CFDEM_LIB_DIR) \
+    -L$(FOAM_USER_LIBBIN) \
     -lfiniteVolume \
     -lincompressibleRASModels \
     -lincompressibleLESModels \
@@ -28,4 +29,10 @@ LIB_LIBS = \
     -L$(CFDEM_LIGGGHTS_SRC_DIR) \
     -Wl,--whole-archive -l$(CFDEM_LIGGGHTS_LIB_NAME) -Wl,--no-whole-archive \
     -L$(CFDEM_M2MLIB_PATH) \
-    -lcouple
+    -lcouple \
+    -L$(CFDEM_Many2ManyLIB_PATH) \
+    -lcoupleMany2Many \
+
+/* add -I$(CFDEM_POEMSLIB_PATH) \ to EXE_INC */
+/*    -L$(CFDEM_POEMSLIB_PATH) \ */
+/*    -lpoems   */
diff --git a/src/lagrangian/cfdemParticle/cfdTools/mathExtra.H b/src/lagrangian/cfdemParticle/cfdTools/mathExtra.H
new file mode 100644
index 00000000..a7c6067d
--- /dev/null
+++ b/src/lagrangian/cfdemParticle/cfdTools/mathExtra.H
@@ -0,0 +1,265 @@
+/* ----------------------------------------------------------------------
+    CFDEMcoupling - Open Source CFD-DEM coupling
+
+    CFDEMcoupling is part of the CFDEMproject
+    www.cfdem.com
+                                Christoph Goniva, christoph.goniva@cfdem.com
+                                Copyright 2009-2012 JKU Linz
+                                Copyright 2012-     DCS Computing GmbH, Linz
+-------------------------------------------------------------------------------
+License
+    This file is part of CFDEMcoupling.
+
+    CFDEMcoupling is free software; you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by the
+    Free Software Foundation; either version 3 of the License, or (at your
+    option) any later version.
+
+    CFDEMcoupling is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with CFDEMcoupling; if not, write to the Free Software Foundation,
+    Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+Description
+    This code is designed to realize coupled CFD-DEM simulations using LIGGGHTS
+    and OpenFOAM(R). Note: this code is not part of OpenFOAM(R) (see DISCLAIMER).
+    
+    Copyright of this contribution:
+    Copyright 2014-     TU Graz, IPPT
+------------------------------------------------------------------------- */
+
+#ifndef CFDEM_MATH_EXTRA_H
+#define CFDEM_MATH_EXTRA_H
+
+#include <vector>
+//#include "math.h"
+#include "stdio.h"
+#include "string.h"
+#include "error.h"
+#include "ctype.h"
+
+#define TOLERANCE_ORTHO 1e-10
+
+namespace MathExtra 
+{
+
+//  inline void     outerProduct(double *vec1, double *vec2, double **m);
+//  inline double  spheroidGeometry(int index, double bi, double ai);
+
+
+
+//--------------------------------------------------------------------
+//   Outer Product of two vectors
+inline void outerProduct(double *vec1, double *vec2, double **m)
+{
+  int i, j;
+  //debug output
+  for( i = 0; i < 3; ++i )
+//    printf("OUTER PRODUCT: Input: vec1 element %d = %g", i, vec1[i]);
+  for( i = 0; i < 3; ++i )
+//    printf("OUTER PRODUCT: Input: vec2 element %d=%g", i, vec2[i]);
+  
+  //calculation
+  for( i = 0; i < 3; ++i )
+    for( j = 0; j < 3; ++j )
+    {
+      m[i][j] = vec1[i] * vec2[j];
+      printf("OUTER PRODUCT: Result: m[%d][%d]=%g", i, j, m[i][j]);
+    }
+    
+    
+}
+//--------------------------------------------------------------------
+// Compute the major, minor axis and eccentricity parameters of a prolate spheroid
+inline bool spheroidGeometry(double radius, double aspectRatio,  //inputs
+                                 double& ai,    double& bi,              //outputs
+                                 double& ei,    double& Le              //outputs
+                                ) //
+{
+    //INPUT
+    //  radius      ...volume-equivalent radius of the spheroid
+    //  aspectRatio ...major/minor aspect ratio
+    
+    //OUTPUT
+    //  ai  ... 
+    //  bi  ... 
+    //  ei  ... 
+    //  Le  ... 
+  
+    if(radius<=0.0)     //avoid troubles in case radius is 0 or negative
+        return false;
+
+ 	ai = radius * std::pow(aspectRatio*aspectRatio,0.33333333333333333333333);
+    bi = ai / aspectRatio;
+ 	ei = std::sqrt( 
+ 	                1.0
+                   - 1.0 / (aspectRatio*aspectRatio)
+ 	              );
+ 	Le = std::log(
+ 		           (1.0+ei)
+ 		          /(1.0-ei)
+                 );
+
+	return true; 		    
+}
+
+//--------------------------------------------------------------------
+// Compute the major, minor axis and eccentricity parameters of a prolate spheroid
+inline double Pi()
+{
+    return 3.1415926535897932384626433832795;
+}
+
+
+//--------------------------------------------------------------------
+// Compute the major, minor axis and eccentricity parameters of a prolate spheroid
+inline bool spheroidGeometry2(double radius, double aspectRatio,  //inputs
+                                 double& ai,     double& bi,           //outputs
+                                 double& XAe,    double& YAe,         //outputs
+                                 double& XCe,    double& YCe,         //outputs
+                                 double& YHe
+                                ) //
+{
+    //INPUT
+    //  radius      ...volume-equivalent radius of the spheroid
+    //  aspectRatio ...major/minor aspect ratio
+    
+    //OUTPUT
+    //  XAe  ...Eccentricity dependet parameter
+    //  YAe  ...Eccentricity dependet parameter
+    //  XCe  ...Eccentricity dependet parameter
+    //  XCe  ...Eccentricity dependet parameter
+    //  YCe  ...Eccentricity dependet parameter
+    //  YHe  ...Eccentricity dependet parameter
+  
+
+    double ei(0.0), Le(0.0);  
+    bool   result = 
+          spheroidGeometry(radius, aspectRatio,  //inputs
+                                 ai,     bi,              //outputs
+                                 ei,    Le              //outputs
+                              );
+    if(!result)
+        return false;
+        
+    XAe= 2.6666666666666666666666667
+         *ei*ei*ei
+         /(-2.0*ei+(1.0+ei*ei)*Le); 
+    YAe= 5.333333333333333333333333333
+        *ei*ei*ei
+        /(2.0*ei+(3*ei*ei-1.0)*Le);
+    XCe= 1.333333333333333333333333333
+        *ei*ei*ei
+        *(1.0-ei*ei)
+        /(2.0*ei-(1.0-ei*ei)*Le);
+    YCe= 1.3333333333333333333333
+        *ei*ei*ei
+        *(2.0-ei*ei)
+        /(-2.0*ei+(1.0+ei*ei)*Le);
+    YHe= 1.3333333333333333333333
+        *ei*ei*ei*ei*ei
+        /(-2.0*ei+(1.0+ei*ei)*Le);
+
+	return true; 		    
+
+}
+
+//--------------------------------------------------------------------
+// zeroize a 3x3x3 tensor
+inline void zeroize333(double tensor[3][3][3] )
+{
+    for(int iX=0; iX<3; iX++)
+        for(int iY=0; iY<3; iY++)
+            for(int iZ=0; iZ<3; iZ++)
+                tensor[iX][iY][iZ] = 0.0;
+}
+
+//--------------------------------------------------------------------
+// zeroize a 3x3 tensor
+inline void zeroize33(double tensor[3][3] )
+{
+    for(int iX=0; iX<3; iX++)
+        for(int iY=0; iY<3; iY++)
+                tensor[iX][iY] = 0.0;
+}
+
+//--------------------------------------------------------------------
+// multiply a 3x3x3 tensor with a scalar
+inline void multiply333(double scalar, double tensor[3][3][3] )
+{
+    for(int iX=0; iX<3; iX++)
+        for(int iY=0; iY<3; iY++)
+            for(int iZ=0; iZ<3; iZ++)
+                tensor[iX][iY][iZ] *= scalar;
+}
+//--------------------------------------------------------------------
+// Compute a dot and dyadic product of with a vector
+inline void permutationTensor(double tensor[3][3][3] )
+{
+    zeroize333(tensor);
+    tensor[0][1][2] = 1.0; 
+    tensor[1][2][0] = 1.0; 
+    tensor[2][0][1] = 1.0; 
+    tensor[0][2][1] =-1.0; 
+    tensor[2][1][0] =-1.0; 
+    tensor[1][0][2] =-1.0; 
+}
+
+
+
+//--------------------------------------------------------------------
+// Compute a dot product of the permutation tensor and 
+// then a dyadic product of with a vector
+inline bool permutationDotDyadic( 
+                                 double vector[3], 
+                                 double  tensor[3][3][3]      
+                                )
+{
+    //Generate permutation tensor
+    double permutationT[3][3][3];
+    permutationTensor(permutationT);
+    
+    //Step 1: compute dot prodcut of permutation tensor
+    double permutationDotProd[3][3];
+    zeroize33(permutationDotProd);
+
+    for(int iX=0; iX<3; iX++)
+        for(int iY=0; iY<3; iY++)
+            for(int iZ=0; iZ<3; iZ++)
+                permutationDotProd[iX][iY] += permutationT[iX][iY][iZ]
+                                            * vector[iZ];
+
+    for(int iX=0; iX<3; iX++)
+        for(int iY=0; iY<3; iY++)
+            for(int iZ=0; iZ<3; iZ++)
+                tensor[iX][iY][iZ] = permutationDotProd[iX][iY]
+                                   * vector[iZ]; 
+    return true;
+
+}
+
+//--------------------------------------------------------------------
+// Compute a dot and dyadic product of with a vector
+inline bool doubleDotTensor333Tensor33(double tensor333[3][3][3], 
+                                       double tensor33[3][3],
+                                       double result[3]
+                                      )
+{
+    result[0]=0.0;result[1]=0.0;result[2]=0.0;
+    
+    for(int iX=0; iX<3; iX++)
+        for(int iY=0; iY<3; iY++)
+            for(int iZ=0; iZ<3; iZ++)
+                result[iX] += tensor333[iX][iY][iZ] * tensor33[iY][iZ];
+
+    return true;
+}
+
+
+}; //end of namespace
+
+#endif
diff --git a/src/lagrangian/cfdemParticle/etc/library-liggghts-list.txt b/src/lagrangian/cfdemParticle/etc/library-liggghts-list.txt
deleted file mode 100644
index 90e6eb87..00000000
--- a/src/lagrangian/cfdemParticle/etc/library-liggghts-list.txt
+++ /dev/null
@@ -1,4 +0,0 @@
-#syntax: makefileName/dir
-#note: dir is not a path, just a keyword here
-###############################################
-M2M/dir
diff --git a/src/lagrangian/cfdemParticle/etc/library-list.txt b/src/lagrangian/cfdemParticle/etc/library-list.txt
index 27e798d4..e34c2541 100644
--- a/src/lagrangian/cfdemParticle/etc/library-list.txt
+++ b/src/lagrangian/cfdemParticle/etc/library-list.txt
@@ -1 +1,11 @@
 lagrangian/cfdemParticle/dir
+
+#====================================================='
+#- RADL
+fvOptions/dir
+cylPorousMedia/dir
+
+#====================================================='
+#- other
+finiteVolume/dir
+
diff --git a/src/lagrangian/cfdemParticle/etc/solver-list.txt b/src/lagrangian/cfdemParticle/etc/solver-list.txt
index 925e9976..960eff2f 100644
--- a/src/lagrangian/cfdemParticle/etc/solver-list.txt
+++ b/src/lagrangian/cfdemParticle/etc/solver-list.txt
@@ -1,4 +1,14 @@
-cfdemSolverPisoMS/dir
 cfdemSolverPiso/dir
 cfdemSolverIB/dir
 cfdemSolverPisoScalar/dir
+cfdemSolverPimpleImEx/dir
+cfdemSolverIBInterLubrication/dir
+cfdemSolverIBScalar/dir
+cfdemSolverInterDyM/dir
+cfdemSolverInterDyMPC/dir
+cfdemSolverBubble/dir
+cfdemSolverPisoMS/dir
+cfdemSolverPimpleDyM_22x/dir
+cfdemSolverPimpleDyMMS_22x/dir
+cfdemSolverPimpleDyMScalar_22x/dir
+myPimpleDyMFoam/dir
diff --git a/src/lagrangian/cfdemParticle/etc/tutorial-list.txt b/src/lagrangian/cfdemParticle/etc/tutorial-list.txt
index d893faf1..92456e25 100644
--- a/src/lagrangian/cfdemParticle/etc/tutorial-list.txt
+++ b/src/lagrangian/cfdemParticle/etc/tutorial-list.txt
@@ -7,15 +7,51 @@
 #===================================================================#
 
 cfdemSolverPiso/settlingTestMPI/dir
-
 cfdemSolverPiso/ErgunTestMPI/dir
-
 cfdemSolverPiso/ErgunTestMPI_cgs/dir
-
 cfdemSolverPiso/ErgunTestMPI_restart/dir
-
 cfdemSolverIB/twoSpheresGlowinskiMPI/dir
-
 cfdemSolverPisoScalar/packedBedTemp/dir
 
+#===================================================================#
+# RADL
+cfdemSolverPimpleImEx/settlingTestMPI/dir
+cfdemSolverPimpleImEx/ErgunTestMPI/dir
+#cfdemSolverPimpleImEx/crossFlow/dir
+#cfdemSolverIB/periodicCase/dir
+#cfdemSolverIB/cfdemIBPeriodicCubicalBox_fullyPeriodic/dir
+#cfdemSolverIBInterLubrication/twoCoatedParticlesRelMotion_smallTest/dir
+#cfdemSolverIBScalar/cfdemIBPeriodicCubicalBoxScalar/dir
+
+#===================================================================#
+# NesteJacobs
+#Projects/Neste/cfdemSolverBubble/3pFBreactor/dir
+#Projects/Neste/cfdemSolverInterDyM/3pFBreactor/dir
+
+#===================================================================#
+# not in release:
+
+#cfdemSolverPiso/settlingTestBigParticleMPI/dir
 cfdemSolverPiso/ErgunTestCG/dir
+cfdemSolverPiso/ErgunTestM2M/dir
+#cfdemSolverPiso/HopperEmptying/dir
+
+cfdemSolverPimpleDyM/ErgunTestMPI/dir
+
+#cfdemSolverPisoMS/settlingTestMPI/dir
+#cfdemSolverPisoMS/ErgunTestMPI/dir
+
+#cfdemSolverInterDyM/twoPhaseSettlingTest/dir
+#cfdemSolverInterDyM/ErgunTestMPI/dir
+#cfdemSolverInterDyM/granularPiston/dir
+#cfdemSolverInterDyM/sugarNcoffee/dir
+
+#cfdemSolverBubble/ErgunTestMPI_pureLiquid/dir
+
+#- these examples are already designed for 2.3.x
+#cfdemSolverInterDyMPC/sugarNcoffee/dir
+#cfdemSolverInterDyMPC/granularPiston/dir
+#cfdemSolverInterDyMPC/meltingPot/dir
+
+
+
diff --git a/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/library/Makefile.fedora_fpic b/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/Makefile.fedora_fpic
similarity index 95%
rename from src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/library/Makefile.fedora_fpic
rename to src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/Makefile.fedora_fpic
index cbb77dca..146e483b 100644
--- a/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/library/Makefile.fedora_fpic
+++ b/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/Makefile.fedora_fpic
@@ -16,7 +16,7 @@ SIZE =		size
 
 # Files
 
-LIB = 	libcouple.a
+LIB = 	libcoupleMany2Many.a
 SRC =	$(wildcard *.cpp)
 INC =	$(wildcard *.h)
 OBJ = 	$(SRC:.cpp=.o)
diff --git a/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/library/Makefile.g++ b/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/Makefile.g++
similarity index 100%
rename from src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/library/Makefile.g++
rename to src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/Makefile.g++
diff --git a/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/library/README b/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/README
similarity index 100%
rename from src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/library/README
rename to src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/README
diff --git a/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/library/error.cpp b/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/error.cpp
similarity index 100%
rename from src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/library/error.cpp
rename to src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/error.cpp
diff --git a/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/library/error.d b/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/error.d
similarity index 100%
rename from src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/library/error.d
rename to src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/error.d
diff --git a/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/library/error.h b/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/error.h
similarity index 100%
rename from src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/library/error.h
rename to src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/error.h
diff --git a/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/library/files.cpp b/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/files.cpp
similarity index 100%
rename from src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/library/files.cpp
rename to src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/files.cpp
diff --git a/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/library/files.d b/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/files.d
similarity index 100%
rename from src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/library/files.d
rename to src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/files.d
diff --git a/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/library/files.h b/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/files.h
similarity index 100%
rename from src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/library/files.h
rename to src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/files.h
diff --git a/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/library/irregular.cpp b/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/irregular.cpp
similarity index 100%
rename from src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/library/irregular.cpp
rename to src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/irregular.cpp
diff --git a/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/library/irregular.d b/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/irregular.d
similarity index 100%
rename from src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/library/irregular.d
rename to src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/irregular.d
diff --git a/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/library/irregular.h b/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/irregular.h
similarity index 100%
rename from src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/library/irregular.h
rename to src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/irregular.h
diff --git a/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/library/lammps_data_write.cpp b/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/lammps_data_write.cpp
similarity index 100%
rename from src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/library/lammps_data_write.cpp
rename to src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/lammps_data_write.cpp
diff --git a/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/library/lammps_data_write.d b/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/lammps_data_write.d
similarity index 100%
rename from src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/library/lammps_data_write.d
rename to src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/lammps_data_write.d
diff --git a/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/library/lammps_data_write.h b/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/lammps_data_write.h
similarity index 100%
rename from src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/library/lammps_data_write.h
rename to src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/lammps_data_write.h
diff --git a/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/libcouple.a b/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/libcouple.a
new file mode 100644
index 00000000..37d6722a
Binary files /dev/null and b/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/libcouple.a differ
diff --git a/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/library/libcouple.a b/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/libcoupleMany2Many.a
similarity index 92%
rename from src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/library/libcouple.a
rename to src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/libcoupleMany2Many.a
index 545482fb..fc190a46 100644
Binary files a/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/library/libcouple.a and b/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/libcoupleMany2Many.a differ
diff --git a/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/library/many2many.cpp b/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/many2many.cpp
similarity index 100%
rename from src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/library/many2many.cpp
rename to src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/many2many.cpp
diff --git a/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/library/many2many.d b/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/many2many.d
similarity index 100%
rename from src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/library/many2many.d
rename to src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/many2many.d
diff --git a/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/library/many2many.h b/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/many2many.h
similarity index 100%
rename from src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/library/many2many.h
rename to src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/many2many.h
diff --git a/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/library/many2one.cpp b/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/many2one.cpp
similarity index 100%
rename from src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/library/many2one.cpp
rename to src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/many2one.cpp
diff --git a/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/library/many2one.d b/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/many2one.d
similarity index 100%
rename from src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/library/many2one.d
rename to src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/many2one.d
diff --git a/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/library/many2one.h b/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/many2one.h
similarity index 100%
rename from src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/library/many2one.h
rename to src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/many2one.h
diff --git a/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/library/memory.cpp b/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/memory.cpp
similarity index 100%
rename from src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/library/memory.cpp
rename to src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/memory.cpp
diff --git a/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/library/memory.d b/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/memory.d
similarity index 100%
rename from src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/library/memory.d
rename to src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/memory.d
diff --git a/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/library/memory.h b/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/memory.h
similarity index 100%
rename from src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/library/memory.h
rename to src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/memory.h
diff --git a/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/library/one2many.cpp b/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/one2many.cpp
similarity index 100%
rename from src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/library/one2many.cpp
rename to src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/one2many.cpp
diff --git a/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/library/one2many.d b/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/one2many.d
similarity index 100%
rename from src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/library/one2many.d
rename to src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/one2many.d
diff --git a/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/library/one2many.h b/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/one2many.h
similarity index 100%
rename from src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/library/one2many.h
rename to src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/one2many.h
diff --git a/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/library/send2one.cpp b/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/send2one.cpp
similarity index 100%
rename from src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/library/send2one.cpp
rename to src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/send2one.cpp
diff --git a/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/library/send2one.d b/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/send2one.d
similarity index 100%
rename from src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/library/send2one.d
rename to src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/send2one.d
diff --git a/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/library/send2one.h b/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/send2one.h
similarity index 100%
rename from src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/library/send2one.h
rename to src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/library/send2one.h
diff --git a/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/twoWayM2M.C b/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/twoWayMany2Many.C
similarity index 70%
rename from src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/twoWayM2M.C
rename to src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/twoWayMany2Many.C
index df671bbe..843245bb 100644
--- a/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/twoWayM2M.C
+++ b/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/twoWayMany2Many.C
@@ -30,7 +30,7 @@ Description
 \*---------------------------------------------------------------------------*/
 
 #include "error.H"
-#include "twoWayM2M.H"
+#include "twoWayMany2Many.H"
 #include "addToRunTimeSelectionTable.H"
 #include "clockModel.H"
 #include "memory.h"
@@ -41,12 +41,12 @@ namespace Foam
 
 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
 
-defineTypeNameAndDebug(twoWayM2M, 0);
+defineTypeNameAndDebug(twoWayMany2Many, 0);
 
 addToRunTimeSelectionTable
 (
     dataExchangeModel,
-    twoWayM2M,
+    twoWayMany2Many,
     dictionary
 );
 
@@ -54,7 +54,7 @@ addToRunTimeSelectionTable
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 
 // Construct from components
-twoWayM2M::twoWayM2M
+twoWayMany2Many::twoWayMany2Many
 (
     const dictionary& dict,
     cfdemCloud& sm
@@ -82,35 +82,45 @@ twoWayM2M::twoWayM2M
     MPI_Comm_split(MPI_COMM_WORLD,liggghts,0,&comm_liggghts);
 
     // open LIGGGHTS input script
-    char *liggghtsPathChar = new char[256];
-    int n = 0;
+    FILE *fp=NULL;
     if (me == 0)
     {
       // read path from dictionary
       const fileName liggghtsPath(propsDict_.lookup("liggghtsPath"));
-      strcpy(liggghtsPathChar, liggghtsPath.c_str());
-      n = strlen(liggghtsPathChar) + 1;
+      char * liggghtsPathChar = (char*)liggghtsPath.c_str();
 
       Info<<"Executing input script '"<< liggghtsPath.c_str() <<"'"<<endl;
+
+      fp = fopen(liggghtsPathChar,"r");
+
+      if (fp == NULL) {
+        printf("ERROR: Could not open LIGGGHTS input script\n");
+        MPI_Abort(MPI_COMM_WORLD,1);
+      }
     }
 
     if (liggghts == 1) lmp = new LAMMPS_NS::LAMMPS(0,NULL,comm_liggghts);
 
-    MPI_Bcast(&n,1,MPI_INT,0,MPI_COMM_WORLD);
-    if (n > 0) {
-        MPI_Bcast(liggghtsPathChar,n,MPI_CHAR,0,MPI_COMM_WORLD);
-        if (liggghts == 1) lmp->input->file(liggghtsPathChar);
+    int n;
+    char line[1024];
+    while (1) {
+      if (me == 0) {
+        if (fgets(line,1024,fp) == NULL) n = 0;
+        else n = strlen(line) + 1;
+        if (n == 0) fclose(fp);
+      }
+      MPI_Bcast(&n,1,MPI_INT,0,MPI_COMM_WORLD);
+      if (n == 0) break;
+      MPI_Bcast(line,n,MPI_CHAR,0,MPI_COMM_WORLD);
+      if (liggghts == 1) lmp->input->one(line);
     }
 
-    delete [] liggghtsPathChar;
-
     // get DEM time step size
     DEMts_ = lmp->update->dt;
     checkTSsize();
 
     // m2m stuff
     firstRun_=true;
-	safeRun_=false;
     particleLost_=false;
     Npart_=-1;
     lmp2foam_ = NULL;
@@ -133,13 +143,12 @@ twoWayM2M::twoWayM2M
     lost_posAll = NULL;
     cellID_foam_ = NULL;
     pos_foam_ = NULL;
-	if (propsDict_.found("safeRun")) safeRun_=true;
 }
 
 
 // * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
 
-twoWayM2M::~twoWayM2M()
+twoWayMany2Many::~twoWayMany2Many()
 {
     free(id_lammps_);
     free(id_lammpsVec_);
@@ -163,7 +172,7 @@ twoWayM2M::~twoWayM2M()
 
 
 // * * * * * * * * * * * * * * * private Member Functions  * * * * * * * * * * * * * //
-char* twoWayM2M::wordToChar(word& inWord) const
+char* twoWayMany2Many::wordToChar(word& inWord) const
 {
     string HH = string(inWord);
     return const_cast<char*>(HH.c_str());
@@ -171,7 +180,7 @@ char* twoWayM2M::wordToChar(word& inWord) const
 
 
 // * * * * * * * * * * * * * * * public Member Functions  * * * * * * * * * * * * * //
-void twoWayM2M::getData
+void twoWayMany2Many::getData
 (
     word name,
     word type,
@@ -183,32 +192,15 @@ void twoWayM2M::getData
     if ( type == "vector-atom" && name != "x")
     {
         double **tmp_ = (double **) lammps_extract_atom(lmp,charName);
-        //for (int i = 0; i < nlocal_lammps_; i++)
-        //    for(int j=0;j<3; j++)
-        //       Pout << couplingStep_ << "st tmp_[" << i << "][j]=" << tmp_[i][j] << "  - name="<< name <<endl;
-
         lmp2foam_vec_->exchange(tmp_ ? tmp_[0] : NULL, field[0]);
-
-        //for (int i = 0; i < nlocal_foam_; i++)
-        //    for(int j=0;j<3; j++)
-        //       Pout << couplingStep_ << "st field[" << i << "][j]=" << field[i][j] << "  - name="<< name <<endl;
-
     }else if (name != "x")
     {
         double *tmp_ = (double *) lammps_extract_atom(lmp,charName);
-
-        //for (int i = 0; i < nlocal_lammps_; i++)
-        //    Pout << couplingStep_ << "st tmp_[" << i << "]=" << tmp_[i] << "  -  name=" << name <<endl;
-
         lmp2foam_->exchange(tmp_, field[0]);
-
-        //for (int i = 0; i < nlocal_foam_; i++)
-        //    Pout << couplingStep_ << "st field[0][" << i << "]=" << field[0][i] << "  -  name=" << name <<endl;
     }
-    //Info << "getData done for :" << name << endl;
 }
 
-void twoWayM2M::getData
+void twoWayMany2Many::getData
 (
     word name,
     word type,
@@ -221,7 +213,7 @@ void twoWayM2M::getData
     data_liggghts_to_of(charName,charType, lmp, (void*&) field,"int");
 }
 
-void twoWayM2M::giveData
+void twoWayMany2Many::giveData
 (
     word name,
     word type,
@@ -230,7 +222,7 @@ void twoWayM2M::giveData
 ) const
 {
     char* charName = wordToChar(name);
-    if ( type == "vector-atom" )
+    if ( type == "vector-atom")
     {
         double **tmp_=NULL;
         LAMMPS_NS::Fix *fix = NULL;
@@ -238,17 +230,9 @@ void twoWayM2M::giveData
         if(fix)
             tmp_ = (double **) static_cast<LAMMPS_NS::FixPropertyAtom*>(fix)->array_atom;
         else
-            FatalError << "coupling fix not found!"<< abort(FatalError);
+            Warning << "coupling fix not found!"<<endl;
 
         foam2lmp_vec_->exchange(field[0],tmp_ ? tmp_[0] : NULL);
-
-        //==================
-        //for(int index = 0;index <  nlocal_lammps_; ++index){
-        //    vector forceField(field[index][0],field[index][1],field[index][2]); 
-        //    vector tmpField(tmp_[index][0],tmp_[index][1],tmp_[index][2]); 
-        //    Pout << "particle=" << index << " ,forceField=" << forceField<< " ,tmpField=" << tmpField <<  endl;
-        //}
-        //==================
     }else if( type == "scalar-atom" )
     {
         Warning << "LIGGGHTS not ready for use of impleDEM and Many2Many" << endl;
@@ -264,22 +248,14 @@ void twoWayM2M::giveData
             allocateArray(tmp_,0,nlocal_lammps_);
 
         foam2lmp_->exchange(field[0],tmp_ ? tmp_ : NULL);
-
-        //==================
-        //for(int index = 0;index <  nlocal_lammps_; ++index){
-        //    scalar forceField(field[index][0]); 
-        //    scalar tmpField(tmp_[index]); 
-        //    Pout << "particle=" << index << " ,forceField=" << forceField<< " ,tmpField=" << tmpField <<  endl;
-        //}
-        //==================
     }else{
-        FatalError << "twoWayM2M::giveData requested type "<< type <<" not implemented! \n"<< abort(FatalError);
+        FatalError << "twoWayMany2Many::giveData requested type not implemented! \n"<< abort(FatalError);
     }
 }
 
 //============
 // double **
-void Foam::twoWayM2M::allocateArray
+void Foam::twoWayMany2Many::allocateArray
 (
     double**& array,
     double initVal,
@@ -294,7 +270,7 @@ void Foam::twoWayM2M::allocateArray
             array[i][j] = initVal;
 }
 
-void Foam::twoWayM2M::allocateArray
+void Foam::twoWayMany2Many::allocateArray
 (
     double**& array,
     double initVal,
@@ -308,14 +284,13 @@ void Foam::twoWayM2M::allocateArray
         for (int j = 0; j < width; j++)
             array[i][j] = initVal;
 }
-void inline Foam::twoWayM2M::destroy(double** array,int len) const
+void inline Foam::twoWayMany2Many::destroy(double** array,int len) const
 {
-    // destroy array[i] first?
     lmp->memory->destroy(array);
 }
 //============
 // int **
-void Foam::twoWayM2M::allocateArray
+void Foam::twoWayMany2Many::allocateArray
 (
     int**& array,
     int initVal,
@@ -330,7 +305,7 @@ void Foam::twoWayM2M::allocateArray
             array[i][j] = initVal;
 }
 
-void Foam::twoWayM2M::allocateArray
+void Foam::twoWayMany2Many::allocateArray
 (
     int**& array,
     int initVal,
@@ -344,41 +319,40 @@ void Foam::twoWayM2M::allocateArray
         for (int j = 0; j < width; j++)
             array[i][j] = initVal;
 }
-void inline Foam::twoWayM2M::destroy(int** array,int len) const
+void inline Foam::twoWayMany2Many::destroy(int** array,int len) const
 {
-    // destroy array[i] first?
     lmp->memory->destroy(array);
 }
 //============
 // double *
-void Foam::twoWayM2M::allocateArray(double*& array, double initVal, int length) const
+void Foam::twoWayMany2Many::allocateArray(double*& array, double initVal, int length) const
 {
     int len = max(length,1);
     lmp->memory->grow(array, len, "m2m:dbl*");
     for (int i = 0; i < len; i++)
         array[i] = initVal;
 }
-void inline Foam::twoWayM2M::destroy(double* array) const
+void inline Foam::twoWayMany2Many::destroy(double* array) const
 {
     lmp->memory->destroy(array);
 }
 //==============
 // int *
-void Foam::twoWayM2M::allocateArray(int*& array, int initVal, int length) const
+void Foam::twoWayMany2Many::allocateArray(int*& array, int initVal, int length) const
 {
     int len = max(length,1);
     lmp->memory->grow(array, len, "m2m:int*");
     for (int i = 0; i < len; i++)
         array[i] = initVal;
 }
-void inline Foam::twoWayM2M::destroy(int* array) const
+void inline Foam::twoWayMany2Many::destroy(int* array) const
 {
     lmp->memory->destroy(array);
 }
 //==============
 
 
-bool Foam::twoWayM2M::couple() const
+bool Foam::twoWayMany2Many::couple() const
 {
     bool coupleNow = false;
     if (doCoupleNow())
@@ -406,59 +380,49 @@ bool Foam::twoWayM2M::couple() const
             Info<<"LIGGGHTS finished"<<endl;
         }
 
-        // give nr of particles to cloud
         double newNpart = liggghts_get_maxtag(lmp);
         setNumberOfParticles(newNpart);
 
-        // check if nr changed
-        // this is not bullet proof - fails if same # particles leave & enter
         if (Npart_ != newNpart)
         {
             Npart_ = newNpart;
             firstRun_ = true;
         }
 
-        // m2m stuff
         particleCloud_.clockM().start(4,"CoupleSyncIDs()");
         syncIDs();
         firstRun_=false;
         particleCloud_.clockM().stop("CoupleSyncIDs()");
 
-        // give nr of particles to cloud
         setNumberOfParticles(nlocal_foam_);
 
         // re-allocate arrays of cloud      
         particleCloud_.reAllocArrays();
 
-        // give existing position and cellID data to cloud
         setPositions(nlocal_foam_,pos_foam_);
         setCellIDs(nlocal_foam_,cellID_foam_);
-        //free(cellID_foam_); // cannot free here?
-        //free(pos_foam_); // cannot free here?
 
-        //MPI_Barrier(MPI_COMM_WORLD);
-        Info <<"Foam::twoWayM2M::couple() done." << endl;
+        Info <<"Foam::twoWayMany2Many::couple() done." << endl;
     }
     return coupleNow;
 }
 
-int Foam::twoWayM2M::getNumberOfParticles() const
+int Foam::twoWayMany2Many::getNumberOfParticles() const
 {
     return liggghts_get_maxtag(lmp);
 }
 
-int Foam::twoWayM2M::getNumberOfClumps() const
+int Foam::twoWayMany2Many::getNumberOfClumps() const
 {
-    Warning << "Foam::twoWayM2M::getNumberOfClumps() - changes necessary here" << endl;
+    Warning << "Foam::twoWayMany2Many::getNumberOfClumps() - changes necessary here" << endl;
     //return liggghts_get_maxtag_ms(lmp);
     return 1;
 }
 
-void Foam::twoWayM2M::syncIDs() const
+void Foam::twoWayMany2Many::syncIDs() const
 {
     particleCloud_.clockM().start(5,"recv_DEM_ids");
 
-    // update communication 
     delete lmp2foam_;
     delete lmp2foam_vec_;
     delete foam2lmp_vec_;
@@ -468,7 +432,6 @@ void Foam::twoWayM2M::syncIDs() const
     foam2lmp_vec_ = new Many2Many(MPI_COMM_WORLD);
     foam2lmp_ = new Many2Many(MPI_COMM_WORLD);
 
-    // get data from lammps
     nlocal_lammps_ = *((int *) lammps_extract_global(lmp,"nlocal"));
 
     int*  id_lammpsSync=NULL;
@@ -478,13 +441,10 @@ void Foam::twoWayM2M::syncIDs() const
 
     if(firstRun_ || particleLost_)
     {
-        // get access to id array
         id_lammps_ = NULL;
         id_lammps_ = (int *) lammps_extract_atom(lmp,"id");
 
-        // genereate vector IDs
         allocateArray(id_lammpsVec_,0,nlocal_lammps_*3);
-        id_lammpsSync = NULL;
         allocateArray(id_lammpsSync,0,nlocal_lammps_);
         for (int i = 0; i < nlocal_lammps_; i++)
         {
@@ -501,27 +461,21 @@ void Foam::twoWayM2M::syncIDs() const
     }
     else
     {
-        // re-arrange data using map
-       
-        // get access to id array
         id_lammps_ = (int *) lammps_extract_atom(lmp,"id");
         allocateArray(id_lammpsSync,0,nlocal_lammps_);
         for (int i = 0; i < nlocal_lammps_; i++)
             id_lammpsSync[i] = id_lammps_[i];
 
-        // generate vector IDs
         allocateArray(id_lammpsVec_,0,nlocal_lammps_*3);
         for (int i = 0; i < nlocal_lammps_; i++)
             for (int j=0;j<3;j++)
                 id_lammpsVec_[i*3+j] = id_lammpsSync[i]*3+j;
 
-        // make setup of m2m
         lmp2foam_->setup(nlocal_lammps_,id_lammpsSync,nlocal_foam_,id_foam_);
         lmp2foam_vec_->setup(nlocal_lammps_*3,id_lammpsVec_,nlocal_foam_*3,id_foamVec_);
         foam2lmp_vec_->setup(nlocal_foam_*3,id_foamVec_,nlocal_lammps_*3,id_lammpsVec_);
         foam2lmp_->setup(nlocal_foam_,id_foam_,nlocal_lammps_,id_lammpsSync);
 
-        // map data according to last TS
         id_lammps_=NULL;
         allocateArray(id_lammps_,-1.,nlocal_foam_);
 		id_lammps_alloc_flag=true;
@@ -531,7 +485,6 @@ void Foam::twoWayM2M::syncIDs() const
         for(int i=0;i<nlocal_foam_;i++)
             id_lammps_[i]=tmpI_[i];
 
-        // get access to "x"
         pos_lammpsSync = (double **) lammps_extract_atom(lmp,"x");
 
         allocateArray(tmp_,-1.,nlocal_foam_*3);
@@ -543,8 +496,6 @@ void Foam::twoWayM2M::syncIDs() const
             for(int j=0;j<3;j++)
                 pos_lammps_[i][j]=tmp_[i*3+j];
 
-        //for (int i = 0; i < nlocal_foam_; i++)
-        //    Pout << couplingStep_ << "st pos exchanged:" <<"=" << pos_lammps_[i][0]<<","<<pos_lammps_[i][1]<<","<<pos_lammps_[i][2] <<endl;
     }
     particleCloud_.clockM().stop("recv_DEM_ids");
 
@@ -553,46 +504,8 @@ void Foam::twoWayM2M::syncIDs() const
 	id_lammps_alloc_flag=true;
     particleCloud_.clockM().stop("locateParticle()");
 
-//MPI_Barrier(MPI_COMM_WORLD);
-//Pout << couplingStep_ << "st == syncIDs  " << endl;
-//if(couplingStep_==30){
-//FatalError<<"stop!!!"<< abort(FatalError);
-//}
-
-        /*//=======================================================================
-        // output
-        Info << "LAMMPS " << endl;
-        for (int i = 0; i < nlocal_lammps_; i++)
-        {
-            if(firstRun_ || particleLost_)
-            {
-                Pout << couplingStep_ << "st id_lammps_[" << i << "]=" << id_lammps_[i] << "  -  "<<endl;
-            }else{
-                Pout << couplingStep_ << "st id_lammpsSync[" << i << "]=" << id_lammpsSync[i] << "  -  "<<endl;
-            }
-        }
-        for (int i = 0; i < nlocal_lammps_*3; i++)
-        {
-            Pout << couplingStep_ << "st id_lammpsVec_[" << i << "]=" << id_lammpsVec_[i] << "  -  "<<endl;
-        }
-        Info << "FOAM "<< endl;
-        for (int i = 0; i < nlocal_foam_; i++)
-        {
-            Pout << couplingStep_ << "st id_foam_[" << i << "]=" << id_foam_[i] << "  -  "<<endl;
-        }
-        for (int i = 0; i < nlocal_foam_*3; i++)
-        {
-            Pout << couplingStep_ << "st id_foamVec_[" << i << "]=" << id_foamVec_[i] << "  -  "<<endl;
-        }
-        MPI_Barrier(MPI_COMM_WORLD);
-        Pout << couplingStep_ << "st nlocal_lammps_=" << nlocal_lammps_ << endl;
-        Pout << couplingStep_ << "st nlocal_foam_=" << nlocal_foam_ << endl;*/
-        //=======================================================================
-
-    // correct mapping
     particleCloud_.clockM().start(11,"setup_Comm");
 
-    // update communication
     delete lmp2foam_;
     delete lmp2foam_vec_;
     delete foam2lmp_vec_;
@@ -624,17 +537,14 @@ void Foam::twoWayM2M::syncIDs() const
 
     particleCloud_.clockM().stop("setup_Comm");
 
-//MPI_Barrier(MPI_COMM_WORLD);
-//Info << "update communication - done." << endl;
 }
 
-void Foam::twoWayM2M::locateParticle(int* id_lammpsSync, bool id_lammps_alloc_flag) const
+void Foam::twoWayMany2Many::locateParticle(int* id_lammpsSync, bool id_lammps_alloc_flag) const
 {
 #if defined(version21)
 
     int nop = particleCloud_.numberOfParticles();
 
-    // realloc array of lost particles     // these arrays will be too long, but we do not know their length a priori???
     allocateArray(id_foamLost_,0,nop);
     allocateArray(lost_pos_,0.,nop*3);
     allocateArray(id_foam_,0,nop);
@@ -643,7 +553,6 @@ void Foam::twoWayM2M::locateParticle(int* id_lammpsSync, bool id_lammps_alloc_fl
     if(firstRun_)
 		allocateArray(pos_foam_,0,nop*3);
 
-    // stage 1 - look on proc or send or prepare for all-to-all
     particleCloud_.clockM().start(7,"locate_Stage1");
     int iterate;
     if(firstRun_ || particleLost_) iterate=nlocal_lammps_;
@@ -663,13 +572,10 @@ void Foam::twoWayM2M::locateParticle(int* id_lammpsSync, bool id_lammps_alloc_fl
         cellID = particleCloud_.locateM().findSingleCell(pos,searchCellID);
         point oldPos(pos_foam_[nlocal_foam_*3+0],pos_foam_[nlocal_foam_*3+1],pos_foam_[nlocal_foam_*3+2]);
 
-        // found particle on cfd proc
         if (cellID >= 0)
         {
-            // IDs for scalars
             id_foam_[nlocal_foam_] = id_lammps_[i];
 
-            // IDs for vectors
             for (int j=0;j<3;j++)
             {
                 id_foamVec_[nlocal_foam_*3+j] = id_lammps_[i]*3+j;                
@@ -678,16 +584,12 @@ void Foam::twoWayM2M::locateParticle(int* id_lammpsSync, bool id_lammps_alloc_fl
             cellID_foam_[nlocal_foam_] = cellID;
 
             nlocal_foam_ += 1;
-            //Pout << couplingStep_ << "st stage1 found particle at pos=" << pos << " ,id_lammps_[i]=" << id_lammps_[i] << endl;
-
         }
         else
         {
-            // find out where particle has migrated (must have passed a CFD proc border)
             bool commPart=false;
             point newPos=pos;
             label nearestFace = particleCloud_.locateM().intersection(oldPos,newPos);
-            //Pout << couplingStep_ << "st nearestFace="<< nearestFace << " oldPos="<< oldPos <<" at pos=" << pos << endl;
 
             if (nearestFace >= particleCloud_.mesh().nInternalFaces())
             {
@@ -705,31 +607,24 @@ void Foam::twoWayM2M::locateParticle(int* id_lammpsSync, bool id_lammps_alloc_fl
                     particleTransferID[n].append(id_lammps_[i]);
                     particleTransferPos[n].append(pos);
                     commPart=true;
-					//Pout << couplingStep_ << "st communicating particle " << id_lammps_[i] << ", to proc# " << n << endl;
                 }
             }
             if (!commPart)
             {
-                // prepare for all to all comm
                 id_foamLost_[nlocal_foam_lost_] = id_lammps_[i];
           
                 for (int j=0; j<3; j++)
                     lost_pos_[nlocal_foam_lost_*3+j] = pos[j];
 
                 nlocal_foam_lost_ += 1;
-                //Pout << couplingStep_ << "st cellID="<< cellID << " lost particle id="<< id_lammps_[i] <<" at pos=" << pos << endl;
             }
         }
     }
     particleCloud_.clockM().stop("locate_Stage1");
-
-    // stage 2 - recv particle, locate or all-to-all
     particleCloud_.clockM().start(8,"locate_Stage2");
  
-    // Allocate transfer buffers
     PstreamBuffers pBufs(Pstream::nonBlocking);
 
-    // Stream into send buffers
     forAll(particleTransferID, i)
     {
         if (particleTransferID[i].size())
@@ -743,8 +638,6 @@ void Foam::twoWayM2M::locateParticle(int* id_lammpsSync, bool id_lammps_alloc_fl
         }
     }
 
-    // Set up transfers when in non-blocking mode. Returns sizes (in bytes)
-    // to be sent/received.
     labelListList allNTrans(Pstream::nProcs());
 
     pBufs.finishedSends(allNTrans);
@@ -760,7 +653,6 @@ void Foam::twoWayM2M::locateParticle(int* id_lammpsSync, bool id_lammps_alloc_fl
         }
     }
 
-    // Retrieve from receive buffers
     label neighbProci;
     label nRec;
     forAll(neighbourProcs_, i)
@@ -777,17 +669,13 @@ void Foam::twoWayM2M::locateParticle(int* id_lammpsSync, bool id_lammps_alloc_fl
 
             forAll(recvParticleTransferID,i)
             {
-                //Pout << couplingStep_ << " stage 2 received id="<<recvParticleTransferID[i] << "received pos="<<recvParticleTransferPos[i] << endl;
                 pos = recvParticleTransferPos[i];
                 cellID = particleCloud_.locateM().findSingleCell(pos,searchCellID);
         
-                // found particle on cfd proc
                 if (cellID >= 0)
                 {
-                    // IDs for scalars
                     id_foam_[nlocal_foam_] = recvParticleTransferID[i];
 
-                    // IDs for vectors
                     for (int j=0;j<3;j++)
                     {
                         id_foamVec_[nlocal_foam_*3+j] = recvParticleTransferID[i]*3+j;
@@ -795,15 +683,10 @@ void Foam::twoWayM2M::locateParticle(int* id_lammpsSync, bool id_lammps_alloc_fl
                     }
                     cellID_foam_[nlocal_foam_] = cellID;
 
-                    // mark that ID was finally found
-                    //id_foamLostAll[i]=-1;
-
                     nlocal_foam_ += 1;
-                    //Pout << couplingStep_ << "st stage2 found particle at pos=" << pos << " ,id_foam_[i]=" << id_foam_[i] << endl;
                 }
-                else // might have been comm to wrong proc
+                else
                 {
-                    // prepare for all to all comm
                     id_foamLost_[nlocal_foam_lost_] = recvParticleTransferID[i];
           
                     for (int j=0; j<3; j++)
@@ -817,12 +700,10 @@ void Foam::twoWayM2M::locateParticle(int* id_lammpsSync, bool id_lammps_alloc_fl
     }
     particleCloud_.clockM().stop("locate_Stage2");
 
-    // stage 3 - all-to-all
     particleCloud_.clockM().start(9,"locate_Stage3");
 
-    // check if all-to-all is necessary
     int nlocal_foam_lostAll(-1);
-	if (firstRun_ || safeRun_)
+	if (firstRun_)
     {
         particleCloud_.clockM().start(10,"locate_Stage3_1");
 	    MPI_Allreduce(&nlocal_foam_lost_, &nlocal_foam_lostAll, 1, MPI_INT, MPI_SUM, MPI_COMM_WORLD);
@@ -842,24 +723,19 @@ void Foam::twoWayM2M::locateParticle(int* id_lammpsSync, bool id_lammps_alloc_fl
            delete[] id_foamLostAll;
            id_foamLostAll = NULL;
         }
-        int nlocal_foam_lostAll = LAMMPS_NS::MPI_Allgather_Vector(lost_pos_, nlocal_foam_lost_*3, lost_posAll, MPI_COMM_WORLD)/3; // new[] für lost_posAll!!!
+        int nlocal_foam_lostAll = LAMMPS_NS::MPI_Allgather_Vector(lost_pos_, nlocal_foam_lost_*3, lost_posAll, MPI_COMM_WORLD)/3;
         LAMMPS_NS::MPI_Allgather_Vector(id_foamLost_, nlocal_foam_lost_, id_foamLostAll, MPI_COMM_WORLD);
         Info << couplingStep_ << "st nlocal_foam_lostAll=" << nlocal_foam_lostAll << endl;
 
-        // locate lost particles
         for (int i = 0; i < nlocal_foam_lostAll; i++)
         {
             pos = vector(lost_posAll[i*3+0],lost_posAll[i*3+1],lost_posAll[i*3+2]);
-            //Pout << "stage3 look for particle at pos=" << pos << endl;
             cellID = particleCloud_.locateM().findSingleCell(pos,searchCellID);
         
-            // found particle on cfd proc
             if (cellID >= 0)
             {
-                // IDs for scalars
                 id_foam_[nlocal_foam_] = id_foamLostAll[i];
 
-                // IDs for vectors
                 for (int j=0;j<3;j++)
                 {
                     id_foamVec_[nlocal_foam_*3+j] = id_foamLostAll[i]*3+j;
@@ -867,22 +743,17 @@ void Foam::twoWayM2M::locateParticle(int* id_lammpsSync, bool id_lammps_alloc_fl
                 }
                 cellID_foam_[nlocal_foam_] = cellID;
 
-                // mark that ID was finally found
                 id_foamLostAll[i]=-1;
 
                 nlocal_foam_ += 1;
-                //Pout << "stage3 found particle at pos=" << pos << " ,id="<< id_foamLostAll[i] << endl;
             }
         }
     }
     particleCloud_.clockM().stop("locate_Stage3");
 
-    // check if really all particles were found
-    // safeRun_ needs to switched on, if particles might get lost
     particleLost_=false;
-    if (firstRun_ || safeRun_)
+    if (firstRun_)
     {
-//    particleCloud_.clockM().start(10,"locate_Stage3");
         int* id_foam_nowhere_all;
         Foam::dataExchangeModel::allocateArray(id_foam_nowhere_all,1,nlocal_foam_lostAll);
         MPI_Allreduce(id_foamLostAll, id_foam_nowhere_all, nlocal_foam_lostAll, MPI_INT, MPI_MIN, MPI_COMM_WORLD);
@@ -890,23 +761,16 @@ void Foam::twoWayM2M::locateParticle(int* id_lammpsSync, bool id_lammps_alloc_fl
         int i=0;
         while (i < nlocal_foam_lostAll)
         {
-            // these particles where found nowhere
             if (id_foam_nowhere_all[i] > 0)
             {
-                // a particle was found nowhere
                 particleLost_=true;
 
-                // delete liggghts particles which were not found from comm
                 for (int j=0;j<nlocal_lammps_;j++)
                 {
                     if (id_lammpsSync[j]==id_foam_nowhere_all[i])
                     {
-
-                        //Pout << "lost ID: id_foam_nowhere_all[i]=" << id_foam_nowhere_all[i] << endl;
-                        // re-arrange IDs
                         id_lammpsSync[j] = id_lammpsSync[nlocal_lammps_-1];
 
-                        // re-arrange IDs for vectors
                         for (int k=0;k<3;k++)
                             id_lammpsVec_[j*3+k] = id_lammpsVec_[(nlocal_lammps_-1)*3+k];
               
@@ -920,23 +784,11 @@ void Foam::twoWayM2M::locateParticle(int* id_lammpsSync, bool id_lammps_alloc_fl
 		Foam::dataExchangeModel::destroy(id_foam_nowhere_all);
 		id_foam_nowhere_all=NULL;
 		if (id_lammps_alloc_flag) destroy(id_lammps_);
-
-        // make cpy
         id_lammps_=NULL;
         allocateArray(id_lammps_,-1.,nlocal_lammps_);
         for (int i = 0; i < nlocal_lammps_; i++)
             id_lammps_[i]=id_lammpsSync[i];
-
-//    particleCloud_.clockM().stop("locate_Stage3");
-//Pout << "particleLost_=" << particleLost_ << endl;
     }
-
-//int gaga;
-//MPI_Allreduce(&nlocal_foam_, &gaga, 1, MPI_INT, MPI_SUM, MPI_COMM_WORLD);
-//int gugu;
-//MPI_Allreduce(&nlocal_lammps_, &gugu, 1, MPI_INT, MPI_SUM, MPI_COMM_WORLD);
-//Info << "nlocal_foam_ALL=" << gaga << " ,nlocal_lammps_ALL=" << gugu << endl;
-
 #elif defined(version16ext)
     Info << "M2M does not work with 1.6.x" << endl;
 #endif
diff --git a/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/twoWayM2M.H b/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/twoWayMany2Many.H
similarity index 94%
rename from src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/twoWayM2M.H
rename to src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/twoWayMany2Many.H
index 4553a307..180c3a0b 100644
--- a/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayM2M/twoWayM2M.H
+++ b/src/lagrangian/cfdemParticle/subModels/dataExchangeModel/twoWayMany2Many/twoWayMany2Many.H
@@ -31,15 +31,15 @@ Description
     two way DEM-CFD coupling via many 2 many scheme
 
 Class
-    twoWayM2M
+    twoWayMany2Many
 
 SourceFiles
-    twoWayM2M.C
+    twoWayMany2Many.C
 
 \*---------------------------------------------------------------------------*/
 
-#ifndef twoWayM2M_H
-#define twoWayM2M_H
+#ifndef twoWayMany2Many_H
+#define twoWayMany2Many_H
 
 #include "dataExchangeModel.H"
 #include "liggghtsCommandModel.H"
@@ -49,20 +49,18 @@ SourceFiles
 #include "force.h"
 #include "forceModel.H"
 
-//=================================//
-//LAMMPS/LIGGGHTS
 #include <stdio.h>
 #include <stdlib.h>
 #include <string.h>
 #include "mpi.h"
-#include <lammps.h>         // these are LAMMPS include files
+#include <lammps.h>
 #include <input.h>
 #include <atom.h>
 #include <library.h>
 #include <error.h>
 #include <library_cfd_coupling.h>
 #include <update.h>
-// m2m stuff
+
 #include <many2many.h>
 #include <mpi_liggghts.h>
 #include "locateModel.H"
@@ -71,7 +69,6 @@ SourceFiles
 #include <modify.h>
 #include <fix_property_atom.h>
 #include "meshSearch.H"
-//=================================//
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
@@ -82,7 +79,7 @@ namespace Foam
                            Class noDrag Declaration
 \*---------------------------------------------------------------------------*/
 
-class twoWayM2M
+class twoWayMany2Many
 :
     public dataExchangeModel
 {
@@ -93,7 +90,6 @@ private:
 
     // new vars
     mutable bool firstRun_;
-    bool safeRun_;
     mutable bool particleLost_;
     mutable double Npart_;
     mutable Many2Many * lmp2foam_;
@@ -142,13 +138,13 @@ private:
 public:
 
     //- Runtime type information
-    TypeName("twoWayM2M");
+    TypeName("twoWayMany2Many");
 
 
     // Constructors
 
         //- Construct from components
-        twoWayM2M
+        twoWayMany2Many
         (
             const dictionary& dict,
             cfdemCloud& sm
@@ -156,7 +152,7 @@ public:
 
     // Destructor
 
-        ~twoWayM2M();
+        ~twoWayMany2Many();
 
 
     // Member Functions
diff --git a/src/lagrangian/cfdemParticle/subModels/locateModel/turboEngineSearchM2M/turboEngineSearchM2M.C b/src/lagrangian/cfdemParticle/subModels/locateModel/engineSearchMany2Many/engineSearchMany2Many.C
similarity index 76%
rename from src/lagrangian/cfdemParticle/subModels/locateModel/turboEngineSearchM2M/turboEngineSearchM2M.C
rename to src/lagrangian/cfdemParticle/subModels/locateModel/engineSearchMany2Many/engineSearchMany2Many.C
index 1417d198..0df1959b 100755
--- a/src/lagrangian/cfdemParticle/subModels/locateModel/turboEngineSearchM2M/turboEngineSearchM2M.C
+++ b/src/lagrangian/cfdemParticle/subModels/locateModel/engineSearchMany2Many/engineSearchMany2Many.C
@@ -31,7 +31,7 @@ Description
 
 #include "error.H"
 
-#include "turboEngineSearchM2M.H"
+#include "engineSearchMany2Many.H"
 #include "addToRunTimeSelectionTable.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -41,12 +41,12 @@ namespace Foam
 
 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
 
-defineTypeNameAndDebug(turboEngineSearchM2M, 0);
+defineTypeNameAndDebug(engineSearchMany2Many, 0);
 
 addToRunTimeSelectionTable
 (
     locateModel,
-    turboEngineSearchM2M,
+    engineSearchMany2Many,
     dictionary
 );
 
@@ -54,38 +54,25 @@ addToRunTimeSelectionTable
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 
 // Construct from components
-turboEngineSearchM2M::turboEngineSearchM2M
+engineSearchMany2Many::engineSearchMany2Many
 (
     const dictionary& dict,
     cfdemCloud& sm
 )
 :
-    turboEngineSearch(dict.subDict(typeName + "Props"),sm)
+    engineSearch(dict.subDict(typeName + "Props"),sm)
 {}
 
 
 // * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
 
-turboEngineSearchM2M::~turboEngineSearchM2M()
+engineSearchMany2Many::~engineSearchMany2Many()
 {}
 
 
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
-label turboEngineSearchM2M::findCell
-(
-    double** const& mask,
-    double**& positions,
-    double**& cellIDs,
-    int size
-) const
-{
-    // search should already be done by M2M
-
-    return 1;
-}
-
-label turboEngineSearchM2M::intersection
+label engineSearchMany2Many::intersection
 (
     const point& pStart,
     const point& pEnd
@@ -93,24 +80,17 @@ label turboEngineSearchM2M::intersection
 {
     // find intersection with boundary
     label face = searchEngine_.findNearestBoundaryFace(pEnd);
-    //Pout << "a1- found face=" << face << " with findNearestBoundaryFace" << endl;
 
     // try alternative
     if (face==-1)
     {
-        // try alternative
         face = searchEngine_.intersection(pStart,pEnd).index();
-        //Pout << "a2- found face=" << face << " with intersection" << endl;
 
-        // might have been first search
         if (face==-1 && mag(pStart-point(0,0,0))<SMALL)
         {
             point pStart2 = pEnd+0.0001*(pStart-pEnd)/mag(pStart-pEnd);
             face = searchEngine_.intersection(pStart2,pEnd).index();
-            //Pout << "a2- found face=" << face << " with intersection"
-            //     << " with pStart2="<< pStart2  << endl;
         }
-
     }
     return face;
 }
diff --git a/src/lagrangian/cfdemParticle/subModels/locateModel/turboEngineSearchM2M/turboEngineSearchM2M.H b/src/lagrangian/cfdemParticle/subModels/locateModel/engineSearchMany2Many/engineSearchMany2Many.H
similarity index 84%
rename from src/lagrangian/cfdemParticle/subModels/locateModel/turboEngineSearchM2M/turboEngineSearchM2M.H
rename to src/lagrangian/cfdemParticle/subModels/locateModel/engineSearchMany2Many/engineSearchMany2Many.H
index 98de363c..5fe4ffce 100755
--- a/src/lagrangian/cfdemParticle/subModels/locateModel/turboEngineSearchM2M/turboEngineSearchM2M.H
+++ b/src/lagrangian/cfdemParticle/subModels/locateModel/engineSearchMany2Many/engineSearchMany2Many.H
@@ -29,17 +29,17 @@ Description
     and OpenFOAM(R). Note: this code is not part of OpenFOAM(R) (see DISCLAIMER).
 
 Class
-    turboEngineSearchM2M
+    engineSearchMany2Many
 
 SourceFiles
-    turboEngineSearchM2M.C
+    engineSearchMany2Many.C
 
 \*---------------------------------------------------------------------------*/
 
-#ifndef turboEngineSearchM2M_H
-#define turboEngineSearchM2M_H
+#ifndef engineSearchMany2Many_H
+#define engineSearchMany2Many_H
 
-#include "turboEngineSearch.H"
+#include "engineSearch.H"
 
 #include "meshSearch.H"
 
@@ -52,9 +52,9 @@ namespace Foam
                            Class noDrag Declaration
 \*---------------------------------------------------------------------------*/
 
-class turboEngineSearchM2M
+class engineSearchMany2Many
 :
-    public turboEngineSearch
+    public engineSearch
 {
 
 private:
@@ -63,13 +63,13 @@ private:
 public:
 
     //- Runtime type information
-    TypeName("turboEngineM2M");
+    TypeName("engineSearchMany2Many");
 
 
     // Constructors
 
         //- Construct from components
-        turboEngineSearchM2M
+        engineSearchMany2Many
         (
             const dictionary& dict,
             cfdemCloud& sm
@@ -77,18 +77,10 @@ public:
 
     // Destructor
 
-        ~turboEngineSearchM2M();
+        ~engineSearchMany2Many();
 
 
     // Member Functions
-        label findCell
-        (
-            double** const& mask,
-            double**& positions,
-            double**& cellIDs,
-            int size
-        ) const;
-
         label intersection
         (
                 const point& pStart,
diff --git a/tutorials/cfdemSolverPisoScalar/packedBedTemp/DEM/in.packing b/tutorials/cfdemSolverPisoScalar/packedBedTemp/DEM/in.packing
deleted file mode 100644
index a6649c46..00000000
--- a/tutorials/cfdemSolverPisoScalar/packedBedTemp/DEM/in.packing
+++ /dev/null
@@ -1,73 +0,0 @@
-#Particle packing by insertion and successive growing of particles
-
-atom_style	granular
-atom_modify	map array
-boundary	m m m
-newton		off
-echo		both
-
-communicate	single vel yes
-
-units		si
-processors	1 1 2
-
-region		reg block 0. 0.1 0. 0.1 0. 1.1 units box
-create_box	1 reg
-
-neighbor	0.002 bin
-neigh_modify	delay 0
-
-
-#Material properties required for new pair styles
-
-fix 		m1 all property/global youngsModulus peratomtype 5.e6
-fix 		m2 all property/global poissonsRatio peratomtype 0.45
-fix 		m3 all property/global coefficientRestitution peratomtypepair 1 0.3
-fix 		m4 all property/global coefficientFriction peratomtypepair 1 0.5
-
-#New pair style
-pair_style gran model hertz tangential history #Hertzian without cohesion
-pair_coeff	* *
-
-timestep	0.00001
-
-#walls
-fix xwalls1 all wall/gran model hertz tangential history primitive type 1  xplane 0.
-fix xwalls2 all wall/gran model hertz tangential history primitive type 1  xplane 0.1
-fix ywalls1 all wall/gran model hertz tangential history primitive type 1  yplane 0.
-fix ywalls2 all wall/gran model hertz tangential history primitive type 1  yplane 0.1
-fix zwalls1 all wall/gran model hertz tangential history primitive type 1  zplane 0.
-fix zwalls2 all wall/gran model hertz tangential history primitive type 1  zplane 1.1
-
-fix		gravi all gravity 9.81 vector 0.0 0.0 -1.0
-
-#heat transfer
-fix 		ftco all property/global thermalConductivity peratomtype 5.
-fix 		ftca all property/global thermalCapacity peratomtype 10.
-fix         heattransfer all heat/gran initial_temperature 263. 
-
-#particle distributions and insertion
-region		bc block 0. 0.1 0. 0.1 0. 1.1 units box
-fix		    pts1 all particletemplate/sphere 1 atom_type 1 density constant 2500 radius constant 0.011
-fix		    pdd1 all particledistribution/discrete 1.  1 pts1 1.0
-
-fix		ins all insert/pack seed 100001 distributiontemplate pdd1 vel constant 0. 0. -3. insert_every 10000 overlapcheck yes all_in yes particles_in_region 1005 region bc
-
-#apply nve integration to all particles that are inserted as single particles
-fix		integr all nve/sphere
-
-#output settings, include total thermal energy
-compute		1 all erotate/sphere
-thermo_style	custom step atoms ke c_1 f_heattransfer vol
-thermo		1000
-thermo_modify	lost ignore norm no
-compute_modify	thermo_temp dynamic yes
-
-#insert the first particles
-run		1
-dump		dmp all custom 1000 post/dump.packing id type type x y z ix iy iz vx vy vz fx fy fz omegax omegay omegaz radius f_heattransfer[0] f_heatFlux[0]
-
-run		150000
-
-write_restart 	liggghts.restart
-