Pythonでの統合

Question

私はPythonで統合することを意味してきましたが、私はScipy、Numpy、またはPythonに統合できる他のプログラムは使用していません。コーディングに関してはかなり初心者ですが、統合の手助けが必要です。私は小さなコードをコピーしましたが、それでも改善する必要があります。

def LeftEndSum(startingx, endingx, numberofRectangles) : 

    width = (float(endingx) - float(startingx))/numberofRectangles 

    runningSum = 0 

    for i in range(numberofRectangles) : 
     height = f(startingx + i*width) 
     area = height * width 
     runningSum += area 
    return runningSum 

私は統合しようとしているが、私は、私はその後、統合

の終わりでのプロットにグラフ私は間隔を定義する考えを持っていたことができ、データポイントのリストを取得したいです[a、b]とdelta n =（ポイント間のボックスの変化）で、ループを止めてポイントを得るためにコンバージェンステストを行うことができます。収束テストが

行くかのI（N（古い値）+デルタ（N）） - I（N（古い値））/ I（N（旧））<イプシロン ここで、ε= 1×10^-6コードは

Answer 1

は、あなたが、yは0.0〜10.0のx * xを=統合したいとしましょう破る

積算値であればここで、我々は答えは333.33333である知っています。ここではそれを行うにはいくつかのコードです：私たちは最後に反復回数に対する答えをプロットし、反復回数が増えると、それは非常にゆっくりと本当の答えに近づく

def y_equals_x_squared(x): 
    y = x*x 
    return y 

def LeftEndSum(startingx, endingx, numberofRectangles) : 
    width = (float(endingx) - float(startingx))/numberofRectangles 
    #print "width = " + str(width) 
    runningSum = 0 
    i = 1 
    while i <= numberofRectangles: 
     x = (endingx - startingx)/(numberofRectangles) * (i - 1) + startingx 
     height = y_equals_x_squared(x) 
     area = height * width 
     #print "i, x , height, area = " + str(i) + ", " + str(x) + ", " + str(height) + ", " + str(area) 
     runningSum += area 
     i += 1 
    return runningSum 
#----------------------------------------------------------------------------- 
startingx = 0.0 
endingx = 10.0 
# 
numberofRectangles = 3 
old_answer = LeftEndSum(startingx, endingx, numberofRectangles) 
# 
numberofRectangles = 4 
new_answer = LeftEndSum(startingx, endingx, numberofRectangles) 
# 
delta_answer = abs(new_answer - old_answer) 
# 
tolerance = 0.0001 
max_iterations = 500 
iteration_count = 0 
iterations = [] 
answers = [] 
while delta_answer > tolerance: 
    numberofRectangles += 100 
    new_answer = LeftEndSum(startingx, endingx, numberofRectangles) 
    delta_answer = abs(new_answer - old_answer) 
    old_answer = new_answer 
    iteration_count += 1 
    iterations.append(iteration_count) 
    answers.append(new_answer) 
    print "iteration_count, new_answer = " + str(iteration_count) + ", " + str(new_answer) 
    if(iteration_count > max_iterations): 
     print "reached max_iterations, breaking" 
     break 
# 
OutputFile = "Integration_Results.txt" 
with open(OutputFile, 'a') as the_file: 
    for index in range(len(answers)): 
     the_file.write(str(index) + " " + str(answers[index]) + "\n") 
# 
import matplotlib.pyplot as plt 
# 
fig, ax = plt.subplots() 
ax.plot(iterations, answers, 'r-', label = "Increasing # Rectangles") 
title_temp = "Simple Integration" 
plt.title(title_temp, fontsize=12, fontweight='bold', color='green') 
ax.legend(loc='best', ncol=1, fancybox=True, shadow=True) 
plt.xlabel('Number of Iterations') 
plt.ylabel('Answer') 
ax.grid(True) 
plt.show(block=True) 

注意してください。台形ルールなどの単純な四角形よりも優れた統合方法があります。あなたがwhileループを入れて公差をチェックするときは、常にmax_iterationsチェックを入れて、無限ループに陥ることはありません。

あなたはここにあなたの答えを確認することができます。 http://www.integral-calculator.com/ 我々は答えはそれが関数にパラメータとしてfを含めるために、より理にかなっている333.3333

Answer 2

である知っている方法ですこと。なぜ固定式のハードワイヤですか？さらに、コードが大幅に発電機に適用される組み込み関数sumを使用することによって単純化することができます。

def riemannSum(f,a,b,n,sample = 'L'): 
    """computes Riemann sum of f over [a,b] using n rectangles and left ('L'), right ('R') or midpoints ('M')""" 
    h = (b-a)/float(n) 
    if sample.upper() == 'L': 
     s = a #s = first sample point 
    elif sample.upper() == 'R': 
     s = a + h 
    else: 
     s = a + h/2.0 
    return h*sum(f(s+i*h) for i in range(n)) 

あなたが明示的に関数を定義し、それらを統合することができます。

>>> def reciprocal(x): return 1/x 

>>> riemannSum(reciprocal,1,2,100) 
0.6956534304818242 

（正確な値があります

：あなたは匿名関数（ラムダ式）を使用することができ、約0.693147）

であるか、2の自然対数、

それとも、あなたはmathモジュールですでに機能を使用することができます。これらのメソッドの

>>> riemannSum(math.sin,0,math.pi,10) 
1.9835235375094546 

どれも非常に正確ではありません。

>>> simpsonsRule(math.sin,0,math.pi,10) 
2.0001095173150043 

これは、はるかに正確な10枚の長方形とリーマン和よりもあるが（真の値は次のとおりです。たとえば

def simpsonsRule(f,a,b,n): 
    if n%2 == 1: 
     return "Not applicable" 
    else: 
     h = (b-a)/float(n) 
     s = f(a) + sum((4 if i%2 == 1 else 2)*f(a+i*h) for i in range(1,n)) + f(b) 
     return s*h/3.0 

：より正確には、Pythonで行うことも非常に簡単であるSimpson's Ruleです2）。