関数で関数を微積しよう
平均変化率
定積分と短冊
1.微分積分のもとは細分
このワークシートはMath by Codeの一部です。
<微分と積分の意味づけの確認>
微分は曲線y=f(x)を細分することで、曲線を線分のあつまりとみなす。
すると、x=aでの線分の傾き(平均変化率)がΔf/Δx=m(a)={f(a+Δx)-f(a)}/Δxとなった。
その細分を微小化して瞬間的な傾きにして、それをxの各点で出せる関数y=m(x)にすることだった。
それを数式操作とすると、細分Δをさらに微小にしたdであらわすことで、
微分は曲線の関数式y=f(x)から、瞬間傾きを表す数式y=df(x)/dxを出すことになったね。
積分は数式操作では、微分の逆操作だったけど、
定積分は曲線y=f(x)とx軸y=0と積分区間x=a,z=bの2直線で囲む領域の面積を出すことだった。
そのために、領域をキュウリのように同じ厚さに薄切りにする。
横からみると、中央での高さと厚みをかけた、短冊の面積をはじからはじまで合計しても出せたね。
その細分を微小化して極限としての面積にすることだった。
<微分をコードにする>
haskellで微分してみよう。
t秒後の車の位置がcarPos(t)=cos tで表せるとき、
t秒後の車の位置の平均変化率つまり、速さを微分で求める関数を作り、
1秒後の速さを求めてみる。もちろん
数式微分ではd(cos x)/dx=-sin x だから、-sin(1)になるはずだ。
導関数derivativeはdt x tの3変数関数にみえるけど
dt=0.01, x=carPosとして、カリー化したから、
carV=derivative 0.01 carPosは、
時間を入れると速度がでる1変数関数になりました。
[IN]Haskell
---====================================================
type R= Double
type Derivative = (R -> R) -> R -> R
derivative :: R -> Derivative
derivative dt x t = (x (t + dt/2) - x (t - dt/2)) / dt
carPos :: R -> R
carPos = cos
carV :: R -> R
carV = derivative 0.01 carPos
carV1 = carV 1
carVf = -sin 1
--====================================================
[OUT]
ghci> carV1
-0.8414674786831666
ghci> carVf
-0.8414709848078965
2.微分を作ろう。
質問:juliaやgeogebraで数値微分をどうやって作りますか。
<juliaで数値微分>
t秒後の車の位置がcarPos(t)=cos tで表せるとき、微分してcarV関数を作る。
そして、1秒後の速さを求めよう。
juliaでは関数は第1級オブジェクトだから、値として渡したり、返したりができる。
なので、Haskellと同様なコードがかけるね。Haskellよりも型の書き方が簡単でいいね。
#[IN]julia
#======================================
function derivative(f::Function)::Function
h=0.001
df(x)= (f(x +h) - f(x)) / h
return df
end
f(x)=cos(x)
carV = derivative(f)
x=1.0
println(carV(x), ":",-sin(x))
#======================================
[OUT]
-0.8414674786831666:-0.8414709848078965
また、xを点ではなく点集合としてリストにし、それを導関数でブロードキャストしてy座標を出す。
これをzipで組にすれば関数グラフのもとになるタプル(x,y)が作れるね。
さらにplots.jlなどの視覚化パッケージをインストールすれば、グラフを書いたりできる。
[IN]julia
#======================================
xs= -1.0:0.1:1.0
ys= carV.(x)
dot=[(x,y) for (x,y) in zip(xs,ys)]
println(dot)
#======================================
[OUT]
[(-1.0, 0.841200693432298), (-0.9, 0.7830159741147869), (-0.8, 0.7170076180145202), (-0.7, 0.643835158806283), (-0.6, 0.5642297115148187), (-0.5, 0.47898666745549434), (-0.4, 0.3889577469469163), (-0.3, 0.2950424892033121), (-0.2, 0.19817926443543854), (-0.1, 0.0993358979667347), (0.0, -0.0004999999583255033), (0.1, -0.10033090204919493), (0.2, -0.19915933093161975), (0.3, -0.29599782561273713), (0.4, -0.3898788078642301), (0.5, -0.479864249944395), (0.6, -0.5650550470610582), (0.7, -0.644600000929918), (0.8, -0.7177043246657799), (0.9, -0.783637584031216), (1.0, -0.8417409956931188)]
さて、geogebraにはderivative(関数)という関数があるので、苦も無く導関数は作れるだろう。
位置と速度
2.定積分を作ろう。
定積分は短冊の面積を足し込めば得られる。
t=aからt=bの積分区間を細分した1つの幅をΔt=0.5とすると、
t=a+Δt/2, a+Δt+Δt/2,.....,b-Δt/2が短冊の幅の中点の時刻になるね。
そのときの高さf(t)と幅Δtの積和が領域の面積になるから、
integral(f, a, b)=sum( [ f(t)Δt | t = [ a+Δt/2, a+3Δt/2,.....,b-Δt/2 ] ] )
この細分Δtを微小部分dtにしていけば積分になる。
インテグラルを定義してから、0から1までのy=x2を定積分すると、
integral(x2,0,1)=[1/3 x3] from 0 to 1= 1/3となるはず。
haskellでやってみよう。
[IN]Haskell
#=======================================
type Integration = (R->R)
->R --Lower limit
->R --Upper limit
->R --result
integral :: R -> Integration
integral dt f a b
= sum [f t * dt | t <- [a+dt/2, a+3*dt/2 .. b - dt/2]]
#=======================================
[comand prompt]
ghci> integral 0.01(\x->x**2) 0 1
0.33332499999999987
質問:juliaやgeogebraでは定積分はどうやりますか。
integraは、関数を与えると数値定積分したあとの関数を返す
だからinetgralの引数の関数を変えると、それに合わせて
定積分します。
数学的な定義をかくだけでプログラミングができてしまう。
これって、とても楽だし気分がいいです。
#[IN]julia
#================================================
#数値積分で定積分を作る
function integral(f::Function,a,b)
dt=0.01
return sum([f(t) * dt for t in a+dt/2 : dt : b-dt/2])
end
g(x)=x^2
integral(g,0,1)
#================================================
[OUT]
0.33332500000000004
さて、geogebraではintegral(fx,a,b1)を作らなくても、
integral(関数,開始,終了)だけで定積分ができるね。