Project Eulerをシェル芸で解いてみる(Problem 27) #シェル芸

概要

Project Euler¹は、数学的な問題解決能力を要求する、プログラミング問題集・オンラインプログラミングジャッジである。
一般的なオンラインプログラミングジャッジとは異なり、解答するプログラミング言語は、特に限定されていない。
フォームに入力した値だけで正誤を判別するため、適切なアルゴリズムさえ考案できれば、あらゆるプログラミング言語で解答することが可能である。
そのため、Web上では様々な手法による解答が公開されている。

しかし、シェル芸²による解答は、私が探した限りでは、 @eban さんの記事しか見つからなかった。
そこで、私自身の勉強も兼ねて、Project Eulerをシェル芸で解いてみることにした。

本記事では、Project EulerのProblem 27を、シェル芸で解いた際の過程を述べる。

実行環境

• Arch Linux (x86_64, 4.5.1-1-ARCH)
• bash (4.3.42)
• GNU coreutils (8.25)
• gawk (4.1.3)
• bsd-games (2.17)

問題文

[原文 (https://projecteuler.net/problem=27)]

Quadratic primes

Euler discovered the remarkable quadratic formula:

n^2 + n + 41

It turns out that the formula will produce $40$ primes for the consecutive values $n = 0$ to $39$. However, when $n = 40$, $40^2 + 40 + 41 = 40(40 + 1) + 41$ is divisible by $41$, and certainly when $n = 41$, $41^2 + 41 + 41$ is clearly divisible by $41$.

The incredible formula $n^2 − 79n + 1601$ was discovered, which produces $80$ primes for the consecutive values $n = 0$ to $79$. The product of the coefficients, $−79$ and $1601$, is $−126479$.

Considering quadratics of the form:

$\hspace{3em}$$n^2 + an + b$, where $|a| \lt 1000$ and $|b| \lt 1000$

$\hspace{3em}$where $|n|$ is the modulus/absolute value of $n$
$\hspace{3em}$e.g. $|11| = 11$ and $|−4| = 4$

Find the product of the coefficients, $a$ and $b$, for the quadratic expression that produces the maximum number of primes for consecutive values of $n$, starting with $n = 0$.

[日本語訳 (http://odz.sakura.ne.jp/projecteuler/index.php?cmd=read&page=Problem%2027)]

二次式素数

オイラーは以下の二次式を考案している:

n2 + n + 41

この式は, $n$ を$0$から$39$までの連続する整数としたときに$40$個の素数を生成する. しかし, $n = 40$ のとき $40^2 + 40 + 41 = 40(40 + 1) + 41$ となり$41$で割り切れる. また, $n = 41$ のときは $41^2 + 41 + 41$ であり明らかに$41$で割り切れる.

計算機を用いて, 二次式 $n^2 - 79n + 1601$ という式が発見できた. これは $n = 0$ から $79$ の連続する整数で$80$個の素数を生成する. 係数の積は, $-79 \times 1601$ で $-126479$である.

さて, $|a| \lt 1000$, $|b| \lt 1000$ として以下の二次式を考える (ここで $|a|$ は絶対値): 例えば $|11| = 11$， $|-4| = 4$である.

n^2 + an + b

$n = 0$ から始めて連続する整数で素数を生成したときに最長の長さとなる上の二次式の, 係数 $a$, $b$ の積を答えよ.

正答は、注釈³に記載する。

解法

1	#!/bin/bash
2	
3	cat <( \
4	        primes 2 99999 | xargs
5	        # 素数表の生成(2~99999まで)
6	    ) \
7	    <( \
8	        primes 2 999 | awk '{for(a=-$1+1; a<1000; a++){print a,$1}}'
9	        # a,bの候補値の生成
10	    ) \
11	| # 素数表およびa,bの候補値の出力
12	awk 'NR==1{for(i=1;i<=NF;i++){prime[$i]=$i}; next} {ORS=""; print $0" "; n=0; while(1){eq=n^2+$1*n+$2; if(prime[eq]==0){print "\n"; next}; print eq" "; n++}}' \
13	| # n=0から始めて、素数を生成できなくなるまで末尾フィールドに出力
14	awk '$0=NF" "$0' \
15	| # フィールド数を第一フィールドとして付加
16	sort -k 1,1nr \
17	| # 第一フィールド(フィールド数)をキーとして、逆順数値ソート
18	head -n 1 \
19	| # 生成する素数の数が最大となる組み合わせの抽出
20	awk '$0=$2*$3'
21	  # 積a*bの算出

以下、コマンドライン上での実行用コードである。

cat <(primes 2 99999 | xargs) <(primes 2 999 | awk '{for(a=-$1+1; a<1000; a++){print a,$1}}') | awk 'NR==1{for(i=1; i<=NF; i++){prime[$i]=$i}; next} {ORS=""; print $0" "; n=0; while(1){eq=n^2+$1*n+$2; if(prime[eq]==0){print "\n"; next}; print eq" "; n++}}' | awk '$0=NF" "$0' | sort -k 1,1nr | head -n 1 | awk '$0=$2*$3'

解説

まず、解答方針を述べる。

今回の問題は、以下の手順で解けば良い。

1. $a$と$b$の候補を列挙する。
2. $a$, $b$の候補に対して、素数が生成されなくなるまで、$n=0$から連続する整数を順に$n^2+an+b$へ代入する。
3. 生成された素数の数が最大になる$a$, $b$の組み合わせを抽出する。
4. 抽出した$a$と$b$の積を算出する。

ただし、単純に$a$と$b$の全ての組み合わせを列挙すると、探索空間が非常に広くなってしまう。
そのため、$a$と$b$の候補を絞り、探索数を減らしておく。

まず、$n=0$の場合を考える。
式に当てはめると、$0^2+a*0+b = b$となる。
ここで、式$n^2+an+b$は素数となるため、$b$は自ずと素数であることがわかる。

次に、$n=1$の場合を考える。
式に当てはめると、$1^2+a*1+b = 1+a+b$となる。
ここで、先ほどと同様に、式 $n^2+an+b$ は素数となる。
$1$以下の素数は存在しないため、先ほど求めた$1+a+b$は$1$以上、$1+a+b \gt 1 \Rightarrow a \gt -b$であることがわかる。

以上の操作により、$b$は$1000$未満の素数、$a$は$a \gt -b$であることが判明した。
この限定された範囲の$a$, $b$の候補を列挙し、設問に該当する組み合わせを探せば良い。

---

次に、ソースコードの解説を述べる。

まず、手順1にしたがい、$a$, $b$の候補を列挙する。
このとき、後の素数判定に用いるために、1レコード目に素数表を出力しておく。
これらの処理は、3~10行目において、catとプロセス置換を用いておこなう。
今回の解答では、コードを簡潔にするために、素数の生成にbsd-gamesのprimesコマンドを利用した。

次に、手順2にしたがい、素数が生成されなくなるまで、各レコードに対して$n$を$0$から順に代入していく。
この処理は、12行目のawkを用いておこなう。
ここで、$n$に整数を代入した結果が素数であるかを判定するために、1レコード目に出力した素数表を利用する。
1レコード目には、$2 \sim 99999$の範囲にある素数が、各フィールドに格納されている。
そこで、各フィールドの値を添字として、連想配列にフィールド値を格納しておく。
その後、2レコード目以降において、1,2フィールド目の$a$, $b$の値を式に代入し、素数が生成されなくなるまで$n$に$0$から順に整数を代入していく。
この素数の判定処理では、1レコード目の処理で作成した連想配列の値を利用する。
以上の処理によって、各レコードの出力は以下のようになる。

# a値 b値 素数(n=0) 素数(n=1) 素数(n=2) ...
-1 2 2 2
0 2 2 3
1 2 2
2 2 2 5
3 2 2
4 2 2 7
︙

そして、手順3にしたがい、生成された素数の数が最大となる$a$と$b$の組み合わせを抽出する。
この処理は、14~18のawk、sort、headを用いておこなう。
まず、最大となる組を抽出する下準備として、各レコードの素数の数を求める。
この操作は、単純にawkの特殊変数NFを利用すればよい。
その後、出力されたレコード数をキーとして逆順数値ソートし、最大値となるもの(1レコード目)を抽出すればよい。

最後に、手順4にしたがい、$a$, $b$の積を算出する。
この処理は、20行目のawkを用いておこなう。

以上が、$n = 0$ から始めて連続する整数で素数を生成したとき、最長の長さとなる上の二次式の係数$a$, $b$の積を求める方法である。

参考にしたもの

1. "Project Euler Problem 27 #シェル芸 - jarp," http://jarp.does.notwork.org/diary/201602a.html#201602091

雑記

• 「こんな別解があるよ」や、「こうすればもっと効率化できるぜ」、「こう書けばエレガントですわよ」等の意見がありましたら、気軽に教えていただけますと幸いです🐺
• 現在の進捗状況 -> https://github.com/gin135/Project_Euler/tree/master/sh

• @eban さんが、ご自身の日記にて別解を投稿していらっしゃいます。
  • "Project Euler Problem 27 #シェル芸 - jarp," http://jarp.does.notwork.org/diary/201602a.html#201602091
  • 本投稿はこの日記を大いに参考にしたため、素数判定の処理以外は、ほぼ同じ考え方になっています。
  • また、 @eban さんの解答では、awkのみで実装されています。

• 今回の解法は、ほぼ @eban さんのやり方をパクったお借りしたような感じになってしまったorz

  • いちおう、素数判定の方法を変えることで、オリジナリティは出した...つもり。
  • ただ、素数表の出力範囲が決め打ちになっている点が、ちょっと残念...

• bashのプロセス置換、シェル芸では無茶苦茶便利なのだけれど...

  • 俺々シェル芸書式との相性が良くないのが難点。すっごく見辛い...

---

1. About - Project Euler (https://projecteuler.net) ↩

2. シェル芸の定義 >> シェル芸 - 上田ブログ (https://blog.ueda.asia/?page_id=1434) ↩

3. -59231 ↩