Вычисление неправильного значения интеграла с использованием правила Симпсонов

Question

Я написал следующий фрагмент кода для интеграции правил Симпсона, чтобы приблизить грех. Уравнение в этом приложении. Я написал отдельные циклы для четных и нечетных терминов, так как они показаны сгруппированными в приложении.

import math 

x1=0
x2=math.pi
N=6
delta_x=(x2-x1)/N
f=math.sin
sum=f(0)

#odd summation
for i in range (1,N+1):
    sum=sum+f(x1+(2*i+1)*delta_x)

sum=4*sum
print(sum)

#even summation

for i in range(2,N+1):
    sumeven=0
    sumeven=sumeven+f(x1+(2*i)*delta_x)
sumeven=2*sumeven

sumeven=sumeven+f(N)
print(sumeven)

integral=(delta_x/3)*(sum+sumeven) 
print(integral)   

Но когда я печатаю значения, это дает мне очень маленькие отрицательные числа.

Кто-нибудь может увидеть, что не так с моим кодом?

Answer 1

Прежде всего, изучите основы отладки. Вы сделали только два простых print строки для этого кода, ни одна из которых не отслеживает более подробные проблемы. Смотрите этот прекрасный блог отладки для помощи. Я добавил некоторые инструменты, немного очистил стиль и снова запустил код, несколько лучше отслеживая логику и поток данных. Диагноз ниже.

sum_odd = f(0)

# odd summation
for i in range (1, N+1):
    x_val = x1 + (2*i+1)*delta_x
    sum_odd = sum_odd + f(x_val)
    print (x_val/math.pi, "pi", f(x_val), sum_odd)

sum_odd = 4*sum_odd

print(sum_odd)

# even summation
for i in range(2, N+1):
    sum_even = 0
    x_val = x1 + (2*i)*delta_x
    sum_even = sum_even + f(x_val)
    print (x_val/math.pi, "pi", f(x_val), sum_even)

Выход:

0.5 pi 1.0 1.0
0.8333333333333333 pi 0.5000000000000003 1.5000000000000004
1.1666666666666665 pi -0.4999999999999997 1.0000000000000007
1.5 pi -1.0 6.661338147750939e-16
1.8333333333333333 pi -0.5000000000000004 -0.4999999999999998
2.1666666666666665 pi 0.4999999999999993 -4.996003610813204e-16
-1.9984014443252818e-15
0.6666666666666666 pi 0.8660254037844387 0.8660254037844387
1.0 pi 1.2246467991473532e-16 1.2246467991473532e-16
1.3333333333333333 pi -0.8660254037844384 -0.8660254037844384
1.6666666666666665 pi -0.866025403784439 -0.866025403784439
2.0 pi -2.4492935982947064e-16 -2.4492935982947064e-16
-0.27941549819892636
-0.04876720424671586

ДИАГНОЗ

У вас есть две насущные проблемы:

1. Вы инициализируете sum_even внутри цикла, который отбрасывает предыдущие вычисления для этой половины серии. Переместите его до цикла.
2. Вы интегрировали более чем в два раза указанный диапазон; это само по себе даст вам крошечное значение (в идеале это будет 0, интеграл для полного цикла функции).

Исправьте инициализацию, установите ограничения, и вы должны увидеть хорошие результаты, например:

0.16666666666666666 pi 0.49999999999999994 0.49999999999999994
0.5 pi 1.0 1.5
0.8333333333333333 pi 0.5000000000000003 2.0000000000000004
8.000000000000002
0.3333333333333333 pi 0.8660254037844386 0.8660254037844386
0.6666666666666666 pi 0.8660254037844387 1.7320508075688772
1.0 pi 1.2246467991473532e-16 1.7320508075688774
3.184686116938829
1.9520959854268212

Кроме того, я настоятельно рекомендую вам изучить руководства по отладке и стилю кодирования; это поможет вам в дальнейшей работе. Я знаю из опыта.:-)