From b97596f625253be23e3e6c3203daa16d7f8e32cf Mon Sep 17 00:00:00 2001
From: Anthony Geay <anthony.geay@edf.fr>
Date: Tue, 28 Aug 2018 16:44:46 +0200
Subject: [PATCH] Addition of seg2b reference coordinates for Gauss Points

---
 .../GaussPoints/InterpKernelGaussCoords.cxx   | 40 +++++++++++++++----
 .../GaussPoints/InterpKernelGaussCoords.hxx   |  6 ++-
 .../MEDCouplingBasicsTest6.py                 | 23 +++++++++++
 3 files changed, 59 insertions(+), 10 deletions(-)

diff --git a/src/INTERP_KERNEL/GaussPoints/InterpKernelGaussCoords.cxx b/src/INTERP_KERNEL/GaussPoints/InterpKernelGaussCoords.cxx
index ef4988b86..1fbdfd9a5 100644
--- a/src/INTERP_KERNEL/GaussPoints/InterpKernelGaussCoords.cxx
+++ b/src/INTERP_KERNEL/GaussPoints/InterpKernelGaussCoords.cxx
@@ -28,7 +28,9 @@
 
 using namespace INTERP_KERNEL;
 
-const double GaussInfo::SEG2_REF[2]={-1.0, 1.0};
+const double GaussInfo::SEG2A_REF[2]={-1.0, 1.0};
+
+const double GaussInfo::SEG2B_REF[2]={0., 1.0};
 
 const double GaussInfo::SEG3_REF[3]={-1.0, 1.0, 0.0};
 
@@ -226,7 +228,7 @@ GaussInfo GaussInfo::convertToLinear() const
         std::vector<double> a(SEG3_REF,SEG3_REF+3);
         if(IsSatisfy(a,_my_reference_coord))
           {
-            std::vector<double> c(SEG2_REF,SEG2_REF+2);
+            std::vector<double> c(SEG2A_REF,SEG2A_REF+2);
             return GaussInfo(NORM_SEG2,_my_gauss_coord,getNbGauss(),c,2);
           }
         throw INTERP_KERNEL::Exception("GaussInfo::convertToLinear : not recognized pattern for SEG3 !");
@@ -454,9 +456,14 @@ void GaussInfo::initLocalInfo()
     case NORM_SEG2:
       _my_local_ref_dim = 1;
       _my_local_nb_ref  = 2;
-      seg2Init();
+      seg2aInit();
       aSatify = isSatisfy();
-      CHECK_MACRO;
+      if(!aSatify)
+        {
+          seg2bInit();
+          aSatify = isSatisfy();
+          CHECK_MACRO;
+        }
       break;
 
     case NORM_SEG3:
@@ -700,14 +707,14 @@ void GaussInfo::point1Init()
 /*!
  * Init Segment 2 Reference coordinates ans Shape function.
  */
-void GaussInfo::seg2Init() 
+void GaussInfo::seg2aInit() 
 {
   LOCAL_COORD_MACRO_BEGIN;
   case 0:
-    coords[0] = SEG2_REF[0];
+    coords[0] = SEG2A_REF[0];
     break;
   case 1:
-    coords[0] = SEG2_REF[1];
+    coords[0] = SEG2A_REF[1];
     break;
   LOCAL_COORD_MACRO_END;
   
@@ -717,6 +724,23 @@ void GaussInfo::seg2Init()
   SHAPE_FUN_MACRO_END;
 }
 
+void GaussInfo::seg2bInit() 
+{
+  LOCAL_COORD_MACRO_BEGIN;
+  case 0:
+    coords[0] = SEG2B_REF[0];
+    break;
+  case 1:
+    coords[0] = SEG2B_REF[1];
+    break;
+  LOCAL_COORD_MACRO_END;
+  
+  SHAPE_FUN_MACRO_BEGIN;
+  funValue[0] = 1.0 - gc[0];
+  funValue[1] = gc[0];
+  SHAPE_FUN_MACRO_END;
+}
+
 /*!
  * Init Segment 3 Reference coordinates ans Shape function.
  */
@@ -1023,7 +1047,7 @@ void GaussInfo::quad4cInit()
 }
 
 /*!
- * This shapefunc map is same as degenerated seg2Init
+ * This shapefunc map is same as degenerated seg2aInit
  */
 void GaussInfo::quad4DegSeg2Init()
 {
diff --git a/src/INTERP_KERNEL/GaussPoints/InterpKernelGaussCoords.hxx b/src/INTERP_KERNEL/GaussPoints/InterpKernelGaussCoords.hxx
index 7f010686c..409cbc5f2 100644
--- a/src/INTERP_KERNEL/GaussPoints/InterpKernelGaussCoords.hxx
+++ b/src/INTERP_KERNEL/GaussPoints/InterpKernelGaussCoords.hxx
@@ -63,7 +63,8 @@ namespace INTERP_KERNEL
     INTERPKERNEL_EXPORT static std::vector<double> NormalizeCoordinatesIfNecessary(NormalizedCellType ct, int inputDim, const std::vector<double>& inputArray);
 
   public:
-    static const double SEG2_REF[2];
+    static const double SEG2A_REF[2];
+    static const double SEG2B_REF[2];
     static const double SEG3_REF[3];
     static const double TRIA3A_REF[6];
     static const double TRIA3B_REF[6];
@@ -99,7 +100,8 @@ namespace INTERP_KERNEL
     void point1Init();
     
     //1D
-    void seg2Init();
+    void seg2aInit();
+    void seg2bInit();
     void seg3Init();
 
     //2D
diff --git a/src/MEDCoupling_Swig/MEDCouplingBasicsTest6.py b/src/MEDCoupling_Swig/MEDCouplingBasicsTest6.py
index 4e92707eb..978cfd221 100644
--- a/src/MEDCoupling_Swig/MEDCouplingBasicsTest6.py
+++ b/src/MEDCoupling_Swig/MEDCouplingBasicsTest6.py
@@ -212,6 +212,8 @@ class MEDCouplingBasicsTest6(unittest.TestCase):
         pass
 
     def testRenumberNodesInConnOpt(self):
+        """ Test of MEDCouplingPointSet.renumberNodesInConn with map as input coming from DataArrayInt.invertArrayN2O2O2NOptimized
+        """
         m=MEDCouplingUMesh("mesh",2)
         m.allocateCells()
         m.insertNextCell(NORM_QUAD4,[10000,10002,10001,10003])
@@ -244,6 +246,27 @@ class MEDCouplingBasicsTest6(unittest.TestCase):
         self.assertTrue(m.getNodalConnectivity().isEqual(DataArrayInt([0,2,1,3])))
         self.assertTrue(m.getNodalConnectivityIndex().isEqual(DataArrayInt([0,4])))
         pass
+
+    def testSeg2bGP(self):
+        """Test of Gauss points on SEG2 using SEG2B style as ref coords
+        """
+        coo=DataArrayDouble([[0.,0.,0.],[1.,1.,1.]])
+        m=MEDCouplingUMesh("mesh",1) ; m.setCoords(coo)
+        m.allocateCells()
+        # the cell description is exactly those described in the description of HEXA27 in MED file 3.0.7 documentation
+        m.insertNextCell(NORM_SEG2,[0,1])
+        refCoo=[0.,1.]
+        weights=[0.8,0.1,0.1]
+        gCoords=[0.2,0.5,0.9]
+        fGauss=MEDCouplingFieldDouble(ON_GAUSS_PT) ; fGauss.setName("fGauss")
+        fGauss.setMesh(m)
+        fGauss.setGaussLocalizationOnType(NORM_SEG2,refCoo,gCoords,weights)
+        arr=DataArrayDouble(fGauss.getNumberOfTuplesExpected()) ; arr.iota()
+        fGauss.setArray(arr)
+        fGauss.checkConsistencyLight()
+        arrOfDisc=fGauss.getLocalizationOfDiscr()
+        self.assertTrue(arrOfDisc.isEqual(DataArrayDouble([0.2,0.2,0.2,0.5,0.5,0.5,0.9,0.9,0.9],3,3),1e-12))
+        pass
     
     pass
 
-- 
2.39.2