diff --git a/src/main.cpp b/src/main.cpp
index 1a1b1d8..a21aa53 100644
--- a/src/main.cpp
+++ b/src/main.cpp
@@ -19,6 +19,22 @@
 #include "matrix.hpp"
 #include "jacobi.hpp"
 
+// void analytical_eigenvec(arma::mat &analytic, int N) 
+// {   
+//     arma::vec analytic_vec = arma::vec(N);
+//     for (int i=0; i < N; i++) {
+//         for (int j=0; j < N; j++) {
+//             analytic_vec(j) = std::sin(((j+1.)*(i+1.)*M_PI) / (N+1.));
+//         }
+//         analytic_vec = arma::normalise(analytic_vec);
+
+//         // Flip the sign of the analytic vector if they are different
+//         if (analytic_vec(0)*v(0) < 0.) {
+//             analytic_vec *= -1;
+//         }
+//     }
+// }
+
 
 void write_transformation_dense(int N)
 {   
@@ -87,6 +103,7 @@ void write_eigenvec(int N)
 
     // Create tridiagonal matrix
     arma::mat A = create_symmetric_tridiagonal(N, a, d);
+    arma::mat analytic = arma::mat(N, N);
     arma::vec eigval;
     arma::mat eigvec;
 
@@ -96,15 +113,38 @@ void write_eigenvec(int N)
     // Solve using Jacobi rotation method
     jacobi_eigensolver(A, 10e-14, eigval, eigvec, 100000, iters, converged);
 
+    // Build analytic eigenvectors 
+    arma::vec v, analytic_vec = arma::vec(N);
+    for (int i=0; i < N; i++) {
+        v = eigvec.col(i);
+
+        for (int j=0; j < N; j++) {
+            analytic_vec(j) = std::sin(((j+1.)*(i+1.)*M_PI) / (N+1.));
+        }
+
+        analytic_vec = arma::normalise(analytic_vec);
+
+        // Flip the sign of the analytic vector if they are different
+        if (analytic_vec(0)*v(0) < 0.) {
+            analytic_vec *= -1;
+        }
+        analytic.col(i) = analytic_vec;
+    }
+
     std::ofstream ofile;
     // Create file based on matrix size, and write header line to file 
     ofile.open("../latex/output/eigenvector_" + std::to_string(N) + ".csv");
-    ofile << "x,Vector 1,Vector 2,Vector 3" << std::endl;
+    ofile << "x," 
+        << "Vector 1,Vector 2,Vector 3," 
+        << "Analytic 1,Analytic 2,Analytic 3" << std::endl;
 
-    // Add boundary value for x=0
+    // Add boundary value for x=0 
     ofile << scientific_format(0., 16) << "," 
             << scientific_format(0., 16) << "," 
             << scientific_format(0., 16) << "," 
+            << scientific_format(0., 16) << "," 
+            << scientific_format(0., 16) << "," 
+            << scientific_format(0., 16) << ","
             << scientific_format(0., 16) << std::endl;
 
     // Add x-value and element i of each eigenvector to same line
@@ -113,12 +153,18 @@ void write_eigenvec(int N)
         ofile << scientific_format(x, 16)<< "," 
             << scientific_format(eigvec(i,0), 16) << "," 
             << scientific_format(eigvec(i,1), 16) << "," 
-            << scientific_format(eigvec(i,2), 16) << std::endl;
+            << scientific_format(eigvec(i,2), 16) << ","
+            << scientific_format(analytic(i,0), 16) << ","
+            << scientific_format(analytic(i,1), 16) << ","
+            << scientific_format(analytic(i,2), 16) << std::endl;
     }
     // Add boundary value for x=1
     ofile << scientific_format(1., 16) << "," 
             << scientific_format(0., 16) << "," 
             << scientific_format(0., 16) << "," 
+            << scientific_format(0., 16) << ","
+            << scientific_format(0., 16) << ","
+            << scientific_format(0., 16) << ","
             << scientific_format(0., 16) << std::endl;
     ofile.close();
 }
@@ -126,8 +172,8 @@ void write_eigenvec(int N)
 
 int main()
 {
-    write_transformation_tridiag(100);
-    write_transformation_dense(100);
+    // write_transformation_tridiag(100);
+    // write_transformation_dense(100);
     write_eigenvec(6);
     write_eigenvec(100);
     return 0;