diff --git a/include/environment.h b/include/environment.h
index 9bfcb79..558b3fc 100644
--- a/include/environment.h
+++ b/include/environment.h
@@ -32,6 +32,3 @@ class GravitationalEnvironment{
 
 // Helper functions
 int getLargestLabelNumber(const std::vector<std::string>& filenames, const std::string logFilePrefix);
-
-
-typedef GravitationalEnvironment<Particle> ParticleEnvironment;
diff --git a/include/octree.h b/include/octree.h
new file mode 100644
index 0000000..1a59aa3
--- /dev/null
+++ b/include/octree.h
@@ -0,0 +1,39 @@
+#pragma once
+
+#include <vector>
+#include "body.h"
+
+template <typename T>
+class Octree {
+    public:
+
+        // Octree constructor
+        Octree(std::array<float, 2>& xCoords, std::array<float, 2>& yCoords, std::array<float, 2>& zCoords);
+
+        // Member functions
+        void clear();
+        void insert(std::shared_ptr<T> objPtr);
+        void build();
+        void updateCenterOfMass();
+
+        // Members
+        std::vector<Body*> objPtrs;
+        std::array<float, 3> centerOfMass;
+        float* totalMass;
+
+        // Dimensions of the current octant
+        std::array<float, 2>* xCoords;
+        std::array<float, 2>* yCoords;
+        std::array<float, 2>* zCoords;
+
+        // Octree children --> 0-7 based on 2D convention in postive z, and then 2D convention in negative z, observing from above        
+        Octree<T>* child0;
+        Octree<T>* child1;
+        Octree<T>* child2;
+        Octree<T>* child3;
+        Octree<T>* child4;
+        Octree<T>* child5;
+        Octree<T>* child6;
+        Octree<T>* child7;
+
+};
diff --git a/src/octree.cpp b/src/octree.cpp
new file mode 100644
index 0000000..bf636f2
--- /dev/null
+++ b/src/octree.cpp
@@ -0,0 +1,124 @@
+#include "./../include/octree.h"
+#include <array>
+
+template <typename T>
+Octree<T>::Octree(std::array<float, 2>& xCoords, std::array<float, 2>& yCoords, std::array<float, 2>& zCoords)
+    : xCoords(xCoords), yCoords(yCoords), zCoords(zCoords), totalMass(0) {};
+
+// Recursively set every child to null in the tree, but preserving the tree
+template <typename T>
+void Octree<T>::clear() {    
+    
+    if (child0 != nullptr) {
+        child0.clear();
+        child0 = nullptr;
+    }
+
+    if (child1 != nullptr) {
+        child1.clear();
+        child1 = nullptr;
+    }
+
+    if (child2 != nullptr) {
+        child2.clear();
+        child2 = nullptr;
+    }
+
+    if (child3 != nullptr) {
+        child3.clear();
+        child3 = nullptr;
+    }
+
+    if (child4 != nullptr) {
+        child4.clear();
+        child4 = nullptr;
+    }
+
+    if (child5 != nullptr) {
+        child5.clear();
+        child5 = nullptr;
+    }
+
+    if (child6 != nullptr) {
+        child6.clear();
+        child6 = nullptr;
+    }
+
+    if (child7 != nullptr) {
+        child7.clear();
+        child7 = nullptr;
+    }
+}
+
+template <typename T>
+void Octree<T>::insert(std::shared_ptr<T> objPtr) {
+
+    // Append objPtr to the vector of objects
+    objPtrs.push_back(objPtr);
+
+    // Instantiate the center of mass if it doesn't exist; update it otherwise
+    float newTotalMass = totalMass + objPtr->mass;
+    if (totalMass == 0) {
+        centerOfMass = objPtr->position;
+    } else {
+        for (int i = 0; i < 3; i++) {
+            centerOfMass[i] = (((centerOfMass[i] * totalMass)) + ((objPtr->position[i]) * (objPtr->mass))) / (newTotalMass);
+        }
+    }
+
+    // Update the total mass
+    totalMass = newTotalMass;
+
+    // Midpoints of coordinates
+    float mX = (xCoords[0] + xCoords[1]) / 2.;
+    float mY = (yCoords[0] + yCoords[1]) / 2.;
+    float mZ = (zCoords[0] + zCoords[1]) / 2.;
+
+
+    // Check for recursive insertion --> if greater than 1 pointer, then need to recursive insert
+    if (objPtrs.size() > 1) {
+        
+        // Initialize flags for octant location
+        bool xFlag = objPtr->position[0] > mX;
+        bool yFlag = objPtr->position[1] > mY; 
+        bool zFlag = objPtr->position[2] > mZ;
+
+        std::array<float, 2> xCoordsNew;
+        std::array<float, 2> yCoordsNew;
+        std::array<float, 2> zCoordsNew;
+
+
+        // Get the X coordinates
+        if (xFlag) {
+            xCoordsNew[0] = mX;
+            xCoordsNew[1] = xCoords[1];
+        } else {
+            xCoordsNew[0] = xCoords[0];
+            xCoordsNew[1] = mX;
+        }
+
+        // Get the Y coordinates
+        if (yFlag) {
+            yCoordsNew[0] = mY;
+            yCoordsNew[1] = yCoords[1];
+        } else {
+            yCoordsNew[0] = yCoords[0];
+            yCoordsNew[1] = mY;
+        }
+
+        // Get the Z coordinates
+        if (zFlag) {
+            zCoordsNew[0] = mZ;
+            zCoordsNew[1] = zCoords[1];
+        } else {
+            zCoordsNew[0] = zCoords[0];
+            zCoordsNew[1] = mZ;
+        }
+
+        // TO DO: instantiate a new octree with the calculated coordinates and the correct child
+        
+
+
+
+    }
+}
\ No newline at end of file