计算方法C语言编程

1、 计算方法C语言编程1、已知方程在区间1，2内有一根，试问用二分法求根，使其具有5位有效数字至少应二分多少次？【程序设计】%-二分法-%clc;cleara=1;b=2;n=0;x=1.0;x0=0;while(abs(x-x0)>0.00005); x0=x; x=(a+b)/2 if f1(x)>0 b=x; else a=x; end n=n+1;endN=n运行结果x =1.3788N =152、用迭代法求的正根，要求准确到小数点后第5位。【程序设计】%-迭代法-%clc;clearx0=0.00000;x=1.50000;k=0;y=(log(x+0.20000)/5.0

2、000;x0=x;x=exp(y);while(abs(x-x0)>0.00005);y=(log(x+0.20000)/5.0000;x0=x;x=exp(y);k=k+1;end N=kX=vpa(x,6)%准确到小数点后第五位运行结果3、用牛顿法求方程在x0=2附近的根，要求准确到小数点后第3位。【程序设计】function f=f3(x)f=x3-3*x-1;%-牛顿法-%clc;clearx0=1.0;x=2.0;k=0;if diff(f3(x0)=0 breakendwhile(abs(x-x0)>0.0005); x0=x; x=x-(x3-3*x-1)/(3*x2

3、-3); k=k+1;endX=vpa(x,4),%精确到小数点后第三位N=k运行结果 4、分别用单点和双点弦截法求方程在1,1.5内的根。要求|xn+1-xn|0.000005【程序设计】%-单点弦截法-%clc;clearx1=1.0;x2=2.0;k=0;while abs(x2-x1)>0.000005; t=x2-(x22-x2-1)*(x2-x1)/(x23-x2-1)-(x13-x1-1); x1=x1; x2=t; k=k+1;endX=k;Y=x2 %-双点弦截法-%clc;clearx1=1.0;x2=2.0;k=0;while abs(x2-x1)>0.000

4、005; t=x2-(x22-x2-1)*(x2-x1)/(x23-x2-1)-(x13-x1-1); x1=x2; x2=t; k=k+1; endX=k,Y=x2 运行结果X = 28Y = 1.61805、分别用列主元素消去法求解下列方程组。(计算取4位小数)【程序设计】%-列主元素消去法-%clc;clearA=1.1161 0.1254 0.1397 0.1490; 0.1582 1.1675 0.1768 0.1871;0.1968 0.2071 1.2168 0.2271;0.2368 0.2471 0.2568 1.26717;b=1.5471;1.6471;1.7471;1.

5、8471;B=A,b; n=length(b);RA=rank(A);RB=rank(B);d=RB-RA;if d>0 disp('此方程组无解')endif RA=RB&RA=n disp('此方程组有无穷解')endif RA=RB&RA=n disp('此方程组有唯一解') n,n=size(A); x=zeros(n,1);for k = 1:n-1 piv,r = max(abs(B(k:n,k); %找列主元所在子矩阵的行r r = r + k - 1; % 列主元所在大矩阵的行 if r>k temp=

6、B(k,:); B(k,:)=B(r,:); B(r,:)=temp; endend if B(k,k)=0, error('对角元出现0'), end % 把增广矩阵消元成为上三角 for p = k+1:n B(p,:)=B(p,:)-B(k,:)*B(p,k)/B(k,k); end end % 解上三角方程组 A =B(:,1:n); b = B(:,n+1); x(n) = b(n)/A(n,n); for k = n-1:-1:1 x(k)=b(k); for p=n:-1:k+1 x(k) = x(k)-A(k,p)*x(p); end x(k)=x(k)/A(k

7、,k); endx运行结果此方程组有唯一解x= 0.9568 1.03331.2095 1.21256、设有方程组 ，取初始向量,用雅可比迭代法求解,要求时迭代终止。【程序设计】%-雅克比迭代法-%clc;clearA=5 2 1;-1 1 1;2 -3 10;b=-12; 20; 3;x0=-3;1;1;ep=3;D=diag(diag(A); %对角阵L=-tril(A,-1); %下三角U=-triu(A,1); %上三角B=D(L+U); f=Db;x=B*x0+f;n=1;while norm(x-x0)>=ep x0=x; x=B*x0+f; n=n+1;endn,x运行结果

8、n = 10x = -6.4468 10.8406 3.9510 7、设有方程组(1) 证明解此方程组的雅可比迭代法收敛,而相应的赛德尔迭代法发散.(2) 取初始向量,用雅可比迭代法求解,要求迭代三次.【程序设计】%-雅克比迭代法 迭代三次-%clc;clearx1=zeros(1,10);x2=zeros(1,10);x3=zeros(1,10);X1=0;X2=0;X3=0;for k=1:3; x1(k+1)=1-2*x2(k)+2*x3(k); x2(k+1)=3-x1(k)-x3(k); x3(k+1)=5-2*x1(k)-2*x2(k);endx=x1(4) x2(4) x3(4)

9、运行结果x = 1 1 18、设有方程组 ，其等价形式为(1) 证明解等价方程组的简单迭代法发散,而赛德尔法收敛(2) 取初始向量,用赛德尔迭代法求解,要求迭代四次【程序设计】%-赛德尔迭代法-迭代四次-%clc;clearx1=zeros(1,10);x2=zeros(1,10);x3=zeros(1,10);X1=0.0;X2=0.0;X3=0.0;for k=1:4; x1(k+1)=(x2(k)+x3(k)/2; x2(k+1)=(x1(k+1)+1)/2; x3(k+1)=x1(k+1);endx=x1(5) x2(5) x3(5)运行结果 x =0.5781 0.7891 0.57

10、81x1.12751.15031.17351.1972Y=f(x)0.11910.139540.159320.179039、已知函数表：应用拉格朗日插值公式计算f(1.1300)的近似值。（计算取4位小数）【程序设计】%-拉格朗日插值-%clc;clearn=4;x=1.1300;X=1.1275;1.1503;1.1735;1.1972;Y=0.1191;0.1395;0.1593;0.1790;y=0;p=1;for k=1:n for j=1:n if j=k p=p.*(x-X(j)/(X(k)-X(j); end end y=y+p*Y(k);endx,y=vpa(y,4)%取四位小

11、数运行结果x =1.1300y =0.124710、利用函数表造出差商表，并利用牛顿插值公式计算f(x)在x=1.682,1.813处的近似值.（计算取5位小数）.【程序设计】%-牛顿插值法-%clc;clearx=1.615,1.634,1.072,1.828,1.92;y=2.41450,2.46459 ,2.65271,3.03035,3.34066;x1=1.68200;x2=1.81300;A=zeros(5,4);A=y',A;for j=2:5 for i=j:5 A(i,j)=(A(i,j-1)-A(i-1,j-1)/(x(i)-x(i+1-j); endenddisp

12、('chashangbiao:');A=x',A;disp(A); %计算f(x1),f(x2)的近似值N=y(1);yd=0;New=1;for i=1:4 for j=i+1:5 yd(j)=(y(j)-y(i)/(x(j)-x(i); endc(i)=yd(i+1);New1=New*(x1-x(i);N1=N+c(i)*New1;New2=New*(x2-x(i);N2=N+c(i)*New2;endN1=vpa(N1,6),N2=vpa(N2,6)%取五位小数运行结果11、已知函数f(x)的下列数据:用柯特斯公式计算积分.(计算取5位小数)xk1.82.02.

13、22.42.6f(xk)3.120144.426596.042418.0301410.46675【程序设计】 %-科特斯公式求积分-%function y=KTS(X,Y)clc;cleara=1.8;b=2.6;X=1.8;2.0;2.2;2.4;2.6;Y=3.12014;4.42659;6.04241;8.03014;10.46675;y=(1.0/90)*(b-a)*(7*Y(1)+32*Y(2)+12*Y(3)+32*Y(4)+7*Y(5);Y=vpa(y,6)%取五位小数运行结果xk00.10.20.30.40.5f(xk)11.0049711.0195361.0426681.07

14、27071.107432xk0.60.70.80.91.0f(xk)1.1441571.7598591.2113071.2352111.24837512、已知函数f(x)的下列数据:用复化梯形公式和复化辛浦生公式计算积分.(计算取6位小数)注意:本题f(x)=cosx+sin²x,积分真值I=1.11414677.【程序设计】function y=fly(x)y=sin(x)*sin(x)+cos(x)%-复化梯形公式、复化辛普森公式-%clc;clearh=0.1;n=11;m=0;l=0;X=0;0.1;0.2;0.3;0.4;0.5;0.6;0.7;0.8;0.9;1.0;f=

15、1;1.004791;1.019536;1.042668;1.072707;1.107432;1.44157;1.759859;1.211307;1.235211;1.248375;for k=1:n-1 m=m+f(k);l=l+fly(0.1*k+0.5*h); y1=0.5*h*(f(1)+2*m+f(n);%复化梯形公式 y2=(1.0/6)*h*(f(1)+4*l+2*m+f(10);%复化辛普森公式endy1=vpa(y1,7),y2=vpa(y2,7)%取六位小数运行结果13、用龙贝格法计算积分I=,要求误差不超过(计算取6位小数).【程序设计】function y=fx13(x

16、)y=exp(-x);lc;clear;format long e;%精确显示y=18; T=zeros(y,y);%定义长度为y的矩阵Tfor i=0:y-1 T(i+1)=i;%对矩阵T的第一行赋初值enda=0;b=1;T(1,2)=(b-a)*(fx13(a)+fx13(b)/2;%算出T1for l=2:y sum=0; for j=1:2(l-1) sum=sum+fx13(a+(2*j-1)*(b-a)/2l); end T(l,2)=T(l-1,2)/2+(b-a)/2l)*sum; %算出梯形序列 if (T(l,2)-T(l-1,2)<0.0000005; %判断是否

17、符合精度 break; endendh=1;for m=1:y-2 for k=h:(y+h-m-1) T(k+1,m+2)=(4m)*T(k+1,m+1)-T(k,m+1)/(4m-1); %按公式算出接下去的序列 if (T(k+1,m+2)-T(k,m+2)<0.0000005; break; end end h=h+1;endI=T(k+1,m+2)*2/sqrt(pi)I=vpa(I,7)%取六位小数运行结果14、用梯形法与欧拉预估-校正法求解初值问题: 取步长h=0.1,并与准确解相比较.(计算取6位小数).【程序设计】function f=f(x,y)f=x+y;%-梯形法-%clc;clearh=0.1;a=0;b=1;n=(b-a)/h;y0=1;x=a:h:b;y(1)=y0;for i=1:n y(i+1)=