diff --git a/Cargo.toml b/Cargo.toml
index 29e63e2..c1c64c0 100644
--- a/Cargo.toml
+++ b/Cargo.toml
@@ -1,6 +1,6 @@
 [package]
 name = "rstats"
-version = "2.1.9"
+version = "2.1.10"
 authors = ["Libor Spacek"]
 edition = "2021"
 description = "Statistics, Information Measures, Data Analysis, Linear Algebra, Clifford Algebra, Machine Learning, Geometric Median, Matrix Decompositions, Mahalanobis Distance, Hulls, Multithreading.."
diff --git a/README.md b/README.md
index a5aa7cf..44796b4 100644
--- a/README.md
+++ b/README.md
@@ -299,7 +299,7 @@ datavec.somemethod::<f64>(arg)
 
 Methods implemented by this trait:
 
-* Vector additions, subtractions and products (scalar, kronecker, outer),
+* Vector additions, subtractions and products (scalar, Kronecker, outer),
 * Other relationships and measures of difference,
 * Pearson's, Spearman's and Kendall's correlations,
 * Joint pdf, joint entropy, statistical independence (based on mutual information).
@@ -341,6 +341,8 @@ Methods which take an additional generic vector argument, such as a vector of we
 
 ## Appendix: Recent Releases
 
+* **Version 2.1.10** - Added `project` of a `TriangMat` to a subspace given by a subspace index.
+
 * **Version 2.1.9** - Added multiplications and more tests for `TriangMat`.
 
 * **Version 2.1.8** - Improved `TriangMat::diagonal()`, restored `TriangMat::determinant()`, tidied up `triangmat` test.
diff --git a/src/triangmat.rs b/src/triangmat.rs
index 262d8bf..8b3d682 100644
--- a/src/triangmat.rs
+++ b/src/triangmat.rs
@@ -89,25 +89,25 @@ impl TriangMat {
         TriangMat { kind: 2, data }
     }
 
-    /// Project symmetric/antisymmetric triangmat to a smaller one of the same kind,
-    /// into a subspace specified by an ascending index of dimensions.
-    /// Deletes all rows and columns of the missing dimensions.
+    /// Projects to a smaller TriangMat of the same kind,
+    /// in a subspace given by a subspace index. 
+    /// Deletes all the rows and columns of the other dimensions.
+    /// The kept ones retain their original order.
     pub fn project(&self, index: &[usize]) -> Self {
         let mut res = Vec::with_capacity(sumn(index.len()));
-        for &row_idx in index {
-            let row = self.row(row_idx);
-            for &column_idx in index {
-                if column_idx >= row.len() {
-                    break;
-                };
-                res.push(row[column_idx]);
-            }
-        }
+        for &rownum in index { 
+            let row = self.row(rownum);
+            for &colnum in index {
+                if colnum > rownum { break; };
+                res.push(row[colnum]);
+            };
+        }; 
         TriangMat {
             kind: self.kind,
             data: res,
         }
     }
+
     /// Copy one raw data row from TriangMat
     /// To interpret the kind (plain, symmetric, assymetric, transposed),
     /// use `realrow,realcolumn,to_full`
diff --git a/tests/tests.rs b/tests/tests.rs
index ef474cc..336f660 100644
--- a/tests/tests.rs
+++ b/tests/tests.rs
@@ -246,14 +246,14 @@ fn triangmat() -> Result<(), RE> {
     println!("Symmetric positive definite matrix A\n{}",cov.gr());
     println!("Full form of A:\n{}",cov.to_full().gr());
     let mut chol = cov.cholesky()?;
-    println!("Cholesky L matrix, such that A=LL'\n{}",chol.gr());
+    println!("Cholesky L matrix, such that A=LL'\n{}",chol.gr());  
     println!("Diagonal of L: {}",chol.diagonal().gr());
     println!("Determinant det(A): {}",chol.determinant().gr());
     let full = chol.to_full();
     println!("Full L matrix\n{}",full.gr()); 
     let tchol = &chol.clone_transpose();
-    println!("A reconstructed from full(L)full(L'):\n{}",full.matmult(&tchol.to_full())?.gr());
-    println!("A reconstructed by more efficient triangular multiplication LL'\n{}",chol.mult(tchol).gr());
+    println!("A reconstructed from full(L)full(L')\n{}",full.matmult(&tchol.to_full())?.gr());
+    println!("A reconstructed by direct triangular multiplication LL'\n{}",chol.mult(tchol).gr());
     let d = 4_usize;
     let v = ranv_f64(d)?; 
     println!("Random test vector:\n{}", v.gr());
@@ -281,11 +281,13 @@ fn triangmat() -> Result<(), RE> {
     println!("Row[2] of C\n{}",cov.realrow(2).gr());
     chol = cov.cholesky()?;
     println!("Cholesky of C:\n{GR}{chol}{UN} ");  
+    let small_chol = chol.project(&[0,2,4,5]);
+    println!("Projected chol:\n{GR}{small_chol}{UN} ");    
     println!("Determinant of C: {}",chol.determinant().gr());
     println!("Row[2] of Cholesky\n{}",chol.realrow(2).gr()); 
     println!("Column[2] of Cholesky\n{}",chol.realcolumn(2).gr()); 
     println!("C reconstructed by triangular multiplication LL'\n{}",
-    chol.mult(&TriangMat{kind:chol.kind+3,data:chol.data.clone()}).gr());
+    chol.mult(&chol.clone_transpose()).gr());
     Ok(())
 }