From: Anthony Geay <anthony.geay@edf.fr>
Date: Fri, 25 Jul 2014 11:30:53 +0000 (+0200)
Subject: Boost of expression evaluator DataArrayDouble::applyFunc* + DataArrayDouble::applyFun... 
X-Git-Tag: V7_5_0a1~2^2~13
X-Git-Url: http://git.salome-platform.org/gitweb/?a=commitdiff_plain;h=5a8dfa8819c5f0ed754b0d0dc51dc51617b3e710;p=tools%2Fmedcoupling.git

Boost of expression evaluator DataArrayDouble::applyFunc* + DataArrayDouble::applyFuncOnThis
---

diff --git a/src/INTERP_KERNEL/Bases/InterpKernelAutoPtr.hxx b/src/INTERP_KERNEL/Bases/InterpKernelAutoPtr.hxx
index 1c805553a..312633862 100644
--- a/src/INTERP_KERNEL/Bases/InterpKernelAutoPtr.hxx
+++ b/src/INTERP_KERNEL/Bases/InterpKernelAutoPtr.hxx
@@ -29,7 +29,7 @@ namespace INTERP_KERNEL
   public:
     AutoPtr(T *ptr=0):_ptr(ptr) {  }
     ~AutoPtr() { destroyPtr(); }
-    AutoPtr &operator=(T *ptr) { destroyPtr(); _ptr=ptr; return *this; }
+    AutoPtr &operator=(T *ptr) { if(_ptr!=ptr) { destroyPtr(); _ptr=ptr; } return *this; }
     T *operator->() { return _ptr ; }
     const T *operator->() const { return _ptr; }
     T& operator*() { return *_ptr; }
@@ -42,13 +42,32 @@ namespace INTERP_KERNEL
     T *_ptr;
   };
 
+  template<class T>
+  class AutoCppPtr
+  {
+  public:
+    AutoCppPtr(T *ptr=0):_ptr(ptr) {  }
+    ~AutoCppPtr() { destroyPtr(); }
+    AutoCppPtr &operator=(T *ptr) { if(_ptr!=ptr) { destroyPtr(); _ptr=ptr; } return *this; }
+    T *operator->() { return _ptr ; }
+    const T *operator->() const { return _ptr; }
+    T& operator*() { return *_ptr; }
+    const T& operator*() const { return *_ptr; }
+    operator T *() { return _ptr; }
+    operator const T *() const { return _ptr; }
+  private:
+    void destroyPtr() { delete _ptr; }
+  private:
+    T *_ptr;
+  };
+
   template<class T>
   class AutoCPtr
   {
   public:
     AutoCPtr(T *ptr=0):_ptr(ptr) {  }
     ~AutoCPtr() { destroyPtr(); }
-    AutoCPtr &operator=(T *ptr) { destroyPtr(); _ptr=ptr; return *this; }
+    AutoCPtr &operator=(T *ptr) { if(_ptr!=ptr) { destroyPtr(); _ptr=ptr; } return *this; }
     T *operator->() { return _ptr ; }
     const T *operator->() const { return _ptr; }
     T& operator*() { return *_ptr; }
diff --git a/src/INTERP_KERNEL/ExprEval/InterpKernelExprParser.cxx b/src/INTERP_KERNEL/ExprEval/InterpKernelExprParser.cxx
index 960d39108..c82cdfe3f 100644
--- a/src/INTERP_KERNEL/ExprEval/InterpKernelExprParser.cxx
+++ b/src/INTERP_KERNEL/ExprEval/InterpKernelExprParser.cxx
@@ -21,6 +21,7 @@
 #include "InterpKernelExprParser.hxx"
 #include "InterpKernelValue.hxx"
 #include "InterpKernelAsmX86.hxx"
+#include "InterpKernelAutoPtr.hxx"
 
 #include <cctype>
 #include <sstream>
@@ -84,13 +85,22 @@ void LeafExprVal::replaceValues(const std::vector<double>& valuesInExpr)
   _value=valuesInExpr[pos];
 }
 
-LeafExprVar::LeafExprVar(const std::string& var):_fast_pos(-1),_var_name(var)
+LeafExprVal *LeafExprVal::deepCpy() const
+{
+  return new LeafExprVal(*this);
+}
+
+LeafExprVar::LeafExprVar(const std::string& var):_fast_pos(-1),_var_name(var),_val(0)
 {
 }
 
 void LeafExprVar::fillValue(Value *val) const
 {
-  val->setVarname(_fast_pos,_var_name);
+  if(_val)
+    val->setDouble(_val[_fast_pos]);
+  else
+    val->setVarname(_fast_pos,_var_name);
+
 }
 
 void LeafExprVar::prepareExprEvaluation(const std::vector<std::string>& vars, int nbOfCompo, int targetNbOfCompo) const
@@ -125,6 +135,22 @@ void LeafExprVar::prepareExprEvaluation(const std::vector<std::string>& vars, in
     }
 }
 
+/*!
+ * \param [in] vars - the sorted list of vars
+ * \param [in] nbOfCompo - the size of the input tuples (it is used to scan if no problem occurs)
+ * \param [in] targetNbOfCompo - the size of the output tuple (it is used to check that no problem occurs)
+ * \param [in] refPos - is an integer in [0,targetNbOfCompo), that tell the id of \a this. It is for multi interpreters.
+ * \sa evaluateDouble
+ */
+void LeafExprVar::prepareExprEvaluationDouble(const std::vector<std::string>& vars, int nbOfCompo, int targetNbOfCompo, int refPos, const double *ptOfInputStart, const double *ptOfInputEnd) const
+{
+  if((int)vars.size()!=std::distance(ptOfInputStart,ptOfInputEnd))
+    throw INTERP_KERNEL::Exception("LeafExprVar::prepareExprEvaluationDouble : size of input vector must be equal to the input vector !");
+  prepareExprEvaluation(vars,nbOfCompo,targetNbOfCompo);
+  _ref_pos=refPos;
+  _val=ptOfInputStart;
+}
+
 void LeafExprVar::prepareExprEvaluationVec() const
 {
   if(!isRecognizedKeyVar(_var_name,_fast_pos))
@@ -145,6 +171,11 @@ bool LeafExprVar::isRecognizedKeyVar(const std::string& var, int& pos)
   return true;
 }
 
+LeafExprVar *LeafExprVar::deepCpy() const
+{
+  return new LeafExprVar(*this);
+}
+
 /*!
  * Nothing to do it is not a bug.
  */
@@ -156,6 +187,26 @@ LeafExprVar::~LeafExprVar()
 {
 }
 
+void ExprParserOfEval::clearSortedMemory()
+{
+  delete _leaf;
+  for(std::vector<ExprParserOfEval>::iterator it=_sub_parts.begin();it!=_sub_parts.end();it++)
+    (*it).clearSortedMemory();
+  for(std::vector<Function *>::iterator it=_funcs.begin();it!=_funcs.end();it++)
+    delete *it;
+}
+
+void ExprParserOfEval::sortMemory()
+{
+  for(std::vector<ExprParserOfEval>::iterator it=_sub_parts.begin();it!=_sub_parts.end();it++)
+    (*it).sortMemory();
+  if(_leaf)
+    _leaf=_leaf->deepCpy();
+  for(std::vector<Function *>::iterator it=_funcs.begin();it!=_funcs.end();it++)
+    if(*it)
+      *it=(*it)->deepCpy();
+}
+
 ExprParser::ExprParser(const std::string& expr, ExprParser *father):_father(father),_is_parsed(false),_leaf(0),_is_parsing_ok(false),_expr(expr)
 {
   _expr=deleteWhiteSpaces(_expr);
@@ -171,6 +222,7 @@ ExprParser::ExprParser(const char *expr, int lgth, ExprParser *father):_father(f
 ExprParser::~ExprParser()
 {
   delete _leaf;
+  _for_eval.clearSortedMemory();
   releaseFunctions();
 }
 
@@ -242,17 +294,15 @@ void ExprParser::parse()
       replaceValues(valuesInExpr);
       _expr=tmp;
     }
+  reverseThis();
   _is_parsing_ok=true;
 }
 
 double ExprParser::evaluate() const
 {
-  Value *gen=new ValueDouble;
-  ValueDouble *res=(ValueDouble *)evaluateLowLev(gen);
-  delete gen;
-  double ret=res->getData();
-  delete res;
-  return ret;
+  AutoCppPtr<Value> gen(new ValueDouble);
+  AutoCppPtr<ValueDouble> res(static_cast<ValueDouble *>(evaluateLowLev(gen)));
+  return res->getData();
 }
 
 DecompositionInUnitBase ExprParser::evaluateUnit() const
@@ -276,20 +326,9 @@ DecompositionInUnitBase ExprParser::evaluateUnit() const
 
 void ExprParser::evaluateExpr(int szOfOutParam, const double *inParam, double *outParam) const
 {
-  Value *gen=new ValueDoubleExpr(szOfOutParam,inParam);
-  ValueDoubleExpr *res=0;
-  try
-    {
-      res=(ValueDoubleExpr *)evaluateLowLev(gen);
-    }
-  catch(INTERP_KERNEL::Exception& e)
-    {
-      delete gen;
-      throw e;
-    }
-  delete gen;
+  AutoCppPtr<Value> gen(new ValueDoubleExpr(szOfOutParam,inParam));
+  AutoCppPtr<ValueDoubleExpr> res(static_cast<ValueDoubleExpr *>(evaluateLowLev(gen)));
   std::copy(res->getData(),res->getData()+szOfOutParam,outParam);
-  delete res;
 }
 
 void ExprParser::prepareExprEvaluation(const std::vector<std::string>& vars, int nbOfCompo, int targetNbOfCompo) const
@@ -301,10 +340,58 @@ void ExprParser::prepareExprEvaluation(const std::vector<std::string>& vars, int
         leafC->prepareExprEvaluation(vars,nbOfCompo,targetNbOfCompo);
     }
   else
-    for(std::list<ExprParser>::const_iterator iter=_sub_expr.begin();iter!=_sub_expr.end();iter++)
+    for(std::vector<ExprParser>::const_iterator iter=_sub_expr.begin();iter!=_sub_expr.end();iter++)
       (*iter).prepareExprEvaluation(vars,nbOfCompo,targetNbOfCompo);
 }
 
+/*!
+ * \param [in] vars - the sorted list of vars
+ * \param [in] nbOfCompo - the size of the input tuples (it is used to scan if no problem occurs)
+ * \param [in] targetNbOfCompo - the size of the output tuple (it is used to check that no problem occurs)
+ * \param [in] refPos - is an integer in [0,targetNbOfCompo), that tell the id of \a this. It is for multi interpreters.
+ * \sa evaluateDouble
+ */
+void ExprParser::prepareExprEvaluationDouble(const std::vector<std::string>& vars, int nbOfCompo, int targetNbOfCompo, int refPos, const double *ptOfInputStart, const double *ptOfInputEnd) const
+{
+  if((int)vars.size()!=std::distance(ptOfInputStart,ptOfInputEnd))
+    throw INTERP_KERNEL::Exception("ExprParser::prepareExprEvaluationDouble : size of input vector must be equal to the input vector !");
+  if(_leaf)
+    {
+      LeafExprVar *leafC=dynamic_cast<LeafExprVar *>(_leaf);
+      if(leafC)
+        leafC->prepareExprEvaluationDouble(vars,nbOfCompo,targetNbOfCompo,refPos,ptOfInputStart,ptOfInputEnd);
+    }
+  else
+    for(std::vector<ExprParser>::const_iterator iter=_sub_expr.begin();iter!=_sub_expr.end();iter++)
+      (*iter).prepareExprEvaluationDouble(vars,nbOfCompo,targetNbOfCompo,refPos,ptOfInputStart,ptOfInputEnd);
+}
+
+void ExprParser::prepareFastEvaluator() const
+{
+  _for_eval.clearSortedMemory();
+  _for_eval=convertMeTo();
+  _for_eval.sortMemory();
+}
+
+/*!
+ * \sa prepareExprEvaluationDouble
+ */
+double ExprParser::evaluateDouble() const
+{
+  checkForEvaluation();
+  std::vector<double> stackOfVal;
+  evaluateDoubleInternal(stackOfVal);
+  return stackOfVal.back();
+}
+
+void ExprParser::checkForEvaluation() const
+{
+  if(!_is_parsing_ok)
+    throw INTERP_KERNEL::Exception("checkForEvaluation : Parsing fails ! Invalid expression !");
+  if(_sub_expr.empty() && !_leaf)
+    throw INTERP_KERNEL::Exception("checkForEvaluation : Empty expression !");
+}
+
 void ExprParser::prepareExprEvaluationVec() const
 {
   std::set<std::string> trueVars;
@@ -328,16 +415,13 @@ void ExprParser::prepareExprEvaluationVecLowLev() const
         leafC->prepareExprEvaluationVec();
     }
   else
-    for(std::list<ExprParser>::const_iterator iter=_sub_expr.begin();iter!=_sub_expr.end();iter++)
+    for(std::vector<ExprParser>::const_iterator iter=_sub_expr.begin();iter!=_sub_expr.end();iter++)
       (*iter).prepareExprEvaluationVecLowLev();
 }
 
 Value *ExprParser::evaluateLowLev(Value *valGen) const
 {
-  if(!_is_parsing_ok)
-    throw INTERP_KERNEL::Exception("Parsing fails ! Invalid expression !");
-  if(_sub_expr.empty() && !_leaf)
-    throw INTERP_KERNEL::Exception("Empty expression !");
+  checkForEvaluation();
   std::vector<Value *> stackOfVal;
   try
     {
@@ -359,12 +443,11 @@ Value *ExprParser::evaluateLowLev(Value *valGen) const
       else
         {
           stackOfVal.resize(_sub_expr.size());
-          std::vector<Value *>::reverse_iterator iter2=stackOfVal.rbegin();
-          for(std::list<ExprParser>::const_iterator iter=_sub_expr.begin();iter!=_sub_expr.end();iter++,iter2++)
+          std::vector<Value *>::iterator iter2=stackOfVal.begin();
+          for(std::vector<ExprParser>::const_iterator iter=_sub_expr.begin();iter!=_sub_expr.end();iter++,iter2++)
             *iter2=(*iter).evaluateLowLev(valGen);
         }
-      std::list<Function *>::const_iterator iter3;
-      for(iter3=_func_btw_sub_expr.begin();iter3!=_func_btw_sub_expr.end();iter3++)
+      for(std::vector<Function *>::const_iterator iter3=_func_btw_sub_expr.begin();iter3!=_func_btw_sub_expr.end();iter3++)
         (*iter3)->operate(stackOfVal);
     }
   catch(INTERP_KERNEL::Exception& e)
@@ -376,6 +459,33 @@ Value *ExprParser::evaluateLowLev(Value *valGen) const
   return stackOfVal.back();
 }
 
+void ExprParser::reverseThis()
+{
+  if(_leaf)
+    return ;
+  for(std::vector<ExprParser>::iterator iter=_sub_expr.begin();iter!=_sub_expr.end();iter++)
+    (*iter).reverseThis();
+  AutoPtr<char> buf(new char[sizeof(ExprParser)]);
+  char *loc(reinterpret_cast<char *>(&_sub_expr[0])),*bufPtr(buf);
+  std::size_t sz(_sub_expr.size());
+  std::size_t nbOfTurn(sz/2);
+  for(std::size_t i=0;i<nbOfTurn;i++)
+    {
+      std::copy(loc+i*sizeof(ExprParser),loc+(i+1)*sizeof(ExprParser),bufPtr);
+      std::copy(loc+(sz-i-1)*sizeof(ExprParser),loc+(sz-i)*sizeof(ExprParser),loc+i*sizeof(ExprParser));
+      std::copy(bufPtr,bufPtr+sizeof(ExprParser),loc+(sz-i-1)*sizeof(ExprParser));
+    }
+}
+
+ExprParserOfEval ExprParser::convertMeTo() const
+{
+  std::size_t sz(_sub_expr.size());
+  std::vector<ExprParserOfEval> subExpr(sz);
+  for(std::size_t i=0;i<sz;i++)
+    subExpr[i]=_sub_expr[i].convertMeTo();
+  return ExprParserOfEval(_leaf,subExpr,_func_btw_sub_expr);
+}
+
 void ExprParser::getSetOfVars(std::set<std::string>& vars) const
 {
   if(_leaf)
@@ -385,7 +495,7 @@ void ExprParser::getSetOfVars(std::set<std::string>& vars) const
         vars.insert(leafC->getVar());
     }
   else
-    for(std::list<ExprParser>::const_iterator iter=_sub_expr.begin();iter!=_sub_expr.end();iter++)
+    for(std::vector<ExprParser>::const_iterator iter=_sub_expr.begin();iter!=_sub_expr.end();iter++)
       (*iter).getSetOfVars(vars);
 }
 
@@ -404,7 +514,7 @@ void ExprParser::getTrueSetOfVars(std::set<std::string>& trueVars) const
 
 void ExprParser::parseDeeper()
 {
-  for(std::list<ExprParser>::iterator iter=_sub_expr.begin();iter!=_sub_expr.end();iter++)
+  for(std::vector<ExprParser>::iterator iter=_sub_expr.begin();iter!=_sub_expr.end();iter++)
     if(!(*iter).simplify())
       (*iter).parseDeeper();
 }
@@ -727,7 +837,7 @@ void ExprParser::parseForPow()
 
 void ExprParser::releaseFunctions()
 {
-  for(std::list<Function *>::iterator iter=_func_btw_sub_expr.begin();iter!=_func_btw_sub_expr.end();iter++)
+  for(std::vector<Function *>::iterator iter=_func_btw_sub_expr.begin();iter!=_func_btw_sub_expr.end();iter++)
     delete *iter;
   _func_btw_sub_expr.clear();
 }
@@ -924,7 +1034,7 @@ void ExprParser::replaceValues(const std::vector<double>& valuesInExpr)
     _leaf->replaceValues(valuesInExpr);
   else
     {
-      for(std::list<ExprParser>::iterator iter=_sub_expr.begin();iter!=_sub_expr.end();iter++)
+      for(std::vector<ExprParser>::iterator iter=_sub_expr.begin();iter!=_sub_expr.end();iter++)
         (*iter).replaceValues(valuesInExpr);
     }
 }
@@ -986,10 +1096,10 @@ void ExprParser::compileX86LowLev(std::vector<std::string>& ass) const
     _leaf->compileX86(ass);
   else
     {
-      for(std::list<ExprParser>::const_iterator iter=_sub_expr.begin();iter!=_sub_expr.end();iter++)
+      for(std::vector<ExprParser>::const_iterator iter=_sub_expr.begin();iter!=_sub_expr.end();iter++)
         (*iter).compileX86LowLev(ass);
     }
-  for(std::list<Function *>::const_iterator iter2=_func_btw_sub_expr.begin();iter2!=_func_btw_sub_expr.end();iter2++)
+  for(std::vector<Function *>::const_iterator iter2=_func_btw_sub_expr.begin();iter2!=_func_btw_sub_expr.end();iter2++)
     (*iter2)->operateX86(ass);
 }
 
@@ -999,13 +1109,18 @@ void ExprParser::compileX86_64LowLev(std::vector<std::string>& ass) const
     _leaf->compileX86_64(ass);
   else
     {
-      for(std::list<ExprParser>::const_iterator iter=_sub_expr.begin();iter!=_sub_expr.end();iter++)
+      for(std::vector<ExprParser>::const_iterator iter=_sub_expr.begin();iter!=_sub_expr.end();iter++)
         (*iter).compileX86_64LowLev(ass);
     }
-  for(std::list<Function *>::const_iterator iter2=_func_btw_sub_expr.begin();iter2!=_func_btw_sub_expr.end();iter2++)
+  for(std::vector<Function *>::const_iterator iter2=_func_btw_sub_expr.begin();iter2!=_func_btw_sub_expr.end();iter2++)
     (*iter2)->operateX86(ass);
 }
 
+double LeafExprVal::getDoubleValue() const
+{
+  return _value;
+}
+
 void LeafExprVal::compileX86(std::vector<std::string>& ass) const
 {
   ass.push_back("sub esp,8");
@@ -1038,6 +1153,17 @@ void LeafExprVal::compileX86_64(std::vector<std::string>& ass) const
   ass.push_back("add rsp,8");
 }
 
+double LeafExprVar::getDoubleValue() const
+{
+  if(_fast_pos>=0)
+    return _val[_fast_pos];
+  else
+    {
+      int pos(-7-_fast_pos);
+      return pos==_ref_pos?1.:0.;
+    }
+}
+
 void LeafExprVar::compileX86(std::vector<std::string>& ass) const
 {
   ass.push_back("fld qword [ebp+8]");
@@ -1057,7 +1183,7 @@ int ExprParser::getStackSizeToPlayX86(const ExprParser *asker) const
     {
       int sz=_father->getStackSizeToPlayX86(this);
       int i=0;
-      for(std::list<ExprParser>::const_reverse_iterator iter=_sub_expr.rbegin();iter!=_sub_expr.rend();iter++,i++)
+      for(std::vector<ExprParser>::const_reverse_iterator iter=_sub_expr.rbegin();iter!=_sub_expr.rend();iter++,i++)
         {
           const ExprParser& obj=(*iter);
           const ExprParser *pt=&obj;
diff --git a/src/INTERP_KERNEL/ExprEval/InterpKernelExprParser.hxx b/src/INTERP_KERNEL/ExprEval/InterpKernelExprParser.hxx
index 0d1ed3e0c..6a03b32cf 100644
--- a/src/INTERP_KERNEL/ExprEval/InterpKernelExprParser.hxx
+++ b/src/INTERP_KERNEL/ExprEval/InterpKernelExprParser.hxx
@@ -39,10 +39,12 @@ namespace INTERP_KERNEL
   {
   public:
     INTERPKERNEL_EXPORT virtual ~LeafExpr();
+    INTERPKERNEL_EXPORT virtual double getDoubleValue() const = 0;
     INTERPKERNEL_EXPORT virtual void fillValue(Value *val) const = 0;
     INTERPKERNEL_EXPORT virtual void compileX86(std::vector<std::string>& ass) const = 0;
     INTERPKERNEL_EXPORT virtual void compileX86_64(std::vector<std::string>& ass) const = 0;
     INTERPKERNEL_EXPORT virtual void replaceValues(const std::vector<double>& valuesInExpr) = 0;
+    INTERPKERNEL_EXPORT virtual LeafExpr *deepCpy() const = 0;
     INTERPKERNEL_EXPORT static LeafExpr *buildInstanceFrom(const std::string& expr);
   };
 
@@ -51,10 +53,12 @@ namespace INTERP_KERNEL
   public:
     INTERPKERNEL_EXPORT LeafExprVal(double value);
     INTERPKERNEL_EXPORT ~LeafExprVal();
+    INTERPKERNEL_EXPORT double getDoubleValue() const;
     INTERPKERNEL_EXPORT void compileX86(std::vector<std::string>& ass) const;
     INTERPKERNEL_EXPORT void compileX86_64(std::vector<std::string>& ass) const;
     INTERPKERNEL_EXPORT void fillValue(Value *val) const;
     INTERPKERNEL_EXPORT void replaceValues(const std::vector<double>& valuesInExpr);
+    INTERPKERNEL_EXPORT LeafExprVal *deepCpy() const;
   private:
     double _value;
   };
@@ -62,21 +66,60 @@ namespace INTERP_KERNEL
   class LeafExprVar : public LeafExpr
   {
   public:
+    INTERPKERNEL_EXPORT LeafExprVar(const LeafExprVar& other):_fast_pos(other._fast_pos),_ref_pos(other._ref_pos),_var_name(other._var_name),_val(other._val) { }
     INTERPKERNEL_EXPORT LeafExprVar(const std::string& var);
     INTERPKERNEL_EXPORT ~LeafExprVar();
+    INTERPKERNEL_EXPORT double getDoubleValue() const;
     INTERPKERNEL_EXPORT void compileX86(std::vector<std::string>& ass) const;
     INTERPKERNEL_EXPORT void compileX86_64(std::vector<std::string>& ass) const;
     INTERPKERNEL_EXPORT void fillValue(Value *val) const;
     INTERPKERNEL_EXPORT std::string getVar() const { return _var_name; }
     INTERPKERNEL_EXPORT void prepareExprEvaluation(const std::vector<std::string>& vars, int nbOfCompo, int targetNbOfCompo) const;
+    INTERPKERNEL_EXPORT void prepareExprEvaluationDouble(const std::vector<std::string>& vars, int nbOfCompo, int targetNbOfCompo, int refPos, const double *ptOfInputStart, const double *ptOfInputEnd) const;
     INTERPKERNEL_EXPORT void prepareExprEvaluationVec() const;
     INTERPKERNEL_EXPORT void replaceValues(const std::vector<double>& valuesInExpr);
     INTERPKERNEL_EXPORT static bool isRecognizedKeyVar(const std::string& var, int& pos);
+    INTERPKERNEL_EXPORT LeafExprVar *deepCpy() const;
   public:
     static const char END_OF_RECOGNIZED_VAR[];
   private:
     mutable int _fast_pos;
+    mutable int _ref_pos;
     std::string _var_name;
+    mutable const double *_val;
+  };
+
+  class ExprParserOfEval
+  {
+  public:
+    ExprParserOfEval():_leaf(0) { }
+    ExprParserOfEval(LeafExpr *leaf, const std::vector<ExprParserOfEval>& subParts, const std::vector<Function *>& funcs):_leaf(leaf),_sub_parts(subParts),_funcs(funcs) { }
+    void evaluateDoubleInternal(std::vector<double>& stck) const
+    {
+      if(_leaf)
+        stck.push_back(_leaf->getDoubleValue());
+      else
+        for(std::vector<ExprParserOfEval>::const_iterator iter=_sub_parts.begin();iter!=_sub_parts.end();iter++)
+          (*iter).evaluateDoubleInternal(stck);
+      for(std::vector<Function *>::const_iterator iter3=_funcs.begin();iter3!=_funcs.end();iter3++)
+        (*iter3)->operateStackOfDouble(stck);
+    }
+    void evaluateDoubleInternalSafe(std::vector<double>& stck) const
+    {
+      if(_leaf)
+        stck.push_back(_leaf->getDoubleValue());
+      else
+        for(std::vector<ExprParserOfEval>::const_iterator iter=_sub_parts.begin();iter!=_sub_parts.end();iter++)
+          (*iter).evaluateDoubleInternalSafe(stck);
+      for(std::vector<Function *>::const_iterator iter3=_funcs.begin();iter3!=_funcs.end();iter3++)
+        (*iter3)->operateStackOfDoubleSafe(stck);
+    }
+    void clearSortedMemory();
+    void sortMemory();
+  private:
+    LeafExpr *_leaf;
+    std::vector<ExprParserOfEval> _sub_parts;
+    std::vector<Function *> _funcs;
   };
 
   class ExprParser
@@ -90,8 +133,14 @@ namespace INTERP_KERNEL
     INTERPKERNEL_EXPORT double evaluate() const;
     INTERPKERNEL_EXPORT DecompositionInUnitBase evaluateUnit() const;
     INTERPKERNEL_EXPORT void prepareExprEvaluation(const std::vector<std::string>& vars, int nbOfCompo, int targetNbOfCompo) const;
-    INTERPKERNEL_EXPORT void evaluateExpr(int szOfOutParam, const double *inParam, double *outParam) const;
+    INTERPKERNEL_EXPORT void prepareExprEvaluationDouble(const std::vector<std::string>& vars, int nbOfCompo, int targetNbOfCompo, int refPos, const double *ptOfInputStart, const double *ptOfInputEnd) const;
+    INTERPKERNEL_EXPORT void prepareFastEvaluator() const;
     INTERPKERNEL_EXPORT void prepareExprEvaluationVec() const;
+    INTERPKERNEL_EXPORT double evaluateDouble() const;
+    INTERPKERNEL_EXPORT void evaluateDoubleInternal(std::vector<double>& stck) const { _for_eval.evaluateDoubleInternal(stck); }
+    INTERPKERNEL_EXPORT void evaluateDoubleInternalSafe(std::vector<double>& stck) const { _for_eval.evaluateDoubleInternalSafe(stck); }
+    INTERPKERNEL_EXPORT void checkForEvaluation() const;
+    INTERPKERNEL_EXPORT void evaluateExpr(int szOfOutParam, const double *inParam, double *outParam) const;
     INTERPKERNEL_EXPORT void getSetOfVars(std::set<std::string>& vars) const;
     INTERPKERNEL_EXPORT void getTrueSetOfVars(std::set<std::string>& vars) const;
     //
@@ -105,6 +154,8 @@ namespace INTERP_KERNEL
     INTERPKERNEL_EXPORT static std::string deleteWhiteSpaces(const std::string& expr);
   private:
     Value *evaluateLowLev(Value *valGen) const;
+    void reverseThis();
+    ExprParserOfEval convertMeTo() const;
   private:
     void prepareExprEvaluationVecLowLev() const;
     bool tryToInterpALeaf();
@@ -128,8 +179,9 @@ namespace INTERP_KERNEL
     LeafExpr *_leaf;
     bool _is_parsing_ok;
     std::string _expr;
-    std::list<ExprParser> _sub_expr;
-    std::list<Function *> _func_btw_sub_expr;
+    mutable ExprParserOfEval _for_eval;
+    std::vector<ExprParser> _sub_expr;
+    std::vector<Function *> _func_btw_sub_expr;
   private:
     static const int MAX_X86_FP_ST=8;
     static const char WHITE_SPACES[];
diff --git a/src/INTERP_KERNEL/ExprEval/InterpKernelFunction.cxx b/src/INTERP_KERNEL/ExprEval/InterpKernelFunction.cxx
index 5c8e395ef..0c036bb30 100644
--- a/src/INTERP_KERNEL/ExprEval/InterpKernelFunction.cxx
+++ b/src/INTERP_KERNEL/ExprEval/InterpKernelFunction.cxx
@@ -22,6 +22,7 @@
 #include "InterpKernelValue.hxx"
 
 #include <cmath>
+#include <limits>
 
 using namespace INTERP_KERNEL;
 
@@ -191,7 +192,7 @@ IdentityFunction::~IdentityFunction()
 {
 }
 
-void IdentityFunction::operate(std::vector<Value *>& stack) const
+void IdentityFunction::operate(std::vector<Value *>& stck) const
 {
 }
 
@@ -199,6 +200,10 @@ void IdentityFunction::operateX86(std::vector<std::string>& asmb) const
 {
 }
 
+void IdentityFunction::operateStackOfDouble(std::vector<double>& stck) const
+{
+}
+
 const char *IdentityFunction::getRepr() const
 {
   return REPR;
@@ -218,7 +223,7 @@ int UnaryFunction::getNbInputParams() const
   return 1;
 }
 
-void PositiveFunction::operate(std::vector<Value *>& stack) const
+void PositiveFunction::operate(std::vector<Value *>& stck) const
 {
 }
 
@@ -226,6 +231,10 @@ void PositiveFunction::operateX86(std::vector<std::string>& asmb) const
 {
 }
 
+void PositiveFunction::operateStackOfDouble(std::vector<double>& stck) const
+{
+}
+
 const char *PositiveFunction::getRepr() const
 {
   return REPR;
@@ -240,9 +249,9 @@ NegateFunction::~NegateFunction()
 {
 }
 
-void NegateFunction::operate(std::vector<Value *>& stack) const
+void NegateFunction::operate(std::vector<Value *>& stck) const
 {
-  Value *val=stack.back();
+  Value *val=stck.back();
   val->negate();
 }
 
@@ -251,6 +260,12 @@ void NegateFunction::operateX86(std::vector<std::string>& asmb) const
   asmb.push_back("fchs");
 }
 
+void NegateFunction::operateStackOfDouble(std::vector<double>& stck) const
+{
+  double v(stck.back());
+  stck.back()=-v;
+}
+
 const char *NegateFunction::getRepr() const
 {
   return REPR;
@@ -265,9 +280,9 @@ CosFunction::~CosFunction()
 {
 }
 
-void CosFunction::operate(std::vector<Value *>& stack) const
+void CosFunction::operate(std::vector<Value *>& stck) const
 {
-  Value *val=stack.back();
+  Value *val=stck.back();
   val->cos();
 }
 
@@ -276,6 +291,12 @@ void CosFunction::operateX86(std::vector<std::string>& asmb) const
   asmb.push_back("fcos");
 }
 
+void CosFunction::operateStackOfDouble(std::vector<double>& stck) const
+{
+  double v(stck.back());
+  stck.back()=cos(v);
+}
+
 const char *CosFunction::getRepr() const
 {
   return REPR;
@@ -290,9 +311,9 @@ SinFunction::~SinFunction()
 {
 }
 
-void SinFunction::operate(std::vector<Value *>& stack) const
+void SinFunction::operate(std::vector<Value *>& stck) const
 {
-  Value *val=stack.back();
+  Value *val=stck.back();
   val->sin();
 }
 
@@ -301,6 +322,12 @@ void SinFunction::operateX86(std::vector<std::string>& asmb) const
   asmb.push_back("fsin");
 }
 
+void SinFunction::operateStackOfDouble(std::vector<double>& stck) const
+{
+  double v(stck.back());
+  stck.back()=sin(v);
+}
+
 const char *SinFunction::getRepr() const
 {
   return REPR;
@@ -315,9 +342,9 @@ TanFunction::~TanFunction()
 {
 }
 
-void TanFunction::operate(std::vector<Value *>& stack) const
+void TanFunction::operate(std::vector<Value *>& stck) const
 {
-  Value *val=stack.back();
+  Value *val=stck.back();
   val->tan();
 }
 
@@ -326,6 +353,12 @@ void TanFunction::operateX86(std::vector<std::string>& asmb) const
   throw INTERP_KERNEL::Exception("Assembly Not implemented yet !");
 }
 
+void TanFunction::operateStackOfDouble(std::vector<double>& stck) const
+{
+  double v(stck.back());
+  stck.back()=tan(v);
+}
+
 const char *TanFunction::getRepr() const
 {
   return REPR;
@@ -340,9 +373,9 @@ ACosFunction::~ACosFunction()
 {
 }
 
-void ACosFunction::operate(std::vector<Value *>& stack) const
+void ACosFunction::operate(std::vector<Value *>& stck) const
 {
-  Value *val=stack.back();
+  Value *val=stck.back();
   val->acos();
 }
 
@@ -351,6 +384,20 @@ void ACosFunction::operateX86(std::vector<std::string>& asmb) const
   throw INTERP_KERNEL::Exception("Assembly Not implemented yet !");
 }
 
+void ACosFunction::operateStackOfDouble(std::vector<double>& stck) const
+{
+  double v(stck.back());
+  stck.back()=acos(v);
+}
+
+void ACosFunction::operateStackOfDoubleSafe(std::vector<double>& stck) const
+{
+  double v(stck.back());
+  if(fabs(v)>1.)
+    throw INTERP_KERNEL::Exception("acos on a value which absolute is > 1 !");
+  stck.back()=acos(v);
+}
+
 const char *ACosFunction::getRepr() const
 {
   return REPR;
@@ -365,9 +412,9 @@ ASinFunction::~ASinFunction()
 {
 }
 
-void ASinFunction::operate(std::vector<Value *>& stack) const
+void ASinFunction::operate(std::vector<Value *>& stck) const
 {
-  Value *val=stack.back();
+  Value *val=stck.back();
   val->asin();
 }
 
@@ -376,6 +423,20 @@ void ASinFunction::operateX86(std::vector<std::string>& asmb) const
   throw INTERP_KERNEL::Exception("Assembly Not implemented yet !");
 }
 
+void ASinFunction::operateStackOfDouble(std::vector<double>& stck) const
+{
+  double v(stck.back());
+  stck.back()=asin(v);
+}
+
+void ASinFunction::operateStackOfDoubleSafe(std::vector<double>& stck) const
+{
+  double v(stck.back());
+  if(fabs(v)>1.)
+    throw INTERP_KERNEL::Exception("asin on a value which absolute is > 1 !");
+  stck.back()=asin(v);
+}
+
 const char *ASinFunction::getRepr() const
 {
   return REPR;
@@ -390,9 +451,9 @@ ATanFunction::~ATanFunction()
 {
 }
 
-void ATanFunction::operate(std::vector<Value *>& stack) const
+void ATanFunction::operate(std::vector<Value *>& stck) const
 {
-  Value *val=stack.back();
+  Value *val=stck.back();
   val->atan();
 }
 
@@ -401,6 +462,12 @@ void ATanFunction::operateX86(std::vector<std::string>& asmb) const
   throw INTERP_KERNEL::Exception("Assembly Not implemented yet !");
 }
 
+void ATanFunction::operateStackOfDouble(std::vector<double>& stck) const
+{
+  double v(stck.back());
+  stck.back()=atan(v);
+}
+
 const char *ATanFunction::getRepr() const
 {
   return REPR;
@@ -415,9 +482,9 @@ CoshFunction::~CoshFunction()
 {
 }
 
-void CoshFunction::operate(std::vector<Value *>& stack) const
+void CoshFunction::operate(std::vector<Value *>& stck) const
 {
-  Value *val=stack.back();
+  Value *val=stck.back();
   val->cosh();
 }
 
@@ -426,6 +493,12 @@ void CoshFunction::operateX86(std::vector<std::string>& asmb) const
   throw INTERP_KERNEL::Exception("Assembly Not implemented yet !");
 }
 
+void CoshFunction::operateStackOfDouble(std::vector<double>& stck) const
+{
+  double v(stck.back());
+  stck.back()=cosh(v);
+}
+
 const char *CoshFunction::getRepr() const
 {
   return REPR;
@@ -440,9 +513,9 @@ SinhFunction::~SinhFunction()
 {
 }
 
-void SinhFunction::operate(std::vector<Value *>& stack) const
+void SinhFunction::operate(std::vector<Value *>& stck) const
 {
-  Value *val=stack.back();
+  Value *val=stck.back();
   val->sinh();
 }
 
@@ -451,6 +524,12 @@ void SinhFunction::operateX86(std::vector<std::string>& asmb) const
   throw INTERP_KERNEL::Exception("Assembly Not implemented yet !");
 }
 
+void SinhFunction::operateStackOfDouble(std::vector<double>& stck) const
+{
+  double v(stck.back());
+  stck.back()=sinh(v);
+}
+
 const char *SinhFunction::getRepr() const
 {
   return REPR;
@@ -465,9 +544,9 @@ TanhFunction::~TanhFunction()
 {
 }
 
-void TanhFunction::operate(std::vector<Value *>& stack) const
+void TanhFunction::operate(std::vector<Value *>& stck) const
 {
-  Value *val=stack.back();
+  Value *val=stck.back();
   val->tanh();
 }
 
@@ -476,6 +555,12 @@ void TanhFunction::operateX86(std::vector<std::string>& asmb) const
   throw INTERP_KERNEL::Exception("Assembly Not implemented yet !");
 }
 
+void TanhFunction::operateStackOfDouble(std::vector<double>& stck) const
+{
+  double v(stck.back());
+  stck.back()=tanh(v);
+}
+
 const char *TanhFunction::getRepr() const
 {
   return REPR;
@@ -490,9 +575,9 @@ SqrtFunction::~SqrtFunction()
 {
 }
 
-void SqrtFunction::operate(std::vector<Value *>& stack) const
+void SqrtFunction::operate(std::vector<Value *>& stck) const
 {
-  Value *val=stack.back();
+  Value *val=stck.back();
   val->sqrt();
 }
 
@@ -501,6 +586,20 @@ void SqrtFunction::operateX86(std::vector<std::string>& asmb) const
   asmb.push_back("fsqrt");
 }
 
+void SqrtFunction::operateStackOfDouble(std::vector<double>& stck) const
+{
+  double v(stck.back());
+  stck.back()=sqrt(v);
+}
+
+void SqrtFunction::operateStackOfDoubleSafe(std::vector<double>& stck) const
+{
+  double v(stck.back());
+  if(v<0.)
+    throw INTERP_KERNEL::Exception("sqrt on a value < 0. !");
+  stck.back()=sqrt(v);
+}
+
 const char *SqrtFunction::getRepr() const
 {
   return REPR;
@@ -515,9 +614,9 @@ AbsFunction::~AbsFunction()
 {
 }
 
-void AbsFunction::operate(std::vector<Value *>& stack) const
+void AbsFunction::operate(std::vector<Value *>& stck) const
 {
-  Value *val=stack.back();
+  Value *val=stck.back();
   val->abs();
 }
 
@@ -526,6 +625,12 @@ void AbsFunction::operateX86(std::vector<std::string>& asmb) const
   asmb.push_back("fabs");
 }
 
+void AbsFunction::operateStackOfDouble(std::vector<double>& stck) const
+{
+  double v(stck.back());
+  stck.back()=fabs(v);
+}
+
 const char *AbsFunction::getRepr() const
 {
   return REPR;
@@ -536,9 +641,9 @@ bool AbsFunction::isACall() const
   return false;
 }
 
-void ExpFunction::operate(std::vector<Value *>& stack) const
+void ExpFunction::operate(std::vector<Value *>& stck) const
 {
-  Value *val=stack.back();
+  Value *val=stck.back();
   val->exp();
 }
 
@@ -547,6 +652,12 @@ void ExpFunction::operateX86(std::vector<std::string>& asmb) const
   throw INTERP_KERNEL::Exception("Assembly Not implemented yet !");
 }
 
+void ExpFunction::operateStackOfDouble(std::vector<double>& stck) const
+{
+  double v(stck.back());
+  stck.back()=std::exp(v);
+}
+
 const char *ExpFunction::getRepr() const
 {
   return REPR;
@@ -561,9 +672,9 @@ LnFunction::~LnFunction()
 {
 }
 
-void LnFunction::operate(std::vector<Value *>& stack) const
+void LnFunction::operate(std::vector<Value *>& stck) const
 {
-  Value *val=stack.back();
+  Value *val=stck.back();
   val->ln();
 }
 
@@ -572,6 +683,20 @@ void LnFunction::operateX86(std::vector<std::string>& asmb) const
   throw INTERP_KERNEL::Exception("Assembly Not implemented yet !");
 }
 
+void LnFunction::operateStackOfDouble(std::vector<double>& stck) const
+{
+  double v(stck.back());
+  stck.back()=std::log(v);
+}
+
+void LnFunction::operateStackOfDoubleSafe(std::vector<double>& stck) const
+{
+  double v(stck.back());
+  if(v<0.)
+    throw INTERP_KERNEL::Exception("ln on a value < 0. !");
+  stck.back()=std::log(v);
+}
+
 const char *LnFunction::getRepr() const
 {
   return REPR;
@@ -586,9 +711,9 @@ LogFunction::~LogFunction()
 {
 }
 
-void LogFunction::operate(std::vector<Value *>& stack) const
+void LogFunction::operate(std::vector<Value *>& stck) const
 {
-  Value *val=stack.back();
+  Value *val=stck.back();
   val->ln();
 }
 
@@ -597,6 +722,20 @@ void LogFunction::operateX86(std::vector<std::string>& asmb) const
   throw INTERP_KERNEL::Exception("Assembly for log Not implemented yet !");
 }
 
+void LogFunction::operateStackOfDouble(std::vector<double>& stck) const
+{
+  double v(stck.back());
+  stck.back()=std::log(v);
+}
+
+void LogFunction::operateStackOfDoubleSafe(std::vector<double>& stck) const
+{
+  double v(stck.back());
+  if(v<0.)
+    throw INTERP_KERNEL::Exception("log on a value < 0. !");
+  stck.back()=std::log(v);
+}
+
 const char *LogFunction::getRepr() const
 {
   return REPR;
@@ -611,9 +750,9 @@ Log10Function::~Log10Function()
 {
 }
 
-void Log10Function::operate(std::vector<Value *>& stack) const
+void Log10Function::operate(std::vector<Value *>& stck) const
 {
-  Value *val=stack.back();
+  Value *val=stck.back();
   val->log10();
 }
 
@@ -622,6 +761,20 @@ void Log10Function::operateX86(std::vector<std::string>& asmb) const
   throw INTERP_KERNEL::Exception("Assembly for log Not implemented yet !");
 }
 
+void Log10Function::operateStackOfDouble(std::vector<double>& stck) const
+{
+  double v(stck.back());
+  stck.back()=std::log10(v);
+}
+
+void Log10Function::operateStackOfDoubleSafe(std::vector<double>& stck) const
+{
+  double v(stck.back());
+  if(v<0.)
+    throw INTERP_KERNEL::Exception("log10 on a value < 0. !");
+  stck.back()=std::log10(v);
+}
+
 const char *Log10Function::getRepr() const
 {
   return REPR;
@@ -641,11 +794,11 @@ PlusFunction::~PlusFunction()
 {
 }
 
-void PlusFunction::operate(std::vector<Value *>& stack) const
+void PlusFunction::operate(std::vector<Value *>& stck) const
 {
-  Value *val1=stack.back();
-  stack.pop_back();
-  Value *& val2=stack.back();
+  Value *val1=stck.back();
+  stck.pop_back();
+  Value *& val2=stck.back();
   Value *val3;
   try
     {
@@ -666,6 +819,13 @@ void PlusFunction::operateX86(std::vector<std::string>& asmb) const
   asmb.push_back("faddp st1");
 }
 
+void PlusFunction::operateStackOfDouble(std::vector<double>& stck) const
+{
+  double a(stck.back());
+  stck.pop_back();
+  stck.back()=a+stck.back();
+}
+
 const char *PlusFunction::getRepr() const
 {
   return REPR;
@@ -680,11 +840,11 @@ MinusFunction::~MinusFunction()
 {
 }
 
-void MinusFunction::operate(std::vector<Value *>& stack) const
+void MinusFunction::operate(std::vector<Value *>& stck) const
 {
-  Value *val1=stack.back();
-  stack.pop_back();
-  Value *& val2=stack.back();
+  Value *val1=stck.back();
+  stck.pop_back();
+  Value *& val2=stck.back();
   Value *val3;
   try
     {
@@ -705,6 +865,13 @@ void MinusFunction::operateX86(std::vector<std::string>& asmb) const
   asmb.push_back("fsubp st1");
 }
 
+void MinusFunction::operateStackOfDouble(std::vector<double>& stck) const
+{
+  double a(stck.back());
+  stck.pop_back();
+  stck.back()=a-stck.back();
+}
+
 const char *MinusFunction::getRepr() const
 {
   return REPR;
@@ -719,11 +886,11 @@ MultFunction::~MultFunction()
 {
 }
 
-void MultFunction::operate(std::vector<Value *>& stack) const
+void MultFunction::operate(std::vector<Value *>& stck) const
 {
-  Value *val1=stack.back();
-  stack.pop_back();
-  Value *& val2=stack.back();
+  Value *val1=stck.back();
+  stck.pop_back();
+  Value *& val2=stck.back();
   Value *val3=val1->mult(val2);
   delete val1;
   delete val2;
@@ -735,6 +902,13 @@ void MultFunction::operateX86(std::vector<std::string>& asmb) const
   asmb.push_back("fmulp st1");
 }
 
+void MultFunction::operateStackOfDouble(std::vector<double>& stck) const
+{
+  double a(stck.back());
+  stck.pop_back();
+  stck.back()=a*stck.back();
+}
+
 const char *MultFunction::getRepr() const
 {
   return REPR;
@@ -749,11 +923,11 @@ DivFunction::~DivFunction()
 {
 }
 
-void DivFunction::operate(std::vector<Value *>& stack) const
+void DivFunction::operate(std::vector<Value *>& stck) const
 {
-  Value *val1=stack.back();
-  stack.pop_back();
-  Value *& val2=stack.back();
+  Value *val1=stck.back();
+  stck.pop_back();
+  Value *& val2=stck.back();
   Value *val3;
   try
     {
@@ -774,6 +948,22 @@ void DivFunction::operateX86(std::vector<std::string>& asmb) const
   asmb.push_back("fdivp st1");
 }
 
+void DivFunction::operateStackOfDouble(std::vector<double>& stck) const
+{
+  double a(stck.back());
+  stck.pop_back();
+  stck.back()=a/stck.back();
+}
+
+void DivFunction::operateStackOfDoubleSafe(std::vector<double>& stck) const
+{
+  double a(stck.back());
+  stck.pop_back();
+  if(stck.back()==0.)
+    throw INTERP_KERNEL::Exception("division by 0. !");
+  stck.back()=a/stck.back();
+}
+
 const char *DivFunction::getRepr() const
 {
   return REPR;
@@ -788,11 +978,11 @@ PowFunction::~PowFunction()
 {
 }
 
-void PowFunction::operate(std::vector<Value *>& stack) const
+void PowFunction::operate(std::vector<Value *>& stck) const
 {
-  Value *val1=stack.back();
-  stack.pop_back();
-  Value *& val2=stack.back();
+  Value *val1=stck.back();
+  stck.pop_back();
+  Value *& val2=stck.back();
   Value *val3;
   try
     {
@@ -813,6 +1003,23 @@ void PowFunction::operateX86(std::vector<std::string>& asmb) const
   throw INTERP_KERNEL::Exception("Assembly Not implemented yet !");
 }
 
+void PowFunction::operateStackOfDouble(std::vector<double>& stck) const
+{
+  double a(stck.back());
+  stck.pop_back();
+  stck.back()=std::pow(a,stck.back());
+}
+
+void PowFunction::operateStackOfDoubleSafe(std::vector<double>& stck) const
+{
+  double a(stck.back());
+  stck.pop_back();
+  double b(stck.back());
+  if(a<0.)
+    throw INTERP_KERNEL::Exception("pow with val < 0. !");
+  stck.back()=std::pow(a,b);
+}
+
 const char *PowFunction::getRepr() const
 {
   return REPR;
@@ -831,11 +1038,11 @@ MaxFunction::~MaxFunction()
 {
 }
 
-void MaxFunction::operate(std::vector<Value *>& stack) const
+void MaxFunction::operate(std::vector<Value *>& stck) const
 {
-  Value *val1=stack.back();
-  stack.pop_back();
-  Value *& val2=stack.back();
+  Value *val1=stck.back();
+  stck.pop_back();
+  Value *& val2=stck.back();
   Value *val3;
   try
     {
@@ -856,6 +1063,13 @@ void MaxFunction::operateX86(std::vector<std::string>& asmb) const
   throw INTERP_KERNEL::Exception("Assembly Not implemented yet !");
 }
 
+void MaxFunction::operateStackOfDouble(std::vector<double>& stck) const
+{
+  double a(stck.back());
+  stck.pop_back();
+  stck.back()=std::max(stck.back(),a);
+}
+
 const char *MaxFunction::getRepr() const
 {
   return REPR;
@@ -870,11 +1084,11 @@ MinFunction::~MinFunction()
 {
 }
 
-void MinFunction::operate(std::vector<Value *>& stack) const
+void MinFunction::operate(std::vector<Value *>& stck) const
 {
-  Value *val1=stack.back();
-  stack.pop_back();
-  Value *& val2=stack.back();
+  Value *val1=stck.back();
+  stck.pop_back();
+  Value *& val2=stck.back();
   Value *val3;
   try
     {
@@ -895,6 +1109,13 @@ void MinFunction::operateX86(std::vector<std::string>& asmb) const
   throw INTERP_KERNEL::Exception("Assembly Not implemented yet !");
 }
 
+void MinFunction::operateStackOfDouble(std::vector<double>& stck) const
+{
+  double a(stck.back());
+  stck.pop_back();
+  stck.back()=std::min(stck.back(),a);
+}
+
 const char *MinFunction::getRepr() const
 {
   return REPR;
@@ -909,11 +1130,11 @@ GreaterThanFunction::~GreaterThanFunction()
 {
 }
 
-void GreaterThanFunction::operate(std::vector<Value *>& stack) const
+void GreaterThanFunction::operate(std::vector<Value *>& stck) const
 {
-  Value *val1=stack.back();
-  stack.pop_back();
-  Value *& val2=stack.back();
+  Value *val1=stck.back();
+  stck.pop_back();
+  Value *& val2=stck.back();
   Value *val3;
   try
     {
@@ -934,6 +1155,14 @@ void GreaterThanFunction::operateX86(std::vector<std::string>& asmb) const
   throw INTERP_KERNEL::Exception("Assembly Not implemented yet !");
 }
 
+void GreaterThanFunction::operateStackOfDouble(std::vector<double>& stck) const
+{
+  double a(stck.back());
+  stck.pop_back();
+  double b(stck.back());
+  stck.back()=a>b?std::numeric_limits<double>::max():-std::numeric_limits<double>::max();
+}
+
 const char *GreaterThanFunction::getRepr() const
 {
   return REPR;
@@ -948,11 +1177,11 @@ LowerThanFunction::~LowerThanFunction()
 {
 }
 
-void LowerThanFunction::operate(std::vector<Value *>& stack) const
+void LowerThanFunction::operate(std::vector<Value *>& stck) const
 {
-  Value *val1=stack.back();
-  stack.pop_back();
-  Value *& val2=stack.back();
+  Value *val1=stck.back();
+  stck.pop_back();
+  Value *& val2=stck.back();
   Value *val3;
   try
     {
@@ -973,6 +1202,14 @@ void LowerThanFunction::operateX86(std::vector<std::string>& asmb) const
   throw INTERP_KERNEL::Exception("Assembly Not implemented yet !");
 }
 
+void LowerThanFunction::operateStackOfDouble(std::vector<double>& stck) const
+{
+  double a(stck.back());
+  stck.pop_back();
+  double b(stck.back());
+  stck.back()=a<b?std::numeric_limits<double>::max():-std::numeric_limits<double>::max();
+}
+
 const char *LowerThanFunction::getRepr() const
 {
   return REPR;
@@ -992,13 +1229,13 @@ IfFunction::~IfFunction()
 {
 }
 
-void IfFunction::operate(std::vector<Value *>& stack) const
+void IfFunction::operate(std::vector<Value *>& stck) const
 {
-  Value *val1=stack.back();
-  stack.pop_back();
-  Value *val2=stack.back();
-  stack.pop_back();
-  Value *&val3=stack.back();
+  Value *val1=stck.back();
+  stck.pop_back();
+  Value *val2=stck.back();
+  stck.pop_back();
+  Value *&val3=stck.back();
   Value *val4;
   try
     {
@@ -1021,6 +1258,28 @@ void IfFunction::operateX86(std::vector<std::string>& asmb) const
   throw INTERP_KERNEL::Exception("Assembly Not implemented yet !");
 }
 
+void IfFunction::operateStackOfDouble(std::vector<double>& stck) const
+{
+  double cond(stck.back());
+  stck.pop_back();
+  double the(stck.back());
+  stck.pop_back();
+  if(cond==std::numeric_limits<double>::max())
+    stck.back()=the;
+}
+
+void IfFunction::operateStackOfDoubleSafe(std::vector<double>& stck) const
+{
+  double cond(stck.back());
+  stck.pop_back();
+  double the(stck.back());
+  stck.pop_back();
+  if(cond!=std::numeric_limits<double>::max() && cond!=-std::numeric_limits<double>::max())
+    throw INTERP_KERNEL::Exception("ifFunc : first parameter of ternary func is NOT a consequence of a boolean op !");
+  if(cond==std::numeric_limits<double>::max())
+    stck.back()=the;
+}
+
 const char *IfFunction::getRepr() const
 {
   return REPR;
diff --git a/src/INTERP_KERNEL/ExprEval/InterpKernelFunction.hxx b/src/INTERP_KERNEL/ExprEval/InterpKernelFunction.hxx
index 95ccbbd99..682e0ebc4 100644
--- a/src/INTERP_KERNEL/ExprEval/InterpKernelFunction.hxx
+++ b/src/INTERP_KERNEL/ExprEval/InterpKernelFunction.hxx
@@ -47,10 +47,13 @@ namespace INTERP_KERNEL
   public:
     virtual ~Function();
     virtual int getNbInputParams() const = 0;
-    virtual void operate(std::vector<Value *>& stack) const = 0;
+    virtual void operate(std::vector<Value *>& stck) const = 0;
     virtual void operateX86(std::vector<std::string>& asmb) const = 0;
+    virtual void operateStackOfDouble(std::vector<double>& stck) const = 0;
+    virtual void operateStackOfDoubleSafe(std::vector<double>& stck) const { operateStackOfDouble(stck); }
     virtual const char *getRepr() const = 0;
     virtual bool isACall() const = 0;
+    virtual Function *deepCpy() const = 0;
   };
 
   class INTERPKERNEL_EXPORT UnaryFunction : public Function
@@ -63,10 +66,12 @@ namespace INTERP_KERNEL
   {
   public:
     ~IdentityFunction();
-    void operate(std::vector<Value *>& stack) const;
+    void operate(std::vector<Value *>& stck) const;
     void operateX86(std::vector<std::string>& asmb) const;
+    void operateStackOfDouble(std::vector<double>& stck) const;
     const char *getRepr() const;
     bool isACall() const;
+    IdentityFunction *deepCpy() const { return new IdentityFunction; }
   public:
     static const char REPR[];
   };
@@ -75,10 +80,12 @@ namespace INTERP_KERNEL
   {
   public:
     ~PositiveFunction();
-    void operate(std::vector<Value *>& stack) const;
+    void operate(std::vector<Value *>& stck) const;
     void operateX86(std::vector<std::string>& asmb) const;
+    void operateStackOfDouble(std::vector<double>& stck) const;
     const char *getRepr() const;
     bool isACall() const;
+    PositiveFunction *deepCpy() const { return new PositiveFunction; }
   public:
     static const char REPR[];
   };
@@ -87,10 +94,12 @@ namespace INTERP_KERNEL
   {
   public:
     ~NegateFunction();
-    void operate(std::vector<Value *>& stack) const;
+    void operate(std::vector<Value *>& stck) const;
     void operateX86(std::vector<std::string>& asmb) const;
+    void operateStackOfDouble(std::vector<double>& stck) const;
     const char *getRepr() const;
     bool isACall() const;
+    NegateFunction *deepCpy() const { return new NegateFunction; }
   public:
     static const char REPR[];
   };
@@ -99,10 +108,12 @@ namespace INTERP_KERNEL
   {
   public:
     ~CosFunction();
-    void operate(std::vector<Value *>& stack) const;
+    void operate(std::vector<Value *>& stck) const;
     void operateX86(std::vector<std::string>& asmb) const;
+    void operateStackOfDouble(std::vector<double>& stck) const;
     const char *getRepr() const;
     bool isACall() const;
+    CosFunction *deepCpy() const { return new CosFunction; }
   public:
     static const char REPR[];
   };
@@ -111,10 +122,12 @@ namespace INTERP_KERNEL
   {
   public:
     ~SinFunction();
-    void operate(std::vector<Value *>& stack) const;
+    void operate(std::vector<Value *>& stck) const;
     void operateX86(std::vector<std::string>& asmb) const;
+    void operateStackOfDouble(std::vector<double>& stck) const;
     const char *getRepr() const;
     bool isACall() const;
+    SinFunction *deepCpy() const { return new SinFunction; }
   public:
     static const char REPR[];
   };
@@ -123,10 +136,12 @@ namespace INTERP_KERNEL
   {
   public:
     ~TanFunction();
-    void operate(std::vector<Value *>& stack) const;
+    void operate(std::vector<Value *>& stck) const;
     void operateX86(std::vector<std::string>& asmb) const;
+    void operateStackOfDouble(std::vector<double>& stck) const;
     const char *getRepr() const;
     bool isACall() const;
+    TanFunction *deepCpy() const { return new TanFunction; }
   public:
     static const char REPR[];
   };
@@ -135,10 +150,13 @@ namespace INTERP_KERNEL
   {
   public:
     ~ACosFunction();
-    void operate(std::vector<Value *>& stack) const;
+    void operate(std::vector<Value *>& stck) const;
     void operateX86(std::vector<std::string>& asmb) const;
+    void operateStackOfDouble(std::vector<double>& stck) const;
+    void operateStackOfDoubleSafe(std::vector<double>& stck) const;
     const char *getRepr() const;
     bool isACall() const;
+    ACosFunction *deepCpy() const { return new ACosFunction; }
   public:
     static const char REPR[];
   };
@@ -147,10 +165,13 @@ namespace INTERP_KERNEL
   {
   public:
     ~ASinFunction();
-    void operate(std::vector<Value *>& stack) const;
+    void operate(std::vector<Value *>& stck) const;
     void operateX86(std::vector<std::string>& asmb) const;
+    void operateStackOfDouble(std::vector<double>& stck) const;
+    void operateStackOfDoubleSafe(std::vector<double>& stck) const;
     const char *getRepr() const;
     bool isACall() const;
+    ASinFunction *deepCpy() const { return new ASinFunction; }
   public:
     static const char REPR[];
   };
@@ -159,10 +180,12 @@ namespace INTERP_KERNEL
   {
   public:
     ~ATanFunction();
-    void operate(std::vector<Value *>& stack) const;
+    void operate(std::vector<Value *>& stck) const;
     void operateX86(std::vector<std::string>& asmb) const;
+    void operateStackOfDouble(std::vector<double>& stck) const;
     const char *getRepr() const;
     bool isACall() const;
+    ATanFunction *deepCpy() const { return new ATanFunction; }
   public:
     static const char REPR[];
   };
@@ -171,10 +194,12 @@ namespace INTERP_KERNEL
   {
   public:
     ~CoshFunction();
-    void operate(std::vector<Value *>& stack) const;
+    void operate(std::vector<Value *>& stck) const;
     void operateX86(std::vector<std::string>& asmb) const;
+    void operateStackOfDouble(std::vector<double>& stck) const;
     const char *getRepr() const;
     bool isACall() const;
+    CoshFunction *deepCpy() const { return new CoshFunction; }
   public:
     static const char REPR[];
   };
@@ -183,10 +208,12 @@ namespace INTERP_KERNEL
   {
   public:
     ~SinhFunction();
-    void operate(std::vector<Value *>& stack) const;
+    void operate(std::vector<Value *>& stck) const;
     void operateX86(std::vector<std::string>& asmb) const;
+    void operateStackOfDouble(std::vector<double>& stck) const;
     const char *getRepr() const;
     bool isACall() const;
+    SinhFunction *deepCpy() const { return new SinhFunction; }
   public:
     static const char REPR[];
   };
@@ -195,10 +222,12 @@ namespace INTERP_KERNEL
   {
   public:
     ~TanhFunction();
-    void operate(std::vector<Value *>& stack) const;
+    void operate(std::vector<Value *>& stck) const;
     void operateX86(std::vector<std::string>& asmb) const;
+    void operateStackOfDouble(std::vector<double>& stck) const;
     const char *getRepr() const;
     bool isACall() const;
+    TanhFunction *deepCpy() const { return new TanhFunction; }
   public:
     static const char REPR[];
   };
@@ -207,10 +236,13 @@ namespace INTERP_KERNEL
   {
   public:
     ~SqrtFunction();
-    void operate(std::vector<Value *>& stack) const;
     void operateX86(std::vector<std::string>& asmb) const;
+    void operate(std::vector<Value *>& stck) const;
+    void operateStackOfDouble(std::vector<double>& stck) const;
+    void operateStackOfDoubleSafe(std::vector<double>& stck) const;
     const char *getRepr() const;
     bool isACall() const;
+    SqrtFunction *deepCpy() const { return new SqrtFunction; }
   public:
     static const char REPR[];
   };
@@ -219,10 +251,12 @@ namespace INTERP_KERNEL
   {
   public:
     ~AbsFunction();
-    void operate(std::vector<Value *>& stack) const;
+    void operate(std::vector<Value *>& stck) const;
     void operateX86(std::vector<std::string>& asmb) const;
+    void operateStackOfDouble(std::vector<double>& stck) const;
     const char *getRepr() const;
     bool isACall() const;
+    AbsFunction *deepCpy() const { return new AbsFunction; }
   public:
     static const char REPR[];
   };
@@ -231,10 +265,12 @@ namespace INTERP_KERNEL
   {
   public:
     ~ExpFunction();
-    void operate(std::vector<Value *>& stack) const;
+    void operate(std::vector<Value *>& stck) const;
     void operateX86(std::vector<std::string>& asmb) const;
+    void operateStackOfDouble(std::vector<double>& stck) const;
     const char *getRepr() const;
     bool isACall() const;
+    ExpFunction *deepCpy() const { return new ExpFunction; }
   public:
     static const char REPR[];
   };
@@ -243,10 +279,13 @@ namespace INTERP_KERNEL
   {
   public:
     ~LnFunction();
-    void operate(std::vector<Value *>& stack) const;
+    void operate(std::vector<Value *>& stck) const;
     void operateX86(std::vector<std::string>& asmb) const;
+    void operateStackOfDouble(std::vector<double>& stck) const;
+    void operateStackOfDoubleSafe(std::vector<double>& stck) const;
     const char *getRepr() const;
     bool isACall() const;
+    LnFunction *deepCpy() const { return new LnFunction; }
   public:
     static const char REPR[];
   };
@@ -255,10 +294,13 @@ namespace INTERP_KERNEL
   {
   public:
     ~LogFunction();
-    void operate(std::vector<Value *>& stack) const;
+    void operate(std::vector<Value *>& stck) const;
     void operateX86(std::vector<std::string>& asmb) const;
+    void operateStackOfDouble(std::vector<double>& stck) const;
+    void operateStackOfDoubleSafe(std::vector<double>& stck) const;
     const char *getRepr() const;
     bool isACall() const;
+    LogFunction *deepCpy() const { return new LogFunction; }
   public:
     static const char REPR[];
   };
@@ -267,10 +309,13 @@ namespace INTERP_KERNEL
   {
   public:
     ~Log10Function();
-    void operate(std::vector<Value *>& stack) const;
+    void operate(std::vector<Value *>& stck) const;
     void operateX86(std::vector<std::string>& asmb) const;
+    void operateStackOfDouble(std::vector<double>& stck) const;
+    void operateStackOfDoubleSafe(std::vector<double>& stck) const;
     const char *getRepr() const;
     bool isACall() const;
+    Log10Function *deepCpy() const { return new Log10Function; }
   public:
     static const char REPR[];
   };
@@ -285,10 +330,12 @@ namespace INTERP_KERNEL
   {
   public:
     ~PlusFunction();
-    void operate(std::vector<Value *>& stack) const;
+    void operate(std::vector<Value *>& stck) const;
     void operateX86(std::vector<std::string>& asmb) const;
+    void operateStackOfDouble(std::vector<double>& stck) const;
     const char *getRepr() const;
     bool isACall() const;
+    PlusFunction *deepCpy() const { return new PlusFunction; }
   public:
     static const char REPR[];
   };
@@ -297,10 +344,12 @@ namespace INTERP_KERNEL
   {
   public:
     ~MinusFunction();
-    void operate(std::vector<Value *>& stack) const;
+    void operate(std::vector<Value *>& stck) const;
     void operateX86(std::vector<std::string>& asmb) const;
+    void operateStackOfDouble(std::vector<double>& stck) const;
     const char *getRepr() const;
     bool isACall() const;
+    MinusFunction *deepCpy() const { return new MinusFunction; }
   public:
     static const char REPR[];
   };
@@ -309,10 +358,12 @@ namespace INTERP_KERNEL
   {
   public:
     ~MultFunction();
-    void operate(std::vector<Value *>& stack) const;
+    void operate(std::vector<Value *>& stck) const;
     void operateX86(std::vector<std::string>& asmb) const;
+    void operateStackOfDouble(std::vector<double>& stck) const;
     const char *getRepr() const;
     bool isACall() const;
+    MultFunction *deepCpy() const { return new MultFunction; }
   public:
     static const char REPR[];
   };
@@ -321,10 +372,13 @@ namespace INTERP_KERNEL
   {
   public:
     ~DivFunction();
-    void operate(std::vector<Value *>& stack) const;
+    void operate(std::vector<Value *>& stck) const;
     void operateX86(std::vector<std::string>& asmb) const;
+    void operateStackOfDouble(std::vector<double>& stck) const;
+    void operateStackOfDoubleSafe(std::vector<double>& stck) const;
     const char *getRepr() const;
     bool isACall() const;
+    DivFunction *deepCpy() const { return new DivFunction; }
   public:
     static const char REPR[];
   };
@@ -333,10 +387,13 @@ namespace INTERP_KERNEL
   {
   public:
     ~PowFunction();
-    void operate(std::vector<Value *>& stack) const;
+    void operate(std::vector<Value *>& stck) const;
     void operateX86(std::vector<std::string>& asmb) const;
+    void operateStackOfDouble(std::vector<double>& stck) const;
+    void operateStackOfDoubleSafe(std::vector<double>& stck) const;
     const char *getRepr() const;
     bool isACall() const;
+    PowFunction *deepCpy() const { return new PowFunction; }
   public:
     static const char REPR[];
   };
@@ -345,10 +402,12 @@ namespace INTERP_KERNEL
   {
   public:
     ~MaxFunction();
-    void operate(std::vector<Value *>& stack) const;
+    void operate(std::vector<Value *>& stck) const;
     void operateX86(std::vector<std::string>& asmb) const;
+    void operateStackOfDouble(std::vector<double>& stck) const;
     const char *getRepr() const;
     bool isACall() const;
+    MaxFunction *deepCpy() const { return new MaxFunction; }
   public:
     static const char REPR[];
   };
@@ -357,10 +416,12 @@ namespace INTERP_KERNEL
   {
   public:
     ~MinFunction();
-    void operate(std::vector<Value *>& stack) const;
+    void operate(std::vector<Value *>& stck) const;
     void operateX86(std::vector<std::string>& asmb) const;
+    void operateStackOfDouble(std::vector<double>& stck) const;
     const char *getRepr() const;
     bool isACall() const;
+    MinFunction *deepCpy() const { return new MinFunction; }
   public:
     static const char REPR[];
   };
@@ -369,10 +430,12 @@ namespace INTERP_KERNEL
   {
   public:
     ~GreaterThanFunction();
-    void operate(std::vector<Value *>& stack) const;
+    void operate(std::vector<Value *>& stck) const;
     void operateX86(std::vector<std::string>& asmb) const;
+    void operateStackOfDouble(std::vector<double>& stck) const;
     const char *getRepr() const;
     bool isACall() const;
+    GreaterThanFunction *deepCpy() const { return new GreaterThanFunction; }
   public:
     static const char REPR[];
   };
@@ -381,10 +444,12 @@ namespace INTERP_KERNEL
   {
   public:
     ~LowerThanFunction();
-    void operate(std::vector<Value *>& stack) const;
+    void operate(std::vector<Value *>& stck) const;
     void operateX86(std::vector<std::string>& asmb) const;
+    void operateStackOfDouble(std::vector<double>& stck) const;
     const char *getRepr() const;
     bool isACall() const;
+    LowerThanFunction *deepCpy() const { return new LowerThanFunction; }
   public:
     static const char REPR[];
   };
@@ -399,10 +464,13 @@ namespace INTERP_KERNEL
   {
   public:
     ~IfFunction();
-    void operate(std::vector<Value *>& stack) const;
+    void operate(std::vector<Value *>& stck) const;
     void operateX86(std::vector<std::string>& asmb) const;
+    void operateStackOfDouble(std::vector<double>& stck) const;
+    void operateStackOfDoubleSafe(std::vector<double>& stck) const;
     const char *getRepr() const;
     bool isACall() const;
+    IfFunction *deepCpy() const { return new IfFunction; }
   public:
     static const char REPR[];
   };
diff --git a/src/INTERP_KERNELTest/ExprEvalInterpTest.cxx b/src/INTERP_KERNELTest/ExprEvalInterpTest.cxx
index ab447534d..bfbe36b6f 100644
--- a/src/INTERP_KERNELTest/ExprEvalInterpTest.cxx
+++ b/src/INTERP_KERNELTest/ExprEvalInterpTest.cxx
@@ -538,3 +538,159 @@ void ExprEvalInterpTest::testInterpreter5()
   expr4.parse();
   CPPUNIT_ASSERT_DOUBLES_EQUAL(6994207.8359543988,expr4.evaluate(),1e-5);
 }
+
+/*!
+ * Test focusing on fast evaluator.
+ */
+void ExprEvalInterpTest::testInterpreter6()
+{
+  std::vector<double> stackOfVal;
+  //
+  stackOfVal.clear();
+  INTERP_KERNEL::ExprParser expr1("1.-2./3.");
+  expr1.parse();
+  expr1.prepareFastEvaluator();
+  expr1.evaluateDoubleInternal(stackOfVal);
+  CPPUNIT_ASSERT_EQUAL(1,(int)stackOfVal.size());
+  CPPUNIT_ASSERT_DOUBLES_EQUAL(0.33333333333333333,stackOfVal.back(),1e-14);
+  //
+  stackOfVal.clear();
+  INTERP_KERNEL::ExprParser expr2("1.-2.^3.");
+  expr2.parse();
+  expr2.prepareFastEvaluator();
+  expr2.evaluateDoubleInternal(stackOfVal);
+  CPPUNIT_ASSERT_EQUAL(1,(int)stackOfVal.size());
+  CPPUNIT_ASSERT_DOUBLES_EQUAL(-7.,stackOfVal.back(),1e-14);
+  //
+  stackOfVal.clear();
+  INTERP_KERNEL::ExprParser expr3("(7.-2.)^3.");
+  expr3.parse();
+  expr3.prepareFastEvaluator();
+  expr3.evaluateDoubleInternal(stackOfVal);
+  CPPUNIT_ASSERT_EQUAL(1,(int)stackOfVal.size());
+  CPPUNIT_ASSERT_DOUBLES_EQUAL(125.,stackOfVal.back(),1e-12);
+  // now playing with one parameter - calling several times
+  stackOfVal.clear();
+  INTERP_KERNEL::ExprParser expr4("1.2/(7.-2.*cos(x/3))");
+  expr4.parse();
+  expr4.prepareFastEvaluator();
+  double z;
+  expr4.prepareExprEvaluationDouble(std::vector<std::string>(1,"x"),1,1,0,&z,&z+1);
+  expr4.prepareFastEvaluator();
+  z=0.77;
+  expr4.evaluateDoubleInternal(stackOfVal);
+  CPPUNIT_ASSERT_EQUAL(1,(int)stackOfVal.size());
+  CPPUNIT_ASSERT_DOUBLES_EQUAL(0.23689586281558844,stackOfVal.back(),1e-12);
+  stackOfVal.pop_back();
+  z=0.55;
+  expr4.evaluateDoubleInternal(stackOfVal);
+  CPPUNIT_ASSERT_EQUAL(1,(int)stackOfVal.size());
+  CPPUNIT_ASSERT_DOUBLES_EQUAL(0.2384018932069258,stackOfVal.back(),1e-12);
+  //
+  stackOfVal.clear();
+  INTERP_KERNEL::ExprParser expr5("x-2*cos(y/3.)");
+  expr5.parse();
+  expr5.prepareFastEvaluator();
+  double *aa(new double[2]);
+  std::vector<std::string> vv(2); vv[0]="x"; vv[1]="y";
+  expr5.prepareExprEvaluationDouble(vv,2,1,0,aa,aa+2);
+  expr5.prepareFastEvaluator();
+  aa[0]=0.3; aa[1]=0.5;
+  expr5.evaluateDoubleInternal(stackOfVal);
+  CPPUNIT_ASSERT_EQUAL(1,(int)stackOfVal.size());
+  CPPUNIT_ASSERT_DOUBLES_EQUAL(-1.67228646312585,stackOfVal.back(),1e-14);
+  stackOfVal.pop_back();
+  aa[0]=0.5; aa[1]=0.3;
+  expr5.evaluateDoubleInternal(stackOfVal);
+  CPPUNIT_ASSERT_EQUAL(1,(int)stackOfVal.size());
+  CPPUNIT_ASSERT_DOUBLES_EQUAL(-1.4900083305560516,stackOfVal.back(),1e-14);
+  stackOfVal.pop_back();
+  delete [] aa;
+  //
+  stackOfVal.clear();
+  INTERP_KERNEL::ExprParser expr6("x*IVec-2*cos(y/3.)*JVec");
+  expr6.parse();
+  aa=new double[2];
+  vv.resize(2); vv[0]="x"; vv[1]="y";
+  expr6.prepareExprEvaluationDouble(vv,2,2,0,aa,aa+2);//emulate 1st interpreter
+  expr6.prepareFastEvaluator();
+  aa[0]=0.3; aa[1]=0.5;
+  expr6.evaluateDoubleInternal(stackOfVal);
+  CPPUNIT_ASSERT_EQUAL(1,(int)stackOfVal.size());
+  CPPUNIT_ASSERT_DOUBLES_EQUAL(0.3,stackOfVal.back(),1e-14);
+  stackOfVal.pop_back();
+  expr6.prepareExprEvaluationDouble(vv,2,2,1,aa,aa+2);//emulate 2nd interpreter
+  expr6.prepareFastEvaluator();
+  expr6.evaluateDoubleInternal(stackOfVal);
+  CPPUNIT_ASSERT_EQUAL(1,(int)stackOfVal.size());
+  CPPUNIT_ASSERT_DOUBLES_EQUAL(-1.97228646312585,stackOfVal.back(),1e-14);
+  stackOfVal.pop_back();
+  delete [] aa;
+  //
+  stackOfVal.clear();
+  INTERP_KERNEL::ExprParser expr7("if(x>3.,-6,7.)");
+  expr7.parse();
+  expr7.prepareExprEvaluationDouble(std::vector<std::string>(1,"x"),1,1,0,&z,&z+1);
+  expr7.prepareFastEvaluator();
+  z=3.1;
+  expr7.evaluateDoubleInternal(stackOfVal);
+  CPPUNIT_ASSERT_EQUAL(1,(int)stackOfVal.size());
+  CPPUNIT_ASSERT_DOUBLES_EQUAL(-6.,stackOfVal.back(),1e-14);
+  stackOfVal.pop_back();
+  z=2.8;
+  expr7.evaluateDoubleInternal(stackOfVal);
+  CPPUNIT_ASSERT_EQUAL(1,(int)stackOfVal.size());
+  CPPUNIT_ASSERT_DOUBLES_EQUAL(7.,stackOfVal.back(),1e-14);
+  stackOfVal.pop_back();
+  //
+  stackOfVal.clear();
+  INTERP_KERNEL::ExprParser expr8("if(x<3.,-6,7.)");
+  expr8.parse();
+  expr8.prepareExprEvaluationDouble(std::vector<std::string>(1,"x"),1,1,0,&z,&z+1);
+  expr8.prepareFastEvaluator();
+  z=3.1;
+  expr8.evaluateDoubleInternal(stackOfVal);
+  CPPUNIT_ASSERT_EQUAL(1,(int)stackOfVal.size());
+  CPPUNIT_ASSERT_DOUBLES_EQUAL(7.,stackOfVal.back(),1e-14);
+  stackOfVal.pop_back();
+  z=2.8;
+  expr8.evaluateDoubleInternal(stackOfVal);
+  CPPUNIT_ASSERT_EQUAL(1,(int)stackOfVal.size());
+  CPPUNIT_ASSERT_DOUBLES_EQUAL(-6.,stackOfVal.back(),1e-14);
+  stackOfVal.pop_back();
+  //
+  stackOfVal.clear();
+  INTERP_KERNEL::ExprParser expr9("x*x/2");//this test is very important to test for iso priority with more than one
+  expr9.parse();
+  expr9.prepareExprEvaluationDouble(std::vector<std::string>(1,"x"),1,1,0,&z,&z+1);
+  expr9.prepareFastEvaluator();
+  z=3.;
+  expr9.evaluateDoubleInternal(stackOfVal);
+  CPPUNIT_ASSERT_EQUAL(1,(int)stackOfVal.size());
+  CPPUNIT_ASSERT_DOUBLES_EQUAL(4.5,stackOfVal.back(),1e-14);
+  stackOfVal.pop_back();
+  //
+  stackOfVal.clear();
+  INTERP_KERNEL::ExprParser expr10("2./3./4./5.");
+  expr10.parse();
+  expr10.prepareFastEvaluator();
+  expr10.evaluateDoubleInternal(stackOfVal);
+  CPPUNIT_ASSERT_EQUAL(1,(int)stackOfVal.size());
+  CPPUNIT_ASSERT_DOUBLES_EQUAL(0.03333333333333333,stackOfVal.back(),1e-13);
+  //
+  stackOfVal.clear();
+  INTERP_KERNEL::ExprParser expr11("2./3./4.*5.");
+  expr11.parse();
+  expr11.prepareFastEvaluator();
+  expr11.evaluateDoubleInternal(stackOfVal);
+  CPPUNIT_ASSERT_EQUAL(1,(int)stackOfVal.size());
+  CPPUNIT_ASSERT_DOUBLES_EQUAL(0.83333333333333333,stackOfVal.back(),1e-14);
+  //
+  stackOfVal.clear();
+  INTERP_KERNEL::ExprParser expr12("2./3.*4.*5.");
+  expr12.parse();
+  expr12.prepareFastEvaluator();
+  expr12.evaluateDoubleInternal(stackOfVal);
+  CPPUNIT_ASSERT_EQUAL(1,(int)stackOfVal.size());
+  CPPUNIT_ASSERT_DOUBLES_EQUAL(13.333333333333333,stackOfVal.back(),1e-14);
+}
diff --git a/src/INTERP_KERNELTest/ExprEvalInterpTest.hxx b/src/INTERP_KERNELTest/ExprEvalInterpTest.hxx
index 0d0496307..f19cb208c 100644
--- a/src/INTERP_KERNELTest/ExprEvalInterpTest.hxx
+++ b/src/INTERP_KERNELTest/ExprEvalInterpTest.hxx
@@ -40,6 +40,7 @@ namespace INTERP_TEST
     CPPUNIT_TEST( testInterpreter3 );
     CPPUNIT_TEST( testInterpreter4 );
     CPPUNIT_TEST( testInterpreter5 );
+    CPPUNIT_TEST( testInterpreter6 );
     CPPUNIT_TEST_SUITE_END();
   public:
     void setUp() { }
@@ -53,6 +54,7 @@ namespace INTERP_TEST
     void testInterpreter3();
     void testInterpreter4();
     void testInterpreter5();
+    void testInterpreter6();
     void testInterpreterUnit0();
     void testInterpreterUnit1();
   };
diff --git a/src/MEDCoupling/MEDCouplingMemArray.cxx b/src/MEDCoupling/MEDCouplingMemArray.cxx
index 6556e4964..d6d65c9e8 100644
--- a/src/MEDCoupling/MEDCouplingMemArray.cxx
+++ b/src/MEDCoupling/MEDCouplingMemArray.cxx
@@ -4307,6 +4307,8 @@ DataArrayDouble *DataArrayDouble::applyFunc(int nbOfComp, FunctionToEvaluate fun
  *  \param [in] nbOfComp - number of components in the result array.
  *  \param [in] func - the expression defining how to transform a tuple of \a this array.
  *              Supported expressions are described \ref MEDCouplingArrayApplyFuncExpr "here".
+ *  \param [in] isSafe - By default true. If true invalid operation (division by 0. acos of value > 1. ...) leads to a throw of an exception.
+ *              If false the computation is carried on without any notification. When false the evaluation is a little faster.
  *  \return DataArrayDouble * - the new instance of DataArrayDouble containing the
  *          same number of tuples as \a this array and \a nbOfComp components.
  *          The caller is to delete this result array using decrRef() as it is no more
@@ -4314,144 +4316,210 @@ DataArrayDouble *DataArrayDouble::applyFunc(int nbOfComp, FunctionToEvaluate fun
  *  \throw If \a this is not allocated.
  *  \throw If computing \a func fails.
  */
-DataArrayDouble *DataArrayDouble::applyFunc(int nbOfComp, const std::string& func) const
+DataArrayDouble *DataArrayDouble::applyFunc(int nbOfComp, const std::string& func, bool isSafe) const
 {
-  checkAllocated();
   INTERP_KERNEL::ExprParser expr(func);
   expr.parse();
   std::set<std::string> vars;
   expr.getTrueSetOfVars(vars);
-  int oldNbOfComp=getNumberOfComponents();
-  if((int)vars.size()>oldNbOfComp)
-    {
-      std::ostringstream oss; oss << "The field has " << oldNbOfComp << " components and there are ";
-      oss << vars.size() << " variables : ";
-      std::copy(vars.begin(),vars.end(),std::ostream_iterator<std::string>(oss," "));
-      throw INTERP_KERNEL::Exception(oss.str().c_str());
-    }
   std::vector<std::string> varsV(vars.begin(),vars.end());
-  expr.prepareExprEvaluation(varsV,oldNbOfComp,nbOfComp);
-  //
-  DataArrayDouble *newArr=DataArrayDouble::New();
-  int nbOfTuples=getNumberOfTuples();
-  newArr->alloc(nbOfTuples,nbOfComp);
-  const double *ptr=getConstPointer();
-  double *ptrToFill=newArr->getPointer();
-  for(int i=0;i<nbOfTuples;i++)
-    {
-      try
-      {
-          expr.evaluateExpr(nbOfComp,ptr+i*oldNbOfComp,ptrToFill+i*nbOfComp);
-      }
-      catch(INTERP_KERNEL::Exception& e)
-      {
-          std::ostringstream oss; oss << "For tuple # " << i << " with value (";
-          std::copy(ptr+oldNbOfComp*i,ptr+oldNbOfComp*(i+1),std::ostream_iterator<double>(oss,", "));
-          oss << ") : Evaluation of function failed !" << e.what();
-          newArr->decrRef();
-          throw INTERP_KERNEL::Exception(oss.str().c_str());
-      }
-    }
-  return newArr;
+  return applyFunc3(nbOfComp,varsV,func,isSafe);
 }
 
 /*!
  * Returns a new DataArrayDouble created from \a this one by applying a function to every
- * tuple of \a this array. Textual data is not copied.
+ * tuple of \a this array. Textual data is not copied. This method works by tuples (whatever its size).
+ * If \a this is a one component array, call applyFuncOnThis instead that performs the same work faster.
+ *
  * For more info see \ref MEDCouplingArrayApplyFunc0.
  *  \param [in] func - the expression defining how to transform a tuple of \a this array.
  *              Supported expressions are described \ref MEDCouplingArrayApplyFuncExpr "here".
+ *  \param [in] isSafe - By default true. If true invalid operation (division by 0. acos of value > 1. ...) leads to a throw of an exception.
+ *                       If false the computation is carried on without any notification. When false the evaluation is a little faster.
  *  \return DataArrayDouble * - the new instance of DataArrayDouble containing the
  *          same number of tuples and components as \a this array.
  *          The caller is to delete this result array using decrRef() as it is no more
  *          needed.
+ *  \sa applyFuncOnThis
  *  \throw If \a this is not allocated.
  *  \throw If computing \a func fails.
  */
-DataArrayDouble *DataArrayDouble::applyFunc(const std::string& func) const
+DataArrayDouble *DataArrayDouble::applyFunc(const std::string& func, bool isSafe) const
 {
+  int nbOfComp(getNumberOfComponents());
+  if(nbOfComp<=0)
+    throw INTERP_KERNEL::Exception("DataArrayDouble::applyFunc : output number of component must be > 0 !");
   checkAllocated();
+  int nbOfTuples(getNumberOfTuples());
+  MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> newArr(DataArrayDouble::New());
+  newArr->alloc(nbOfTuples,nbOfComp);
   INTERP_KERNEL::ExprParser expr(func);
   expr.parse();
-  expr.prepareExprEvaluationVec();
-  //
-  DataArrayDouble *newArr=DataArrayDouble::New();
-  int nbOfTuples=getNumberOfTuples();
-  int nbOfComp=getNumberOfComponents();
-  newArr->alloc(nbOfTuples,nbOfComp);
-  const double *ptr=getConstPointer();
-  double *ptrToFill=newArr->getPointer();
-  for(int i=0;i<nbOfTuples;i++)
+  std::set<std::string> vars;
+  expr.getTrueSetOfVars(vars);
+  if((int)vars.size()>1)
     {
-      try
-      {
-          expr.evaluateExpr(nbOfComp,ptr+i*nbOfComp,ptrToFill+i*nbOfComp);
-      }
-      catch(INTERP_KERNEL::Exception& e)
-      {
-          std::ostringstream oss; oss << "For tuple # " << i << " with value (";
-          std::copy(ptr+nbOfComp*i,ptr+nbOfComp*(i+1),std::ostream_iterator<double>(oss,", "));
-          oss << ") : Evaluation of function failed ! " << e.what();
-          newArr->decrRef();
-          throw INTERP_KERNEL::Exception(oss.str().c_str());
-      }
+      std::ostringstream oss; oss << "DataArrayDouble::applyFunc : this method works only with at most one var func expression ! If you need to map comps on variables please use applyFunc2 or applyFunc3 instead ! Vars in expr are : ";
+      std::copy(vars.begin(),vars.end(),std::ostream_iterator<std::string>(oss," "));
+      throw INTERP_KERNEL::Exception(oss.str().c_str());
     }
-  return newArr;
+  if(vars.empty())
+    {
+      expr.prepareFastEvaluator();
+      newArr->rearrange(1);
+      newArr->fillWithValue(expr.evaluateDouble());
+      newArr->rearrange(nbOfComp);
+      return newArr.retn();
+    }
+  std::vector<std::string> vars2(vars.begin(),vars.end());
+  double buff,*ptrToFill(newArr->getPointer());
+  const double *ptr(begin());
+  std::vector<double> stck;
+  expr.prepareExprEvaluationDouble(vars2,1,1,0,&buff,&buff+1);
+  expr.prepareFastEvaluator();
+  if(!isSafe)
+    {
+      for(int i=0;i<nbOfTuples;i++)
+        {
+          for(int iComp=0;iComp<nbOfComp;iComp++,ptr++,ptrToFill++)
+            {
+              buff=*ptr;
+              expr.evaluateDoubleInternal(stck);
+              *ptrToFill=stck.back();
+              stck.pop_back();
+            }
+        }
+    }
+  else
+    {
+      for(int i=0;i<nbOfTuples;i++)
+        {
+          for(int iComp=0;iComp<nbOfComp;iComp++,ptr++,ptrToFill++)
+            {
+              buff=*ptr;
+              try
+              {
+                  expr.evaluateDoubleInternalSafe(stck);
+              }
+              catch(INTERP_KERNEL::Exception& e)
+              {
+                  std::ostringstream oss; oss << "For tuple # " << i << " component # " << iComp << " with value (";
+                  oss << buff;
+                  oss << ") : Evaluation of function failed !" << e.what();
+                  throw INTERP_KERNEL::Exception(oss.str().c_str());
+              }
+              *ptrToFill=stck.back();
+              stck.pop_back();
+            }
+        }
+    }
+  return newArr.retn();
 }
 
 /*!
- * Returns a new DataArrayDouble created from \a this one by applying a function to every
- * tuple of \a this array. Textual data is not copied.
- * For more info see \ref MEDCouplingArrayApplyFunc2.
- *  \param [in] nbOfComp - number of components in the result array.
+ * This method is a non const method that modify the array in \a this.
+ * This method only works on one component array. It means that function \a func must
+ * contain at most one variable.
+ * This method is a specialization of applyFunc method with one parameter on one component array.
+ *
  *  \param [in] func - the expression defining how to transform a tuple of \a this array.
  *              Supported expressions are described \ref MEDCouplingArrayApplyFuncExpr "here".
- *  \return DataArrayDouble * - the new instance of DataArrayDouble containing the
- *          same number of tuples as \a this array.
- *          The caller is to delete this result array using decrRef() as it is no more
- *          needed.
- *  \throw If \a this is not allocated.
- *  \throw If \a func contains vars that are not in \a this->getInfoOnComponent().
- *  \throw If computing \a func fails.
+ *  \param [in] isSafe - By default true. If true invalid operation (division by 0. acos of value > 1. ...) leads to a throw of an exception.
+ *              If false the computation is carried on without any notification. When false the evaluation is a little faster.
+ *
+ * \sa applyFunc
  */
-DataArrayDouble *DataArrayDouble::applyFunc2(int nbOfComp, const std::string& func) const
+void DataArrayDouble::applyFuncOnThis(const std::string& func, bool isSafe)
 {
+  int nbOfComp(getNumberOfComponents());
+  if(nbOfComp<=0)
+    throw INTERP_KERNEL::Exception("DataArrayDouble::applyFuncOnThis : output number of component must be > 0 !");
   checkAllocated();
+  int nbOfTuples(getNumberOfTuples());
   INTERP_KERNEL::ExprParser expr(func);
   expr.parse();
   std::set<std::string> vars;
   expr.getTrueSetOfVars(vars);
-  int oldNbOfComp=getNumberOfComponents();
-  if((int)vars.size()>oldNbOfComp)
+  if((int)vars.size()>1)
     {
-      std::ostringstream oss; oss << "The field has " << oldNbOfComp << " components and there are ";
-      oss << vars.size() << " variables : ";
+      std::ostringstream oss; oss << "DataArrayDouble::applyFuncOnThis : this method works only with at most one var func expression ! If you need to map comps on variables please use applyFunc2 or applyFunc3 instead ! Vars in expr are : ";
       std::copy(vars.begin(),vars.end(),std::ostream_iterator<std::string>(oss," "));
       throw INTERP_KERNEL::Exception(oss.str().c_str());
     }
-  expr.prepareExprEvaluation(getVarsOnComponent(),oldNbOfComp,nbOfComp);
-  //
-  DataArrayDouble *newArr=DataArrayDouble::New();
-  int nbOfTuples=getNumberOfTuples();
-  newArr->alloc(nbOfTuples,nbOfComp);
-  const double *ptr=getConstPointer();
-  double *ptrToFill=newArr->getPointer();
-  for(int i=0;i<nbOfTuples;i++)
+  if(vars.empty())
     {
-      try
-      {
-          expr.evaluateExpr(nbOfComp,ptr+i*oldNbOfComp,ptrToFill+i*nbOfComp);
-      }
-      catch(INTERP_KERNEL::Exception& e)
-      {
-          std::ostringstream oss; oss << "For tuple # " << i << " with value (";
-          std::copy(ptr+oldNbOfComp*i,ptr+oldNbOfComp*(i+1),std::ostream_iterator<double>(oss,", "));
-          oss << ") : Evaluation of function failed !" << e.what();
-          newArr->decrRef();
-          throw INTERP_KERNEL::Exception(oss.str().c_str());
-      }
+      expr.prepareFastEvaluator();
+      std::vector<std::string> compInfo(getInfoOnComponents());
+      rearrange(1);
+      fillWithValue(expr.evaluateDouble());
+      rearrange(nbOfComp);
+      setInfoOnComponents(compInfo);
+      return ;
+    }
+  std::vector<std::string> vars2(vars.begin(),vars.end());
+  double buff,*ptrToFill(getPointer());
+  const double *ptr(begin());
+  std::vector<double> stck;
+  expr.prepareExprEvaluationDouble(vars2,1,1,0,&buff,&buff+1);
+  expr.prepareFastEvaluator();
+  if(!isSafe)
+    {
+      for(int i=0;i<nbOfTuples;i++)
+        {
+          for(int iComp=0;iComp<nbOfComp;iComp++,ptr++,ptrToFill++)
+            {
+              buff=*ptr;
+              expr.evaluateDoubleInternal(stck);
+              *ptrToFill=stck.back();
+              stck.pop_back();
+            }
+        }
+    }
+  else
+    {
+      for(int i=0;i<nbOfTuples;i++)
+        {
+          for(int iComp=0;iComp<nbOfComp;iComp++,ptr++,ptrToFill++)
+            {
+              buff=*ptr;
+              try
+              {
+                  expr.evaluateDoubleInternalSafe(stck);
+              }
+              catch(INTERP_KERNEL::Exception& e)
+              {
+                  std::ostringstream oss; oss << "For tuple # " << i << " component # " << iComp << " with value (";
+                  oss << buff;
+                  oss << ") : Evaluation of function failed !" << e.what();
+                  throw INTERP_KERNEL::Exception(oss.str().c_str());
+              }
+              *ptrToFill=stck.back();
+              stck.pop_back();
+            }
+        }
     }
-  return newArr;
+}
+
+/*!
+ * Returns a new DataArrayDouble created from \a this one by applying a function to every
+ * tuple of \a this array. Textual data is not copied.
+ * For more info see \ref MEDCouplingArrayApplyFunc2.
+ *  \param [in] nbOfComp - number of components in the result array.
+ *  \param [in] func - the expression defining how to transform a tuple of \a this array.
+ *              Supported expressions are described \ref MEDCouplingArrayApplyFuncExpr "here".
+ *  \param [in] isSafe - By default true. If true invalid operation (division by 0. acos of value > 1. ...) leads to a throw of an exception.
+ *              If false the computation is carried on without any notification. When false the evaluation is a little faster.
+ *  \return DataArrayDouble * - the new instance of DataArrayDouble containing the
+ *          same number of tuples as \a this array.
+ *          The caller is to delete this result array using decrRef() as it is no more
+ *          needed.
+ *  \throw If \a this is not allocated.
+ *  \throw If \a func contains vars that are not in \a this->getInfoOnComponent().
+ *  \throw If computing \a func fails.
+ */
+DataArrayDouble *DataArrayDouble::applyFunc2(int nbOfComp, const std::string& func, bool isSafe) const
+{
+  return applyFunc3(nbOfComp,getVarsOnComponent(),func,isSafe);
 }
 
 /*!
@@ -4462,6 +4530,8 @@ DataArrayDouble *DataArrayDouble::applyFunc2(int nbOfComp, const std::string& fu
  *  \param [in] varsOrder - sequence of vars defining their order.
  *  \param [in] func - the expression defining how to transform a tuple of \a this array.
  *              Supported expressions are described \ref MEDCouplingArrayApplyFuncExpr "here".
+ *  \param [in] isSafe - By default true. If true invalid operation (division by 0. acos of value > 1. ...) leads to a throw of an exception.
+ *              If false the computation is carried on without any notification. When false the evaluation is a little faster.
  *  \return DataArrayDouble * - the new instance of DataArrayDouble containing the
  *          same number of tuples as \a this array.
  *          The caller is to delete this result array using decrRef() as it is no more
@@ -4470,14 +4540,20 @@ DataArrayDouble *DataArrayDouble::applyFunc2(int nbOfComp, const std::string& fu
  *  \throw If \a func contains vars not in \a varsOrder.
  *  \throw If computing \a func fails.
  */
-DataArrayDouble *DataArrayDouble::applyFunc3(int nbOfComp, const std::vector<std::string>& varsOrder, const std::string& func) const
+DataArrayDouble *DataArrayDouble::applyFunc3(int nbOfComp, const std::vector<std::string>& varsOrder, const std::string& func, bool isSafe) const
 {
+  if(nbOfComp<=0)
+    throw INTERP_KERNEL::Exception("DataArrayDouble::applyFunc3 : output number of component must be > 0 !");
+  std::vector<std::string> varsOrder2(varsOrder);
+  int oldNbOfComp(getNumberOfComponents());
+  for(int i=(int)varsOrder.size();i<oldNbOfComp;i++)
+    varsOrder2.push_back(std::string());
   checkAllocated();
+  int nbOfTuples(getNumberOfTuples());
   INTERP_KERNEL::ExprParser expr(func);
   expr.parse();
   std::set<std::string> vars;
   expr.getTrueSetOfVars(vars);
-  int oldNbOfComp=getNumberOfComponents();
   if((int)vars.size()>oldNbOfComp)
     {
       std::ostringstream oss; oss << "The field has " << oldNbOfComp << " components and there are ";
@@ -4485,29 +4561,49 @@ DataArrayDouble *DataArrayDouble::applyFunc3(int nbOfComp, const std::vector<std
       std::copy(vars.begin(),vars.end(),std::ostream_iterator<std::string>(oss," "));
       throw INTERP_KERNEL::Exception(oss.str().c_str());
     }
-  expr.prepareExprEvaluation(varsOrder,oldNbOfComp,nbOfComp);
-  //
-  DataArrayDouble *newArr=DataArrayDouble::New();
-  int nbOfTuples=getNumberOfTuples();
+  MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> newArr(DataArrayDouble::New());
   newArr->alloc(nbOfTuples,nbOfComp);
-  const double *ptr=getConstPointer();
-  double *ptrToFill=newArr->getPointer();
-  for(int i=0;i<nbOfTuples;i++)
-    {
-      try
-      {
-          expr.evaluateExpr(nbOfComp,ptr+i*oldNbOfComp,ptrToFill+i*nbOfComp);
-      }
-      catch(INTERP_KERNEL::Exception& e)
-      {
-          std::ostringstream oss; oss << "For tuple # " << i << " with value (";
-          std::copy(ptr+oldNbOfComp*i,ptr+oldNbOfComp*(i+1),std::ostream_iterator<double>(oss,", "));
-          oss << ") : Evaluation of function failed !" << e.what();
-          newArr->decrRef();
-          throw INTERP_KERNEL::Exception(oss.str().c_str());
-      }
+  INTERP_KERNEL::AutoPtr<double> buff(new double[oldNbOfComp]);
+  double *buffPtr(buff),*ptrToFill;
+  std::vector<double> stck;
+  for(int iComp=0;iComp<nbOfComp;iComp++)
+    {
+      expr.prepareExprEvaluationDouble(varsOrder2,oldNbOfComp,nbOfComp,iComp,buffPtr,buffPtr+oldNbOfComp);
+      expr.prepareFastEvaluator();
+      const double *ptr(getConstPointer());
+      ptrToFill=newArr->getPointer()+iComp;
+      if(!isSafe)
+        {
+          for(int i=0;i<nbOfTuples;i++,ptrToFill+=nbOfComp,ptr+=oldNbOfComp)
+            {
+              std::copy(ptr,ptr+oldNbOfComp,buffPtr);
+              expr.evaluateDoubleInternal(stck);
+              *ptrToFill=stck.back();
+              stck.pop_back();
+            }
+        }
+      else
+        {
+          for(int i=0;i<nbOfTuples;i++,ptrToFill+=nbOfComp,ptr+=oldNbOfComp)
+            {
+              std::copy(ptr,ptr+oldNbOfComp,buffPtr);
+              try
+              {
+                  expr.evaluateDoubleInternalSafe(stck);
+                  *ptrToFill=stck.back();
+                  stck.pop_back();
+              }
+              catch(INTERP_KERNEL::Exception& e)
+              {
+                  std::ostringstream oss; oss << "For tuple # " << i << " with value (";
+                  std::copy(ptr+oldNbOfComp*i,ptr+oldNbOfComp*(i+1),std::ostream_iterator<double>(oss,", "));
+                  oss << ") : Evaluation of function failed !" << e.what();
+                  throw INTERP_KERNEL::Exception(oss.str().c_str());
+              }
+            }
+        }
     }
-  return newArr;
+  return newArr.retn();
 }
 
 void DataArrayDouble::applyFuncFast32(const std::string& func)
diff --git a/src/MEDCoupling/MEDCouplingMemArray.hxx b/src/MEDCoupling/MEDCouplingMemArray.hxx
index a032dc7f6..bfd822845 100644
--- a/src/MEDCoupling/MEDCouplingMemArray.hxx
+++ b/src/MEDCoupling/MEDCouplingMemArray.hxx
@@ -341,10 +341,11 @@ namespace ParaMEDMEM
     MEDCOUPLING_EXPORT void applyRPow(double val);
     MEDCOUPLING_EXPORT DataArrayDouble *negate() const;
     MEDCOUPLING_EXPORT DataArrayDouble *applyFunc(int nbOfComp, FunctionToEvaluate func) const;
-    MEDCOUPLING_EXPORT DataArrayDouble *applyFunc(int nbOfComp, const std::string& func) const;
-    MEDCOUPLING_EXPORT DataArrayDouble *applyFunc(const std::string& func) const;
-    MEDCOUPLING_EXPORT DataArrayDouble *applyFunc2(int nbOfComp, const std::string& func) const;
-    MEDCOUPLING_EXPORT DataArrayDouble *applyFunc3(int nbOfComp, const std::vector<std::string>& varsOrder, const std::string& func) const;
+    MEDCOUPLING_EXPORT DataArrayDouble *applyFunc(int nbOfComp, const std::string& func, bool isSafe=true) const;
+    MEDCOUPLING_EXPORT DataArrayDouble *applyFunc(const std::string& func, bool isSafe=true) const;
+    MEDCOUPLING_EXPORT void applyFuncOnThis(const std::string& func, bool isSafe=true);
+    MEDCOUPLING_EXPORT DataArrayDouble *applyFunc2(int nbOfComp, const std::string& func, bool isSafe=true) const;
+    MEDCOUPLING_EXPORT DataArrayDouble *applyFunc3(int nbOfComp, const std::vector<std::string>& varsOrder, const std::string& func, bool isSafe=true) const;
     MEDCOUPLING_EXPORT void applyFuncFast32(const std::string& func);
     MEDCOUPLING_EXPORT void applyFuncFast64(const std::string& func);
     MEDCOUPLING_EXPORT DataArrayInt *getIdsInRange(double vmin, double vmax) const;
diff --git a/src/MEDCoupling/Test/MEDCouplingExamplesTest.cxx b/src/MEDCoupling/Test/MEDCouplingExamplesTest.cxx
index c466cea48..9bf0c7f4a 100644
--- a/src/MEDCoupling/Test/MEDCouplingExamplesTest.cxx
+++ b/src/MEDCoupling/Test/MEDCouplingExamplesTest.cxx
@@ -189,8 +189,8 @@ void CppExample_MEDCouplingFieldDouble_applyFunc_same_nb_comp()
   MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> field =
     MEDCouplingFieldDouble::New( ParaMEDMEM::ON_CELLS );
   field->setArray( array );
-  const char func[] = "IVec * v + JVec * v*v + 10";
-  field->applyFunc( func );
+  const char func[] = "IVec * v + JVec * w*w + 10";
+  field->applyFunc( 2, func );
   CPPUNIT_ASSERT( field->getNumberOfComponents() == 2 ); // 2 components remains
   //! [CppSnippet_MEDCouplingFieldDouble_applyFunc_same_nb_comp_1]
   //! [CppSnippet_MEDCouplingFieldDouble_applyFunc_same_nb_comp_2]
diff --git a/src/MEDCoupling_Swig/MEDCouplingBasicsTest.py b/src/MEDCoupling_Swig/MEDCouplingBasicsTest.py
index 691a760c1..19301bb7d 100644
--- a/src/MEDCoupling_Swig/MEDCouplingBasicsTest.py
+++ b/src/MEDCoupling_Swig/MEDCouplingBasicsTest.py
@@ -15894,6 +15894,17 @@ class MEDCouplingBasicsTest(unittest.TestCase):
         self.assertEqual(expRes2, res2.getValues())
         pass
 
+    def testDADApplyFuncOnThis1(self):
+        d=DataArrayDouble(5) ; d.iota(0.)
+        d.applyFuncOnThis("2*x+1")
+        self.assertTrue(d.isEqual(DataArrayDouble([1.,3.,5.,7.,9.]),1e-12))
+        d=DataArrayDouble(6) ; d.iota(0.) ; d.rearrange(2)
+        d.applyFuncOnThis("2*x+1")
+        self.assertTrue(d.isEqual(DataArrayDouble([1.,3.,5.,7.,9.,11.],3,2),1e-12))
+        d.applyFuncOnThis("1+2*3")
+        self.assertTrue(d.isEqual(DataArrayDouble([(7.,7.),(7.,7.),(7.,7.)]),1e-12))
+        pass
+
     pass
 
 if __name__ == '__main__':
diff --git a/src/MEDCoupling_Swig/MEDCouplingExamplesTest.py b/src/MEDCoupling_Swig/MEDCouplingExamplesTest.py
index b0e7647f7..3b9249675 100644
--- a/src/MEDCoupling_Swig/MEDCouplingExamplesTest.py
+++ b/src/MEDCoupling_Swig/MEDCouplingExamplesTest.py
@@ -154,8 +154,8 @@ class MEDCouplingBasicsTest(unittest.TestCase):
         array = DataArrayDouble.New( v, 2, 2 ) # 2 tuples per 2 components
         field = MEDCouplingFieldDouble.New( ON_CELLS )
         field.setArray( array )
-        func = "IVec * v + JVec * v*v + 10"
-        field.applyFunc( func )
+        func = "IVec * v + JVec * w*w + 10"
+        field.applyFunc( 2, func )
         self.assertTrue( field.getNumberOfComponents() == 2 ) # 2 components remains
         #! [PySnippet_MEDCouplingFieldDouble_applyFunc_same_nb_comp_1]
         #! [PySnippet_MEDCouplingFieldDouble_applyFunc_same_nb_comp_2]
@@ -2357,7 +2357,7 @@ class MEDCouplingBasicsTest(unittest.TestCase):
         self.assertTrue(d1.isEqual(DataArrayDouble([1.,4.,121.,144.,441.,484.,961.,1681.],4,2),1e-12))
 # ! [PySnippetDataArrayApplyFunc1_2]
 # ! [PySnippetDataArrayApplyFunc1_3]
-        d2=d.applyFunc("smth*IVec+2*smth*JVec")
+        d2=d.applyFunc(2,"smth1*IVec+2*smth2*JVec")
         self.assertTrue(d2.isEqual(DataArrayDouble([1.,4.,11.,24.,21.,44.,31.,82.],4,2),1e-12))
 # ! [PySnippetDataArrayApplyFunc1_3]
 # ! [PySnippetDataArrayApplyFunc1_4]
diff --git a/src/MEDCoupling_Swig/MEDCouplingMemArray.i b/src/MEDCoupling_Swig/MEDCouplingMemArray.i
index 27f0d46cf..5d2fb64f2 100644
--- a/src/MEDCoupling_Swig/MEDCouplingMemArray.i
+++ b/src/MEDCoupling_Swig/MEDCouplingMemArray.i
@@ -564,10 +564,11 @@ namespace ParaMEDMEM
     void applyRPow(double val) throw(INTERP_KERNEL::Exception);
     DataArrayDouble *negate() const throw(INTERP_KERNEL::Exception);
     DataArrayDouble *applyFunc(int nbOfComp, FunctionToEvaluate func) const throw(INTERP_KERNEL::Exception);
-    DataArrayDouble *applyFunc(int nbOfComp, const std::string& func) const throw(INTERP_KERNEL::Exception);
-    DataArrayDouble *applyFunc(const std::string& func) const throw(INTERP_KERNEL::Exception);
-    DataArrayDouble *applyFunc2(int nbOfComp, const std::string& func) const throw(INTERP_KERNEL::Exception);
-    DataArrayDouble *applyFunc3(int nbOfComp, const std::vector<std::string>& varsOrder, const std::string& func) const throw(INTERP_KERNEL::Exception);
+    DataArrayDouble *applyFunc(int nbOfComp, const std::string& func, bool isSafe=true) const throw(INTERP_KERNEL::Exception);
+    DataArrayDouble *applyFunc(const std::string& func, bool isSafe=true) const throw(INTERP_KERNEL::Exception);
+    void applyFuncOnThis(const std::string& func, bool isSafe=true) throw(INTERP_KERNEL::Exception);
+    DataArrayDouble *applyFunc2(int nbOfComp, const std::string& func, bool isSafe=true) const throw(INTERP_KERNEL::Exception);
+    DataArrayDouble *applyFunc3(int nbOfComp, const std::vector<std::string>& varsOrder, const std::string& func, bool isSafe=true) const throw(INTERP_KERNEL::Exception);
     void applyFuncFast32(const std::string& func) throw(INTERP_KERNEL::Exception);
     void applyFuncFast64(const std::string& func) throw(INTERP_KERNEL::Exception);
     DataArrayInt *getIdsInRange(double vmin, double vmax) const throw(INTERP_KERNEL::Exception);