#EDITS: minor fixes/updates

experiment/eigen_ex.cpp

@@ -0,0 +1,41 @@
+#include <iostream>
+#include <Eigen/Dense>
+
+using namespace Eigen;
+using namespace std;
+
+int main() {
+    // Define the matrix A
+    MatrixXd A(3, 3);
+    /*A << 1, 2, 3,
+         4, 5, 6,
+         7, 8, 9;
+    */
+    A << 2.0, -1.0, 0.0,
+        -1.0, 2.0, -1.0,
+        0.0, -1.0, 2.0;
+
+    // Initial guess for the eigenvector
+    VectorXd v = VectorXd::Random(3);
+
+    // Power iteration method
+    int iterations = 1000; // Maximum number of iterations
+    double tolerance = 1e-3; // Tolerance for convergence
+
+    for (int i = 0; i < iterations; ++i) {
+        VectorXd Av = A * v;
+        double lambda = v.dot(Av) / v.dot(v); // Approximation of the dominant eigenvalue
+        VectorXd v_new = Av / Av.norm(); // Normalized eigenvector
+
+        if ((v_new - v).norm() < tolerance) {
+            cout << "Dominant eigenvalue: " << lambda << endl;
+            cout << "Eigenvector: \n" << v_new << endl;
+            break; // Convergence achieved
+        }
+
+        v = v_new;
+    }
+
+    return 0;
+}
+

experiment/eigen_ex.py

@@ -0,0 +1,49 @@
+import numpy as np 
+
+# Define the matrix A 
+A = np.array([
+                [2.0, -1.0, 0.0], 
+                [-1.0, 2.0, -1.0],
+                [0.0, -1.0, 2.0]
+              ]) 
+
+
+
+# Choose the initial vector x 
+x = np.array([[1, 2]]).T 
+
+# Define the tolerance for the eigenvalue 
+# and eigenvector approximations 
+# (i.e. the maximum allowed difference between 
+# the approximations and the actual values) 
+tol = 1e-6
+
+# Define the maximum number of iterations 
+max_iter = 100
+
+# Define the variable lam_prev to store the 
+# previous approximation for the largest eigenvalue 
+lam_prev = 0
+
+# Iteratively improve the approximations 
+# for the largest eigenvalue and eigenvector 
+# using the power method 
+for i in range(max_iter): 
+    # Compute the updated approximation for the eigenvector 
+    x = A @ x / np.linalg.norm(A @ x) 
+
+    # Compute the updated approximation for the largest eigenvalue 
+    lam = (x.T @ A @ x) / (x.T @ x) 
+
+    # Check if the approximations have converged 
+    if np.abs(lam - lam_prev) < tol: 
+        break
+
+    # Store the current approximation for the largest eigenvalue 
+    lam_prev = lam 
+
+# Print the approximations for the 
+# largest eigenvalue and eigenvector 
+print(float(lam)) 
+print(x) 
+

experiment/eigen_ex0.cpp

@@ -0,0 +1,15 @@
+#include <iostream>
+#include <Eigen/Dense>
+
+using Eigen::MatrixXd;
+
+int main()
+{
+  MatrixXd m(2,2);
+  m(0,0) = 3;
+  m(1,0) = 2.5;
+  m(0,1) = -1;
+  m(1,1) = m(1,0) + m(0,1);
+  std::cout << m << std::endl;
+}
+

modules/linalg/CMakeLists.txt

@@ -4,6 +4,7 @@ message("-- Including Linear Algebra Module")
 
 option(BUILD_GPU_OPENCL "Build OpenCL Linear Algebra Extension" OFF)
 
+# BASE SOURCE FILES
 set(SOURCE_FILES
     mtx.cpp
     mtx_routines.f90
@@ -22,7 +23,7 @@ set(SOURCE_FILES
     dgemm_arr.cpp
 )
 
-# Add files depending on the detected SIMD ISA
+# Add files depending on the detected ISA
 if (SIMD_ISA STREQUAL "AVX2")
     list(APPEND SOURCE_FILES
         avx/mtx_avx2_arr_i8.cpp

modules/linalg/mtx.cpp

@@ -21,8 +21,7 @@
  * You may opt to use, copy, modify, merge, publish, distribute
  * and/or sell copies of the Software, and permit persons to whom
  * the Software is furnished to do so, under the terms of the
- * LICENSE file. As this is an Open Source effort, all implementations
- * must be of the same methodology.
+ * LICENSE file.
  *
  *
  *
@@ -34,6 +33,7 @@
 /** Main Matrix operation interface for matrices as arrays & vectors */
 
 #include <cmath>
+#include <vector>
 #include <limits>
 #include <openGPMP/linalg/_dgemm.hpp>
 
@@ -42,3 +42,12 @@
 #include <cblas.h>
 
 #endif
+
+void mtx_mult(std::vector<double> A, std::vector<double> B, std::vector<double> C);
+
+
+//void mtx_mult(char *A, char *B, char *C) {
+//    gpmp::linalg::DGEMM dgemm;
+//}
+
+

samples/cpp/eigen.cpp

@@ -0,0 +1,17 @@
+#include <iostream>
+#include <openGPMP/linalg/eigen.hpp>
+#include <vector>
+
+const double TOLERANCE = 1e-3;
+
+int main () {
+    std::vector<std::vector<double>> squareMat = {{2.0, -1.0, 0.0},
+                                                  {-1.0, 2.0, -1.0},
+                                                  {0.0, -1.0, 2.0}};
+    gpmp::linalg::Eigen eigenSquare(squareMat);
+
+    double eigenvalueSquare = eigenSquare.power_iter(1000, TOLERANCE);
+
+    std::cout << "Eigenvalue square: " << eigenvalueSquare << std::endl;
+
+}