我、京大生ぞ

現役京大生の雑記。メインテーマは大学受験とプログラミングと英語と京大の日常

2018年の京大数学をプログラミング(Python、Sympy)で解いてみた


こんにちは、現役京大生ブロガーのゲーテ(@goethe_kyodai)です。


普通に解いても面白くないのでプログラミング(Python)で2018年の理系の京大数学を解いてみました


プログラミングで解くことに重きを置き、厳密な議論は飛ばしていて模範解答とは程遠い解答となっているので高校生のみんなは参考にしないように。


使用言語:Python2系
使用ライブラリ:Sympy,Numpy,Math,matplotlib


Sympyの使い方はこちらを参考にしました。



はてなブログにコードを埋め込むのはこちらを参考にしました。

 

目次

問題1

f:id:a86223990:20180310142205p:plain

(1)

C_1C_2からyを消去して整理して、

(b - a)x^2 - 2bx + b + c ・・・(1)

を得る。

上式の判別式をdiscriminant()で求めると、

 D = b^2 - (b-a)(b+c)

となる。

これをsolve()でbについて解くと、

b = \frac{ac}{c-a}

となる。

これをsubs()で(1)に代入してsolve()でxについて解くと、

x = \frac{c}{a}

を得る。

これをy = ax^2にsubs()で代入して、

t = \frac{c^2}{a}

を得る。

これらが答えの接点の座標である。

まとめると答えは

(\frac{c}{a}, \frac{c^2}{a})


(2)

まず先ほど求めた接点を(x,y) =(\frac{c}{a}, \frac{c^2}{a}) と置きます。

そしてsolve()でこれらの連立方程式をx、yについて解くと、

(a,c) = (\frac{y}{x^2}, \frac{y}{x})

となります。

これを条件(1)の不等式を等式に置き換えた、1 + c^2 -2aにsubs()で代入してsolve()でyについて解くと、

y = -\sqrt(-x^2 + 1) + 1 , \sqrt(-x^2 + 1) + 1

を得られます。

このyと閉区間 [-1,1] の範囲のxを使って、matplotlibで図示したのが下です。

f:id:a86223990:20180310175749p:plain

本来の条件(1)は不等式なので上図の円に囲まれた部分が答えです。(本来はa,b,cの条件から除外する点があるがめんどくさいので省略)


問題2

f:id:a86223990:20180310142250p:plain


f(n) = n^3 - 7*n + 9と置きます。

因数分解できなさそうなので代入して構造を調べました。

その前に整数nじゃ範囲が広すぎるので、f(n)が素数になり得るnの範囲を調べておきましょう。

g(x) = x^3 - 7*x + 9と置きます。

これをmatplotlibのpyplotで図示すると下のようになります。

f:id:a86223990:20180310153334p:plain

見た感じg(x) > 0となるのはx \geqq -3ですね。これで調べる範囲をn \geqq -3に限定できました。

では実際に代入して素数かどうかの判定をとりあえずn = -3 ~ 20 までしてきます。
プログラムで実際に代入した結果がこちら↓

n = -3,f(n) = 3で、「f(n)が素数である」はTrue
n = -2,f(n) = 15で、「f(n)が素数である」はFalse
n = -1,f(n) = 15で、「f(n)が素数である」はFalse
n = 0,f(n) = 9で、「f(n)が素数である」はFalse
n = 1,f(n) = 3で、「f(n)が素数である」はTrue
n = 2,f(n) = 3で、「f(n)が素数である」はTrue
n = 3,f(n) = 15で、「f(n)が素数である」はFalse
n = 4,f(n) = 45で、「f(n)が素数である」はFalse
n = 5,f(n) = 99で、「f(n)が素数である」はFalse
n = 6,f(n) = 183で、「f(n)が素数である」はFalse
n = 7,f(n) = 303で、「f(n)が素数である」はFalse
n = 8,f(n) = 465で、「f(n)が素数である」はFalse
n = 9,f(n) = 675で、「f(n)が素数である」はFalse
n = 10,f(n) = 939で、「f(n)が素数である」はFalse
n = 11,f(n) = 1263で、「f(n)が素数である」はFalse
n = 12,f(n) = 1653で、「f(n)が素数である」はFalse
n = 13,f(n) = 2115で、「f(n)が素数である」はFalse
n = 14,f(n) = 2655で、「f(n)が素数である」はFalse
n = 15,f(n) = 3279で、「f(n)が素数である」はFalse
n = 16,f(n) = 3993で、「f(n)が素数である」はFalse
n = 17,f(n) = 4803で、「f(n)が素数である」はFalse
n = 18,f(n) = 5715で、「f(n)が素数である」はFalse
n = 19,f(n) = 6735で、「f(n)が素数である」はFalse
n = 20,f(n) = 7869で、「f(n)が素数である」はFalse


どうやら n \geqq 3ではf(n)は素数にならないっぽいので、n=3m,3m+1,3m+2(m は 整数で1以上)と置いてそれを証明してみます。

とりあえず代入してみます。

f(3m) = 3(9m^3 - 7m + 3)
f(3m+1) = 3(9m^3 + 9m + 1)
f(3m+2) = 3(18m^3 + 5m + 1)

いい感じに3 \times(多項式)になってますね。

あとはf(3m)、f(3m+1)、f(3m+1)を3で割った 9m^3 - 7m + 39m^3 + 9m + 118m^3 + 5m + 1が2以上の整数になっていれば、f(3m)、f(3m+1)、f(3m+1)が合成数になるので素数にならないことを証明できます。

それをmatplotlibのpyplotで図示して目視で確認しましょう。

f:id:a86223990:20180310161320p:plain

それぞれ2以上になってますね。

以上より答えはn = -3、1、2で解けました。


問題3

f:id:a86223990:20180310142256p:plain


f:id:a86223990:20180310180621p:plain


変数\theta (0 < \theta < \alpha)を導入して、条件から上図のように角を表します。

外接円の半径が1であることと正弦定理より
 \frac{AB}{sin(\pi  - ( \alpha + \theta)) } = \frac{BC}{sin(\theta)}
 = \frac{CD}{sin(\alpha - \theta)} = \frac{DA}{sin(\theta)} =2

よって
AB = 2sin(\pi  - ( \alpha + \theta))
BC = 2sin(\theta)
CD = 2sin(\alpha - \theta)
DA = 2sin(\theta)
と表せます。

これらを使ってkは

k = 16sin(\theta)^2sin(\pi  - ( \alpha + \theta)sin(\alpha - \theta)

と表せます。(ここまでプログラミング登場なし)

入試問題的には最大値を\alphaで表せばいいのですが、プログラミングを駆使して解きたいので0 < \alpha \leqq \frac{\pi}{2} を満たす任意の実数の入力に対して、線形探索で最大値を出すプログラムを書きました。

そのプログラムのアルゴリズムは
1:\alphaを受け取る
2:0 < \theta < \alpha の範囲で任意個に均等に分けた\thetaを用意する。
3:さっき分けた\thetaを小さい方から順に\alphaと共にkに代入して計算し、kの最大値を線形探索で求める


線形探索なので計算量はO(n)(n は \theta を分割した数)です。

問題4

f:id:a86223990:20180310142303p:plain

解法1

 z_n = 1 となる確率をP( z_n = 1)とおく。

P( z_n = 1 ) = \frac{z_n = 1となるz_nの数}{全ての z_n} です。

[全ての z_n ]は2^nで、[ z_n= 1となるz_nの数] ]は再帰を使った全探索で求められます。

計算量がO(2^n)になってしまいますが、計算時間の制限がなければ任意の入力nに対して一応この解法で解けます。


解法2

解法1は計算量が大きすぎるので速い方法で解いてみます。

a = 1、b = \frac{-1+\sqrt{3} i}{2}、c = \frac{-1-\sqrt{3} i}{2}と置きます。

z_n がa、b、cとなる確率をa_nb_nc_nと置きます。

それらを用いた以下の漸化式が成り立ちます。

a_n + c_n = 2a_{n+1}
a_n + c_n = 2b_{n+1}
b_n = c_{n+1}

A = \left(
    \begin{array}{ccc}
      1 & 0 & 1 \\
      1 & 0 & 1 \\
      0 & 1 & 0
    \end{array}
  \right)
、B = \left(
    \begin{array}{ccc}
      2 & 0 & 0 \\
      0 & 2 & 0 \\
      0 & 0 & 1
    \end{array}
  \right)x_{n} = (a_n,b_n,c_n)とすると上の漸化式は

A x_n = B x_{n+1}

と表せます。

式変形すると

x_{n} = B^{-1}Ax_{n-1}

と表せます。C = B^{-1}Aと置くと

x_{n} = C x_{n-1}

となるので

x_{n} = C^{n-1} x_1

と表せます。

C^{n-1}をプログラミングの計算力を使ってCをn-1回掛けて求めてもいいですが、計算量がO((3^2)^n)になってしまうので対角化しましょう。

D = P^{-1} C PとなるようにCを対角化します。この時、Cの固有値をdiagnolize()で求めて、\lambda_1、\lambda_2、\lambda_3とすると、

C^{n-1} = P D^{n-1} P^{-1} = P \left(
    \begin{array}{ccc}
      \lambda_1^{n-1} & 0 & 0 \\
      0 & \lambda_2^{n-1} & 0\\
      0 & 0 & \lambda_3^{n-1}
    \end{array}
  \right) P^{-1}

と表せて計算量がO(n)ですみます。

このC^{n-1}x_1 = (\frac{1}{2},\frac{1}{2},0)を使ってx_n = (a_n,b_n,c_n)を求められます。

このa_nが求める確率P( z_n = 1)です。


問題5

f:id:a86223990:20180310142307p:plain

(1)

u(t)、v(t)を普通に求めると、

(u(t), v(t))

 = (t + \frac{1}{\sqrt{t^2 + 1}},-\frac{t}{\sqrt{t^2 + 1}} + log(t))

これらをdiff()を使ってtで微分すると答えを出せます。

答えは

(\frac{du}{dt}, \frac{dv}{dt})
 = (1 - \frac{t}{(t^2+1)\sqrt{t^2 + 1}},\frac{1}{t} - \frac{1}{(t^2+1)\sqrt{t^2 + 1}})

です。

(2)

 L_1(r) = \int_1^r \sqrt{1+\frac{1}{t^2}} dtL_2(r) = \int_1^r \sqrt{ ( \frac{du}{dt})^2+( \frac{dv}{dt} )^2 } dtと表せます。

integrate()で積分してL_1(r)L_2(r)を求めてlimit()で\lim_{ r \to \infty} L_1(r) - L_2(r)を求めようと試みましたが、、、

L_1(r)はすぐ求められました、L_2(r)の式がゴツすぎるせいか中々結果を吐き出してくれませんでした。

なので式を簡単化します。

(1)で求めた導関数をよく見ると、\frac{du}{dt} = t \times \frac{dv}{dt}の関係になってます。

これを使うとL_2(r) = \int_1^r \frac{du}{dt} \sqrt{ 1+( \frac{1}{t^2} ) } dtと表せ、さらにこれを使ってL_1(r) - L_2(r)を簡単化して

L_1(r) - L_2(r) = \int_1^r \frac{1}{1+t^2} dt

と表せます。(Sympyのcancel()やsimplify()を使って\frac{du}{dt} = t \times \frac{dv}{dt}からここまで簡単化できると思ったんですが中々できなくて手で計算しました。)

このように式変形してlim()でやっと解けるようになります。


問題6

f:id:a86223990:20180310142313p:plain

(1)

まずD= (0,0,0)、A=(a_1,a_2,a_3)、B=(b_1,b_2,b_3)、C= (c_1,c_2,c_3)と置きます。

その座標を使って計算したPQ•ABをLと置きます。PQとABが垂直になることを示せばいいのでL=0を示せばいいですね。

条件AC = BD、AD = BCをsolve()で解いてa_1、a_2をa_3、b_1、b_2、b_3、c_1、c_2、c_3で表して、それらをLにsubs()してsymplify()とexpand()してもLが0にならなかったので方針を変更しました。

条件AC = BD、AD = BCを二乗したAC^2 = BD^2、AD^2 = BC^2を辺々足してexpand()で展開すると、

2a_1 ^2 - 2a_1 b_1 -2a_1 c_1 + 2a_2 ^2 - 2a_2 b_2

- 2a_2 c_2 + 2a_3 ^2 -2a_3 b_3 - 2a_3 c_3 = 0

になります。左辺をMと置きます。

Lをexpand()で展開すると、

L

 = -0.5a_1 ^2 + 0.5a_1 b_1 +0.5a_1 c_1 - 0.5a_2 ^2

+ 0.5a_2 b_2 + 0.5a_2 c_2 - 0.5a_3 ^2 +0.5a_3 b_3 + 0.5a_3 c_3

となります。

二つの式を目視で比べると、L = - 0.25 Mとなっていることがわかるので、L = 0となり、証明できました。(結局目視)


(2)

分からん。

PQを通る平面の方程式を出して、どの線分に交わるかで場合分けして、分けられた図形の体積を導出して比べたら解けそうだけどめんどくさいのでやってません。

いい感じにプログラミングで解ける方法思いついたらコメントで教えてください。

感想

最初はなんでも計算機の計算能力でゴリ押しすれば京大数学なんて解けるだろ!と思ってたんですが全然そんなことありませんでした。解けるように上手く式変形するのは人間の手でやるしかなさそうです。


どんな数学の問題を特にしても、計算機が得意なところと人間が得意なところを分担できる程度の計算機の理解が必要 だとも思いました。あと計算機が解ける形にするために数学はやっぱ大事だなって改めて思いました。数学勉強頑張ります。

参考書


行列式、固有値、固有ベクトルなど線形代数の基本的な内容を図で直感的に理解できる良本です。
 


京大が誇るハイスペック人材たち!


本当にいるハイスペック京大生 - 我、京大生ぞ