近世における対数関数や指数関数の発見は人類の計算能力を飛躍的に向上させました(乗除算を加減算に置き換える魔法)。それがそれ以上のものだと知るまでにはそれほどの時間は要さなかったのです。
【数理考古学】常用対数表(Table of Common Logarithms)を使った計算
「2」と「3」の挟間
使い勝手を考えると「xy直交座標系においてxの増分が1の時にyの増分も1となる(傾きが1となる)根」を見定めたい。とりあえずざっくりとした目安として物理学で等加速運動の時間あたりの平均速度を計算する式$\frac{b-a}{2}$を用いると$\frac{α^{+1}-α^{-1}}{2}$となり、2を根とした冪算$2^n$の結果が$\frac{2^{+1}-2^{-1}}{2}=\frac{2-\frac{1}{2}}{2}=\frac{3}{4}$となって$2^0=1$を下回り、3を根とした冪算$3^n$の結果が$\frac{3^{+1}-3^{-1}}{2}=\frac{3-\frac{1}{3}}{2}=\frac{4}{3}$となって$3^0=1$を上回る事から探してる根αの値が2より大きく3より小さい事が明らかとなる(逆数関係にあって掛け合わせると1になるのが興味深い)。
実はこの式自体はハイパボリックサインSinhのものなので1でなく$Sinh_e(1)=1.17520119$に収束してしまう。おや?見慣れない添字が? ここでは対数みたいに「底」概念を導入して考えた方が話が早いのである。
$a^{+x}=\frac{a^{+x}+a^{-x}}{2}+\frac{a^{+x}-a^{-x}}{2}=Cosh_a(x)$
$a^{-x}=\frac{a^{+x}+a^{-x}}{2}-\frac{a^{+x}-a^{-x}}{2}=Sinh_a(x)$
正しい答えはヤコブ・ベルヌーイ(Jakob Bernoulli,1654年~1705年)が会計学の複利計算式より導入した$(1-\frac{1}{n})^n$によって得られ、その値2.718281828をネイピア数という。
%matplotlib nbagg
import math as m
import cmath as c
import numpy as np
import matplotlib.pyplot as plt
import matplotlib.animation as animation
# 等差数列を生成
X = np.linspace(-10,10,101,endpoint =True)
Y0=np.exp(X)
plt.style.use('default')
fig = plt.figure()
Tcode=np.linspace(0.1,4,40,endpoint =True)
def exponential(n):
plt.cla()
#値域の決定
Y1=Tcode[n]**X
def tangent(x):return(np.log(Tcode[n])*x+1)
Y2=tangent(X)
# 可視化
ex01="y="+str(round(Tcode[n],3))+"^x"
ex02="y=log("+str(round(Tcode[n],3))+")x+1"
plt.plot(X, X+1,color="green", label="y=x+1")
plt.plot(-X, X+1,color="green", label="y=-x+1")
plt.plot(X, Y0,color="purple", label="y=exp(x)")
plt.plot(-X, Y0,color="purple", label="y=exp(-x)")
plt.plot(X, Y1,color="red", label=ex01)
plt.plot(X, Y2,color="blue", label=ex02)
slopeX1=[0,1,1,0]
slopeY1=[1,1,tangent(1),1]
plt.fill(slopeX1,slopeY1,color="y",alpha=0.5)
slopeX2=[0,-1,-1,0]
slopeY2=[1,1,tangent(-1),1]
plt.fill(slopeX2,slopeY2,color="y",alpha=0.5)
#頂点
plt.plot(0, 1,color="black", marker="o")
plt.plot(1, 1,color="black", marker="o")
plt.plot(-1, 1,color="black", marker="o")
plt.plot(1, tangent(1),color="blue", marker="o")
plt.plot(-1, tangent(-1),color="blue", marker="o")
#座標
plt.hlines(0,-10,10,color="black")
plt.hlines(1,-10,10,color="black",lw=0.5)
plt.hlines(-1,-10,10,color="black",lw=0.5)
plt.vlines(0,-10,10,color="black")
plt.ylim([-4,4])
plt.xlim([-4,4])
plt.xlabel("X")
plt.ylabel("Y")
plt.title("Variation of Exponentiation Root")
ax = fig.add_subplot(111)
ax.set_aspect('equal')
ax.text(0, 1," (0,1)",color="black", size=9)
ax.text(1, 1," (1,1)",color="black", size=9)
ax.text(-1, 1," (-1,1)",color="black", size=9)
tan1="(1,"+str(round(tangent(1),3))+")"
ax.text(1, tangent(1),tan1,color="blue", size=9)
tan2="(-1,"+str(round(tangent(-1),3))+")"
ax.text(-1, tangent(-1),tan1,color="blue", size=9)
ax.legend(loc='lower center')
#plt.show()
ani = animation.FuncAnimation(fig, exponential, interval=50,frames=len(Tcode))
ani.save("nap012.gif", writer="pillow")
この考え方の効用
%matplotlib nbagg
import math as m
import cmath as c
import numpy as num
import matplotlib.pyplot as plt
from functools import reduce
from itertools import accumulate
import matplotlib.animation as animation
def Napier_culc(ABS,n):
RIM=ABS*m.pi*(0+1j)
Tarm0=num.repeat(RIM/n,n)
Tarm=num.concatenate([[(1+0j)], Tarm0])
return accumulate(Tarm,lambda x,y:x+x*y)
#単位円データ作成
c0=num.linspace(0,2*m.pi,61,endpoint = True)
s0=[]
for nm in range(len(c0)):
s0.append(complex(m.cos(c0[nm]),m.sin(c0[nm])))
s1=num.array(s0)
plt.style.use('default')
fig = plt.figure(111)
#関数定義
def Eulerien_Identity(n):
plt.cla()
#ネイピア数算出
ei01=num.array(list(Napier_culc(ABS,Time_Code[n])))
EI0=[]
for nm in range(len(ei01)):
EI0.append(ei01[nm].conjugate())
ei02=num.array(EI0)
#円周描画
plt.plot(s1.real,s1.imag,color="green", label="Circle")
plt.plot(ei01.real,ei01.imag,color="blue", label="(1+1/n)^n")
plt.plot(ei02.real,ei02.imag,color="red", label="(1-1/n)^n")
plt.xlim([-m.pi,m.pi])
plt.ylim([-m.pi,m.pi])
plt.title("Eulerean Identity")
plt.xlabel("Real")
plt.ylabel("Imaginal")
ax = fig.add_subplot()
ax.set_aspect('equal', adjustable='box')
#補助線描画
plt.axvline(0, 0, 1,color="black")
plt.axhline(0, 0, 1,color="black")
plt.plot([-1,0],[0,0],color="red",lw=1)
plt.plot([0,1],[0,0],color="blue",lw=1)
ax.text(0, 0, "0",color="black",size=13)
ax.text(1, 0, "1",color="blue",size=13)
ax.text(-1, 0, "-1",color="red",size=13)
ax.text(0, 1, "π",color="green",size=13)
ax.text(0, -1, "-π",color="green",size=13)
#凡例描画
ax.legend(loc='lower right')
ABS=1
Time_Code=[1,2,4,8,16,32,64,128]
#Eulerien_Identity(0)
#plt.show()
ani = animation.FuncAnimation(fig, Eulerien_Identity, interval=100,frames=len(Time_Code))
ani.save("EI01.gif", writer="pillow")
問題は計算精度で、どんなに高性能のコンピューターを使って無限に周回運動を続ける事は出来ません。そこで全体を半周(πラジアン)なり1周(2πラジアン)なりを無限に繰り返してると考え、誤差を最小限に抑え込もうとするオイラーの定理$e^{θi}=Cos(θ)+Sin(θ)i$の出番となるという…
ABS=6
Time_Code=[1,2,4,8,16,32,64,128,256,512,1024,2048]
ちなみに複素関数論ではこの様に「1を保守する」作業を「正規化」ないしは「規格化」と呼ぶ様です。
複素数値関数の正規化または規格化
それにつけても冪算は「幾何級数(Geometric-series)的に増える」なる表現の由来となってるくらいだから増率が尋常じゃありません。例えばある生物の大きさが昨年${\frac{1}{2}}$で来年2倍になる見込みなら${2^n}$のオーダーで、昨年${\frac{1}{3}}$で来年3倍になる見込みなら${3^n}$のオーダーで成長してると推察されます。計算上、比較年次を「前年(-1)」「当年(0)」「来年(+1)」としたかったのでこの設定。まぁネイピア数はこの挟間で見つかった訳ですが…
その概念の登場こそがマルサス「人口論(An Essay on the Principle of Population,1978年)」)を生み、そのあまりに悲観的な人口爆発予測を覆す為にロジスティック方程式が開発されたという次第。
ロジスティック方程式とその関連範囲
N=\frac{K}{1+(\frac{K}{N_0}-1)e^{-rt}}
なお…
- N=観測終了時間$t_x$時点での生物個体数(マルサス式だと人口爆発)
- $N_0$=観測開始時間$t_0$時点での生物個体数
- K=環境収容力(その環境が維持できる個体数)。ほぼ定数扱い。
- r=内的自然増加率(一個体当たりの増加率=その生物が実現する可能性のある最大増加率)。ほぼ定数扱い。
おや、ここにもネイピア数が? そんな感じで以下続報…