FieldDouble::multiply deals with multiplication betwenn 1-comp n-comp fields.

diff --git a/src/MEDCoupling/MEDCouplingFieldDouble.cxx b/src/MEDCoupling/MEDCouplingFieldDouble.cxx
index 41b9514a160ea684bb58e08e053e61bf4c31dbaf..ed4e7746266401f914febc6f4c3dbff83f208aaa 100644 (file)
--- a/src/MEDCoupling/MEDCouplingFieldDouble.cxx
+++ b/src/MEDCoupling/MEDCouplingFieldDouble.cxx
@@ -74,6 +74,20 @@ bool MEDCouplingFieldDouble::areCompatible(const MEDCouplingField *other) const
   return true;
 }
 
+bool MEDCouplingFieldDouble::areCompatibleForMul(const MEDCouplingField *other) const
+{
+  if(!MEDCouplingField::areCompatible(other))
+    return false;
+  const MEDCouplingFieldDouble *otherC=dynamic_cast<const MEDCouplingFieldDouble *>(other);
+  if(!otherC)
+    return false;
+  if(_nature!=otherC->_nature)
+    return false;
+  if(!_time_discr->areCompatibleForMul(otherC->_time_discr))
+    return false;
+  return true;
+}
+
 TypeOfTimeDiscretization MEDCouplingFieldDouble::getTimeDiscretization() const
 {
   return _time_discr->getEnum();
@@ -414,7 +428,7 @@ MEDCouplingFieldDouble *MEDCouplingFieldDouble::substractFields(const MEDCouplin
 
 MEDCouplingFieldDouble *MEDCouplingFieldDouble::multiplyFields(const MEDCouplingFieldDouble *f1, const MEDCouplingFieldDouble *f2)
 {
-  if(!f1->areCompatible(f2))
+  if(!f1->areCompatibleForMul(f2))
     throw INTERP_KERNEL::Exception("Fields are not compatible ; unable to applymultiplyFields  on them !");
   MEDCouplingTimeDiscretization *td=f1->_time_discr->multiply(f2->_time_discr);
   MEDCouplingFieldDouble *ret=new MEDCouplingFieldDouble(f1->getNature(),td,f1->getTypeOfField());

diff --git a/src/MEDCoupling/MEDCouplingFieldDouble.hxx b/src/MEDCoupling/MEDCouplingFieldDouble.hxx
index 9c3525dcfb1d702e66dd634df79e0cfdc7db5ae0..8e92075ac73b5a0d483b89cf5b6714b9c7bba2ba 100644 (file)
--- a/src/MEDCoupling/MEDCouplingFieldDouble.hxx
+++ b/src/MEDCoupling/MEDCouplingFieldDouble.hxx
@@ -33,6 +33,7 @@ namespace ParaMEDMEM
     static MEDCouplingFieldDouble *New(TypeOfField type, TypeOfTimeDiscretization td=NO_TIME);
     bool isEqual(const MEDCouplingField *other, double meshPrec, double valsPrec) const;
     bool areCompatible(const MEDCouplingField *other) const;
+    bool areCompatibleForMul(const MEDCouplingField *other) const;
     MEDCouplingFieldDouble *clone(bool recDeepCpy) const;
     MEDCouplingFieldDouble *buildNewTimeReprFromThis(TypeOfTimeDiscretization td, bool deepCpy) const;
     TypeOfTimeDiscretization getTimeDiscretization() const;

diff --git a/src/MEDCoupling/MEDCouplingMemArray.cxx b/src/MEDCoupling/MEDCouplingMemArray.cxx
index f3e356477d1a5f7c7478ff77ee3633572ad423b5..71d65c10899c9652c0092d4d2d536cd03fc0ba7e 100644 (file)
--- a/src/MEDCoupling/MEDCouplingMemArray.cxx
+++ b/src/MEDCoupling/MEDCouplingMemArray.cxx
@@ -163,15 +163,46 @@ DataArrayDouble *DataArrayDouble::substract(const DataArrayDouble *a1, const Dat
 
 DataArrayDouble *DataArrayDouble::multiply(const DataArrayDouble *a1, const DataArrayDouble *a2)
 {
-  int nbOfComp=a1->getNumberOfComponents();
-  if(nbOfComp!=a2->getNumberOfComponents())
-    throw INTERP_KERNEL::Exception("Nb of components mismatch for array multiply !");
   int nbOfTuple=a1->getNumberOfTuples();
-  if(nbOfTuple!=a2->getNumberOfTuples())
+  int nbOfTuple2=a2->getNumberOfTuples();
+  int nbOfComp=a1->getNumberOfComponents();
+  int nbOfComp2=a2->getNumberOfComponents();
+  if(nbOfTuple!=nbOfTuple2)
     throw INTERP_KERNEL::Exception("Nb of tuples mismatch for array multiply !");
-  DataArrayDouble *ret=DataArrayDouble::New();
-  ret->alloc(nbOfTuple,nbOfComp);
-  std::transform(a1->getConstPointer(),a1->getConstPointer()+nbOfTuple*nbOfComp,a2->getConstPointer(),ret->getPointer(),std::multiplies<double>());
+  DataArrayDouble *ret=0;
+  if(nbOfComp==nbOfComp2)
+    {
+      ret=DataArrayDouble::New();
+      ret->alloc(nbOfTuple,nbOfComp);
+      std::transform(a1->getConstPointer(),a1->getConstPointer()+nbOfTuple*nbOfComp,a2->getConstPointer(),ret->getPointer(),std::multiplies<double>());
+    }
+  else
+    {
+      int nbOfCompMin,nbOfCompMax;
+      const DataArrayDouble *aMin, *aMax;
+      if(nbOfComp>nbOfComp2)
+        {
+          nbOfCompMin=nbOfComp2; nbOfCompMax=nbOfComp;
+          aMin=a2; aMax=a1;
+        }
+      else
+        {
+          nbOfCompMin=nbOfComp; nbOfCompMax=nbOfComp2;
+          aMin=a1; aMax=a2;
+        }
+      if(nbOfCompMin==1)
+        {
+          ret=DataArrayDouble::New();
+          ret->alloc(nbOfTuple,nbOfCompMax);
+          const double *aMinPtr=aMin->getConstPointer();
+          const double *aMaxPtr=aMax->getConstPointer();
+          double *res=ret->getPointer();
+          for(int i=0;i<nbOfTuple;i++)
+            res=std::transform(aMaxPtr+i*nbOfCompMax,aMaxPtr+(i+1)*nbOfCompMax,res,std::bind2nd(std::multiplies<double>(),aMinPtr[i]));
+        }
+      else
+        throw INTERP_KERNEL::Exception("Nb of components mismatch for array multiply !");
+    }
   ret->copyStringInfoFrom(*a1);
   return ret;
 }

diff --git a/src/MEDCoupling/MEDCouplingTimeDiscretization.cxx b/src/MEDCoupling/MEDCouplingTimeDiscretization.cxx
index b163086bf7baf76ec8d46ec7469e30d96dfaab54..43ecbdca0d318f6c068adb78110595913d303960 100644 (file)
--- a/src/MEDCoupling/MEDCouplingTimeDiscretization.cxx
+++ b/src/MEDCoupling/MEDCouplingTimeDiscretization.cxx
@@ -68,6 +68,22 @@ bool MEDCouplingTimeDiscretization::areCompatible(const MEDCouplingTimeDiscretiz
   return true;
 }
 
+bool MEDCouplingTimeDiscretization::areCompatibleForMul(const MEDCouplingTimeDiscretization *other) const
+{
+  if(std::fabs(_time_tolerance-other->_time_tolerance)>1.e-16)
+    return false;
+  if(_array==0 && other->_array==0)
+    return true;
+  if(_array==0 || other->_array==0)
+    return false;
+  int nbC1=_array->getNumberOfComponents();
+  int nbC2=other->_array->getNumberOfComponents();
+  int nbMin=std::min(nbC1,nbC2);
+  if(nbC1!=nbC2 && nbMin!=1)
+    return false;
+  return true;
+}
+
 bool MEDCouplingTimeDiscretization::isEqual(const MEDCouplingTimeDiscretization *other, double prec) const
 {
   if(!areCompatible(other))
@@ -322,6 +338,14 @@ bool MEDCouplingNoTimeLabel::areCompatible(const MEDCouplingTimeDiscretization *
   return otherC!=0;
 }
 
+bool MEDCouplingNoTimeLabel::areCompatibleForMul(const MEDCouplingTimeDiscretization *other) const
+{
+  if(!MEDCouplingTimeDiscretization::areCompatibleForMul(other))
+    return false;
+  const MEDCouplingNoTimeLabel *otherC=dynamic_cast<const MEDCouplingNoTimeLabel *>(other);
+  return otherC!=0;
+}
+
 bool MEDCouplingNoTimeLabel::isEqual(const MEDCouplingTimeDiscretization *other, double prec) const
 {
   const MEDCouplingNoTimeLabel *otherC=dynamic_cast<const MEDCouplingNoTimeLabel *>(other);
@@ -491,6 +515,16 @@ bool MEDCouplingWithTimeStep::areCompatible(const MEDCouplingTimeDiscretization
   return std::fabs(_time-otherC->_time)<_time_tolerance;
 }
 
+bool MEDCouplingWithTimeStep::areCompatibleForMul(const MEDCouplingTimeDiscretization *other) const
+{
+  if(!MEDCouplingTimeDiscretization::areCompatibleForMul(other))
+    return false;
+  const MEDCouplingWithTimeStep *otherC=dynamic_cast<const MEDCouplingWithTimeStep *>(other);
+  if(!otherC)
+    return false;
+  return std::fabs(_time-otherC->_time)<_time_tolerance;
+}
+
 bool MEDCouplingWithTimeStep::isEqual(const MEDCouplingTimeDiscretization *other, double prec) const
 {
   const MEDCouplingWithTimeStep *otherC=dynamic_cast<const MEDCouplingWithTimeStep *>(other);
@@ -709,6 +743,16 @@ bool MEDCouplingConstOnTimeInterval::areCompatible(const MEDCouplingTimeDiscreti
   return (std::fabs(_start_time-otherC->_start_time)<_time_tolerance && std::fabs(_end_time-otherC->_end_time)<_time_tolerance);
 }
 
+bool MEDCouplingConstOnTimeInterval::areCompatibleForMul(const MEDCouplingTimeDiscretization *other) const
+{
+  if(!MEDCouplingTimeDiscretization::areCompatible(other))
+    return false;
+  const MEDCouplingConstOnTimeInterval *otherC=dynamic_cast<const MEDCouplingConstOnTimeInterval *>(other);
+  if(!otherC)
+    return false;
+  return (std::fabs(_start_time-otherC->_start_time)<_time_tolerance && std::fabs(_end_time-otherC->_end_time)<_time_tolerance);
+}
+
 bool MEDCouplingConstOnTimeInterval::isEqual(const MEDCouplingTimeDiscretization *other, double prec) const
 {
   const MEDCouplingConstOnTimeInterval *otherC=dynamic_cast<const MEDCouplingConstOnTimeInterval *>(other);

diff --git a/src/MEDCoupling/MEDCouplingTimeDiscretization.hxx b/src/MEDCoupling/MEDCouplingTimeDiscretization.hxx
index f97924786c61b4ede6104e59334270336a98e0e7..e13742ee04cd170d6c3ff82a82f2ac9dad9da999 100644 (file)
--- a/src/MEDCoupling/MEDCouplingTimeDiscretization.hxx
+++ b/src/MEDCoupling/MEDCouplingTimeDiscretization.hxx
@@ -40,6 +40,7 @@ namespace ParaMEDMEM
     void updateTime();
     static MEDCouplingTimeDiscretization *New(TypeOfTimeDiscretization type);
     virtual bool areCompatible(const MEDCouplingTimeDiscretization *other) const;
+    virtual bool areCompatibleForMul(const MEDCouplingTimeDiscretization *other) const;
     virtual bool isEqual(const MEDCouplingTimeDiscretization *other, double prec) const;
     virtual MEDCouplingTimeDiscretization *buildNewTimeReprFromThis(const MEDCouplingTimeDiscretization *other,
                                                                     TypeOfTimeDiscretization type, bool deepCpy) const;
@@ -103,6 +104,7 @@ namespace ParaMEDMEM
     MEDCouplingTimeDiscretization *divide(const MEDCouplingTimeDiscretization *other) const;
     bool isEqual(const MEDCouplingTimeDiscretization *other, double prec) const;
     bool areCompatible(const MEDCouplingTimeDiscretization *other) const;
+    bool areCompatibleForMul(const MEDCouplingTimeDiscretization *other) const;
     MEDCouplingTimeDiscretization *performCpy(bool deepCpy) const;
     void checkNoTimePresence() const throw(INTERP_KERNEL::Exception) { }
     void checkTimePresence(double time) const throw(INTERP_KERNEL::Exception);
@@ -136,6 +138,7 @@ namespace ParaMEDMEM
     MEDCouplingTimeDiscretization *divide(const MEDCouplingTimeDiscretization *other) const;
     bool isEqual(const MEDCouplingTimeDiscretization *other, double prec) const;
     bool areCompatible(const MEDCouplingTimeDiscretization *other) const;
+    bool areCompatibleForMul(const MEDCouplingTimeDiscretization *other) const;
     void getTinySerializationIntInformation(std::vector<int>& tinyInfo) const;
     void getTinySerializationDbleInformation(std::vector<double>& tinyInfo) const;
     void finishUnserialization(const std::vector<int>& tinyInfoI, const std::vector<double>& tinyInfoD, const std::vector<std::string>& tinyInfoS);
@@ -172,6 +175,7 @@ namespace ParaMEDMEM
     void finishUnserialization(const std::vector<int>& tinyInfoI, const std::vector<double>& tinyInfoD, const std::vector<std::string>& tinyInfoS);
     MEDCouplingTimeDiscretization *performCpy(bool deepCpy) const;
     bool areCompatible(const MEDCouplingTimeDiscretization *other) const;
+    bool areCompatibleForMul(const MEDCouplingTimeDiscretization *other) const;
     bool isEqual(const MEDCouplingTimeDiscretization *other, double prec) const;
     std::vector< const DataArrayDouble *> getArraysForTime(double time) const throw(INTERP_KERNEL::Exception);
     void getValueForTime(double time, const std::vector<double>& vals, double *res) const;