# linear congruential generator (lcg)
lcg_prng <- function(n, y1, m, a, b) {
    if (n == 1) return(y1)
    y <- (a*y1 + b) %% m
    c(y1, lcg_prng(n-1, y, m, a, b))
}
par(mfrow=c(1,2))
m <- 2^11
x <- lcg_prng(500, y1=1, m=m, a=1017, b=1)
plot(x[1:499]/m, x[2:500]/m, pch=".", cex=4, xlab="u1, ..., u499", ylab="u2, ..., u500", main="lcg_prng()")
# Vergleich mit R-Funktion sample.int(), siehe ?sample
x <- sample.int(m, 500)
plot(x[1:499]/m, x[2:500]/m, pch=".", cex=4, xlab="u1, ..., u499", ylab="u2, ..., u500", main="sample.int()")