数值分析课程设计

1、课程设计（论文）任务书学 院理 学 院专 业信息与计算科学学生姓名班级学号课程名称数值分析课程设计课 程 设 计（论文）题目数值分析算法设计与应用案例设计要求（技术参数）：1、熟练掌握Matlab软件。2、运用Matlab软件进行各种数值算法编程。3、每个同学至少做四个算法编程实验，每种算法最少有一个实例。 4、每种算法在关键语句中要有注释。设计任务：实验一：LU分解法解线性方程组实验二：Lagrange插值法数值求解实验三：龙贝格求积公式求数值积分实验四：用Runge-Kutta方法求常微分方程数值解实验五：线性规划任务分配问题实验六：求电阻问题计划与进度安排：第一周：查阅资料 学习使用Ma

2、tlab软件 完成实验一、实验二第二周：完成实验三 、实验四、实验五 攥写课程设计报告 答辩成绩：指导教师（签字）：年 月 日专业负责人（签字）：年 月 日主管院长（签字） 年 月 日目 录实验一 LU分解法解线性方程组11.1实验目的与要求11.2实验基本原理11.3算法程序21.4实例及结果分析3实验二 Lagrange插值法数值求解42.1实验目的与要求42.2实验基本原理42.3算法程序52.4实例及结果分析5实验三 龙贝格求积公式求数值积分63.1实验目的与要求63.2实验基本原理63.3算法程序63.4实例及结果分析7实验四 用Runge-Kutta方法求常微分方程数值解84.1实

3、验目的与要求84.2实验基本原理84.3算法程序94.4实例及结果分析9实验五 线性规划任务分配问题115.1问题提出115.2模型分析115.3模型建立125.4模型求解12实验六 求电阻问题136.1问题提出136.2问题分析与模型建立135.3模型求解146.4结果分析15实验一 LU分解法解线性方程组1.1实验目的与要求1.了解LU分解法解线性方程组的基本原理；2.熟悉计算方法的技巧和过程，能用LU分解法解实际问题；3.用matlab实现LU分解。1.2实验基本原理1.若一个线性方程组系数矩阵为n阶方阵A且各阶顺序主子式均不为0则A的LU分解存在且唯一。2.在满足1的条件下课推导得出以

4、下公式（1） （2） 3.公式（1）用于求解矩阵L、U，公式（2）用于会带求解y、x。从公式中可以看出：L对角线上元素为1，U第一行与A第一行相同。4LU分解的具体过程和顺序如下：（1）第一步分解：（2）第二步分解： （3）第三步分解： （n）第n步分解：依次计算：、.，. 1.3算法程序function L,U,x=Lu_x(A,d)n,m=size(A);if n=merror(The rows and columns of matrix A must be equal!);return;endfor ii=1:nfor i=1:iifor j=1:iiAA(i,j)=A(i,j);end

5、endif (det(AA)=0)error(The matrix can not be divided by LU!)return;endendA n,n=size(A); L=zeros(n,n); U=zeros(n,n); for i=1:n L(i,i)=1; end for k=1:n for j=k:n U(k,j)=A(k,j)-sum(L(k,1:k-1).*U(1:k-1,j); end for i=k+1:n L(i,k)=(A(i,k)-sum(L(i,1:k-1).*U(1:k-1,k)/U(k,k); end end y(1)=d(1);for i=2:n for

6、j=1:i-1d(i)=d(i)-L(i,j)*y(j);endy(i)=d(i);endx(n)=y(n)/U(n,n);for i=(n-1):-1:1for j=n:-1:i+1y(i)=y(i)-U(i,j)*x(j);endx(i)=y(i)/U(i,i);end1.4实例及结果分析 MATLAB命令窗体输入如下： A=3 5 6 2;9 8 7 5;2 7 1 6;3 12 5 9; b=18 10 16 17; L,U,x=Lu_x(A,b)得到的结果如下：A = 3 5 6 2 9 8 7 5 2 7 1 6 3 12 5 9L = 1.0000 0 0 0 3.0000 1.

7、0000 0 0 0.6667 -0.5238 1.0000 0 1.0000 -1.0000 1.3696 1.0000U = 3.0000 5.0000 6.0000 2.0000 0 -7.0000 -11.0000 -1.0000 0 0 -8.7619 4.1429 0 0 0 0.3261x =4.7500 54.2500 -25.2500 -58.0000实验二 Lagrange插值法数值求解2.1实验目的与要求1.了解拉格朗日插值的基本概念；2.了解插值公式的基本原理，运用插值公式求解实际问题；3.编写matlab程序，实现拉格朗日插值法，观察所得结果的精确性。2.2实验基本原

8、理设已知，.,及=f()(i=0,1,.,n),为不超过n次多项式且满足（i=0,1,.n）.易知其中，均为n次多项式,再由（ji）为n次多项式的n个根知.最后，由,i=0,1,.,n.总之，=，=式为n阶Lagrange插值公式，其中，（i=0,1,.n）称为n阶Lagrange插值的基函数。2.3算法程序function yy=nalagr(x,y,xx)m=length(x);n=length(y);if m=n,error( 向量x与y的长度必须一致);ends=0;for i=1:n t=ones(1,length(xx); for j=1:n if j=i t=t.*(xx-x(j

9、)/(x(i)-x(j); end end s=s+t*y(i);endyy=s;2.4实例及结果分析给出的数值表，如下图所示，用Lagrange插值计算的近似值。0.40.50.60.70.8-0.91629-0.69315-0.51083-0.35777-0.22314MATLAB命令窗体输入如下： x=0.4 0.5 0.6 0.7 0.8; y=-0.91629 -0.69315 -0.51083 -0.35777 -0.22314 xx=0.54; nalagr(x,y,xx)得到的结果如下：ans = -0.6160用Lagrange插值法求得的近似值为-0.6160，与计算机求得

10、的值-0.6162，比较，是极为相近的。实验三 龙贝格求积公式求数值积分3.1实验目的与要求1. 掌握龙贝格求积的基本思路和步骤；2. 培养编程与上机调试能力。3.2实验基本原理（1） 取，,求。令（记为区间的二分次数）。（2） 求梯形值，即按递推公式计算。（3） 求加速值，按公式 （=1,2）逐个求出T表的第行其余各元素（=1,2）。（4） 若（预先给出的精度），则终止计算，并取；否则令转（2继续计算）。 排成三角数表3.3算法程序function R = romberg(f, a, b, e)% 参数介绍:% f - 被积函数f(x)% a - x的左区间.% b - x的右区间.% e

11、- 误差限.% 结果:% R - 返回Romberg 表. n = 1; %区间二分次数while 1 %在此仅代表多次二分，在后面判断循环终止 R = zeros(n + 1, n + 1);%生成(n+1)*(n+1)的0矩阵 R(0+1, 0+1) = (b - a) / 2 * (feval(f, a) + feval(f, b); % 初始值(2点梯形公式). for i = 1 : n % 按照公式计算Romberg 表的第一列. h = (b - a) / 2i; s = 0; for k = 1 : 2(i-1) s = s + feval(f, a + (2*k - 1)*h

12、); end R(i+1, 0+1) = R(i-1+1, 0+1)/2 + h*s; end for j = 1 : n % 计算Romberg 表的其他列. fac = 1 / (4j - 1); for m = j : n R(m+1, j+1) = R(m+1, j-1+1) + fac*(R(m+1, j-1+1) - R(m-1+1, j-1+1); end end if abs(R(n,n) - R(n+1,n+1) fun=inline(4./(1+x2),x) ; romberg(fun,0,1,1e-6)得到的结果如下：t = 3.1416ans = 3.0000 0 0

13、0 0 0 3.1000 3.1333 0 0 0 0 3.1312 3.1416 3.1421 0 0 0 3.1390 3.1416 3.1416 3.1416 0 0 3.1409 3.1416 3.1416 3.1416 3.1416 0 3.1414 3.1416 3.1416 3.1416 3.1416 3.1416实验四 用Runge-Kutta方法求常微分方程数值解4.1实验目的与要求（1）掌握Runge-Kutta方法的基本原理；（2）用matlab编写4阶Runge-Kutta程序；（3）用写好的程序求解常微分方程，并与系统自带的ode45函数进行比较。4.2实验基本原理对

14、于一阶精度的欧拉公式有： 当用点处的斜率近似值与右端点处的斜率的算术平均值作为平均斜率K*的近似值，那么就会得到二阶精度的改进欧拉公式： 依次类推，如果在区间内多预估几个点上的斜率值、.，并用他们的加权平均数作为平均斜率K*的近似值，显然能构造出具有很高精度的高阶计算公式。经数学推导、求解，可以得出四阶龙格库塔公式，也就是在工程中应用广泛的经典龙格库塔算法： 4.3算法程序建立lgkt4j.m文件；function y,x=lgkt4j(x0,xn,y0,h) x=x0:h:xn;n=length(x);y1=x;y1(1)=y0;for i=1:n-1 K1=f(x(i),y1(i); K2

15、=f(x(i)+h/2,y1(i)+h/2*K1); K3=f(x(i)+h/2,y1(i)+h/2*K2); K4=f(x(i)+h,y1(i)+h*K3); y1(i+1)=y1(i)+h/6*(K1+2*K2+2*K3+K4);end y=y1;4.4实例及结果分析求微分方程，在时的数值解。（1） 首先建立f.m文件； function z=f(x,y) z=-2*x*y2;（2） 在MATLAB命令窗体输入如下； x0=0;xn=1.2;y0=1;h=0.1; y,x=lgkt4j(x0,xn,y0,h); n=length(x); fprintf(自编的龙格库塔函数效果n); fpr

16、intf( i x(i) y(i)n); for i=1:n fprintf(%2d %12.4f %12.4fn,i,x(i),y(i); end plot(x,y,ro-); x1,y1=ode45(f,0,1,1);fprintf(ode45效果n); fprintf( i x1(i) y1(i)n); for i=1:n fprintf(%2d %12.4f %12.4fn,i,x(i),y(i); end plot(x1,y1,bv-);运行结果为：自编的龙格库塔函数效果 i x(i) y(i) 1 0.0000 1.0000 2 0.1000 0.9901 3 0.2000 0.9

17、615 4 0.3000 0.9174 5 0.4000 0.8621 6 0.5000 0.8000 7 0.6000 0.7353 8 0.7000 0.6711 9 0.8000 0.609810 0.9000 0.552511 1.0000 0.500012 1.1000 0.452513 1.2000 0.4098ode45效果 i x1(i) y1(i) 1 0.0000 1.0000 2 0.1000 0.9901 3 0.2000 0.9615 4 0.3000 0.9174 5 0.4000 0.8621 6 0.5000 0.8000 7 0.6000 0.7353 8 0

18、.7000 0.6711 9 0.8000 0.609810 0.9000 0.552511 1.0000 0.500012 1.1000 0.452513 1.2000 0.4098上文分别用自编的龙格贝塔函数和matlab自带的ode45函数求常微分函数的数值解，从数值解和图对比看出，用自带的龙格库塔函数求的值是非常接近的。实验五 线性规划任务分配问题5.1问题提出某车间有甲、乙两台机床，可用于加工三种工件.假定这两台车床的可用台时数分别为800和900，三种工件的数量分别为400、600和500，且已知用三种不同车床加工单位数量不同工件所需的台时数和加工费用如下表.问怎样分配车床的加工任

19、务，才能既满足加工工件的要求，又使加工费用最低？车床类 型单位工件所需加工台时数单位工件的加工费用可用台时数工件1工件2工件3工件1工件2工件3甲0.41.11.013910800乙0.51.21.3111289005.2模型分析设在甲车床上加工工件1、2、3的数量分别为、，在乙车床上加工工件1、2、3的数量分别为、。约束条件如下：（1）3种工件的数量分别为400、600、500，则有约束条件、；（2）甲车床的使用时间必须小于800，则；（3）乙车床的使用时间必须小于900，则；（4）非负约束：。此题为线性规划模型，可用到matlab中的优化工具箱解线性规划；模型： s.t 命令：x=linprog（c，A，b，Aeq,beq, VLB，VUB）5.3模型建立5.4模型求解改写为：s.t ，在命令窗口输入命令如下： f = 13 9 10 11 12 8; A = 0.4 1.1 1 0 0 0;0 0 0 0.5 1.2 1.3 ;b = 800; 900 ;Aeq=1 0 0 1 0 0;0 1 0 0 1 0;0 0 1 0 0 1;beq=400 600 500; vlb = zeros(6,1); vub=; x,fval = linprog(f,A,b,Aeq,beq,vlb,vub)得到结果：x = 0.0000 600.0000 0.0