From: Anthony Geay Date: Wed, 13 Nov 2024 07:01:03 +0000 (+0100) Subject: [EDF31315] : Better polyedra management X-Git-Url: http://git.salome-platform.org/gitweb/?a=commitdiff_plain;h=369222ec9b82f4298023ce091c8854e4984c2ddb;p=tools%2Fmedcoupling.git [EDF31315] : Better polyedra management --- diff --git a/src/MEDCoupling/MEDCouplingMemArray.hxx b/src/MEDCoupling/MEDCouplingMemArray.hxx index 058502d2b..5f39a3a82 100755 --- a/src/MEDCoupling/MEDCouplingMemArray.hxx +++ b/src/MEDCoupling/MEDCouplingMemArray.hxx @@ -700,6 +700,8 @@ namespace MEDCoupling const DataArrayType *srcArr, const DataArrayIdType *srcArrIndex); static bool RemoveIdsFromIndexedArrays(const T *idsToRemoveBg, const T *idsToRemoveEnd, DataArrayType *arr, DataArrayIdType *arrIndx, mcIdType offsetForRemoval=0); + static void FromVTKInternalReprOfPolyedra(const DataArrayType *arrIn, const DataArrayIdType *arrIndxIn, + MCAuto &arrOut, MCAuto &arrIndexOut); static DataArrayType *Range(T begin, T end, T step); public: void getTinySerializationIntInformation(std::vector& tinyInfo) const; diff --git a/src/MEDCoupling/MEDCouplingMemArray.txx b/src/MEDCoupling/MEDCouplingMemArray.txx index 46d5ac131..940101fff 100755 --- a/src/MEDCoupling/MEDCouplingMemArray.txx +++ b/src/MEDCoupling/MEDCouplingMemArray.txx @@ -7479,6 +7479,52 @@ struct NotInRange arrOut=arro.retn(); arrIndexOut=arrIo.retn(); } + + /*! + * This method converts from VTK polyhedra nodal connectivity to MED. + * + * \param [in] arrIn arr origin array from which the extraction will be done. + * \param [in] arrIndxIn is the input index array allowing to walk into \b arrIn + * \param [out] arrOut the resulting array + * \param [out] arrIndexOut the index array of the resulting array \b arrOut + */ + template + void DataArrayDiscrete::FromVTKInternalReprOfPolyedra(const DataArrayType *arrIn, const DataArrayIdType *arrIndxIn, + MCAuto &arrOut, MCAuto &arrIndexOut) + { + if(!arrIn || !arrIndxIn) + throw INTERP_KERNEL::Exception("DataArrayInt::FromVTKInternalReprOfPolyedra : input pointer is NULL !"); + arrIn->checkAllocated(); arrIndxIn->checkAllocated(); + arrIn->checkNbOfComps(1,"1st array must have single component"); + arrIndxIn->checkNbOfComps(1,"2nd array must have single component"); + if(arrIndxIn->getNumberOfTuples()<1) + THROW_IK_EXCEPTION("2nd input array must be of size >= 1"); + mcIdType nbCells(arrIndxIn->getNumberOfTuples()-1); + const T *arrInPt(arrIn->begin()); + const mcIdType *arrIndxInPt(arrIndxIn->begin()); + arrIndexOut = DataArrayIdType::New(); arrIndexOut->alloc(arrIndxIn->getNumberOfTuples(),1); + arrOut = DataArrayType::New(); arrOut->alloc(arrIn->getNumberOfTuples() - 2*nbCells,1); + T *arrOutPt(arrOut->getPointer()); + mcIdType *arrOutIdxPt(arrIndexOut->getPointer()); *arrOutIdxPt = 0; + for(auto i = 0 ; i < nbCells ; ++i) + { + T nbFaces = arrInPt[ arrIndxInPt[i] ]; + T *arrOutPtStart(arrOutPt); + const T *facePtr = arrInPt + arrIndxInPt[i] + 1; + for(T iFace = 0 ; iFace < nbFaces ; ++iFace) + { + T nbNodesInFace = *facePtr++; + if(iFace>0) + { + *arrOutPt++ = -1; + } + arrOutPt = std::copy(facePtr,facePtr+nbNodesInFace,arrOutPt); + facePtr += nbNodesInFace; + } + arrOutIdxPt[1] = arrOutIdxPt[0] + std::distance(arrOutPtStart,arrOutPt); + ++arrOutIdxPt; + } + } /*! * This method works on a pair input (\b arrIn, \b arrIndxIn) where \b arrIn indexes is in \b arrIndxIn diff --git a/src/MEDCoupling_Swig/DataArrayInt.i b/src/MEDCoupling_Swig/DataArrayInt.i index 339c847d3..c36a32523 100644 --- a/src/MEDCoupling_Swig/DataArrayInt.i +++ b/src/MEDCoupling_Swig/DataArrayInt.i @@ -2210,6 +2210,17 @@ PyTuple_SetItem(ret,1,SWIG_NewPointerObj(SWIG_as_voidptr(arrIndexOut),SWIGTITraits::TI, SWIG_POINTER_OWN | 0 )); return ret; } + + static PyObject *FromVTKInternalReprOfPolyedra(const ARRAY *arrIn, const DataArrayIdType *arrIndxIn) throw(INTERP_KERNEL::Exception) + { + MCAuto arrOut; + MCAuto arrIndexOut; + ARRAY::FromVTKInternalReprOfPolyedra(arrIn,arrIndxIn,arrOut,arrIndexOut); + PyObject *ret=PyTuple_New(2); + PyTuple_SetItem(ret,0,SWIG_NewPointerObj(SWIG_as_voidptr(arrOut.retn()),SWIGTITraits::TI, SWIG_POINTER_OWN | 0 )); + PyTuple_SetItem(ret,1,SWIG_NewPointerObj(SWIG_as_voidptr(arrIndexOut.retn()),SWIGTITraits::TI, SWIG_POINTER_OWN | 0 )); + return ret; + } static void SetPartOfIndexedArraysSameIdx(PyObject *li, ARRAY *arrIn, const DataArrayIdType *arrIndxIn, const ARRAY *srcArr, const DataArrayIdType *srcArrIndex) throw(INTERP_KERNEL::Exception) diff --git a/src/MEDCoupling_Swig/MEDCouplingBasicsTest7.py b/src/MEDCoupling_Swig/MEDCouplingBasicsTest7.py index f827cbdf3..87e51a18a 100644 --- a/src/MEDCoupling_Swig/MEDCouplingBasicsTest7.py +++ b/src/MEDCoupling_Swig/MEDCouplingBasicsTest7.py @@ -1544,6 +1544,24 @@ class MEDCouplingBasicsTest7(unittest.TestCase): c,ci = arr.findCommonTuples(2) self.assertTrue( ci.isEqual( DataArrayInt([0,2,5]) ) ) self.assertTrue( c.isEqual( DataArrayInt([0,4, 1,6,7]) ) ) + + def testDAIFromVTKInternalReprOfPolyedra(self): + """ + EDF31315 : VTK internal representation of polyedra data structure + """ + faces = DataArrayInt64( [11, 4, 7199, 6757, 2950, 6758, 5, 6455, 1794, 2400, 6620, 7200, 6, 6620, 2400, 5864, 2950, 6757, 7244, 6, 6758, 2950, 5864, 3223, 6818, 7245, 5, 7246, 6818, 3223, 1794, 6455, 4, 1794, 3223, 5864, 2400, 4, 0, 7246, 6455, 7200, 4, 0, 7200, 6620, 7244, 4, 0, 7244, 6757, 7199, 4, 0, 7199, 6758, 7245, 4, 0, 7245, 6818, 7246, + 12, 5, 6408, 978, 1721, 6441, 7203, 4, 7204, 7007, 3987, 7008, 6, 6441, 1721, 5813, 3987, 7007, 7247, 5, 7248, 7130, 4480, 978, 6408, 6, 7008, 3987, 5813, 4480, 7130, 7249, 4, 978, 4480, 5813, 1721, 4, 1, 7248, 6408, 7203, 4, 1, 7203, 6441, 7247, 4, 1, 7247, 7007, 7204, 4, 1, 7204, 7008, 7249, 3, 1, 7249, 7130, 5, 1, 7, 6, 5, 4] ) + facesIndex = DataArrayInt( [0, 62, 129] ) + facesOut, facesIndexOut = DataArrayInt64.FromVTKInternalReprOfPolyedra(faces,facesIndex) + m =MEDCoupling1DGTUMesh("mesh",NORM_POLYHED) + m.setCoords(DataArrayDouble(10000,3)) + m.setNodalConnectivity(facesOut,facesIndexOut) + m = m.buildUnstructured() + self.assertTrue( m.computeNbOfFacesPerCell().isEqual(DataArrayInt([11,12])) ) + connExp = DataArrayInt64( [31, 7199, 6757, 2950, 6758, -1, 6455, 1794, 2400, 6620, 7200, -1, 6620, 2400, 5864, 2950, 6757, 7244, -1, 6758, 2950, 5864, 3223, 6818, 7245, -1, 7246, 6818, 3223, 1794, 6455, -1, 1794, 3223, 5864, 2400, -1, 0, 7246, 6455, 7200, -1, 0, 7200, 6620, 7244, -1, 0, 7244, 6757, 7199, -1, 0, 7199, 6758, 7245, -1, 0, 7245, 6818, 7246, 31, 6408, 978, 1721, 6441, 7203, -1, 7204, 7007, 3987, 7008, -1, 6441, 1721, 5813, 3987, 7007, 7247, -1, 7248, 7130, 4480, 978, 6408, -1, 7008, 3987, 5813, 4480, 7130, 7249, -1, 978, 4480, 5813, 1721, -1, 1, 7248, 6408, 7203, -1, 1, 7203, 6441, 7247, -1, 1, 7247, 7007, 7204, -1, 1, 7204, 7008, 7249, -1, 1, 7249, 7130, -1, 1, 7, 6, 5, 4] ) + connIndexExp = DataArrayInt( [0, 61, 127] ) + self.assertTrue( m.getNodalConnectivity().isEqual( connExp ) ) + self.assertTrue( m.getNodalConnectivityIndex().isEqual( connIndexExp ) ) if __name__ == '__main__': unittest.main() diff --git a/src/MEDCoupling_Swig/vtk2medcoupling.py b/src/MEDCoupling_Swig/vtk2medcoupling.py index 887db5729..022e30105 100644 --- a/src/MEDCoupling_Swig/vtk2medcoupling.py +++ b/src/MEDCoupling_Swig/vtk2medcoupling.py @@ -25,6 +25,22 @@ from vtk.util import numpy_support import medcoupling as mc import numpy as np +def patchForPolyedra(polyhedCellIds, ug, mesh): + """ + Method in charge to change the connectivity of polyedra contained in mesh using ug vtkUnstructuredGrid. + + :param in polyhedCellIds: mc.DataArrayInt of cells ids in mesh to be patched + :param in ug: vtkUnstructuredGrid containing polyhedra + :param in-out mesh: mc.MEDCouplingUMesh. 3D Mesh whose polyedra cells connectivity will be modified + """ + facesLoc = mc.DataArrayInt( numpy_support.vtk_to_numpy( ug.GetFaceLocations() ) ) + faces = mc.DataArrayInt( numpy_support.vtk_to_numpy( ug.GetFaces() ) ) + facesLoc = mc.DataArrayInt.Aggregate( [ facesLoc, mc.DataArrayInt([ len(faces) ]) ] ) + connForPoly, facesLoc = mc.DataArrayInt.FromVTKInternalReprOfPolyedra(faces,facesLoc) + meshPoly = mc.MEDCoupling1DGTUMesh(mesh.getName(),mc.NORM_POLYHED) ; meshPoly.setCoords( mesh.getCoords() ) ; meshPoly.setNodalConnectivity(connForPoly,facesLoc) + mesh[polyhedCellIds] = meshPoly.buildUnstructured() + pass + def mesh_convertor(fileName): #vtk.vtkDataSetReader() reader = vtk.vtkXMLUnstructuredGridReader() @@ -69,4 +85,9 @@ def mesh_convertor_mem(ug): m.setConnectivity(c,ci,True) m.checkConsistencyLight() # + if m.getMeshDimension() == 3: + polyhedCellIds = ct.findIdsEqual(mc.NORM_POLYHED) + if not polyhedCellIds.empty(): + patchForPolyedra(polyhedCellIds, ug, m) + # return m