From: ageay Date: Fri, 17 Sep 2010 11:02:40 +0000 (+0000) Subject: Modification of normL1 and L2 API and their implementation. X-Git-Tag: V5_1_main_FINAL~41 X-Git-Url: http://git.salome-platform.org/gitweb/?a=commitdiff_plain;h=4822ae828becb825c7e56517d2478c93b56e0d6e;p=tools%2Fmedcoupling.git Modification of normL1 and L2 API and their implementation. --- diff --git a/src/MEDCoupling/MEDCouplingField.cxx b/src/MEDCoupling/MEDCouplingField.cxx index b50f34fe5..d8b9e7d9f 100644 --- a/src/MEDCoupling/MEDCouplingField.cxx +++ b/src/MEDCoupling/MEDCouplingField.cxx @@ -78,6 +78,17 @@ TypeOfField MEDCouplingField::getTypeOfField() const return _type->getEnum(); } +/*! + * This method retrieves the weighting field of 'this'. If no '_mesh' is defined an exception will be thrown. + * Warning the retrieved field life cycle is the responsability of caller. + */ +MEDCouplingFieldDouble *MEDCouplingField::buildWeightingField(bool isAbs) const throw(INTERP_KERNEL::Exception) +{ + if(_mesh==0) + throw INTERP_KERNEL::Exception("MEDCouplingField::getWeightingField : no mesh defined !!!"); + return _type->getWeightingField(_mesh,isAbs); +} + void MEDCouplingField::setMesh(const MEDCouplingMesh *mesh) { if(mesh!=_mesh) diff --git a/src/MEDCoupling/MEDCouplingField.hxx b/src/MEDCoupling/MEDCouplingField.hxx index d7d45f29f..8287acecd 100644 --- a/src/MEDCoupling/MEDCouplingField.hxx +++ b/src/MEDCoupling/MEDCouplingField.hxx @@ -33,6 +33,7 @@ namespace ParaMEDMEM { class DataArrayInt; class MEDCouplingMesh; + class MEDCouplingFieldDouble; class MEDCouplingFieldDiscretization; class MEDCouplingGaussLocalization; @@ -50,6 +51,7 @@ namespace ParaMEDMEM void setDescription(const char *desc) { _desc=desc; } const char *getName() const { return _name.c_str(); } TypeOfField getTypeOfField() const; + MEDCouplingFieldDouble *buildWeightingField(bool isAbs) const throw(INTERP_KERNEL::Exception); MEDCouplingMesh *buildSubMeshData(const int *start, const int *end, DataArrayInt *&di) const; MEDCouplingFieldDiscretization *getDiscretization() const { return _type; } // Gauss point specific methods diff --git a/src/MEDCoupling/MEDCouplingFieldDiscretization.cxx b/src/MEDCoupling/MEDCouplingFieldDiscretization.cxx index f34dbf5da..47630d093 100644 --- a/src/MEDCoupling/MEDCouplingFieldDiscretization.cxx +++ b/src/MEDCoupling/MEDCouplingFieldDiscretization.cxx @@ -97,22 +97,19 @@ void MEDCouplingFieldDiscretization::updateTime() * @param res output parameter expected to be of size arr->getNumberOfComponents(); * @throw when the field discretization fails on getWeighting fields (gauss points for example) */ -void MEDCouplingFieldDiscretization::normL1(const MEDCouplingMesh *mesh, const DataArrayDouble *arr, bool isWAbs, double *res) const throw(INTERP_KERNEL::Exception) +void MEDCouplingFieldDiscretization::normL1(const MEDCouplingMesh *mesh, const DataArrayDouble *arr, double *res) const throw(INTERP_KERNEL::Exception) { - MEDCouplingFieldDouble *vol=getWeightingField(mesh,isWAbs); + MEDCouplingFieldDouble *vol=getWeightingField(mesh,true); int nbOfCompo=arr->getNumberOfComponents(); int nbOfElems=getNumberOfTuples(mesh); std::fill(res,res+nbOfCompo,0.); - double totVol=0.; const double *arrPtr=arr->getConstPointer(); const double *volPtr=vol->getArray()->getConstPointer(); for(int i=0;i(),1/totVol)); vol->decrRef(); } @@ -121,22 +118,19 @@ void MEDCouplingFieldDiscretization::normL1(const MEDCouplingMesh *mesh, const D * @param res output parameter expected to be of size arr->getNumberOfComponents(); * @throw when the field discretization fails on getWeighting fields (gauss points for example) */ -void MEDCouplingFieldDiscretization::normL2(const MEDCouplingMesh *mesh, const DataArrayDouble *arr, bool isWAbs, double *res) const throw(INTERP_KERNEL::Exception) +void MEDCouplingFieldDiscretization::normL2(const MEDCouplingMesh *mesh, const DataArrayDouble *arr, double *res) const throw(INTERP_KERNEL::Exception) { - MEDCouplingFieldDouble *vol=getWeightingField(mesh,isWAbs); + MEDCouplingFieldDouble *vol=getWeightingField(mesh,true); int nbOfCompo=arr->getNumberOfComponents(); int nbOfElems=getNumberOfTuples(mesh); std::fill(res,res+nbOfCompo,0.); - double totVol=0.; const double *arrPtr=arr->getConstPointer(); const double *volPtr=vol->getArray()->getConstPointer(); for(int i=0;i(),1/totVol)); std::transform(res,res+nbOfCompo,res,std::ptr_fun(std::sqrt)); vol->decrRef(); } diff --git a/src/MEDCoupling/MEDCouplingFieldDiscretization.hxx b/src/MEDCoupling/MEDCouplingFieldDiscretization.hxx index 1e24df87d..5e2ab13fc 100644 --- a/src/MEDCoupling/MEDCouplingFieldDiscretization.hxx +++ b/src/MEDCoupling/MEDCouplingFieldDiscretization.hxx @@ -49,8 +49,8 @@ namespace ParaMEDMEM virtual MEDCouplingFieldDiscretization *clone() const = 0; virtual const char *getStringRepr() const = 0; virtual int getNumberOfTuples(const MEDCouplingMesh *mesh) const = 0; - virtual void normL1(const MEDCouplingMesh *mesh, const DataArrayDouble *arr, bool isWAbs, double *res) const throw(INTERP_KERNEL::Exception); - virtual void normL2(const MEDCouplingMesh *mesh, const DataArrayDouble *arr, bool isWAbs, double *res) const throw(INTERP_KERNEL::Exception); + virtual void normL1(const MEDCouplingMesh *mesh, const DataArrayDouble *arr, double *res) const throw(INTERP_KERNEL::Exception); + 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 checkCompatibilityWithNature(NatureOfField nat) const throw(INTERP_KERNEL::Exception) = 0; diff --git a/src/MEDCoupling/MEDCouplingFieldDouble.cxx b/src/MEDCoupling/MEDCouplingFieldDouble.cxx index 3eb999f21..bfefcb24c 100644 --- a/src/MEDCoupling/MEDCouplingFieldDouble.cxx +++ b/src/MEDCoupling/MEDCouplingFieldDouble.cxx @@ -295,7 +295,7 @@ void MEDCouplingFieldDouble::accumulate(double *res) const * \f] * If compId>=nbOfComponent an exception is thrown. */ -double MEDCouplingFieldDouble::normL1(int compId, bool isWAbs) const throw(INTERP_KERNEL::Exception) +double MEDCouplingFieldDouble::normL1(int compId) const throw(INTERP_KERNEL::Exception) { if(!_mesh) throw INTERP_KERNEL::Exception("No mesh underlying this field to perform normL1"); @@ -305,7 +305,7 @@ double MEDCouplingFieldDouble::normL1(int compId, bool isWAbs) const throw(INTER double *res=new double[nbComps]; try { - _type->normL1(_mesh,getArray(),isWAbs,res); + _type->normL1(_mesh,getArray(),res); } catch(INTERP_KERNEL::Exception& e) { @@ -324,11 +324,11 @@ double MEDCouplingFieldDouble::normL1(int compId, bool isWAbs) const throw(INTER * \f] * The res is expected to be of size getNumberOfComponents(). */ -void MEDCouplingFieldDouble::normL1(bool isWAbs, double *res) const throw(INTERP_KERNEL::Exception) +void MEDCouplingFieldDouble::normL1(double *res) const throw(INTERP_KERNEL::Exception) { if(!_mesh) throw INTERP_KERNEL::Exception("No mesh underlying this field to perform normL1"); - _type->normL1(_mesh,getArray(),isWAbs,res); + _type->normL1(_mesh,getArray(),res); } /*! @@ -338,17 +338,17 @@ void MEDCouplingFieldDouble::normL1(bool isWAbs, double *res) const throw(INTERP * \f] * If compId>=nbOfComponent an exception is thrown. */ -double MEDCouplingFieldDouble::normL2(int compId, bool isWAbs) const throw(INTERP_KERNEL::Exception) +double MEDCouplingFieldDouble::normL2(int compId) const throw(INTERP_KERNEL::Exception) { if(!_mesh) - throw INTERP_KERNEL::Exception("No mesh underlying this field to perform normL1"); + throw INTERP_KERNEL::Exception("No mesh underlying this field to perform normL2"); int nbComps=getArray()->getNumberOfComponents(); if(compId>=nbComps) throw INTERP_KERNEL::Exception("Invalid compId specified : No such nb of components !"); double *res=new double[nbComps]; try { - _type->normL2(_mesh,getArray(),isWAbs,res); + _type->normL2(_mesh,getArray(),res); } catch(INTERP_KERNEL::Exception& e) { @@ -367,11 +367,11 @@ double MEDCouplingFieldDouble::normL2(int compId, bool isWAbs) const throw(INTER * \f] * The res is expected to be of size getNumberOfComponents(). */ -void MEDCouplingFieldDouble::normL2(bool isWAbs, double *res) const throw(INTERP_KERNEL::Exception) +void MEDCouplingFieldDouble::normL2(double *res) const throw(INTERP_KERNEL::Exception) { if(!_mesh) - throw INTERP_KERNEL::Exception("No mesh underlying this field to perform normL1"); - _type->normL2(_mesh,getArray(),isWAbs,res); + throw INTERP_KERNEL::Exception("No mesh underlying this field to perform normL2"); + _type->normL2(_mesh,getArray(),res); } /*! diff --git a/src/MEDCoupling/MEDCouplingFieldDouble.hxx b/src/MEDCoupling/MEDCouplingFieldDouble.hxx index b57ae8f40..b6524128a 100644 --- a/src/MEDCoupling/MEDCouplingFieldDouble.hxx +++ b/src/MEDCoupling/MEDCouplingFieldDouble.hxx @@ -63,10 +63,10 @@ namespace ParaMEDMEM double accumulate(int compId) const; void accumulate(double *res) const; double getMaxValue() const throw(INTERP_KERNEL::Exception); - double normL1(int compId, bool isWAbs) const throw(INTERP_KERNEL::Exception); - void normL1(bool isWAbs, double *res) const throw(INTERP_KERNEL::Exception); - double normL2(int compId, bool isWAbs) const throw(INTERP_KERNEL::Exception); - void normL2(bool isWAbs, double *res) const throw(INTERP_KERNEL::Exception); + double normL1(int compId) const throw(INTERP_KERNEL::Exception); + void normL1(double *res) const throw(INTERP_KERNEL::Exception); + double normL2(int compId) const throw(INTERP_KERNEL::Exception); + void normL2(double *res) const throw(INTERP_KERNEL::Exception); double integral(int compId, bool isWAbs) const throw(INTERP_KERNEL::Exception); void integral(bool isWAbs, double *res) const throw(INTERP_KERNEL::Exception); void getValueOnPos(int i, int j, int k, double *res) const throw(INTERP_KERNEL::Exception); diff --git a/src/MEDCoupling/Test/MEDCouplingBasicsTest2.cxx b/src/MEDCoupling/Test/MEDCouplingBasicsTest2.cxx index ef9f175f3..d0ff53675 100644 --- a/src/MEDCoupling/Test/MEDCouplingBasicsTest2.cxx +++ b/src/MEDCoupling/Test/MEDCouplingBasicsTest2.cxx @@ -372,29 +372,28 @@ void MEDCouplingBasicsTest::testNormL12Integ1D() for(int i=0;i<3;i++) CPPUNIT_ASSERT_DOUBLES_EQUAL(expected4[i],res[i],1e-12); //normL1 - f1->normL1(false,res); - double expected5[3]={16.377,24.3225,34.6715}; + f1->normL1(res); + double expected5[3]={11.3068,27.3621,43.7881}; for(int i=0;i<3;i++) CPPUNIT_ASSERT_DOUBLES_EQUAL(expected5[i],res[i],1e-12); - CPPUNIT_ASSERT_DOUBLES_EQUAL(expected5[0],f1->normL1(0,false),1e-12); - CPPUNIT_ASSERT_DOUBLES_EQUAL(expected5[1],f1->normL1(1,false),1e-12); - CPPUNIT_ASSERT_DOUBLES_EQUAL(expected5[2],f1->normL1(2,false),1e-12); - f1->normL1(true,res); - double expected6[3]={6.97950617284,16.8901851851852,27.0296913580247}; - for(int i=0;i<3;i++) - CPPUNIT_ASSERT_DOUBLES_EQUAL(expected6[i],res[i],1e-12); + CPPUNIT_ASSERT_DOUBLES_EQUAL(expected5[0],f1->normL1(0),1e-12); + CPPUNIT_ASSERT_DOUBLES_EQUAL(expected5[1],f1->normL1(1),1e-12); + CPPUNIT_ASSERT_DOUBLES_EQUAL(expected5[2],f1->normL1(2),1e-12); //normL2 - f1->normL2(false,res); - double expected7[3]={13.3073310622,23.1556650736,33.8113075021}; + f1->normL2(res); + double expected7[3]={9.0252562290496776, 21.545259176904789, 34.433193070059595}; for(int i=0;i<3;i++) CPPUNIT_ASSERT_DOUBLES_EQUAL(expected7[i],res[i],1e-9); - CPPUNIT_ASSERT_DOUBLES_EQUAL(expected7[0],f1->normL2(0,false),1e-9); - CPPUNIT_ASSERT_DOUBLES_EQUAL(expected7[1],f1->normL2(1,false),1e-9); - CPPUNIT_ASSERT_DOUBLES_EQUAL(expected7[2],f1->normL2(2,false),1e-9); - f1->normL2(true,res); - double expected8[3]={7.09091097945239,16.9275542960123,27.0532714641609}; - for(int i=0;i<3;i++) - CPPUNIT_ASSERT_DOUBLES_EQUAL(expected8[i],res[i],1e-9); + CPPUNIT_ASSERT_DOUBLES_EQUAL(expected7[0],f1->normL2(0),1e-9); + CPPUNIT_ASSERT_DOUBLES_EQUAL(expected7[1],f1->normL2(1),1e-9); + CPPUNIT_ASSERT_DOUBLES_EQUAL(expected7[2],f1->normL2(2),1e-9); + //buildWeightingField + MEDCouplingFieldDouble *f4=f1->buildWeightingField(false); + CPPUNIT_ASSERT_DOUBLES_EQUAL(-0.2,f4->accumulate(0),1e-12); + f4->decrRef(); + f4=f1->buildWeightingField(true); + CPPUNIT_ASSERT_DOUBLES_EQUAL(1.62,f4->accumulate(0),1e-12); + f4->decrRef(); // f1->decrRef(); m1->decrRef(); @@ -437,18 +436,12 @@ void MEDCouplingBasicsTest::testNormL12Integ1D() f1->integral(true,res); for(int i=0;i<3;i++) CPPUNIT_ASSERT_DOUBLES_EQUAL(sqrt(2.)*expected4[i],res[i],1e-12); - f1->normL1(false,res); - for(int i=0;i<3;i++) - CPPUNIT_ASSERT_DOUBLES_EQUAL(expected6[i],res[i],1e-12); - f1->normL1(true,res); - for(int i=0;i<3;i++) - CPPUNIT_ASSERT_DOUBLES_EQUAL(expected6[i],res[i],1e-12); - f1->normL2(false,res); + f1->normL1(res); for(int i=0;i<3;i++) - CPPUNIT_ASSERT_DOUBLES_EQUAL(expected8[i],res[i],1e-12); - f1->normL2(true,res); + CPPUNIT_ASSERT_DOUBLES_EQUAL(sqrt(2.)*expected5[i],res[i],1e-12); + f1->normL2(res); for(int i=0;i<3;i++) - CPPUNIT_ASSERT_DOUBLES_EQUAL(expected8[i],res[i],1e-12); + CPPUNIT_ASSERT_DOUBLES_EQUAL(sqrt(sqrt(2.))*expected7[i],res[i],1e-12); // f1->decrRef(); m1->decrRef(); diff --git a/src/MEDCoupling_Swig/MEDCouplingBasicsTest.py b/src/MEDCoupling_Swig/MEDCouplingBasicsTest.py index 3c715cdb7..6773e5c21 100644 --- a/src/MEDCoupling_Swig/MEDCouplingBasicsTest.py +++ b/src/MEDCoupling_Swig/MEDCouplingBasicsTest.py @@ -2010,34 +2010,28 @@ class MEDCouplingBasicsTest(unittest.TestCase): self.assertTrue(abs(expected4[i]-res[i])<1e-12); pass #normL1 - res=f1.normL1(False); - expected5=[16.377,24.3225,34.6715] + res=f1.normL1(); + expected5=[11.3068,27.3621,43.7881] for i in xrange(3): self.assertTrue(abs(expected5[i]-res[i])<1e-12); pass - self.assertTrue(abs(expected5[0]-f1.normL1(0,False))<1e-12); - self.assertTrue(abs(expected5[1]-f1.normL1(1,False))<1e-12); - self.assertTrue(abs(expected5[2]-f1.normL1(2,False))<1e-12); - res=f1.normL1(True); - expected6=[6.97950617284,16.8901851851852,27.0296913580247] - for i in xrange(3): - self.assertTrue(abs(expected6[i]-res[i])<1e-12); - pass + self.assertTrue(abs(expected5[0]-f1.normL1(0))<1e-12); + self.assertTrue(abs(expected5[1]-f1.normL1(1))<1e-12); + self.assertTrue(abs(expected5[2]-f1.normL1(2))<1e-12); #normL2 - res=f1.normL2(False); - expected7=[13.3073310622,23.1556650736,33.8113075021] + res=f1.normL2(); + expected7=[9.0252562290496776, 21.545259176904789, 34.433193070059595] for i in xrange(3): self.assertTrue(abs(expected7[i]-res[i])<1e-9); pass - self.assertTrue(abs(expected7[0]-f1.normL2(0,False))<1e-9); - self.assertTrue(abs(expected7[1]-f1.normL2(1,False))<1e-9); - self.assertTrue(abs(expected7[2]-f1.normL2(2,False))<1e-9); - res=f1.normL2(True); - expected8=[7.09091097945239,16.9275542960123,27.0532714641609] - for i in xrange(3): - self.assertTrue(abs(expected8[i]-res[i])<1e-9); - pass - # + self.assertTrue(abs(expected7[0]-f1.normL2(0))<1e-9); + self.assertTrue(abs(expected7[1]-f1.normL2(1))<1e-9); + self.assertTrue(abs(expected7[2]-f1.normL2(2))<1e-9); + #buildWeightingField + f4=f1.buildWeightingField(False); + self.assertTrue(abs(-0.2-f4.accumulate(0))<1e-12); + f4=f1.buildWeightingField(True); + self.assertTrue(abs(1.62-f4.accumulate(0))<1e-12); # Testing with 2D Curve m1=MEDCouplingDataForTest.build2DCurveTargetMesh_3(); f2=m1.getMeasureField(False); @@ -2077,21 +2071,13 @@ class MEDCouplingBasicsTest(unittest.TestCase): for i in xrange(3): self.assertTrue(abs(sqrt(2.)*expected4[i]-res[i])<1e-12); pass - res=f1.normL1(False); - for i in xrange(3): - self.assertTrue(abs(expected6[i]-res[i])<1e-12); - pass - res=f1.normL1(True); - for i in xrange(3): - self.assertTrue(abs(expected6[i]-res[i])<1e-12); - pass - res=f1.normL2(False); + res=f1.normL1(); for i in xrange(3): - self.assertTrue(abs(expected8[i]-res[i])<1e-12); + self.assertTrue(abs(sqrt(2.)*expected5[i]-res[i])<1e-12); pass - res=f1.normL2(True); + res=f1.normL2(); for i in xrange(3): - self.assertTrue(abs(expected8[i]-res[i])<1e-12); + self.assertTrue(abs(sqrt(sqrt(2.))*expected7[i]-res[i])<1e-12); pass pass @@ -2645,7 +2631,7 @@ class MEDCouplingBasicsTest(unittest.TestCase): self.assertAlmostEqual(9.5,f.getMaxValue(),14); pass - def test(self): + def testSubstractInPlaceDM1(self): mesh1=MEDCouplingDataForTest.build2DTargetMesh_3(); mesh2=MEDCouplingDataForTest.build2DTargetMesh_3(); f1=MEDCouplingFieldDouble.New(ON_CELLS,NO_TIME); diff --git a/src/MEDCoupling_Swig/libMEDCoupling_Swig.i b/src/MEDCoupling_Swig/libMEDCoupling_Swig.i index edbdb3770..be50f0f19 100644 --- a/src/MEDCoupling_Swig/libMEDCoupling_Swig.i +++ b/src/MEDCoupling_Swig/libMEDCoupling_Swig.i @@ -56,6 +56,7 @@ using namespace INTERP_KERNEL; %newobject ParaMEDMEM::DataArrayDouble::convertToIntArr; %newobject ParaMEDMEM::DataArrayInt::convertToDblArr; %newobject ParaMEDMEM::MEDCouplingUMesh::New; +%newobject ParaMEDMEM::MEDCouplingField::buildWeightingField; %newobject ParaMEDMEM::MEDCouplingFieldDouble::New; %newobject ParaMEDMEM::MEDCouplingFieldDouble::mergeFields; %newobject ParaMEDMEM::MEDCouplingFieldDouble::operator+; @@ -686,6 +687,7 @@ namespace ParaMEDMEM void setDescription(const char *desc); const char *getName() const; TypeOfField getTypeOfField() const; + MEDCouplingFieldDouble *buildWeightingField(bool isAbs) const throw(INTERP_KERNEL::Exception); MEDCouplingFieldDiscretization *getDiscretization() const; void setGaussLocalizationOnType(INTERP_KERNEL::NormalizedCellType type, const std::vector& refCoo, const std::vector& gsCoo, const std::vector& wg) throw(INTERP_KERNEL::Exception); @@ -772,8 +774,8 @@ namespace ParaMEDMEM double accumulate(int compId) const throw(INTERP_KERNEL::Exception); double getMaxValue() const throw(INTERP_KERNEL::Exception); double integral(int compId, bool isWAbs) const throw(INTERP_KERNEL::Exception); - double normL1(int compId, bool isWAbs) const throw(INTERP_KERNEL::Exception); - double normL2(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); static MEDCouplingFieldDouble *mergeFields(const MEDCouplingFieldDouble *f1, const MEDCouplingFieldDouble *f2) throw(INTERP_KERNEL::Exception); MEDCouplingFieldDouble *operator+(const MEDCouplingFieldDouble& other) const throw(INTERP_KERNEL::Exception); const MEDCouplingFieldDouble &operator+=(const MEDCouplingFieldDouble& other) throw(INTERP_KERNEL::Exception); @@ -889,20 +891,20 @@ namespace ParaMEDMEM delete [] tmp; return ret; } - PyObject *normL1(bool isWAbs) const throw(INTERP_KERNEL::Exception) + PyObject *normL1() const throw(INTERP_KERNEL::Exception) { int sz=self->getNumberOfTuples(); double *tmp=new double[sz]; - self->normL1(isWAbs,tmp); + self->normL1(tmp); PyObject *ret=convertDblArrToPyList(tmp,sz); delete [] tmp; return ret; } - PyObject *normL2(bool isWAbs) const throw(INTERP_KERNEL::Exception) + PyObject *normL2() const throw(INTERP_KERNEL::Exception) { int sz=self->getNumberOfTuples(); double *tmp=new double[sz]; - self->normL2(isWAbs,tmp); + self->normL2(tmp); PyObject *ret=convertDblArrToPyList(tmp,sz); delete [] tmp; return ret;