diff --git a/Geo/MElementOctree.cpp b/Geo/MElementOctree.cpp
index 0b1da70271fe8ef2ffcd1172ab384794dd88f8d8..854347d30b7ba45e7b09a0a0eb727b9eaad908ce 100644
--- a/Geo/MElementOctree.cpp
+++ b/Geo/MElementOctree.cpp
@@ -8,21 +8,45 @@
 #include "MElementOctree.h"
 #include "Octree.h"
 #include "Context.h"
+#include "fullMatrix.h"
+#include "bezierBasis.h"
+#include "BasisFactory.h"
+#include "FuncSpaceData.h"
 
 void MElementBB(void *a, double *min, double *max)
 {
   MElement *e = (MElement*)a;
-  MVertex *v = e->getVertex(0);
-  min[0] = max[0] = v->x();
-  min[1] = max[1] = v->y();
-  min[2] = max[2] = v->z();
-  for(int i = 1; i < e->getNumVertices(); i++){
-    v = e->getVertex(i);
-    min[0] = std::min(min[0], v->x()); max[0] = std::max(max[0], v->x());
-    min[1] = std::min(min[1], v->y()); max[1] = std::max(max[1], v->y());
-    min[2] = std::min(min[2], v->z()); max[2] = std::max(max[2], v->z());
+  if (e->getPolynomialOrder() == 1) {
+    MVertex *v = e->getVertex(0);
+    min[0] = max[0] = v->x();
+    min[1] = max[1] = v->y();
+    min[2] = max[2] = v->z();
+    for(int i = 1; i < e->getNumVertices(); i++){
+      v = e->getVertex(i);
+      min[0] = std::min(min[0], v->x()); max[0] = std::max(max[0], v->x());
+      min[1] = std::min(min[1], v->y()); max[1] = std::max(max[1], v->y());
+      min[2] = std::min(min[2], v->z()); max[2] = std::max(max[2], v->z());
+    }
   }
+  else {
+    fullMatrix<double> nodesXYZ(e->getNumVertices(), 3);
+    e->getNodesCoord(nodesXYZ);
 
+    fullMatrix<double> bezNodes(e->getNumVertices(), 3);
+    bezierBasis *bez = BasisFactory::getBezierBasis(FuncSpaceData(e));
+    bez->lag2Bez(nodesXYZ, bezNodes);
+    min[0] = max[0] = bezNodes(0, 0);
+    min[1] = max[1] = bezNodes(0, 1);
+    min[2] = max[2] = bezNodes(0, 2);
+    for(int i = 1; i < e->getNumVertices(); i++){
+      min[0] = std::min(min[0], bezNodes(0, 0));
+      max[0] = std::max(max[0], bezNodes(0, 0));
+      min[1] = std::min(min[1], bezNodes(0, 1));
+      max[1] = std::max(max[1], bezNodes(0, 1));
+      min[2] = std::min(min[2], bezNodes(0, 2));
+      max[2] = std::max(max[2], bezNodes(0, 2));
+    }
+  }
   // make bounding boxes larger up to (absolute) geometrical tolerance
   double eps = CTX::instance()->geom.tolerance;
   for(int i = 0; i < 3; i++){