楕円周長の数値積分のエクセルＶＢＡのコード

本日、楕円積分で質問をいくつかしてしまい、疲れて、結果エクセル
ＶＢＡで検証したくなりました。
そこでSimpsonの公式を利用した
f(θ) = ∫{0～π/2} √{1 - k^2 sin^2 θ} dθ
（１/４楕円の周長）
のエクセルＶＢＡのコード例をどなたか作成していただけないでしょうか。

同時にd^4(f(θ）)/dθ^4 の計算結果をご教示お願いします。

No.1 ベストアンサー 回答者： age_momo 回答日時： 2007/03/21 18:45

質問として成立しているかどうか怪しいですが、

プログラム自体は簡単ですから作ってみました。
mainのところの長径、短径を1,1で計算してみましたが、
2πとの誤差は10^-10以下でしたので、これそのものはあっていると
思います。

Option Explicit
Private Function daen(sita As Double, k As Double) As Double

daen = (1 - k * (Sin(sita)) ^ 2) ^ 0.5

End Function

Sub test()

Const bunkatu = 512
Const pai = 3.14159265358979
Const chokei = 8
Const tankei = 6
Dim sum As Double
Dim ct As Integer
Dim keisu As Integer
Dim hiritu As Double

hiritu = 1 - (tankei / chokei) ^ 2
For ct = 0 To bunkatu
If ct = 0 Or ct = bunkatu Then
keisu = 1
ElseIf ct Mod 2 = 0 Then
keisu = 2
Else
keisu = 4
End If
sum = sum + daen(pai / 2 / bunkatu * ct, hiritu) * keisu
Next
sum = sum / bunkatu / 3 * 4 * chokei * pai / 2
MsgBox (sum)
End Sub

この回答への補足

うわーー大変感謝です。

ついでにどなたか
f(θ) = ∫{0～π/2} √{1 - k^2 sin^2 θ} dθ
について
d^4(f(θ）)/dθ^4
を　解いていただけないものでしょうか。

現在上のｆ（θ）の級数展開式をエクセルＶＢＡで組んでいます。
これと比較して数値が一致するようでしたら即採用予定です。

補足日時：2007/03/21 19:49

No.2 回答者： age_momo 回答日時： 2007/03/22 00:34

やっとwikiが見れました。

（こちら中国じゃ、串指さないとwikiにアクセスできないんで)

>現在上のｆ（θ）の級数展開式をエクセルＶＢＡで組んでいます。

こちらのプログラムの方が更に簡単ですよ。載せておきます。
(2n-1)!!/2n!!は逐次計算の際にそれぞれをかけていけば簡単に解決できます。

Sub daenshu()

Const chokei = 8#
Const tankei = 6#
Const pai = 3.14159265358979
Dim ct As Integer
Dim risin As Double
Dim dExcl As Double
Dim sum As Double

risin = 1 - (tankei / chokei) ^ 2
dExcl = 1#
sum = 1#
For ct = 1 To 100
dExcl = dExcl * (2 * ct - 1) ^ 2 / 4 / ct ^ 2
sum = sum - risin ^ ct / (2 * ct - 1) * dExcl
Next
sum = sum * pai * 2 * chokei
MsgBox (sum)

End Sub

この回答へのお礼

ほんとうにありがとうございました。別の方法の下記式と結果が
見事に一致しました。おそらく、楕円の周長計算の方法に決定打
が出たと思います。当方の確認式です。
h1=(a-b)
h2=(a+b)

L=PI()*h2*(1+(h1/h2)^2/4+(h1/h2)^4/64+(h1/h2)^6/256
L=44.20698431

お礼日時：2007/03/22 08:01