この下の問題だけど，よく解析してからプログラムを書くと，あっと驚く速さになる。
まずは，探索範囲を狭めること。
> system.time({
+ euclid = function(m, n) {
+ while ((temp = n%%m) != 0) {
+ n = m
+ m = temp
+ }
+ return(m)
+ }
+ N = 1000
+ mx = floor(sqrt((1 + sqrt(2)) * N/2))
+ d = expand.grid(1:mx, 1:mx)
+ d = d[d[, 1] > d[, 2], ]
+ n = d[, 1]
+ m = d[, 2]
+ a = n^2 - m^2
+ b = 2 * n * m
+ c = n^2 + m^2
+ d = cbind(ifelse(a < b, a, b), ifelse(a > b, a, b), c)
+ d = d[d[, 2] < N, ]
+ GCD = apply(d, 1, function(x) Reduce(euclid, x))
+ d = d/GCD
+ d = unique(d)
+ })
ユーザ システム 経過
0.027 0.002 0.036
