From a015f58077df6fbadce6dbd424d35e3de6774fdf Mon Sep 17 00:00:00 2001 From: Anthony Geay Date: Thu, 25 Aug 2016 18:34:58 +0200 Subject: [PATCH] Merge 'agy/br810_2' branch. --- src/MEDCoupling/MEDCouplingUMesh.cxx | 24 +----- src/MEDLoader/MEDFileField.cxx | 7 +- src/MEDLoader/MEDFileMesh.cxx | 102 +++++++++++++++++++++++++ src/MEDLoader/MEDFileMesh.hxx | 2 + src/MEDLoader/Swig/MEDLoaderCommon.i | 16 ++++ src/MEDLoader/Swig/MEDLoaderTest3.py | 109 +++++++++++++++++++++++++++ 6 files changed, 239 insertions(+), 21 deletions(-) diff --git a/src/MEDCoupling/MEDCouplingUMesh.cxx b/src/MEDCoupling/MEDCouplingUMesh.cxx index 38da9fe17..ea521f675 100644 --- a/src/MEDCoupling/MEDCouplingUMesh.cxx +++ b/src/MEDCoupling/MEDCouplingUMesh.cxx @@ -1448,27 +1448,11 @@ void MEDCouplingUMesh::simplifyPolyhedra(double eps) DataArrayInt *MEDCouplingUMesh::computeFetchedNodeIds() const { checkConnectivityFullyDefined(); - int nbOfCells=getNumberOfCells(); - const int *connIndex=_nodal_connec_index->getConstPointer(); - const int *conn=_nodal_connec->getConstPointer(); - const int *maxEltPt=std::max_element(_nodal_connec->begin(),_nodal_connec->end()); - int maxElt=maxEltPt==_nodal_connec->end()?0:std::abs(*maxEltPt)+1; + const int *maxEltPt(std::max_element(_nodal_connec->begin(),_nodal_connec->end())); + int maxElt(maxEltPt==_nodal_connec->end()?0:std::abs(*maxEltPt)+1); std::vector retS(maxElt,false); - for(int i=0;i=0) - retS[conn[j]]=true; - int sz=0; - for(int i=0;ialloc(sz,1); - int *retPtr=ret->getPointer(); - for(int i=0;i& arrOut, const MEDFileFieldNameScope& nasc) const { - static const char MSG0[]="MEDFileAnyTypeField1TSWithoutSDA::fieldOnMesh : the field is too complex to be able to extract a field ! Call getFieldOnMeshAtLevel method instead to deal with complexity !"; + static const char MSG0[]="MEDFileAnyTypeField1TSWithoutSDA::fieldOnMesh : the field is too complex to be able to be extracted with \"field\" method ! Call getFieldOnMeshAtLevel method instead to deal with complexity !"; if(_field_per_mesh.empty()) throw INTERP_KERNEL::Exception("MEDFileAnyTypeField1TSWithoutSDA::fieldOnMesh : the field is empty ! Nothing to extract !"); if(_field_per_mesh.size()>1) @@ -6424,6 +6424,7 @@ DataArrayDouble *MEDFileField1TS::ReturnSafelyDataArrayDouble(MCAuto& * The keys of \a extractDef is level relative to max ext of \a mm mesh. * * \return A new object that the caller is responsible to deallocate. + * \sa MEDFileUMesh::deduceNodeSubPartFromCellSubPart , MEDFileUMesh::extractPart */ MEDFileField1TS *MEDFileField1TS::extractPart(const std::map >& extractDef, MEDFileMesh *mm) const { @@ -6443,6 +6444,8 @@ MEDFileField1TS *MEDFileField1TS::extractPart(const std::map t((*it2).second); + if(t.isNull()) + throw INTERP_KERNEL::Exception("MEDFileField1TS::extractPart : presence of a value with null pointer 1 !"); MCAuto f(getFieldOnMeshAtLevel(ON_CELLS,(*lev),mm)); MCAuto fOut(f->buildSubPart(t)); ret->setFieldNoProfileSBT(fOut); @@ -6455,6 +6458,8 @@ MEDFileField1TS *MEDFileField1TS::extractPart(const std::map t((*it2).second); + if(t.isNull()) + throw INTERP_KERNEL::Exception("MEDFileField1TS::extractPart : presence of a value with null pointer 1 !"); MCAuto f(getFieldOnMeshAtLevel(ON_NODES,0,mm)); MCAuto fOut(f->deepCopy()); DataArrayDouble *arr(f->getArray()); diff --git a/src/MEDLoader/MEDFileMesh.cxx b/src/MEDLoader/MEDFileMesh.cxx index 90c34e17e..d187a455e 100644 --- a/src/MEDLoader/MEDFileMesh.cxx +++ b/src/MEDLoader/MEDFileMesh.cxx @@ -4062,6 +4062,108 @@ DataArrayInt *MEDFileUMesh::zipCoords() return ret.retn(); } +/*! + * This method is a const method. It computes the minimal set of node ids covered by the cell extraction of \a this. + * The extraction of \a this is specified by the extractDef \a input map. + * This map tells for each level of cells, the cells kept in the extraction. + * + * \return - a new reference of DataArrayInt that represents sorted node ids, the extraction is lying on. + * \sa MEDFileField1TS::extractPart, MEDFileUMesh::extractPart + */ +DataArrayInt *MEDFileUMesh::deduceNodeSubPartFromCellSubPart(const std::map >& extractDef) const +{ + std::vector levs(getNonEmptyLevels()); + std::vector fetchedNodes(getNumberOfNodes(),false); + for(std::map >::const_iterator it=extractDef.begin();it!=extractDef.end();it++) + { + if((*it).first>1) + throw INTERP_KERNEL::Exception("MEDFileUMesh::deduceNodeSubPartFromCellSubPart : invalid key ! Must be <=1 !"); + if((*it).second.isNull()) + throw INTERP_KERNEL::Exception("MEDFileUMesh::deduceNodeSubPartFromCellSubPart : presence of a value with null pointer !"); + if((*it).first==1) + continue; + if(std::find(levs.begin(),levs.end(),(*it).first)==levs.end()) + { + std::ostringstream oss; oss << "MEDFileUMesh::deduceNodeSubPartFromCellSubPart : invalid level " << (*it).first << " ! Not present in this !"; + throw INTERP_KERNEL::Exception(oss.str()); + } + MCAuto m(getMeshAtLevel((*it).first)); + MCAuto mPart(m->buildPartOfMySelf((*it).second->begin(),(*it).second->end(),true)); + mPart->computeNodeIdsAlg(fetchedNodes); + } + return DataArrayInt::BuildListOfSwitchedOn(fetchedNodes); +} + +/*! + * This method returns a new MEDFileUMesh that is the result of the extraction of cells/nodes in \a this. + * + * \return - a new reference of MEDFileUMesh + * \sa MEDFileUMesh::deduceNodeSubPartFromCellSubPart, MEDFileFields::extractPart + */ +MEDFileUMesh *MEDFileUMesh::extractPart(const std::map >& extractDef) const +{ + MCAuto ret(MEDFileUMesh::New()); ret->setName(getName()); ret->copyFamGrpMapsFrom(*this); + std::vector levs(getNonEmptyLevels()); + for(std::map >::const_iterator it=extractDef.begin();it!=extractDef.end();it++) + { + if((*it).first>1) + throw INTERP_KERNEL::Exception("MEDFileUMesh::extractPart : invalid key ! Must be <=1 !"); + if((*it).second.isNull()) + throw INTERP_KERNEL::Exception("MEDFileUMesh::extractPart : presence of a value with null pointer !"); + if((*it).first==1) + continue; + if(std::find(levs.begin(),levs.end(),(*it).first)==levs.end()) + { + std::ostringstream oss; oss << "MEDFileUMesh::extractPart : invalid level " << (*it).first << " ! Not present in this !"; + throw INTERP_KERNEL::Exception(oss.str()); + } + MCAuto m(getMeshAtLevel((*it).first)); + MCAuto mPart(m->buildPartOfMySelf((*it).second->begin(),(*it).second->end(),true)); + ret->setMeshAtLevel((*it).first,mPart); + const DataArrayInt *fam(getFamilyFieldAtLevel((*it).first)),*num(getNumberFieldAtLevel((*it).first)); + if(fam) + { + MCAuto famPart(fam->selectByTupleIdSafe((*it).second->begin(),(*it).second->end())); + ret->setFamilyFieldArr((*it).first,famPart); + } + if(num) + { + MCAuto numPart(num->selectByTupleIdSafe((*it).second->begin(),(*it).second->end())); + ret->setFamilyFieldArr((*it).first,numPart); + } + } + std::map >::const_iterator it2(extractDef.find(1)); + if(it2!=extractDef.end()) + { + const DataArrayDouble *coo(ret->getCoords()); + if(!coo) + throw INTERP_KERNEL::Exception("MEDFileUMesh::extractPart : trying to extract nodes whereas there is no nodes !"); + MCAuto o2nNodes(((*it2).second)->invertArrayN2O2O2N(coo->getNumberOfTuples())); + MCAuto cooPart(coo->selectByTupleIdSafe((*it2).second->begin(),(*it2).second->end())); + ret->setCoords(cooPart); + const DataArrayInt *fam(getFamilyFieldAtLevel(1)),*num(getNumberFieldAtLevel(1)); + if(fam) + { + MCAuto famPart(fam->selectByTupleIdSafe((*it2).second->begin(),(*it2).second->end())); + ret->setFamilyFieldArr(1,famPart); + } + if(num) + { + MCAuto numPart(num->selectByTupleIdSafe((*it2).second->begin(),(*it2).second->end())); + ret->setFamilyFieldArr(1,numPart); + } + for(std::map >::const_iterator it3=extractDef.begin();it3!=extractDef.end();it3++) + { + if((*it3).first==1) + continue; + MCAuto m(ret->getMeshAtLevel((*it3).first)); + m->renumberNodesInConn(o2nNodes->begin()); + ret->setMeshAtLevel((*it3).first,m); + } + } + return ret.retn(); +} + /*! * This method performs an extrusion along a path defined by \a m1D. * \a this is expected to be a mesh with max mesh dimension equal to 2. diff --git a/src/MEDLoader/MEDFileMesh.hxx b/src/MEDLoader/MEDFileMesh.hxx index f40acc5c7..9f0231313 100644 --- a/src/MEDLoader/MEDFileMesh.hxx +++ b/src/MEDLoader/MEDFileMesh.hxx @@ -329,6 +329,8 @@ namespace MEDCoupling MEDLOADER_EXPORT void buildInnerBoundaryAlongM1Group(const std::string& grpNameM1, DataArrayInt *&nodesDuplicated, DataArrayInt *&cellsModified, DataArrayInt *&cellsNotModified); MEDLOADER_EXPORT bool unPolyze(std::vector& oldCode, std::vector& newCode, DataArrayInt *& o2nRenumCell); MEDLOADER_EXPORT DataArrayInt *zipCoords(); + MEDLOADER_EXPORT DataArrayInt *deduceNodeSubPartFromCellSubPart(const std::map >& extractDef) const; + MEDLOADER_EXPORT MEDFileUMesh *extractPart(const std::map >& extractDef) const; MEDLOADER_EXPORT MEDFileUMesh *buildExtrudedMesh(const MEDCouplingUMesh *m1D, int policy) const; MEDLOADER_EXPORT MEDFileUMesh *linearToQuadratic(int conversionType=0, double eps=1e-12) const; MEDLOADER_EXPORT MEDFileUMesh *quadraticToLinear(double eps=1e-12) const; diff --git a/src/MEDLoader/Swig/MEDLoaderCommon.i b/src/MEDLoader/Swig/MEDLoaderCommon.i index b50ddecbf..818147457 100644 --- a/src/MEDLoader/Swig/MEDLoaderCommon.i +++ b/src/MEDLoader/Swig/MEDLoaderCommon.i @@ -120,6 +120,8 @@ using namespace MEDCoupling; %newobject MEDCoupling::MEDFileUMesh::extractFamilyFieldOnGeoType; %newobject MEDCoupling::MEDFileUMesh::extractNumberFieldOnGeoType; %newobject MEDCoupling::MEDFileUMesh::zipCoords; +%newobject MEDCoupling::MEDFileUMesh::deduceNodeSubPartFromCellSubPart; +%newobject MEDCoupling::MEDFileUMesh::extractPart; %newobject MEDCoupling::MEDFileUMesh::buildExtrudedMesh; %newobject MEDCoupling::MEDFileUMesh::linearToQuadratic; %newobject MEDCoupling::MEDFileUMesh::quadraticToLinear; @@ -1392,6 +1394,20 @@ namespace MEDCoupling self->removeMeshAtLevel(meshDimRelToMax); } + DataArrayInt *deduceNodeSubPartFromCellSubPart(PyObject *extractDef) const throw(INTERP_KERNEL::Exception) + { + std::map > extractDefCpp; + convertToMapIntDataArrayInt(extractDef,extractDefCpp); + return self->deduceNodeSubPartFromCellSubPart(extractDefCpp); + } + + MEDFileUMesh *extractPart(PyObject *extractDef) const throw(INTERP_KERNEL::Exception) + { + std::map > extractDefCpp; + convertToMapIntDataArrayInt(extractDef,extractDefCpp); + return self->extractPart(extractDefCpp); + } + void setMeshes(PyObject *li, bool renum=false) throw(INTERP_KERNEL::Exception) { std::vector ms; diff --git a/src/MEDLoader/Swig/MEDLoaderTest3.py b/src/MEDLoader/Swig/MEDLoaderTest3.py index 0c6efa452..9ec3797c1 100644 --- a/src/MEDLoader/Swig/MEDLoaderTest3.py +++ b/src/MEDLoader/Swig/MEDLoaderTest3.py @@ -5418,6 +5418,115 @@ class MEDLoaderTest3(unittest.TestCase): mm=MEDFileMesh.New(fname) ; f1ts=MEDFileField1TS(fname,fieldName,6,7) self.assertTrue(f.isEqual(f1ts.field(mm),1e-12,1e-12)) pass + + def testExtractPart1(self): + coo=DataArrayDouble([(0,0),(1,0),(2,0),(3,0),(4,0),(0,1),(1,1),(2,1),(3,1),(4,1),(0,2),(1,2),(2,2),(3,2),(4,2)]) + meshName="mesh" + m0=MEDCouplingUMesh(meshName,2) ; m0.setCoords(coo) ; m0.allocateCells() + m0.insertNextCell(NORM_TRI3,[8,4,3]) + m0.insertNextCell(NORM_TRI3,[8,9,4]) + m0.insertNextCell(NORM_TRI3,[7,13,8]) + m0.insertNextCell(NORM_TRI3,[7,12,13]) + m0.insertNextCell(NORM_TRI3,[0,6,1]) + m0.insertNextCell(NORM_TRI3,[0,5,6]) + m0.insertNextCell(NORM_QUAD4,[1,6,7,2]) + m0.insertNextCell(NORM_QUAD4,[2,7,8,3]) + m0.insertNextCell(NORM_QUAD4,[8,13,14,9]) + m0.insertNextCell(NORM_QUAD4,[6,11,12,7]) + m0.insertNextCell(NORM_QUAD4,[5,10,11,6]) + # + m1=MEDCouplingUMesh(meshName,1) ; m1.setCoords(coo) ; m1.allocateCells() + m1.insertNextCell(NORM_SEG2,[10,5]) + m1.insertNextCell(NORM_SEG2,[5,0]) + m1.insertNextCell(NORM_SEG2,[0,1]) + m1.insertNextCell(NORM_SEG2,[1,2]) + m1.insertNextCell(NORM_SEG2,[2,3]) + m1.insertNextCell(NORM_SEG2,[3,4]) + m1.insertNextCell(NORM_SEG2,[4,9]) + m1.insertNextCell(NORM_SEG2,[9,14]) + m1.insertNextCell(NORM_SEG2,[14,13]) + m1.insertNextCell(NORM_SEG2,[13,12]) + m1.insertNextCell(NORM_SEG2,[12,11]) + m1.insertNextCell(NORM_SEG2,[11,10]) + mm=MEDFileUMesh() + mm[0]=m0 ; mm[-1]=m1 + arr0=DataArrayInt([0,1,2,3,4,6,7,8,12,13]) + tab={} # + tab[0]=DataArrayInt([0,2,3,4,6,7]) + tab[-1]=DataArrayInt([2,3,4,5,9]) + fs=MEDFileFields() + self.assertTrue(mm.deduceNodeSubPartFromCellSubPart(tab).isEqual(arr0)) + tab[1]=arr0 + # + fname0="Field0" + fmts=MEDFileFieldMultiTS() ; fs.pushField(fmts) + t0=(16.5,3,4) + ic=["toto [m]"] + arr0_0=DataArrayDouble([100,101,102,103,104,105,106,107,108,109,110]) ; arr0_0.setInfoOnComponents(ic) + f0=MEDCouplingFieldDouble(ON_CELLS) ; f0.setTime(*t0) ; f0.setArray(arr0_0) + f0.setMesh(m0) ; f0.setName(fname0) + f1=MEDCouplingFieldDouble(ON_CELLS) ; f1.setTime(*t0) ; f1.setArray(DataArrayDouble([200,201,202,203,204,205,206,207,208,209,210,211])) + f1.setMesh(m1) ; f1.setName(fname0) ; f1.getArray().setInfoOnComponents(ic) + f2=MEDCouplingFieldDouble(ON_NODES) ; f2.setTime(*t0) ; f2.setArray(DataArrayDouble([300,301,302,303,304,305,306,307,308,309,310,311,312,313,314])) + f2.setMesh(m0) ; f2.setName(fname0) ; f2.getArray().setInfoOnComponents(ic) + f1ts=MEDFileField1TS() ; f1ts.setFieldNoProfileSBT(f0) ; f1ts.setFieldNoProfileSBT(f1) ; f1ts.setFieldNoProfileSBT(f2) + fmts.pushBackTimeStep(f1ts) + # + mmOut=mm.extractPart(tab) + # + fsPart0=fs.extractPart(tab,mm) + self.assertEqual(len(fsPart0),1) + fmtsP=fsPart0[0] + self.assertEqual(len(fmtsP),1) + f1ts=fmtsP[0] + self.assertRaises(InterpKernelException,f1ts.field,mmOut) + # + self.assertTrue(mmOut[0].computeCellCenterOfMass().isEqual(m0[tab[0]].computeCellCenterOfMass(),1e-12)) + self.assertTrue(mmOut[-1].computeCellCenterOfMass().isEqual(m1[tab[-1]].computeCellCenterOfMass(),1e-12)) + # + m0Part=m0.deepCopy()[tab[0]] ; m0Part.renumberNodes(tab[1].invertArrayN2O2O2N(mm.getNumberOfNodes()),len(tab[1])) ; m0Part.setName(m0.getName()) + self.assertTrue(mmOut[0].isEqual(m0Part,1e-12)) + m1Part=m1.deepCopy()[tab[-1]] ; m1Part.renumberNodes(tab[1].invertArrayN2O2O2N(mm.getNumberOfNodes()),len(tab[1])) ; m1Part.setName(m0.getName()) + self.assertTrue(mmOut[0].isEqual(m0Part,1e-12)) + self.assertTrue(mmOut[-1].isEqual(m1Part,1e-12)) + # + f0Part=f1ts.getFieldOnMeshAtLevel(ON_CELLS,0,mmOut) ; f0Part.checkConsistencyLight() + self.assertEqual(f0Part.getTypeOfField(),ON_CELLS) + self.assertTrue(f0Part.getMesh().isEqual(m0Part,1e-12)) + arr0Exp=DataArrayDouble([100,102,103,104,106,107]) ; arr0Exp.setInfoOnComponents(ic) + self.assertTrue(f0Part.getArray().isEqual(arr0Exp,1e-12)) ; self.assertEqual(f0Part.getTime(),list(t0)) + f1Part=f1ts.getFieldOnMeshAtLevel(ON_CELLS,-1,mmOut) ; f1Part.checkConsistencyLight() + self.assertEqual(f1Part.getTypeOfField(),ON_CELLS) + self.assertTrue(f1Part.getMesh().isEqual(m1Part,1e-12)) + arr1Exp=DataArrayDouble([202,203,204,205,209]) ; arr1Exp.setInfoOnComponents(ic) + self.assertTrue(f1Part.getArray().isEqual(arr1Exp,1e-12)) ; self.assertEqual(f1Part.getTime(),list(t0)) + # + f2Part=f1ts.getFieldOnMeshAtLevel(ON_NODES,0,mmOut) ; f2Part.checkConsistencyLight() + arr2Exp=DataArrayDouble([300,301,302,303,304,306,307,308,312,313]) ; arr2Exp.setInfoOnComponents(ic) + self.assertTrue(f2Part.getArray().isEqual(arr2Exp,1e-12)) ; self.assertEqual(f2Part.getTime(),list(t0)) + # multisteps + fs=MEDFileFields() ; fmts=MEDFileFieldMultiTS() ; fs.pushField(fmts) + tss=[(16.5,3,4),(17.5,4,5),(18.5,5,6)] + for i,tt in enumerate(tss): + f0=MEDCouplingFieldDouble(ON_CELLS) ; f0.setTime(*tt) + myarr=arr0_0+i*1000. + f0.setArray(myarr) + f0.setMesh(m0) ; f0.setName(fname0) ; f0.getArray().setInfoOnComponents(ic) + f1ts=MEDFileField1TS() ; f1ts.setFieldNoProfileSBT(f0) ; fmts.pushBackTimeStep(f1ts) + pass + fsPart1=fs.extractPart(tab,mm) + self.assertEqual(len(fsPart1),1) + fmtsP=fsPart1[0] + self.assertEqual(len(fmtsP),len(tss)) + for i,(f1tsP,tt) in enumerate(zip(fmtsP,tss)): + fPart=f1tsP.field(mmOut) ; fPart.checkConsistencyLight() + self.assertEqual(fPart.getTypeOfField(),ON_CELLS) + arr0Exp=DataArrayDouble([100,102,103,104,106,107]) ; arr0Exp.setInfoOnComponents(ic) ; arr0Exp+=i*1000. + self.assertTrue(fPart.getMesh().isEqual(m0Part,1e-12)) + self.assertTrue(fPart.getArray().isEqual(arr0Exp,1e-12)) + self.assertEqual(fPart.getTime(),list(tt)) + pass + pass pass -- 2.39.2