From d558cd0c72a5b600c7bd184cd23a593369faa8d0 Mon Sep 17 00:00:00 2001
From: Anthony Geay <anthony.geay@edf.fr>
Date: Thu, 18 Aug 2022 11:08:14 +0200
Subject: [PATCH] Implementation of
 MEDCouplingGaussLocalization.getFunctionValues

---
 .../MEDCouplingGaussLocalization.cxx          | 35 ++++++++++++++++---
 .../MEDCouplingGaussLocalization.hxx          |  1 +
 src/MEDCoupling_Swig/MEDCouplingCommon.i      |  7 ++++
 3 files changed, 39 insertions(+), 4 deletions(-)

diff --git a/src/MEDCoupling/MEDCouplingGaussLocalization.cxx b/src/MEDCoupling/MEDCouplingGaussLocalization.cxx
index 3026d8ac1..6337caa15 100644
--- a/src/MEDCoupling/MEDCouplingGaussLocalization.cxx
+++ b/src/MEDCoupling/MEDCouplingGaussLocalization.cxx
@@ -76,15 +76,20 @@ void MEDCouplingGaussLocalization::checkConsistencyLight() const
 int MEDCouplingGaussLocalization::getDimension() const
 {
   if(_weight.empty())
-    return -1;
+    THROW_IK_EXCEPTION("getDimension : weight is empty !");
   return (int)_gauss_coord.size()/(int)_weight.size();
 }
 
 int MEDCouplingGaussLocalization::getNumberOfPtsInRefCell() const
 {
-  int dim=getDimension();
-  if(dim==0)
-    return -1;
+  if(_gauss_coord.empty())
+  {
+    if(!_weight.empty())
+      THROW_IK_EXCEPTION("getNumberOfPtsInRefCell : gauss_coords are empty whereas weights are not empty !");
+    const INTERP_KERNEL::CellModel& cm = INTERP_KERNEL::CellModel::GetCellModel(_type);
+    return ((int)_ref_coord.size()) / ((int)cm.getDimension());
+  }
+  int dim( getDimension() );
   return (int)_ref_coord.size()/dim;
 }
 
@@ -241,6 +246,28 @@ MCAuto<MEDCouplingUMesh> MEDCouplingGaussLocalization::buildRefCell() const
   return MCAuto<MEDCouplingUMesh>(ret->buildUnstructured());
 }
 
+/*!
+ * This method returns an array containing for each Gauss Points in \a this, function values relative to the points of the
+ * reference cell. Number of components of returned array is equal to the number of points of the reference cell.
+ */
+MCAuto<DataArrayDouble> MEDCouplingGaussLocalization::getFunctionValues() const
+{
+  MCAuto<DataArrayDouble> ret(DataArrayDouble::New());
+  int nbGaussPt(getNumberOfGaussPt()),nbPtsRefCell(getNumberOfPtsInRefCell()),dim(getDimension());
+  ret->alloc(nbGaussPt,nbPtsRefCell);
+  double *retPtr(ret->getPointer());
+  for(int iGaussPt = 0 ; iGaussPt < nbGaussPt ; ++iGaussPt)
+  {
+    std::vector<double> curGaussPt(_gauss_coord.begin()+iGaussPt*dim,_gauss_coord.begin()+(iGaussPt+1)*dim);
+    INTERP_KERNEL::GaussInfo gi(_type,curGaussPt,1,_ref_coord,nbPtsRefCell);
+    gi.initLocalInfo();
+    const double *funcVal( gi.getFunctionValues(0) );
+    std::copy(funcVal,funcVal+nbPtsRefCell,retPtr);
+    retPtr += nbPtsRefCell;
+  }
+  return ret;
+}
+
 /*!
  * This method sets the comp_th component of ptIdInCell_th point coordinate of reference element of type this->_type.
  * @throw if not 0<=ptIdInCell<nbOfNodePerCell or if not 0<=comp<dim
diff --git a/src/MEDCoupling/MEDCouplingGaussLocalization.hxx b/src/MEDCoupling/MEDCouplingGaussLocalization.hxx
index 56e061e93..8ed585953 100644
--- a/src/MEDCoupling/MEDCouplingGaussLocalization.hxx
+++ b/src/MEDCoupling/MEDCouplingGaussLocalization.hxx
@@ -54,6 +54,7 @@ namespace MEDCoupling
     //
     MEDCOUPLING_EXPORT MCAuto<DataArrayDouble> localizePtsInRefCooForEachCell(const DataArrayDouble *ptsInRefCoo, const MEDCouplingUMesh *mesh) const;
     MEDCOUPLING_EXPORT MCAuto<MEDCouplingUMesh> buildRefCell() const;
+    MEDCOUPLING_EXPORT MCAuto<DataArrayDouble> getFunctionValues() const;
     //
     MEDCOUPLING_EXPORT const std::vector<double>& getRefCoords() const { return _ref_coord; }
     MEDCOUPLING_EXPORT double getRefCoord(int ptIdInCell, int comp) const;
diff --git a/src/MEDCoupling_Swig/MEDCouplingCommon.i b/src/MEDCoupling_Swig/MEDCouplingCommon.i
index c760c94d2..4fe5a1e41 100644
--- a/src/MEDCoupling_Swig/MEDCouplingCommon.i
+++ b/src/MEDCoupling_Swig/MEDCouplingCommon.i
@@ -478,6 +478,7 @@ typedef long mcPyPtrType;
 %newobject MEDCoupling::DenseMatrix::__mul__;
 %newobject MEDCoupling::MEDCouplingGaussLocalization::localizePtsInRefCooForEachCell;
 %newobject MEDCoupling::MEDCouplingGaussLocalization::buildRefCell;
+%newobject MEDCoupling::MEDCouplingGaussLocalization::getFunctionValues;
 %newobject MEDCoupling::MEDCouplingSkyLineArray::BuildFromPolyhedronConn;
 %newobject MEDCoupling::MEDCouplingSkyLineArray::getSuperIndexArray;
 %newobject MEDCoupling::MEDCouplingSkyLineArray::getIndexArray;
@@ -1304,6 +1305,12 @@ namespace MEDCoupling
         MCAuto<MEDCouplingUMesh> ret(self->buildRefCell());
         return ret.retn();
       }
+      
+      DataArrayDouble *getFunctionValues() const
+      {
+        MCAuto<DataArrayDouble> ret(self->getFunctionValues());
+        return ret.retn();
+      }
     }
   };
 
-- 
2.39.2