From 4f531e95b1c0c6881db0573c719f4c099b8ab6ae Mon Sep 17 00:00:00 2001 From: ageay Date: Tue, 26 Mar 2013 09:58:50 +0000 Subject: [PATCH] pow on DataArrayDouble --- src/MEDCoupling/MEDCouplingMemArray.cxx | 142 +++++++++++++++++- src/MEDCoupling/MEDCouplingMemArray.hxx | 4 + src/MEDCoupling_Swig/MEDCoupling.i | 3 + src/MEDCoupling_Swig/MEDCouplingBasicsTest.py | 45 ++++++ src/MEDCoupling_Swig/MEDCouplingCommon.i | 113 ++++++++++++++ src/MEDCoupling_Swig/MEDCouplingFinalize.i | 1 + 6 files changed, 307 insertions(+), 1 deletion(-) diff --git a/src/MEDCoupling/MEDCouplingMemArray.cxx b/src/MEDCoupling/MEDCouplingMemArray.cxx index d2a322a0d..81cdfe77c 100644 --- a/src/MEDCoupling/MEDCouplingMemArray.cxx +++ b/src/MEDCoupling/MEDCouplingMemArray.cxx @@ -3619,6 +3619,66 @@ DataArrayDouble *DataArrayDouble::negate() const throw(INTERP_KERNEL::Exception) return newArr; } +/*! + * Modify all elements of \a this array, so that + * an element _x_ becomes val ^ x . Contrary to DataArrayInt::applyPow + * all values in \a this have to be >= 0 if val is \b not integer. + * \param [in] val - the value used to apply pow on all array elements. + * \throw If \a this is not allocated. + * \warning If an exception is thrown because of presence of 0 element in \a this + * array and \a val is \b not integer, all elements processed before detection of the zero element remain + * modified. + */ +void DataArrayDouble::applyPow(double val) throw(INTERP_KERNEL::Exception) +{ + checkAllocated(); + double *ptr=getPointer(); + std::size_t nbOfElems=getNbOfElems(); + int val2=(int)val; + bool isInt=((double)val2)==val; + if(!isInt) + { + for(std::size_t i=0;i=0) + *ptr=pow(*ptr,val); + else + { + std::ostringstream oss; oss << "DataArrayDouble::applyPow (double) : At elem # " << i << " value is " << *ptr << " ! must be >=0. !"; + throw INTERP_KERNEL::Exception(oss.str().c_str()); + } + } + } + else + { + for(std::size_t i=0;i= 0 !"); + double *ptr=getPointer(); + std::size_t nbOfElems=getNbOfElems(); + for(std::size_t i=0;igetNumberOfTuples() != \a a2->getNumberOfTuples() + * \throw If \a a1->getNumberOfComponents() != 1 or \a a2->getNumberOfComponents() != 1. + * \throw If there is a negative value in \a a1. + */ +DataArrayDouble *DataArrayDouble::Pow(const DataArrayDouble *a1, const DataArrayDouble *a2) throw(INTERP_KERNEL::Exception) +{ + if(!a1 || !a2) + throw INTERP_KERNEL::Exception("DataArrayDouble::Pow : at least one of input instances is null !"); + int nbOfTuple=a1->getNumberOfTuples(); + int nbOfTuple2=a2->getNumberOfTuples(); + int nbOfComp=a1->getNumberOfComponents(); + int nbOfComp2=a2->getNumberOfComponents(); + if(nbOfTuple!=nbOfTuple2) + throw INTERP_KERNEL::Exception("DataArrayDouble::Pow : number of tuples mismatches !"); + if(nbOfComp!=1 || nbOfComp2!=1) + throw INTERP_KERNEL::Exception("DataArrayDouble::Pow : number of components of both arrays must be equal to 1 !"); + MEDCouplingAutoRefCountObjectPtr ret=DataArrayDouble::New(); ret->alloc(nbOfTuple,1); + const double *ptr1(a1->begin()),*ptr2(a2->begin()); + double *ptr=ret->getPointer(); + for(int i=0;i=0) + { + *ptr=pow(*ptr1,*ptr2); + } + else + { + std::ostringstream oss; oss << "DataArrayDouble::Pow : on tuple #" << i << " of a1 value is < 0 (" << *ptr1 << ") !"; + throw INTERP_KERNEL::Exception(oss.str().c_str()); + } + } + return ret.retn(); +} + +/*! + * Apply pow on values of another DataArrayDouble to values of \a this one. + * + * \param [in] other - an array to pow to \a this one. + * \throw If \a other is NULL. + * \throw If \a this->getNumberOfTuples() != \a other->getNumberOfTuples() + * \throw If \a this->getNumberOfComponents() != 1 or \a other->getNumberOfComponents() != 1 + * \throw If there is a negative value in \a this. + */ +void DataArrayDouble::powEqual(const DataArrayDouble *other) throw(INTERP_KERNEL::Exception) +{ + if(!other) + throw INTERP_KERNEL::Exception("DataArrayDouble::powEqual : input instance is null !"); + int nbOfTuple=getNumberOfTuples(); + int nbOfTuple2=other->getNumberOfTuples(); + int nbOfComp=getNumberOfComponents(); + int nbOfComp2=other->getNumberOfComponents(); + if(nbOfTuple!=nbOfTuple2) + throw INTERP_KERNEL::Exception("DataArrayDouble::powEqual : number of tuples mismatches !"); + if(nbOfComp!=1 || nbOfComp2!=1) + throw INTERP_KERNEL::Exception("DataArrayDouble::powEqual : number of components of both arrays must be equal to 1 !"); + double *ptr=getPointer(); + const double *ptrc=other->begin(); + for(int i=0;i=0) + *ptr=pow(*ptr,*ptrc); + else + { + std::ostringstream oss; oss << "DataArrayDouble::powEqual : on tuple #" << i << " of this value is < 0 (" << *ptr << ") !"; + throw INTERP_KERNEL::Exception(oss.str().c_str()); + } + } + declareAsNew(); +} + /*! * Useless method for end user. Only for MPI/Corba/File serialsation for multi arrays class. * Server side. @@ -8094,7 +8234,7 @@ void DataArrayInt::applyRPow(int val) throw(INTERP_KERNEL::Exception) int tmp=1; for(int j=0;j<*ptr;j++) tmp*=val; - *ptr=val; + *ptr=tmp; } else { diff --git a/src/MEDCoupling/MEDCouplingMemArray.hxx b/src/MEDCoupling/MEDCouplingMemArray.hxx index 6bb2c48ee..cfebca035 100644 --- a/src/MEDCoupling/MEDCouplingMemArray.hxx +++ b/src/MEDCoupling/MEDCouplingMemArray.hxx @@ -287,6 +287,8 @@ namespace ParaMEDMEM MEDCOUPLING_EXPORT void applyLin(double a, double b, int compoId) throw(INTERP_KERNEL::Exception); MEDCOUPLING_EXPORT void applyLin(double a, double b) throw(INTERP_KERNEL::Exception); MEDCOUPLING_EXPORT void applyInv(double numerator) throw(INTERP_KERNEL::Exception); + MEDCOUPLING_EXPORT void applyPow(double val) throw(INTERP_KERNEL::Exception); + MEDCOUPLING_EXPORT void applyRPow(double val) throw(INTERP_KERNEL::Exception); MEDCOUPLING_EXPORT DataArrayDouble *negate() const throw(INTERP_KERNEL::Exception); MEDCOUPLING_EXPORT DataArrayDouble *applyFunc(int nbOfComp, FunctionToEvaluate func) const throw(INTERP_KERNEL::Exception); MEDCOUPLING_EXPORT DataArrayDouble *applyFunc(int nbOfComp, const char *func) const throw(INTERP_KERNEL::Exception); @@ -312,6 +314,8 @@ namespace ParaMEDMEM MEDCOUPLING_EXPORT void multiplyEqual(const DataArrayDouble *other) throw(INTERP_KERNEL::Exception); MEDCOUPLING_EXPORT static DataArrayDouble *Divide(const DataArrayDouble *a1, const DataArrayDouble *a2) throw(INTERP_KERNEL::Exception); MEDCOUPLING_EXPORT void divideEqual(const DataArrayDouble *other) throw(INTERP_KERNEL::Exception); + MEDCOUPLING_EXPORT static DataArrayDouble *Pow(const DataArrayDouble *a1, const DataArrayDouble *a2) throw(INTERP_KERNEL::Exception); + MEDCOUPLING_EXPORT void powEqual(const DataArrayDouble *other) throw(INTERP_KERNEL::Exception); MEDCOUPLING_EXPORT void updateTime() const { } public: MEDCOUPLING_EXPORT void getTinySerializationIntInformation(std::vector& tinyInfo) const; diff --git a/src/MEDCoupling_Swig/MEDCoupling.i b/src/MEDCoupling_Swig/MEDCoupling.i index de843c7b6..754e78b4a 100644 --- a/src/MEDCoupling_Swig/MEDCoupling.i +++ b/src/MEDCoupling_Swig/MEDCoupling.i @@ -32,6 +32,9 @@ def ParaMEDMEMDataArrayDoubleImul(self,*args): def ParaMEDMEMDataArrayDoubleIdiv(self,*args): import _MEDCoupling return _MEDCoupling.DataArrayDouble____idiv___(self, self, *args) +def ParaMEDMEMDataArrayDoubleIpow(self,*args): + import _MEDCoupling + return _MEDCoupling.DataArrayDouble____ipow___(self, self, *args) def ParaMEDMEMMEDCouplingFieldDoubleIadd(self,*args): import _MEDCoupling return _MEDCoupling.MEDCouplingFieldDouble____iadd___(self, self, *args) diff --git a/src/MEDCoupling_Swig/MEDCouplingBasicsTest.py b/src/MEDCoupling_Swig/MEDCouplingBasicsTest.py index 82d5e46a9..e20306b2a 100644 --- a/src/MEDCoupling_Swig/MEDCouplingBasicsTest.py +++ b/src/MEDCoupling_Swig/MEDCouplingBasicsTest.py @@ -11714,6 +11714,51 @@ class MEDCouplingBasicsTest(unittest.TestCase): self.assertTrue(d.isEqual(DataArrayInt([6,16,5,15,4,14,3,13,2,12,1,11,0,10],7,2))) pass + def testSwig2DAPow1(self): + d=DataArrayInt(10) + d.iota(0) + d1=d.deepCpy() + d.setIJ(2,0,-2) + self.assertTrue((d**2).isEqual(DataArrayInt([0,1,4,9,16,25,36,49,64,81]))) + self.assertTrue((d**3).isEqual(DataArrayInt([0,1,-8,27,64,125,216,343,512,729]))) + for elt in [d]: + elt**=2 + pass + self.assertTrue(d.isEqual(DataArrayInt([0,1,4,9,16,25,36,49,64,81]))) + self.assertTrue((d1[:4]**d1[:4]).isEqual(DataArrayInt([1,1,4,27]))) + self.assertTrue((3**d1[:4]).isEqual(DataArrayInt([1,3,9,27]))) + d2=d1[:4] + d2**=d2 + self.assertTrue(d2.isEqual(DataArrayInt([1,1,4,27]))) + self.assertRaises(InterpKernelException,d2.__pow__,-1)#non supporting negative pow in DataArrayInt.__pow__ + self.assertRaises(InterpKernelException,d2.__ipow__,-1)#non supporting negative pow in DataArrayInt.__pow__ + # + d=DataArrayDouble(10) + d.iota(0) + d1=d.deepCpy() + d.setIJ(2,0,-2.) + self.assertTrue((d**2).isEqual(DataArrayDouble([0,1,4,9,16,25,36,49,64,81]),1e-12)) + self.assertTrue((d**3).isEqual(DataArrayDouble([0,1,-8,27,64,125,216,343,512,729]),1e-12)) + self.assertRaises(InterpKernelException,d.__pow__,3.1)#3.1 is double not integer -> not supporting negative values in d + for elt in [d]: + elt**=2 + pass + self.assertTrue(d.isEqual(DataArrayDouble([0,1,4,9,16,25,36,49,64,81]),1e-12)) + self.assertTrue((d1[:4]**d1[:4]).isEqual(DataArrayDouble([1,1,4,27]),1e-12)) + self.assertTrue((3**d1[:4]).isEqual(DataArrayDouble([1,3,9,27]),1e-12)) + d2=d1[:4] + d2**=d2 + self.assertTrue(d2.isEqual(DataArrayDouble([1,1,4,27]),1e-12)) + d2**=-0.5 + self.assertTrue(d2.isEqual(DataArrayDouble([1,1,1./2,1./sqrt(27.)]),1e-14)) + d3=-1./d1[1:5] + self.assertTrue((3**d3).isEqual(DataArrayDouble([0.3333333333333333,0.5773502691896257,0.6933612743506348,0.7598356856515925]),1e-14)) + d4=d3.deepCpy() ; d4.abs() + self.assertTrue((d4**d3).isEqual(DataArrayDouble([1.,sqrt(2.),1.4422495703074083,sqrt(2.)]),1e-14)) + d4**=d3 + self.assertTrue(d4.isEqual(DataArrayDouble([1.,sqrt(2.),1.4422495703074083,sqrt(2.)]),1e-14)) + pass + def setUp(self): pass pass diff --git a/src/MEDCoupling_Swig/MEDCouplingCommon.i b/src/MEDCoupling_Swig/MEDCouplingCommon.i index b60d7000e..f1adae767 100644 --- a/src/MEDCoupling_Swig/MEDCouplingCommon.i +++ b/src/MEDCoupling_Swig/MEDCouplingCommon.i @@ -244,6 +244,7 @@ using namespace INTERP_KERNEL; %newobject ParaMEDMEM::DataArrayDouble::Substract; %newobject ParaMEDMEM::DataArrayDouble::Multiply; %newobject ParaMEDMEM::DataArrayDouble::Divide; +%newobject ParaMEDMEM::DataArrayDouble::Pow; %newobject ParaMEDMEM::DataArrayDouble::substr; %newobject ParaMEDMEM::DataArrayDouble::changeNbOfComponents; %newobject ParaMEDMEM::DataArrayDouble::keepSelectedComponents; @@ -288,6 +289,8 @@ using namespace INTERP_KERNEL; %newobject ParaMEDMEM::DataArrayDouble::__rmul__; %newobject ParaMEDMEM::DataArrayDouble::__div__; %newobject ParaMEDMEM::DataArrayDouble::__rdiv__; +%newobject ParaMEDMEM::DataArrayDouble::__pow__; +%newobject ParaMEDMEM::DataArrayDouble::__rpow__; %newobject ParaMEDMEM::DataArrayDoubleTuple::buildDADouble; %newobject ParaMEDMEM::MEDCouplingMesh::deepCpy; %newobject ParaMEDMEM::MEDCouplingMesh::checkDeepEquivalOnSameNodesWith; @@ -5437,6 +5440,116 @@ namespace ParaMEDMEM } } + DataArrayDouble *__pow__(PyObject *obj) throw(INTERP_KERNEL::Exception) + { + const char msg[]="Unexpected situation in __pow__ !"; + double val; + DataArrayDouble *a; + DataArrayDoubleTuple *aa; + std::vector bb; + int sw; + convertObjToPossibleCpp5(obj,sw,val,a,aa,bb); + switch(sw) + { + case 1: + { + MEDCouplingAutoRefCountObjectPtr ret=self->deepCpy(); + ret->applyPow(val); + return ret.retn(); + } + case 2: + { + return DataArrayDouble::Pow(self,a); + } + case 3: + { + MEDCouplingAutoRefCountObjectPtr aaa=aa->buildDADouble(1,self->getNumberOfComponents()); + return DataArrayDouble::Pow(self,aaa); + } + case 4: + { + MEDCouplingAutoRefCountObjectPtr aaa=DataArrayDouble::New(); aaa->useArray(&bb[0],false,CPP_DEALLOC,1,(int)bb.size()); + return DataArrayDouble::Pow(self,aaa); + } + default: + throw INTERP_KERNEL::Exception(msg); + } + } + + DataArrayDouble *__rpow__(PyObject *obj) throw(INTERP_KERNEL::Exception) + { + const char msg[]="Unexpected situation in __rpow__ !"; + double val; + DataArrayDouble *a; + DataArrayDoubleTuple *aa; + std::vector bb; + int sw; + convertObjToPossibleCpp5(obj,sw,val,a,aa,bb); + switch(sw) + { + case 1: + { + MEDCouplingAutoRefCountObjectPtr ret=self->deepCpy(); + ret->applyRPow(val); + return ret.retn(); + } + case 3: + { + MEDCouplingAutoRefCountObjectPtr aaa=aa->buildDADouble(1,self->getNumberOfComponents()); + return DataArrayDouble::Pow(aaa,self); + } + case 4: + { + MEDCouplingAutoRefCountObjectPtr aaa=DataArrayDouble::New(); aaa->useArray(&bb[0],false,CPP_DEALLOC,1,(int)bb.size()); + return DataArrayDouble::Pow(aaa,self); + } + default: + throw INTERP_KERNEL::Exception(msg); + } + } + + PyObject *___ipow___(PyObject *trueSelf, PyObject *obj) throw(INTERP_KERNEL::Exception) + { + const char msg[]="Unexpected situation in __ipow__ !"; + double val; + DataArrayDouble *a; + DataArrayDoubleTuple *aa; + std::vector bb; + int sw; + convertObjToPossibleCpp5(obj,sw,val,a,aa,bb); + switch(sw) + { + case 1: + { + self->applyPow(val); + Py_XINCREF(trueSelf); + return trueSelf; + } + case 2: + { + self->powEqual(a); + Py_XINCREF(trueSelf); + return trueSelf; + } + case 3: + { + MEDCouplingAutoRefCountObjectPtr aaa=aa->buildDADouble(1,self->getNumberOfComponents()); + self->powEqual(aaa); + Py_XINCREF(trueSelf); + return trueSelf; + } + case 4: + { + MEDCouplingAutoRefCountObjectPtr aaa=DataArrayDouble::New(); aaa->useArray(&bb[0],false,CPP_DEALLOC,1,(int)bb.size()); + self->powEqual(aaa); + Py_XINCREF(trueSelf); + return trueSelf; + } + default: + throw INTERP_KERNEL::Exception(msg); + } + } + PyObject *computeTupleIdsNearTuples(const DataArrayDouble *other, double eps) { DataArrayInt *c=0,*cI=0; diff --git a/src/MEDCoupling_Swig/MEDCouplingFinalize.i b/src/MEDCoupling_Swig/MEDCouplingFinalize.i index d686aa607..f6c73e1dc 100644 --- a/src/MEDCoupling_Swig/MEDCouplingFinalize.i +++ b/src/MEDCoupling_Swig/MEDCouplingFinalize.i @@ -22,6 +22,7 @@ DataArrayDouble.__iadd__=ParaMEDMEMDataArrayDoubleIadd DataArrayDouble.__isub__=ParaMEDMEMDataArrayDoubleIsub DataArrayDouble.__imul__=ParaMEDMEMDataArrayDoubleImul DataArrayDouble.__idiv__=ParaMEDMEMDataArrayDoubleIdiv +DataArrayDouble.__ipow__=ParaMEDMEMDataArrayDoubleIpow DataArrayInt.__iadd__=ParaMEDMEMDataArrayIntIadd DataArrayInt.__isub__=ParaMEDMEMDataArrayIntIsub -- 2.39.2