#Trapezoidal rule # โจทย์ f(x) = e^x*sin(x) ,n=10 from math import e,sin def f(x) : return e**x*sin(x) #n คือ จำนวนกรแบ่งช่วง #h คือ ค่าของช่วงที่แบ่ง from (b-a)/n def trap(f,a,h,n) : return sum([ (1/2)*h*(f(a+i*h)+f(a+(i+1)*h)) for i in range(n) ]) def area(f,a,b,n): return sum([ (1/2)*((b-a)/n)*(f(a+i*((b-a)/n))+f(a+(i+1)*((b-a)/n))) for i in range(n) ]) trap(f,0.1,0.05,10),area(f,0.1,0.6,10) ............................................................................................. #หา integrate จาก Sympy x,y=sympy.Symbol('x'),sympy.Symbol('y') y=3*x**2+9*x+2 z=sympy.integrate(y,(x,0.1,0.6)) from math import e zz = e**-x zzz = sympy.integrate(zz,(x,0,1)) y,z,zzz #Simpson's rule import sympy from math import e f =e**(-x) sympy.integrate(f,(x,0,1)) #print 0.632120558828558 def simpson(f, a, b, n): h = (b-a)/n x0 = a # x_i = x_0 + i*h #f(x0) + 4*f(x1) + 2*f(x2) + .... + 4*f(x_n-1) + f(x_n) s = f(x0) # odd - x1, x3, x5, ... #s += 4*f(x1) + 4*f(x3) + 4*f(x5) + ... for i in range(1,n,2): s += 4*f(x0+i*h) # even - x2, x4, x6, ... #s += 2*f(x2) + 2*f(x4) + 2*f(x6) + ... for i in range(2,n-1,2): s += 2*f(x0+i*h) s += f(x0+n*h) return s*h/3 from math import e from sympy import * def f(x): return e**(-x) a,b = 0,1 r=0.632120558828558 simpson(f, a, b, 4) for v in [ (n,simpson(f,a,b,n),area(f,a,b,n)) for n in range(4, 22, 2) ] : print ('{:5} {:20} {:20}'.format(v[0],abs(v[1]-r),v[2]))