〜2016年1月中旬〜
13個連続する数字文字列の最大積。この問題も2008年に解いたときは5個連続だった。
サイズ13の窓を1つずつずらすという感じになる。Rubyならeach_cons(13)だな。
% echo 73..50 | #長いので省略
grep -o . |
awk '{m=1; for (i=1; i<=12; i++){ m *= a[i] = a[i+1] }}$0 = m *= a[i] = $0' |
sort -n |
tail -1
23514624000
grepで1000文字1行を1文字1000行に変換し、a[i]へa[13]から順に入れてずらしていく。
0が含まれていれば$0は0になるので出力されない。
最初の12行は初期化してないa[i]が入っているので$0は0になり出力されない。
まあ無駄に計算するが気にしない。
awkでは初期化してない変数は、文字列として扱えば空文字、数値として扱えば0となる。
Rubyだとこんな感じだ。
% ruby -e 'p "73..50".chars.map(&:to_i).each_cons(13).map{|a|a.reduce(:*)}.max'
23514624000
cf: ProjectEuler - Project Eulerをシェル芸で解いてみる(Problem 8) - Qiita
行番号を13で割った余りをインデックスにすれば配列をシフトする必要がないのであった。
% echo 73..50 | #長いので省略
grep -o . |
awk '{m=1; a[NR%13] = $0; for (i=0; i<13; i++){ m *= a[i] }}$0 = m' |
sort -n |
tail -1
23514624000
元のawkの部分。
awk '{m=1; for (i=1; i<=12; i++){ m *= a[i] = a[i+1] }}$0 = m *= a[i] = $0' |
ピタゴラス数。ひねりなし。awkやbcは出てるのでbashで。
% cat p9.sh
for a in {1..332}; do
for ((b=a+1; b<500; b++)); do
((c=1000-a-b))
if ((a*a+b*b == c*c)); then
echo $[a*b*c]
break 2
fi
done
done
% time bash p9.sh
31875000
bash p9.sh 0.92s user 0.00s system 99% cpu 0.926 total
ちょっと高速化を考える。
c=1000-a-bなのでc*cの部分に代入すると両辺のa*aとb*bが相殺される。
最終的には
((c=1000-a-b))
if ((a*a+b*b == c*c)); then
echo $[a*b*c]
は
if ((1000*(a+b)-a*b == 500000)); then
echo $[a*b*(1000-a-b)]
でいい。こうすると計算量が減るのでちょっと速くなる。
% time bash p9.sh 31875000 bash p9.sh 0.68s user 0.00s system 99% cpu 0.685 total
変形してやるときれいな式になる。さらにちょっとだけ速くなる。
if (( (1000-a)*(1000-b) == 500000 )); then
cf: ProjectEuler - Project Eulerをシェル芸で解いてみる(Problem 9) - Qiita
素数の和。200万と数がでかいだけに時間がかかる。
ひねりなし。
% time seq 2000000 | factor | awk '!$3{s+=$2}END{print s}'
142913828922
seq 2000000 0.11s user 0.05s system 17% cpu 0.907 total
factor 0.86s user 0.03s system 98% cpu 0.910 total
awk '!$3{s+=$2}END{print s}' 0.76s user 0.00s system 83% cpu 0.910 total
Ruby版も。
% time ruby -rprime -e 'p Prime.each(2e6).reduce(:+)' 142913828922 ruby -rprime -e 'p Prime.each(2e6).reduce(:+)' 0.89s user 0.01s system 90% cpu 0.992 total
prime.rbはpure rubyなライブラリな割りには速い。
cf: ProjectEuler - Project Eulerをシェル芸で解いてみる(Problem 10) - Qiita
Lucy_Hedgehogという人が書いたPython版が興味深い。とてつもなく速い。gawk(64bit版)に移植した。
mawkだと表示が1.42914e+11のように丸められてしまう。
function P10(n) {
r = int(sqrt(n))
for (i=1; i<=r; i++) {
V[i] = v = int(n/i)
S[v] = v*(v+1)/2-1
}
for (j=V[i-1]-1; j > 0 ; j--) {
V[i++] = j
S[j] = j*(j+1)/2-1
}
l=i-1
for (p=2; p <= r; p++) {
if (S[p] > sp = S[p-1]) {
p2 = p*p
for (i=1; i<=l; i++) {
v = V[i]
if (v < p2)
break
S[v] -= p*(S[int(v/p)]-sp)
}
}
}
return S[n]
}
BEGIN {
print P10(2000000)
}
v*(v+1)/2-1で(Σv)-1(1は素数じゃない)がわかるので、
そこから素数以外を引けば素数の和が求まるという寸法だ。
このときif文が成立するときはpは素数なので同時に素数も求まっているが、
p<=rつまりsqrt(2000000)以下なので、ちょっと使いにくい。
なぜかsが消えていた。chmod +sすれば動くがcapabilityで対処してみた。
% cp =ping . % ./ping 127.0.0.1 ping: icmp open socket: Operation not permitted % getcap ./ping % sudo setcap cap_net_raw=ep ./ping % getcap ./ping ./ping = cap_net_raw+ep % ./ping 127.0.0.1 PING 127.0.0.1 (127.0.0.1) 56(84) bytes of data. 64 bytes from 127.0.0.1: icmp_seq=1 ttl=64 time=0.071 ms 64 bytes from 127.0.0.1: icmp_seq=2 ttl=64 time=0.040 ms ^C --- 127.0.0.1 ping statistics --- 2 packets transmitted, 2 received, 0% packet loss, time 999ms rtt min/avg/max/mdev = 0.040/0.055/0.071/0.017 ms
てな感じでパンピーでも実行できる。
Largest product in a grid.
縦横斜めと同じような処理をどうするかだが、
2次元じゃなく1次元にしてしまうと実は楽になる。
基点からオフセットを考えると
横: [0, 1, 2, 3] 縦: [0,20,40,60] ／: [0,19,38,57] ＼: [0,21,42,63]
となる。つまり[1,20,19,21]*[0,1,2,3]的な感じになる。
Rubyなら
[1, 20, 19, 21].map{|x|4.times.map{|y|x*y}}
だ。ただ、1次元だと右側を超えて左側につながっていることになるので、
0を各行の先頭に入れて掛ければ0になるようにする。
それを考慮すると[1, 21, 20, 22]になる。
% cat p11.sh
echo '
08 02 22 97 38 15 00 40 00 75 04 05 07 78 52 12 50 77 91 08
49 49 99 40 17 81 18 57 60 87 17 40 98 43 69 48 04 56 62 00
81 49 31 73 55 79 14 29 93 71 40 67 53 88 30 03 49 13 36 65
52 70 95 23 04 60 11 42 69 24 68 56 01 32 56 71 37 02 36 91
22 31 16 71 51 67 63 89 41 92 36 54 22 40 40 28 66 33 13 80
24 47 32 60 99 03 45 02 44 75 33 53 78 36 84 20 35 17 12 50
32 98 81 28 64 23 67 10 26 38 40 67 59 54 70 66 18 38 64 70
67 26 20 68 02 62 12 20 95 63 94 39 63 08 40 91 66 49 94 21
24 55 58 05 66 73 99 26 97 17 78 78 96 83 14 88 34 89 63 72
21 36 23 09 75 00 76 44 20 45 35 14 00 61 33 97 34 31 33 95
78 17 53 28 22 75 31 67 15 94 03 80 04 62 16 14 09 53 56 92
16 39 05 42 96 35 31 47 55 58 88 24 00 17 54 24 36 29 85 57
86 56 00 48 35 71 89 07 05 44 44 37 44 60 21 58 51 54 17 58
19 80 81 68 05 94 47 69 28 73 92 13 86 52 17 77 04 89 55 40
04 52 08 83 97 35 99 16 07 97 57 32 16 26 26 79 33 27 98 66
88 36 68 87 57 62 20 72 03 46 33 67 46 55 12 32 63 93 53 69
04 42 16 73 38 25 39 11 24 94 72 18 08 46 29 32 40 62 76 36
20 69 36 41 72 30 23 88 34 62 99 69 82 67 59 85 74 04 36 16
20 73 35 29 78 31 90 01 74 31 49 71 48 86 81 16 23 57 05 54
01 70 54 71 83 51 54 69 16 92 33 48 61 43 52 01 89 19 67 48' |
sed '1d;s/^/0 /' |
awk '
{ for (i = 1; i <= NF; i++) data[++n] = $i }
END {
b[0] = 1; b[1] = 21; b[2] = 20; b[3] = 22
for (i = 1; i <= n; i++) {
for (j = 0; j < 4; j++) {
m = 1
for (k = 0; k < 4; k++) m *= data[i+b[j]*k]
if (max < m) max = m
}
}
print max
}
'
Ruby版。
data = <<EOM.gsub(/^/,"0 ").split.map{|i|i.to_i(10)}
08 02 22 97 38 15 00 40 00 75 04 05 07 78 52 12 50 77 91 08
49 49 99 40 17 81 18 57 60 87 17 40 98 43 69 48 04 56 62 00
81 49 31 73 55 79 14 29 93 71 40 67 53 88 30 03 49 13 36 65
52 70 95 23 04 60 11 42 69 24 68 56 01 32 56 71 37 02 36 91
22 31 16 71 51 67 63 89 41 92 36 54 22 40 40 28 66 33 13 80
24 47 32 60 99 03 45 02 44 75 33 53 78 36 84 20 35 17 12 50
32 98 81 28 64 23 67 10 26 38 40 67 59 54 70 66 18 38 64 70
67 26 20 68 02 62 12 20 95 63 94 39 63 08 40 91 66 49 94 21
24 55 58 05 66 73 99 26 97 17 78 78 96 83 14 88 34 89 63 72
21 36 23 09 75 00 76 44 20 45 35 14 00 61 33 97 34 31 33 95
78 17 53 28 22 75 31 67 15 94 03 80 04 62 16 14 09 53 56 92
16 39 05 42 96 35 31 47 55 58 88 24 00 17 54 24 36 29 85 57
86 56 00 48 35 71 89 07 05 44 44 37 44 60 21 58 51 54 17 58
19 80 81 68 05 94 47 69 28 73 92 13 86 52 17 77 04 89 55 40
04 52 08 83 97 35 99 16 07 97 57 32 16 26 26 79 33 27 98 66
88 36 68 87 57 62 20 72 03 46 33 67 46 55 12 32 63 93 53 69
04 42 16 73 38 25 39 11 24 94 72 18 08 46 29 32 40 62 76 36
20 69 36 41 72 30 23 88 34 62 99 69 82 67 59 85 74 04 36 16
20 73 35 29 78 31 90 01 74 31 49 71 48 86 81 16 23 57 05 54
01 70 54 71 83 51 54 69 16 92 33 48 61 43 52 01 89 19 67 48
EOM
d = [1, 21, 20, 22].map{|x|4.times.map{|y|x*y}}
p data.size.times.map{|i|
d.map{|a|
a.map{|j|data[i+j]||0}.reduce(:*)
}.max
}.max
cf: ProjectEuler - Project Eulerをシェル芸で解いてみる(Problem 11) - Qiita
Highly divisible triangular number
最初に約数の数が500を超える三角数。
素因数分解が肝なのでfactorの出力をなんとかするわけだが、
すでにソートされてるのでuniq -cすれば個数はわかる。
% factor 10000
10000: 2 2 2 2 5 5 5 5
% factor 10000 | tr ' ' '\n'
10000:
2
2
2
2
5
5
5
5
% factor 10000 | tr ' ' '\n' | uniq -c
1 10000:
4 2
4 5
10000なら2**4*5**4とわかる。
うまい具合に:つきで三角数もわかるのでこれを区切りにできる。
% factor 10000 100000 | tr ' ' '\n' | uniq -c
1 10000:
4 2
4 5
1 100000:
5 2
5 5
というわけであとはawkだ。
% yes|awk '$0=NR*++NR/2'|factor|tr ' ' '\n'|uniq -c|awk '/:/{if(m>500){print +t;exit};m=1;t=$2;next}{m*=$1+1}'
cf: ProjectEuler - Project Eulerをシェル芸で解いてみる(Problem 12) - Qiita
Large sum.
50桁の数値100個の和の上10桁。
全然ひねりなし。
% paste -sd+ p13.txt | bc | cut -c-10
64ビットgawkでも可能。
% gawk '{s+=$0}END{print substr(s,1,10)}' p13.txt
mawkは32ビットなので残念な結果に。
% mawk '{s+=$0}END{print substr(s,1,10)}' p13.txt
5.53738e+5
しかし上10桁ということなら何とかなる。
cut -c-11 p13.txtを利用すれば。
まあ、でも長くなって面白くないので省略。
最後にRuby版。
% ruby -0 -ane 'printf "%.10s\n", eval($F*"+")' p13.txt
cf: ProjectEuler - Project Eulerをシェル芸で解いてみる(Problem 13) - Qiita
今日は趣向を変えて、
「ProjectEuler - Project Eulerをシェル芸で解いてみる(Problem 14) - Qiita」
を高速化してみる。
% /usr/bin/time -f %E seq 1 999999 |gawk 'ORS=""; {n=$0; while(n!=1){print n " "; n=(n%2 == 0) ? n/2 : 3*n + 1}; print n "\n"}'| gawk '{print $1, gensub(/[0-9]* /, "s", "g")}' | gawk '{print $1, length($2)}' |sort -k 2n,2 |tail -n 1 |cut -d ' ' -f 1
3:01.54
837799
うちのPCだと3分かかる。まずはprintが呼ばれる回数を減らす。
% /usr/bin/time -f %e seq 1 999999 |gawk 'ORS=""; {n=$0; s=""; while(n!=1){s=s n " "; n=(n%2 == 0) ? n/2 : 3*n + 1}; print s n "\n"}'| gawk '{print $1, gensub(/[0-9]* /, "s", "g")}' | gawk '{print $1, length($2)}' |sort -k 2n,2 |tail -n 1 |cut -d ' ' -f 1
1:57.60
837799
結構縮む。2分を切った。
そもそも最初のgawkの出力が300MBを超えているのが遅い原因の1つ。
あとで数値を"s"に変換してるのでnじゃなくて"s"でもいい。
これでかなりのサイズダウンが期待できる。
% /usr/bin/time -f %E seq 1 999999 |gawk 'ORS=""; {n=$0; s=n" "; while(n!=1){n=(n%2 == 0) ? n/2 : 3*n + 1;s=s "s"}; print s "1\n"}'| gawk '{print $1, length($2)}' |sort -k 2n,2 |tail -n 1 |cut -d ' ' -f 1
0:26.26
837799
30秒を切った。sは文字列じゃなくて長さを求めればいい。あとすでにORSは不要なので整理。
% /usr/bin/time -f %E seq 1 999999 |gawk '{n=$0; s=1; while(n!=1){n=(n%2 == 0) ? n/2 : 3*n + 1; s++}; print $0,s}'|sort -k 2n,2 |tail -n 1 |cut -d ' ' -f 1
0:23.40
837799
ここで3*n+1に注目する。このときのnは奇数なので結果は必ず偶数になる。
だったらそれを次回に回さずにこの場で2で割ってしまえばいいわけだ。
% /usr/bin/time -f %E seq 1 999999 |gawk '{n=$0; s=1; while(n!=1){if(n%2){n=3*n+1;s++}; n/=2; s++}; print $0,s}'|sort -k 2n,2 |tail -n 1 |cut -d ' ' -f 1
0:18.32
837799
さらなる高速化はまた明日。
昨日の最適化で18秒になった。ここでキャッシュを使って高速化してみる。
いくつか方法があるが、まずは問題文にもある13を例に考える。
13 → 40 → 20 → 10 → 5 → 16 → 8 → 4 → 2 → 1
ここから言えるのは40,20,10,5,16,8,4,2,1は13よりも短いシーケンスだということ。
だからもう計算する必要ない。これをキャッシュさせて覚えておく。
効率よくやるには長そうなものから始めるべきで、つまり逆順にするのがいい。
% /usr/bin/time -f %E seq 999999 -1 1 |gawk '!h[$0]{n=$0; s=1; while(n!=1){if(n%2){n=3*n+1;s++}; n/=2; s++; h[n]=1}; print $0,s}'|sort -k 2n,2 |tail -n 1 |cut -d ' ' -f 1
0:10.70
837799
ブロックの先頭に!h[$0]を追加、nを計算するごとにh[n]=1を追加しただけこれだけ速くなった。
もう一つの考え方として、長さ自体をキャッシュする方法。
この場合は1から順番に長さをh[$0]へ入れる。
while loopの終了条件として$0<=nが使える。これがかなり効く。
13の場合なら13未満の10が出てきた時点でh[10]が使えるので、そこでloopをやめられる。
% /usr/bin/time -f %E seq 1 999999 |gawk '{h[1]=1; n=$0; s=1; while(n!=1&&$0<=n){if(n%2){n=3*n+1;s++}; n/=2; s++}; h[$0]=s+h[n]; print $0,h[$0]}'|sort -k 2n,2 |tail -n 1 |cut -d ' ' -f 1
0:03.90
837799
4秒を切った。
cf: ProjectEuler - Project Eulerをシェル芸で解いてみる(Problem 14) - Qiita