2
\$\begingroup\$

My own Givens-based QR decomposition function in R (pseudocode from Politecnico of Turin's math department) is the following:

>Givens.fn<-function(V) {

   tol<-1.0*10^-14
   m<-dim(V)[1]
   n<-dim(V)[2]

   spId<-bandSparse(m,m,0,list(rep(1, m+1))) # computed just once 
   Q<-spId                                   
   R<-V
   sapply(1:n, function(j){                      
      if (j<m) {    
         vapply((j+1):m, function(i){          # vectorized internal loop  
         if(abs(R[i,j])>tol) {                             
            G<-spId                  
            x<-R[j,j]
            y<-R[i,j]
            norm<-sqrt(x^2+y^2)
            c<-x/norm
            s<-y/norm
            G[j,j]=c
            G[i,i]=c
            G[j,i]=s
            G[i,j]=-s                         
            Q<<-G%*%Q
            R<<-G%*%R
         } 
         return(1)                   # --> saves 15% execution time!
      },FUN.VALUE=1.0)
   }
})

Q<-t(Q)

return(list(Q,R))
}

It works fine, but the givens() function provided in R by pracma is 134x faster on an total elapsed time around 10" (on my intel i5-3570 4core, 8GB desktop):

>m<-20; n<-20
>set.seed(1)
>X <- as.matrix(replicate(n, runif(m)))
>library(rbenchmark)
>library(pracma)
>benchmark(Givens.fn(X),
           givens(X),
           order = "elapsed",
           replications = 10)

          test replications elapsed relative user.self sys.self user.child sys.child
2    givens(X)           10    0.08    1.000      0.08     0.00         NA        NA
1 Givens.fn(X)           10   10.77  134.625     10.69     0.08         NA        NA

Unfortunately, the pracma version accepts just square matrices, which is not my case. Any further code improvements? Any other Givens-based QR functions (even in C++)? Thanks.

\$\endgroup\$
  • \$\begingroup\$ What is QR? Is it related to QR codes? Please edit your question to include this information. \$\endgroup\$ – Roland Illig May 22 at 22:13
  • \$\begingroup\$ All the 3 ways to compute QR decomposition are well described here \$\endgroup\$ – Antonio Piemontese May 26 at 20:59
  • \$\begingroup\$ As I said: please edit the question to include this information. A comment is the wrong place for it. \$\endgroup\$ – Roland Illig May 26 at 21:35
  • \$\begingroup\$ Could you perhaps remove the leading > from your code? This about syntax errors when pasting the code. \$\endgroup\$ – Roland Illig May 26 at 21:40
  • \$\begingroup\$ You could also format your code according to the Tidyverse convention. This way, the reviewers can focus on the interesting part of the code instead of complainingthatthecodeisunreadablebecauseofmissingspaces. \$\endgroup\$ – Roland Illig May 26 at 21:42

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Browse other questions tagged or ask your own question.