diff --git a/preprocessor/BlockTriangular.cc b/preprocessor/BlockTriangular.cc
index b9348f6d6..5f9b02c3f 100644
--- a/preprocessor/BlockTriangular.cc
+++ b/preprocessor/BlockTriangular.cc
@@ -108,19 +108,18 @@ BlockTriangular::Allocate_Block(int size, int *count_Equ, int count_Block, Block
   bool *tmp_variable_evaluated;
   bool *Cur_IM;
   bool *IM, OK;
-  ModelBlock->Periods = periods;
   int Lag_Endo, Lead_Endo, Lag_Exo, Lead_Exo, Lag_Other_Endo, Lead_Other_Endo;
 
+  ModelBlock->Periods = periods;
   ModelBlock->Block_List[count_Block].is_linear=true;
   ModelBlock->Block_List[count_Block].Size = size;
   ModelBlock->Block_List[count_Block].Type = type;
-  ModelBlock->Block_List[count_Block].Temporary_terms=new temporary_terms_type ();
-  ModelBlock->Block_List[count_Block].Temporary_terms->clear();
   ModelBlock->Block_List[count_Block].Temporary_InUse=new temporary_terms_inuse_type ();
   ModelBlock->Block_List[count_Block].Temporary_InUse->clear();
   ModelBlock->Block_List[count_Block].Simulation_Type = SimType;
   ModelBlock->Block_List[count_Block].Equation = (int*)malloc(ModelBlock->Block_List[count_Block].Size * sizeof(int));
   ModelBlock->Block_List[count_Block].Variable = (int*)malloc(ModelBlock->Block_List[count_Block].Size * sizeof(int));
+  ModelBlock->Block_List[count_Block].Temporary_Terms_in_Equation = (temporary_terms_type**)malloc(ModelBlock->Block_List[count_Block].Size * sizeof(temporary_terms_type));
   ModelBlock->Block_List[count_Block].Own_Derivative = (int*)malloc(ModelBlock->Block_List[count_Block].Size * sizeof(int));
   Lead = Lag = 0;
   first_count_equ = *count_Equ;
@@ -144,6 +143,8 @@ BlockTriangular::Allocate_Block(int size, int *count_Equ, int count_Block, Block
   memset(tmp_variable_evaluated, 0, symbol_table.endo_nbr*sizeof(bool));
   for (i = 0;i < size;i++)
     {
+      ModelBlock->Block_List[count_Block].Temporary_Terms_in_Equation[i]=new temporary_terms_type ();
+      ModelBlock->Block_List[count_Block].Temporary_Terms_in_Equation[i]->clear();
       ModelBlock->Block_List[count_Block].Equation[i] = Index_Equ_IM[*count_Equ].index;
       ModelBlock->Block_List[count_Block].Variable[i] = Index_Var_IM[*count_Equ].index;
       i_1 = Index_Var_IM[*count_Equ].index;
@@ -453,7 +454,9 @@ BlockTriangular::Free_Block(Model_Block* ModelBlock) const
               }
           }
         free(ModelBlock->Block_List[blk].IM_lead_lag);
-        delete(ModelBlock->Block_List[blk].Temporary_terms);
+        for(i=0; i<ModelBlock->Block_List[blk].Size; i++)
+          delete ModelBlock->Block_List[blk].Temporary_Terms_in_Equation[i];
+        free(ModelBlock->Block_List[blk].Temporary_Terms_in_Equation);
         delete(ModelBlock->Block_List[blk].Temporary_InUse);
       }
     free(ModelBlock->Block_List);
diff --git a/preprocessor/ExprNode.cc b/preprocessor/ExprNode.cc
index e9433b37a..79b1b9ac5 100644
--- a/preprocessor/ExprNode.cc
+++ b/preprocessor/ExprNode.cc
@@ -84,9 +84,10 @@ ExprNode::computeTemporaryTerms(map<NodeID, int> &reference_count,
 void
 ExprNode::computeTemporaryTerms(map<NodeID, int> &reference_count,
                                 temporary_terms_type &temporary_terms,
-                                map<NodeID, int> &first_occurence,
+                                map<NodeID, pair<int, int> > &first_occurence,
                                 int Curr_block,
                                 Model_Block *ModelBlock,
+                                int equation,
                                 map_idx_type &map_idx) const
 {
   // Nothing to do for a terminal node
@@ -148,23 +149,12 @@ NumConstNode::eval(const eval_context_type &eval_context) const throw (EvalExcep
 void
 NumConstNode::compile(ofstream &CompileCode, bool lhs_rhs, ExprNodeOutputType output_type, const temporary_terms_type &temporary_terms, map_idx_type map_idx) const
 {
-  //CompileCode.write(reinterpret_cast<char *>(&FLDT), sizeof(FLDT));
-  /*temporary_terms_type::const_iterator it = temporary_terms.find(const_cast<NumConstNode *>(this));
-  if (it != temporary_terms.end())
-    {
-      CompileCode.write(&FLDT, sizeof(FLDT));
-      idl=
-      CompileCode.write(reinterpret_cast<char *>(&idl),sizeof(idl));
-    }
-  else
-    {*/
-      CompileCode.write(&FLDC, sizeof(FLDC));
-      double vard=atof(datatree.num_constants.get(id).c_str());
+  CompileCode.write(&FLDC, sizeof(FLDC));
+  double vard=atof(datatree.num_constants.get(id).c_str());
 #ifdef DEBUGC
-      cout << "FLDC " << vard << "\n";
+  cout << "FLDC " << vard << "\n";
 #endif
-      CompileCode.write(reinterpret_cast<char *>(&vard),sizeof(vard));
-    /*}*/
+  CompileCode.write(reinterpret_cast<char *>(&vard),sizeof(vard));
 }
 
 
@@ -427,15 +417,9 @@ VariableNode::writeOutput(ostream &output, ExprNodeOutputType output_type,
 double
 VariableNode::eval(const eval_context_type &eval_context) const throw (EvalException)
 {
-  // ModelTree::evaluateJacobian need to have the initval values applied to lead/lagged variables also
-  /*if (lag != 0)
-    throw EvalException();*/
-  /*if(type==eModelLocalVariable)
-    cout << "eModelLocalVariable = " << symb_id << "\n";*/
   eval_context_type::const_iterator it = eval_context.find(make_pair(symb_id, type));
   if (it == eval_context.end())
     {
-      //cout << "unknonw variable type = " << type << "  simb_id = " << symb_id << "\n";
       throw EvalException();
     }
 
@@ -445,15 +429,6 @@ VariableNode::eval(const eval_context_type &eval_context) const throw (EvalExcep
 void
 VariableNode::compile(ofstream &CompileCode, bool lhs_rhs, ExprNodeOutputType output_type, const temporary_terms_type &temporary_terms, map_idx_type map_idx) const
 {
-  // If node is a temporary term
-  /*temporary_terms_type::const_iterator it = temporary_terms.find(const_cast<VariableNode *>(this));
-  if (it != temporary_terms.end())
-    {
-      CompileCode.write(&FLDT, sizeof(FLDT));
-      int var=temporary_terms.count(const_cast<VariableNode *>(this))-1;
-      CompileCode.write(reinterpret_cast<char *>(&var), sizeof(var));
-      return;
-    }*/
   int i, lagl;
 #ifdef DEBUGC
   cout << "output_type=" << output_type << "\n";
@@ -501,6 +476,23 @@ VariableNode::compile(ofstream &CompileCode, bool lhs_rhs, ExprNodeOutputType ou
     }
 }
 
+
+void
+VariableNode::computeTemporaryTerms(map<NodeID, int> &reference_count,
+                                   temporary_terms_type &temporary_terms,
+                                   map<NodeID, pair<int, int> > &first_occurence,
+                                   int Curr_block,
+                                   Model_Block *ModelBlock,
+                                   int equation,
+                                   map_idx_type &map_idx) const
+{
+  if(type== eModelLocalVariable)
+    datatree.local_variables_table[symb_id]->computeTemporaryTerms(reference_count, temporary_terms, first_occurence, Curr_block, ModelBlock, equation, map_idx);
+}
+
+
+
+
 void
 VariableNode::collectEndogenous(set<pair<int, int> > &result) const
 {
@@ -711,9 +703,10 @@ UnaryOpNode::computeTemporaryTerms(map<NodeID, int> &reference_count,
 void
 UnaryOpNode::computeTemporaryTerms(map<NodeID, int> &reference_count,
                                    temporary_terms_type &temporary_terms,
-                                   map<NodeID, int> &first_occurence,
+                                   map<NodeID, pair<int, int> > &first_occurence,
                                    int Curr_block,
                                    Model_Block *ModelBlock,
+                                   int equation,
                                    map_idx_type &map_idx) const
 {
   NodeID this2 = const_cast<UnaryOpNode *>(this);
@@ -721,8 +714,8 @@ UnaryOpNode::computeTemporaryTerms(map<NodeID, int> &reference_count,
   if (it == reference_count.end())
     {
       reference_count[this2] = 1;
-      first_occurence[this2] = Curr_block;
-      arg->computeTemporaryTerms(reference_count, temporary_terms, first_occurence, Curr_block, ModelBlock, map_idx);
+      first_occurence[this2] = make_pair(Curr_block,equation);
+      arg->computeTemporaryTerms(reference_count, temporary_terms, first_occurence, Curr_block, ModelBlock, equation, map_idx);
     }
   else
     {
@@ -730,7 +723,7 @@ UnaryOpNode::computeTemporaryTerms(map<NodeID, int> &reference_count,
       if (reference_count[this2] * cost(temporary_terms, false) > MIN_COST_C)
         {
           temporary_terms.insert(this2);
-          ModelBlock->Block_List[first_occurence[this2]].Temporary_terms->insert(this2);
+          ModelBlock->Block_List[first_occurence[this2].first].Temporary_Terms_in_Equation[first_occurence[this2].second]->insert(this2);
         }
     }
 }
@@ -1151,9 +1144,10 @@ BinaryOpNode::computeTemporaryTerms(map<NodeID, int> &reference_count,
 void
 BinaryOpNode::computeTemporaryTerms(map<NodeID, int> &reference_count,
                                     temporary_terms_type &temporary_terms,
-                                    map<NodeID, int> &first_occurence,
+                                    map<NodeID, pair<int, int> > &first_occurence,
                                     int Curr_block,
                                     Model_Block *ModelBlock,
+                                    int equation,
                                     map_idx_type &map_idx) const
 {
   NodeID this2 = const_cast<BinaryOpNode *>(this);
@@ -1161,17 +1155,19 @@ BinaryOpNode::computeTemporaryTerms(map<NodeID, int> &reference_count,
   if (it == reference_count.end())
     {
       reference_count[this2] = 1;
-      first_occurence[this2] = Curr_block;
-      arg1->computeTemporaryTerms(reference_count, temporary_terms, first_occurence, Curr_block, ModelBlock, map_idx);
-      arg2->computeTemporaryTerms(reference_count, temporary_terms, first_occurence, Curr_block, ModelBlock, map_idx);
+      first_occurence[this2] = make_pair(Curr_block, equation);
+      arg1->computeTemporaryTerms(reference_count, temporary_terms, first_occurence, Curr_block, ModelBlock, equation, map_idx);
+      arg2->computeTemporaryTerms(reference_count, temporary_terms, first_occurence, Curr_block, ModelBlock, equation, map_idx);
     }
   else
     {
       reference_count[this2]++;
       if (reference_count[this2] * cost(temporary_terms, false) > MIN_COST_C)
         {
+          if(this2->idx==2280)
+            cout << "==>Curr_block= " << Curr_block << " equation= " << equation << " first_occurence[this2].first=" << first_occurence[this2].first << " first_occurence[this2].second=" << first_occurence[this2].second << "\n";
           temporary_terms.insert(this2);
-          ModelBlock->Block_List[first_occurence[this2]].Temporary_terms->insert(this2);
+          ModelBlock->Block_List[first_occurence[this2].first].Temporary_Terms_in_Equation[first_occurence[this2].second]->insert(this2);
         }
     }
 }
@@ -1563,9 +1559,10 @@ TrinaryOpNode::computeTemporaryTerms(map<NodeID, int> &reference_count,
 void
 TrinaryOpNode::computeTemporaryTerms(map<NodeID, int> &reference_count,
                                     temporary_terms_type &temporary_terms,
-                                    map<NodeID, int> &first_occurence,
+                                    map<NodeID, pair<int, int> > &first_occurence,
                                     int Curr_block,
                                     Model_Block *ModelBlock,
+                                    int equation,
                                     map_idx_type &map_idx) const
 {
   NodeID this2 = const_cast<TrinaryOpNode *>(this);
@@ -1573,18 +1570,20 @@ TrinaryOpNode::computeTemporaryTerms(map<NodeID, int> &reference_count,
   if (it == reference_count.end())
     {
       reference_count[this2] = 1;
-      first_occurence[this2] = Curr_block;
-      arg1->computeTemporaryTerms(reference_count, temporary_terms, first_occurence, Curr_block, ModelBlock, map_idx);
-      arg2->computeTemporaryTerms(reference_count, temporary_terms, first_occurence, Curr_block, ModelBlock, map_idx);
-      arg3->computeTemporaryTerms(reference_count, temporary_terms, first_occurence, Curr_block, ModelBlock, map_idx);
+      first_occurence[this2] = make_pair(Curr_block,equation);
+      arg1->computeTemporaryTerms(reference_count, temporary_terms, first_occurence, Curr_block, ModelBlock, equation, map_idx);
+      arg2->computeTemporaryTerms(reference_count, temporary_terms, first_occurence, Curr_block, ModelBlock, equation, map_idx);
+      arg3->computeTemporaryTerms(reference_count, temporary_terms, first_occurence, Curr_block, ModelBlock, equation, map_idx);
     }
   else
     {
       reference_count[this2]++;
       if (reference_count[this2] * cost(temporary_terms, false) > MIN_COST_C)
         {
+          if(this2->idx==2280)
+            cout << "==>Curr_block= " << Curr_block << " equation= " << equation << " first_occurence[this2].first=" << first_occurence[this2].first << " first_occurence[this2].second=" << first_occurence[this2].second << "\n";
           temporary_terms.insert(this2);
-          ModelBlock->Block_List[first_occurence[this2]].Temporary_terms->insert(this2);
+          ModelBlock->Block_List[first_occurence[this2].first].Temporary_Terms_in_Equation[first_occurence[this2].second]->insert(this2);
         }
     }
 }
@@ -1743,9 +1742,10 @@ void UnknownFunctionNode::writeOutput(ostream &output, ExprNodeOutputType output
 void
 UnknownFunctionNode::computeTemporaryTerms(map<NodeID, int> &reference_count,
                                            temporary_terms_type &temporary_terms,
-                                           map<NodeID, int> &first_occurence,
+                                           map<NodeID, pair<int, int> > &first_occurence,
                                            int Curr_block,
                                            Model_Block *ModelBlock,
+                                           int equation,
                                            map_idx_type &map_idx) const
 {
   cerr << "UnknownFunctionNode::computeTemporaryTerms: not implemented" << endl;
diff --git a/preprocessor/MatlabFile.cc b/preprocessor/MatlabFile.cc
new file mode 100644
index 000000000..7fff2b498
--- /dev/null
+++ b/preprocessor/MatlabFile.cc
@@ -0,0 +1,1076 @@
+/*
+ * Copyright (C) 2006-2008 Dynare Team
+ *
+ * This file is part of Dynare.
+ *
+ * Dynare is free software: you can redistribute it and/or modify
+ * it under the terms of the GNU General Public License as published by
+ * the Free Software Foundation, either version 3 of the License, or
+ * (at your option) any later version.
+ *
+ * Dynare is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+ * GNU General Public License for more details.
+ *
+ * You should have received a copy of the GNU General Public License
+ * along with Dynare.  If not, see <http://www.gnu.org/licenses/>.
+ */
+/*
+Usefull documentation: Matlab 7 Mat-File Format
+-----------------------------------------------
+revision: October 2008 PDF only Rereleased for Version 7.7 (Release 2008b)
+available at: http://www.mathworks.com/access/helpdesk/help/pdf_doc/matlab/matfile_format.pdf
+*/
+#include "MatlabFile.hh"
+
+
+MatlabFile::MatlabFile()
+{
+}
+
+MatlabFile::~MatlabFile()
+{
+}
+
+ArrayElem::ArrayElem()
+{
+}
+
+ArrayElem::~ArrayElem()
+{
+}
+
+
+string
+UTF8::ReadAlph(char* InBuff, int* pBuff, int Size) const
+{
+  char tmp_c[Size+1];
+  string str;
+  memcpy(&tmp_c, InBuff+*pBuff, Size);
+  tmp_c[Size]=0;
+  str.assign(tmp_c);
+  if(Size<=4)
+    *pBuff += 4;
+  else if(Size % 8)
+    *pBuff += 8*ceil(double(Size) / 8);
+  else
+    *pBuff += Size;
+  return(str);
+}
+
+
+string
+UTF16::ReadAlph(char* InBuff, int* pBuff, int Size) const
+{
+  string str("");
+  for(int i=0;i<Size;i++)
+    str.append(&(InBuff[i*2+*pBuff]));
+  if(Size*2<=4)
+    *pBuff += 4;
+  else if((Size*2) % 8)
+    *pBuff += 8*ceil(double(Size*2) / 8);
+  else
+    *pBuff += Size*2;
+  return(str);
+}
+
+string
+UTF32::ReadAlph(char* InBuff, int* pBuff, int Size) const
+{
+  string str("");
+  for(int i=0;i<Size;i++)
+    str.append(&(InBuff[i*4+*pBuff]));
+  if((Size*4) % 8)
+    *pBuff += 8*ceil(double(Size*4) / 8);
+  else
+    *pBuff += Size*4;
+  return(str);
+}
+
+double
+INT8::ReadNum(char* InBuff, int* pBuff) const
+{
+  char val;
+  val = InBuff[*pBuff];
+  *pBuff += sizeof(val);
+  return(val);
+}
+
+
+double
+INT16::ReadNum(char* InBuff, int* pBuff) const
+{
+  Int val;
+  memcpy(&val, InBuff+*pBuff, sizeof(val));
+  *pBuff += sizeof(val);
+  return(val);
+}
+
+double
+INT32::ReadNum(char* InBuff, int* pBuff) const
+{
+  LongInt val;
+  memcpy(&val, InBuff+*pBuff, sizeof(val));
+  *pBuff += sizeof(val);
+  return(val);
+}
+
+
+double
+INT64::ReadNum(char* InBuff, int* pBuff) const
+{
+  LongLongInt val;
+  memcpy(&val, InBuff+*pBuff, sizeof(val));
+  *pBuff += sizeof(val);
+  return(val);
+}
+
+double
+Single::ReadNum(char* InBuff, int* pBuff) const
+{
+  float val;
+  memcpy(&val, InBuff+*pBuff, sizeof(val));
+  *pBuff += sizeof(val);
+  return(val);
+}
+
+double
+Double::ReadNum(char* InBuff, int* pBuff) const
+{
+  double val;
+  memcpy(&val, InBuff+*pBuff, sizeof(val));
+  *pBuff += sizeof(val);
+  return(val);
+}
+
+
+
+Data_Header_t
+SimpleElem::ReadDataHeader(char* InBuff, int* pBuff) const
+{
+  Data_Header_t data_header;
+  memcpy(&data_header.DataType, InBuff+*pBuff, sizeof(data_header.DataType));
+  *pBuff += 2;
+  memcpy(&data_header.S_Number_of_Bytes, InBuff+*pBuff, sizeof(data_header.S_Number_of_Bytes));
+  *pBuff += 2;
+  if(data_header.S_Number_of_Bytes!=0)
+    data_header.Number_of_Bytes=data_header.S_Number_of_Bytes;
+  else
+    {
+      memcpy(&data_header.Number_of_Bytes, InBuff+*pBuff, sizeof(data_header.Number_of_Bytes));
+      *pBuff += 4;
+    }
+  return(data_header);
+}
+
+
+SimpleElem::SimpleElem()
+{
+  verbose=true;
+  array_elem=NULL;
+  Type=-1;
+}
+
+SimpleElem::~SimpleElem()
+{
+}
+
+returned_ReadData_t
+SimpleElem::Get_Data_Class(Data_Header_t data_header) const
+{
+  returned_ReadData_t ret;
+  switch(data_header.DataType)
+    {
+      case   miINT8:       //8 bit, signed
+      case   miUINT8:      //8 bit, unsigned
+        ret.Type = Numerical;
+        ret.Simple = new INT8;
+        return(ret);
+        break;
+      case   miINT16:      //16-bit, signed
+      case   miUINT16:     //16-bit, unsigned
+        ret.Type = Numerical;
+        ret.Simple = new INT16;
+        return(ret);
+        break;
+      case   miINT32:      //32-bit, signed
+      case   miUINT32:     //32-bit, unsigned
+        ret.Type = Numerical;
+        ret.Simple = new INT32;
+        return(ret);
+        break;
+      case   miINT64:      //64-bit, signed
+      case   miUINT64:     //64-bit, unsigned
+        ret.Type = Numerical;
+        ret.Simple = new INT64;
+        return(ret);
+        break;
+      case   miSINGLE:     //IEEEŽ 754 single format
+        ret.Type = Numerical;
+        ret.Simple = new Single;
+        return(ret);
+        break;
+      case   miDOUBLE:     //IEEE 754 double format
+        ret.Type = Numerical;
+        ret.Simple = new Double;
+        return(ret);
+        break;
+      case   miMATRIX:     //MATLAB array
+        ret.Type = Matrix;
+        ret.Simple = NULL;
+        return(ret);
+        break;
+      case   miCOMPRESSED: //Compressed Data
+        ret.Type = Compressed;
+        ret.Simple = NULL;
+        return(ret);
+        break;
+      case   miUTF8:       //Unicode UTF-8 Encoded Character Data
+        ret.Type = AlphaNumeric;
+        ret.Simple = new UTF8;
+        return(ret);
+        break;
+      case   miUTF16:      //Unicode UTF-16 Encoded Character Data
+        ret.Type = AlphaNumeric;
+        ret.Simple = new UTF16;
+        return(ret);
+        break;
+      case   miUTF32:      //Unicode UTF-32 Encoded Character Data
+        ret.Type = AlphaNumeric;
+        ret.Simple = new UTF32;
+        return(ret);
+        break;
+      default:
+        ret.Type = Unknown;
+        ret.Simple = NULL;
+        return(ret);
+    }
+}
+
+void
+SimpleElem::DataProceed(Data_Header data_header, char* InBuff, int* pBuff, FlagStructure_t flag)
+{
+  ArrayElem matrix;
+  returned_ReadData ret;
+  ret = Get_Data_Class(data_header);
+  Type = ret.Type;
+  double tmpv;
+  switch(ret.Type)
+    {
+      case Numerical:
+        if(data_header.Number_of_Bytes/ret.Simple->size())
+          for(unsigned int i=0;i<data_header.Number_of_Bytes/ret.Simple->size();i++)
+            {
+              tmpv=ret.Simple->ReadNum(InBuff, pBuff);
+              VNumeric.push_back(tmpv);
+            }
+        //to align pBuff on a 4 Bytes
+        if(*pBuff % 4)
+          *pBuff += 4-(*pBuff % 4);
+        delete ret.Simple;
+        break;
+      case AlphaNumeric:
+        for(unsigned int i=0;i<data_header.Number_of_Bytes/ret.Simple->size();i++)
+          Vstr.push_back(ret.Simple->ReadAlph(InBuff, pBuff, data_header.Number_of_Bytes));
+        //to align pBuff on a 4 Bytes
+        if(*pBuff % 4)
+          *pBuff += 4-(*pBuff % 4);
+        delete ret.Simple;
+        break;
+      case Matrix:
+        array_elem = matrix.ReadArray_class(InBuff, pBuff);
+        array_elem->ReadArray(InBuff, pBuff, flag);
+        break;
+      case Compressed:
+        cerr << "Error: Compressed data in Mat-file not implemnted yet!\n set option -v6 when saving to Mat-file\n";
+        exit(EXIT_FAILURE);
+      case Unknown:
+        cerr << "Error: Mat-file format use incomptible format with Matlab 7 specification!\n set option -v6 when saving to Mat-file\n";
+        exit(EXIT_FAILURE);
+    }
+}
+
+returned_ReadData_t
+SimpleElem::ReadData(char* InBuff, int* pBuff) const
+{
+  Data_Header_t data_header;
+  data_header = ReadDataHeader(InBuff, pBuff);
+  return(Get_Data_Class(data_header));
+}
+
+void
+SimpleElem::Collect(const string &name, bool found, CollectStruct &collect_struct) const
+{
+  if(VNumeric.size())
+    {
+      if(found)
+        collect_struct.variable_double_name[collect_struct.tmp_name]=VNumeric;
+    }
+  else if(Vstr.size())
+    {
+      if(found)
+        collect_struct.variable_string_name[collect_struct.tmp_name]=Vstr;
+    }
+  else
+    {
+      if(array_elem)
+        {
+          array_elem->Collect(name, found, collect_struct);
+        }
+   }
+}
+
+
+void
+SimpleElem::Print() const
+{
+  if(VNumeric.size())
+    {
+      for(vector<double>::const_iterator it=VNumeric.begin(); it!=VNumeric.end(); it++)
+        cout << " " << *it;
+    }
+  else if(Vstr.size())
+    {
+      for(vector<string>::const_iterator it=Vstr.begin(); it!=Vstr.end(); it++)
+        cout << " " << *it;
+    }
+  else
+    {
+      if(array_elem)
+        array_elem->Print();
+    }
+}
+
+void
+SimpleElem::Delete() const
+{
+  if(array_elem)
+    {
+      array_elem->Delete();
+      delete array_elem;
+    }
+}
+
+LongInt
+ArrayElem::ReadINT32(char* InBuff, int* pBuff) const
+{
+  LongInt val;
+  memcpy(&val, InBuff+*pBuff, sizeof(val));
+  *pBuff += sizeof(val);
+  return(val);
+}
+
+
+Array_Flag_t
+ArrayElem::ReadArrayFlag(char* InBuff, int* pBuff) /*const*/
+{
+  Array_Flag_t array_flag;
+  memcpy(&array_flag, InBuff+*pBuff, sizeof(array_flag));
+  *pBuff += sizeof(array_flag);
+  array_complex = (bool)(array_flag.flag & 16);
+  array_global = array_flag.flag & 32;
+  array_logical = array_flag.flag & 64;
+  type = (Array_Type)array_flag.classe;
+  if(type==Sparse_array)
+    array_nzmax = array_flag.nzmax;
+  return(array_flag);
+}
+
+void
+ArrayElem::ReadArrayDimension(char* InBuff, int* pBuff) /*const*/
+{
+  Data_Header_t data_header;
+  data_header = ReadDataHeader(InBuff, pBuff);
+  number_of_dimensions = data_header.Number_of_Bytes/4;
+  for(int i=0; i<number_of_dimensions; i++)
+    {
+      double tmp_d=ReadINT32(InBuff, pBuff);
+      dimension.push_back(tmp_d);
+    }
+  if(number_of_dimensions % 2)
+    *pBuff += sizeof(dimension[0]);
+}
+
+void
+ArrayElem::ReadStructureNames(char* InBuff, int* pBuff)
+{
+  Data_Header_t data_header, data_header_2;
+  LongInt Field_name_length;
+  data_header=ReadDataHeader(InBuff, pBuff);
+  Field_name_length=ReadINT32(InBuff, pBuff);
+  data_header_2=ReadDataHeader(InBuff, pBuff);
+  Structure_number = data_header_2.Number_of_Bytes/Field_name_length;
+  char tmp_c[Field_name_length];
+  for(int i=0; i<Structure_number;i++)
+    {
+      memcpy(tmp_c, InBuff+*pBuff, Field_name_length);
+      *pBuff += Field_name_length;
+      string variable_name(tmp_c);
+      Structure_Elem_name.push_back(variable_name);
+    }
+}
+
+void
+ArrayElem::ReadArrayName(char* InBuff, int* pBuff) /*const*/
+{
+  Data_Header_t data_header;
+  data_header = ReadDataHeader(InBuff, pBuff);
+  char tmp_c[data_header.Number_of_Bytes+1];
+  memcpy(&tmp_c, InBuff+*pBuff, data_header.Number_of_Bytes);
+  tmp_c[data_header.Number_of_Bytes]=0;
+  variable_name.assign(tmp_c);
+  if(data_header.Number_of_Bytes<=4)
+    *pBuff += 4;
+  else if(data_header.Number_of_Bytes % 8)
+    *pBuff += 8*ceil(double(data_header.Number_of_Bytes) / 8);
+  else
+    *pBuff += data_header.Number_of_Bytes;
+}
+
+
+ArrayElem*
+ArrayElem::ReadArray_class(char* InBuff, int* pBuff) const
+{
+  Array_Flag_t array_flag;
+  ArrayElem array_elem;
+  array_flag=array_elem.ReadArrayFlag(InBuff, pBuff);
+  switch(array_flag.classe)
+    {
+       case Cell_array:
+         return(new CellArray);
+         break;
+       case Structure_:
+         return(new Structure);
+         break;
+       case Object_:
+         return(new Object);
+         break;
+       case Character_array:
+         return(new CharacterArray);
+         break;
+       case Sparse_array:
+         return(new SparseArray);
+         break;
+       case Double_precision_array:
+         return(new DoublePrecisionArray);
+         break;
+       case Single_precision_array:
+         return(new SinglePrecisionArray);
+         break;
+       case Signed_integer_8_bit:
+         return(new Bit8SignedInteger);
+         break;
+       case Unsigned_integer_8_bit:
+         return(new Bit8UnsignedInteger);
+         break;
+       case Signed_integer_16_bit:
+         return(new Bit16SignedInteger);
+         break;
+       case Unsigned_integer_16_bit:
+         return(new Bit16UnsignedInteger);
+         break;
+       case Signed_integer_32_bit:
+         return(new Bit32SignedInteger);
+         break;
+       case Unsigned_integer_32_bit:
+         return(new Bit32UnsignedInteger);
+         break;
+       default:
+         return(NULL);
+    }
+  return(NULL);
+}
+
+void
+ArrayElem::Collect(const string &name, bool found, CollectStruct &collect_struct) const
+{
+}
+
+void
+ArrayElem::Print() const
+{
+}
+
+void
+ArrayElem::Delete() const
+{
+}
+
+void
+CellArray::Collect(const string &name, bool found, CollectStruct &collect_struct) const
+{
+}
+
+void
+CellArray::Print() const
+{
+  //cout << "CellArray: "<< variable_name << "\n";
+}
+
+void
+CellArray::Delete() const
+{
+  //cout << "CellArray: "<< variable_name << "\n";
+}
+
+void
+CellArray::ReadArray(char* InBuff, int* pBuff, FlagStructure_t flag)
+{
+  SimpleElem* simple;
+  ReadArrayDimension(InBuff, pBuff);
+  Cell_number = 1;
+  for(unsigned int i=0;i<dimension.size();i++)
+    Cell_number *= dimension[i];
+  ReadArrayName(InBuff, pBuff);
+  flag.character=true;
+  for(int i=0;i<Cell_number;i++)
+    {
+      simple=new SimpleElem;
+      Data_Header_t data_header=simple->ReadDataHeader(InBuff, pBuff);
+      VCell.push_back(simple);
+      simple->DataProceed(data_header, InBuff, pBuff, flag);
+    }
+}
+
+
+void
+Structure::Collect(const string &name, bool found, CollectStruct &collect_struct) const
+{
+  if(name==variable_name || found)
+    {
+      found = true;
+      vector<string>::const_iterator it2=Structure_Elem_name.begin();
+      for(vector<PSimpleElem>::const_iterator it=VCell.begin(); it!=VCell.end(); it++)
+        {
+          collect_struct.tmp_name = *it2;
+          it2++;
+          (*it)->Collect(name, found, collect_struct);
+        }
+    }
+}
+
+void
+Structure::Print() const
+{
+  cout << "Structure: " << variable_name << "\n";
+  vector<string>::const_iterator it2=Structure_Elem_name.begin();
+  int i=0;
+  for(vector<PSimpleElem>::const_iterator it=VCell.begin(); it!=VCell.end(); it++)
+    {
+      cout << ++i << " -> " << *it2 << " :";
+      it2++;
+      (*it)->Print();
+      cout << "\n";
+    }
+}
+
+void
+Structure::Delete() const
+{
+  vector<string>::const_iterator it2=Structure_Elem_name.begin();
+  for(vector<PSimpleElem>::const_iterator it=VCell.begin(); it!=VCell.end(); it++)
+    {
+      (*it)->Delete();
+      delete *it;
+    }
+}
+
+void
+Structure::ReadArray(char* InBuff, int* pBuff, FlagStructure_t flag)
+{
+  SimpleElem* simple;
+  int i;
+  ReadArrayDimension(InBuff, pBuff);
+  ReadArrayName(InBuff, pBuff);
+  ReadStructureNames(InBuff, pBuff);
+  flag.no_name=true;
+  for(i=0;i<Structure_number;i++)
+    {
+      simple=new SimpleElem;
+      Data_Header_t data_header=simple->ReadDataHeader(InBuff, pBuff);
+      simple->DataProceed(data_header, InBuff, pBuff, flag);
+      data_header=simple->ReadDataHeader(InBuff, pBuff);
+      VCell.push_back(simple);
+      simple->DataProceed(data_header, InBuff, pBuff, flag);
+    }
+}
+
+void
+Object::Collect(const string &name, bool found, CollectStruct &collect_struct) const
+{
+}
+
+
+void
+Object::Print() const
+{
+}
+
+void
+Object::Delete() const
+{
+}
+
+
+void
+Object::ReadArray(char* InBuff, int* pBuff, FlagStructure_t flag)
+{
+  cerr << "Error: Object not implemented\n";
+  exit(EXIT_FAILURE);
+}
+
+void
+CharacterArray::Collect(const string &name, bool found, CollectStruct &collect_struct) const
+{
+  if(name==variable_name || found)
+    {
+      found = true;
+      collect_struct.tmp_name=variable_name;
+      vector<PSimpleElem>::const_iterator it=VCell.begin();
+      (*it)->Collect(name, found, collect_struct);
+    }
+}
+
+void
+CharacterArray::Print() const
+{
+  cout << "CharacterArray: " << variable_name << "\n";
+  vector<PSimpleElem>::const_iterator it=VCell.begin();
+  (*it)->Print();
+  cout << "\n";
+}
+
+void
+CharacterArray::Delete() const
+{
+  //cout << "CharacterArray: " << variable_name << "\n";
+  vector<PSimpleElem>::const_iterator it=VCell.begin();
+  (*it)->Delete();
+  delete *it;
+  //cout << "\n";
+}
+
+
+void
+CharacterArray::ReadArray(char* InBuff, int* pBuff, FlagStructure_t flag)
+{
+  Data_Header_t data_header;
+  SimpleElem* simple;
+  ReadArrayDimension(InBuff, pBuff);
+  Matrix_Elem_number = 1;
+  for(unsigned int i=0;i<dimension.size();i++)
+    {
+      Matrix_Elem_number *= dimension[i];
+    }
+  simple=new SimpleElem;
+  VCell.push_back(simple);
+  data_header=simple->ReadDataHeader(InBuff, pBuff);
+  simple->DataProceed(data_header, InBuff, pBuff, flag);
+}
+
+void
+SparseArray::Collect(const string &name, bool found, CollectStruct &collect_struct) const
+{
+}
+
+void
+SparseArray::Print() const
+{
+  cout << "Sparse Array: " << variable_name << "\n";
+}
+
+void
+SparseArray::Delete() const
+{
+  //cout << "Sparse Array: " << variable_name << "\n";
+}
+
+
+void
+SparseArray::ReadArray(char* InBuff, int* pBuff, FlagStructure_t flag)
+{
+  cerr << "Error: Sparse Array not implemented\n";
+  exit(EXIT_FAILURE);
+}
+
+void
+DoublePrecisionArray::Collect(const string &name, bool found, CollectStruct &collect_struct) const
+{
+  if(name==variable_name || found)
+    {
+      found = true;
+      collect_struct.tmp_name=variable_name;
+      vector<PSimpleElem>::const_iterator it=VCell.begin();
+      (*it)->Collect(name, found, collect_struct);
+    }
+}
+
+
+void
+DoublePrecisionArray::Print() const
+{
+  cout << "DoublePrecisionArray: " << variable_name << "\n";
+  vector<PSimpleElem>::const_iterator it=VCell.begin();
+  (*it)->Print();
+  cout << "\n";
+}
+
+void
+DoublePrecisionArray::Delete() const
+{
+  vector<PSimpleElem>::const_iterator it=VCell.begin();
+  (*it)->Delete();
+  delete *it;
+}
+
+void
+DoublePrecisionArray::ReadArray(char* InBuff, int* pBuff, FlagStructure_t flag)
+{
+  Data_Header_t data_header;
+  SimpleElem* simple;
+  ReadArrayDimension(InBuff, pBuff);
+  Matrix_Elem_number = 1;
+  for(unsigned int i=0;i<dimension.size();i++)
+    {
+      Matrix_Elem_number *= dimension[i];
+    }
+  if(!flag.no_name)
+    ReadArrayName(InBuff, pBuff);
+  simple=new SimpleElem;
+  VCell.push_back(simple);
+  data_header=simple->ReadDataHeader(InBuff, pBuff);
+  simple->DataProceed(data_header, InBuff, pBuff, flag);
+}
+
+void
+SinglePrecisionArray::Collect(const string &name, bool found, CollectStruct &collect_struct) const
+{
+  if(name==variable_name || found)
+    {
+      found = true;
+      collect_struct.tmp_name=variable_name;
+      vector<PSimpleElem>::const_iterator it=VCell.begin();
+      (*it)->Collect(name, found, collect_struct);
+    }
+}
+
+
+void
+SinglePrecisionArray::Print() const
+{
+  cout << "SinglePrecisionArray: " << variable_name << "\n";
+  vector<PSimpleElem>::const_iterator it=VCell.begin();
+  (*it)->Print();
+  cout << "\n";
+}
+
+void
+SinglePrecisionArray::Delete() const
+{
+  vector<PSimpleElem>::const_iterator it=VCell.begin();
+  (*it)->Delete();
+  delete *it;
+}
+
+void
+SinglePrecisionArray::ReadArray(char* InBuff, int* pBuff, FlagStructure_t flag)
+{
+  Data_Header_t data_header;
+  SimpleElem* simple;
+  ReadArrayDimension(InBuff, pBuff);
+  Matrix_Elem_number = 1;
+  for(unsigned int i=0;i<dimension.size();i++)
+    {
+      Matrix_Elem_number *= dimension[i];
+    }
+  if(!flag.no_name)
+    ReadArrayName(InBuff, pBuff);
+  simple=new SimpleElem;
+  VCell.push_back(simple);
+  data_header=simple->ReadDataHeader(InBuff, pBuff);
+  simple->DataProceed(data_header, InBuff, pBuff, flag);
+}
+
+void
+Bit8SignedInteger::Collect(const string &name, bool found, CollectStruct &collect_struct) const
+{
+}
+
+
+
+void
+Bit8SignedInteger::Print() const
+{
+  //cout << "Bit8SignedInteger: \n";
+}
+
+void
+Bit8SignedInteger::Delete() const
+{
+  //cout << "Bit8SignedInteger: \n";
+}
+
+void
+Bit8SignedInteger::ReadArray(char* InBuff, int* pBuff, FlagStructure_t flag)
+{
+  Data_Header_t data_header;
+  SimpleElem* simple;
+  simple=new SimpleElem;
+  VCell.push_back(simple);
+  data_header=simple->ReadDataHeader(InBuff, pBuff);
+  simple->DataProceed(data_header, InBuff, pBuff, flag);
+}
+
+void
+Bit8UnsignedInteger::Collect(const string &name, bool found, CollectStruct &collect_struct) const
+{
+}
+
+
+void
+Bit8UnsignedInteger::Print() const
+{
+  //cout << "Bit8UnsignedInteger: \n";
+}
+
+void
+Bit8UnsignedInteger::Delete() const
+{
+  //cout << "Bit8UnsignedInteger: \n";
+}
+
+void
+Bit8UnsignedInteger::ReadArray(char* InBuff, int* pBuff, FlagStructure_t flag)
+{
+  Data_Header_t data_header;
+  SimpleElem* simple;
+  simple=new SimpleElem;
+  VCell.push_back(simple);
+  data_header=simple->ReadDataHeader(InBuff, pBuff);
+  simple->DataProceed(data_header, InBuff, pBuff, flag);
+}
+
+void
+Bit16SignedInteger::Collect(const string &name, bool found, CollectStruct &collect_struct) const
+{
+}
+
+void
+Bit16SignedInteger::Print() const
+{
+  //cout << "Bit16SignedInteger: \n";
+}
+
+void
+Bit16SignedInteger::Delete() const
+{
+  //cout << "Bit16SignedInteger: \n";
+}
+
+void
+Bit16SignedInteger::ReadArray(char* InBuff, int* pBuff, FlagStructure_t flag)
+{
+  Data_Header_t data_header;
+  SimpleElem* simple;
+  simple=new SimpleElem;
+  VCell.push_back(simple);
+  data_header=simple->ReadDataHeader(InBuff, pBuff);
+  simple->DataProceed(data_header, InBuff, pBuff, flag);
+}
+
+void
+Bit16UnsignedInteger::Collect(const string &name, bool found, CollectStruct &collect_struct) const
+{
+}
+
+void
+Bit16UnsignedInteger::Print() const
+{
+  //cout << "Bit16UnsignedInteger: \n";
+}
+
+void
+Bit16UnsignedInteger::Delete() const
+{
+  //cout << "Bit16UnsignedInteger: \n";
+}
+
+void
+Bit16UnsignedInteger::ReadArray(char* InBuff, int* pBuff, FlagStructure_t flag)
+{
+  Data_Header_t data_header;
+  SimpleElem* simple;
+  simple=new SimpleElem;
+  VCell.push_back(simple);
+  data_header=simple->ReadDataHeader(InBuff, pBuff);
+  simple->DataProceed(data_header, InBuff, pBuff, flag);
+}
+
+void
+Bit32SignedInteger::Collect(const string &name, bool found, CollectStruct &collect_struct) const
+{
+}
+
+void
+Bit32SignedInteger::Print() const
+{
+  //cout << "Bit32SignedInteger: \n";
+}
+
+void
+Bit32SignedInteger::Delete() const
+{
+  //cout << "Bit32SignedInteger: \n";
+}
+
+void
+Bit32SignedInteger::ReadArray(char* InBuff, int* pBuff, FlagStructure_t flag)
+{
+  Data_Header_t data_header;
+  SimpleElem* simple;
+  simple=new SimpleElem;
+  VCell.push_back(simple);
+  data_header=simple->ReadDataHeader(InBuff, pBuff);
+  simple->DataProceed(data_header, InBuff, pBuff, flag);
+}
+
+void
+Bit32UnsignedInteger::Collect(const string &name, bool found, CollectStruct &collect_struct) const
+{
+}
+
+void
+Bit32UnsignedInteger::Print() const
+{
+  //cout << "Bit32UnsignedInteger: \n";
+}
+
+void
+Bit32UnsignedInteger::Delete() const
+{
+  //cout << "Bit32UnsignedInteger: \n";
+}
+
+void
+Bit32UnsignedInteger::ReadArray(char* InBuff, int* pBuff, FlagStructure_t flag)
+{
+  Data_Header_t data_header;
+  SimpleElem* simple;
+  simple=new SimpleElem;
+  VCell.push_back(simple);
+  data_header=simple->ReadDataHeader(InBuff, pBuff);
+  simple->DataProceed(data_header, InBuff, pBuff, flag);
+}
+
+Data_Header_t
+MatlabFile::ReadDataHeader(ifstream &MatFile)
+{
+  Data_Header_t data_header;
+  MatFile.read(reinterpret_cast<char*>(&data_header.DataType),sizeof(data_header.DataType));
+  MatFile.read(reinterpret_cast<char*>(&data_header.S_Number_of_Bytes),sizeof(data_header.S_Number_of_Bytes));
+  if(data_header.S_Number_of_Bytes!=0)
+    data_header.Number_of_Bytes=data_header.S_Number_of_Bytes;
+  else
+    {
+      MatFile.read(reinterpret_cast<char*>(&data_header.Number_of_Bytes),sizeof(data_header.Number_of_Bytes));
+    }
+  if(data_header.Number_of_Bytes<8)
+    data_header.Number_of_Bytes = 8;
+  return(data_header);
+}
+
+
+
+
+
+void
+MatlabFile::MatFileRead(string filename)
+{
+  ifstream MatFile;
+  Data_Header_t data_header;
+  SimpleElem *simpl;
+  int pBuff;
+  FlagStructure_t flag;
+  //ArrayElem elem;
+  MatFile.open(filename.c_str(),std::ios::in | std::ios::binary);
+  if (!MatFile.is_open())
+    {
+      cerr << filename.c_str() << " Cannot be opened\n";
+      exit(EXIT_FAILURE);
+    }
+  // Read the Header of the Mat-File
+  MatFile.read(reinterpret_cast<char*>(&header),sizeof(header));
+  do
+    {
+      data_header=ReadDataHeader(MatFile);
+      char* InBuff;
+      InBuff = (char*)malloc(data_header.Number_of_Bytes+1);
+      MatFile.read(InBuff,data_header.Number_of_Bytes+1);
+      pBuff = 0;
+      simpl = new SimpleElem;
+      VSimpl.push_back(simpl);
+      flag.no_name=false;
+      flag.character=false;
+      simpl->DataProceed(data_header, InBuff, &pBuff, flag);
+      free(InBuff);
+    }
+  while(!MatFile.eof());
+  MatFile.close();
+}
+
+bool
+MatlabFile::Collect(const string &name, CollectStruct &collect_struct) const
+{
+  for(vector<PSimpleElem>::const_iterator it=VSimpl.begin(); it!=VSimpl.end(); it++)
+      (*it)->Collect(name, false, collect_struct);
+  return(!(collect_struct.variable_double_name.empty() and collect_struct.variable_string_name.empty()));
+}
+
+
+void
+MatlabFile::MatFilePrint()
+{
+  for(vector<PSimpleElem>::iterator it=VSimpl.begin(); it!=VSimpl.end(); it++)
+      (*it)->Print();
+}
+
+
+void
+MatlabFile::Delete()
+{
+  for(vector<PSimpleElem>::iterator it=VSimpl.begin(); it!=VSimpl.end(); it++)
+    {
+      (*it)->Delete();
+      delete *it;
+    }
+}
+
+/*
+int
+main(int argc, char** argv)
+{
+  CollectStruct collect_struct;
+  MatlabFile matlab_file;
+  matlab_file.MatFileRead("gimf_steady.mat");
+  //matlab_file.MatFileRead("essai.mat");
+  matlab_file.MatFilePrint();
+  bool tmp_b=false;
+  //string tmp_s("stored_values");
+  tmp_b=matlab_file.Collect("stored_values", collect_struct);
+  if(tmp_b)
+    {
+      int i=0;
+      for(map<string,vector<double> >::iterator it2=collect_struct.variable_double_name.begin();it2!=collect_struct.variable_double_name.end();it2++)
+        {
+          i++;
+          cout << i << " " << it2->first.c_str() << " : ";
+          for(vector<double>::iterator it=it2->second.begin();it!=it2->second.end();it++)
+            cout << *it;
+          cout << "\n";
+        }
+    }
+}
+
+
+*/
diff --git a/preprocessor/ModFile.cc b/preprocessor/ModFile.cc
index 554f80d76..c2081049d 100644
--- a/preprocessor/ModFile.cc
+++ b/preprocessor/ModFile.cc
@@ -131,7 +131,6 @@ ModFile::evalAllExpressions()
         {
           if(global_eval_context.find(make_pair(j, eEndogenous))==global_eval_context.end())
             {
-              //it2=mat_file.variable.find(symbol_table.getNameByID(eEndogenous, j));
               map<string,vector<double> >::iterator it2=collect_struct.variable_double_name.find(symbol_table.getNameByID(eEndogenous, j));
               if(it2!=collect_struct.variable_double_name.end())
                 {
diff --git a/preprocessor/ModelTree.cc b/preprocessor/ModelTree.cc
index 127d70578..9c972fba2 100644
--- a/preprocessor/ModelTree.cc
+++ b/preprocessor/ModelTree.cc
@@ -277,7 +277,8 @@ ModelTree::writeModelEquations(ostream &output, ExprNodeOutputType output_type)
 void
 ModelTree::computeTemporaryTermsOrdered(int order, Model_Block *ModelBlock)
 {
-  map<NodeID, int> reference_count, first_occurence;
+  map<NodeID, pair<int, int> > first_occurence;
+  map<NodeID, int> reference_count;
   int i, j, m, eq, var, lag;
   temporary_terms_type vect;
   ostringstream tmp_output;
@@ -293,7 +294,7 @@ ModelTree::computeTemporaryTermsOrdered(int order, Model_Block *ModelBlock)
       for (i = 0;i < ModelBlock->Block_List[j].Size;i++)
         {
           eq_node = equations[ModelBlock->Block_List[j].Equation[i]];
-          eq_node->computeTemporaryTerms(reference_count, temporary_terms, first_occurence, j, ModelBlock, map_idx);
+          eq_node->computeTemporaryTerms(reference_count, temporary_terms, first_occurence, j, ModelBlock, i, map_idx);
         }
       for (m=0;m<=ModelBlock->Block_List[j].Max_Lead+ModelBlock->Block_List[j].Max_Lag;m++)
         {
@@ -303,7 +304,7 @@ ModelTree::computeTemporaryTermsOrdered(int order, Model_Block *ModelBlock)
               eq=ModelBlock->Block_List[j].IM_lead_lag[m].Equ_Index[i];
               var=ModelBlock->Block_List[j].IM_lead_lag[m].Var_Index[i];
               it=first_derivatives.find(make_pair(eq,variable_table.getID(eEndogenous, var,lag)));
-              it->second->computeTemporaryTerms(reference_count, temporary_terms, first_occurence, j, ModelBlock, map_idx);
+              it->second->computeTemporaryTerms(reference_count, temporary_terms, first_occurence, j, ModelBlock, ModelBlock->Block_List[j].Size-1, map_idx);
             }
         }
       for (m=0;m<=ModelBlock->Block_List[j].Max_Lead+ModelBlock->Block_List[j].Max_Lag;m++)
@@ -314,7 +315,7 @@ ModelTree::computeTemporaryTermsOrdered(int order, Model_Block *ModelBlock)
               eq=ModelBlock->Block_List[j].IM_lead_lag[m].Equ_X_Index[i];
               var=ModelBlock->Block_List[j].IM_lead_lag[m].Exogenous_Index[i];
               it=first_derivatives.find(make_pair(eq,variable_table.getID(eExogenous, var,lag)));
-              it->second->computeTemporaryTerms(reference_count, temporary_terms, first_occurence, j, ModelBlock, map_idx);
+              it->second->computeTemporaryTerms(reference_count, temporary_terms, first_occurence, j, ModelBlock, ModelBlock->Block_List[j].Size-1, map_idx);
             }
         }
       //jacobian_max_exo_col=(variable_table.max_exo_lag+variable_table.max_exo_lead+1)*symbol_table.exo_nbr;
@@ -328,7 +329,7 @@ ModelTree::computeTemporaryTermsOrdered(int order, Model_Block *ModelBlock)
                   eq=ModelBlock->Block_List[j].IM_lead_lag[m].Equ_Index_other_endo[i];
                   var=ModelBlock->Block_List[j].IM_lead_lag[m].Var_Index_other_endo[i];
                   it=first_derivatives.find(make_pair(eq,variable_table.getID(eEndogenous, var,lag)));
-                  it->second->computeTemporaryTerms(reference_count, temporary_terms, first_occurence, j, ModelBlock, map_idx);
+                  it->second->computeTemporaryTerms(reference_count, temporary_terms, first_occurence, j, ModelBlock, ModelBlock->Block_List[j].Size-1, map_idx);
                 }
             }
         }
@@ -470,7 +471,6 @@ ModelTree::writeModelEquationsOrdered_M( Model_Block *ModelBlock, const string &
           }
 
         output << "  g2=0;g3=0;\n";
-        temporary_terms_type tt2;
         if(ModelBlock->Block_List[j].Temporary_InUse->size())
           {
             tmp_output.str("");
@@ -502,24 +502,24 @@ ModelTree::writeModelEquationsOrdered_M( Model_Block *ModelBlock, const string &
             sps = "  ";
           else
             sps="";
-        if (ModelBlock->Block_List[j].Temporary_terms->size())
-          output << "  " << sps << "% //Temporary variables" << endl;
-        i=0;
-        for (temporary_terms_type::const_iterator it = ModelBlock->Block_List[j].Temporary_terms->begin();
-             it != ModelBlock->Block_List[j].Temporary_terms->end(); it++)
-          {
-            output << "  " <<  sps;
-            (*it)->writeOutput(output, oMatlabDynamicModelSparse, temporary_terms);
-            output << " = ";
-            (*it)->writeOutput(output, oMatlabDynamicModelSparse, tt2);
-            // Insert current node into tt2
-            tt2.insert(*it);
-            output << ";" << endl;
-            i++;
-          }
         // The equations
         for (i = 0;i < ModelBlock->Block_List[j].Size;i++)
           {
+            temporary_terms_type tt2;
+            tt2.clear();
+            if (ModelBlock->Block_List[j].Temporary_Terms_in_Equation[i]->size())
+              output << "  " << sps << "% //Temporary variables" << endl;
+            for (temporary_terms_type::const_iterator it = ModelBlock->Block_List[j].Temporary_Terms_in_Equation[i]->begin();
+                 it != ModelBlock->Block_List[j].Temporary_Terms_in_Equation[i]->end(); it++)
+              {
+                output << "  " <<  sps;
+                (*it)->writeOutput(output, oMatlabDynamicModelSparse, temporary_terms);
+                output << " = ";
+                (*it)->writeOutput(output, oMatlabDynamicModelSparse, tt2);
+                // Insert current node into tt2
+                tt2.insert(*it);
+                output << ";" << endl;
+              }
             string sModel = symbol_table.getNameByID(eEndogenous, ModelBlock->Block_List[j].Variable[i]) ;
             eq_node = equations[ModelBlock->Block_List[j].Equation[i]];
             lhs = eq_node->arg1;
@@ -904,7 +904,6 @@ ModelTree::writeModelStaticEquationsOrdered_M(Model_Block *ModelBlock, const str
         << "  % ////////////////////////////////////////////////////////////////////////" << endl;
         //The Temporary terms
         //output << global_output.str();
-        temporary_terms_type tt2;
         if(ModelBlock->Block_List[j].Temporary_InUse->size())
           {
             tmp_output.str("");
@@ -946,24 +945,24 @@ ModelTree::writeModelStaticEquationsOrdered_M(Model_Block *ModelBlock, const str
             output << "  residual=zeros(" << ModelBlock->Block_List[j].Size << ",1);\n";
           }
         sps="";
-        if (ModelBlock->Block_List[j].Temporary_terms->size())
-          output << "  " << sps << "% //Temporary variables" << endl;
-        i=0;
-        for (temporary_terms_type::const_iterator it = ModelBlock->Block_List[j].Temporary_terms->begin();
-             it != ModelBlock->Block_List[j].Temporary_terms->end(); it++)
-          {
-            output << "  " <<  sps;
-            (*it)->writeOutput(output, oMatlabStaticModelSparse, temporary_terms);
-            output << " = ";
-            (*it)->writeOutput(output, oMatlabStaticModelSparse, tt2);
-            // Insert current node into tt2
-            tt2.insert(*it);
-            output << ";" << endl;
-            i++;
-          }
         // The equations
         for (i = 0;i < ModelBlock->Block_List[j].Size;i++)
           {
+            temporary_terms_type tt2;
+            tt2.clear();
+            if (ModelBlock->Block_List[j].Temporary_Terms_in_Equation[i]->size())
+              output << "  " << sps << "% //Temporary variables" << endl;
+            for (temporary_terms_type::const_iterator it = ModelBlock->Block_List[j].Temporary_Terms_in_Equation[i]->begin();
+              it != ModelBlock->Block_List[j].Temporary_Terms_in_Equation[i]->end(); it++)
+              {
+                output << "  " <<  sps;
+                (*it)->writeOutput(output, oMatlabStaticModelSparse, temporary_terms);
+                output << " = ";
+                (*it)->writeOutput(output, oMatlabStaticModelSparse, tt2);
+                // Insert current node into tt2
+                tt2.insert(*it);
+                output << ";" << endl;
+              }
             string sModel = symbol_table.getNameByID(eEndogenous, ModelBlock->Block_List[j].Variable[i]) ;
             output << sps << "  % equation " << ModelBlock->Block_List[j].Equation[i]+1 << " variable : "
             << sModel << " (" << ModelBlock->Block_List[j].Variable[i]+1 << ")" << endl;
@@ -1202,38 +1201,41 @@ ModelTree::writeModelEquationsCodeOrdered(const string file_name, const Model_Bl
               }
             else
               lhs_rhs_done=false;
-            //The Temporary terms
-            temporary_terms_type tt2;
-            i=0;
-            for (temporary_terms_type::const_iterator it = ModelBlock->Block_List[j].Temporary_terms->begin();
-                 it != ModelBlock->Block_List[j].Temporary_terms->end(); it++)
-              {
-                (*it)->compile(code_file,false, output_type, tt2, map_idx);
-                code_file.write(&FSTPT, sizeof(FSTPT));
-                map_idx_type::const_iterator ii=map_idx.find((*it)->idx);
-                v=(int)ii->second;
-                code_file.write(reinterpret_cast<char *>(&v), sizeof(v));
-                // Insert current node into tt2
-                tt2.insert(*it);
-#ifdef DEBUGC
-                cout << "FSTPT " << v << "\n";
-                code_file.write(&FOK, sizeof(FOK));
-                code_file.write(reinterpret_cast<char *>(&i), sizeof(i));
-#endif
-                i++;
-              }
-            for (temporary_terms_type::const_iterator it = ModelBlock->Block_List[j].Temporary_terms->begin();
-                 it != ModelBlock->Block_List[j].Temporary_terms->end(); it++)
-              {
-                map_idx_type::const_iterator ii=map_idx.find((*it)->idx);
-#ifdef DEBUGC
-                cout << "map_idx[" << (*it)->idx <<"]=" << ii->second << "\n";
-#endif
-              }
             // The equations
             for (i = 0;i < ModelBlock->Block_List[j].Size;i++)
               {
                 //ModelBlock->Block_List[j].Variable_Sorted[i] = variable_table.getID(eEndogenous, ModelBlock->Block_List[j].Variable[i], 0);
+                //The Temporary terms
+                temporary_terms_type tt2;
+#ifdef DEBUGC
+                k=0;
+#endif
+                for (temporary_terms_type::const_iterator it = ModelBlock->Block_List[j].Temporary_Terms_in_Equation[i]->begin();
+                     it != ModelBlock->Block_List[j].Temporary_Terms_in_Equation[i]->end(); it++)
+                  {
+                    (*it)->compile(code_file,false, output_type, tt2, map_idx);
+                    code_file.write(&FSTPT, sizeof(FSTPT));
+                    map_idx_type::const_iterator ii=map_idx.find((*it)->idx);
+                    v=(int)ii->second;
+                    code_file.write(reinterpret_cast<char *>(&v), sizeof(v));
+                    // Insert current node into tt2
+                    tt2.insert(*it);
+#ifdef DEBUGC
+                    cout << "FSTPT " << v << "\n";
+                    code_file.write(&FOK, sizeof(FOK));
+                    code_file.write(reinterpret_cast<char *>(&k), sizeof(k));
+                    ki++;
+#endif
+
+                  }
+#ifdef DEBUGC
+                for (temporary_terms_type::const_iterator it = ModelBlock->Block_List[j].Temporary_terms->begin();
+                     it != ModelBlock->Block_List[j].Temporary_terms->end(); it++)
+                  {
+                    map_idx_type::const_iterator ii=map_idx.find((*it)->idx);
+                    cout << "map_idx[" << (*it)->idx <<"]=" << ii->second << "\n";
+                  }
+#endif
                 if (!lhs_rhs_done)
                   {
                     eq_node = equations[ModelBlock->Block_List[j].Equation[i]];
diff --git a/preprocessor/include/ExprNode.hh b/preprocessor/include/ExprNode.hh
index 06d9228e0..9fb52b9df 100644
--- a/preprocessor/include/ExprNode.hh
+++ b/preprocessor/include/ExprNode.hh
@@ -148,9 +148,10 @@ public:
   virtual void collectTemporary_terms(const temporary_terms_type &temporary_terms, Model_Block *ModelBlock, int Curr_Block) const = 0;
   virtual void computeTemporaryTerms(map<NodeID, int> &reference_count,
                                      temporary_terms_type &temporary_terms,
-                                     map<NodeID, int> &first_occurence,
+                                     map<NodeID, pair<int, int> > &first_occurence,
                                      int Curr_block,
                                      Model_Block *ModelBlock,
+                                     int equation,
                                      map_idx_type &map_idx) const;
 
   class EvalException
@@ -204,6 +205,13 @@ public:
   virtual void writeOutput(ostream &output, ExprNodeOutputType output_type, const temporary_terms_type &temporary_terms = temporary_terms_type()) const;
   virtual void collectEndogenous(set<pair<int, int> > &result) const;
   virtual void collectExogenous(set<pair<int, int> > &result) const;
+  virtual void computeTemporaryTerms(map<NodeID, int> &reference_count,
+                                   temporary_terms_type &temporary_terms,
+                                   map<NodeID, pair<int, int> > &first_occurence,
+                                   int Curr_block,
+                                   Model_Block *ModelBlock,
+                                   int equation,
+                                   map_idx_type &map_idx) const;
   virtual void collectTemporary_terms(const temporary_terms_type &temporary_terms, Model_Block *ModelBlock, int Curr_Block) const;
   virtual double eval(const eval_context_type &eval_context) const throw (EvalException);
   virtual void compile(ofstream &CompileCode, bool lhs_rhs, ExprNodeOutputType output_type, const temporary_terms_type &temporary_terms, map_idx_type map_idx) const;
@@ -225,9 +233,10 @@ public:
   virtual void writeOutput(ostream &output, ExprNodeOutputType output_type, const temporary_terms_type &temporary_terms) const;
   virtual void computeTemporaryTerms(map<NodeID, int> &reference_count,
                                      temporary_terms_type &temporary_terms,
-                                     map<NodeID, int> &first_occurence,
+                                     map<NodeID, pair<int, int> > &first_occurence,
                                      int Curr_block,
                                      Model_Block *ModelBlock,
+                                     int equation,
                                      map_idx_type &map_idx) const;
   virtual void collectEndogenous(set<pair<int, int> > &result) const;
   virtual void collectExogenous(set<pair<int, int> > &result) const;
@@ -254,9 +263,10 @@ public:
   virtual void writeOutput(ostream &output, ExprNodeOutputType output_type, const temporary_terms_type &temporary_terms) const;
   virtual void computeTemporaryTerms(map<NodeID, int> &reference_count,
                                      temporary_terms_type &temporary_terms,
-                                     map<NodeID, int> &first_occurence,
+                                     map<NodeID, pair<int, int> > &first_occurence,
                                      int Curr_block,
                                      Model_Block *ModelBlock,
+                                     int equation,
                                      map_idx_type &map_idx) const;
   virtual void collectEndogenous(set<pair<int, int> > &result) const;
   virtual void collectExogenous(set<pair<int, int> > &result) const;
@@ -290,9 +300,10 @@ public:
   virtual void writeOutput(ostream &output, ExprNodeOutputType output_type, const temporary_terms_type &temporary_terms) const;
   virtual void computeTemporaryTerms(map<NodeID, int> &reference_count,
                                      temporary_terms_type &temporary_terms,
-                                     map<NodeID, int> &first_occurence,
+                                     map<NodeID, pair<int, int> > &first_occurence,
                                      int Curr_block,
                                      Model_Block *ModelBlock,
+                                     int equation,
                                      map_idx_type &map_idx) const;
   virtual void collectEndogenous(set<pair<int, int> > &result) const;
   virtual void collectExogenous(set<pair<int, int> > &result) const;
@@ -317,9 +328,10 @@ public:
   virtual void writeOutput(ostream &output, ExprNodeOutputType output_type, const temporary_terms_type &temporary_terms) const;
   virtual void computeTemporaryTerms(map<NodeID, int> &reference_count,
                                      temporary_terms_type &temporary_terms,
-                                     map<NodeID, int> &first_occurence,
+                                     map<NodeID, pair<int, int> > &first_occurence,
                                      int Curr_block,
                                      Model_Block *ModelBlock,
+                                     int equation,
                                      map_idx_type &map_idx) const;
   virtual void collectEndogenous(set<pair<int, int> > &result) const;
   virtual void collectExogenous(set<pair<int, int> > &result) const;
@@ -349,7 +361,8 @@ struct Block
   bool is_linear;
   int *Equation, *Own_Derivative;
   int *Variable, *Other_Endogenous, *Exogenous;
-  temporary_terms_type *Temporary_terms;
+  temporary_terms_type **Temporary_Terms_in_Equation;
+  //temporary_terms_type *Temporary_terms;
   temporary_terms_inuse_type *Temporary_InUse;
   IM_compact *IM_lead_lag;
   int Code_Start, Code_Length;
diff --git a/preprocessor/include/MatlabFile.hh b/preprocessor/include/MatlabFile.hh
new file mode 100644
index 000000000..807b17439
--- /dev/null
+++ b/preprocessor/include/MatlabFile.hh
@@ -0,0 +1,426 @@
+/*
+ * Copyright (C) 2006-2008 Dynare Team
+ *
+ * This file is part of Dynare.
+ *
+ * Dynare is free software: you can redistribute it and/or modify
+ * it under the terms of the GNU General Public License as published by
+ * the Free Software Foundation, either version 3 of the License, or
+ * (at your option) any later version.
+ *
+ * Dynare is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+ * GNU General Public License for more details.
+ *
+ * You should have received a copy of the GNU General Public License
+ * along with Dynare.  If not, see <http://www.gnu.org/licenses/>.
+ */
+/*
+Usefull documentation: Matlab 7 Mat-File Format
+-----------------------------------------------
+revision: October 2008 PDF only Rereleased for Version 7.7 (Release 2008b)
+available at: http://www.mathworks.com/access/helpdesk/help/pdf_doc/matlab/matfile_format.pdf
+*/
+
+#ifndef _MAT_FILE_HH
+#define _MAT_FILE_HH
+
+#include <cstdio>
+#include <cstdlib>
+#include <string>
+#include <fstream>
+#include <cstring>
+#include <iostream>
+#include <map>
+#include <vector>
+#include <math.h>
+//! zlib needed to uncompress the mat-file. It is available with GCC 4.3.2 but it needs a dll !!
+//! => to avoid compress MatFile, save is used with option '-v6' in save_params_and_strady_state.m
+//#include "zlib.h"
+using namespace std;
+
+
+enum Data_Type
+{
+  miINT8       = 1,  //8 bit, signed
+  miUINT8      = 2,  //8 bit, unsigned
+  miINT16      = 3,  //16-bit, signed
+  miUINT16     = 4,  //16-bit, unsigned
+  miINT32      = 5,  //32-bit, signed
+  miUINT32     = 6,  //32-bit, unsigned
+  miSINGLE     = 7,  //IEEEŽ 754 single format
+  miDOUBLE     = 9,  //IEEE 754 double format
+  miINT64      = 12, //64-bit, signed
+  miUINT64     = 13, //64-bit, unsigned
+  miMATRIX     = 14, //MATLAB array
+  miCOMPRESSED = 15, //Compressed Data
+  miUTF8       = 16, //Unicode UTF-8 Encoded Character Data
+  miUTF16      = 17, //Unicode UTF-16 Encoded Character Data
+  miUTF32      = 18  //Unicode UTF-32 Encoded Character Data
+};
+
+enum Array_Type
+{
+  Cell_array              = 1,
+  Structure_              = 2,
+  Object_                 = 3,
+  Character_array         = 4,
+  Sparse_array            = 5,
+  Double_precision_array  = 6,
+  Single_precision_array  = 7,
+  Signed_integer_8_bit    = 8,
+  Unsigned_integer_8_bit  = 9,
+  Signed_integer_16_bit   = 10,
+  Unsigned_integer_16_bit = 11,
+  Signed_integer_32_bit   = 12,
+  Unsigned_integer_32_bit = 13
+};
+
+
+enum Function_Returned_Type
+{
+  Numerical    =1,
+  AlphaNumeric =2,
+  Matrix       =3,
+  Compressed   =4,
+  Unknown      =5
+};
+
+class ArrayElem;
+class SimpleElem;
+
+typedef long long LongLongInt;
+typedef long int LongInt;
+typedef short int Int;
+//typedef char ShortInt;
+typedef unsigned long int uLongInt ;
+typedef unsigned short int uShortInt ;
+typedef short int ShortInt ;
+typedef class SimpleElem *PSimpleElem;
+
+//!Header of MatFile
+typedef struct Header
+{
+  char Theader[124];
+  short int Version;
+  char Edian_Indicator[2];
+}
+Header_t;
+
+typedef struct Data_Header
+{
+  ShortInt S_Number_of_Bytes;
+  ShortInt DataType;
+  uLongInt Number_of_Bytes;
+}
+Data_Header_t;
+
+typedef struct Array_Flag
+{
+  Data_Header_t tag;
+  unsigned char classe;
+  unsigned char flag;
+  char undef1[2];
+  uLongInt nzmax;
+}
+Array_Flag_t;
+
+
+typedef struct returned_ReadData
+{
+  SimpleElem* Simple;
+  Function_Returned_Type Type;
+} returned_ReadData_t;
+
+
+typedef struct FlagStructure
+{
+  bool no_name;
+  bool character;
+} FlagStructure_t;
+
+typedef struct CollectStruct
+{
+  /*vector<string> variable_name;
+  vector< vector<double> > variable_double;
+  vector< vector<string> > variable_string;*/
+  string tmp_name;
+  map<string,vector<string> > variable_string_name;
+  map<string,vector<double> > variable_double_name;
+}
+CollectStruct;
+
+typedef vector<int> Array_Dimensions_t;
+
+
+
+
+
+//! Base class for simple elements in Mat-File
+class SimpleElem
+{
+public:
+   bool verbose;
+   vector<double> VNumeric;
+   vector<string> Vstr;
+   ArrayElem *array_elem;
+   int Type;
+   SimpleElem();
+   virtual ~SimpleElem();
+   virtual double ReadNum(char* InBuff, int* pBuff) const{return(0);};
+   virtual string ReadAlph(char* InBuff, int* pBuff, int Size) const{return(NULL);};
+   virtual Data_Header_t ReadDataHeader(char* InBuff, int* pBuff) const;
+   virtual int size() const{cout << "oups\n";return(0);};
+   void Print() const;
+   void Delete() const;
+   void Collect(const string &name, bool found, CollectStruct &collect_struct) const;
+   returned_ReadData_t ReadData(char* InBuff, int* pBuff) const;
+   returned_ReadData_t Get_Data_Class(Data_Header_t data_header) const;
+   void DataProceed(Data_Header data_header, char* InBuff, int* pBuff, FlagStructure flag);
+};
+
+
+class UTF8 : public SimpleElem
+{
+public:
+  virtual int size() const {return(1);};
+  virtual double ReadNum(char* InBuff, int* pBuff) const {return(NULL);};
+  virtual string ReadAlph(char* InBuff, int* pBuff, int Size) const;
+};
+
+class UTF16 : public SimpleElem
+{
+public:
+  virtual int size() const {return(2);};
+  virtual double ReadNum(char* InBuff, int* pBuff) const {return(NULL);};
+  virtual string ReadAlph(char* InBuff, int* pBuff, int Size) const;
+};
+
+class UTF32 : public SimpleElem
+{
+public:
+  virtual int size() const {return(4);};
+  virtual double ReadNum(char* InBuff, int* pBuff) const {return(NULL);};
+  virtual string ReadAlph(char* InBuff, int* pBuff, int Size) const;
+};
+
+class INT8 : public SimpleElem
+{
+public:
+   virtual int size() const {return(1);};
+   virtual double ReadNum(char* InBuff, int* pBuff) const;
+   virtual string ReadAlph(char* InBuff, int* pBuff, int Size) const {return(NULL);};
+};
+
+class INT16 : public SimpleElem
+{
+public:
+   virtual int size() const {return(2);};
+   virtual double ReadNum(char* InBuff, int* pBuff) const;
+   virtual string ReadAlph(char* InBuff, int* pBuff, int Size) const {return(NULL);};
+};
+
+
+class INT32 : public SimpleElem
+{
+public:
+   virtual int size() const {return(4);};
+   virtual double ReadNum(char* InBuff, int* pBuff) const;
+   virtual string ReadAlph(char* InBuff, int* pBuff, int Size) const {return(NULL);};
+};
+
+
+class INT64 : public SimpleElem
+{
+public:
+   virtual int size() const {return(8);};
+   virtual double ReadNum(char* InBuff, int* pBuff) const;
+   virtual string ReadAlph(char* InBuff, int* pBuff, int Size) const {return(NULL);};
+};
+
+
+class Single : public SimpleElem
+{
+public:
+   virtual int size() const {return(4);};
+   virtual double ReadNum(char* InBuff, int* pBuff) const;
+   virtual string ReadAlph(char* InBuff, int* pBuff, int Size) const {return(NULL);};
+};
+
+
+class Double : public SimpleElem
+{
+public:
+   virtual int size() const {return(8);};
+   virtual double ReadNum(char* InBuff, int* pBuff) const;
+   virtual string ReadAlph(char* InBuff, int* pBuff, int Size) const {return(NULL);};
+};
+
+//! Base class for Array Element in Mat-File
+class ArrayElem : public SimpleElem
+{
+private:
+
+protected:
+
+public:
+  int Cell_number, Structure_number, Matrix_Elem_number;
+  Array_Type Type;
+  vector<PSimpleElem> VCell;
+  vector<string> Structure_Elem_name;
+  bool array_complex, array_global, array_logical;
+  string variable_name;
+  int array_nzmax;
+  int number_of_dimensions;
+  vector<int> dimension;
+  vector<double> Double_value;
+  vector<string> String_value;
+
+  Array_Type type;
+  ArrayElem();
+  virtual LongInt ReadINT32(char* InBuff, int* pBuff) const;
+  virtual Array_Flag_t ReadArrayFlag(char* InBuff, int* pBuff) /*const*/;
+  virtual void ReadArrayDimension(char* InBuff, int* pBuff) /*const*/;
+  virtual void ReadArrayName(char* InBuff, int* pBuff) /*const*/;
+  virtual void ReadStructureNames(char* InBuff, int* pBuff);
+  virtual ArrayElem* ReadArray_class(char* InBuff, int* pBuff) const;
+  virtual void ReadArray(char* InBuff, int* pBuff, FlagStructure_t flag) /*const*/ {cout << "oups..\n";};
+  virtual void Collect(const string &name, bool found, CollectStruct &collect_struct) const;
+  virtual void Print() const;
+  virtual void Delete() const;
+  virtual ~ArrayElem();
+};
+
+class CellArray : public ArrayElem
+{
+  public:
+  virtual void Collect(const string &name, bool found, CollectStruct &collect_struct) const;
+  virtual void Print() const;
+  virtual void Delete() const;
+  virtual void ReadArray(char* InBuff, int* pBuff, FlagStructure_t flag);
+};
+
+class Structure : public ArrayElem
+{
+  public:
+  virtual void Collect(const string &name, bool found, CollectStruct &collect_struct) const;
+  virtual void Print() const;
+  virtual void Delete() const;
+  virtual void ReadArray(char* InBuff, int* pBuff, FlagStructure_t flag);
+};
+
+class Object : public ArrayElem
+{
+  public:
+  virtual void Collect(const string &name, bool found, CollectStruct &collect_struct) const;
+  virtual void Print() const;
+  virtual void Delete() const;
+  virtual void ReadArray(char* InBuff, int* pBuff, FlagStructure_t flag);
+};
+
+class CharacterArray : public ArrayElem
+{
+  public:
+  virtual void Collect(const string &name, bool found, CollectStruct &collect_struct) const;
+  virtual void Print() const;
+  virtual void Delete() const;
+  virtual void ReadArray(char* InBuff, int* pBuff, FlagStructure_t flag);
+};
+
+class SparseArray : public ArrayElem
+{
+  public:
+  virtual void Collect(const string &name, bool found, CollectStruct &collect_struct) const;
+  virtual void Print() const;
+  virtual void Delete() const;
+  virtual void ReadArray(char* InBuff, int* pBuff, FlagStructure_t flag);
+};
+
+class DoublePrecisionArray : public ArrayElem
+{
+  public:
+  virtual void Collect(const string &name, bool found, CollectStruct &collect_struct) const;
+  virtual void Print() const;
+  virtual void Delete() const;
+  virtual void ReadArray(char* InBuff, int* pBuff, FlagStructure_t flag);
+};
+
+class SinglePrecisionArray : public ArrayElem
+{
+  public:
+  virtual void Collect(const string &name, bool found, CollectStruct &collect_struct) const;
+  virtual void Print() const;
+  virtual void Delete() const;
+  virtual void ReadArray(char* InBuff, int* pBuff, FlagStructure_t flag);
+};
+
+class Bit8SignedInteger : public ArrayElem
+{
+  public:
+  virtual void Collect(const string &name, bool found, CollectStruct &collect_struct) const;
+  virtual void Print() const;
+  virtual void Delete() const;
+  virtual void ReadArray(char* InBuff, int* pBuff, FlagStructure_t flag);
+};
+
+class Bit8UnsignedInteger : public ArrayElem
+{
+  public:
+  virtual void Collect(const string &name, bool found, CollectStruct &collect_struct) const;
+  virtual void Print() const;
+  virtual void Delete() const;
+  virtual void ReadArray(char* InBuff, int* pBuff, FlagStructure_t flag);
+};
+
+class Bit16SignedInteger : public ArrayElem
+{
+  public:
+  virtual void Collect(const string &name, bool found, CollectStruct &collect_struct) const;
+  virtual void Print() const;
+  virtual void Delete() const;
+  virtual void ReadArray(char* InBuff, int* pBuff, FlagStructure_t flag);
+};
+
+class Bit16UnsignedInteger : public ArrayElem
+{
+  public:
+  virtual void Collect(const string &name, bool found, CollectStruct &collect_struct) const;
+  virtual void Print() const;
+  virtual void Delete() const;
+  virtual void ReadArray(char* InBuff, int* pBuff, FlagStructure_t flag);
+};
+
+class Bit32SignedInteger : public ArrayElem
+{
+  public:
+  virtual void Collect(const string &name, bool found, CollectStruct &collect_struct) const;
+  virtual void Print() const;
+  virtual void Delete() const;
+  virtual void ReadArray(char* InBuff, int* pBuff, FlagStructure_t flag);
+};
+
+class Bit32UnsignedInteger : public ArrayElem
+{
+  public:
+  virtual void Collect(const string &name, bool found, CollectStruct &collect_struct) const;
+  virtual void Print() const;
+  virtual void Delete() const;
+  virtual void ReadArray(char* InBuff, int* pBuff, FlagStructure_t flag);
+};
+
+
+class MatlabFile
+{
+  public:
+  Header_t header;
+  vector<PSimpleElem> VSimpl;
+  MatlabFile();
+  ~MatlabFile();
+  void MatFileRead(string filename);
+  void MatFilePrint();
+  void Delete();
+  bool Collect(const string &name, CollectStruct &collect_struct) const;
+  Data_Header_t ReadDataHeader(ifstream &MatFile);
+};
+#endif