diff --git a/FunctionSpace/FunctionSpace.cpp b/FunctionSpace/FunctionSpace.cpp
index 521aaf36e06fb0d77dcb60fb84e2c102444ecc72..871ae5056549a53a04853502253bddd1e97843e2 100644
--- a/FunctionSpace/FunctionSpace.cpp
+++ b/FunctionSpace/FunctionSpace.cpp
@@ -37,7 +37,6 @@ void FunctionSpace::build(const GroupOfElement& goe, string family){
   basis.resize(nGeoType, NULL);
 
   // Generate Bases //
-  const vector<const MElement*>& element = goe.getAll();
 
   // Get geomtrical type statistics
   const vector<size_t>& geoTypeStat = goe.getTypeStats();
@@ -49,27 +48,39 @@ void FunctionSpace::build(const GroupOfElement& goe, string family){
       basis[i] = BasisGenerator::generate(i, form, order, family);
 
   // Get Number of Function per Entity //
-  // Same for all basis since we have a uniform order
-  MElement& myElement = const_cast<MElement&>(*element[0]);
-
-  int nVertex = myElement.getNumPrimaryVertices();
-  int nEdge   = myElement.getNumEdges();
-  int nFace   = myElement.getNumFaces();
-  int type    = myElement.getType();
-
-  fPerVertex = basis[type]->getNVertexBased() / nVertex;
-
-  if(nEdge)
-    fPerEdge = this->basis[type]->getNEdgeBased() / nEdge;
-  else
-    fPerEdge = 0;
-
-  if(nFace)
-    fPerFace = this->basis[type]->getNFaceBased() / nFace;
-  else
-    fPerFace = 0;
+  fPerVertex.resize(nGeoType);
+  fPerEdge.resize(nGeoType);
+  fPerFace.resize(nGeoType);
+  fPerCell.resize(nGeoType);
+
+  for(size_t i = 0; i < nGeoType; i++){
+    if(geoTypeStat[i]){
+      int nVertex = ReferenceSpaceManager::getNVertex(i);
+      int nEdge   = ReferenceSpaceManager::getNEdge(i);
+      int nFace   = ReferenceSpaceManager::getNFace(i);
+
+      fPerVertex[i] = basis[i]->getNVertexBased() / nVertex;
+
+      if(nEdge)
+        fPerEdge[i] = this->basis[i]->getNEdgeBased() / nEdge;
+      else
+        fPerEdge[i] = 0;
+
+      if(nFace)
+        fPerFace[i] = this->basis[i]->getNFaceBased() / nFace;
+      else
+        fPerFace[i] = 0;
+
+      fPerCell[i] = this->basis[i]->getNCellBased();
+    }
 
-  fPerCell = this->basis[type]->getNCellBased();
+    else{
+      fPerVertex[i] = 0;
+      fPerEdge[i]   = 0;
+      fPerFace[i]   = 0;
+      fPerCell[i]   = 0;
+    }
+  }
 
   // Build Dof //
   buildDof();
@@ -107,13 +118,6 @@ vector<Dof> FunctionSpace::getUnorderedKeys(const MElement& elem) const{
   const size_t nEdge   = element.getNumEdges();
   const size_t nFace   = element.getNumFaces();
 
-  // Number of cells //
-  size_t nCell;
-  if(element.getDim() == 3) // Need to be 3D to have one Cell
-    nCell = 1;
-  else
-    nCell = 0;
-
   vector<MVertex*> vertex(nVertex);
   vector<MEdge>    edge(nEdge);
   vector<MFace>    face(nFace);
@@ -128,11 +132,17 @@ vector<Dof> FunctionSpace::getUnorderedKeys(const MElement& elem) const{
     face[i] = element.getFace(i);
 
   // Create Dof //
+  const size_t type       = element.getType();
+  const size_t fPerVertex = this->fPerVertex[type];
+  const size_t fPerEdge   = this->fPerEdge[type];
+  const size_t fPerFace   = this->fPerFace[type];
+  const size_t fPerCell   = this->fPerCell[type];
+
   size_t nDof =
     fPerVertex * nVertex +
     fPerEdge   * nEdge   +
     fPerFace   * nFace   +
-    fPerCell   * nCell;
+    fPerCell;
 
   vector<Dof> myDof(nDof);
 
@@ -163,11 +173,9 @@ vector<Dof> FunctionSpace::getUnorderedKeys(const MElement& elem) const{
   }
 
   // Add Cell Based Dof //
-  for(size_t i = 0; i < nCell; i++){
-    for(size_t j = 0; j < fPerCell; j++){
-      myDof[it].setDof(mesh->getGlobalId(element), j);
-      it++;
-    }
+  for(size_t j = 0; j < fPerCell; j++){
+    myDof[it].setDof(mesh->getGlobalId(element), j);
+    it++;
   }
 
   return myDof;
@@ -221,3 +229,17 @@ void FunctionSpace::getKeys(const GroupOfElement& goe,
       dof.insert(myDof[d]);
   }
 }
+
+void FunctionSpace::getKeys(const GroupOfElement& goe,
+                            std::vector<std::vector<Dof> >& dof) const{
+  // Get Elements //
+  const size_t nElement = goe.getNumber();
+  const vector<const MElement*>& element = goe.getAll();
+
+  // Init Struct //
+  dof.resize(nElement);
+
+  // Create Dofs //
+  for(size_t i = 0; i < nElement; i++)
+    dof[i] = getKeys(*(element[i]));
+}
diff --git a/FunctionSpace/FunctionSpace.h b/FunctionSpace/FunctionSpace.h
index decf948b1a61769d27f65b140b27b2c5db48d492..9c3db8ccbf35898d3f47782f2e675c582b98ef09 100644
--- a/FunctionSpace/FunctionSpace.h
+++ b/FunctionSpace/FunctionSpace.h
@@ -42,10 +42,10 @@ class FunctionSpace{
   // Basis //
   std::vector<const Basis*> basis;
 
-  size_t fPerVertex;
-  size_t fPerEdge;
-  size_t fPerFace;
-  size_t fPerCell;
+  std::vector<size_t> fPerVertex;
+  std::vector<size_t> fPerEdge;
+  std::vector<size_t> fPerFace;
+  std::vector<size_t> fPerCell;
 
   // Differential From & Order //
   bool   scalar;
@@ -72,6 +72,8 @@ class FunctionSpace{
   std::vector<Dof> getKeys(const MElement& element) const;
 
   void getKeys(const GroupOfElement& goe, std::set<Dof>& dof) const;
+  void getKeys(const GroupOfElement& goe,
+               std::vector<std::vector<Dof> >& dof) const;
 
   const std::set<Dof>&                  getAllDofs(void)   const;
   const std::vector<std::vector<Dof> >& getAllGroups(void) const;
diff --git a/FunctionSpace/ReferenceSpace.h b/FunctionSpace/ReferenceSpace.h
index b926f8f38f9eecd5249a10d7f63f51cf8ce1897a..c875e582a0bb79a84f463ab31059baa5a7373058 100644
--- a/FunctionSpace/ReferenceSpace.h
+++ b/FunctionSpace/ReferenceSpace.h
@@ -92,6 +92,10 @@ class ReferenceSpace{
   ReferenceSpace(const std::string& path);
   virtual ~ReferenceSpace(void);
 
+  size_t getNVertex(void) const;
+  size_t getNEdge(void)   const;
+  size_t getNFace(void)   const;
+
   size_t getNOrientation(void) const;
   size_t getOrientation(const MElement& element) const;
 
@@ -348,6 +352,18 @@ class ReferenceSpace{
 // Inline Functions //
 //////////////////////
 
+inline size_t ReferenceSpace::getNVertex(void) const{
+  return nVertex;
+}
+
+inline size_t ReferenceSpace::getNEdge(void) const{
+  return refEdgeNodeIdx.size();
+}
+
+inline size_t ReferenceSpace::getNFace(void) const{
+  return refFaceNodeIdx.size();
+}
+
 inline size_t ReferenceSpace::getNOrientation(void) const{
   return refSpaceNodeId.size();
 }
diff --git a/FunctionSpace/ReferenceSpaceManager.h b/FunctionSpace/ReferenceSpaceManager.h
index 8ec773bd0be2b56c08bdba6883764039cd50833a..c36dc3a77a62187f8a2ebec4ce447c64b4c3b0f6 100644
--- a/FunctionSpace/ReferenceSpaceManager.h
+++ b/FunctionSpace/ReferenceSpaceManager.h
@@ -30,6 +30,10 @@ class ReferenceSpaceManager{
   static void clear(void);
   static const ReferenceSpace& getReferenceSpace(int elementType);
 
+  static size_t getNVertex(int elementType);
+  static size_t getNEdge(int elementType);
+  static size_t getNFace(int elementType);
+
   static size_t getNOrientation(int elementType);
   static size_t getOrientation(const MElement& element);
 
@@ -65,6 +69,27 @@ class ReferenceSpaceManager{
 // Inline Methods //
 ////////////////////
 
+inline size_t ReferenceSpaceManager::getNVertex(int elementType){
+  if(!refSpace[elementType])
+    init(elementType);
+
+  return refSpace[elementType]->getNVertex();
+}
+
+inline size_t ReferenceSpaceManager::getNEdge(int elementType){
+  if(!refSpace[elementType])
+    init(elementType);
+
+  return refSpace[elementType]->getNEdge();
+}
+
+inline size_t ReferenceSpaceManager::getNFace(int elementType){
+  if(!refSpace[elementType])
+    init(elementType);
+
+  return refSpace[elementType]->getNFace();
+}
+
 inline
 const ReferenceSpace& ReferenceSpaceManager::getReferenceSpace(int elementType){
   if(!refSpace[elementType])