diff --git a/cpp/Makefile b/cpp/Makefile
index bcffd1cbe550055ad9bdc9b98855dafa23afd236..8ffa5bd74bcdec1dcdf1f8286bf75dada688fdfd 100644
--- a/cpp/Makefile
+++ b/cpp/Makefile
@@ -1,9 +1,14 @@
 #CXX=clang++
 
-INCLUDES += -I.
-
-CXXFLAGS += -std=c++14 -O3 -Wall -Wextra -Wfatal-errors -Werror -pedantic -Wno-deprecated-declarations $(INCLUDES)
-
+INCLUDES += -I. \
+			-I/usr/include/eigen3
+
+CXXFLAGS += -std=c++14 -O3
+CXXFLAGS += -Wall -Wextra -Wfatal-errors -Werror -pedantic -Wno-deprecated-declarations
+CXXFLAGS += $(INCLUDES)
+# Disable Eigen's debug info and disable a warning generated by Eigen
+CXXFLAGS += -DNDEBUG
+CXXFLAGS += -Wno-int-in-bool-context
 
 all: switchchain switchchain_exponent switchchain_initialtris switchchain_dsp
 
diff --git a/cpp/graph.hpp b/cpp/graph.hpp
index ec318ad7193fbffdbac51fc7167fe825ac4dd653..6585dbd12d79407616158c9bd3d7108ed08eef51 100644
--- a/cpp/graph.hpp
+++ b/cpp/graph.hpp
@@ -58,6 +58,8 @@ class Graph {
 
     const auto& getAdj() const { return adj; }
 
+    const auto& getBooleanAdj() const { return badj; }
+
     // When the degree sequence is not graphics, the Graph can be
     // in any state, it is not neccesarily empty
     bool createFromDegreeSequence(const DegreeSequence &d) {
diff --git a/cpp/graph_spectrum.hpp b/cpp/graph_spectrum.hpp
new file mode 100644
index 0000000000000000000000000000000000000000..4be14e52ec06aff130a54c02533c228d011857a6
--- /dev/null
+++ b/cpp/graph_spectrum.hpp
@@ -0,0 +1,65 @@
+#include "graph.hpp"
+#include <Eigen/Dense>
+#include <Eigen/Eigenvalues>
+
+using MatrixType =
+    Eigen::Matrix<float, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor>;
+
+// A: Adjacency matrix
+//    lambda_max <= d_max
+//
+// L: Laplacian
+//    L = D - A
+///
+// P: Random walk matrix
+//    lambda_max = 1
+
+class GraphSpectrum {
+  public:
+    GraphSpectrum(const Graph& g) : graph(g) {}
+    ~GraphSpectrum() {}
+
+    std::vector<float> computeAdjacencySpectrum() const {
+        // matrix stored as std::vector<std::vector<bool>>
+        auto& badj = graph.getBooleanAdj();
+
+        // Convert it to MatrixType
+        auto n = badj.size();
+        MatrixType m(n, n);
+        for (auto i = 0u; i < n; ++i)
+            for (auto j = 0u; j < n; ++j)
+                m(i, j) = badj[i][j] ? 1.0f : 0.0f;
+
+        return getEigenvalues_(m);
+    }
+
+    std::vector<float> computeLaplacianSpectrum() const {
+        // matrix stored as std::vector<std::vector<bool>>
+        auto& badj = graph.getBooleanAdj();
+        auto& adj = graph.getAdj();
+
+        // - A
+        auto n = badj.size();
+        MatrixType m(n, n);
+        for (auto i = 0u; i < n; ++i)
+            for (auto j = 0u; j < n; ++j)
+                m(i, j) = badj[i][j] ? -1.0f : 0.0f;
+
+        // + D
+        for (auto i = 0u; i < n; ++i)
+            m(i, i) = float(adj[i].size());
+
+        return getEigenvalues_(m);
+    }
+
+  private:
+    const Graph& graph;
+
+    std::vector<float> getEigenvalues_(const MatrixType& m) const {
+        Eigen::SelfAdjointEigenSolver<MatrixType> es(
+            m, Eigen::DecompositionOptions::EigenvaluesOnly);
+        auto ev = es.eigenvalues();
+        return std::vector<float>(ev.data(), ev.data() + ev.rows() * ev.cols());
+    }
+};
+
diff --git a/cpp/switchchain.cpp b/cpp/switchchain.cpp
index 46c0ec8fee28118b9e38a9a7f3e5bba08d96ba1b..c28de798386c5f684874bda92e5ab545330f8104 100644
--- a/cpp/switchchain.cpp
+++ b/cpp/switchchain.cpp
@@ -1,6 +1,7 @@
 #include "exports.hpp"
 #include "graph.hpp"
 #include "powerlaw.hpp"
+#include "graph_spectrum.hpp"
 #include <algorithm>
 #include <array>
 #include <fstream>