The R Project SVN R

Rev

Rev 29013 | Rev 41316 | Go to most recent revision | Show entire file | Ignore whitespace | Details | Blame | Last modification | View Log | RSS feed

Rev 29013 Rev 30396
Line 49... Line 49...
49 
 
 49 
 
50 
 
 50 
 
51 
 
 51 
 
52 
##  -------  tests of random real and complex matrices ------
 52 
##  -------  tests of random real and complex matrices ------
53 
 
 53 
 
54 
# 100 may cause failures here.
 54 
##			       100  may cause failures here.
55 
eigenok <- function(A, E, Eps=1000*.Machine$double.eps)
 55 
eigenok <- function(A, E, Eps=1000*.Machine$double.eps)
56 
{
 56 
{
57 
    V <- E$vect; lam <- E$values
 57 
    V <- E$vect; lam <- E$values
58 
    stopifnot(abs(A %*% V - V %*% diag(lam)) < Eps,
 58 
    stopifnot(abs(A %*% V - V %*% diag(lam)) < Eps,
- 
 
 59 
              abs(lam[length(lam)]/lam[1]) < Eps || # this one not for singular A :
59 
              abs(A - V %*% diag(lam) %*% t(V)) < Eps)
 60 
              abs(A - V %*% diag(lam) %*% t(V)) < Eps)
60 
}
 61 
}
61 
 
 62 
 
62 
Ceigenok <- function(A, E, Eps=1000*.Machine$double.eps)
 63 
Ceigenok <- function(A, E, Eps=1000*.Machine$double.eps)
63 
{
 64 
{
64 
    V <- E$vect; lam <- E$values
 65 
    V <- E$vect; lam <- E$values
65 
    stopifnot(Mod(A %*% V - V %*% diag(lam)) < Eps,
 66 
    stopifnot(Mod(A %*% V - V %*% diag(lam)) < Eps,
66 
              Mod(A - V %*% diag(lam) %*% Conj(t(V))) < Eps)
 67 
              Mod(A - V %*% diag(lam) %*% Conj(t(V))) < Eps)
67 
}
 68 
}
68 
 
 69 
 
- 
 
 70 
## failed for some 64bit-Lapack-gcc combinations:
- 
 
 71 
sm <- cbind(1, 3:1, 1:3)
- 
 
 72 
eigenok(sm, eigen(sm))
- 
 
 73 
eigenok(sm, eigen(sm, sym=FALSE))
- 
 
 74 
 
69 
set.seed(123)
 75 
set.seed(123)
70 
sm <- matrix(rnorm(25), 5, 5)
 76 
sm <- matrix(rnorm(25), 5, 5)
71 
sm <- 0.5 * (sm + t(sm))
 77 
sm <- 0.5 * (sm + t(sm))
72 
eigenok(sm, eigen(sm))
 78 
eigenok(sm, eigen(sm))
73 
eigenok(sm, eigen(sm, sym=FALSE))
 79 
eigenok(sm, eigen(sm, sym=FALSE))