Wie man mit Simpsons-Regel in Python integrieren?

Question

Ich versuche, diese Funktion zu integrieren:

def f(x): 
    return (x**4)-(2*x)+1 

N=10 
a=0.0 
b=2.0 
h=(b-a)/N 

s=f(a)+f(b) 

for k in range(1,N/2): 
    s+=4*f(a+(2*k-1)*h) 

for k in range(1,N/(2-1)): 
    s1 +=f(a+(2*k*h) 

M=(s)+(2*s1) 
print((1/3.0)*h)*(3) 

Aber ich habe diesen Fehler:

File "<ipython-input-29-6107592420b6>", line 17 
    M=(s)+(2*s1): 
^
SyntaxError: invalid syntax 

Ich habe versucht, das Schreiben es in verschiedenen x^4 - 2x + 1 von 0 bis 2

ich dieses Programm geschrieben Formulare, aber ich bekomme immer einen Fehler in M

Answer 1

Sie haben eine schließende Klammern in Ihrem zweitenvergessenSchleife hier: s1 += f(a+(2*k*h). Es sollte sein:

Answer 2

Für zukünftige Referenz könnten Sie auch über die Verwendung von scipy.integrate denken. Schauen Sie hier for some methods, die je nach Art und Auflösung Ihres Datensatzes eine bessere Genauigkeit haben können.

Ein Code könnte wie folgt aussehen:

import scipy.integrate as int 
x = [ ii/10. for ii in range(21)] 
y = [ xi**4 - 2*xi + 1 for xi in x] 
tahdah = int.simps(y,x,even='avg') 
print(tahdah) 

Welche ergibt und Antwort von 4.4, die Sie mit Bleistift und Papier bestätigen können.

Answer 3

Haben Sie die code example von Simpsons Rule auf Wikipedia (in Python geschrieben) gesehen? Ich werde es hier zugunsten zukünftiger Leser erneut veröffentlichen.

#!/usr/bin/env python3 
from __future__ import division # Python 2 compatibility 

def simpson(f, a, b, n): 
    """Approximates the definite integral of f from a to b by the 
    composite Simpson's rule, using n subintervals (with n even)""" 

    if n % 2: 
     raise ValueError("n must be even (received n=%d)" % n) 

    h = (b - a)/n 
    s = f(a) + f(b) 

    for i in range(1, n, 2): 
     s += 4 * f(a + i * h) 
    for i in range(2, n-1, 2): 
     s += 2 * f(a + i * h) 

    return s * h/3 

# Demonstrate that the method is exact for polynomials up to 3rd order 
print(simpson(lambda x:x**3, 0.0, 10.0, 2))  # 2500.0 
print(simpson(lambda x:x**3, 0.0, 10.0, 100000)) # 2500.0 

print(simpson(lambda x:x**4, 0.0, 10.0, 2))  # 20833.3333333 
print(simpson(lambda x:x**4, 0.0, 10.0, 100000)) # 20000.0