matlab数值分析程序--高等数学,数值代数的matlab实现-文字版, matlab电子书, 和matlab 有关的电子书:

9.3.8 利用拉格朗日插值等方法求表格型数值积分的MATLAB方法

9.3.8  利用拉格朗日插值等方法求表格型数值积分的MATLAB方法

例9.3.19  测得定积分d中被积函数在积分变量的某几个特定点X =0,0.157 1, 0.471 2, 0.549 8, 0.628 3, 0.863 9, 1.021 0, 1.099 6, 1.570 8处的值分别为Y= 0, 0.156 5, 0.454 0, 0.522 5, 0.587 8, 0.760 4, 0.852 6, 0.891 0, 1.000 0,试根据这些值首先分别用分段二、三阶拉格朗日插值多项式和三次样条插值函数构造被积函数,然后分别用quad函数和quad8函数计算,取精度为,并将计算结果与精确值1进行比较.

解 方法1   (1) 拉格朗日插值多项式的MATLAB主程序

function L=lfun(X,Y)

m=length(X); n=M1; L=ones(m,m);

for k=1:n+1

    V=1;

    for i=1:n+1

     if k~=i

        V=conv(V,poly(X(i)))/(X(k)-X(i));

     end

end

l(k,:)=poly2sym (V);

end

L=Y* l;

(2) 用二阶拉格朗日插值多项式构造被积函数. 输入程序

>>  X =[0,0.1571, 0.4712,0.5498,0.6283,0.8639, 1.0210,1.0996,1.5708];

Y=[0,0.1565,0.4540,0.5225,0.5878,0.7604,0.8526,0.8910,1.0000];

L1=lfun(X(1:3),Y(1:3)), L2=lfun(X(3:5),Y(3:5)),

L3=lfun(X(5:7),Y(5:7)), L4=lfun(X(7:9),Y(7:9)),

运行后输出多项式L1,L2,L3,L4.

(3)输入程序

>> syms x

L1=inline('-14644281526214207/140737488355328000*x.^2+2280009553133670687/2251799813685248000*x');

L2=inline('-88810055001957943/351843720888320000*x.^2+15892292661979563/14073748835532800*x-15511422992738729/703687441776640000');

L3=inline('-261101353037957127/703687441776640000*x.^2+362054074498830351/281474976710656000*x-207985426949036897/2814749767106560000');

L4=inline('-131688789655946847/281474976710656000*x.^2+104193434945799191/70368744177664000*x-9309740592764125768759/54295819320043765760000');

 [Q41,FCNT41] = quad(L1,0, 0.4712,1.e-4);

[Q42,FCNT42] =quad(L2, 0.4712, 0.6283,1.e-4);

[Q43,FCNT43] = quad(L3, 0.6283, 1.0210,1.e-4);

[Q44,FCNT44] = quad(L4, 1.0210, pi/2,1.e-4);

Q= Q41+ Q42+ Q43+ Q44,

FCNT= FCNT41+ FCNT42+ FCNT43+ FCNT44

[Q841,FCNT841] = quad8(L1,0, 0.4712,1.e-4);

[Q842,FCNT842] =quad8(L2, 0.4712, 0.6283,1.e-4);

[Q843,FCNT843] = quad8(L3, 0.6283, 1.0210,1.e-4);

[Q844,FCNT844] = quad8(L4, 1.0210, pi/2,1.e-4);

Q8= Q841+ Q842+ Q843+ Q844,

FCNT8= FCNT841+ FCNT842+ FCNT843+ FCNT844

运行后屏幕显示分别用quad函数和quad8函数计算的近似值Q、Q8和迭代次数FCNT、FCNT8,取精度为如下

Q =                  FCNT =        Q8 =            FCNT8 =

0.99957595925413           52       0.99957595925413          132

方法2  (1) 用三阶拉格朗日插值多项式构造被积函数.输入程序

>>X=[0,0.1571,0.4712,0.5498,0.6283,0.8639,1.0210,1.0996,1.5708];

Y=[0,0.1565,0.4540,0.5225,0.5878,0.7604,0.8526,0.8910,1.0000];

L1=lfun(X(1:4),Y(1:4)), L2=lfun(X(4:7),Y(4:7)),

L3=lfun(X(7:9),Y(7:9)),

运行后输出多项式L1,L2,L3.

 (2)输入程序

>> syms x

L1=inline('-9070382563990323/56294995342131200*x.^3-793956227651389/281474976710656000*x.^2+2253151961697736691/2251799813685248000*x');

L2=inline('-41118170146065463/351843720888320000*x.^3-1700780706742359/21990232555520000*x.^2+735208215754549839/703687441776640000*x-25688208496779769/2814749767106560000');

L3 =

inline('-131688789655946847/281474976710656000*x.^2+104193434945799191/70368744177664000*x-9309740592764125768759/54295819320043765760000');

 [Q41,FCNT41] = quad(L1,0, 0.5498,1.e-4);

[Q42,FCNT42] =quad(L2, 0.5498, 1.0210,1.e-4);

[Q43,FCNT43] = quad(L3, 1.0210, pi/2,1.e-4);

Q= Q41+ Q42+ Q43,

FCNT= FCNT41+ FCNT42+ FCNT43

[Q841,FCNT841] = quad8(L1,0, 0.5498,1.e-4,2);

[Q842,FCNT842] =quad8(L2, 0.5498, 1.0210,1.e-4,2);

[Q843,FCNT843] = quad8(L3, 1.0210, pi/2,1.e-4,2);

Q8= Q841+ Q842+ Q843,

FCNT8= FCNT841+ FCNT842+ FCNT843

运行后屏幕显示分别用quad函数和quad8函数计算的近似值Q、Q8和迭代次数FCNT、FCNT8,取精度为如下

Q =                    FCNT=                         

0.99975250938584         39          

Q84=                   FCNT8=

0.99975250938584          99

方法3  三次样条插值多项式计算数值积分.输入程序

>> a=0,b=1.5708,

X=[0,0.1571,0.4712,0.5498,0.6283,0.8639,1.0210,1.0996,1.5708];

Y=[0,0.1565,0.4540,0.5225,0.5878,0.7604,0.8526,0.8910,1.0000];

pp = spline(X,Y)

[Q8,FCNT8]= quad8(@ppval,a,b, ,1.e-4,[],pp)

[Q,FCNT]= quad(@ppval,a,b,[],[],pp)

运行后屏幕显示      

Q8 =                   FCNT8 =     Q =                 FCNT =

        1.00012965800727        33       1.00012956660843          17

欢迎转载,转载请注明来自一手册:http://yishouce.com/book/3/3090308.html
友情链接It题库(ittiku.com)| 版权归yishouce.com所有| 友链等可联系 admin#yishouce.com|粤ICP备16001685号-1