From 67b4d2a8cf293f94bc6d74c9eb6479fd9ef8a116 Mon Sep 17 00:00:00 2001
From: Anthony Geay <anthony.geay@edf.fr>
Date: Wed, 24 May 2017 16:04:49 +0200
Subject: [PATCH] Huge bug in matrix inversion revealed by
 MEDCouplingUMesh.computePlaneEquationOf3DFaces.

---
 src/INTERP_KERNEL/InterpKernelMatrixTools.cxx  |  6 +++---
 src/MEDCoupling_Swig/MEDCouplingBasicsTest5.py | 11 +++++++++++
 2 files changed, 14 insertions(+), 3 deletions(-)

diff --git a/src/INTERP_KERNEL/InterpKernelMatrixTools.cxx b/src/INTERP_KERNEL/InterpKernelMatrixTools.cxx
index 80d5f1d02..cc1c08130 100644
--- a/src/INTERP_KERNEL/InterpKernelMatrixTools.cxx
+++ b/src/INTERP_KERNEL/InterpKernelMatrixTools.cxx
@@ -301,12 +301,12 @@ namespace INTERP_KERNEL
     for(int k=0;k<n-1;k++)
       {
         //  Find L=pivot index.
-        l=idamax(n-k-1,a+k+k*lda,1)+k;
+        l=idamax(n-k,a+k+k*lda,1)+k;
         ipvt[k]=l;
         // Zero pivot implies this column already triangularized.
         if(a[l+k*lda]==0.0)
           {
-            info=k+1;
+            info=k;
             continue;
           }
         //Interchange if necessary.
@@ -323,7 +323,7 @@ namespace INTERP_KERNEL
         for(int j=k+1;j<n;j++)
           {
             t=a[l+j*lda];
-            if(l!=k-1)
+            if(l!=k)
               {
                 a[l+j*lda]=a[k+j*lda];
                 a[k+j*lda]=t;
diff --git a/src/MEDCoupling_Swig/MEDCouplingBasicsTest5.py b/src/MEDCoupling_Swig/MEDCouplingBasicsTest5.py
index 344dd7dd3..0bc7f7898 100644
--- a/src/MEDCoupling_Swig/MEDCouplingBasicsTest5.py
+++ b/src/MEDCoupling_Swig/MEDCouplingBasicsTest5.py
@@ -4664,6 +4664,17 @@ class MEDCouplingBasicsTest5(unittest.TestCase):
         d=DataArrayInt([3,4,3])
         self.assertRaises(InterpKernelException,d.checkUniformAndGuess)# non uniform
         pass
+
+    def testUMComputePlaneEquationOf3DFaces1(self):
+        """ Consequence of an invalid traduction of matrix inversion transposition."""
+        m=MEDCoupling1SGTUMesh("msh",NORM_QUAD4)
+        m.setCoords(DataArrayDouble([(0,0,0),(1,0,0),(2,0,0),(0,2,0),(1,2,0),(2,2,0),(0,4,0),(1,4,0),(2,4,0),(0,0,3),(1,0,3),(2,0,3),(0,2,3),(1,2,3),(2,2,3),(0,4,3),(1,4,3),(2,4,3)]))
+        m.setNodalConnectivity(DataArrayInt([0,1,4,3,9,12,13,10,0,9,10,1,1,10,13,4,4,13,12,3,3,12,9,0,1,2,5,4,10,13,14,11,1,10,11,2,2,11,14,5,5,14,13,4,3,4,7,6,12,15,16,13,4,13,16,7,7,16,15,6,6,15,12,3,4,5,8,7,13,16,17,14,5,14,17,8,8,17,16,7]))
+        m=m.buildUnstructured()
+        ref=DataArrayDouble([(0,0,1,0),(0,0,1,-3),(0,1,0,0),(1,0,0,-1),(0,1,0,-2),(1,0,0,0),(0,0,1,0),(0,0,1,-3),(0,1,0,0),(1,0,0,-2),(0,1,0,-2),(0,0,1,0),(0,0,1,-3),(1,0,0,-1),(0,1,0,-4),(1,0,0,0),(0,0,1,0),(0,0,1,-3),(1,0,0,-2),(0,1,0,-4)])
+        res=m.computePlaneEquationOf3DFaces()
+        self.assertTrue(res.isEqual(ref,1e-12))
+        pass
     
     pass
 
-- 
2.39.2