1. 分別用三點和四點gauss-chebyshev公式計算積分
並與準確積分值2arctan4比較誤差。若用同樣的三點和四點gauss-legendre公式計算,也給出誤差比較結果。
2*atan(4)
ans =
2.6516
gauss-chebyshev:
function i = gausscheby(f,a,b,n)
syms t;
t= findsym(sym(f));
ta = (b-a)/2;
tb = (a+b)/2;
switch n
case 3
i=pi/n*ta*(subs(sym(f),t,ta*cos(pi/(2*n))+tb)*sqrt(1-cos(pi/(2*n))^2)+...
subs(sym(f),t,ta*cos(3*pi/(2*n))+tb)*sqrt(1-cos(3*pi/(2*n))^2)+...
subs(sym(f),t,ta*cos(5*pi/(2*n))+tb)*sqrt(1-cos(5*pi/(2*n))^2));
case 4
i=pi/n*ta*(subs(sym(f),t,ta*cos(pi/(2*n))+tb)*sqrt(1-cos(pi/(2*n))^2)+...
subs(sym(f),t,ta*cos(3*pi/(2*n))+tb)*sqrt(1-cos(3*pi/(2*n))^2)+...
subs(sym(f),t,ta*cos(5*pi/(2*n))+tb)*sqrt(1-cos(5*pi/(2*n))^2)+...
subs(sym(f),t,ta*cos(7*pi/(2*n))+tb)*sqrt(1-cos(7*pi/(2*n))^2));
endi=simplify(i);
i=vpa(i,6);
syms x
f=1/(1+x^2);
a=-4;b=4;
n=3;
y=gausscheby(f,a,b,n)
y =4.511
n=4:
y =1.90041
gauss-legendre:
function i = intgausslegen(f,a,b,n)
syms t;
t= findsym(sym(f));
ta = (b-a)/2;
tb = (a+b)/2;
switch n
case 0,
i=2*ta*subs(sym(f),t,tb);
case 1,
i=ta*(subs(sym(f),t,ta*0.5773503+tb)+...
subs(sym(f),t,-ta*0.5773503+tb));
case 2,
i=ta*(0.55555556*subs(sym(f),t,ta*0.7745967+tb)+...
0.55555556*subs(sym(f),t,-ta*0.7745967+tb)+...
0.88888889*subs(sym(f),t,tb));
case 3,
i=ta*(0.3478548*subs(sym(f),t,ta*0.8611363+tb)+...
0.3478548*subs(sym(f),t,-ta*0.8611363+tb)+...
0.6521452*subs(sym(f),t,ta*0.3398810+tb) +...
0.6521452*subs(sym(f),t,-ta*0.3398810+tb));
case 4,
i=ta*(0.2369269*subs(sym(f),t,ta*0.9061793+tb)+...
0.2369269*subs(sym(f),t,-ta*0.9061793+tb)+...
0.4786287*subs(sym(f),t,ta*0.5384693+tb) +...
0.4786287*subs(sym(f),t,-ta*0.5384693+tb)+...
0.5688889*subs(sym(f),t,tb));
case 5,
i=ta*(0.1713245*subs(sym(f),t,ta*0.9324695+tb)+...
0.1713245*subs(sym(f),t,-ta*0.9324695+tb)+...
0.3607616*subs(sym(f),t,ta*0.6612094+tb)+...
0.3607616*subs(sym(f),t,-ta*0.6612094+tb)+...
0.4679139*subs(sym(f),t,ta*0.2386292+tb)+...
0.4679139*subs(sym(f),t,-ta*0.2386292+tb));
endi=simplify(i);
i=vpa(i,6);
y =2.04798
n=4:
y =3.08862
2. 分別用三點和四點gauss-lagurre公式計算積分
function i = gausslagurre(f,n)
syms t;
t= findsym(sym(f));
switch n
case 2
i=0.7110930*subs(sym(f),t,0.4157746)+...
0.2785177*subs(sym(f),t,2.2942804)+...
0.0103893*subs(sym(f),t,6.2899451);
case 3
i=0.6031541*subs(sym(f),t,0.3225477)+...
0.3574187*subs(sym(f),t,1.7457611)+...
0.0388879*subs(sym(f),t,4.5366203) +...
0.0005393*subs(sym(f),t,9.3950710);
endi=simplify(i);
i=vpa(i,6);
syms x
f=exp(-10*x)*sin(x);
f=f./exp(-x);
a=0;b=inf;
n=2;
y= gausslagurre(f,n)
y =0.00680897
n=4:
y =0.0104892
3. 設,分別取,,用以下三個公式計算,
列表比較三個公式的計算誤差,從誤差可以得出什麼結論?
function [df1,df2,df3,w1,w2,w3]=midpoint(func,a)
if (nargin == 3 && h == 0.0)
disp('h不能為0');
return;
endfor k=1:6
h=1/10^k;
y0=subs(sym(func), findsym(sym(func)),a);
y1 = subs(sym(func), findsym(sym(func)),a+h);
y2 = subs(sym(func), findsym(sym(func)),a-h);
df1(k) = (y1-y0)/h;
df2(k) = (y1-y2)/(2*h);
y3=subs(sym(func), findsym(sym(func)),a+2*h);
y4=subs(sym(func), findsym(sym(func)),a-2*h);
df3(k)=(y4-8*y2+8*y1-y3)/(12*h);
w1(k)=1/a-df1(k);
w2(k)=1/a-df2(k);
w3(k)=1/a-df3(k);
enddf1=simplify(df1); df1=vpa(df1,6);
df2=simplify(df2); df2=vpa(df2,6);
df3=simplify(df3); df3=vpa(df3,6);
w1=simplify(w1); w1=vpa(w1,6);
w2=simplify(w2); w2=vpa(w2,6);
w3=simplify(w3); w3=vpa(w3,6);
syms x
f=log(x);
a=0.7;
[y1,y2,y3,w1,w2,w3]=midpoint(f,a)
y1 =
[ 1.33531, 1.41846, 1.42755, 1.42847, 1.42856, 1.42857]
y2 =
[ 1.43841, 1.42867, 1.42857, 1.42857, 1.42857, 1.42857]
y3 =
[ 1.42806, 1.42857, 1.42857, 1.42857, 1.42857, 1.42857]
w1 =
[ 0.0932575, 0.0101079, 0.00101944, 0.000102031, 0.0000102043, 0.00000101738]
w2 =
[ -0.00983893, -0.0000971936, -9.71814e-7, -9.72193e-9, 1.92131e-10, -3.02526e-9]
w3 =
[ 0.000513166, 4.76342e-8, 8.04334e-12, -1.10191e-10, 3.37991e-10, -4.5859e-9]
Matlab Romberg求積公式求積分
用romberg求積法計算積分 1 111 100x2d x int frac d x 11 1 100 x21 dx的近似值,要求誤差不超過0.5 1 0 70.5 times10 0.5 10 7。程式1 romberg通用程式 function result romberg a,b,f,eps...
高斯型求積公式 勒讓德 拉蓋爾 切比雪夫
gauss型求積公式 若機械求積公式具有階代數精度,則稱為gauss型求積公式,而在上關於權函式的次正交多項式的零點就是gauss型求積公式的gauss點。在gauss型求積公式中,若權函式,區間為,則公式為 特別的稱為gauss legendre公式。下表列出gauss legendre公式的結點...
高斯那些公式
貝葉斯基本公式與理解 的姊妹篇。多元正態分佈 multivariate normal distribution fx x1,xn 1 2 n 2 1 2 exp 1 2 x t 1 x 一維 f x 1 2 exp x 22 2 二維 f x,y 12 2exp x 2 y2 2 2 單高斯 非高斯...