diff --git a/Geo/CellComplex.cpp b/Geo/CellComplex.cpp
index 64224caf67db307e58e1c1e69ac5554450901fef..ecf672d2f9e7737156915d9983cf4e5e6290954c 100644
--- a/Geo/CellComplex.cpp
+++ b/Geo/CellComplex.cpp
@@ -421,14 +421,16 @@ int CellComplex::coreduceComplex(bool omitLowdim){
   
   int count = 0;
   
-  removeSubdomain();
+  CellComplex::removeSubdomain();
   
   
   for(int dim = 0; dim < 4; dim++){
     citer cit = firstCell(dim);
-    if(cit != lastCell(dim)){
+    while(cit != lastCell(dim)){
       Cell* cell = *cit;
       count = count + coreduction(cell);
+      if(count != 0) break;
+      cit++;
     }
   } 
   
diff --git a/Geo/Homology.cpp b/Geo/Homology.cpp
index d39404115bbe9df0ddde19d86cad7feb37e22c12..7174b65e0f2540f3bc5fa5e847c794931a45e33f 100644
--- a/Geo/Homology.cpp
+++ b/Geo/Homology.cpp
@@ -86,7 +86,9 @@ void Homology::findGenerators(std::string fileName){
   chains->computeHomology();
   Msg::Info("Homology Computation complete.");
   
+  int HRank[4];
   for(int j = 0; j < 4; j++){
+    HRank[j] = 0;
     for(int i = 1; i <= chains->getBasisSize(j); i++){
       
       std::string generator;
@@ -97,16 +99,17 @@ void Homology::findGenerators(std::string fileName){
       std::string name = dimension + "D Generator " + generator;
       Chain* chain = new Chain(_cellComplex->getCells(j), chains->getCoeffVector(j,i), _cellComplex, name);
       chain->writeChainMSH(fileName);
+      if(chain->getSize() != 0) HRank[j] = HRank[j] + 1;
       delete chain;
     }
   }
   
   Msg::Info("Ranks of homology groups for primal cell complex:");
-  Msg::Info("H0 = %d", chains->getBasisSize(0));
-  Msg::Info("H1 = %d", chains->getBasisSize(1));
-  Msg::Info("H2 = %d", chains->getBasisSize(2));
-  Msg::Info("H3 = %d", chains->getBasisSize(3));
-  if(omitted != 0) Msg::Info("%d 0D generators are not shown completely.", omitted);
+  Msg::Info("H0 = %d", HRank[0]);
+  Msg::Info("H1 = %d", HRank[1]);
+  Msg::Info("H2 = %d", HRank[2]);
+  Msg::Info("H3 = %d", HRank[3]);
+  if(omitted != 0) Msg::Info("Computation of %dD generators was omitted.", _cellComplex->getDim());
   
   
   Msg::Info("Wrote results to %s.", fileName.c_str());
@@ -118,6 +121,8 @@ void Homology::findGenerators(std::string fileName){
 
 void Homology::findThickCuts(std::string fileName){
   
+  _cellComplex->removeSubdomain();
+  
   Msg::Info("Reducing Cell Complex...");
   int omitted = _cellComplex->coreduceComplex(true);
   _cellComplex->cocombine(0);
@@ -135,28 +140,33 @@ void Homology::findThickCuts(std::string fileName){
   chains->computeHomology(true);
   Msg::Info("Homology Computation complete.");
   
+  int dim = _cellComplex->getDim();
+  
+  int HRank[4];
   for(int j = 3; j > -1; j--){
+    HRank[j] = 0;
     for(int i = 1; i <= chains->getBasisSize(j); i++){
       
       std::string generator;
       std::string dimension;
       convert(i, generator);
-      convert(3-j, dimension);
+      convert(dim-j, dimension);
       
       std::string name = dimension + "D Thick cut " + generator;
       Chain* chain = new Chain(_cellComplex->getCells(j), chains->getCoeffVector(j,i), _cellComplex, name);
       chain->writeChainMSH(fileName);
+      if(chain->getSize() != 0) HRank[j] = HRank[j] + 1;
       delete chain;
             
     }
   }
   
   Msg::Info("Ranks of homology groups for dual cell complex:");
-  Msg::Info("H0 = %d", chains->getBasisSize(3));
-  Msg::Info("H1 = %d", chains->getBasisSize(2));
-  Msg::Info("H2 = %d", chains->getBasisSize(1));
-  Msg::Info("H3 = %d", chains->getBasisSize(0));
-  if(omitted != 0) Msg::Info("%d 3D thick cuts are not shown completely.", omitted);
+  Msg::Info("H0 = %d", HRank[0]);
+  Msg::Info("H1 = %d", HRank[1]);
+  Msg::Info("H2 = %d", HRank[2]);
+  Msg::Info("H3 = %d", HRank[3]);
+  if(omitted != 0) Msg::Info("Computation of %dD thick cuts was omitted.", dim);
   
   
   Msg::Info("Wrote results to %s.", fileName.c_str());
diff --git a/benchmarks/3d/induction1.geo b/benchmarks/3d/induction1.geo
index 5b3aa2eb3c89de45ca7ec2f2d93d62c04d655b07..1692eba488e7467bec97a142aeaac32ecb88c3ef 100644
--- a/benchmarks/3d/induction1.geo
+++ b/benchmarks/3d/induction1.geo
@@ -313,6 +313,7 @@ Physical Surface(1001) = {64,66,68,105,54,58};  //gh0
 Physical Surface(700) = {75,1034}; //elec0
 Physical Surface(701) = {79,1032}; //elec1
 
+Physical Surface(702) = {75,1034,79,1032}; //elec1+2