From ec42999812b9a4a311795e54618e84e3fa86cdc7 Mon Sep 17 00:00:00 2001
From: Jiri Pittner <jiri@pittnerovi.com>
Date: Fri, 28 Jul 2023 15:38:24 +0200
Subject: [PATCH] bugfix in permutation of smatrix

---
 smat.cc | 28 +++++++++++++++++++++++++---
 t.cc    |  4 +---
 2 files changed, 26 insertions(+), 6 deletions(-)

diff --git a/smat.cc b/smat.cc
index 9939a62..7a8c39b 100644
--- a/smat.cc
+++ b/smat.cc
@@ -311,13 +311,35 @@ const NRSMat<T> NRSMat<T>::permuted(const NRPerm<int> &p, const bool inverse) co
 if(!p.is_valid()) laerror("invalid permutation of smatrix");
 #endif
 int n=p.size();
-if(n!=(*this).size()) laerror("incompatible permutation and smatrix");
+if(n!=(*this).nrows()) laerror("incompatible permutation and smatrix");
 #ifdef CUDALA
         if(this->getlocation() != cpu || p.getlocation() != cpu ) laerror("permutations can be done only in CPU memory");
 #endif
 NRSMat<T> r(n);
-if(inverse) for(int i=1; i<=n; ++i) {int pi = p[i]-1; r(i-1,i-1) = (*this)(pi,pi);}
-else for(int i=1; i<=n; ++i) {int pi = p[i]-1; r(pi,pi) = (*this)(i-1,i-1);}
+if(inverse) 
+	{
+	for(int i=1; i<=n; ++i) 
+                {
+                int pi = p[i]-1; 
+                for(int j=1; j<=i; ++j)
+                        {
+                        int pj = p[j] - 1;
+                        r(i-1,j-1) = (*this)(pi,pj);
+                        }
+                }
+	}
+else 
+	{
+	for(int i=1; i<=n; ++i) 
+		{
+		int pi = p[i]-1; 
+		for(int j=1; j<=i; ++j)
+			{
+			int pj = p[j] - 1;
+			r(pi,pj) = (*this)(i-1,j-1);
+			}
+		}
+	}
 return r;
 }
 
diff --git a/t.cc b/t.cc
index 6f23cdf..7e7a640 100644
--- a/t.cc
+++ b/t.cc
@@ -2760,9 +2760,7 @@ cout <<"Clique cover is "<<cover<<endl;
 NRPerm<int> p(cover.size());
 cover.sort(0,p);
 cout<<"permutation to disentabgle the cliques = "<<p<<endl;
-NRVec<int> perm1(p.size());
-for(int i=0; i<p.size(); ++i) perm1[i] = p[i+1]-1;
-NRSMat<char> adjperm = adj.submatrix(perm1);
+NRSMat<char> adjperm = adj.permuted(p);
 cout <<"resorted graph = "<<adjperm<<endl;
 }