Jorge Ivan Velez wrote:
>
> I think the function could be better but try this:
>
> # Function: M is your matrix and n MUST be an integer>0
> mat.pow<-function(M,n) {
> result<-M
> if(n>1){
> for ( iter in 2:n) result<-M%*%result
> result
> }
> else {result}
> result
> }
>
There are much more efficient ways to compute a power,
just check:
http://en.wikipedia.org/wiki/Exponentiation_by_squaring
Or, translating the pseudo-code to R:
mat.pow <- function(M, n) {
result <- diag(1, ncol(M))
while (n > 0) {
if (n %% 2 == 1) {
result <- result %*% M
n <- n - 1
}
M <- M %*% M
n <- n / 2
}
return(result)
}
# The matrix
m <- rbind(c(1,1,0), c(0,1,1), c(0,0,1))
# Goal m^6 = m x m x m x m x m x m
goal= m %*% m %*% m %*% m %*% m %*% m
# matpow
res=mat.pow(m,6)
# Check point
all.equal(goal,res)
This algorithm would be fast, unless n is a _very_ big number.
Alberto Monteiro
______________________________________________
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
