diff --git a/src/lagrangian/basic/InteractionLists/InteractionLists.C b/src/lagrangian/basic/InteractionLists/InteractionLists.C
index 8f14c8b7f8d857b33563e4555e8f30f9d1cd300e..80d49bb4316c058c6ae22314aa0375091df5afa6 100644
--- a/src/lagrangian/basic/InteractionLists/InteractionLists.C
+++ b/src/lagrangian/basic/InteractionLists/InteractionLists.C
@@ -881,6 +881,33 @@ Foam::InteractionLists<ParticleType>::InteractionLists
     // Store wallFaceIndexAndTransformToDistribute
     wallFaceIndexAndTransformToDistribute_.transfer(wallFaceIAndTToExchange);
 
+    // Determine inverse addressing for referred wallFaces
+
+    rwfilInverse_.setSize(mesh_.nCells());
+
+    // Temporary Dynamic lists for accumulation
+    List<DynamicList<label> > rwfilInverseTemp(rwfilInverse_.size());
+
+    // Loop over all referred wall faces
+    forAll(rwfil_, refWallFaceI)
+    {
+        const labelList& realCells = rwfil_[refWallFaceI];
+
+        // Loop over all real cells in that the referred wall face is
+        // to supply interactions to and record the index of this
+        // referred wall face in the real cells entry in rwfilInverse
+
+        forAll(realCells, realCellI)
+        {
+            rwfilInverseTemp[realCells[realCellI]].append(refWallFaceI);
+        }
+    }
+
+    forAll(rwfilInverse_, cellI)
+    {
+        rwfilInverse_[cellI].transfer(rwfilInverseTemp[cellI]);
+    }
+
     // Refer wall faces to the appropriate processor
     referredWallFaces_.setSize(wallFaceIndexAndTransformToDistribute_.size());
 
@@ -895,18 +922,17 @@ Foam::InteractionLists<ParticleType>::InteractionLists
             globalTransforms_.transformIndex(wfiat)
         );
 
-        Tuple2<face, pointField>& referredWallFace = referredWallFaces_[rwfI];
+        Tuple2<face, pointField>& rwf = referredWallFaces_[rwfI];
 
         const labelList& facePts = mesh_.faces()[wallFaceIndex];
 
-        referredWallFace.first() = face(identity(facePts.size()));
+        rwf.first() = face(identity(facePts.size()));
 
-        referredWallFace.second().setSize(facePts.size());
+        rwf.second().setSize(facePts.size());
 
         forAll(facePts, fPtI)
         {
-            referredWallFace.second()[fPtI] =
-            mesh_.points()[facePts[fPtI]] - transform;
+            rwf.second()[fPtI] = mesh_.points()[facePts[fPtI]] - transform;
         }
     }
 
diff --git a/src/lagrangian/basic/InteractionLists/InteractionLists.H b/src/lagrangian/basic/InteractionLists/InteractionLists.H
index dc99558e44bb44b636129e6101188392bc25e661..7c7e405e13660c3653660e71d75aa32e5a76ed39 100644
--- a/src/lagrangian/basic/InteractionLists/InteractionLists.H
+++ b/src/lagrangian/basic/InteractionLists/InteractionLists.H
@@ -128,6 +128,11 @@ class InteractionLists
         //  wall face interaction list)
         labelListList rwfil_;
 
+        //- Inverse addressing for referred wall faces, specifies
+        //  which referred wall faces interact with the real cells
+        //  indexed in this container.
+        labelListList rwfilInverse_;
+
         //- Which cells are to be sent via the cellMap, and an index
         //  specifying how they should be transformed.
         List<labelPair> cellIndexAndTransformToDistribute_;
@@ -252,6 +257,10 @@ public:
             //- Return access to the referred wall face interaction list
             inline const labelListList& rwfil() const;
 
+            //- Return access to the inverse referred wall face
+            //  interaction list
+            inline const labelListList& rwfilInverse() const;
+
             //- Return access to the cellIndexAndTransformToDistribute list
             inline const List<labelPair>&
             cellIndexAndTransformToDistribute() const;
diff --git a/src/lagrangian/basic/InteractionLists/InteractionListsI.H b/src/lagrangian/basic/InteractionLists/InteractionListsI.H
index 26a49eede30f3ffe65cf7a2e4c2fb2ba4b83b3ac..5175ecdb465cfd64213ce3375f7383cdb63b9196 100644
--- a/src/lagrangian/basic/InteractionLists/InteractionListsI.H
+++ b/src/lagrangian/basic/InteractionLists/InteractionListsI.H
@@ -80,6 +80,14 @@ const Foam::labelListList& Foam::InteractionLists<ParticleType>::ril() const
 }
 
 
+template<class ParticleType>
+const Foam::labelListList&
+Foam::InteractionLists<ParticleType>::rilInverse() const
+{
+    return rilInverse_;
+}
+
+
 template<class ParticleType>
 const Foam::labelListList& Foam::InteractionLists<ParticleType>::rwfil() const
 {
@@ -87,6 +95,14 @@ const Foam::labelListList& Foam::InteractionLists<ParticleType>::rwfil() const
 }
 
 
+template<class ParticleType>
+const Foam::labelListList&
+Foam::InteractionLists<ParticleType>::rwfilInverse() const
+{
+    return rwfilInverse_;
+}
+
+
 template<class ParticleType>
 const Foam::List<Foam::labelPair>&
 Foam::InteractionLists<ParticleType>::cellIndexAndTransformToDistribute() const
@@ -104,14 +120,6 @@ wallFaceIndexAndTransformToDistribute() const
 }
 
 
-template<class ParticleType>
-const Foam::labelListList&
-Foam::InteractionLists<ParticleType>::rilInverse() const
-{
-    return rilInverse_;
-}
-
-
 template<class ParticleType>
 const Foam::List<Foam::Tuple2<Foam::face, Foam::pointField> >&
 Foam::InteractionLists<ParticleType>::referredWallFaces() const
diff --git a/src/lagrangian/intermediate/submodels/Kinematic/CollisionModel/PairCollision/PairCollision.C b/src/lagrangian/intermediate/submodels/Kinematic/CollisionModel/PairCollision/PairCollision.C
index c208a9c3efb7e2db225c20889d392f0077baf563..76ca569acd32c28908e6ce0f53250e62fbafec63 100644
--- a/src/lagrangian/intermediate/submodels/Kinematic/CollisionModel/PairCollision/PairCollision.C
+++ b/src/lagrangian/intermediate/submodels/Kinematic/CollisionModel/PairCollision/PairCollision.C
@@ -175,8 +175,6 @@ void Foam::PairCollision<CloudType>::wallInteraction()
 
     const labelListList directWallFaces = il_.dwfil();
 
-    // const ReferredCellList<typename CloudType::parcelType>& ril = il_.ril();
-
     // Storage for the wall interaction sites
     DynamicList<point> flatSites;
     DynamicList<scalar> flatSiteExclusionDistancesSqr;
@@ -264,74 +262,66 @@ void Foam::PairCollision<CloudType>::wallInteraction()
                 }
             }
 
-            // // referred wallFace interactions
-
-            // // The labels of referred cells in range of this real cell
-            // const labelList& referredCellsInRange =
-            // dil.referredCellsForInteraction()[realCellI];
-
-            // forAll(referredCellsInRange, refCellInRangeI)
-            // {
-            //     const ReferredCell<typename CloudType::parcelType>& refCell =
-            //     ril[referredCellsInRange[refCellInRangeI]];
-
-            //     const labelList& refWallFaces = refCell.wallFaces();
-
-            //     forAll(refWallFaces, refWallFaceI)
-            //     {
-            //         label refFaceI = refWallFaces[refWallFaceI];
-
-            //         pointHit nearest = refCell.faces()[refFaceI].nearestPoint
-            //         (
-            //             pos,
-            //             refCell.points()
-            //         );
-
-            //         if (nearest.distance() < r)
-            //         {
-            //             vector normal = refCell.faceAreas()[refFaceI];
-
-            //             normal /= mag(normal);
-
-            //             const vector& nearPt = nearest.rawPoint();
-
-            //             vector pW = nearPt - pos;
-
-            //             scalar normalAlignment = normal & pW/mag(pW);
-
-            //             if (normalAlignment > cosPhiMinFlatWall)
-            //             {
-            //                 // Guard against a flat interaction being
-            //                 // present on the boundary of two or more
-            //                 // faces, which would create duplicate contact
-            //                 // points. Duplicates are discarded.
-            //                 if
-            //                 (
-            //                     !duplicatePointInList
-            //                     (
-            //                         flatSites,
-            //                         nearPt,
-            //                         sqr(r*flatWallDuplicateExclusion)
-            //                     )
-            //                 )
-            //                 {
-            //                     flatSites.append(nearPt);
-
-            //                     flatSiteExclusionDistancesSqr.append
-            //                     (
-            //                         sqr(r) - sqr(nearest.distance())
-            //                     );
-            //                 }
-            //             }
-            //             else
-            //             {
-            //                 otherSites.append(nearPt);
-
-            //                 otherSiteDistances.append(nearest.distance());
-            //             }
-            //         }
-            //     }
-            // }
+            // referred wallFace interactions
+
+            // The labels of referred wall faces in range of this real cell
+            const labelList& cellRefWallFaces = il_.rwfilInverse()[realCellI];
+
+            forAll(cellRefWallFaces, rWFI)
+            {
+                const Tuple2<face, pointField>& rwf =
+                il_.referredWallFaces()[cellRefWallFaces[rWFI]];
+
+                const face& f = rwf.first();
+
+                const pointField& pts = rwf.second();
+
+                pointHit nearest = f.nearestPoint(pos, pts);
+
+                if (nearest.distance() < r)
+                {
+                    vector normal = f.normal(pts);
+
+                    normal /= mag(normal);
+
+                    const vector& nearPt = nearest.rawPoint();
+
+                    vector pW = nearPt - pos;
+
+                    scalar normalAlignment = normal & pW/mag(pW);
+
+                    if (normalAlignment > cosPhiMinFlatWall)
+                    {
+                        // Guard against a flat interaction being
+                        // present on the boundary of two or more
+                        // faces, which would create duplicate contact
+                        // points. Duplicates are discarded.
+                        if
+                        (
+                            !duplicatePointInList
+                            (
+                                flatSites,
+                                nearPt,
+                                sqr(r*flatWallDuplicateExclusion)
+                            )
+                        )
+                        {
+                            flatSites.append(nearPt);
+
+                            flatSiteExclusionDistancesSqr.append
+                            (
+                                sqr(r) - sqr(nearest.distance())
+                            );
+                        }
+                    }
+                    else
+                    {
+                        otherSites.append(nearPt);
+
+                        otherSiteDistances.append(nearest.distance());
+                    }
+                }
+            }
 
             // All flat interaction sites found, now classify the
             // other sites as being in range of a flat interaction, or