0

I have a matrix: 

P <- matrix(c(0, 0.5, 0, 0, 0.5, 0, 
              0.5, 0, 0.5, 0, 0, 0, 
              0, 0.5, 0, 0.5, 0, 0, 
              0, 0, 0.5, 0, 0, 0.5, 
              0, 0, 0, 0, 1, 0, 
              0, 0, 0, 0, 0, 1), 
            nrow = 6, ncol = 6, byrow = TRUE)

Can anyone help me with the code to find the n that makes P^n = P^n+1 with an accuracy of 4 decimal places.

The code under here, is how I´ve attempted the problem. Although, I'm not so good at coding so the last part is not correct. Is there something that i can add to it so it works? Or is there another way to solve it? If so please explain the code you write.

mpow <- function(P, n) {
if (n == 0) {
return(diag(nrow(P)))} 
else if (n == 1){ 
return(P)
} else {
return(P %*% mpow(P, n - 1))
}
}

matrices_equal <- function(A, B, d = 4) { 
A_new <- trunc(A * 10^d) 
B_new <-trunc(B * 10^d) 
if (all(A_new == B_new)) { 
return(TRUE) } 
else { 
return(FALSE) } }

when_converged <- function(P) 
for (n in 1:50) {{
A <- mpow(P, n)
B <- mpow(P, n+1)
} 
{ if matrices_equal(A, B, 4)
return(n)
}}

Thank you.

user20650
  • 24,654
  • 5
  • 56
  • 91
PeterNiklas
  • 75
  • 1
  • 9
  • For you criterion: "P^n = P^n+1 with an accuracy of 4 decimal places," do you mean that each element of P^n and P^n+1 are within 4 decimal places? – lmo Jul 06 '16 at 17:53
  • You also have some invalid syntax on the last 3 lines just FYI. It would help if you gave more description about what you're trying to do. – Hack-R Jul 06 '16 at 17:54
  • What i meant was that an the difference between an element in P^n and the same element P^n+1 is less than 0.0005. To Hack-R i added my code for matrices_equal, maybe now you can see what i was trying to do. The code is still not working tho. – PeterNiklas Jul 06 '16 at 18:07
  • To clarify abit extra, the question is for what n is P^n+1 = P^n – PeterNiklas Jul 06 '16 at 18:10
  • 2
    `when_converged <- function(P, tol=0.0001) {; n = 1; diff = 1;; while(diff > tol){; A <- mpow(P, n); B <- mpow(P, n+1); diff <- max(abs(A - B)); n <- n + 1; } ; return(n); }` – user20650 Jul 06 '16 at 18:11
  • User20650, Thank you. It looks good, although the function is not returning a "n". Also, kan you explain diff <- max(abs(A-B) why are we using maximum & absolute values? – PeterNiklas Jul 06 '16 at 18:16
  • when_converged(P, 0.0001)? – PeterNiklas Jul 06 '16 at 18:23
  • 1
    It returns the n for me. As we want P^n = P^n+1 up to tolerance, of say, 0.0001, we can think of this as P^n - P^n+1 < 0.0001. I used `abs`solute to treat positive or negative differences in the matrices the same. I used `max` as this will be the largest absolute difference between the matrices. if the largest difference is less than the tolerance then all differences will be. [ps: im sure there will be a clever way to do this] – user20650 Jul 06 '16 at 18:25

0 Answers0