From 2d899b8551c8f73c9aab32ca98419efc7418b41e Mon Sep 17 00:00:00 2001
From: Hrvoje Jasak <h.jasak@wikki.co.uk>
Date: Tue, 7 Apr 2015 08:53:35 +0100
Subject: [PATCH] Parallel bug fix, processor boundaries comms.  Zeljko Tukovic

---
 .../coupledLduMatrix/coupledLduMatrix.C       | 114 +++++++++++++++---
 1 file changed, 94 insertions(+), 20 deletions(-)

diff --git a/src/coupledMatrix/coupledLduMatrix/coupledLduMatrix.C b/src/coupledMatrix/coupledLduMatrix/coupledLduMatrix.C
index 8a5319e3c..2fddadef1 100644
--- a/src/coupledMatrix/coupledLduMatrix/coupledLduMatrix.C
+++ b/src/coupledMatrix/coupledLduMatrix/coupledLduMatrix.C
@@ -30,9 +30,10 @@ Description
 Author
     Hrvoje Jasak, Wikki Ltd.  All rights reserved
 
-\*----------------------------------------------------------------------------*/
+\*---------------------------------------------------------------------------*/
 
 #include "coupledLduMatrix.H"
+#include "processorLduInterfaceField.H"
 
 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
 
@@ -42,9 +43,6 @@ namespace Foam
 }
 
 
-// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
-
-
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 
 // Construct given size
@@ -202,16 +200,101 @@ void Foam::coupledLduMatrix::initMatrixInterfaces
 {
     const PtrList<lduMatrix>& matrices = *this;
 
+    // Note.  The comms design requires all non-processor interfaces
+    // to be updated first, followed by the update of processor
+    // interfaces.  The reason is that non-processor coupled
+    // interfaces require a complex comms pattern involving more than
+    // pairwise communications.
+    // Under normal circumstances this is achieved naturall, since
+    // processor interfaces come last on the list and other coupled
+    // interfaces execute complex comms at init() level.
+    // For coupled matrices, the update loop needs to be split over
+    // all matrices by hand
+    // Bug fix: Zeljko Tukovic, 7/Apr/2015
+
+    // Init update all non-processor coupled interfaces
     forAll (matrices, rowI)
     {
-        matrices[rowI].initMatrixInterfaces
+        if
         (
-            coupleCoeffs[rowI],
-            interfaces[rowI],
-            x[rowI],
-            result[rowI],
-            cmpt
-        );
+            Pstream::defaultCommsType == Pstream::blocking
+         || Pstream::defaultCommsType == Pstream::nonBlocking
+        )
+        {
+            forAll (interfaces[rowI], interfaceI)
+            {
+                if (interfaces[rowI].set(interfaceI))
+                {
+                    if
+                    (
+                       !isA<processorLduInterfaceField>
+                        (
+                            interfaces[rowI][interfaceI]
+                        )
+                    )
+                    {
+                        interfaces[rowI][interfaceI].initInterfaceMatrixUpdate
+                        (
+                            x[rowI],
+                            result[rowI],
+                            matrices[rowI],
+                            coupleCoeffs[rowI][interfaceI],
+                            cmpt,
+                            Pstream::defaultCommsType,
+                            false
+                        );
+                    }
+                }
+            }
+        }
+        else
+        {
+            matrices[rowI].initMatrixInterfaces
+            (
+                coupleCoeffs[rowI],
+                interfaces[rowI],
+                x[rowI],
+                result[rowI],
+                cmpt
+            );
+        }
+    }
+
+    // Init update for all processor interfaces
+    forAll (matrices, rowI)
+    {
+        if
+        (
+            Pstream::defaultCommsType == Pstream::blocking
+         || Pstream::defaultCommsType == Pstream::nonBlocking
+        )
+        {
+            forAll (interfaces[rowI], interfaceI)
+            {
+                if (interfaces[rowI].set(interfaceI))
+                {
+                    if
+                    (
+                        isA<processorLduInterfaceField>
+                        (
+                            interfaces[rowI][interfaceI]
+                        )
+                    )
+                    {
+                        interfaces[rowI][interfaceI].initInterfaceMatrixUpdate
+                        (
+                            x[rowI],
+                            result[rowI],
+                            matrices[rowI],
+                            coupleCoeffs[rowI][interfaceI],
+                            cmpt,
+                            Pstream::defaultCommsType,
+                            false
+                        );
+                    }
+                }
+            }
+        }
     }
 }
 
@@ -241,13 +324,4 @@ void Foam::coupledLduMatrix::updateMatrixInterfaces
 }
 
 
-// * * * * * * * * * * * * * * * Member Operators  * * * * * * * * * * * * * //
-
-
-// * * * * * * * * * * * * * * * Friend Functions  * * * * * * * * * * * * * //
-
-
-// * * * * * * * * * * * * * * * Friend Operators  * * * * * * * * * * * * * //
-
-
 // ************************************************************************* //