From 03f22ba168cadcc2bf77086517b0fd82140aa269 Mon Sep 17 00:00:00 2001 From: ageay Date: Mon, 27 Sep 2010 05:50:48 +0000 Subject: [PATCH] Addition of functionalities of extraction of a fielddouble Modification of semantic of returned DataArray of zipCoordsTraducer. --- src/MEDCoupling/MEDCouplingCMesh.cxx | 12 ++ src/MEDCoupling/MEDCouplingCMesh.hxx | 2 + src/MEDCoupling/MEDCouplingExtrudedMesh.cxx | 12 ++ src/MEDCoupling/MEDCouplingExtrudedMesh.hxx | 2 + src/MEDCoupling/MEDCouplingField.cxx | 2 +- .../MEDCouplingFieldDiscretization.cxx | 77 +++++---- .../MEDCouplingFieldDiscretization.hxx | 21 ++- src/MEDCoupling/MEDCouplingFieldDouble.cxx | 63 ++++++++ src/MEDCoupling/MEDCouplingFieldDouble.hxx | 3 + src/MEDCoupling/MEDCouplingMemArray.cxx | 74 +++++++++ src/MEDCoupling/MEDCouplingMemArray.hxx | 4 + src/MEDCoupling/MEDCouplingMesh.hxx | 2 + src/MEDCoupling/MEDCouplingPointSet.cxx | 27 ++++ src/MEDCoupling/MEDCouplingPointSet.hxx | 2 + src/MEDCoupling/MEDCouplingUMesh.cxx | 23 ++- .../Test/MEDCouplingBasicsTest.hxx | 4 + .../Test/MEDCouplingBasicsTest1.cxx | 4 +- .../Test/MEDCouplingBasicsTest2.cxx | 146 ++++++++++++++++++ src/MEDCoupling_Swig/MEDCouplingBasicsTest.py | 134 +++++++++++++++- src/MEDCoupling_Swig/libMEDCoupling_Swig.i | 27 ++++ 20 files changed, 590 insertions(+), 51 deletions(-) diff --git a/src/MEDCoupling/MEDCouplingCMesh.cxx b/src/MEDCoupling/MEDCouplingCMesh.cxx index 45f1e8e2e..a69362f47 100644 --- a/src/MEDCoupling/MEDCouplingCMesh.cxx +++ b/src/MEDCoupling/MEDCouplingCMesh.cxx @@ -403,6 +403,18 @@ void MEDCouplingCMesh::setCoords(DataArrayDouble *coordsX, DataArrayDouble *coor declareAsNew(); } +MEDCouplingMesh *MEDCouplingCMesh::buildPart(const int *start, const int *end) const +{ + //not implemented yet ! + return 0; +} + +MEDCouplingMesh *MEDCouplingCMesh::buildPartAndReduceNodes(const int *start, const int *end, DataArrayInt*& arr) const +{ + //not implemented yet ! + return 0; +} + void MEDCouplingCMesh::getBoundingBox(double *bbox) const { //not implemented yet ! diff --git a/src/MEDCoupling/MEDCouplingCMesh.hxx b/src/MEDCoupling/MEDCouplingCMesh.hxx index 9d35de6d8..754f5b54d 100644 --- a/src/MEDCoupling/MEDCouplingCMesh.hxx +++ b/src/MEDCoupling/MEDCouplingCMesh.hxx @@ -59,6 +59,8 @@ namespace ParaMEDMEM DataArrayDouble *coordsY=0, DataArrayDouble *coordsZ=0); // tools + MEDCouplingMesh *buildPart(const int *start, const int *end) const; + MEDCouplingMesh *buildPartAndReduceNodes(const int *start, const int *end, DataArrayInt*& arr) const; void getBoundingBox(double *bbox) const; MEDCouplingFieldDouble *getMeasureField(bool isAbs) const; MEDCouplingFieldDouble *getMeasureFieldOnNode(bool isAbs) const; diff --git a/src/MEDCoupling/MEDCouplingExtrudedMesh.cxx b/src/MEDCoupling/MEDCouplingExtrudedMesh.cxx index a06532e50..a95d38e4e 100644 --- a/src/MEDCoupling/MEDCouplingExtrudedMesh.cxx +++ b/src/MEDCoupling/MEDCouplingExtrudedMesh.cxx @@ -540,6 +540,18 @@ void MEDCouplingExtrudedMesh::translate(const double *vector) _mesh1D->translate(vector); } +MEDCouplingMesh *MEDCouplingExtrudedMesh::buildPart(const int *start, const int *end) const +{ + // not implemented yet ! + return 0; +} + +MEDCouplingMesh *MEDCouplingExtrudedMesh::buildPartAndReduceNodes(const int *start, const int *end, DataArrayInt*& arr) const +{ + // not implemented yet ! + return 0; +} + MEDCouplingMesh *MEDCouplingExtrudedMesh::mergeMyselfWith(const MEDCouplingMesh *other) const { // not implemented yet ! diff --git a/src/MEDCoupling/MEDCouplingExtrudedMesh.hxx b/src/MEDCoupling/MEDCouplingExtrudedMesh.hxx index 1481b58d2..916beae9b 100644 --- a/src/MEDCoupling/MEDCouplingExtrudedMesh.hxx +++ b/src/MEDCoupling/MEDCouplingExtrudedMesh.hxx @@ -72,6 +72,8 @@ namespace ParaMEDMEM MEDCouplingUMesh *&m1r, MEDCouplingUMesh *&m2r, double *v) throw(INTERP_KERNEL::Exception); void rotate(const double *center, const double *vector, double angle); void translate(const double *vector); + MEDCouplingMesh *buildPart(const int *start, const int *end) const; + MEDCouplingMesh *buildPartAndReduceNodes(const int *start, const int *end, DataArrayInt*& arr) const; MEDCouplingMesh *mergeMyselfWith(const MEDCouplingMesh *other) const; DataArrayDouble *getCoordinatesAndOwner() const; DataArrayDouble *getBarycenterAndOwner() const; diff --git a/src/MEDCoupling/MEDCouplingField.cxx b/src/MEDCoupling/MEDCouplingField.cxx index 6425a21a1..8982d4760 100644 --- a/src/MEDCoupling/MEDCouplingField.cxx +++ b/src/MEDCoupling/MEDCouplingField.cxx @@ -251,5 +251,5 @@ MEDCouplingField::MEDCouplingField(const MEDCouplingField& other):_name(other._n */ MEDCouplingMesh *MEDCouplingField::buildSubMeshData(const int *start, const int *end, DataArrayInt *&di) const { - return _type->buildSubMeshData(start,end,_mesh,di); + return _type->buildSubMeshData(_mesh,start,end,di); } diff --git a/src/MEDCoupling/MEDCouplingFieldDiscretization.cxx b/src/MEDCoupling/MEDCouplingFieldDiscretization.cxx index f6027769c..77a958052 100644 --- a/src/MEDCoupling/MEDCouplingFieldDiscretization.cxx +++ b/src/MEDCoupling/MEDCouplingFieldDiscretization.cxx @@ -18,7 +18,6 @@ // #include "MEDCouplingFieldDiscretization.hxx" -#include "MEDCouplingPointSet.hxx" #include "MEDCouplingCMesh.hxx" #include "MEDCouplingFieldDouble.hxx" #include "CellModel.hxx" @@ -28,7 +27,7 @@ #include #include #include -#include +#include using namespace ParaMEDMEM; @@ -300,6 +299,14 @@ DataArrayDouble *MEDCouplingFieldDiscretizationP0::getLocalizationOfDiscValues(c return mesh->getBarycenterAndOwner(); } +void MEDCouplingFieldDiscretizationP0::computeMeshRestrictionFromTupleIds(const MEDCouplingMesh *mesh, const int *partBg, const int *partEnd, + DataArrayInt *&cellRest) +{ + cellRest=DataArrayInt::New(); + cellRest->alloc(std::distance(partBg,partEnd),1); + std::copy(partBg,partEnd,cellRest->getPointer()); +} + void MEDCouplingFieldDiscretizationP0::checkCompatibilityWithNature(NatureOfField nat) const throw(INTERP_KERNEL::Exception) { } @@ -349,9 +356,9 @@ void MEDCouplingFieldDiscretizationP0::renumberValuesOnNodes(const int *, DataAr * @ 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' */ -MEDCouplingMesh *MEDCouplingFieldDiscretizationP0::buildSubMeshData(const int *start, const int *end, const MEDCouplingMesh *mesh, DataArrayInt *&di) const +MEDCouplingMesh *MEDCouplingFieldDiscretizationP0::buildSubMeshData(const MEDCouplingMesh *mesh, const int *start, const int *end, DataArrayInt *&di) const { - MEDCouplingPointSet* ret=((const MEDCouplingPointSet *) mesh)->buildPartOfMySelf(start,end,false); + MEDCouplingMesh *ret=mesh->buildPart(start,end); di=DataArrayInt::New(); di->alloc(std::distance(start,end),1); int *pt=di->getPointer(); @@ -398,6 +405,29 @@ DataArrayDouble *MEDCouplingFieldDiscretizationP1::getLocalizationOfDiscValues(c return mesh->getCoordinatesAndOwner(); } +void MEDCouplingFieldDiscretizationP1::computeMeshRestrictionFromTupleIds(const MEDCouplingMesh *mesh, const int *partBg, const int *partEnd, + DataArrayInt *&cellRest) +{ + std::vector crest; + std::vector conn; + std::set p(partBg,partEnd); + int nbOfCells=mesh->getNumberOfCells(); + for(int i=0;i conn; + mesh->getNodeIdsOfCell(i,conn); + bool cont=true; + for(std::vector::const_iterator iter=conn.begin();iter!=conn.end() && cont;iter++) + if(p.find(*iter)==p.end()) + cont=false; + if(cont) + crest.push_back(i); + } + cellRest=DataArrayInt::New(); + cellRest->alloc(crest.size(),1); + std::copy(crest.begin(),crest.end(),cellRest->getPointer()); +} + void MEDCouplingFieldDiscretizationP1::checkCompatibilityWithNature(NatureOfField nat) const throw(INTERP_KERNEL::Exception) { if(nat!=ConservativeVolumic) @@ -493,33 +523,14 @@ void MEDCouplingFieldDiscretizationP1::renumberValuesOnNodes(const int *old2NewP arrCpy->decrRef(); } -/*! - * This method invert array 'di' that is a conversion map from Old to New node numbering to New to Old node numbering. - */ -DataArrayInt *MEDCouplingFieldDiscretizationP1::invertArrayO2N2N2O(const MEDCouplingMesh *mesh, const DataArrayInt *di) -{ - DataArrayInt *ret=DataArrayInt::New(); - ret->alloc(mesh->getNumberOfNodes(),1); - int nbOfOldNodes=di->getNumberOfTuples(); - const int *old2New=di->getConstPointer(); - int *pt=ret->getPointer(); - for(int i=0;i!=nbOfOldNodes;i++) - if(old2New[i]!=-1) - pt[old2New[i]]=i; - 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 nodes ids) in mesh 'mesh' of entity in returned submesh. * Example : The first node id of returned mesh has the (*di)[0] id in 'mesh' */ -MEDCouplingMesh *MEDCouplingFieldDiscretizationP1::buildSubMeshData(const int *start, const int *end, const MEDCouplingMesh *mesh, DataArrayInt *&di) const +MEDCouplingMesh *MEDCouplingFieldDiscretizationP1::buildSubMeshData(const MEDCouplingMesh *mesh, const int *start, const int *end, DataArrayInt *&di) const { - MEDCouplingPointSet* ret=((const MEDCouplingPointSet *) mesh)->buildPartOfMySelf(start,end,true); - DataArrayInt *diInv=ret->zipCoordsTraducer(); - di=invertArrayO2N2N2O(ret,diInv); - diInv->decrRef(); + MEDCouplingMesh *ret=mesh->buildPartAndReduceNodes(start,end,di); return ret; } @@ -682,6 +693,12 @@ DataArrayDouble *MEDCouplingFieldDiscretizationGauss::getLocalizationOfDiscValue throw INTERP_KERNEL::Exception("Not implemented yet !"); } +void MEDCouplingFieldDiscretizationGauss::computeMeshRestrictionFromTupleIds(const MEDCouplingMesh *mesh, const int *partBg, const int *partEnd, + DataArrayInt *&cellRest) +{ + throw INTERP_KERNEL::Exception("Not implemented yet !"); +} + /*! * Empty : not a bug */ @@ -808,7 +825,7 @@ void MEDCouplingFieldDiscretizationGauss::getValueOnPos(const DataArrayDouble *a throw INTERP_KERNEL::Exception("getValueOnPos(i,j,k) : Not applyable for Gauss points !"); } -MEDCouplingMesh *MEDCouplingFieldDiscretizationGauss::buildSubMeshData(const int *start, const int *end, const MEDCouplingMesh *mesh, DataArrayInt *&di) const +MEDCouplingMesh *MEDCouplingFieldDiscretizationGauss::buildSubMeshData(const MEDCouplingMesh *mesh, const int *start, const int *end, DataArrayInt *&di) const { throw INTERP_KERNEL::Exception("Not implemented yet !"); } @@ -1054,6 +1071,12 @@ DataArrayDouble *MEDCouplingFieldDiscretizationGaussNE::getLocalizationOfDiscVal throw INTERP_KERNEL::Exception("Not implemented yet !"); } +void MEDCouplingFieldDiscretizationGaussNE::computeMeshRestrictionFromTupleIds(const MEDCouplingMesh *mesh, const int *partBg, const int *partEnd, + DataArrayInt *&cellRest) +{ + throw INTERP_KERNEL::Exception("Not implemented yet !"); +} + void MEDCouplingFieldDiscretizationGaussNE::checkCompatibilityWithNature(NatureOfField nat) const throw(INTERP_KERNEL::Exception) { } @@ -1096,7 +1119,7 @@ void MEDCouplingFieldDiscretizationGaussNE::getValueOnPos(const DataArrayDouble throw INTERP_KERNEL::Exception("getValueOnPos(i,j,k) : Not applyable for Gauss points !"); } -MEDCouplingMesh *MEDCouplingFieldDiscretizationGaussNE::buildSubMeshData(const int *start, const int *end, const MEDCouplingMesh *mesh, DataArrayInt *&di) const +MEDCouplingMesh *MEDCouplingFieldDiscretizationGaussNE::buildSubMeshData(const MEDCouplingMesh *mesh, const int *start, const int *end, DataArrayInt *&di) const { throw INTERP_KERNEL::Exception("Not implemented yet !"); } diff --git a/src/MEDCoupling/MEDCouplingFieldDiscretization.hxx b/src/MEDCoupling/MEDCouplingFieldDiscretization.hxx index 5e2ab13fc..89d3e7c0f 100644 --- a/src/MEDCoupling/MEDCouplingFieldDiscretization.hxx +++ b/src/MEDCoupling/MEDCouplingFieldDiscretization.hxx @@ -53,6 +53,8 @@ namespace ParaMEDMEM virtual void normL2(const MEDCouplingMesh *mesh, const DataArrayDouble *arr, double *res) const throw(INTERP_KERNEL::Exception); virtual void integral(const MEDCouplingMesh *mesh, const DataArrayDouble *arr, bool isWAbs, double *res) const throw(INTERP_KERNEL::Exception); virtual DataArrayDouble *getLocalizationOfDiscValues(const MEDCouplingMesh *mesh) const = 0; + virtual void computeMeshRestrictionFromTupleIds(const MEDCouplingMesh *mesh, const int *partBg, const int *partEnd, + DataArrayInt *&cellRest) = 0; virtual void checkCompatibilityWithNature(NatureOfField nat) const throw(INTERP_KERNEL::Exception) = 0; virtual void renumberCells(const int *old2NewBg, bool check) throw(INTERP_KERNEL::Exception); virtual void renumberArraysForCell(const MEDCouplingMesh *mesh, const std::vector& arrays, @@ -62,7 +64,7 @@ namespace ParaMEDMEM virtual MEDCouplingFieldDouble *getWeightingField(const MEDCouplingMesh *mesh, bool isAbs) const = 0; 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 MEDCouplingMesh *buildSubMeshData(const int *start, const int *end, const MEDCouplingMesh *mesh, DataArrayInt *&di) const = 0; + virtual MEDCouplingMesh *buildSubMeshData(const MEDCouplingMesh *mesh, const int *start, const int *end, DataArrayInt *&di) const = 0; virtual void renumberValuesOnNodes(const int *old2New, DataArrayDouble *arr) const = 0; virtual void getSerializationIntArray(DataArrayInt *& arr) const; virtual void getTinySerializationIntInformation(std::vector& tinyInfo) const; @@ -100,12 +102,14 @@ namespace ParaMEDMEM const int *old2NewBg, bool check) throw(INTERP_KERNEL::Exception); DataArrayDouble *getLocalizationOfDiscValues(const MEDCouplingMesh *mesh) const; void checkCompatibilityWithNature(NatureOfField nat) const throw(INTERP_KERNEL::Exception); + void computeMeshRestrictionFromTupleIds(const MEDCouplingMesh *mesh, const int *partBg, const int *partEnd, + DataArrayInt *&cellRest); void checkCoherencyBetween(const MEDCouplingMesh *mesh, const DataArrayDouble *da) const throw(INTERP_KERNEL::Exception); MEDCouplingFieldDouble *getWeightingField(const MEDCouplingMesh *mesh, bool isAbs) const; void getValueOn(const DataArrayDouble *arr, const MEDCouplingMesh *mesh, const double *loc, double *res) const; void getValueOnPos(const DataArrayDouble *arr, const MEDCouplingMesh *mesh, int i, int j, int k, double *res) const; void renumberValuesOnNodes(const int *old2New, DataArrayDouble *arr) const; - MEDCouplingMesh *buildSubMeshData(const int *start, const int *end, const MEDCouplingMesh *mesh, DataArrayInt *&di) const; + MEDCouplingMesh *buildSubMeshData(const MEDCouplingMesh *mesh, const int *start, const int *end, DataArrayInt *&di) const; public: static const char REPR[]; static const TypeOfField TYPE; @@ -122,14 +126,15 @@ namespace ParaMEDMEM void renumberArraysForCell(const MEDCouplingMesh *mesh, const std::vector& arrays, const int *old2NewBg, bool check) throw(INTERP_KERNEL::Exception); DataArrayDouble *getLocalizationOfDiscValues(const MEDCouplingMesh *mesh) const; + void computeMeshRestrictionFromTupleIds(const MEDCouplingMesh *mesh, const int *partBg, const int *partEnd, + DataArrayInt *&cellRest); void checkCompatibilityWithNature(NatureOfField nat) const throw(INTERP_KERNEL::Exception); void checkCoherencyBetween(const MEDCouplingMesh *mesh, const DataArrayDouble *da) const throw(INTERP_KERNEL::Exception); MEDCouplingFieldDouble *getWeightingField(const MEDCouplingMesh *mesh, bool isAbs) const; void getValueOn(const DataArrayDouble *arr, const MEDCouplingMesh *mesh, const double *loc, double *res) const; void getValueOnPos(const DataArrayDouble *arr, const MEDCouplingMesh *mesh, int i, int j, int k, double *res) const; - MEDCouplingMesh *buildSubMeshData(const int *start, const int *end, const MEDCouplingMesh *mesh, DataArrayInt *&di) const; + MEDCouplingMesh *buildSubMeshData(const MEDCouplingMesh *mesh, const int *start, const int *end, DataArrayInt *&di) const; void renumberValuesOnNodes(const int *old2New, DataArrayDouble *arr) const; - static DataArrayInt *invertArrayO2N2N2O(const MEDCouplingMesh *mesh, const DataArrayInt *di); public: static const char REPR[]; static const TypeOfField TYPE; @@ -167,6 +172,8 @@ namespace ParaMEDMEM void renumberArraysForCell(const MEDCouplingMesh *mesh, const std::vector& arrays, const int *old2NewBg, bool check) throw(INTERP_KERNEL::Exception); DataArrayDouble *getLocalizationOfDiscValues(const MEDCouplingMesh *mesh) const; + void computeMeshRestrictionFromTupleIds(const MEDCouplingMesh *mesh, const int *partBg, const int *partEnd, + DataArrayInt *&cellRest); void checkCompatibilityWithNature(NatureOfField nat) const throw(INTERP_KERNEL::Exception); void getTinySerializationIntInformation(std::vector& tinyInfo) const; void getTinySerializationDbleInformation(std::vector& tinyInfo) const; @@ -178,7 +185,7 @@ namespace ParaMEDMEM MEDCouplingFieldDouble *getWeightingField(const MEDCouplingMesh *mesh, bool isAbs) const; void getValueOn(const DataArrayDouble *arr, const MEDCouplingMesh *mesh, const double *loc, double *res) const; void getValueOnPos(const DataArrayDouble *arr, const MEDCouplingMesh *mesh, int i, int j, int k, double *res) const; - MEDCouplingMesh *buildSubMeshData(const int *start, const int *end, const MEDCouplingMesh *mesh, DataArrayInt *&di) const; + MEDCouplingMesh *buildSubMeshData(const MEDCouplingMesh *mesh, const int *start, const int *end, DataArrayInt *&di) const; void renumberValuesOnNodes(const int *old2New, DataArrayDouble *arr) const; void setGaussLocalizationOnType(const MEDCouplingMesh *m, INTERP_KERNEL::NormalizedCellType type, const std::vector& refCoo, const std::vector& gsCoo, const std::vector& wg) throw(INTERP_KERNEL::Exception); @@ -218,13 +225,15 @@ namespace ParaMEDMEM void renumberArraysForCell(const MEDCouplingMesh *mesh, const std::vector& arrays, const int *old2NewBg, bool check) throw(INTERP_KERNEL::Exception); DataArrayDouble *getLocalizationOfDiscValues(const MEDCouplingMesh *mesh) const; + void computeMeshRestrictionFromTupleIds(const MEDCouplingMesh *mesh, const int *partBg, const int *partEnd, + DataArrayInt *&cellRest); void checkCompatibilityWithNature(NatureOfField nat) const throw(INTERP_KERNEL::Exception); double getIJK(const MEDCouplingMesh *mesh, const DataArrayDouble *da, int cellId, int nodeIdInCell, int compoId) const throw(INTERP_KERNEL::Exception); void checkCoherencyBetween(const MEDCouplingMesh *mesh, const DataArrayDouble *da) const throw(INTERP_KERNEL::Exception); MEDCouplingFieldDouble *getWeightingField(const MEDCouplingMesh *mesh, bool isAbs) const; void getValueOn(const DataArrayDouble *arr, const MEDCouplingMesh *mesh, const double *loc, double *res) const; void getValueOnPos(const DataArrayDouble *arr, const MEDCouplingMesh *mesh, int i, int j, int k, double *res) const; - MEDCouplingMesh *buildSubMeshData(const int *start, const int *end, const MEDCouplingMesh *mesh, DataArrayInt *&di) const; + MEDCouplingMesh *buildSubMeshData(const MEDCouplingMesh *mesh, const int *start, const int *end, DataArrayInt *&di) const; void renumberValuesOnNodes(const int *old2New, DataArrayDouble *arr) const; protected: MEDCouplingFieldDiscretizationGaussNE(const MEDCouplingFieldDiscretizationGaussNE& other); diff --git a/src/MEDCoupling/MEDCouplingFieldDouble.cxx b/src/MEDCoupling/MEDCouplingFieldDouble.cxx index 80b38b919..e7f431fdf 100644 --- a/src/MEDCoupling/MEDCouplingFieldDouble.cxx +++ b/src/MEDCoupling/MEDCouplingFieldDouble.cxx @@ -250,6 +250,69 @@ void MEDCouplingFieldDouble::renumberNodesWithoutMesh(const int *old2NewBg) thro _type->renumberValuesOnNodes(old2NewBg,*iter); } +/*! + * This method makes the assumption that the default array is set. If not an exception will be thrown. + * This method is usable only if the default array has exactly one component. If not an exception will be thrown too. + * This method returns all tuples ids that fit the range [vmin,vmax]. + * The caller has the responsability of the returned DataArrayInt. + */ +DataArrayInt *MEDCouplingFieldDouble::getIdsInRange(double vmin, double vmax) const throw(INTERP_KERNEL::Exception) +{ + if(getArray()==0) + throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::getIdsInRange : no default array set !"); + return getArray()->getIdsInRange(vmin,vmax); +} + +/*! + * Builds a newly created field, that the caller will have the responsability. + * This method makes the assumption that the field is correctly defined when this method is called, no check of this will be done. + * This method returns a restriction of 'this' so that only tuples id specified in 'part' will be contained in returned field. + */ +MEDCouplingFieldDouble *MEDCouplingFieldDouble::buildSubPart(const DataArrayInt *part) const throw(INTERP_KERNEL::Exception) +{ + if(part==0) + throw INTERP_KERNEL::Exception("MEDCouplingFieldDouble::buildSubPart : not empty array must be passed to this method !"); + const int *start=part->getConstPointer(); + const int *end=start+part->getNbOfElems(); + return buildSubPart(start,end); +} + +/*! + * Builds a newly created field, that the caller will have the responsability. + * This method makes the assumption that the field is correctly defined when this method is called, no check of this will be done. + * This method returns a restriction of 'this' so that only tuples id specified in ['partBg';'partEnd') will be contained in returned field. + */ +MEDCouplingFieldDouble *MEDCouplingFieldDouble::buildSubPart(const int *partBg, const int *partEnd) const throw(INTERP_KERNEL::Exception) +{ + DataArrayInt *cellRest; + _type->computeMeshRestrictionFromTupleIds(_mesh,partBg,partEnd,cellRest); + DataArrayInt *arrSelect; + MEDCouplingMesh *m=_type->buildSubMeshData(_mesh,cellRest->getConstPointer(),cellRest->getConstPointer()+cellRest->getNbOfElems(),arrSelect); + if(cellRest) + cellRest->decrRef(); + MEDCouplingFieldDouble *ret=clone(false);//quick shallow copy. + ret->setMesh(m); + m->decrRef(); + std::vector arrays; + _time_discr->getArrays(arrays); + std::vector arrs; + const int *arrSelBg=arrSelect->getConstPointer(); + const int *arrSelEnd=arrSelBg+arrSelect->getNbOfElems(); + for(std::vector::const_iterator iter=arrays.begin();iter!=arrays.end();iter++) + { + DataArrayDouble *arr=0; + if(*iter) + arr=(*iter)->selectByTupleId(arrSelBg,arrSelEnd); + arrs.push_back(arr); + } + ret->_time_discr->setArrays(arrs,0); + for(std::vector::const_iterator iter=arrs.begin();iter!=arrs.end();iter++) + if(*iter) + (*iter)->decrRef(); + arrSelect->decrRef(); + return ret; +} + TypeOfTimeDiscretization MEDCouplingFieldDouble::getTimeDiscretization() const { return _time_discr->getEnum(); diff --git a/src/MEDCoupling/MEDCouplingFieldDouble.hxx b/src/MEDCoupling/MEDCouplingFieldDouble.hxx index a3c40d86c..d29f8c6ca 100644 --- a/src/MEDCoupling/MEDCouplingFieldDouble.hxx +++ b/src/MEDCoupling/MEDCouplingFieldDouble.hxx @@ -43,6 +43,9 @@ namespace ParaMEDMEM void renumberCellsWithoutMesh(const int *old2NewBg, bool check) throw(INTERP_KERNEL::Exception); void renumberNodes(const int *old2NewBg) throw(INTERP_KERNEL::Exception); void renumberNodesWithoutMesh(const int *old2NewBg) throw(INTERP_KERNEL::Exception); + DataArrayInt *getIdsInRange(double vmin, double vmax) const throw(INTERP_KERNEL::Exception); + MEDCouplingFieldDouble *buildSubPart(const DataArrayInt *part) const throw(INTERP_KERNEL::Exception); + MEDCouplingFieldDouble *buildSubPart(const int *partBg, const int *partEnd) const throw(INTERP_KERNEL::Exception); MEDCouplingFieldDouble *clone(bool recDeepCpy) const; MEDCouplingFieldDouble *cloneWithMesh(bool recDeepCpy) const; MEDCouplingFieldDouble *buildNewTimeReprFromThis(TypeOfTimeDiscretization td, bool deepCpy) const; diff --git a/src/MEDCoupling/MEDCouplingMemArray.cxx b/src/MEDCoupling/MEDCouplingMemArray.cxx index 3d2e5c16b..c778fd8c9 100644 --- a/src/MEDCoupling/MEDCouplingMemArray.cxx +++ b/src/MEDCoupling/MEDCouplingMemArray.cxx @@ -156,6 +156,10 @@ DataArrayInt *DataArrayDouble::convertToIntArr() const return ret; } +/*! + * This method does \b not change the number of tuples after this call. + * Only a permutation is done. If a permutation reduction is needed substr, or selectByTupleId should be used. + */ void DataArrayDouble::renumberInPlace(const int *old2New) { int nbTuples=getNumberOfTuples(); @@ -169,6 +173,10 @@ void DataArrayDouble::renumberInPlace(const int *old2New) declareAsNew(); } +/*! + * This method does \b not change the number of tuples after this call. + * Only a permutation is done. If a permutation reduction is needed substr, or selectByTupleId should be used. + */ DataArrayDouble *DataArrayDouble::renumber(const int *old2New) const { int nbTuples=getNumberOfTuples(); @@ -211,6 +219,23 @@ DataArrayDouble *DataArrayDouble::substr(int tupleIdBg, int tupleIdEnd) const th return ret; } +/*! + * This method is a generalization of DataArrayDouble::substr method because a not contigous range can be specified here. + */ +DataArrayDouble *DataArrayDouble::selectByTupleId(const int *start, const int *end) const +{ + DataArrayDouble *ret=DataArrayDouble::New(); + int nbComp=getNumberOfComponents(); + ret->alloc(std::distance(start,end),nbComp); + ret->copyStringInfoFrom(*this); + double *pt=ret->getPointer(); + const double *srcPt=getConstPointer(); + int i=0; + for(const int *w=start;w!=end;w++,i++) + std::copy(srcPt+(*w)*nbComp,srcPt+((*w)+1)*nbComp,pt+i*nbComp); + return ret; +} + void DataArrayDouble::setArrayIn(DataArrayDouble *newArray, DataArrayDouble* &arrayToSet) { if(newArray!=arrayToSet) @@ -301,6 +326,22 @@ double DataArrayDouble::accumulate(int compId) const return ret; } +DataArrayInt *DataArrayDouble::getIdsInRange(double vmin, double vmax) const throw(INTERP_KERNEL::Exception) +{ + if(getNumberOfComponents()!=1) + throw INTERP_KERNEL::Exception("DataArrayDouble::getIdsInRange : the default array must have only one component !"); + const double *cptr=getConstPointer(); + std::vector res; + int nbOfTuples=getNumberOfTuples(); + for(int i=0;i=vmin && *cptr<=vmax) + res.push_back(i); + DataArrayInt *ret=DataArrayInt::New(); + ret->alloc(res.size(),1); + std::copy(res.begin(),res.end(),ret->getPointer()); + return ret; +} + DataArrayDouble *DataArrayDouble::aggregate(const DataArrayDouble *a1, const DataArrayDouble *a2) throw(INTERP_KERNEL::Exception) { int nbOfComp=a1->getNumberOfComponents(); @@ -643,6 +684,22 @@ void DataArrayInt::transformWithIndArr(const int *indArr) pt[i]=indArr[pt[i]]; } +/*! + * This method invert array 'di' that is a conversion map from Old to New node numbering to New to Old node numbering. + */ +DataArrayInt *DataArrayInt::invertArrayO2N2N2O(int newNbOfElem) const +{ + DataArrayInt *ret=DataArrayInt::New(); + ret->alloc(newNbOfElem,1); + int nbOfOldNodes=getNumberOfTuples(); + const int *old2New=getConstPointer(); + int *pt=ret->getPointer(); + for(int i=0;i!=nbOfOldNodes;i++) + if(old2New[i]!=-1) + pt[old2New[i]]=i; + return ret; +} + bool DataArrayInt::isEqual(const DataArrayInt& other) const { if(!areInfoEquals(other)) @@ -742,6 +799,23 @@ DataArrayInt *DataArrayInt::substr(int tupleIdBg, int tupleIdEnd) const throw(IN return ret; } +/*! + * This method is a generalization of DataArrayDouble::substr method because a not contigous range can be specified here. + */ +DataArrayInt *DataArrayInt::selectByTupleId(const int *start, const int *end) const +{ + DataArrayInt *ret=DataArrayInt::New(); + int nbComp=getNumberOfComponents(); + ret->alloc(std::distance(start,end),nbComp); + ret->copyStringInfoFrom(*this); + int *pt=ret->getPointer(); + const int *srcPt=getConstPointer(); + int i=0; + for(const int *w=start;w!=end;w++,i++) + std::copy(srcPt+(*w)*nbComp,srcPt+((*w)+1)*nbComp,pt+i*nbComp); + return ret; +} + void DataArrayInt::reAlloc(int nbOfTuples) { _mem.reAlloc(_info_on_compo.size()*nbOfTuples); diff --git a/src/MEDCoupling/MEDCouplingMemArray.hxx b/src/MEDCoupling/MEDCouplingMemArray.hxx index 640fe8c8f..0405f3da8 100644 --- a/src/MEDCoupling/MEDCouplingMemArray.hxx +++ b/src/MEDCoupling/MEDCouplingMemArray.hxx @@ -128,6 +128,7 @@ namespace ParaMEDMEM MEDCOUPLING_EXPORT void renumberInPlace(const int *old2New); MEDCOUPLING_EXPORT DataArrayDouble *renumber(const int *old2New) const; MEDCOUPLING_EXPORT DataArrayDouble *substr(int tupleIdBg, int tupleIdEnd=-1) const throw(INTERP_KERNEL::Exception); + MEDCOUPLING_EXPORT DataArrayDouble *selectByTupleId(const int *start, const int *end) const; MEDCOUPLING_EXPORT void getTuple(int tupleId, double *res) const { std::copy(_mem.getConstPointerLoc(tupleId*_info_on_compo.size()),_mem.getConstPointerLoc((tupleId+1)*_info_on_compo.size()),res); } MEDCOUPLING_EXPORT double getIJ(int tupleId, int compoId) const { return _mem[tupleId*_info_on_compo.size()+compoId]; } MEDCOUPLING_EXPORT void setIJ(int tupleId, int compoId, double newVal) { _mem[tupleId*_info_on_compo.size()+compoId]=newVal; declareAsNew(); } @@ -142,6 +143,7 @@ namespace ParaMEDMEM MEDCOUPLING_EXPORT double getAverageValue() const throw(INTERP_KERNEL::Exception); MEDCOUPLING_EXPORT void accumulate(double *res) const; MEDCOUPLING_EXPORT double accumulate(int compId) const; + MEDCOUPLING_EXPORT DataArrayInt *getIdsInRange(double vmin, double vmax) const throw(INTERP_KERNEL::Exception); MEDCOUPLING_EXPORT static DataArrayDouble *aggregate(const DataArrayDouble *a1, const DataArrayDouble *a2) throw(INTERP_KERNEL::Exception); MEDCOUPLING_EXPORT static DataArrayDouble *dot(const DataArrayDouble *a1, const DataArrayDouble *a2) throw(INTERP_KERNEL::Exception); MEDCOUPLING_EXPORT static DataArrayDouble *crossProduct(const DataArrayDouble *a1, const DataArrayDouble *a2) throw(INTERP_KERNEL::Exception); @@ -179,6 +181,7 @@ namespace ParaMEDMEM MEDCOUPLING_EXPORT void reprWithoutNameStream(std::ostream& stream) const; MEDCOUPLING_EXPORT void reprZipWithoutNameStream(std::ostream& stream) const; MEDCOUPLING_EXPORT void transformWithIndArr(const int *indArr); + MEDCOUPLING_EXPORT DataArrayInt *invertArrayO2N2N2O(int newNbOfElem) const; //!alloc or useArray should have been called before. MEDCOUPLING_EXPORT void reAlloc(int nbOfTuples); MEDCOUPLING_EXPORT DataArrayDouble *convertToDblArr() const; @@ -186,6 +189,7 @@ namespace ParaMEDMEM MEDCOUPLING_EXPORT DataArrayInt *renumber(const int *old2New) const; MEDCOUPLING_EXPORT bool isIdentity() const; MEDCOUPLING_EXPORT DataArrayInt *substr(int tupleIdBg, int tupleIdEnd=-1) const throw(INTERP_KERNEL::Exception); + MEDCOUPLING_EXPORT DataArrayInt *selectByTupleId(const int *start, const int *end) const; MEDCOUPLING_EXPORT void getTuple(int tupleId, int *res) const { std::copy(_mem.getConstPointerLoc(tupleId*_info_on_compo.size()),_mem.getConstPointerLoc((tupleId+1)*_info_on_compo.size()),res); } MEDCOUPLING_EXPORT int getIJ(int tupleId, int compoId) const { return _mem[tupleId*_info_on_compo.size()+compoId]; } MEDCOUPLING_EXPORT void setIJ(int tupleId, int compoId, int newVal) { _mem[tupleId*_info_on_compo.size()+compoId]=newVal; } diff --git a/src/MEDCoupling/MEDCouplingMesh.hxx b/src/MEDCoupling/MEDCouplingMesh.hxx index bb1d9adfc..7ca13c6f3 100644 --- a/src/MEDCoupling/MEDCouplingMesh.hxx +++ b/src/MEDCoupling/MEDCouplingMesh.hxx @@ -84,6 +84,8 @@ namespace ParaMEDMEM virtual void translate(const double *vector) = 0; virtual void renumberCells(const int *old2NewBg, bool check) throw(INTERP_KERNEL::Exception) = 0; virtual MEDCouplingMesh *mergeMyselfWith(const MEDCouplingMesh *other) const = 0; + virtual MEDCouplingMesh *buildPart(const int *start, const int *end) const = 0; + virtual MEDCouplingMesh *buildPartAndReduceNodes(const int *start, const int *end, DataArrayInt*& arr) const = 0; virtual bool areCompatibleForMerge(const MEDCouplingMesh *other) const; static MEDCouplingMesh *mergeMeshes(const MEDCouplingMesh *mesh1, const MEDCouplingMesh *mesh2); //serialisation-unserialization diff --git a/src/MEDCoupling/MEDCouplingPointSet.cxx b/src/MEDCoupling/MEDCouplingPointSet.cxx index 142b0e112..709ae8ce2 100644 --- a/src/MEDCoupling/MEDCouplingPointSet.cxx +++ b/src/MEDCoupling/MEDCouplingPointSet.cxx @@ -674,6 +674,33 @@ void MEDCouplingPointSet::rotate3DAlg(const double *center, const double *vect, } } +/*! + * This method implements pure virtual method MEDCouplingMesh::buildPart. + * This method build a part of 'this' by simply keeping cells whose ids are in ['start','end'). + * The coords are kept unchanged contrary to pure virtual method MEDCouplingMesh::buildPartAndReduceNodes. + * The returned mesh has to be managed by the caller. + */ +MEDCouplingMesh *MEDCouplingPointSet::buildPart(const int *start, const int *end) const +{ + return buildPartOfMySelf(start,end,false); +} + +/*! + * This method implements pure virtual method MEDCouplingMesh::buildPartAndReduceNodes. + * This method build a part of 'this' by simply keeping cells whose ids are in ['start','end') \b and potentially reduces the nodes set + * behind returned mesh. This cause an overhead but it is more little in memory. + * This method returns an array too. This array allows to the caller to know the mapping between nodeids in 'this' and nodeids in + * returned mesh. This is quite usefull for MEDCouplingFieldDouble on nodes for example... + * The returned mesh has to be managed by the caller. + */ +MEDCouplingMesh *MEDCouplingPointSet::buildPartAndReduceNodes(const int *start, const int *end, DataArrayInt*& arr) const +{ + MEDCouplingPointSet *ret=buildPartOfMySelf(start,end,true); + arr=ret->zipCoordsTraducer(); + return ret; +} + + /*! * 'This' is expected to be of spaceDim==2. Idem for 'center' and 'vect' */ diff --git a/src/MEDCoupling/MEDCouplingPointSet.hxx b/src/MEDCoupling/MEDCouplingPointSet.hxx index ac69fee23..21f56317e 100644 --- a/src/MEDCoupling/MEDCouplingPointSet.hxx +++ b/src/MEDCoupling/MEDCouplingPointSet.hxx @@ -70,6 +70,8 @@ namespace ParaMEDMEM static MEDCouplingPointSet *buildInstanceFromMeshType(MEDCouplingMeshType type); static void rotate2DAlg(const double *center, double angle, int nbNodes, double *coords); static void rotate3DAlg(const double *center, const double *vect, double angle, int nbNodes, double *coords); + MEDCouplingMesh *buildPart(const int *start, const int *end) const; + MEDCouplingMesh *buildPartAndReduceNodes(const int *start, const int *end, DataArrayInt*& arr) const; virtual MEDCouplingPointSet *buildPartOfMySelf(const int *start, const int *end, bool keepCoords) const = 0; virtual MEDCouplingPointSet *buildPartOfMySelfNode(const int *start, const int *end, bool fullyIn) const = 0; virtual MEDCouplingPointSet *buildFacePartOfMySelfNode(const int *start, const int *end, bool fullyIn) const = 0; diff --git a/src/MEDCoupling/MEDCouplingUMesh.cxx b/src/MEDCoupling/MEDCouplingUMesh.cxx index b45ca0568..32bec6151 100644 --- a/src/MEDCoupling/MEDCouplingUMesh.cxx +++ b/src/MEDCoupling/MEDCouplingUMesh.cxx @@ -490,19 +490,17 @@ void MEDCouplingUMesh::convertToPolyTypes(const std::vector& cellIdsToConve } /*! - * Array returned is the correspondance old to new. + * Array returned is the correspondance new to old. * The maximum value stored in returned array is the number of nodes of 'this' minus 1 after call of this method. * The size of returned array is the number of nodes of the old (previous to the call of this method) number of nodes. * -1 values in returned array means that the corresponding old node is no more used. */ DataArrayInt *MEDCouplingUMesh::zipCoordsTraducer() { - DataArrayInt *ret=DataArrayInt::New(); int nbOfNodes=getNumberOfNodes(); int spaceDim=getSpaceDimension(); int *traducer=new int[nbOfNodes]; std::fill(traducer,traducer+nbOfNodes,-1); - ret->useArray(traducer,true,CPP_DEALLOC,nbOfNodes,1); int nbOfCells=getNumberOfCells(); const int *connIndex=_nodal_connec_index->getConstPointer(); int *conn=_nodal_connec->getPointer(); @@ -516,17 +514,14 @@ DataArrayInt *MEDCouplingUMesh::zipCoordsTraducer() for(int j=connIndex[i]+1;j=0) conn[j]=traducer[conn[j]]; - DataArrayDouble *newCoords=DataArrayDouble::New(); - double *newCoordsPtr=new double[newNbOfNodes*spaceDim]; - const double *oldCoordsPtr=_coords->getConstPointer(); - newCoords->useArray(newCoordsPtr,true,CPP_DEALLOC,newNbOfNodes,spaceDim); - int *work=std::find_if(traducer,traducer+nbOfNodes,std::bind2nd(std::not_equal_to(),-1)); - for(;work!=traducer+nbOfNodes;work=std::find_if(work,traducer+nbOfNodes,std::bind2nd(std::not_equal_to(),-1))) - { - newCoordsPtr=std::copy(oldCoordsPtr+spaceDim*(work-traducer),oldCoordsPtr+spaceDim*(work-traducer+1),newCoordsPtr); - work++; - } - newCoords->copyStringInfoFrom(*_coords); + DataArrayInt *ret=DataArrayInt::New(); + ret->alloc(newNbOfNodes,1); + int *retPtr=ret->getPointer(); + for(int i=0;iselectByTupleId(retPtr,retPtr+newNbOfNodes); setCoords(newCoords); newCoords->decrRef(); return ret; diff --git a/src/MEDCoupling/Test/MEDCouplingBasicsTest.hxx b/src/MEDCoupling/Test/MEDCouplingBasicsTest.hxx index ec4ae8989..23deef32e 100644 --- a/src/MEDCoupling/Test/MEDCouplingBasicsTest.hxx +++ b/src/MEDCoupling/Test/MEDCouplingBasicsTest.hxx @@ -109,6 +109,8 @@ namespace ParaMEDMEM CPPUNIT_TEST( testDotCrossProduct1 ); CPPUNIT_TEST( testMinMaxFields1 ); CPPUNIT_TEST( testApplyLin1 ); + CPPUNIT_TEST( testGetIdsInRange1 ); + CPPUNIT_TEST( testBuildSubPart1 ); //MEDCouplingBasicsTestInterp.cxx CPPUNIT_TEST( test2DInterpP0P0_1 ); CPPUNIT_TEST( test2DInterpP0P0PL_1 ); @@ -246,6 +248,8 @@ namespace ParaMEDMEM void testDotCrossProduct1(); void testMinMaxFields1(); void testApplyLin1(); + void testGetIdsInRange1(); + void testBuildSubPart1(); //MEDCouplingBasicsTestInterp.cxx void test2DInterpP0P0_1(); void test2DInterpP0P0PL_1(); diff --git a/src/MEDCoupling/Test/MEDCouplingBasicsTest1.cxx b/src/MEDCoupling/Test/MEDCouplingBasicsTest1.cxx index 420650895..fb96af12d 100644 --- a/src/MEDCoupling/Test/MEDCouplingBasicsTest1.cxx +++ b/src/MEDCoupling/Test/MEDCouplingBasicsTest1.cxx @@ -691,8 +691,8 @@ void MEDCouplingBasicsTest::testZipCoords() MEDCouplingUMesh *subMesh=dynamic_cast(subMeshPtSet); CPPUNIT_ASSERT(subMesh); DataArrayInt *traducer=subMesh->zipCoordsTraducer(); - const int expectedTraducer[9]={0,1,-1,2,3,4,-1,5,6}; - CPPUNIT_ASSERT(std::equal(expectedTraducer,expectedTraducer+9,traducer->getPointer())); + const int expectedTraducer[7]={0,1,3,4,5,7,8}; + CPPUNIT_ASSERT(std::equal(expectedTraducer,expectedTraducer+7,traducer->getPointer())); traducer->decrRef(); CPPUNIT_ASSERT_EQUAL(INTERP_KERNEL::NORM_QUAD4,*subMesh->getAllTypes().begin()); CPPUNIT_ASSERT_EQUAL(2,subMesh->getNumberOfCells()); diff --git a/src/MEDCoupling/Test/MEDCouplingBasicsTest2.cxx b/src/MEDCoupling/Test/MEDCouplingBasicsTest2.cxx index 479a8d25d..eca5b836e 100644 --- a/src/MEDCoupling/Test/MEDCouplingBasicsTest2.cxx +++ b/src/MEDCoupling/Test/MEDCouplingBasicsTest2.cxx @@ -1113,6 +1113,10 @@ void MEDCouplingBasicsTest::testSubstractInPlaceDM1() f1->setArray(array); array->decrRef(); // + CPPUNIT_ASSERT_EQUAL(10,f1->getNumberOfTuples()); + CPPUNIT_ASSERT_EQUAL(2,f1->getNumberOfComponents()); + CPPUNIT_ASSERT_EQUAL(20,f1->getNumberOfValues()); + // const int renum[]={0,2,1,3,4,5,6,8,7,9}; mesh2->renumberCells(renum,false); // @@ -1249,3 +1253,145 @@ void MEDCouplingBasicsTest::testApplyLin1() mesh1->decrRef(); f1->decrRef(); } + +void MEDCouplingBasicsTest::testGetIdsInRange1() +{ + MEDCouplingUMesh *mesh1=build2DTargetMesh_3(); + MEDCouplingFieldDouble *f1=MEDCouplingFieldDouble::New(ON_CELLS,ONE_TIME); + f1->setTime(2.3,5,6); + f1->setMesh(mesh1); + DataArrayDouble *array=DataArrayDouble::New(); + array->alloc(mesh1->getNumberOfCells(),1); + const double arr1[10]={2.,8.,6.,5.,11.,7.,9.,3.,10.,4.}; + std::copy(arr1,arr1+10,array->getPointer()); + f1->setArray(array); + array->decrRef(); + // + f1->checkCoherency(); + DataArrayInt *da=f1->getIdsInRange(2.9,7.1); + CPPUNIT_ASSERT_EQUAL(5,da->getNbOfElems()); + const int expected1[5]={2,3,5,7,9}; + CPPUNIT_ASSERT(std::equal(expected1,expected1+5,da->getConstPointer())); + da->decrRef(); + da=f1->getIdsInRange(8.,12.); + CPPUNIT_ASSERT_EQUAL(4,da->getNbOfElems()); + const int expected2[4]={1,4,6,8}; + CPPUNIT_ASSERT(std::equal(expected2,expected2+4,da->getConstPointer())); + da->decrRef(); + // + f1->decrRef(); + mesh1->decrRef(); +} + +void MEDCouplingBasicsTest::testBuildSubPart1() +{ + MEDCouplingUMesh *mesh1=build2DTargetMesh_1(); + MEDCouplingFieldDouble *f1=MEDCouplingFieldDouble::New(ON_CELLS,ONE_TIME); + f1->setTime(2.3,5,6); + f1->setMesh(mesh1); + DataArrayDouble *array=DataArrayDouble::New(); + array->alloc(mesh1->getNumberOfCells(),2); + const double arr1[10]={3.,103.,4.,104.,5.,105.,6.,106.,7.,107.}; + std::copy(arr1,arr1+10,array->getPointer()); + f1->setArray(array); + array->decrRef(); + // + const int part1[3]={2,1,4}; + MEDCouplingFieldDouble *f2=f1->buildSubPart(part1,part1+3); + CPPUNIT_ASSERT_EQUAL(3,f2->getNumberOfTuples()); + CPPUNIT_ASSERT_EQUAL(2,f2->getNumberOfComponents()); + const double expected1[6]={5.,105.,4.,104.,7.,107.}; + for(int i=0;i<6;i++) + CPPUNIT_ASSERT_DOUBLES_EQUAL(f2->getIJ(0,i),expected1[i],1e-12); + CPPUNIT_ASSERT_EQUAL(3,f2->getMesh()->getNumberOfCells()); + CPPUNIT_ASSERT_EQUAL(6,f2->getMesh()->getNumberOfNodes()); + CPPUNIT_ASSERT_EQUAL(2,f2->getMesh()->getSpaceDimension()); + CPPUNIT_ASSERT_EQUAL(2,f2->getMesh()->getMeshDimension()); + MEDCouplingUMesh *m2C=dynamic_cast(const_cast(f2->getMesh())); + CPPUNIT_ASSERT_EQUAL(13,m2C->getMeshLength()); + const double expected2[12]={0.2, -0.3, 0.7, -0.3, 0.2, 0.2, 0.7, 0.2, 0.2, 0.7, 0.7, 0.7}; + for(int i=0;i<12;i++) + CPPUNIT_ASSERT_DOUBLES_EQUAL(expected2[i],m2C->getCoords()->getIJ(0,i),1.e-12); + const double expected3[13]={3,2,3,1,3,0,2,1,4,4,5,3,2}; + CPPUNIT_ASSERT(std::equal(expected3,expected3+13,m2C->getNodalConnectivity()->getConstPointer())); + const double expected4[4]={0,4,8,13}; + CPPUNIT_ASSERT(std::equal(expected4,expected4+4,m2C->getNodalConnectivityIndex()->getConstPointer())); + f2->decrRef(); + f1->decrRef(); + // Test with field on nodes. + f1=MEDCouplingFieldDouble::New(ON_NODES,ONE_TIME); + f1->setTime(2.3,5,6); + f1->setMesh(mesh1); + array=DataArrayDouble::New(); + array->alloc(mesh1->getNumberOfNodes(),2); + const double arr2[18]={3.,103.,4.,104.,5.,105.,6.,106.,7.,107.,8.,108.,9.,109.,10.,110.,11.,111.}; + std::copy(arr2,arr2+18,array->getPointer()); + f1->setArray(array); + array->decrRef(); + const int part2[4]={1,4,2,5}; + f2=f1->buildSubPart(part2,part2+4); + CPPUNIT_ASSERT_EQUAL(4,f2->getNumberOfTuples()); + CPPUNIT_ASSERT_EQUAL(2,f2->getNumberOfComponents()); + const double expected5[8]={4.,104.,5.,105.,7.,107.,8.,108.}; + for(int i=0;i<8;i++) + CPPUNIT_ASSERT_DOUBLES_EQUAL(f2->getIJ(0,i),expected5[i],1e-12); + CPPUNIT_ASSERT_EQUAL(2,f2->getMesh()->getNumberOfCells()); + CPPUNIT_ASSERT_EQUAL(4,f2->getMesh()->getNumberOfNodes()); + CPPUNIT_ASSERT_EQUAL(2,f2->getMesh()->getSpaceDimension()); + CPPUNIT_ASSERT_EQUAL(2,f2->getMesh()->getMeshDimension()); + m2C=dynamic_cast(const_cast(f2->getMesh())); + CPPUNIT_ASSERT_EQUAL(8,m2C->getMeshLength()); + for(int i=0;i<8;i++)//8 is not an error + CPPUNIT_ASSERT_DOUBLES_EQUAL(expected2[i],m2C->getCoords()->getIJ(0,i),1.e-12); + CPPUNIT_ASSERT(std::equal(expected3,expected3+4,m2C->getNodalConnectivity()->getConstPointer()+4)); + CPPUNIT_ASSERT(std::equal(expected3+4,expected3+8,m2C->getNodalConnectivity()->getConstPointer())); + CPPUNIT_ASSERT(std::equal(expected4,expected4+3,m2C->getNodalConnectivityIndex()->getConstPointer())); + f2->decrRef(); + //idem previous because nodes of cell#4 are not fully present in part3 + const int part3[5]={1,4,2,5,7}; + DataArrayInt *arrr=DataArrayInt::New(); + arrr->alloc(5,1); + std::copy(part3,part3+5,arrr->getPointer()); + f2=f1->buildSubPart(arrr); + arrr->decrRef(); + CPPUNIT_ASSERT_EQUAL(4,f2->getNumberOfTuples()); + CPPUNIT_ASSERT_EQUAL(2,f2->getNumberOfComponents()); + for(int i=0;i<8;i++) + CPPUNIT_ASSERT_DOUBLES_EQUAL(f2->getIJ(0,i),expected5[i],1e-12); + CPPUNIT_ASSERT_EQUAL(2,f2->getMesh()->getNumberOfCells()); + CPPUNIT_ASSERT_EQUAL(4,f2->getMesh()->getNumberOfNodes()); + CPPUNIT_ASSERT_EQUAL(2,f2->getMesh()->getSpaceDimension()); + CPPUNIT_ASSERT_EQUAL(2,f2->getMesh()->getMeshDimension()); + m2C=dynamic_cast(const_cast(f2->getMesh())); + CPPUNIT_ASSERT_EQUAL(8,m2C->getMeshLength()); + for(int i=0;i<8;i++)//8 is not an error + CPPUNIT_ASSERT_DOUBLES_EQUAL(expected2[i],m2C->getCoords()->getIJ(0,i),1.e-12); + CPPUNIT_ASSERT(std::equal(expected3,expected3+4,m2C->getNodalConnectivity()->getConstPointer()+4)); + CPPUNIT_ASSERT(std::equal(expected3+4,expected3+8,m2C->getNodalConnectivity()->getConstPointer())); + CPPUNIT_ASSERT(std::equal(expected4,expected4+3,m2C->getNodalConnectivityIndex()->getConstPointer())); + f2->decrRef(); + // + const int part4[6]={1,4,2,5,7,8}; + f2=f1->buildSubPart(part4,part4+6); + CPPUNIT_ASSERT_EQUAL(6,f2->getNumberOfTuples()); + CPPUNIT_ASSERT_EQUAL(2,f2->getNumberOfComponents()); + const double expected6[12]={4.,104.,5.,105.,7.,107.,8.,108.,10.,110.,11.,111.}; + for(int i=0;i<12;i++) + CPPUNIT_ASSERT_DOUBLES_EQUAL(f2->getIJ(0,i),expected6[i],1e-12); + CPPUNIT_ASSERT_EQUAL(3,f2->getMesh()->getNumberOfCells()); + CPPUNIT_ASSERT_EQUAL(6,f2->getMesh()->getNumberOfNodes()); + CPPUNIT_ASSERT_EQUAL(2,f2->getMesh()->getSpaceDimension()); + CPPUNIT_ASSERT_EQUAL(2,f2->getMesh()->getMeshDimension()); + m2C=dynamic_cast(const_cast(f2->getMesh())); + CPPUNIT_ASSERT_EQUAL(13,m2C->getMeshLength()); + for(int i=0;i<12;i++) + CPPUNIT_ASSERT_DOUBLES_EQUAL(expected2[i],m2C->getCoords()->getIJ(0,i),1.e-12); + CPPUNIT_ASSERT(std::equal(expected3,expected3+4,m2C->getNodalConnectivity()->getConstPointer()+4)); + CPPUNIT_ASSERT(std::equal(expected3+4,expected3+8,m2C->getNodalConnectivity()->getConstPointer())); + CPPUNIT_ASSERT(std::equal(expected3+8,expected3+13,m2C->getNodalConnectivity()->getConstPointer()+8)); + CPPUNIT_ASSERT(std::equal(expected4,expected4+4,m2C->getNodalConnectivityIndex()->getConstPointer())); + f2->decrRef(); + // + f1->decrRef(); + mesh1->decrRef(); +} diff --git a/src/MEDCoupling_Swig/MEDCouplingBasicsTest.py b/src/MEDCoupling_Swig/MEDCouplingBasicsTest.py index 4c001de19..a764cbfcb 100644 --- a/src/MEDCoupling_Swig/MEDCouplingBasicsTest.py +++ b/src/MEDCoupling_Swig/MEDCouplingBasicsTest.py @@ -446,7 +446,7 @@ class MEDCouplingBasicsTest(unittest.TestCase): subMesh=mesh.buildPartOfMySelf(tab1,True); self.assertTrue(isinstance(subMesh,MEDCouplingUMesh)) traducer=subMesh.zipCoordsTraducer(); - expectedTraducer=[0,1,-1,2,3,4,-1,5,6] + expectedTraducer=[0,1,3,4,5,7,8] self.assertEqual(expectedTraducer,traducer.getValues()); self.assertEqual(NORM_QUAD4,subMesh.getAllTypes()[0]); self.assertEqual(2,subMesh.getNumberOfCells()); @@ -2648,6 +2648,10 @@ class MEDCouplingBasicsTest(unittest.TestCase): array.setValues(arr,mesh1.getNumberOfCells(),2); f1.setArray(array); # + self.assertEqual(10,f1.getNumberOfTuples()); + self.assertEqual(2,f1.getNumberOfComponents()); + self.assertEqual(20,f1.getNumberOfValues()); + # renum=[0,2,1,3,4,5,6,8,7,9] mesh2.renumberCells(renum,False); # @@ -2756,6 +2760,134 @@ class MEDCouplingBasicsTest(unittest.TestCase): pass # pass + + def testGetIdsInRange1(self): + mesh1=MEDCouplingDataForTest.build2DTargetMesh_3(); + f1=MEDCouplingFieldDouble.New(ON_CELLS,ONE_TIME); + f1.setTime(2.3,5,6); + f1.setMesh(mesh1); + array=DataArrayDouble.New(); + arr1=[2.,8.,6.,5.,11.,7.,9.,3.,10.,4.] + array.setValues(arr1,mesh1.getNumberOfCells(),1); + f1.setArray(array); + # + f1.checkCoherency(); + da=f1.getIdsInRange(2.9,7.1); + self.failUnlessEqual(5,da.getNbOfElems()); + expected1=[2,3,5,7,9] + self.failUnlessEqual(expected1,da.getValues()); + da=f1.getIdsInRange(8.,12.); + self.failUnlessEqual(4,da.getNbOfElems()); + expected2=[1,4,6,8] + self.failUnlessEqual(expected2,da.getValues()); + # + pass + + def testBuildSubPart1(self): + mesh1=MEDCouplingDataForTest.build2DTargetMesh_1(); + f1=MEDCouplingFieldDouble.New(ON_CELLS,ONE_TIME); + f1.setTime(2.3,5,6); + f1.setMesh(mesh1); + array=DataArrayDouble.New(); + arr1=[3.,103.,4.,104.,5.,105.,6.,106.,7.,107.] + array.setValues(arr1,mesh1.getNumberOfCells(),2); + f1.setArray(array); + # + part1=[2,1,4] + f2=f1.buildSubPart(part1); + self.failUnlessEqual(3,f2.getNumberOfTuples()); + self.failUnlessEqual(2,f2.getNumberOfComponents()); + expected1=[5.,105.,4.,104.,7.,107.] + for i in xrange(6): + self.assertAlmostEqual(f2.getIJ(0,i),expected1[i],12); + pass + self.failUnlessEqual(3,f2.getMesh().getNumberOfCells()); + self.failUnlessEqual(6,f2.getMesh().getNumberOfNodes()); + self.failUnlessEqual(2,f2.getMesh().getSpaceDimension()); + self.failUnlessEqual(2,f2.getMesh().getMeshDimension()); + m2C=f2.getMesh(); + self.failUnlessEqual(13,m2C.getMeshLength()); + expected2=[0.2, -0.3, 0.7, -0.3, 0.2, 0.2, 0.7, 0.2, 0.2, 0.7, 0.7, 0.7] + for i in xrange(12): + self.assertAlmostEqual(expected2[i],m2C.getCoords().getIJ(0,i),12); + pass + expected3=[3,2,3,1,3,0,2,1,4,4,5,3,2] + self.failUnlessEqual(expected3,m2C.getNodalConnectivity().getValues()); + expected4=[0,4,8,13] + self.failUnlessEqual(expected4,m2C.getNodalConnectivityIndex().getValues()); + # Test with field on nodes. + f1=MEDCouplingFieldDouble.New(ON_NODES,ONE_TIME); + f1.setTime(2.3,5,6); + f1.setMesh(mesh1); + array=DataArrayDouble.New(); + arr2=[3.,103.,4.,104.,5.,105.,6.,106.,7.,107.,8.,108.,9.,109.,10.,110.,11.,111.] + array.setValues(arr2,mesh1.getNumberOfNodes(),2); + f1.setArray(array); + part2=[1,4,2,5] + f2=f1.buildSubPart(part2); + self.failUnlessEqual(4,f2.getNumberOfTuples()); + self.failUnlessEqual(2,f2.getNumberOfComponents()); + expected5=[4.,104.,5.,105.,7.,107.,8.,108.] + for i in xrange(8): + self.assertAlmostEqual(f2.getIJ(0,i),expected5[i],12); + pass + self.failUnlessEqual(2,f2.getMesh().getNumberOfCells()); + self.failUnlessEqual(4,f2.getMesh().getNumberOfNodes()); + self.failUnlessEqual(2,f2.getMesh().getSpaceDimension()); + self.failUnlessEqual(2,f2.getMesh().getMeshDimension()); + m2C=f2.getMesh(); + self.failUnlessEqual(8,m2C.getMeshLength()); + for i in xrange(8):#8 is not an error + self.assertAlmostEqual(expected2[i],m2C.getCoords().getIJ(0,i),12); + pass + self.failUnlessEqual(expected3[:4],m2C.getNodalConnectivity().getValues()[4:]); + self.failUnlessEqual(expected3[4:8],m2C.getNodalConnectivity().getValues()[:4]); + self.failUnlessEqual(expected4[:3],m2C.getNodalConnectivityIndex().getValues()); + #idem previous because nodes of cell#4 are not fully present in part3 + part3=[1,4,2,5,7] + arrr=DataArrayInt.New(); + arrr.setValues(part3,5,1); + f2=f1.buildSubPart(arrr); + self.failUnlessEqual(4,f2.getNumberOfTuples()); + self.failUnlessEqual(2,f2.getNumberOfComponents()); + for i in xrange(8): + self.assertAlmostEqual(f2.getIJ(0,i),expected5[i],12); + pass + self.failUnlessEqual(2,f2.getMesh().getNumberOfCells()); + self.failUnlessEqual(4,f2.getMesh().getNumberOfNodes()); + self.failUnlessEqual(2,f2.getMesh().getSpaceDimension()); + self.failUnlessEqual(2,f2.getMesh().getMeshDimension()); + m2C=f2.getMesh(); + self.failUnlessEqual(8,m2C.getMeshLength()); + for i in xrange(8):#8 is not an error + self.assertAlmostEqual(expected2[i],m2C.getCoords().getIJ(0,i),12); + pass + self.failUnlessEqual(expected3[:4],m2C.getNodalConnectivity().getValues()[4:8]); + self.failUnlessEqual(expected3[4:8],m2C.getNodalConnectivity().getValues()[:4]); + self.failUnlessEqual(expected4[:3],m2C.getNodalConnectivityIndex().getValues()); + # + part4=[1,4,2,5,7,8] + f2=f1.buildSubPart(part4); + self.failUnlessEqual(6,f2.getNumberOfTuples()); + self.failUnlessEqual(2,f2.getNumberOfComponents()); + expected6=[4.,104.,5.,105.,7.,107.,8.,108.,10.,110.,11.,111.] + for i in xrange(12): + self.assertAlmostEqual(f2.getIJ(0,i),expected6[i],12); + pass + self.failUnlessEqual(3,f2.getMesh().getNumberOfCells()); + self.failUnlessEqual(6,f2.getMesh().getNumberOfNodes()); + self.failUnlessEqual(2,f2.getMesh().getSpaceDimension()); + self.failUnlessEqual(2,f2.getMesh().getMeshDimension()); + m2C=f2.getMesh(); + self.failUnlessEqual(13,m2C.getMeshLength()); + for i in xrange(12): + self.assertAlmostEqual(expected2[i],m2C.getCoords().getIJ(0,i),12); + pass + self.failUnlessEqual(expected3[0:4],m2C.getNodalConnectivity().getValues()[4:8]); + self.failUnlessEqual(expected3[4:8],m2C.getNodalConnectivity().getValues()[0:4]); + self.failUnlessEqual(expected3[8:13],m2C.getNodalConnectivity().getValues()[8:13]); + self.failUnlessEqual(expected4,m2C.getNodalConnectivityIndex().getValues()); + pass def setUp(self): pass diff --git a/src/MEDCoupling_Swig/libMEDCoupling_Swig.i b/src/MEDCoupling_Swig/libMEDCoupling_Swig.i index a1ee44c36..777765c8e 100644 --- a/src/MEDCoupling_Swig/libMEDCoupling_Swig.i +++ b/src/MEDCoupling_Swig/libMEDCoupling_Swig.i @@ -68,6 +68,8 @@ using namespace INTERP_KERNEL; %newobject ParaMEDMEM::MEDCouplingFieldDouble::max; %newobject ParaMEDMEM::MEDCouplingFieldDouble::minFields; %newobject ParaMEDMEM::MEDCouplingFieldDouble::min; +%newobject ParaMEDMEM::MEDCouplingFieldDouble::getIdsInRange; +%newobject ParaMEDMEM::MEDCouplingFieldDouble::buildSubPart; %newobject ParaMEDMEM::MEDCouplingFieldDouble::operator+; %newobject ParaMEDMEM::MEDCouplingFieldDouble::operator-; %newobject ParaMEDMEM::MEDCouplingFieldDouble::operator*; @@ -78,6 +80,7 @@ using namespace INTERP_KERNEL; %newobject ParaMEDMEM::DataArrayInt::deepCopy; %newobject ParaMEDMEM::DataArrayInt::performCpy; %newobject ParaMEDMEM::DataArrayInt::substr; +%newobject ParaMEDMEM::DataArrayInt::selectByTupleId; %newobject ParaMEDMEM::DataArrayDouble::aggregate; %newobject ParaMEDMEM::DataArrayDouble::dot; %newobject ParaMEDMEM::DataArrayDouble::crossProduct; @@ -86,12 +89,15 @@ using namespace INTERP_KERNEL; %newobject ParaMEDMEM::DataArrayDouble::multiply; %newobject ParaMEDMEM::DataArrayDouble::divide; %newobject ParaMEDMEM::DataArrayDouble::substr; +%newobject ParaMEDMEM::DataArrayDouble::getIdsInRange; +%newobject ParaMEDMEM::DataArrayDouble::selectByTupleId; %newobject ParaMEDMEM::MEDCouplingFieldDouble::clone; %newobject ParaMEDMEM::MEDCouplingFieldDouble::cloneWithMesh; %newobject ParaMEDMEM::MEDCouplingFieldDouble::buildNewTimeReprFromThis; %newobject ParaMEDMEM::MEDCouplingMesh::getCoordinatesAndOwner; %newobject ParaMEDMEM::MEDCouplingMesh::getBarycenterAndOwner; %newobject ParaMEDMEM::MEDCouplingMesh::buildOrthogonalField; +%newobject ParaMEDMEM::MEDCouplingMesh::buildPart; %newobject ParaMEDMEM::MEDCouplingMesh::mergeMyselfWith; %newobject ParaMEDMEM::MEDCouplingMesh::fillFromAnalytic; %newobject ParaMEDMEM::MEDCouplingPointSet::zipCoordsTraducer; @@ -187,6 +193,7 @@ namespace ParaMEDMEM virtual MEDCouplingFieldDouble *buildOrthogonalField() const = 0; virtual void rotate(const double *center, const double *vector, double angle) = 0; virtual void translate(const double *vector) = 0; + virtual MEDCouplingMesh *buildPart(const int *start, const int *end) const = 0; virtual MEDCouplingMesh *mergeMyselfWith(const MEDCouplingMesh *other) const throw(INTERP_KERNEL::Exception) = 0; virtual bool areCompatibleForMerge(const MEDCouplingMesh *other) const; static MEDCouplingMesh *mergeMeshes(const MEDCouplingMesh *mesh1, const MEDCouplingMesh *mesh2); @@ -799,6 +806,8 @@ namespace ParaMEDMEM double integral(int compId, bool isWAbs) const throw(INTERP_KERNEL::Exception); double normL1(int compId) const throw(INTERP_KERNEL::Exception); double normL2(int compId) const throw(INTERP_KERNEL::Exception); + DataArrayInt *getIdsInRange(double vmin, double vmax) const throw(INTERP_KERNEL::Exception); + MEDCouplingFieldDouble *buildSubPart(const DataArrayInt *part) const throw(INTERP_KERNEL::Exception); static MEDCouplingFieldDouble *mergeFields(const MEDCouplingFieldDouble *f1, const MEDCouplingFieldDouble *f2) throw(INTERP_KERNEL::Exception); static MEDCouplingFieldDouble *dotFields(const MEDCouplingFieldDouble *f1, const MEDCouplingFieldDouble *f2); MEDCouplingFieldDouble *dot(const MEDCouplingFieldDouble& other) const; @@ -971,6 +980,24 @@ namespace ParaMEDMEM } delete [] tmp; } + + MEDCouplingFieldDouble *buildSubPart(PyObject *li) const throw(INTERP_KERNEL::Exception) + { + int size; + int *tmp=convertPyToNewIntArr2(li,&size); + MEDCouplingFieldDouble *ret=0; + try + { + ret=self->buildSubPart(tmp,tmp+size); + } + catch(INTERP_KERNEL::Exception& e) + { + delete [] tmp; + throw e; + } + delete [] tmp; + return ret; + } } }; } -- 2.39.2