diff --git a/src/lagrangian/DSMC/clouds/Templates/DSMCCloud/DSMCCloud.C b/src/lagrangian/DSMC/clouds/Templates/DSMCCloud/DSMCCloud.C
index 250c3f8e740014f1b872c77042587f05b5b90c3d..7e128716da169fe1cadf7431b1d99f2864f287b1 100644
--- a/src/lagrangian/DSMC/clouds/Templates/DSMCCloud/DSMCCloud.C
+++ b/src/lagrangian/DSMC/clouds/Templates/DSMCCloud/DSMCCloud.C
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           |
+    \\  /    A nd           | Copyright (C) 2019 OpenCFD Ltd.
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
                             | Copyright (C) 2011-2017 OpenFOAM Foundation
@@ -50,7 +50,7 @@ void Foam::DSMCCloud<ParcelType>::buildConstProps()
 
     forAll(typeIdList_, i)
     {
-        const word& id(typeIdList_[i]);
+        const word& id = typeIdList_[i];
 
         Info<< "    " << id << endl;
 
@@ -64,14 +64,14 @@ void Foam::DSMCCloud<ParcelType>::buildConstProps()
 template<class ParcelType>
 void Foam::DSMCCloud<ParcelType>::buildCellOccupancy()
 {
-    forAll(cellOccupancy_, cO)
+    for (auto& list : cellOccupancy_)
     {
-        cellOccupancy_[cO].clear();
+        list.clear();
     }
 
-    forAllIter(typename DSMCCloud<ParcelType>, *this, iter)
+    for (ParcelType& p : *this)
     {
-        cellOccupancy_[iter().cell()].append(&iter());
+        cellOccupancy_[p.cell()].append(&p);
     }
 }
 
@@ -400,9 +400,8 @@ void Foam::DSMCCloud<ParcelType>::calculateFields()
     scalarField& iDof = iDof_.primitiveFieldRef();
     vectorField& momentum = momentum_.primitiveFieldRef();
 
-    forAllConstIter(typename DSMCCloud<ParcelType>, *this, iter)
+    for (const ParcelType& p : *this)
     {
-        const ParcelType& p = iter();
         const label celli = p.cell();
 
         rhoN[celli]++;
@@ -692,7 +691,7 @@ Foam::DSMCCloud<ParcelType>::DSMCCloud
     const fvMesh& mesh,
     const IOdictionary& dsmcInitialiseDict
 )
-    :
+:
     Cloud<ParcelType>(mesh, cloudName, false),
     DSMCBaseCloud(),
     cloudName_(cloudName),
@@ -1038,13 +1037,11 @@ void Foam::DSMCCloud<ParcelType>::dumpParticlePositions() const
       + this->db().time().timeName() + ".obj"
     );
 
-    forAllConstIter(typename DSMCCloud<ParcelType>, *this, iter)
+    for (const ParcelType& p : *this)
     {
-        const ParcelType& p = iter();
-
         pObj<< "v " << p.position().x()
-            << " "  << p.position().y()
-            << " "  << p.position().z()
+            << ' '  << p.position().y()
+            << ' '  << p.position().z()
             << nl;
     }
 
diff --git a/src/lagrangian/DSMC/clouds/Templates/DSMCCloud/DSMCCloudI.H b/src/lagrangian/DSMC/clouds/Templates/DSMCCloud/DSMCCloudI.H
index d0e96515af56f2d70ac08133764fcca8929648b2..91182dfd01adb30bf8fcaabc83ac20d14e3b9ba4 100644
--- a/src/lagrangian/DSMC/clouds/Templates/DSMCCloud/DSMCCloudI.H
+++ b/src/lagrangian/DSMC/clouds/Templates/DSMCCloud/DSMCCloudI.H
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           |
+    \\  /    A nd           | Copyright (C) 2019 OpenCFD Ltd.
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
                             | Copyright (C) 2011-2016 OpenFOAM Foundation
@@ -259,10 +259,8 @@ inline Foam::scalar Foam::DSMCCloud<ParcelType>::massInSystem() const
 {
     scalar sysMass = 0.0;
 
-    forAllConstIter(typename DSMCCloud<ParcelType>, *this, iter)
+    for (const ParcelType& p : *this)
     {
-        const ParcelType& p = iter();
-
         const typename ParcelType::constantProperties& cP = constProps
         (
             p.typeId()
@@ -280,10 +278,8 @@ inline Foam::vector Foam::DSMCCloud<ParcelType>::linearMomentumOfSystem() const
 {
     vector linearMomentum(Zero);
 
-    forAllConstIter(typename DSMCCloud<ParcelType>, *this, iter)
+    for (const ParcelType& p : *this)
     {
-        const ParcelType& p = iter();
-
         const typename ParcelType::constantProperties& cP = constProps
         (
             p.typeId()
@@ -302,10 +298,8 @@ Foam::DSMCCloud<ParcelType>::linearKineticEnergyOfSystem() const
 {
     scalar linearKineticEnergy = 0.0;
 
-    forAllConstIter(typename DSMCCloud<ParcelType>, *this, iter)
+    for (const ParcelType& p : *this)
     {
-        const ParcelType& p = iter();
-
         const typename ParcelType::constantProperties& cP = constProps
         (
             p.typeId()
@@ -324,10 +318,8 @@ Foam::DSMCCloud<ParcelType>::internalEnergyOfSystem() const
 {
     scalar internalEnergy = 0.0;
 
-    forAllConstIter(typename DSMCCloud<ParcelType>, *this, iter)
+    for (const ParcelType& p : *this)
     {
-        const ParcelType& p = iter();
-
         internalEnergy += p.Ei();
     }
 
diff --git a/src/lagrangian/DSMC/parcels/Templates/DSMCParcel/DSMCParcelIO.C b/src/lagrangian/DSMC/parcels/Templates/DSMCParcel/DSMCParcelIO.C
index c066219452775bf9be2e3d0364c5ef2e2d951d83..556b0a8f285f87411fae453b53ff280361a88683 100644
--- a/src/lagrangian/DSMC/parcels/Templates/DSMCParcel/DSMCParcelIO.C
+++ b/src/lagrangian/DSMC/parcels/Templates/DSMCParcel/DSMCParcelIO.C
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2016 OpenCFD Ltd.
+    \\  /    A nd           | Copyright (C) 2016-2019 OpenCFD Ltd.
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
                             | Copyright (C) 2011-2017 OpenFOAM Foundation
@@ -94,14 +94,12 @@ void Foam::DSMCParcel<ParcelType>::readFields(Cloud<DSMCParcel<ParcelType>>& c)
     c.checkFieldIOobject(c, typeId);
 
     label i = 0;
-    forAllIter(typename Cloud<DSMCParcel<ParcelType>>, c, iter)
+    for (DSMCParcel<ParcelType>& p : c)
     {
-        DSMCParcel<ParcelType>& p = iter();
-
         p.U_ = U[i];
         p.Ei_ = Ei[i];
         p.typeId_ = typeId[i];
-        i++;
+        ++i;
     }
 }
 
@@ -121,14 +119,12 @@ void Foam::DSMCParcel<ParcelType>::writeFields
     IOField<label> typeId(c.fieldIOobject("typeId", IOobject::NO_READ), np);
 
     label i = 0;
-    forAllConstIter(typename Cloud<DSMCParcel<ParcelType>>, c, iter)
+    for (const DSMCParcel<ParcelType>& p : c)
     {
-        const DSMCParcel<ParcelType>& p = iter();
-
         U[i] = p.U();
         Ei[i] = p.Ei();
         typeId[i] = p.typeId();
-        i++;
+        ++i;
     }
 
     U.write(np > 0);
diff --git a/src/lagrangian/molecularDynamics/molecularMeasurements/distribution/distribution.C b/src/lagrangian/molecularDynamics/molecularMeasurements/distribution/distribution.C
index a8ac9eec4aa296cb0874a5612375c91de47d7f3e..82d25034fe0c77511538ffd0a8150a443175aa7c 100644
--- a/src/lagrangian/molecularDynamics/molecularMeasurements/distribution/distribution.C
+++ b/src/lagrangian/molecularDynamics/molecularMeasurements/distribution/distribution.C
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           |
+    \\  /    A nd           | Copyright (C) 2019 OpenCFD Ltd.
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
                             | Copyright (C) 2011-2016 OpenFOAM Foundation
@@ -87,9 +87,9 @@ Foam::label Foam::distribution::totalEntries() const
 {
     label sumOfEntries = 0;
 
-    forAllConstIter(Map<label>, *this, iter)
+    forAllConstIters(*this, iter)
     {
-        sumOfEntries += iter();
+        sumOfEntries += iter.val();
 
         if (sumOfEntries < 0)
         {
@@ -115,9 +115,9 @@ Foam::scalar Foam::distribution::approxTotalEntries() const
 {
     scalar sumOfEntries = 0;
 
-    forAllConstIter(Map<label>, *this, iter)
+    forAllConstIters(*this, iter)
     {
-        sumOfEntries += scalar(iter());
+        sumOfEntries += scalar(iter.val());
     }
 
     return sumOfEntries;
diff --git a/src/lagrangian/molecularDynamics/molecularMeasurements/distribution/distribution.H b/src/lagrangian/molecularDynamics/molecularMeasurements/distribution/distribution.H
index fc739019c96f4dcf2012ca1b7d166c99afcea00e..2514313d1f67acab4187ba6296471314f3e93e72 100644
--- a/src/lagrangian/molecularDynamics/molecularMeasurements/distribution/distribution.H
+++ b/src/lagrangian/molecularDynamics/molecularMeasurements/distribution/distribution.H
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           |
+    \\  /    A nd           | Copyright (C) 2019 OpenCFD Ltd.
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
                             | Copyright (C) 2011-2016 OpenFOAM Foundation
diff --git a/src/lagrangian/molecularDynamics/molecule/mdTools/calculateAutoCorrelationFunctions.H b/src/lagrangian/molecularDynamics/molecule/mdTools/calculateAutoCorrelationFunctions.H
index 2418e86fdb6d0dce348bbd340ed11d02800674fb..d41bca24c34b06ddace81c86ce2c1bd344944d7b 100644
--- a/src/lagrangian/molecularDynamics/molecule/mdTools/calculateAutoCorrelationFunctions.H
+++ b/src/lagrangian/molecularDynamics/molecule/mdTools/calculateAutoCorrelationFunctions.H
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           |
+    \\  /    A nd           | Copyright (C) 2019 OpenCFD Ltd.
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
                             | Copyright (C) 2011-2016 OpenFOAM Foundation
@@ -31,9 +31,11 @@ if (mesh.time().timeIndex() % vacf.sampleSteps() == 0)
 
     label uV = 0;
 
-    forAllConstIter(IDLList<molecule>, molecules, mol)
+    for (const molecule& mol : molecules)
     {
-        uVals[uV++] = mol().U();
+        uVals[uV] = mol.U();
+
+        ++uV;
     }
 
     vacf.calculateCorrelationFunction(uVals);
@@ -43,19 +45,19 @@ if (mesh.time().timeIndex() % pacf.sampleSteps() == 0)
 {
     vector p = Zero;
 
-    forAllConstIter(IDLList<molecule>, molecules, mol)
+    for (const molecule& mol : molecules)
     {
         p.x() +=
-            mol().mass() * mol().U().y() * mol().U().z()
-          + 0.5*mol().rf().yz();
+            mol.mass() * mol.U().y() * mol.U().z()
+          + 0.5*mol.rf().yz();
 
         p.y() +=
-            mol().mass() * mol().U().z() * mol().U().x()
-          + 0.5*mol().rf().zx();
+            mol.mass() * mol.U().z() * mol.U().x()
+          + 0.5*mol.rf().zx();
 
         p.z() +=
-            mol().mass() * mol().U().x() * mol().U().y()
-          + 0.5*mol().rf().xy();
+            mol.mass() * mol.U().x() * mol.U().y()
+          + 0.5*mol.rf().xy();
     }
 
     pacf.calculateCorrelationFunction(p);
@@ -65,14 +67,14 @@ if (mesh.time().timeIndex() % hfacf.sampleSteps() == 0)
 {
     vector s = Zero;
 
-    forAllConstIter(IDLList<molecule>, molecules, mol)
+    for (const molecule& mol : molecules)
     {
         s +=
         (
-            0.5*mol().mass()*magSqr(mol().U())
-          + mol().potentialEnergy()
-        )*mol().U()
-      + 0.5*(mol().rf() & mol().U());
+            0.5*mol.mass()*magSqr(mol.U())
+          + mol.potentialEnergy()
+        )*mol.U()
+      + 0.5*(mol.rf() & mol.U());
     }
 
     hfacf.calculateCorrelationFunction(s);
diff --git a/src/lagrangian/molecularDynamics/molecule/mdTools/meanMomentumEnergyAndNMols.H b/src/lagrangian/molecularDynamics/molecule/mdTools/meanMomentumEnergyAndNMols.H
index 34aab16ee9707e4672ae212786e236fd294ee393..08a8492f87fd77fecb3ddfbf60e7e36a376bdd6c 100644
--- a/src/lagrangian/molecularDynamics/molecule/mdTools/meanMomentumEnergyAndNMols.H
+++ b/src/lagrangian/molecularDynamics/molecule/mdTools/meanMomentumEnergyAndNMols.H
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           |
+    \\  /    A nd           | Copyright (C) 2019 OpenCFD Ltd.
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
                             | Copyright (C) 2011-2016 OpenFOAM Foundation
@@ -56,15 +56,15 @@ label singleStepNMols = molecules.size();
 label singleStepDOFs = 0;
 
 {
-    forAllConstIter(IDLList<molecule>, molecules, mol)
+    for (const molecule& mol : molecules)
     {
-        const label molId = mol().id();
+        const label molId = mol.id();
 
         scalar molMass(molecules.constProps(molId).mass());
 
         singleStepTotalMass += molMass;
 
-        //singleStepCentreOfMass += mol().position()*molMass;
+        //singleStepCentreOfMass += mol.position()*molMass;
     }
 
     // if (singleStepNMols)
@@ -72,9 +72,9 @@ label singleStepDOFs = 0;
     //     singleStepCentreOfMass /= singleStepTotalMass;
     // }
 
-    forAllConstIter(IDLList<molecule>, molecules, mol)
+    for (const molecule& mol : molecules)
     {
-        const label molId = mol().id();
+        const label molId = mol.id();
 
         const molecule::constantProperties cP(molecules.constProps(molId));
 
@@ -82,16 +82,16 @@ label singleStepDOFs = 0;
 
         const diagTensor& molMoI(cP.momentOfInertia());
 
-        const vector& molV(mol().v());
+        const vector& molV(mol.v());
 
-        const vector molOmega(inv(molMoI) & mol().pi());
+        const vector molOmega(inv(molMoI) & mol.pi());
 
-        vector molPiGlobal = mol().Q() & mol().pi();
+        vector molPiGlobal = mol.Q() & mol.pi();
 
         singleStepTotalLinearMomentum += molV * molMass;
 
         singleStepTotalAngularMomentum += molPiGlobal;
-        //+((mol().position() - singleStepCentreOfMass) ^ (molV * molMass));
+        //+((mol.position() - singleStepCentreOfMass) ^ (molV * molMass));
 
         if (mag(molV) > singleStepMaxVelocityMag)
         {
@@ -102,9 +102,9 @@ label singleStepDOFs = 0;
 
         singleStepTotalAngularKE += 0.5*(molOmega & molMoI & molOmega);
 
-        singleStepTotalPE += mol().potentialEnergy();
+        singleStepTotalPE += mol.potentialEnergy();
 
-        singleStepTotalrDotf += tr(mol().rf());
+        singleStepTotalrDotf += tr(mol.rf());
 
         singleStepDOFs += cP.degreesOfFreedom();
     }
diff --git a/src/lagrangian/molecularDynamics/molecule/molecule/moleculeIO.C b/src/lagrangian/molecularDynamics/molecule/molecule/moleculeIO.C
index d15a321bd979e12923ab6cf70c2498518d487bd5..a210a2c55ebd633d96b92352a965535f13ed439d 100644
--- a/src/lagrangian/molecularDynamics/molecule/molecule/moleculeIO.C
+++ b/src/lagrangian/molecularDynamics/molecule/molecule/moleculeIO.C
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2016 OpenCFD Ltd.
+    \\  /    A nd           | Copyright (C) 2016-2019 OpenCFD Ltd.
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
                             | Copyright (C) 2011-2017 OpenFOAM Foundation
@@ -128,10 +128,8 @@ void Foam::molecule::readFields(Cloud<molecule>& mC)
     mC.checkFieldIOobject(mC, id);
 
     label i = 0;
-    forAllIter(moleculeCloud, mC, iter)
+    for (molecule& mol : mC)
     {
-        molecule& mol = iter();
-
         mol.Q_ = Q[i];
         mol.v_ = v[i];
         mol.a_ = a[i];
@@ -140,7 +138,8 @@ void Foam::molecule::readFields(Cloud<molecule>& mC)
         mol.specialPosition_ = specialPosition[i];
         mol.special_ = special[i];
         mol.id_ = id[i];
-        i++;
+
+        ++i;
     }
 }
 
@@ -197,10 +196,8 @@ void Foam::molecule::writeFields(const Cloud<molecule>& mC)
     );
 
     label i = 0;
-    forAllConstIter(moleculeCloud, mC, iter)
+    for (const molecule& mol : mC)
     {
-        const molecule& mol = iter();
-
         Q[i] = mol.Q_;
         v[i] = mol.v_;
         a[i] = mol.a_;
@@ -217,7 +214,7 @@ void Foam::molecule::writeFields(const Cloud<molecule>& mC)
         orientation2[i] = mol.Q_ & vector(0,1,0);
         orientation3[i] = mol.Q_ & vector(0,0,1);
 
-        i++;
+        ++i;
     }
 
     const bool valid = np > 0;
diff --git a/src/lagrangian/molecularDynamics/molecule/moleculeCloud/moleculeCloud.C b/src/lagrangian/molecularDynamics/molecule/moleculeCloud/moleculeCloud.C
index 8a20490eb725e58549683e38976ba340e9e4716c..e9ae8c838ba99dd28c8fdc51f8d33cdac8f77611 100644
--- a/src/lagrangian/molecularDynamics/molecule/moleculeCloud/moleculeCloud.C
+++ b/src/lagrangian/molecularDynamics/molecule/moleculeCloud/moleculeCloud.C
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           |
+    \\  /    A nd           | Copyright (C) 2019 OpenCFD Ltd.
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
                             | Copyright (C) 2011-2017 OpenFOAM Foundation
@@ -97,32 +97,32 @@ void Foam::moleculeCloud::buildConstProps()
 
 void Foam::moleculeCloud::setSiteSizesAndPositions()
 {
-    forAllIter(moleculeCloud, *this, mol)
+    for (molecule& mol : *this)
     {
-        const molecule::constantProperties& cP = constProps(mol().id());
+        const molecule::constantProperties& cP = constProps(mol.id());
 
-        mol().setSiteSizes(cP.nSites());
+        mol.setSiteSizes(cP.nSites());
 
-        mol().setSitePositions(cP);
+        mol.setSitePositions(cP);
     }
 }
 
 
 void Foam::moleculeCloud::buildCellOccupancy()
 {
-    forAll(cellOccupancy_, cO)
+    for (auto& list : cellOccupancy_)
     {
-        cellOccupancy_[cO].clear();
+        list.clear();
     }
 
-    forAllIter(moleculeCloud, *this, mol)
+    for (molecule& mol : *this)
     {
-        cellOccupancy_[mol().cell()].append(&mol());
+        cellOccupancy_[mol.cell()].append(&mol);
     }
 
-    forAll(cellOccupancy_, cO)
+    for (auto& list : cellOccupancy_)
     {
-        cellOccupancy_[cO].shrink();
+        list.shrink();
     }
 }
 
@@ -191,12 +191,7 @@ void Foam::moleculeCloud::calculatePairForce()
 
             IDLList<molecule>& refMols = referredMols[r];
 
-            forAllIter
-            (
-                IDLList<molecule>,
-                refMols,
-                refMol
-            )
+            for (molecule& refMol : refMols)
             {
                 forAll(realCells, rC)
                 {
@@ -206,7 +201,7 @@ void Foam::moleculeCloud::calculatePairForce()
                     {
                         molI = celli[cellIMols];
 
-                        evaluatePair(*molI, refMol());
+                        evaluatePair(*molI, refMol);
                     }
                 }
             }
@@ -219,24 +214,24 @@ void Foam::moleculeCloud::calculateTetherForce()
 {
     const tetherPotentialList& tetherPot(pot_.tetherPotentials());
 
-    forAllIter(moleculeCloud, *this, mol)
+    for (molecule& mol : *this)
     {
-        if (mol().tethered())
+        if (mol.tethered())
         {
-            vector rIT = mol().position() - mol().specialPosition();
+            vector rIT = mol.position() - mol.specialPosition();
 
-            label idI = mol().id();
+            label idI = mol.id();
 
             scalar massI = constProps(idI).mass();
 
             vector fIT = tetherPot.force(idI, rIT);
 
-            mol().a() += fIT/massI;
+            mol.a() += fIT/massI;
 
-            mol().potentialEnergy() += tetherPot.energy(idI, rIT);
+            mol.potentialEnergy() += tetherPot.energy(idI, rIT);
 
             // What to do here?
-            // mol().rf() += rIT*fIT;
+            // mol.rf() += rIT*fIT;
         }
     }
 }
@@ -244,9 +239,9 @@ void Foam::moleculeCloud::calculateTetherForce()
 
 void Foam::moleculeCloud::calculateExternalForce()
 {
-    forAllIter(moleculeCloud, *this, mol)
+    for (molecule& mol : *this)
     {
-        mol().a() += pot_.gravity();
+        mol.a() += pot_.gravity();
     }
 }
 
@@ -382,14 +377,9 @@ void Foam::moleculeCloud::removeHighEnergyOverlaps()
         {
             IDLList<molecule>& refMols = referredMols[r];
 
-            forAllIter
-            (
-                IDLList<molecule>,
-                refMols,
-                refMol
-            )
+            for (molecule& refMol : refMols)
             {
-                molJ = &refMol();
+                molJ = &refMol;
 
                 const List<label>& realCells = ril[r];
 
@@ -489,10 +479,8 @@ void Foam::moleculeCloud::initialiseMolecules
             << abort(FatalError);
     }
 
-    forAll(cellZones, z)
+    for (const cellZone& zone : cellZones)
     {
-        const cellZone& zone(cellZones[z]);
-
         if (zone.size())
         {
             if (!mdInitialiseDict.found(zone.name()))
@@ -1041,9 +1029,9 @@ Foam::label Foam::moleculeCloud::nSites() const
 {
     label n = 0;
 
-    forAllConstIter(moleculeCloud, *this, mol)
+    for (const molecule& mol : *this)
     {
-        n += constProps(mol().id()).nSites();
+        n += constProps(mol.id()).nSites();
     }
 
     return n;
@@ -1135,13 +1123,13 @@ void Foam::moleculeCloud::calculateForce()
     buildCellOccupancy();
 
     // Set accumulated quantities to zero
-    forAllIter(moleculeCloud, *this, mol)
+    for (molecule& mol : *this)
     {
-        mol().siteForces() = Zero;
+        mol.siteForces() = Zero;
 
-        mol().potentialEnergy() = 0.0;
+        mol.potentialEnergy() = 0.0;
 
-        mol().rf() = Zero;
+        mol.rf() = Zero;
     }
 
     calculatePairForce();
@@ -1172,11 +1160,11 @@ void Foam::moleculeCloud::applyConstraintsAndThermostats
         << "----------------------------------------"
         << endl;
 
-    forAllIter(moleculeCloud, *this, mol)
+    for (molecule& mol : *this)
     {
-        mol().v() *= temperatureCorrectionFactor;
+        mol.v() *= temperatureCorrectionFactor;
 
-        mol().pi() *= temperatureCorrectionFactor;
+        mol.pi() *= temperatureCorrectionFactor;
     }
 }
 
@@ -1187,13 +1175,13 @@ void Foam::moleculeCloud::writeXYZ(const fileName& fName) const
 
     os  << nSites() << nl << "moleculeCloud site positions in angstroms" << nl;
 
-    forAllConstIter(moleculeCloud, *this, mol)
+    for (const molecule& mol : *this)
     {
-        const molecule::constantProperties& cP = constProps(mol().id());
+        const molecule::constantProperties& cP = constProps(mol.id());
 
-        forAll(mol().sitePositions(), i)
+        forAll(mol.sitePositions(), i)
         {
-            const point& sP = mol().sitePositions()[i];
+            const point& sP = mol.sitePositions()[i];
 
             os  << pot_.siteIdList()[cP.siteIds()[i]]
                 << ' ' << sP.x()*1e10
diff --git a/src/lagrangian/molecularDynamics/molecule/moleculeCloud/moleculeCloud.H b/src/lagrangian/molecularDynamics/molecule/moleculeCloud/moleculeCloud.H
index 4543d127d7f83c6cc67de3d1abfbec0192d7f0a0..d5bbf7e7e1ede5592a35fe0d6e2cd4be8d18f058 100644
--- a/src/lagrangian/molecularDynamics/molecule/moleculeCloud/moleculeCloud.H
+++ b/src/lagrangian/molecularDynamics/molecule/moleculeCloud/moleculeCloud.H
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           |
+    \\  /    A nd           | Copyright (C) 2019 OpenCFD Ltd.
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
                             | Copyright (C) 2011-2017 OpenFOAM Foundation
@@ -32,7 +32,6 @@ SourceFiles
     moleculeCloudI.H
     moleculeCloud.C
 
-
 \*---------------------------------------------------------------------------*/
 
 #ifndef moleculeCloud_H
@@ -60,10 +59,7 @@ class moleculeCloud
 :
     public Cloud<molecule>
 {
-
-private:
-
-    // Private data
+    // Private Data
 
         const polyMesh& mesh_;
 
@@ -179,22 +175,21 @@ public:
         );
 
 
-        // Access
+    // Access
 
-            inline const polyMesh& mesh() const;
+        inline const polyMesh& mesh() const;
 
-            inline const potential& pot() const;
+        inline const potential& pot() const;
 
-            inline const List<DynamicList<molecule*>>& cellOccupancy() const;
+        inline const List<DynamicList<molecule*>>& cellOccupancy() const;
 
-            inline const InteractionLists<molecule>& il() const;
+        inline const InteractionLists<molecule>& il() const;
 
-            inline const List<molecule::constantProperties> constProps() const;
+        inline const List<molecule::constantProperties> constProps() const;
 
-            inline const molecule::constantProperties&
-                constProps(label id) const;
+        inline const molecule::constantProperties& constProps(label id) const;
 
-            inline Random& rndGen();
+        inline Random& rndGen();
 
 
     // Member Operators