https://github.com/cran/pracma
Tip revision: 0e3ae6b4108dcce6052c1e366d66aa83f8502596 authored by HwB on 10 January 2013, 00:00:00 UTC
version 1.3.8
version 1.3.8
Tip revision: 0e3ae6b
agmean.R
##
## a g m e a n . R Arithmetic-geometric Mean
##
agmean <- function(a, b) {
eps <- .Machine$double.eps
stopifnot(is.numeric(a) || is.complex(a),
is.numeric(b) || is.complex(b))
if (is.numeric(a) && any(a < 0) || is.numeric(b) && any(b < 0)) {
a <- as.complex(a)
b <- as.complex(b)
}
if (length(a) == 1) {
n <- length(b)
a <- rep(a, n)
} else if (length(b) == 1) {
n <- length(a)
b <- rep(b, n)
} else
if (length(a) != length(b))
stop("Arguments must have the same length or one has length 1.")
while ( any(abs(a-b) >= eps) ) {
a1 <- (a + b) / 2
b1 <- sqrt(a * b)
if (max(abs(a-a1)) < eps && max(abs(b-b1)) < eps) break
a <- a1
b <- b1
}
return( (a+b)/2 )
}