Project Euler 6~10をpythonで解く

Problem 6

問題（英語）

The sum of the squares of the first ten natural numbers is,
1^2+2^2+...+10^2=385
The square of the sum of the first ten natural numbers is,
(1+2+...+10)^2=552=3025
Hence the difference between the sum of the squares of the first ten natural numbers and the square of the sum is 3025−385=2640.

Find the difference between the sum of the squares of the first one hundred natural numbers and the square of the sum.
「 https://projecteuler.net/problem=6 」

問題（日本語）

最初の10個の自然数について, その二乗の和は,

1^2 + 2^2 + ... + 10^2 = 385

最初の10個の自然数について, その和の二乗は,

(1 + 2 + ... + 10)^2 = 3025

これらの数の差は 3025 - 385 = 2640 となる.

同様にして, 最初の100個の自然数について二乗の和と和の二乗の差を求めよ.
「 http://odz.sakura.ne.jp/projecteuler/index.php?cmd=read&page=Problem%206 」

コード

#aが累乗の和、bが和の累乗 
a = 0
b = 0
for i in range(100):
a += (i + 1) * (i + 1)
b += (i + 1)
b = b * b
print(b - a)

#aが累乗の和、bが和の累乗

a = 0

b = 0

for i in range(100):

a += (i + 1) * (i + 1)

b += (i + 1)

b = b * b

print(b - a)

考え方

そのまま。

どことなくコードからもやる気のなさを感じる。

これは問題とは関係ないが、tex関連の記載をされるとコピーが上手くできなくて困るw

最初、1^2 + 2^2 + ... + 10^2 = 385 が

12 + 22 + ... + 102 = 385 となっていて直すのが少々面倒だった(^^;)

Problem 7

問題（英語）

By listing the first six prime numbers: 2, 3, 5, 7, 11, and 13, we can see that the 6th prime is 13.
What is the 10 001st prime number?
「 https://projecteuler.net/problem=7 」

問題（日本語）

素数を小さい方から6つ並べると 2, 3, 5, 7, 11, 13 であり, 6番目の素数は 13 である.

10 001 番目の素数を求めよ.
「 http://odz.sakura.ne.jp/projecteuler/index.php?cmd=read&page=Problem%207 」

コード

#10001番目の素数
import math
n = 3
#counterは素数の数を記録するもの。最初から「2」を数えている。
counter = 1
while counter < 10001:
    #nが素数かどうか。ループの度にTrueに初期化。
    n_is_prime = True
    sqrt_n = int(math.sqrt(n))
    for i in range(sqrt_n):
        if n % (i+2) == 0:
            n_is_prime = False
            break
    if n_is_prime:
        counter += 1
    n += 1

print(n-1)

#10001番目の素数

import math

n = 3

#counterは素数の数を記録するもの。最初から「2」を数えている。

counter = 1

while counter < 10001:

#nが素数かどうか。ループの度にTrueに初期化。

n_is_prime = True

sqrt_n = int(math.sqrt(n))

for i in range(sqrt_n):

if n % (i+2) == 0:

n_is_prime = False

break

if n_is_prime:

counter += 1

n += 1

print(n-1)

考え方

問題はザ・シンプル。

前回素因数は出てきだが、素数が出てくるのはこの問題が初めて。

全体的な流れは自然数を一つずつ素数かどうかチェックしていき、countが10001になったら結果を出すというもの。

問題は素数の判定である。

素数をどう判定するか

これはいろんなサイトに書いてあることでこのコードを書く際にもいくつかのサイトを参考にしたと思う。

基本的には判定したい数を自然数で割っていきその余りが0になるかを調べていく。

一度でも0になれば何かで割り切れることを意味するため、その数は素数ではない。

肝心なのはどう高速化するかである。

どう高速化するか

今、自然数nが素数か調べたいとする。

そもそも素数か判定するには、

割る数は約√nまででいい。
また、割る数は素数だけで十分である。

これが高速化の鍵である。

今回のコードでは、箇条書きの一つ目だけ以下の部分で採用した。

sqrt_n = int(math.sqrt(n)) 
for i in range(sqrt_n):

1 2	sqrt_n = int(math.sqrt(n)) for i in range(sqrt_n):

自分の関数作っちゃうのがいいよね

project eulerには素数問題が非常にたくさんあるので、よく使う関数は自分で作ってしまうのが手っ取り早い。

僕の場合は、問題の解答のコードを保存するファイルと別に、自作の関数のコード保存するファイルを作って管理している。

僕が現在採用してるのは、以下の素数判定の関数である。

#参考「https://oku.edu.mie-u.ac.jp/~okumura/python/prime.html」
import math

def is_prime(n):
    if n < 2:
        return False
    if n == 2:
        return True
    if n % 2 == 0:
        return False
    m = math.floor(math.sqrt(n)) + 1
    for p in range(3, m, 2):
        if n % p == 0:
            return False
    return True

#参考「https://oku.edu.mie-u.ac.jp/~okumura/python/prime.html」

import math

def is_prime(n):

if n < 2:

return False

if n == 2:

return True

if n % 2 == 0:

return False

m = math.floor(math.sqrt(n)) + 1

for p in range(3, m, 2):

if n % p == 0:

return False

return True

参考「https://oku.edu.mie-u.ac.jp/~okumura/python/prime.html」

割る数から2の倍数を抜いている。（正確には割る数を2ずつ増やしている。）

割る数が半分になるので単純計算で2倍の早さになる。

2も素数なので割らないとダメだが、先に2以外の偶数は「False」と判定しているので大丈夫。

そもそも素数のリストを作ってしまえばいいのでは

と思う人がいるかもしれないが、その通りで「エラトネスのふるい」という素数をリスト化する方法がある。

この方法も非常に強力。

単純に素数判定
エラトネスのふるいで素数のリストを作り、そこから探す

どちらの方法を使うかはケースバイケースである。

初めての素数、初めての汎用的な関数ということで長くなりました。

Problem 8

問題（英語）

The four adjacent digits in the 1000-digit number that have the greatest product are 9 × 9 × 8 × 9 = 5832.

73167176531330624919225119674426574742355349194934
96983520312774506326239578318016984801869478851843
85861560789112949495459501737958331952853208805511
12540698747158523863050715693290963295227443043557
66896648950445244523161731856403098711121722383113
62229893423380308135336276614282806444486645238749
30358907296290491560440772390713810515859307960866
70172427121883998797908792274921901699720888093776
65727333001053367881220235421809751254540594752243
52584907711670556013604839586446706324415722155397
53697817977846174064955149290862569321978468622482
83972241375657056057490261407972968652414535100474
82166370484403199890008895243450658541227588666881
16427171479924442928230863465674813919123162824586
17866458359124566529476545682848912883142607690042
24219022671055626321111109370544217506941658960408
07198403850962455444362981230987879927244284909188
84580156166097919133875499200524063689912560717606
05886116467109405077541002256983155200055935729725
71636269561882670428252483600823257530420752963450

Find the thirteen adjacent digits in the 1000-digit number that have the greatest product. What is the value of this product?
「 https://projecteuler.net/problem=8 」

問題（日本語）

次の1000桁の数字のうち, 隣接する4つの数字の総乗の中で, 最大となる値は, 9 × 9 × 8 × 9 = 5832である.

73167176531330624919225119674426574742355349194934
96983520312774506326239578318016984801869478851843
85861560789112949495459501737958331952853208805511
12540698747158523863050715693290963295227443043557
66896648950445244523161731856403098711121722383113
62229893423380308135336276614282806444486645238749
30358907296290491560440772390713810515859307960866
70172427121883998797908792274921901699720888093776
65727333001053367881220235421809751254540594752243
52584907711670556013604839586446706324415722155397
53697817977846174064955149290862569321978468622482
83972241375657056057490261407972968652414535100474
82166370484403199890008895243450658541227588666881
16427171479924442928230863465674813919123162824586
17866458359124566529476545682848912883142607690042
24219022671055626321111109370544217506941658960408
07198403850962455444362981230987879927244284909188
84580156166097919133875499200524063689912560717606
05886116467109405077541002256983155200055935729725
71636269561882670428252483600823257530420752963450

この1000桁の数字から13個の連続する数字を取り出して, それらの総乗を計算する. では、それら総乗のうち、最大となる値はいくらか.

EX 6桁の数123789から5個の連続する数字を取り出す場合, 1*2*3*7*8と2*3*7*8*9の二通りとなり, 後者の2*3*7*8*9=3024が最大の総乗となる.
「 http://odz.sakura.ne.jp/projecteuler/index.php?cmd=read&page=Problem%208 」

コード

a = 7316717653133062491922511967442657474235534919493496983520312774506326239578318016984801869478851843858615607891129494954595017379583319528532088055111254069874715852386305071569329096329522744304355766896648950445244523161731856403098711121722383113622298934233803081353362766142828064444866452387493035890729629049156044077239071381051585930796086670172427121883998797908792274921901699720888093776657273330010533678812202354218097512545405947522435258490771167055601360483958644670632441572215539753697817977846174064955149290862569321978468622482839722413756570560574902614079729686524145351004748216637048440319989000889524345065854122758866688116427171479924442928230863465674813919123162824586178664583591245665294765456828489128831426076900422421902267105562632111110937054421750694165896040807198403850962455444362981230987879927244284909188845801561660979191338754992005240636899125607176060588611646710940507754100225698315520005593572972571636269561882670428252483600823257530420752963450
a_list = [int(c) for c in str(a)]
max_pro = 0
Digits = 13
for i in range(1000 - Digits + 1):
num = 1
for j in range(Digits):
num *= a_list[i+j]
max_pro = max(max_pro,num)
print(max_pro)

a = 7316717653133062491922511967442657474235534919493496983520312774506326239578318016984801869478851843858615607891129494954595017379583319528532088055111254069874715852386305071569329096329522744304355766896648950445244523161731856403098711121722383113622298934233803081353362766142828064444866452387493035890729629049156044077239071381051585930796086670172427121883998797908792274921901699720888093776657273330010533678812202354218097512545405947522435258490771167055601360483958644670632441572215539753697817977846174064955149290862569321978468622482839722413756570560574902614079729686524145351004748216637048440319989000889524345065854122758866688116427171479924442928230863465674813919123162824586178664583591245665294765456828489128831426076900422421902267105562632111110937054421750694165896040807198403850962455444362981230987879927244284909188845801561660979191338754992005240636899125607176060588611646710940507754100225698315520005593572972571636269561882670428252483600823257530420752963450

a_list = [int(c) for c in str(a)]

max_pro = 0

Digits = 13

for i in range(1000 - Digits + 1):

num = 1

for j in range(Digits):

num *= a_list[i+j]

max_pro = max(max_pro,num)

print(max_pro)

考え方

データの扱い方が試される問題。

改行されているデータをどう数字に戻すか

まず、1000桁の数字と言われているが、実際に与えられるのは何行かに改行された数字。

このデータを1000桁の数字に戻すのが最初の課題なのだが、やり方が分からず、メモ帳にペーストして改行を消していくという強引かつ面倒な解決策をとる。

このデータの扱い方を勉強したのはProblem13で、Problem11も似た問題だが、強引に解決した。

スマートな解決法は、まずデータをstr型として取得して、split関数を使ったりするとうまくいく。

気が向いたらこのコードも書き直してみる。（これは一生やらないやつ）

本題

上の桁から13連の組みの積を計算し、下の桁へとずらしていけばいいだけである。

と言っても、具体的にどう書くかは意外と難しかった。

定数（今回のではDigit）は頭文字だけでなく、全部の文字を大文字で書くのが一般的のよう。

Problem 9

問題（英語）

A Pythagorean triplet is a set of three natural numbers, a < b < c, for which,

a² + b² = c²

For example, 3² + 4² = 9 + 16 = 25 = 5².

There exists exactly one Pythagorean triplet for which a + b + c = 1000.
Find the product abc.
「 https://projecteuler.net/problem=9 」

問題（日本語）

ピタゴラス数(ピタゴラスの定理を満たす自然数)とは a < b < c で以下の式を満たす数の組である.

a^2 + b^2 = c^2

例えば, 3^2 + 4^2 = 9 + 16 = 25 = 5^2 である.

a + b + c = 1000 となるピタゴラスの三つ組が一つだけ存在する.
これらの積 abc を計算しなさい.
「 http://odz.sakura.ne.jp/projecteuler/index.php?cmd=read&page=Problem%209 」

コード

import math
end = 0
a = 0
while end == 0:
a += 1
b = 0
while end == 0 and b < 500:
b += 1
c = math.sqrt(a * a + b * b)
c_is_int = (c == int(c))
if c_is_int:
if a + b + c == 1000:
result = a * b * c
end = 1
print(result)
print('a:'+str(a)+' b:'+str(b)+' c:'+str(c))

import math

end = 0

a = 0

while end == 0:

a += 1

b = 0

while end == 0 and b < 500:

b += 1

c = math.sqrt(a * a + b * b)

c_is_int = (c == int(c))

if c_is_int:

if a + b + c == 1000:

result = a * b * c

end = 1

print(result)

print('a:'+str(a)+' b:'+str(b)+' c:'+str(c))

考え方

みんな大好きピタゴラス数。

c = √(a^2+b^2) として、cがint（整数）ならばピタゴラス数として判定する。

endは0,1のビット表現するより、True,Falseで表した方が分かりやすかったと反省。

Problem 10

問題（英語）

The sum of the primes below 10 is 2 + 3 + 5 + 7 = 17.

Find the sum of all the primes below two million.
「 https://projecteuler.net/problem=10 」

問題（日本語）

10以下の素数の和は 2 + 3 + 5 + 7 = 17 である.

200万以下の全ての素数の和を求めよ.
「 http://odz.sakura.ne.jp/projecteuler/index.php?cmd=read&page=Problem%2010 」

コード

import math

def isPrime(n):
  if n == 2:
    return True
  if n % 2 == 0:
    return False
  m = math.floor(math.sqrt(n)) + 1
  for p in range(3, m, 2):
      if n % p == 0:
        return False
  return True

num = 0
Top = 2000000
for i in range(Top):
    if isPrime(i):
        num += i
print(num)

import math

def isPrime(n):

if n == 2:

return True

if n % 2 == 0:

return False

m = math.floor(math.sqrt(n)) + 1

for p in range(3, m, 2):

if n % p == 0:

return False

return True

num = 0

Top = 2000000

for i in range(Top):

if isPrime(i):

num += i

print(num)

考え方

先程のProblem7の際に紹介した素数を判定する自作の関数isPrimeをこの時に作った。

mainでやるのは、1から200万までの数字を素数判定し、素数ならnumに足していくという作業。

これでは時間がかかりすぎるのでは、と先程時間を測ってみた。

プログラムの処理時間を計測する方法

これももう少し経ってから身につけたスキル。

先程のコードに計測用のコードを追加する。

import math

def isPrime(n):
  if n == 2:
    return True
  if n % 2 == 0:
    return False
  m = math.floor(math.sqrt(n)) + 1
  for p in range(3, m, 2):
      if n % p == 0:
        return False
  return True

import time 
t1 = time.time() #開始
num = 0
Top = 2000000
for i in range(Top):
    if isPrime(i):
        num += i
print(num)
t2 = time.time() #終了
print(t2-t1) #処理時間

import math

def isPrime(n):

if n == 2:

return True

if n % 2 == 0:

return False

m = math.floor(math.sqrt(n)) + 1

for p in range(3, m, 2):

if n % p == 0:

return False

return True

import time

t1 = time.time() #開始

num = 0

Top = 2000000

for i in range(Top):

if isPrime(i):

num += i

print(num)

t2 = time.time() #終了

print(t2-t1) #処理時間

このページを参考にしました。（「 https://qiita.com/fantm21/items/3dc7fbf4e935311488bc 」）

実行結果:

142913828923
7.814446926116943

1 2	142913828923 7.814446926116943

下が計測時間（単位は秒）。

約8秒と遅くも早くもない結果に。

エラトネスのふるいを使うと？

エラトネスのふるいを使えばもう少し早くなりそうと思い、先程書いてみた。

import math

def get_sieve_of_eratosthenes(n):
    if not isinstance(n, int):
        raise TypeError('n is int type.')
    if n < 2:
        raise ValueError('n is more than 2')
    prime = [False] * 2 + [True] * (n-1)
    limit = int(math.sqrt(n))
    for i in range(2,limit+1):
        if prime[i]:
            for j in range(i*2,n+1,i):
                prime[j] = False
    return prime

import time
t1 = time.time()
tot = 0
TOP = 2000000
prime_ls = get_sieve_of_eratosthenes(TOP)
for i in range(TOP):
    if prime_ls[i]:
        tot += i
print(tot)
t2 = time.time()
print(t2-t1) #0.6140727996826172

import math

def get_sieve_of_eratosthenes(n):

if not isinstance(n, int):

raise TypeError('n is int type.')

if n < 2:

raise ValueError('n is more than 2')

prime = [False] * 2 + [True] * (n-1)

limit = int(math.sqrt(n))

for i in range(2,limit+1):

if prime[i]:

for j in range(i*2,n+1,i):

prime[j] = False

return prime

import time

t1 = time.time()

tot = 0

TOP = 2000000

prime_ls = get_sieve_of_eratosthenes(TOP)

for i in range(TOP):

if prime_ls[i]:

tot += i

print(tot)

t2 = time.time()

print(t2-t1) #0.6140727996826172

0.6秒と格段に早くなった。

ここで自作のエラトネスのふるいの関数（get_sieve_of_eratosthenes）について少し説明を加えておく。

エラトネスのふるいとは

指定された整数以下のすべての素数を探すアルゴリズム。

pythonのような遅めの言語でも200万ぐらいまでの素数なら一瞬でリスト化できる。

詳しくはwikiとかで（「エラトネスの篩」）

ただリスト化するだけだと...

実はただ素数のリストを作るだけでは実行速度の問題に悩まされる。

例えば次の例。

prime = [2,3,5,7,11,13,17,19]
a = 13
if a in prime:
    print('a is prime')
else:
    print('a is not prime')

prime = [2,3,5,7,11,13,17,19]

a = 13

if a in prime:

print('a is prime')

else:

print('a is not prime')

実行結果:

a is prime

1	a is prime

primeは20までの素数をリスト化したもの、調べたい数をaとした時、

a in prime

1	a in prime

はリストprimeにaと同じ値が入っていた場合True,入ってない場合Falseを返す。

つまり、これで素数判定ができる。

しかし、これではa=13がリストに入っているか探す時に、どうやらリストの最初から一つ一つ調べるため無駄に時間がかかる。

2,3,5,7,11,13と6回も探索しなくてはならない。

リストと探す数が大きくなればなるほどこの無駄が大きくなり、速度が遅くなる。

解決策

結論から言えば、n番目の要素に整数nが素数かどうかの情報を入れたリストを作れば良い。

次のようなものである。

prime = [False, False, True, True, False, True, False, True, False, False, False]
a = 8
if prime[a]:
    print('a is prime')
else:
    print('a is not prime')

prime = [False, False, True, True, False, True, False, True, False, False, False]

a = 8

if prime[a]:

print('a is prime')

else:

print('a is not prime')

実行結果:

a is not prime

1	a is not prime

これならば、

prime[8]

prime[8]

はリストの8番目に入っている値を持ってくるだけなので1回で済む。

この方法（リストのn番目の整数nのbool値を入れる方法）は他の場面でも使われるので大事な考え。

キャッシュを作るのと似た考え？

my function

ということで、自作のエラトネスのふるいは次のようにbool値を入れたリストを作るように書いた。

#参考「https://qiita.com/fantm21/items/5e270dce9f4f1d963c1e」
#参考「https://oku.edu.mie-u.ac.jp/~okumura/python/sieve.html」
import math

def get_sieve_of_eratosthenes(n):
    if not isinstance(n, int):
        raise TypeError('n is int type.')
    if n < 2:
        raise ValueError('n is more than 2')
    prime = [False] * 2 + [True] * (n-1)
    limit = int(math.sqrt(n))
    for i in range(2,limit+1):
        if prime[i]:
            for j in range(i*2,n+1,i):
                prime[j] = False
    return prime

if __name__ == '__main__':
    print(get_sieve_of_eratosthenes(10))
#[False, False, True, True, False, True, False, True, False, False, False]

#参考「https://qiita.com/fantm21/items/5e270dce9f4f1d963c1e」

#参考「https://oku.edu.mie-u.ac.jp/~okumura/python/sieve.html」

import math

def get_sieve_of_eratosthenes(n):

if not isinstance(n, int):

raise TypeError('n is int type.')

if n < 2:

raise ValueError('n is more than 2')

prime = [False] * 2 + [True] * (n-1)

limit = int(math.sqrt(n))

for i in range(2,limit+1):

if prime[i]:

for j in range(i*2,n+1,i):

prime[j] = False

return prime

if __name__ == '__main__':

print(get_sieve_of_eratosthenes(10))

#[False, False, True, True, False, True, False, True, False, False, False]

参考「 https://qiita.com/fantm21/items/5e270dce9f4f1d963c1e 」

参考「 https://oku.edu.mie-u.ac.jp/~okumura/python/sieve.html 」

以上！

Problem10まででもこれだけ勉強になるProject Eulerは偉大。

次の記事→「 https://py-study.net/project-ruler-11-15-python/ 」