removed end from mgs #8

src/algorithms/davidson.rs

@@ -102,7 +102,7 @@ impl EigenDavidson {
         for i in 0..conf.max_iters {
             // 2. Generate subpace matrix problem by projecting into the basis
             let subspace = basis.columns(0, dim_sub);
-            let matrix_proj = subspace.transpose() * &h.matrix_matrix_prod(subspace); // (&h * subspace);
+            let matrix_proj = subspace.transpose() * &h.matrix_matrix_prod(subspace);
 
             // 3. compute the eigenvalues and their corresponding ritz_vectors
             let ord_sort = match conf.spectrum_target {
@@ -137,8 +137,7 @@ impl EigenDavidson {
                 update_subspace(&mut basis, correction, dim_sub, dim_sub + nvalues);
 
                 // 6. Orthogonalize the subspace
-                basis = orthogonalize_subspace(basis, dim_sub, dim_sub + nvalues);
-
+                basis = orthogonalize_subspace(basis, 0);
                 // update counter
                 dim_sub += nvalues;
 
@@ -269,8 +268,8 @@ fn update_subspace(basis: &mut DMatrix<f64>, vectors: DMatrix<f64>, start: usize
 }
 
 /// Orthogonalize the subpsace using the QR method
-fn orthogonalize_subspace(vectors: DMatrix<f64>, start: usize, end: usize) -> DMatrix<f64> {
-    let mgs = MGS::new(vectors, start, end);
+fn orthogonalize_subspace(vectors: DMatrix<f64>, start: usize) -> DMatrix<f64> {
+    let mgs = MGS::new(vectors, start);
     match mgs {
         Ok(result) => result.basis,
         Err(msg) => panic!("Error orthonormalising the basis:{}", msg),

src/main.rs

@@ -5,8 +5,8 @@ use eigenvalues::davidson::EigenDavidson;
 use eigenvalues::SpectrumTarget;
 
 fn main() {
-    let brr = eigenvalues::utils::generate_diagonal_dominant(20, 0.5);
-    let eig = EigenDavidson::new(brr, 6, "GJD", SpectrumTarget::Lowest).unwrap();
+    let brr = eigenvalues::utils::generate_diagonal_dominant(10, 0.05);
+    let eig = EigenDavidson::new(brr, 2, "GJD", SpectrumTarget::Lowest).unwrap();
     println!("eigenvalues:{}", eig.eigenvalues);
     println!("eigenvectors:{}", eig.eigenvectors);
 }

src/modified_gram_schmidt.rs

@@ -20,10 +20,9 @@ impl MGS {
     /// * `vectors` to diagonalize as columns of the matrix
     /// * `start` index of the column to start orthogonalizing
     /// * `end` last index of the column to diagonalize (non-inclusive)
-    pub fn new(vectors: DMatrix<f64>, start: usize, end: usize) -> Result<Self, &'static str> {
-        let mut basis = vectors.clone();
+    pub fn new(mut basis: DMatrix<f64>, start: usize) -> Result<Self, &'static str> {
+	let end = basis.ncols();
         for i in start..end {
-            basis.set_column(i, &vectors.column(i));
             for j in 0..i {
                 let proj = MGS::project(&basis.column(j), &basis.column(i));
                 basis.set_column(i, &(basis.column(i) - proj));
@@ -59,8 +58,7 @@ mod test {
     }
 
     fn fun_test(vectors: DMatrix<f64>, start: usize) {
-        let end = vectors.ncols();
-        let mgs_result = super::MGS::new(vectors, start, end);
+        let mgs_result = super::MGS::new(vectors, start);
         let basis = match mgs_result {
             Ok(ortho) => ortho.basis,
             Err(message) => panic!(message),

tests/tests_davidson.rs

@@ -35,15 +35,16 @@ fn test_davidson_unsorted() {
 #[test]
 fn test_davidson_highest() {
     // Test the compution of the highest eigenvalues
-    let dim = 100;
-    let nvalues = 7;
+    let dim = 20;
+    let nvalues = 2;
     let arr = generate_diagonal_dominant(dim, 0.005);
     let eig = sort_eigenpairs(na::linalg::SymmetricEigen::new(arr.clone()), false);
 
     let target = SpectrumTarget::Highest;
-
+    println!("running DPR");
     let dav_eig = EigenDavidson::new(arr.clone(), nvalues, "DPR", target.clone()).unwrap();
     test_eigenpairs(&eig, dav_eig, nvalues);
+    println!("running GJD");
     let dav_eig = EigenDavidson::new(arr.clone(), nvalues, "GJD", target).unwrap();
     test_eigenpairs(&eig, dav_eig, nvalues);
 }