From 202cb592e8e8254a04285574de68ec9693b2978a Mon Sep 17 00:00:00 2001 From: ageay Date: Mon, 2 Jul 2012 14:46:24 +0000 Subject: [PATCH] Addition of extractSlice3D in fields. Implementation of buildSubPart for Gauss and GaussNE. --- src/MEDCoupling/MEDCouplingField.cxx | 9 ++ src/MEDCoupling/MEDCouplingField.hxx | 1 + .../MEDCouplingFieldDiscretization.cxx | 89 ++++++++++++++++++- .../MEDCouplingFieldDiscretization.hxx | 5 ++ src/MEDCoupling/MEDCouplingFieldDouble.cxx | 42 ++++++++- src/MEDCoupling/MEDCouplingFieldDouble.hxx | 1 + src/MEDCoupling/MEDCouplingUMesh.cxx | 27 ++++++ src/MEDCoupling/MEDCouplingUMesh.hxx | 1 + .../Test/MEDCouplingBasicsTest5.cxx | 74 +++++++++++++++ .../Test/MEDCouplingBasicsTest5.hxx | 4 + src/MEDCoupling_Swig/MEDCouplingBasicsTest.py | 40 +++++++++ src/MEDCoupling_Swig/MEDCouplingCommon.i | 47 ++++++++++ 12 files changed, 335 insertions(+), 5 deletions(-) diff --git a/src/MEDCoupling/MEDCouplingField.cxx b/src/MEDCoupling/MEDCouplingField.cxx index 46f4f5b29..9b0daf593 100644 --- a/src/MEDCoupling/MEDCouplingField.cxx +++ b/src/MEDCoupling/MEDCouplingField.cxx @@ -312,6 +312,15 @@ MEDCouplingMesh *MEDCouplingField::buildSubMeshData(const int *start, const int return _type->buildSubMeshData(_mesh,start,end,di); } +/*! + * This method returns tuples ids implied by the mesh selection of the cell ids contained in array defined as an interval [start;end). + * \return a newly allocated DataArrayInt instance containing tuples ids. + */ +DataArrayInt *MEDCouplingField::computeTupleIdsToSelectFromCellIds(const int *startCellIds, const int *endCellIds) const +{ + return _type->computeTupleIdsToSelectFromCellIds(_mesh,startCellIds,endCellIds); +} + /*! * This method returns number of tuples expected regarding its discretization and its _mesh attribute. * This method expected a not null _mesh instance. If null, an exception will be thrown. diff --git a/src/MEDCoupling/MEDCouplingField.hxx b/src/MEDCoupling/MEDCouplingField.hxx index 74093090d..87e052597 100644 --- a/src/MEDCoupling/MEDCouplingField.hxx +++ b/src/MEDCoupling/MEDCouplingField.hxx @@ -59,6 +59,7 @@ namespace ParaMEDMEM DataArrayDouble *getLocalizationOfDiscr() const throw(INTERP_KERNEL::Exception); MEDCouplingFieldDouble *buildMeasureField(bool isAbs) const throw(INTERP_KERNEL::Exception); MEDCouplingMesh *buildSubMeshData(const int *start, const int *end, DataArrayInt *&di) const; + DataArrayInt *computeTupleIdsToSelectFromCellIds(const int *startCellIds, const int *endCellIds) const; MEDCouplingFieldDiscretization *getDiscretization() const { return _type; } int getNumberOfTuplesExpected() const throw(INTERP_KERNEL::Exception); int getNumberOfMeshPlacesExpected() const throw(INTERP_KERNEL::Exception); diff --git a/src/MEDCoupling/MEDCouplingFieldDiscretization.cxx b/src/MEDCoupling/MEDCouplingFieldDiscretization.cxx index acbd40f3f..3e3a5e082 100644 --- a/src/MEDCoupling/MEDCouplingFieldDiscretization.cxx +++ b/src/MEDCoupling/MEDCouplingFieldDiscretization.cxx @@ -478,7 +478,23 @@ void MEDCouplingFieldDiscretizationP0::renumberValuesOnCellsR(const MEDCouplingM } /*! - * This method returns a submesh of 'mesh' instance constituting cell ids contained in array defined as an interval [start;end]. + * This method returns a tuple ids selection from cell ids selection [start;end). + * This method is called by MEDCouplingFieldDiscretizationP0::buildSubMeshData to return parameter \b di. + * Here for P0 it's very simple ! + * + * \return a newly allocated array containing ids to select into the DataArrayDouble of the field. + * + */ +DataArrayInt *MEDCouplingFieldDiscretizationP0::computeTupleIdsToSelectFromCellIds(const MEDCouplingMesh *mesh, const int *startCellIds, const int *endCellIds) const +{ + MEDCouplingAutoRefCountObjectPtr ret=DataArrayInt::New(); + ret->alloc((int)std::distance(startCellIds,endCellIds),1); + std::copy(startCellIds,endCellIds,ret->getPointer()); + ret->incrRef(); return ret; +} + +/*! + * This method returns a submesh of 'mesh' instance constituting cell ids contained in array defined as an interval [start;end). * @param di is an array returned that specifies entity ids (here cells ids) in mesh 'mesh' of entity in returned submesh. * Example : The first cell id of returned mesh has the (*di)[0] id in 'mesh' */ @@ -683,6 +699,23 @@ MEDCouplingMesh *MEDCouplingFieldDiscretizationP1::buildSubMeshData(const MEDCou return ret; } +/*! + * This method returns a tuple ids selection from cell ids selection [start;end). + * This method is called by MEDCouplingFieldDiscretizationP0::buildSubMeshData to return parameter \b di. + * Here for P1 only nodes fetched by submesh of mesh[startCellIds:endCellIds) is returned ! + * + * \return a newly allocated array containing ids to select into the DataArrayDouble of the field. + * + */ +DataArrayInt *MEDCouplingFieldDiscretizationP1::computeTupleIdsToSelectFromCellIds(const MEDCouplingMesh *mesh, const int *startCellIds, const int *endCellIds) const +{ + if(!mesh) + throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretizationP1::computeTupleIdsToSelectFromCellIds : null mesh !"); + const MEDCouplingAutoRefCountObjectPtr umesh=mesh->buildUnstructured(); + MEDCouplingAutoRefCountObjectPtr umesh2=static_cast(umesh->buildPartOfMySelf(startCellIds,endCellIds,true)); + return umesh2->computeFetchedNodeIds(); +} + MEDCouplingFieldDiscretizationPerCell::MEDCouplingFieldDiscretizationPerCell():_discr_per_cell(0) { } @@ -1101,7 +1134,38 @@ DataArrayDouble *MEDCouplingFieldDiscretizationGauss::getValueOnMulti(const Data MEDCouplingMesh *MEDCouplingFieldDiscretizationGauss::buildSubMeshData(const MEDCouplingMesh *mesh, const int *start, const int *end, DataArrayInt *&di) const { - throw INTERP_KERNEL::Exception("Not implemented yet !"); + di=computeTupleIdsToSelectFromCellIds(mesh,start,end); + return mesh->buildPart(start,end); +} + +/*! + * This method returns a tuple ids selection from cell ids selection [start;end). + * This method is called by MEDCouplingFieldDiscretizationGauss::buildSubMeshData to return parameter \b di. + * + * \return a newly allocated array containing ids to select into the DataArrayDouble of the field. + * + */ +DataArrayInt *MEDCouplingFieldDiscretizationGauss::computeTupleIdsToSelectFromCellIds(const MEDCouplingMesh *mesh, const int *startCellIds, const int *endCellIds) const +{ + if(!mesh) + throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretizationGauss::computeTupleIdsToSelectFromCellIds : null mesh !"); + if(!_discr_per_cell) + throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretizationGauss::computeTupleIdsToSelectFromCellIds : null discretization ids !"); + int nbOfCells=mesh->getNumberOfCells(); + if(_discr_per_cell->getNumberOfTuples()!=nbOfCells) + throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretizationGauss::computeTupleIdsToSelectFromCellIds : mismatch of nb of tuples of cell ids array and number of cells !"); + MEDCouplingAutoRefCountObjectPtr nbOfNodesPerCell=DataArrayInt::New(); nbOfNodesPerCell->alloc(nbOfCells,1); + int *retPtr=nbOfNodesPerCell->getPointer(); + const int *pt=_discr_per_cell->getConstPointer(); + int nbMaxOfLocId=(int)_loc.size(); + for(int i=0;i=0 && *ptcomputeOffsets2(); + MEDCouplingAutoRefCountObjectPtr sel=DataArrayInt::New(); sel->useArray(startCellIds,false,CPP_DEALLOC,(int)std::distance(startCellIds,endCellIds),1); + return sel->buildExplicitArrByRanges(nbOfNodesPerCell); } /*! @@ -1531,7 +1595,26 @@ DataArrayDouble *MEDCouplingFieldDiscretizationGaussNE::getValueOnMulti(const Da MEDCouplingMesh *MEDCouplingFieldDiscretizationGaussNE::buildSubMeshData(const MEDCouplingMesh *mesh, const int *start, const int *end, DataArrayInt *&di) const { - throw INTERP_KERNEL::Exception("Not implemented yet !"); + di=computeTupleIdsToSelectFromCellIds(mesh,start,end); + return mesh->buildPart(start,end); +} + +/*! + * This method returns a tuple ids selection from cell ids selection [start;end). + * This method is called by MEDCouplingFieldDiscretizationGaussNE::buildSubMeshData to return parameter \b di. + * + * \return a newly allocated array containing ids to select into the DataArrayDouble of the field. + * + */ +DataArrayInt *MEDCouplingFieldDiscretizationGaussNE::computeTupleIdsToSelectFromCellIds(const MEDCouplingMesh *mesh, const int *startCellIds, const int *endCellIds) const +{ + if(!mesh) + throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretizationGaussNE::computeTupleIdsToSelectFromCellIds : null mesh !"); + const MEDCouplingAutoRefCountObjectPtr umesh=mesh->buildUnstructured(); + MEDCouplingAutoRefCountObjectPtr nbOfNodesPerCell=umesh->computeNbOfNodesPerCell(); + nbOfNodesPerCell->computeOffsets2(); + MEDCouplingAutoRefCountObjectPtr sel=DataArrayInt::New(); sel->useArray(startCellIds,false,CPP_DEALLOC,(int)std::distance(startCellIds,endCellIds),1); + return sel->buildExplicitArrByRanges(nbOfNodesPerCell); } /*! diff --git a/src/MEDCoupling/MEDCouplingFieldDiscretization.hxx b/src/MEDCoupling/MEDCouplingFieldDiscretization.hxx index af6b296f0..44aa5a9c6 100644 --- a/src/MEDCoupling/MEDCouplingFieldDiscretization.hxx +++ b/src/MEDCoupling/MEDCouplingFieldDiscretization.hxx @@ -69,6 +69,7 @@ namespace ParaMEDMEM virtual void getValueOn(const DataArrayDouble *arr, const MEDCouplingMesh *mesh, const double *loc, double *res) const = 0; virtual void getValueOnPos(const DataArrayDouble *arr, const MEDCouplingMesh *mesh, int i, int j, int k, double *res) const = 0; virtual DataArrayDouble *getValueOnMulti(const DataArrayDouble *arr, const MEDCouplingMesh *mesh, const double *loc, int nbOfPoints) const = 0; + virtual DataArrayInt *computeTupleIdsToSelectFromCellIds(const MEDCouplingMesh *mesh, const int *startCellIds, const int *endCellIds) const = 0; virtual MEDCouplingMesh *buildSubMeshData(const MEDCouplingMesh *mesh, const int *start, const int *end, DataArrayInt *&di) const = 0; virtual void renumberValuesOnNodes(double epsOnVals, const int *old2New, DataArrayDouble *arr) const = 0; virtual void renumberValuesOnCells(double epsOnVals, const MEDCouplingMesh *mesh, const int *old2New, DataArrayDouble *arr) const = 0; @@ -125,6 +126,7 @@ namespace ParaMEDMEM void renumberValuesOnCells(double epsOnVals, const MEDCouplingMesh *mesh, const int *old2New, DataArrayDouble *arr) const; void renumberValuesOnCellsR(const MEDCouplingMesh *mesh, const int *new2old, int newSz, DataArrayDouble *arr) const; MEDCouplingMesh *buildSubMeshData(const MEDCouplingMesh *mesh, const int *start, const int *end, DataArrayInt *&di) const; + DataArrayInt *computeTupleIdsToSelectFromCellIds(const MEDCouplingMesh *mesh, const int *startCellIds, const int *endCellIds) const; public: static const char REPR[]; static const TypeOfField TYPE; @@ -153,6 +155,7 @@ namespace ParaMEDMEM void getValueOnPos(const DataArrayDouble *arr, const MEDCouplingMesh *mesh, int i, int j, int k, double *res) const; DataArrayDouble *getValueOnMulti(const DataArrayDouble *arr, const MEDCouplingMesh *mesh, const double *loc, int nbOfPoints) const; MEDCouplingMesh *buildSubMeshData(const MEDCouplingMesh *mesh, const int *start, const int *end, DataArrayInt *&di) const; + DataArrayInt *computeTupleIdsToSelectFromCellIds(const MEDCouplingMesh *mesh, const int *startCellIds, const int *endCellIds) const; void renumberValuesOnNodes(double epsOnVals, const int *old2New, DataArrayDouble *arr) const; void renumberValuesOnCells(double epsOnVals, const MEDCouplingMesh *mesh, const int *old2New, DataArrayDouble *arr) const; void renumberValuesOnCellsR(const MEDCouplingMesh *mesh, const int *new2old, int newSz, DataArrayDouble *arr) const; @@ -219,6 +222,7 @@ namespace ParaMEDMEM void getValueOnPos(const DataArrayDouble *arr, const MEDCouplingMesh *mesh, int i, int j, int k, double *res) const; DataArrayDouble *getValueOnMulti(const DataArrayDouble *arr, const MEDCouplingMesh *mesh, const double *loc, int nbOfPoints) const; MEDCouplingMesh *buildSubMeshData(const MEDCouplingMesh *mesh, const int *start, const int *end, DataArrayInt *&di) const; + DataArrayInt *computeTupleIdsToSelectFromCellIds(const MEDCouplingMesh *mesh, const int *startCellIds, const int *endCellIds) const; void renumberValuesOnNodes(double epsOnVals, const int *old2New, DataArrayDouble *arr) const; void renumberValuesOnCells(double epsOnVals, const MEDCouplingMesh *mesh, const int *old2New, DataArrayDouble *arr) const; void renumberValuesOnCellsR(const MEDCouplingMesh *mesh, const int *new2old, int newSz, DataArrayDouble *arr) const; @@ -275,6 +279,7 @@ namespace ParaMEDMEM void getValueOnPos(const DataArrayDouble *arr, const MEDCouplingMesh *mesh, int i, int j, int k, double *res) const; DataArrayDouble *getValueOnMulti(const DataArrayDouble *arr, const MEDCouplingMesh *mesh, const double *loc, int nbOfPoints) const; MEDCouplingMesh *buildSubMeshData(const MEDCouplingMesh *mesh, const int *start, const int *end, DataArrayInt *&di) const; + DataArrayInt *computeTupleIdsToSelectFromCellIds(const MEDCouplingMesh *mesh, const int *startCellIds, const int *endCellIds) const; void renumberValuesOnNodes(double epsOnVals, const int *old2New, DataArrayDouble *arr) const; void renumberValuesOnCells(double epsOnVals, const MEDCouplingMesh *mesh, const int *old2New, DataArrayDouble *arr) const; void renumberValuesOnCellsR(const MEDCouplingMesh *mesh, const int *new2old, int newSz, DataArrayDouble *arr) const; diff --git a/src/MEDCoupling/MEDCouplingFieldDouble.cxx b/src/MEDCoupling/MEDCouplingFieldDouble.cxx index dd1ba78c7..17054653d 100644 --- a/src/MEDCoupling/MEDCouplingFieldDouble.cxx +++ b/src/MEDCoupling/MEDCouplingFieldDouble.cxx @@ -423,8 +423,8 @@ MEDCouplingFieldDouble *MEDCouplingFieldDouble::buildSubPart(const int *partBg, std::vector arrays; _time_discr->getArrays(arrays); std::vector arrs; - const int *arrSelBg=arrSelect->getConstPointer(); - const int *arrSelEnd=arrSelBg+arrSelect->getNbOfElems(); + const int *arrSelBg=arrSelect->begin(); + const int *arrSelEnd=arrSelect->end(); for(std::vector::const_iterator iter=arrays.begin();iter!=arrays.end();iter++) { DataArrayDouble *arr=0; @@ -1309,6 +1309,44 @@ bool MEDCouplingFieldDouble::zipConnectivity(int compType, double epsOnVals) thr return false; } +/*! + * This method calls MEDCouplingUMesh::buildSlice3D method. So this method makes the assumption that underlying mesh exists. + * For the moment, this method is implemented for fields on cells. + * + * \return a newly allocated field double containing the result that the user should deallocate. + */ +MEDCouplingFieldDouble *MEDCouplingFieldDouble::extractSlice3D(const double *origin, const double *vec, double eps) const throw(INTERP_KERNEL::Exception) +{ + const MEDCouplingMesh *mesh=getMesh(); + if(!mesh) + throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::extractSlice3D : underlying mesh is null !"); + if(getTypeOfField()!=ON_CELLS) + throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::extractSlice3D : only implemented for fields on cells !"); + const MEDCouplingAutoRefCountObjectPtr umesh(mesh->buildUnstructured()); + MEDCouplingAutoRefCountObjectPtr ret=clone(false); + ret->setMesh(umesh); + DataArrayInt *cellIds=0; + MEDCouplingAutoRefCountObjectPtr mesh2=umesh->buildSlice3D(origin,vec,eps,cellIds); + MEDCouplingAutoRefCountObjectPtr cellIds2=cellIds; + ret->setMesh(mesh2); + MEDCouplingAutoRefCountObjectPtr tupleIds=computeTupleIdsToSelectFromCellIds(cellIds->begin(),cellIds->end()); + std::vector arrays; + _time_discr->getArrays(arrays); + int i=0; + std::vector newArr(arrays.size()); + std::vector< MEDCouplingAutoRefCountObjectPtr > newArr2(arrays.size()); + for(std::vector::const_iterator iter=arrays.begin();iter!=arrays.end();iter++,i++) + { + if(*iter) + { + newArr2[i]=(*iter)->selectByTupleIdSafe(cellIds->begin(),cellIds->end()); + newArr[i]=newArr2[i]; + } + } + ret->setArrays(newArr); + ret->incrRef(); return ret; +} + /*! * This method applyies ParaMEDMEM::MEDCouplingUMesh::simplexize on 'this->_mesh'. * The semantic of 'policy' is given in ParaMEDMEM::MEDCouplingUMesh::simplexize method. diff --git a/src/MEDCoupling/MEDCouplingFieldDouble.hxx b/src/MEDCoupling/MEDCouplingFieldDouble.hxx index cb6204e01..67d9f8d11 100644 --- a/src/MEDCoupling/MEDCouplingFieldDouble.hxx +++ b/src/MEDCoupling/MEDCouplingFieldDouble.hxx @@ -138,6 +138,7 @@ namespace ParaMEDMEM bool mergeNodes2(double eps, double epsOnVals=1e-15) throw(INTERP_KERNEL::Exception); bool zipCoords(double epsOnVals=1e-15) throw(INTERP_KERNEL::Exception); bool zipConnectivity(int compType, double epsOnVals=1e-15) throw(INTERP_KERNEL::Exception); + MEDCouplingFieldDouble *extractSlice3D(const double *origin, const double *vec, double eps) const throw(INTERP_KERNEL::Exception); bool simplexize(int policy) throw(INTERP_KERNEL::Exception); MEDCouplingFieldDouble *doublyContractedProduct() const throw(INTERP_KERNEL::Exception); MEDCouplingFieldDouble *determinant() const throw(INTERP_KERNEL::Exception); diff --git a/src/MEDCoupling/MEDCouplingUMesh.cxx b/src/MEDCoupling/MEDCouplingUMesh.cxx index b11c90712..a8b16f9ec 100644 --- a/src/MEDCoupling/MEDCouplingUMesh.cxx +++ b/src/MEDCoupling/MEDCouplingUMesh.cxx @@ -1081,6 +1081,33 @@ DataArrayInt *MEDCouplingUMesh::getNodeIdsInUse(int& nbrOfNodesInUse) const thro return ret; } +/*! + * This method returns a newly allocated array containing this->getNumberOfCells() tuples and 1 component. + * For each cell in \b this the number of nodes constituting cell is computed. + * Excepted for poyhedrons, the result can be deduced by performing a deltaShiftIndex on the nodal connectivity index in \b this minus 1. + * For polyhedrons, the face separation (-1) are excluded from the couting. + * + * \return a newly allocated array + */ +DataArrayInt *MEDCouplingUMesh::computeNbOfNodesPerCell() const throw(INTERP_KERNEL::Exception) +{ + checkConnectivityFullyDefined(); + int nbOfCells=getNumberOfCells(); + MEDCouplingAutoRefCountObjectPtr ret=DataArrayInt::New(); + ret->alloc(nbOfCells,1); + int *retPtr=ret->getPointer(); + const int *conn=getNodalConnectivity()->getConstPointer(); + const int *connI=getNodalConnectivityIndex()->getConstPointer(); + for(int i=0;iincrRef(); return ret; +} + /*! * Array returned is the correspondance in \b old \b to \b new format. The returned array is newly created and should be dealt by the caller. * The maximum value stored in returned array is the number of nodes of 'this' minus 1 after call of this method. diff --git a/src/MEDCoupling/MEDCouplingUMesh.hxx b/src/MEDCoupling/MEDCouplingUMesh.hxx index 8760441df..f771a4020 100644 --- a/src/MEDCoupling/MEDCouplingUMesh.hxx +++ b/src/MEDCoupling/MEDCouplingUMesh.hxx @@ -104,6 +104,7 @@ namespace ParaMEDMEM MEDCOUPLING_EXPORT std::vector partitionBySpreadZone() const throw(INTERP_KERNEL::Exception); MEDCOUPLING_EXPORT DataArrayInt *computeFetchedNodeIds() const throw(INTERP_KERNEL::Exception); MEDCOUPLING_EXPORT DataArrayInt *getNodeIdsInUse(int& nbrOfNodesInUse) const throw(INTERP_KERNEL::Exception); + MEDCOUPLING_EXPORT DataArrayInt *computeNbOfNodesPerCell() const throw(INTERP_KERNEL::Exception); MEDCOUPLING_EXPORT DataArrayInt *zipCoordsTraducer() throw(INTERP_KERNEL::Exception); MEDCOUPLING_EXPORT DataArrayInt *zipConnectivityTraducer(int compType) throw(INTERP_KERNEL::Exception); MEDCOUPLING_EXPORT bool areCellsIncludedIn(const MEDCouplingUMesh *other, int compType, DataArrayInt *& arr) const throw(INTERP_KERNEL::Exception); diff --git a/src/MEDCoupling/Test/MEDCouplingBasicsTest5.cxx b/src/MEDCoupling/Test/MEDCouplingBasicsTest5.cxx index 7e52bc1bf..f5cb3df6b 100644 --- a/src/MEDCoupling/Test/MEDCouplingBasicsTest5.cxx +++ b/src/MEDCoupling/Test/MEDCouplingBasicsTest5.cxx @@ -1252,3 +1252,77 @@ void MEDCouplingBasicsTest5::testGiveCellsWithType1() // m->decrRef(); } + +void MEDCouplingBasicsTest5::testBuildSlice3D2() +{ + MEDCouplingUMesh *mesh2D=0; + MEDCouplingUMesh *mesh3D=build3DExtrudedUMesh_1(mesh2D); + mesh2D->decrRef(); + // First slice in the middle of 3D cells + const double vec1[3]={-0.07,1.,0.07}; + const double origin1[3]={1.524,1.4552,1.74768}; + DataArrayInt *ids=0; + MEDCouplingUMesh *slice1=mesh3D->buildSlice3D(origin1,vec1,1e-10,ids); + // + MEDCouplingFieldDouble *f=MEDCouplingFieldDouble::New(ON_CELLS,ONE_TIME); + f->setTime(4.5,6,7) ; f->setMesh(mesh3D); + DataArrayDouble *arr=DataArrayDouble::New(); arr->alloc(mesh3D->getNumberOfCells(),2); + arr->rearrange(1); arr->iota(2.); arr->rearrange(2); + f->setArray(arr); + f->checkCoherency(); + const int exp1[9]={1,3,4,7,9,10,13,15,16}; + DataArrayInt *expected1=DataArrayInt::New(); expected1->alloc(9,1); std::copy(exp1,exp1+9,expected1->getPointer()); + CPPUNIT_ASSERT(expected1->isEqual(*ids)); + DataArrayDouble *arr2=arr->selectByTupleIdSafe(expected1->begin(),expected1->end()); + // + MEDCouplingFieldDouble *f2=f->extractSlice3D(origin1,vec1,1e-10); + CPPUNIT_ASSERT(f2->getArray()->isEqual(*arr2,1e-12)); + CPPUNIT_ASSERT(slice1->isEqual(f2->getMesh(),1e-12)); + int a,b; + double c=f2->getTime(a,b); + CPPUNIT_ASSERT_EQUAL(6,a); + CPPUNIT_ASSERT_EQUAL(7,b); + CPPUNIT_ASSERT_DOUBLES_EQUAL(4.5,c,1e-12); + // + ids->decrRef(); + slice1->decrRef(); + arr2->decrRef(); + arr->decrRef(); + f2->decrRef(); + f->decrRef(); + mesh3D->decrRef(); + expected1->decrRef(); +} + +void MEDCouplingBasicsTest5::testComputeTupleIdsToSelectFromCellIds1() +{ + MEDCouplingUMesh *m=build2DTargetMesh_3(); + MEDCouplingFieldDouble *f=MEDCouplingFieldDouble::New(ON_GAUSS_NE,NO_TIME); + f->setMesh(m); + DataArrayDouble *arr=DataArrayDouble::New(); arr->alloc(52,2) ; arr->rearrange(1) ; arr->iota(7.); arr->rearrange(2); + f->setArray(arr); + // + const int subPart1[3]={1,5,9}; + MEDCouplingFieldDouble *f2=f->buildSubPart(subPart1,subPart1+3); + f2->checkCoherency(); + DataArrayInt *cI=m->computeNbOfNodesPerCell(); + cI->computeOffsets2(); + const int sel1[3]={1,5,9}; + DataArrayInt *sel=DataArrayInt::New(); sel->useArray(sel1,false,CPP_DEALLOC,3,1); + DataArrayInt *res=sel->buildExplicitArrByRanges(cI); + DataArrayDouble *arr2=arr->selectByTupleIdSafe(res->begin(),res->end()); + const double expected1[30]={13.,14.,15.,16.,17.,18.,19.,20.,59.,60.,61.,62.,63.,64.,95.,96.,97.,98.,99.,100.,101.,102.,103.,104.,105.,106.,107.,108.,109.,110.}; + DataArrayDouble *arr3=DataArrayDouble::New(); arr3->useArray(expected1,false,CPP_DEALLOC,15,2); + CPPUNIT_ASSERT(arr2->isEqual(*arr3,1e-12)); + CPPUNIT_ASSERT(arr2->isEqual(*f2->getArray(),1e-12)); + // + cI->decrRef(); + arr3->decrRef(); + arr2->decrRef(); + arr->decrRef(); + m->decrRef(); + f->decrRef(); + f2->decrRef(); + sel->decrRef(); + res->decrRef(); +} diff --git a/src/MEDCoupling/Test/MEDCouplingBasicsTest5.hxx b/src/MEDCoupling/Test/MEDCouplingBasicsTest5.hxx index d7a2b8ca4..67b7b7c93 100644 --- a/src/MEDCoupling/Test/MEDCouplingBasicsTest5.hxx +++ b/src/MEDCoupling/Test/MEDCouplingBasicsTest5.hxx @@ -58,6 +58,8 @@ namespace ParaMEDMEM CPPUNIT_TEST( testDataArraySort1 ); CPPUNIT_TEST( testPartitionBySpreadZone1 ); CPPUNIT_TEST( testGiveCellsWithType1 ); + CPPUNIT_TEST( testBuildSlice3D2 ); + CPPUNIT_TEST( testComputeTupleIdsToSelectFromCellIds1 ); CPPUNIT_TEST_SUITE_END(); public: void testUMeshTessellate2D1(); @@ -83,6 +85,8 @@ namespace ParaMEDMEM void testDataArraySort1(); void testPartitionBySpreadZone1(); void testGiveCellsWithType1(); + void testBuildSlice3D2(); + void testComputeTupleIdsToSelectFromCellIds1(); }; } diff --git a/src/MEDCoupling_Swig/MEDCouplingBasicsTest.py b/src/MEDCoupling_Swig/MEDCouplingBasicsTest.py index a853cf7f6..e43b76719 100644 --- a/src/MEDCoupling_Swig/MEDCouplingBasicsTest.py +++ b/src/MEDCoupling_Swig/MEDCouplingBasicsTest.py @@ -9980,6 +9980,46 @@ class MEDCouplingBasicsTest(unittest.TestCase): pass pass + def testBuildSlice3D2(self): + mesh3D,mesh2D=MEDCouplingDataForTest.build3DExtrudedUMesh_1(); + vec1=[-0.07,1.,0.07] + origin1=[1.524,1.4552,1.74768] + slice1,ids=mesh3D.buildSlice3D(origin1,vec1,1e-10); + f=MEDCouplingFieldDouble(ON_CELLS,ONE_TIME) + f.setTime(4.5,6,7) ; f.setMesh(mesh3D) + arr=DataArrayDouble(mesh3D.getNumberOfCells(),2) + arr.rearrange(1) ; arr.iota(2.) ; arr.rearrange(2) + f.setArray(arr) + f.checkCoherency() + expected1=DataArrayInt([1,3,4,7,9,10,13,15,16]) + self.assertTrue(expected1.isEqual(ids)) + arr2=arr[expected1] + # + f2=f.extractSlice3D(origin1,vec1,1e-10) + self.assertTrue(f2.getArray().isEqual(arr2,1e-12)); + self.assertTrue(slice1.isEqual(f2.getMesh(),1e-12)) + self.assertEqual(6,f2.getTime()[1]) ; self.assertEqual(7,f2.getTime()[2]) + self.assertAlmostEqual(4.5,f2.getTime()[0],12); + pass + + def testComputeTupleIdsToSelectFromCellIds1(self): + m=MEDCouplingDataForTest.build2DTargetMesh_3() + f=MEDCouplingFieldDouble.New(ON_GAUSS_NE,NO_TIME); + f.setMesh(m); + arr=DataArrayDouble(52,2) ; arr.rearrange(1) ; arr.iota(7.) ; arr.rearrange(2) + f.setArray(arr) + # + f2=f.buildSubPart([1,5,9]) + f2.checkCoherency() + cI=m.computeNbOfNodesPerCell() + cI.computeOffsets2() + sel=DataArrayInt([1,5,9]) + res=sel.buildExplicitArrByRanges(cI) + arr2=arr[res] + self.assertTrue(arr2.isEqual(DataArrayDouble([13,14,15,16,17,18,19,20,59,60,61,62,63,64,95,96,97,98,99,100,101,102,103,104,105,106,107,108,109,110],15,2),1e-12)) + self.assertTrue(arr2.isEqual(f2.getArray(),1e-12)) + pass + def setUp(self): pass pass diff --git a/src/MEDCoupling_Swig/MEDCouplingCommon.i b/src/MEDCoupling_Swig/MEDCouplingCommon.i index 455e25331..90b155fac 100644 --- a/src/MEDCoupling_Swig/MEDCouplingCommon.i +++ b/src/MEDCoupling_Swig/MEDCouplingCommon.i @@ -75,6 +75,7 @@ using namespace INTERP_KERNEL; %newobject ParaMEDMEM::MEDCouplingFieldDiscretization::getOffsetArr; %newobject ParaMEDMEM::MEDCouplingField::buildMeasureField; %newobject ParaMEDMEM::MEDCouplingField::getLocalizationOfDiscr; +%newobject ParaMEDMEM::MEDCouplingField::computeTupleIdsToSelectFromCellIds; %newobject ParaMEDMEM::MEDCouplingFieldDouble::New; %newobject ParaMEDMEM::MEDCouplingFieldDouble::getArray; %newobject ParaMEDMEM::MEDCouplingFieldDouble::getEndArray; @@ -90,6 +91,7 @@ using namespace INTERP_KERNEL; %newobject ParaMEDMEM::MEDCouplingFieldDouble::magnitude; %newobject ParaMEDMEM::MEDCouplingFieldDouble::maxPerTuple; %newobject ParaMEDMEM::MEDCouplingFieldDouble::keepSelectedComponents; +%newobject ParaMEDMEM::MEDCouplingFieldDouble::extractSlice3D; %newobject ParaMEDMEM::MEDCouplingFieldDouble::DotFields; %newobject ParaMEDMEM::MEDCouplingFieldDouble::dot; %newobject ParaMEDMEM::MEDCouplingFieldDouble::CrossProductFields; @@ -255,6 +257,7 @@ using namespace INTERP_KERNEL; %newobject ParaMEDMEM::MEDCouplingUMesh::buildDescendingConnectivity2; %newobject ParaMEDMEM::MEDCouplingUMesh::buildExtrudedMesh; %newobject ParaMEDMEM::MEDCouplingUMesh::buildSpreadZonesWithPoly; +%newobject ParaMEDMEM::MEDCouplingUMesh::computeNbOfNodesPerCell; %newobject ParaMEDMEM::MEDCouplingUMesh::MergeUMeshes; %newobject ParaMEDMEM::MEDCouplingUMesh::MergeUMeshesOnSameCoords; %newobject ParaMEDMEM::MEDCouplingUMesh::ComputeSpreadZoneGradually; @@ -1250,6 +1253,7 @@ namespace ParaMEDMEM DataArrayInt *convertCellArrayPerGeoType(const DataArrayInt *da) const throw(INTERP_KERNEL::Exception); DataArrayInt *computeFetchedNodeIds() const throw(INTERP_KERNEL::Exception); DataArrayInt *zipConnectivityTraducer(int compType) throw(INTERP_KERNEL::Exception); + DataArrayInt *computeNbOfNodesPerCell() const throw(INTERP_KERNEL::Exception); MEDCouplingUMesh *buildDescendingConnectivity(DataArrayInt *desc, DataArrayInt *descIndx, DataArrayInt *revDesc, DataArrayInt *revDescIndx) const throw(INTERP_KERNEL::Exception); MEDCouplingUMesh *buildDescendingConnectivity2(DataArrayInt *desc, DataArrayInt *descIndx, DataArrayInt *revDesc, DataArrayInt *revDescIndx) const throw(INTERP_KERNEL::Exception); void orientCorrectlyPolyhedrons() throw(INTERP_KERNEL::Exception); @@ -6106,6 +6110,33 @@ namespace ParaMEDMEM return res; } + DataArrayInt *computeTupleIdsToSelectFromCellIds(PyObject *li) const + { + int sw; + int pos1; + std::vector pos2; + DataArrayInt *pos3=0; + DataArrayIntTuple *pos4=0; + convertObjToPossibleCpp1(li,sw,pos1,pos2,pos3,pos4); + switch(sw) + { + case 1: + { + return self->computeTupleIdsToSelectFromCellIds(&pos1,&pos1+1); + } + case 2: + { + return self->computeTupleIdsToSelectFromCellIds(&pos2[0],&pos2[0]+pos2.size()); + } + case 3: + { + return self->computeTupleIdsToSelectFromCellIds(pos3->begin(),pos3->end()); + } + default: + throw INTERP_KERNEL::Exception("MEDCouplingField::computeTupleIdsToSelectFromCellIds : unexpected input array type recognized !"); + } + } + void setGaussLocalizationOnCells(PyObject *li, const std::vector& refCoo, const std::vector& gsCoo, const std::vector& wg) throw(INTERP_KERNEL::Exception) { @@ -6553,6 +6584,22 @@ namespace ParaMEDMEM self->setSelectedComponents(f,tmp); } + MEDCouplingFieldDouble *extractSlice3D(PyObject *origin, PyObject *vec, double eps) const throw(INTERP_KERNEL::Exception) + { + double val,val2; + DataArrayDouble *a,*a2; + DataArrayDoubleTuple *aa,*aa2; + std::vector bb,bb2; + int sw; + int spaceDim=3; + const char msg[]="Python wrap of MEDCouplingFieldDouble::extractSlice3D : 1st paramater for origin."; + const char msg2[]="Python wrap of MEDCouplingFieldDouble::extractSlice3D : 2nd paramater for vector."; + const double *orig=convertObjToPossibleCpp5_Safe(origin,sw,val,a,aa,bb,msg,1,spaceDim,true); + const double *vect=convertObjToPossibleCpp5_Safe(vec,sw,val2,a2,aa2,bb2,msg2,1,spaceDim,true); + // + return self->extractSlice3D(orig,vect,eps); + } + PyObject *___iadd___(PyObject *trueSelf, const MEDCouplingFieldDouble& other) throw(INTERP_KERNEL::Exception) { *self+=other; -- 2.39.2