7

数値メソッド クラスの場合、Simpson の合成規則を使用して定積分を評価するプログラムを作成する必要があります。私はすでにここまでたどり着きましたが (以下を参照)、私の答えは正しくありません。f(x)=x を 0 から 1 に統合してプログラムをテストしています。結果は 0.5 になるはずです。私は 0.78746... などを取得します。Scipy で利用できるシンプソンの規則があることは知っていますが、実際には自分で書く必要があります。

2 つのループに問題があると思われます。以前に「for i in range(1, n, 2)」と「for i in range(2, n-1, 2)」を試したところ、結果は 0.41668333... などになりました。 x += h」と「x += i*h」を試してみました。最初のオプションは 0.3954 で、2 番目のオプションは 7.9218 でした。

# Write a program to evaluate a definite integral using Simpson's rule with
# n subdivisions

from math import *
from pylab import *

def simpson(f, a, b, n):
    h=(b-a)/n
    k=0.0
    x=a
    for i in range(1,n/2):
        x += 2*h
        k += 4*f(x)
    for i in range(2,(n/2)-1):
        x += 2*h
        k += 2*f(x)
    return (h/3)*(f(a)+f(b)+k)

def function(x): return x

print simpson(function, 0.0, 1.0, 100)
4

6 に答える 6

9

2 番目のループの前に初期化するのを忘れている可能性がxあります。また、開始条件と反復回数がオフになっています。正しい方法は次のとおりです。

def simpson(f, a, b, n):
    h=(b-a)/n
    k=0.0
    x=a + h
    for i in range(1,n/2 + 1):
        k += 4*f(x)
        x += 2*h

    x = a + 2*h
    for i in range(1,n/2):
        k += 2*f(x)
        x += 2*h
    return (h/3)*(f(a)+f(b)+k)

あなたの間違いは、ループ不変条件の概念に関連しています。あまり詳しくは説明しませんが、サイクルの最初ではなく、最後に進むサイクルを理解してデバッグする方が一般的に簡単です。ここでは、x += 2 * h行を最後に移動しました。x = a - hあなたの実装では、ループの最初の行として追加するためだけに、最初のループに奇数を割り当てる必要が2 * hあります。

于 2013-04-14T16:27:06.060 に答える