Project EulerのProblem 66を解く際に平方根の連分数を求める関数を書きました。これについてちょっとメモを。 それと当然ながらProject EulerのProblem 64〜66あたりのネタバレになっちゃいます。注意して下さい。 ## 連分数の循環節を掃きだす rep.frac <- function(n){ if(sqrt(n)%%1==0) return(0) om0 <- om <- sqrt(n) rep.set <- numeric(0) ## a_0 a <- rep.set[1] <- floor(om) p <- n-a^2 l <- a om <- (om0 + l)/p ## 循環節のひとつめ(ここに戻ったら終了) a <- rep.set[2] <- floor(om) l1 <- l <- p*a - l p1 <- p <- (n-l^2)/p