diff --git a/Common/Context.h b/Common/Context.h
index 720f5b3ccfe4994e050cabe7375a27d3df304a4a..91738ca637b1ff884355211103352dc482db7e76 100644
--- a/Common/Context.h
+++ b/Common/Context.h
@@ -20,7 +20,7 @@ struct contextMeshOptions {
   int pointsNum, linesNum, surfacesNum, volumesNum;
   int optimize, optimizeNetgen, refineSteps, qualityType, labelType, remove4triangles;
   double normals, tangents, explode, angleSmoothNormals, allowSwapEdgeAngle;
-  double mshFileVersion, labelFrequency, pointSize, lineWidth;
+  double mshFileVersion, mshFilePartitioned, labelFrequency, pointSize, lineWidth;
   double qualityInf, qualitySup, radiusInf, radiusSup;
   double scalingFactor, lcFactor, randFactor, lcIntegrationPrecision;
   double lcMin, lcMax, toleranceEdgeLength;
diff --git a/Common/CreateFile.cpp b/Common/CreateFile.cpp
index c41a265c848fa0fd7d8ca5f2825af209385a48b9..a58a7998916469a52001669b1c925ac488b53d91 100644
--- a/Common/CreateFile.cpp
+++ b/Common/CreateFile.cpp
@@ -171,10 +171,18 @@ void CreateOutputFile(std::string fileName, int format)
     break;
 
   case FORMAT_MSH:
-    GModel::current()->writeMSH
-      (fileName, CTX::instance()->mesh.mshFileVersion,
-       CTX::instance()->mesh.binary, CTX::instance()->mesh.saveAll,
-       CTX::instance()->mesh.saveParametric, CTX::instance()->mesh.scalingFactor);
+    if(GModel::current()->getMeshPartitions().size() && 
+       CTX::instance()->mesh.mshFilePartitioned){
+      GModel::current()->writePartitionedMSH
+        (fileName, CTX::instance()->mesh.binary, CTX::instance()->mesh.saveAll,
+         CTX::instance()->mesh.saveParametric, CTX::instance()->mesh.scalingFactor);
+    }
+    else{
+      GModel::current()->writeMSH
+        (fileName, CTX::instance()->mesh.mshFileVersion,
+         CTX::instance()->mesh.binary, CTX::instance()->mesh.saveAll,
+         CTX::instance()->mesh.saveParametric, CTX::instance()->mesh.scalingFactor);
+    }
     if (CTX::instance()->mesh.saveDistance){
       GModel::current()->writeDistanceMSH
 	(fileName, CTX::instance()->mesh.mshFileVersion,
diff --git a/Common/DefaultOptions.h b/Common/DefaultOptions.h
index fca933aaf72b6ab082fc8680f25fdf9860688d40..55e38ad1eb7f4317e1784d296a0b059ecb0601ef 100644
--- a/Common/DefaultOptions.h
+++ b/Common/DefaultOptions.h
@@ -1102,8 +1102,10 @@ StringXNumber MeshOptions_Number[] = {
     "Minimum number of points used to mesh a circle" },
   { F|O, "MinimumCurvePoints" , opt_mesh_min_curv_points, 3. ,
     "Minimum number of points used to mesh a (non-straight) curve" },
-  { F|O, "MshFileVersion" , opt_mesh_msh_file_version , 2.1 , 
+  { F|O, "MshFileVersion" , opt_mesh_msh_file_version , 2.2 , 
     "Version of the MSH file format to use" },
+  { F|O, "MshFilePartitioned" , opt_mesh_msh_file_partitioned , 0. , 
+    "Split MSH file by mesh partition" },
 
   { F, "NbHexahedra" , opt_mesh_nb_hexahedra , 0. , 
     "Number of hexahedra in the current mesh (read-only)" },
diff --git a/Common/Gmsh.cpp b/Common/Gmsh.cpp
index 25d5fe73b46005fc0ef35ba886e66a4a987abd1e..5cfb1f00d25735fbf11f5d01345c64260393f248 100644
--- a/Common/Gmsh.cpp
+++ b/Common/Gmsh.cpp
@@ -170,8 +170,10 @@ int GmshBatch()
       RefineMesh(GModel::current(), CTX::instance()->mesh.secondOrderLinear);
 #if defined(HAVE_CHACO) || defined(HAVE_METIS)
     if(CTX::instance()->batchAfterMesh == 1){
-      if (CTX::instance()->partitionOptions.num_partitions > 1)PartitionMesh(GModel::current(), CTX::instance()->partitionOptions);
-      if (CTX::instance()->partitionOptions.renumber)RenumberMesh(GModel::current(), CTX::instance()->partitionOptions);
+      if (CTX::instance()->partitionOptions.num_partitions > 1)
+        PartitionMesh(GModel::current(), CTX::instance()->partitionOptions);
+      if (CTX::instance()->partitionOptions.renumber)
+        RenumberMesh(GModel::current(), CTX::instance()->partitionOptions);
     }
 #endif
 #endif
diff --git a/Common/Options.cpp b/Common/Options.cpp
index 4a259a0810f5353a959a9cb5e9d86781751e3f3d..e707c1cce4df7ff67b381f08c0a6e6b46950efe3 100644
--- a/Common/Options.cpp
+++ b/Common/Options.cpp
@@ -5364,6 +5364,13 @@ double opt_mesh_msh_file_version(OPT_ARGS_NUM)
   return CTX::instance()->mesh.mshFileVersion;
 }
 
+double opt_mesh_msh_file_partitioned(OPT_ARGS_NUM)
+{
+  if(action & GMSH_SET)
+    CTX::instance()->mesh.mshFilePartitioned = val;
+  return CTX::instance()->mesh.mshFilePartitioned;
+}
+
 double opt_mesh_binary(OPT_ARGS_NUM)
 {
   if(action & GMSH_SET)
diff --git a/Common/Options.h b/Common/Options.h
index 31b03f7ae7622e22e0c55509ed7f499870c8cd45..7fd5c9bf3c66a52bb28e5acf651ac0404fe2811c 100644
--- a/Common/Options.h
+++ b/Common/Options.h
@@ -495,6 +495,7 @@ double opt_mesh_light_lines(OPT_ARGS_NUM);
 double opt_mesh_light_two_side(OPT_ARGS_NUM);
 double opt_mesh_format(OPT_ARGS_NUM);
 double opt_mesh_msh_file_version(OPT_ARGS_NUM);
+double opt_mesh_msh_file_partitioned(OPT_ARGS_NUM);
 double opt_mesh_binary(OPT_ARGS_NUM);
 double opt_mesh_bdf_field_format(OPT_ARGS_NUM);
 double opt_mesh_nb_smoothing(OPT_ARGS_NUM);
diff --git a/Fltk/fileDialogs.cpp b/Fltk/fileDialogs.cpp
index 266916a6261a9a9565143e792c12307070035e55..77f4ebdff51cc7184c394da2f615b68ebc7825d3 100644
--- a/Fltk/fileDialogs.cpp
+++ b/Fltk/fileDialogs.cpp
@@ -721,9 +721,7 @@ int mshFileDialog(const char *name)
 {
   struct _mshFileDialog{
     Fl_Window *window;
-    Fl_Check_Button *b;
-    Fl_Check_Button *p;
-    Fl_Check_Button *d;
+    Fl_Check_Button *b[3];
     Fl_Choice *c;
     Fl_Button *ok, *cancel;
   };
@@ -740,19 +738,23 @@ int mshFileDialog(const char *name)
 
   if(!dialog){
     dialog = new _mshFileDialog;
-    int h = 3 * WB + 4 * BH, w = 2 * BBB + 3 * WB, y = WB;
+    int h = 3 * WB + 5 * BH, w = 2 * BBB + 3 * WB, y = WB;
     dialog->window = new Fl_Double_Window(w, h, "MSH Options");
     dialog->window->box(GMSH_WINDOW_BOX);
     dialog->window->set_modal();
     dialog->c = new Fl_Choice(WB, y, BBB + BBB / 2, BH, "Format"); y += BH;
     dialog->c->menu(formatmenu);
     dialog->c->align(FL_ALIGN_RIGHT);
-    dialog->b = new Fl_Check_Button
+    dialog->b[0] = new Fl_Check_Button
       (WB, y, 2 * BBB + WB, BH, "Save all (ignore physical groups)"); y += BH;
-    dialog->b->type(FL_TOGGLE_BUTTON);
-    dialog->p = new Fl_Check_Button
-      (WB, y, 2 * BBB + WB, BH, "Save Parametric Coordinates"); y += BH;
-    dialog->p->type(FL_TOGGLE_BUTTON);
+    dialog->b[0]->type(FL_TOGGLE_BUTTON);
+    dialog->b[1] = new Fl_Check_Button
+      (WB, y, 2 * BBB + WB, BH, "Save parametric coordinates"); y += BH;
+    dialog->b[1]->type(FL_TOGGLE_BUTTON);
+    dialog->b[2] = new Fl_Check_Button
+      (WB, y, 2 * BBB + WB, BH, "Save one file per partition"); y += BH;
+    dialog->b[2]->type(FL_TOGGLE_BUTTON);
+
     dialog->ok = new Fl_Return_Button(WB, y + WB, BBB, BH, "OK");
     dialog->cancel = new Fl_Button(2 * WB + BBB, y + WB, BBB, BH, "Cancel");
     dialog->window->end();
@@ -761,9 +763,9 @@ int mshFileDialog(const char *name)
   
   dialog->c->value((CTX::instance()->mesh.mshFileVersion == 1.0) ? 0 : 
                    CTX::instance()->mesh.binary ? 2 : 1);
-  dialog->b->value(CTX::instance()->mesh.saveAll ? 1 : 0);
-  dialog->p->value(CTX::instance()->mesh.saveParametric ? 1 : 0);
-  dialog->d->value(CTX::instance()->mesh.saveDistance ? 1 : 0);
+  dialog->b[0]->value(CTX::instance()->mesh.saveAll ? 1 : 0);
+  dialog->b[1]->value(CTX::instance()->mesh.saveParametric ? 1 : 0);
+  dialog->b[2]->value(CTX::instance()->mesh.mshFilePartitioned ? 1 : 0);
   dialog->window->show();
 
   while(dialog->window->shown()){
@@ -775,9 +777,9 @@ int mshFileDialog(const char *name)
         opt_mesh_msh_file_version(0, GMSH_SET | GMSH_GUI, 
                                   (dialog->c->value() == 0) ? 1.0 : 2.2);
         opt_mesh_binary(0, GMSH_SET | GMSH_GUI, (dialog->c->value() == 2) ? 1 : 0);
-        opt_mesh_save_all(0, GMSH_SET | GMSH_GUI, dialog->b->value() ? 1 : 0);
-        opt_mesh_save_parametric(0, GMSH_SET | GMSH_GUI, dialog->p->value() ? 1 : 0);
-	opt_mesh_save_distance(0, GMSH_SET | GMSH_GUI, dialog->d->value() ? 1 : 0);
+        opt_mesh_save_all(0, GMSH_SET | GMSH_GUI, dialog->b[0]->value() ? 1 : 0);
+        opt_mesh_save_parametric(0, GMSH_SET | GMSH_GUI, dialog->b[1]->value() ? 1 : 0);
+	opt_mesh_msh_file_partitioned(0, GMSH_SET | GMSH_GUI, dialog->b[2]->value() ? 1 : 0);
         CreateOutputFile(name, FORMAT_MSH);
         dialog->window->hide();
         return 1;
diff --git a/Geo/GModel.h b/Geo/GModel.h
index 070c3b22a84a0c5a56fbd91dc721158d78978bb8..43caee9ce4ade13cb003b85cac625085d5dd7d76 100644
--- a/Geo/GModel.h
+++ b/Geo/GModel.h
@@ -377,7 +377,12 @@ class GModel
   // Gmsh mesh file format
   int readMSH(const std::string &name);
   int writeMSH(const std::string &name, double version=1.0, bool binary=false,
-               bool saveAll=false, bool saveParametric=false, double scalingFactor=1.0);
+               bool saveAll=false, bool saveParametric=false, 
+               double scalingFactor=1.0, int elementStartNum=0, 
+               int saveSinglePartition=0);
+  int writePartitionedMSH(const std::string &baseName, bool binary=false,
+                          bool saveAll=false, bool saveParametric=false, 
+                          double scalingFactor=1.0);
   int writeDistanceMSH(const std::string &name, double version=1.0, bool binary=false,
 		       bool saveAll=false, bool saveParametric=false, 
                        double scalingFactor=1.0);
diff --git a/Geo/GModelIO_Mesh.cpp b/Geo/GModelIO_Mesh.cpp
index 0a9565d5eb8a41d74b3b98d3b0a61e76f825fa80..290ca40eb84f79f2b18a8f1b97c96507081f7144 100644
--- a/Geo/GModelIO_Mesh.cpp
+++ b/Geo/GModelIO_Mesh.cpp
@@ -8,6 +8,8 @@
 #include <string.h>
 #include <map>
 #include <string>
+#include <sstream>
+#include <iomanip>
 #include "GModel.h"
 #include "GmshDefines.h"
 #include "MPoint.h"
@@ -344,10 +346,12 @@ int GModel::readMSH(const std::string &name)
       Msg::ResetProgressMeter();
       if(!binary){
         for(int i = 0; i < numElements; i++) {
-          int num, type, physical = 0, elementary = 0, partition = 0, parent = 0, numVertices;
+          int num, type, physical = 0, elementary = 0, partition = 0, parent = 0;
+          int numVertices;
           std::vector<short> ghosts;
           if(version <= 1.0){
-            if(fscanf(fp, "%d %d %d %d %d", &num, &type, &physical, &elementary, &numVertices) != 5)
+            if(fscanf(fp, "%d %d %d %d %d", &num, &type, &physical, &elementary, 
+                      &numVertices) != 5)
               return 0;
             if(numVertices != MElement::getInfoMSH(type)) return 0;
           }
@@ -517,31 +521,64 @@ static void writeElementMSH(FILE *fp, GModel *model, T *ele, bool saveAll,
 
 template<class T>
 static void writeElementsMSH(FILE *fp, GModel *model, std::vector<T*> &ele, 
-                             bool saveAll, double version, bool binary, int &num, 
-                             int elementary, std::vector<int> &physicals)
+                             bool saveAll, int saveSinglePartition, double version,
+                             bool binary, int &num, int elementary,
+                             std::vector<int> &physicals,
+                             std::map<MElement*, int> *parentsNum=0)
 {
-  for(unsigned int i = 0; i < ele.size(); i++)
+  for(unsigned int i = 0; i < ele.size(); i++){
+    if(saveSinglePartition && ele[i]->getPartition() != saveSinglePartition)
+      continue;
+    int parentNum = 0;
+    if(parentsNum){
+      MElement *parent = ele[i]->getParent();
+      if(parent){
+        std::map<MElement*, int>::iterator it = parentsNum->find(parent);
+        if(it != parentsNum->end()) parentNum = it->second;
+      }
+    }
     writeElementMSH(fp, model, ele[i], saveAll, version, binary, num,
-                    elementary, physicals, 0);
+                    elementary, physicals, parentNum);
+  }
 }
 
-template<class T>
-static void writeElementsMSH(FILE *fp, GModel *model, std::vector<T*> &ele,
-                             bool saveAll, double version, bool binary, int &num,
-                             int elementary, std::vector<int> &physicals,
-                             std::map<MElement*, int> &parentsNum)
+static int getNumElementsMSH(GEntity *ge, bool saveAll, int saveSinglePartition)
 {
-  for(unsigned int i = 0; i < ele.size(); i++) {
-    int parentNum = (ele[i]->getParent() &&
-                     parentsNum.find(ele[i]->getParent()) != parentsNum.end()) ?
-      parentsNum.find(ele[i]->getParent())->second : 0;
-    writeElementMSH(fp, model, ele[i], saveAll, version, binary, num,
-                    elementary, physicals, parentNum);
-  }
+  int n = 0, p = saveAll ? 1 : ge->physicals.size();
+  if(!saveSinglePartition)
+    n = p * ge->getNumMeshElements();
+  else
+    for(unsigned int i = 0; i < ge->getNumMeshElements(); i++)
+      if(ge->getMeshElement(i)->getPartition() == saveSinglePartition)
+        n += p;
+  return n;
+}
+
+static int getNumElementsMSH(GModel *m, bool saveAll, int saveSinglePartition)
+{
+  int n = 0;
+  for(GModel::viter it = m->firstVertex(); it != m->lastVertex(); ++it)
+    n += getNumElementsMSH(*it, saveAll, saveSinglePartition);
+  for(GModel::eiter it = m->firstEdge(); it != m->lastEdge(); ++it)
+    n += getNumElementsMSH(*it, saveAll, saveSinglePartition);
+  for(GModel::fiter it = m->firstFace(); it != m->lastFace(); ++it){
+    n += getNumElementsMSH(*it, saveAll, saveSinglePartition);
+    for(unsigned int i = 0; i < (*it)->polygons.size(); i++)
+      if((*it)->polygons[i]->ownsParent())
+        n += (saveAll ? 1 : (*it)->physicals.size());
+  }
+  for(GModel::riter it = m->firstRegion(); it != m->lastRegion(); ++it){
+    n += getNumElementsMSH(*it, saveAll, saveSinglePartition);
+    for(unsigned int i = 0; i < (*it)->polyhedra.size(); i++)
+      if((*it)->polyhedra[i]->ownsParent())
+        n += (saveAll ? 1 : (*it)->physicals.size());
+  }
+  return n;
 }
 
 int GModel::writeMSH(const std::string &name, double version, bool binary,
-                     bool saveAll, bool saveParametric, double scalingFactor)
+                     bool saveAll, bool saveParametric, double scalingFactor,
+                     int elementStartNum, int saveSinglePartition)
 {
   FILE *fp = fopen(name.c_str(), binary ? "wb" : "w");
   if(!fp){
@@ -555,6 +592,9 @@ int GModel::writeMSH(const std::string &name, double version, bool binary,
   // get the number of vertices and index the vertices in a continuous
   // sequence
   int numVertices = indexMeshVertices(saveAll);
+
+  // FIXME if saveSinglePartition, re-tag some nodes with '0' index
+  // and recompute numVertices
   
   // binary format exists only in version 2
   if(version > 1 || binary) 
@@ -563,32 +603,24 @@ int GModel::writeMSH(const std::string &name, double version, bool binary,
     version = 1.0;
   
   // get the number of elements we need to save
-  int numElements = 0;
+  int numElements = getNumElementsMSH(this, saveAll, saveSinglePartition);
+
+  // get parent elements for polygons and polyhedra
   std::map<MElement*, GEntity*> parents[2];
-  for(viter it = firstVertex(); it != lastVertex(); ++it)
-    numElements += (saveAll ? 1 : (*it)->physicals.size()) * (*it)->points.size();
-  for(eiter it = firstEdge(); it != lastEdge(); ++it)
-    numElements += (saveAll ? 1 : (*it)->physicals.size()) * (*it)->lines.size();
-  for(fiter it = firstFace(); it != lastFace(); ++it){
-    numElements += (saveAll ? 1 : (*it)->physicals.size()) * 
-      ((*it)->triangles.size() + (*it)->quadrangles.size() + (*it)->polygons.size());
-    for(unsigned int i = 0; i < (*it)->polygons.size(); i++){
-      if((*it)->polygons[i]->ownsParent()) {
+  for(fiter it = firstFace(); it != lastFace(); ++it)
+    for(unsigned int i = 0; i < (*it)->polygons.size(); i++)
+      if((*it)->polygons[i]->ownsParent())
         parents[0][(*it)->polygons[i]->getParent()] = (*it);
-        numElements += (saveAll ? 1 : (*it)->physicals.size());
-      }
-    }
-  }
-  for(riter it = firstRegion(); it != lastRegion(); ++it){
-    numElements += (saveAll ? 1 : (*it)->physicals.size()) *
-      ((*it)->tetrahedra.size() + (*it)->hexahedra.size() + (*it)->prisms.size() +
-       (*it)->pyramids.size() + (*it)->polyhedra.size());
-    for(unsigned int i = 0; i < (*it)->polyhedra.size(); i++){
-      if((*it)->polyhedra[i]->ownsParent()) {
+  for(riter it = firstRegion(); it != lastRegion(); ++it)
+    for(unsigned int i = 0; i < (*it)->polyhedra.size(); i++)
+      if((*it)->polyhedra[i]->ownsParent())
         parents[1][(*it)->polyhedra[i]->getParent()] = (*it);
-        numElements += (saveAll ? 1 : (*it)->physicals.size());
-      }
-    }
+
+  // get ghost cells
+  std::multimap<MElement*, int> ghostCells;
+  if(saveSinglePartition && getGhostCells().size()){
+    // XXXX
+    numElements += ghostCells.size();    
   }
 
   if(version >= 2.0){
@@ -642,18 +674,18 @@ int GModel::writeMSH(const std::string &name, double version, bool binary,
   }
 
   fprintf(fp, "%d\n", numElements);
-  int num = 0;
+  int num = elementStartNum;
   std::map<MElement*, int> parentsNum;
 
   // points
   for(viter it = firstVertex(); it != lastVertex(); ++it)
-    writeElementsMSH(fp, this, (*it)->points, saveAll, version, binary, num,
-                     (*it)->tag(), (*it)->physicals);
+    writeElementsMSH(fp, this, (*it)->points, saveAll, saveSinglePartition,
+                     version, binary, num, (*it)->tag(), (*it)->physicals);
 
   // lines
   for(eiter it = firstEdge(); it != lastEdge(); ++it)
-    writeElementsMSH(fp, this, (*it)->lines, saveAll, version, binary, num,
-                     (*it)->tag(), (*it)->physicals);
+    writeElementsMSH(fp, this, (*it)->lines, saveAll, saveSinglePartition,
+                     version, binary, num, (*it)->tag(), (*it)->physicals);
   
   // triangles
   for(std::map<MElement*, GEntity*>::const_iterator it = parents[0].begin(); 
@@ -664,8 +696,8 @@ int GModel::writeMSH(const std::string &name, double version, bool binary,
       parentsNum[it->first] = num;
     }
   for(fiter it = firstFace(); it != lastFace(); ++it)
-    writeElementsMSH(fp, this, (*it)->triangles, saveAll, version, binary, num,
-                     (*it)->tag(), (*it)->physicals);
+    writeElementsMSH(fp, this, (*it)->triangles, saveAll, saveSinglePartition,
+                     version, binary, num, (*it)->tag(), (*it)->physicals);
 
   // quads
   for(std::map<MElement*, GEntity*>::const_iterator it = parents[0].begin(); 
@@ -676,8 +708,8 @@ int GModel::writeMSH(const std::string &name, double version, bool binary,
       parentsNum[it->first] = num;
     }
   for(fiter it = firstFace(); it != lastFace(); ++it)
-    writeElementsMSH(fp, this, (*it)->quadrangles, saveAll, version, binary, num,
-                     (*it)->tag(), (*it)->physicals);
+    writeElementsMSH(fp, this, (*it)->quadrangles, saveAll, saveSinglePartition,
+                     version, binary, num, (*it)->tag(), (*it)->physicals);
 
   // polygons
   for(std::map<MElement*, GEntity*>::const_iterator it = parents[0].begin();
@@ -688,8 +720,9 @@ int GModel::writeMSH(const std::string &name, double version, bool binary,
       parentsNum[it->first] = num;
     }
   for(fiter it = firstFace(); it != lastFace(); it++)
-    writeElementsMSH(fp, this, (*it)->polygons, saveAll, version, binary, num,
-                     (*it)->tag(), (*it)->physicals, parentsNum);
+    writeElementsMSH(fp, this, (*it)->polygons, saveAll, saveSinglePartition,
+                     version, binary, num, (*it)->tag(), (*it)->physicals, 
+                     &parentsNum);
 
   // tets
   for(std::map<MElement*, GEntity*>::const_iterator it = parents[1].begin();
@@ -700,8 +733,8 @@ int GModel::writeMSH(const std::string &name, double version, bool binary,
       parentsNum[it->first] = num;
     }
   for(riter it = firstRegion(); it != lastRegion(); ++it)
-    writeElementsMSH(fp, this, (*it)->tetrahedra, saveAll, version, binary, num,
-                     (*it)->tag(), (*it)->physicals);
+    writeElementsMSH(fp, this, (*it)->tetrahedra, saveAll, saveSinglePartition,
+                     version, binary, num, (*it)->tag(), (*it)->physicals);
 
   // hexas
   for(std::map<MElement*, GEntity*>::const_iterator it = parents[1].begin();
@@ -712,8 +745,8 @@ int GModel::writeMSH(const std::string &name, double version, bool binary,
       parentsNum[it->first] = num;
     }
   for(riter it = firstRegion(); it != lastRegion(); ++it)
-    writeElementsMSH(fp, this, (*it)->hexahedra, saveAll, version, binary, num,
-                     (*it)->tag(), (*it)->physicals);
+    writeElementsMSH(fp, this, (*it)->hexahedra, saveAll, saveSinglePartition,
+                     version, binary, num, (*it)->tag(), (*it)->physicals);
 
   // prisms
   for(std::map<MElement*, GEntity*>::const_iterator it = parents[1].begin();
@@ -724,8 +757,8 @@ int GModel::writeMSH(const std::string &name, double version, bool binary,
       parentsNum[it->first] = num;
     }
   for(riter it = firstRegion(); it != lastRegion(); ++it)
-    writeElementsMSH(fp, this, (*it)->prisms, saveAll, version, binary, num,
-                     (*it)->tag(), (*it)->physicals);
+    writeElementsMSH(fp, this, (*it)->prisms, saveAll, saveSinglePartition,
+                     version, binary, num, (*it)->tag(), (*it)->physicals);
   
   // pyramids
   for(std::map<MElement*, GEntity*>::const_iterator it = parents[1].begin();
@@ -736,8 +769,8 @@ int GModel::writeMSH(const std::string &name, double version, bool binary,
       parentsNum[it->first] = num;
     }
   for(riter it = firstRegion(); it != lastRegion(); ++it)
-    writeElementsMSH(fp, this, (*it)->pyramids, saveAll, version, binary, num,
-                     (*it)->tag(), (*it)->physicals);
+    writeElementsMSH(fp, this, (*it)->pyramids, saveAll, saveSinglePartition,
+                     version, binary, num, (*it)->tag(), (*it)->physicals);
 
   // polyhedra
   for(std::map<MElement*, GEntity*>::const_iterator it = parents[1].begin();
@@ -748,8 +781,9 @@ int GModel::writeMSH(const std::string &name, double version, bool binary,
       parentsNum[it->first] = num;
     }
   for(riter it = firstRegion(); it != lastRegion(); ++it)
-    writeElementsMSH(fp, this, (*it)->polyhedra, saveAll, version, binary, num,
-                     (*it)->tag(), (*it)->physicals, parentsNum);
+    writeElementsMSH(fp, this, (*it)->polyhedra, saveAll, saveSinglePartition,
+                     version, binary, num, (*it)->tag(), (*it)->physicals,
+                     &parentsNum);
   
   if(binary) fprintf(fp, "\n");
 
@@ -778,6 +812,26 @@ int GModel::writeMSH(const std::string &name, double version, bool binary,
   return 1;
 }
 
+int GModel::writePartitionedMSH(const std::string &baseName, bool binary,
+                                bool saveAll, bool saveParametric,
+                                double scalingFactor)
+{
+  int index = 0;
+  for(std::set<int>::iterator it = meshPartitions.begin(); 
+      it != meshPartitions.end(); it++){
+    int partition = *it;
+    
+    std::ostringstream sstream;
+    sstream << baseName << "_" << std::setw(3) << std::setfill('0') << partition;
+
+    int startNum = index ? getNumElementsMSH(this, saveAll, partition) : 0;
+    writeMSH(sstream.str(), 2.2, binary, saveAll, saveParametric, 
+             scalingFactor, startNum, partition);
+    index++;
+  }
+  return 1;
+}
+
 int GModel::writeDistanceMSH(const std::string &name, double version, bool binary,
 			     bool saveAll, bool saveParametric, double scalingFactor)
 {
diff --git a/Geo/GRegion.cpp b/Geo/GRegion.cpp
index 72844a4012c4627c1a92ed442b62bc0b13fcfb99..97a2c00014f99f62d1743ebd1027547fecad6085 100644
--- a/Geo/GRegion.cpp
+++ b/Geo/GRegion.cpp
@@ -50,7 +50,8 @@ void GRegion::deleteMesh()
 
 unsigned int GRegion::getNumMeshElements()
 { 
-  return tetrahedra.size() + hexahedra.size() + prisms.size() + pyramids.size() + polyhedra.size();
+  return tetrahedra.size() + hexahedra.size() + prisms.size() + pyramids.size() +
+    polyhedra.size();
 }
 
 void GRegion::getNumMeshElements(unsigned *const c) const
diff --git a/Mesh/meshPartitionOptions.h b/Mesh/meshPartitionOptions.h
index 05406a7715b44f5f4d547cdd2a0e7558da7c29de..78d4bb796ef75c2b29edc5531c0deed92fe71c8b 100644
--- a/Mesh/meshPartitionOptions.h
+++ b/Mesh/meshPartitionOptions.h
@@ -75,9 +75,9 @@ struct meshPartitionOptions
 //--Constructor
 
   meshPartitionOptions()
-     :
-     goal(0)
-  { }
+  {
+    setDefaults();
+  }
 
 //--Default values
 
diff --git a/Solver/TESTCASES/Advection3D.lua b/Solver/TESTCASES/Advection3D.lua
index c00bb6d9d609db046abb0469be95c50786e08df1..3c35593d581ecfc7520db04e39038858b4198d56 100644
--- a/Solver/TESTCASES/Advection3D.lua
+++ b/Solver/TESTCASES/Advection3D.lua
@@ -1,21 +1,21 @@
 
-model = GModel  ()
-model:load ('box.geo')
-model:load ('box.msh')
-dg = dgSystemOfEquations (model)
+model = GModel()
+model:load('box.geo')
+model:load('box.msh')
+dg = dgSystemOfEquations(model)
 dg:setOrder(1)
 
 -- initial condition
 function initial_condition( xyz , f )
-  for i=0,xyz:size1()-1 do
-    x = xyz:get(i,0)
-    y = xyz:get(i,1)
-    z = xyz:get(i,2)
-    --f:set (i, 0, math.exp(-100*((x-0.5)^2)))
-    --f:set (i, 0, x )
-    f:set (i, 0, math.sqrt((x-0.3)*(x-0.3)+(y-0.3)*(y-0.3)+z*z)-0.2)	
-    --f:set(i, 0, math.exp((x-.3)*(x-.3)+(y-.3)*(y-.3)+z*z)-.2))
-  end
+   for i=0,xyz:size1()-1 do
+      x = xyz:get(i,0)
+      y = xyz:get(i,1)
+      z = xyz:get(i,2)
+      --f:set (i, 0, math.exp(-100*((x-0.5)^2)))
+      --f:set (i, 0, x )
+      f:set (i, 0, math.sqrt((x-0.3)*(x-0.3)+(y-0.3)*(y-0.3)+z*z)-0.2)	
+      --f:set(i, 0, math.exp((x-.3)*(x-.3)+(y-.3)*(y-.3)+z*z)-.2))
+   end
 end
 
 -- conservation law