1. (25分)计算积分
, n=0,1,2,…,20
若用下列两种算法
(A)
(B)
试依据积分In的性质及数值结果说明何种算法更合理。
代码实现(C语言):
#include <stdio.h>
#include <math.h> void main(){ double i=0.18232155679395;//ln6-ln5 printf("算法A:\n"); printf("I(0)=0.18232155679395\n"); for(int n=1;n<=20;n++) { printf("I(%d)=%lf\n",n,i); i=-5*i+pow(n,-1); } double j=0.007997523028232166; printf("\n算法B:\n"); for(int n=20;n>=0;n--) { printf("I(%d)=%lf\n",n,j); j=0.2*(pow(n,-1)-j); } }
运行结果:
分析:算法A每向前推进一步,其计算值的舍入误差便增长5倍;而算法B每向后推进一步,其舍入误差便减少5倍。因此算法B更合理
2. (25分)求解方程f(x)=0有如下牛顿迭代公式
, n≥1,x0给定
(1) 编制上述算法的通用程序,并以(ε为预定精度)作为终止迭代的准则。
(2) 利用牛顿迭代求解方程
设预定精度ε=10-10。
代码实现(MATLAB):
newton.m:
% newton: Newton迭代法
function [x1,k] = newton(f,f_de,varible,x0,eplison) dx = -subs(f_de,varible,x0)\subs(f,varible,x0); x1=x0+dx; k=1; while norm(x1-x0)>=eplison x0=x1; dx=-subs(f_de,varible,x0)\subs(f,varible,x0); x1=x0+dx; k=k+1; end
end
上机作业2.m:
format short
syms x y
f1=5*cos(x)-x*y-3;
f2=y^2+8*x*y-7;
F=[f1;f2];
F_de=[diff(F,'x'),diff(F,'y')];
eplison=10e-10;
variable=[x;y];
X0=[1;1]; [x,k]=newton(F,F_de,variable,X0,eplison);
fprintf('root:\n', k);
for i=1:length(x) fprintf('%1.6f\n',x(i));
end
运行结果:
>> 上机作业2
root:
0.724712
1.025858
3. (25分) 已知,
(1) 利用插值节点x0=1.00,x1=1.02,x2=1.04,x3=1.06,构造三次Lagrange插值公式,由此计算f(1.03)的近似值,并给出其实际误差;
(2) 利用插值节点x0=1,x1=1.05构造三次Hermite插值公式,由此计算f(1.03)的近似值,并给出其实际误差。
(1)代码实现(C语言):
#include <stdio.h>
#include <math.h> double f(double x)
{ double result =((x-1.02)*(x-1.04)*(x-1.06))/((1.00-1.02) *(1.00-1.04)*(1.00-1.06))*0.765789386+ ((x-1.00)*(x-1.04)*(x-1.06))/((1.02-1.00) *(1.02-1.04)*(1.02-1.06))*0.795366778+ ((x-1.00)*(x-1.02)*(x-1.06))/((1.04-1.00) *(1.04-1.02)*(1.04-1.06))*0.82268817+ ((x-1.00)*(x-1.02)*(x-1.04))/((1.06-1.00) *(1.06-1.02)*(1.06-1.04))*0.84752258; return result;
} void main()
{ double x; printf("please input x: "); scanf("%lf",&x); printf("f(%lf)=%lf\n",x,f(x));
}
运行结果:
实际误差:0.8093236189017053-0.809324=-3.8109829469945566e-7
(2) 代码实现(C语言):
#include <stdio.h>
#include <math.h> double L1(double x){ return (x-1.02)*(x-1.05)/((1.00-1.02)*(1-1.05));
}
double L2(double x){ return (x-1.00)*(x-1.05)/((1.02-1.00)*(1.02-1.05));
}
double L3(double x){ return (x-1.00)*(x-1.02)/((1.05-1.00)*(1.05-1.02));
}
double a1(double x){ return (1-2.0*(x-1.00)*(-70.0))*pow(L1(x),2.0);
}
double a2(double x){ return (1-2.0*(x-1.02)*(50.0/3.0))*pow(L2(x),2.0);
}
double a3(double x){ return (1-2.0*(x-1.05)*(160/3))*pow(L3(x),2.0);
}
double belta1(double x){ return (x-1.00)*pow(L1(x),2.0);
}
double belta2(double x){ return (x-1.02)*pow(L2(x),2.0);
}
double belta3(double x){ return (x-1.05)*pow(L3(x),2.0);
}
double H3(double x){ double y1,y2,y3,derivative_y1,derivative_y2,derivative_y3; y1=0.7657893864464859; y2=0.7953667788517541; y3=0.8354311093313167; derivative_y1=1.5315787728929717; derivative_y2=1.4243418718656515; derivative_y3=1.242214550953158; double result = y1*a1(x)+y2*a2(x)+y3*a3(x)+derivative_y1*belta1(x)+derivative_y2*belta2(x)+derivative_y3*belta3(x); return result;
} void main(){ double x; printf("please input x: "); scanf("%lf",&x); printf("f(%lf)=%lf\n",x,H3(x));
}
运行结果:
实际误差:0.8093236189017053-0.808878=0.000445618901705358
4. (25分) 利用复合求积算法计算积分
取步长为10-4
代码实现(MATLAB):
algorithm_1.m:
%% algorithm_1: 复合梯形公式
function F = algorithm_1(f,a,b,h) N=(b-a)/h; f_sum=0; for n=1:N-1 x_n=a+n*h; f_sum=f_sum+f(x_n); end F=h/2*(f(a)+2*f_sum+f(b));
end
algorithm_2.m:
%% algorithm_2: 复合辛普森公式
function F = algorithm_2(f,a,b,h) N=(b-a)/h; f_sum_half=0; f_sum=0; x_n=a; for n=1:N-1 x_n_half=x_n+h/2; x_n=a+n*h; f_sum_half=f_sum_half+f(x_n_half); f_sum=f_sum+f(x_n); end F=h/6*(f(a)+4*f_sum_half+2*f_sum+f(b));
end
上机作业4.m:
syms x
g(x)=sqrt(1+cos(x)^2);
a=0;
b=48;
% h=10e-4;
h=0.1;
F_algorithm_1=algorithm_1(g,a,b,h);
fprintf('algorithm_1: %1.6f\n', F_algorithm_1);
F_algorithm_2=algorithm_2(g,a,b,h);
fprintf('algorithm_2: %1.6f\n', F_algorithm_2);
运行结果:
>> 上机作业4
algorithm_1: 58. 470466
algorithm_2: 58. 462540