Purpose
To code the 4 conditional numbers relevant to solving Ax = b

> A <- matrix(c(1, 1, 1, 1.0001, 1, 1.0001), byrow = T, nrow = 3,
+     ncol = 2)
> b <- matrix(c(2, 1e-04, 4.0001), nrow = 3, ncol = 1)
> x <- solve(t(A) %*% A) %*% t(A) %*% b
> A_plus <- solve(t(A) %*% A) %*% t(A)
> P <- A %*% solve(t(A) %*% A) %*% t(A)
> Pb <- P %*% b
> x <- A_plus %*% b
> y <- Pb

x and y are the exact solutions for Ax ~= b

k(A)

> kappa_A <- max(svd(A)$d)/min(svd(A)$d)

theta

> norm2 <- function(x) sqrt(sum(x^2))
> normA <- function(x) return(max(svd(x)$d)/min(svd(x)$d))
> theta <- acos(norm2(y)/norm2(b))

eta

> eta <- normA(A) * norm2(x)/norm2(A %*% x)

The parameters are

> print(kappa_A)
[1] 42429.24
> print(theta)
[1] 0.6847028
> print(eta)
[1] 17321.09

Perturbation table

> pert <- diag(2)
> pert[1, 1] <- 1/cos(theta)
> pert[1, 2] <- kappa_A/(eta * cos(theta))
> pert[2, 1] <- kappa_A/cos(theta)
> pert[2, 2] <- kappa_A + (kappa_A)^2 * tan(theta)/eta
> print(round(pert, 2))
         [,1]      [,2]
[1,]     1.29      3.16
[2,] 54775.18 127287.70

As can be seen that changes in A can produce large perturbations in y and x Hence it is ill conditioned.