pythonでセンター試験解いてみた～その3～

前回は平成31年の数１Aにある第一問の（１）解きました。その続きになります。

今回は平成31年の数１Aにある第四問をやっていきたいと思います。いわゆる「整数問題」ですが、恐れず行ってみましょう！！！

準備
センター試験を解いてみよう
まとめ
参考

準備

ぼくがセンター試験を解くに当たり準備したものは

センター試験過去問（平成31年版）
jupyternotebook

の二つです。

jupyternotebookの導入や使い方についてはこちらを参照ください。

センター試験を解いてみよう

解く問題

(1)不定方程式$$49x-23y =1 $$の解となる自然数$x,y$の中で、$x$の値が最小のものは$$x=[ア],y=[イウ]$$であり、すべての整数解は、$k$を整数として$$x=[エオ]k+[ア],y=[カキ]k+[イウ]$$と表せる。
(2)49の倍数である自然数$A$と23の倍数である自然数$B$の組($A,B$)を考える。$A$と$B$の差の絶対値が1となる組($A,B$)の中で、$A$が最小になるのは $$(A,B) =(49\times[ク],23\times[ケコ])$$ である。また、$A$と$B$の差の絶対値が2となる組($A,B$)の中で、$A$が最小になるのは$$(A,B) =(49\times[サ],23\times[シス])$$ である。
(3)連続する三つの自然数$a,a+1,a+2$を考える。 $$ \begin{cases} aとa+1の最大公約数は1\\ a+1とa+2の最大公約数は1\\ aとa+2の最大公約数は1または[セ]\\ \end{cases} $$ である。
また、次の条件がすべての自然数$a$で成り立つような自然数$m$のうち、最大のものは$m=[ソ]$である。$$\\条件：a(a+1)(a+2)はmの倍数である。$$
(4)6762を素因数分解すると$$6762 = 2\times [タ]\times 7^{[チ]}\times [ツテ]$$ である。
$b$を$b(b+1)(b+2)$が6762の倍数となる最小の自然数とする。このとき$b=[トナニ]$である。
平成31年センター試験第四問（途中中略あり）

不定方程式

まず、不定方程式を解いていきます。これを解くにはdiophantineという関数を使います。

import sympy
from sympy import symbols
from sympy.solvers.diophantine import diophantine

x,y = symbols("x y")
l = sympy.symbols("l", integer=True)
diophantine(49*x - 23*y - 1)

出力は

{(23*t_0 + 8, 49*t_0 + 17)}

となり、$t_0=0$のときにxが自然数で最小の値は8,そのときのyは17ということがわかります。

よって、 $$x=[8],y=[17]$$ さらにすべての整数解は、 $k$を整数として$$x=[23]k+[8],y=[49]k+[17]$$と表せます。

絶対値の計算

$A$と$B$の差の絶対値が1となる組をもとめるために絶対値の方程式を解いていきます。そのときの解は[ク]が0~9,[ケコ]が10~99と回答欄から読み取れるので、すべての場合を計算していきます。（ゴリ押し）

 import numpy as np from numpy import abs a,b = symbols("a,b",real=True) for p in range(10):     for q in range(10,99):         f = np.abs(49*a-23*b).subs(a,p).subs(b,q)         if f == 1:             print("a = {} ".format(p))             print("b = {} ".format(q))

出力は

a = 8 
b = 17

で $$(A,B) =(49\times[8],23\times[17])$$ が正解です。

絶対値が2になる場合も同様で、

for p in range(10):
    for q in range(10,99):
        f = np.abs(49*a-23*b).subs(a,p).subs(b,q)
        if f == 2:
            print("a = {} ".format(p))
            print("b = {} ".format(q))

出力は

a = 7 
b = 15

で $$(A,B) =(49\times[7],23\times[15])$$ が正解です。

連続する3つの自然数

連続する3つの自然数の最大公約数を求める問題です。

回答の候補が1ともう一つしかないので、とりあえず自然数aに10000（適当）までの数字を代入し、条件に合うものの最大値(正解は1より大きいはずなので)を取ればいいのではないかと考えました。

import math

ans = []
for p in range(10000):
    if math.gcd(p,p+1)==1 and math.gcd(p+1,p+2)==1 :
        q = math.gcd(p,p+2)
        ans = ans+[q]
print(max(ans))

出力は

で $$a+1とa+2の最大公約数は1\\ aとa+2の最大公約数は1または[2]$$ となりました。

また、$$a(a+1)(a+2)はmの倍数である$$ときにこれを満たす$m$の最大値を求めます。

今回もとりあえず自然数aに10000（適当）までの数字を代入し、条件に合うものの個数をカウントしていきます。

import collections

ans = []
for p in range(10000):
    for m in range(1,9):
        if p*(p+1)*(p+2)%m == 0:
            ans = ans+[m]
            
c = collections.Counter(ans)

print(c)

出力

Counter({1: 10000, 2: 10000, 3: 10000, 6: 10000, 4: 7500, 8: 6250, 5: 6000, 7: 4285})

で6までが10000まで自然数aにおけるa(a+1)(a+2)を割ることができています。

aを100000にしようと1000000にしようとこれは成り立つので、mの最大値は$m=[6]$となります。

素因数分解

素因数分解をします。

from sympy.ntheory import factorint
a = 6762
print(factorint(a))

出力は

{2: 1, 3: 1, 7: 2, 23: 1}

となり、 $$6762 = 2\times [3]\times 7^{[2]}\times [23]$$ です。

次にb(b+1)(b+2)が6762の倍数であるときの最小の自然数は、解答欄の形からｂは100～999のいずれかであると考えられるので、これをすべて計算させると

ans=[]
for b in range(100,999):
    if b*(b+1)*(b+2)%6762 == 0:
        ans = ans +[b]
        
print(min(ans))

出力は

です。

まとめ

ここまでのコードをまとめておきます。

import sympy
from sympy import symbols
from sympy.solvers.diophantine import diophantine
import numpy as np
from numpy import abs
import math
import collections
from sympy.ntheory import factorint

x,y = symbols("x y")
l = sympy.symbols("l", integer=True)
diophantine(49*x - 23*y - 1)


a,b = symbols("a,b",real=True)#real=Trueこれがないとエラー

for p in range(10):
    for q in range(10,99):
        f = np.abs(49*a-23*b).subs(a,p).subs(b,q)
        if f == 1:
            print("a = {} ".format(p))
            print("b = {} ".format(q))
            
for p in range(10):
    for q in range(10,99):
        f = np.abs(49*a-23*b).subs(a,p).subs(b,q)
        if f == 2:
            print("a = {} ".format(p))
            print("b = {} ".format(q))


ans = []
for p in range(10000):
    if math.gcd(p,p+1)==1 and math.gcd(p+1,p+2)==1 :
        q = math.gcd(p,p+2)
        ans = ans+[q]
print(max(ans))


ans = []
for p in range(10000):
    for m in range(1,9):
        if p*(p+1)*(p+2)%m == 0:
            ans = ans+[m]
            
c = collections.Counter(ans)

print(c)

a = 6762
print(factorint(a)) 
ans=[]
for b in range(100,999):
    if b*(b+1)*(b+2)%6762 == 0:
        ans = ans +[b]
        
print(min(ans))