上海交通大学计算方法大作业

计算方法A上机报告学院（系）：电气工程学院学生姓名：陶然学号：授课老师：完成日期：2019年12月03日西安交通大学Xi'an Jiaotong University目录1 QR分解法求解线性方程组 (2)1.1 算法原理 (2)1.1.1 基于吉文斯变换的QR分解 (2)1.1.2 基于豪斯霍尔德变换的QR分解 (3)1.2 程序流程图 (4)1.2.1 基于吉文斯变换的QR分解流程图 (4)1.2.2 基于豪斯霍尔德变换的QR分解流程图 (5)1.3 程序使用说明 (5)1.3.1 基于吉文斯变换的QR分解程序说明 (5)1.3.2 基于豪斯霍尔德变换的QR分解程序说明 (7)1.4 算例计算结果 (8)2 共轭梯度法求解线性方程组 (10)2.1 算法原理 (10)2.2 程序流程图 (10)2.3 程序使用说明 (11)2.4 算例计算结果 (12)3 三次样条插值 (14)3.1 算法原理 (14)3.2 程序流程图 (16)3.3 程序使用说明 (17)3.4 算例计算结果 (19)4 龙贝格积分 (21)4.1 算法原理 (21)4.2 程序流程图 (22)4.3 程序使用说明 (23)4.4 算例计算结果 (24)结论 (26)1 QR 分解法求解线性方程组1.1 算法原理矩阵的QR 分解是指，可以将矩阵A 分解为一个正交矩阵Q 和一个上三角矩阵R 的乘积，实际中，QR 分解经常被用来解决线性最小二乘问题，分解情况如图1.1所示。

=⨯图1.1 QR 分解示意图本次上机学习主要进行了两个最基本的正交变换—吉文斯（Givens ）变换和豪斯霍尔德（Householder ）变换，并由此导出矩阵的QR 分解以及求解线性方程组的的方法。

1.1.1 基于吉文斯变换的QR 分解吉文斯矩阵也称初等旋转阵，如式（1.1）所示，它把n 阶单位矩阵I 的第,i j 行的对角元改为c ，将第i 行第j 列的元素改为s ，第j 行第i 列的元素改为s −后形成的矩阵。

计算方法大作业学院：电子工程姓名：班级：学号：大作业选题：分析方程求根问题中牛顿法的性能，包括收敛性等，并用该方法求解一个问题，给出过程和结果。

一、牛顿迭代法介绍：用迭代法求方程0)(=x f 的根时，首先要构造一个迭代函数，迭代函数构造的好坏，不仅影响收敛速度，而且有可能使迭代序列发散，构造迭代函数的一条重要途径，是用近似方程代替原方程去求根，因此如果能将非线性方程0)(=x f 用线性方程来近似代替，那么求近似根问题就容易得到解决，而且十分方便。

牛顿法就是把非线性方程线性化的一种方法。

二、牛顿迭代法原理设已知方程0)(=x f 的近似根0x ,则在0x 附近)(x f 可用一阶泰勒多项式))((')()(000x x x f x f x p -+=近似代替.因此, 方程0)(=x f 可近似地表示为0)(=x p .用1x 表示0)(=x p 的根,它与0)(=x f 的根差异不大.设0)('0≠x f ,由于1x 满足,0))((')(0100=-+x x x f x f 解得)(')(0001x f x f x x -= 重复这一过程,得到迭代公式)(')(1n n n n x f x f x x -=+ 这就是著名的牛顿迭代公式,它相应的不动点方程为)(')()(x f x f x x g -=.用牛顿迭代公式求方程根的方法称为牛顿迭代法，简称牛顿法。

三、牛顿迭代法的几何解析在0x 处作曲线的切线，切线方程为))((')(000x x x f x f y -+=。

令0=y ，可得切线与x 轴的交点坐标)(')(0001x f x f x x -=，这就是牛顿法的迭代公式。

因此，牛顿法又称“切线法”，其几何意义即为0x 点处的切线方程。

四、牛顿迭代法的收敛性 计算可得2)]('[)(")()('x f x f x f x g -=,设*x 是0)(=x f 的单根,有0)(',0)(**≠=x f x f ,则0)]('[)(")()('2****=-=x f x f x f x g , 故在*x 附近,有1)('<x g .根据不动点原理知牛顿迭代法对单根收敛.同理可知当*x 是0)(=x f 的重根时也收敛，则可分析出牛顿法不论对单根还是重根均是局部收敛的，只要初值足够靠近*x ，牛顿迭代序列均收敛于*x 。

《四位二进制数可控加减法》实验报告实验名称: 四位二进制数可控加减法姓名：学号：班级：目录一、实验方案 (3)二、设计思路................................................................................ 错误!未定义书签。

三、程序代码................................................................................ 错误!未定义书签。

四、调试问题 (6)五、心得感想 (7)一、实验方案1）基本功能实现两个四位二进制数的加减法运算,能够在led灯和数码管显示出结果。

2）清零功能利用一个微动开关，当微动开关按下时结果清零显示。

3）数码管显示将结果转换为七段显示器显示。

将运算结果输送到数码管中。

利用到人的视觉误差和短暂延时显示四位运算结果。

4）溢出问题若有溢出，则数码管显示“E”。

二、设计思路基本功能中分为连个模块，主模块用来运算加减法以及记录溢出和结果，子模块用来进行七段数码管的显示。

扩展功能中数码管显示要利用暂留现象，因此利用时钟clk来进行设计。

三、程序代码module show_sub(input [1:0]num,output reg [6:0] a_to_g );always @(*)case(num)2'b00: a_to_g=7'b1000000;2'b01: a_to_g=7'b1111001;2'b10: a_to_g=7'b1111111;2'b11: a_to_g=7'b0000110;default: a_to_g=7'b0000110;endcaseendmodulemodule show_top(input clk,clr,input wire [7:0] sw,input plus,sub,output wire [6:0] a_to_g,output reg [3:0] an,output reg [3:0] led );reg [15:0] clk_cnt;wire [1:0]s;reg [3:0] result; //运算结果reg [1:0] res;reg flag; //溢出标志wire [3:0] data1;wire [3:0] data2;assign data1=sw[7:4];assign data2=sw[3:0];assign s=clk_cnt[15:14];always @(posedge clk)beginclk_cnt=clk_cnt+1;endalways@(posedge plus or posedge sub or posedge clr)。

上海交通大学研究生（非数学专业）数学基础课程《计算方法》教学大纲(2007修改讨论稿)一．概况1.开课学院（系）和学科：理学院数学系计算数学教研室2.课程编码：3.课程名称：计算方法4.学时/学分：54学时/3学分5.预修课程：线性代数，高等数学，程序设计语言6.课程主干内容: 数值代数，数值逼近，非线性方程数值解，常微分方程数值解。

7.适应专业学科：全校的机、电、材、管理、生命和物理、力学诸大学科类，以及人文学科需要的专业。

8.教材/教学参考书：(1)李庆扬、王能超、易大义，数值分析（第4版），华中理工大学出版社， 2003(2)孙志忠，袁慰平，闻震初，数值分析，东南大学出版社，2002(3)J.Stoer and R. Bulirsch, Introduction to Numerical Analysis (secondedition), Springer-Verlag, Berlin-New York, 1993.(4)Atkinson K E，An Introduction to Numerical Analysis，John Wiley & Sons. 1989.二．课程的性质和任务本课程属于数值计算课程的基础部分。

数值计算课程是非数学类研究生数学公共基础课程，该组课程列入计算数学系列，目前按照“分级”的原则，设置《计算方法》（基础部分）、《微分方程数值方法》(扩展部分) 和《高等计算方法》（提高部分）三门课程。

本课程讨论用计算机求解数学问题的几类基本的数值方法及其相关的数学理论。

计算机是对近代科学研究、工程技术和人类社会生活影响最深远的高新技术之一，它对科学技术最深刻的改变，莫过于使科学计算平行于理论分析和实验研究，成为人类探索未知和进行大型工程设计的第三种方法和手段。

计算机的飞速发展正把计算的方法的创新、改进、提高推向人类科技活动的前沿。

人类现代计算能力的巨大更取决于计算方法的效率。

计算方法大作业班级XXXXXX 学号XXXXXXX 任课老师贺力平姓名XXX2013年12月实验一幂法与矩阵特征值1.幂法求主特征值思路幂法的主要思想就是对假设的任意初始列向量作用n次A矩阵（左乘A矩阵）后，初始向量就接近A矩阵的主特征值对应的特征向量。

由于左乘n次A矩阵有可能会造成计算量溢出，所以每次都对列向量作归一化处理。

2.程序代码function [ ] = mifa( A,v )%UNTITLED Summary of this function goes here % Detailed explanation goes hereA =[ 1 21 2 334 5 2 154 2 6 42 2 59 0];v=[1 1 1 1]';u(:,1)=v(:,1);fori=1:100v(:,i+1)=A*u(:,i);if abs(max(v(:,i+1))-max(v(:,i)))<10^-4breakendu(:,i+1)=v(:,i+1)/max(v(:,i+1));enddisp(u(:,i));disp(max(v(:,i)));k=u(:,i)'*A*u(:,i)/(u(:,i)'*u(:,i));disp(k);[x,c]=eig(A);disp(c);disp(x);disp(i);end3.结果比较和结论初始矩阵A =[ 1 21 2 334 5 2 154 2 6 42 2 59 0];初始向量v=[1 1 1 1]';幂法求得特征向量12.514714.914025.693539.6330归一化后特征向量0.31580.37630.64831.0000列向量最大值近似主特征值39.6330Rayleigh商求出主特征值39.6330用eig()函数算出的特征值和特征向量39.6331 0 0 00 -18.2401 + 7.4985i 0 00 0 -18.2401 - 7.4985i 00 0 0 8.8471-0.2450 -0.0471 + 0.0965i -0.0471 - 0.0965i -0.0574 -0.2919 0.1147 - 0.0939i 0.1147 + 0.0939i -0.1747 -0.5029 0.2862 - 0.1187i 0.2862 + 0.1187i 0.1535 -0.7758 -0.9330 -0.9330 0.9709 达到精度要求所需次数：17结论：可以看出初始列向量经过多次迭代后，用幂法求出的特征值和用eig()函数求出的A的特征值，满足计算精度在，并且特征向量也具有数乘关系。

1、多媒体计算机系统同一般计算机相比，特别需要利用计算机的数字化技术和______________ 。

选择一项：広a.通信技术D b.安全技术匕c.人机交互技术7匚d.存储技术2、计算机应用中，英文缩略语CAI所表示的计算机术语是_________ 。

* a.计算机辅助教学D b.计算机辅助设计匕c.计算机辅助制造匚d.计算机辅助工程3、当前气象预报已广泛采用数值预报方法，这主要涉及计算机应用中的____________ 。

匕a.数据处理和辅助设计□ b.科学计算与辅助设计-c.科学计算和数据处理d.科学计算和过程控制4、与二进制数110101110等值的十六进制数是屛b鼠标d.打印机/扫描仪6、假设给定一个十进制整数D,转换成对应的二进制整数B,那么就这两个数字的位数而言，B与D 相比， ________ 。

a. Db. Bc. DA d. B的位数大于等于的位数大于的位数大于的位数大于等于7、为解决某一特定的问题而设计的指令序列称为广a.b.语言文档c.系统卜d.程序丫8、对微机配置描述为中内存为 _________ 。

ra. 1.7G多媒体”，其b. 128MBI-c. 60GBa. 35ED b. 1AE fc. D70匚d. EA25、以下为计算机输出设备的是______________ 匸a.键盘d. 40X9、目前制造计算机所用的电子元件是广a.晶体管A b.集成电路电子管广c.超大规模集成电路10、有一类计算机应用领域称为 "CAD"，也就是b. 逢十六进 a. 计算机辅助教学c. b. 计算机辅助制造d. 逢八进15、八进制数631转成二进制数是c. 计算机辅助设计 计算机辅助管理d. a.10011101111、计算机内部采用二进制数进行运算、存储和控 制的主要原因是 广 b. 101011001a. 现/二进制只有0、1 两种状态，技术上容易实 c. 110011001 d.110100001二进制数数码少, b. 快 16、计算机的主机的组成部分是在计算机网络中传送速度 a. CPU 、外存储器、外部设备c. 二进制数数码少, 比十进制数容易读懂和记 广b. CPU 和存储器系统忆 二进制数数码少, d. c.主机箱、键盘、显示器存储起来占用的存储容量 小 12、在计算机中，所有信息的存储都采用 屛d. CPU 和内存储器17、下面哪一项不是计算机采用二进制的主要原因a. 卜六进制 二进制/b. a. 二进制只有0和1两个状态，技术上容易实c. 八进制 广bb. 二进制数的o 和1与逻辑代数的"真”和"假"十进制 d. 相吻合, 适合于计算机进行逻辑运算 13、计算机操作系统作为一个接口， 连接着广c. 二进制运算规则简单a. 用户与计算机7 金dd.二进制可与十进制直接进行算术运算b. 主机与外设 18、高级语言所编写的程序又称为源程序，此类程 序c. 系统软件与应用软件 a. 不大于100行的程序可以被机器直接执行 用户与软件d. A b.在更高级的大型计算机中能被机器直接执14、十进制数的运算法则是 a.逢二进 c. 不能被机器直接执行d.能被机器直接执行19、光盘是一种已广泛使用的外存储器，英文缩写CD-ROM旨的是_______ 。

一、大作业：本课程共包括3次大作业，旨在培养学生分析实际问题和解决实际问题能力。

要求学生自己实践与尝试，自己去调查、分析和计算，可以进行分组，进行学习小组交流、讨论，形成小组意见，课堂上安排小组代表作简要介绍，任课教师点评和总结。

1、设计“一个”新的金融产品。

2、计算一个具体的投资组合风险（例如VaR）以及解决风险的方法。

3、选择一个具体的金融产品定价（例如权证或者银行的理财产品）。

二、课后习题第１章金融工程概述1、请论述学习金融工程的三个基本目标，并举例说明。

2、根据已有的金融工程几个代表性定义，请阐述你对这几个定义的理解和看法。

3、请论述中国开展金融衍生产品交易的意义及其面临的问题。

第 2 章无套利定价原理1、假设市场的无风险借贷利率为 8 ％，另外存在两种风险证券 A 和 B ，其价格变化情况如图 2-11，不考虑交易成本。

图 2-11 两种风险证券的价格变化情况问题：（１）证券 B 的合理价格为多少呢？（２）如果 B 的市场价格为１１０元，是否存在套利机会？如果有，如何套利？（３）如果存在交易成本，例如，每次卖或买费用均为１元，结论又如何？2、假设无风险借贷半年利率 r ＝ 4 ％（单时期），两种资产的两时期价格变动情况如图2-12 ：图 2-12 两种资产的两时期价格变动情况问题：（１）利用动态组合复制定价技术给证券 B 定价；（２）如果证券 B 的市场价格为１００元，是否存在套利机会？如果有，如何构造套利策略？3 、试分析金融市场套利与商业贸易中的价差盈利的关系？为何金融市场中套利概念如此重要？第 3 章金融产品创新原理1 、如何设计一个更加合理的全流通方案？2 、如何设计一个金融新产品？第 4 章金融风险管理原理1 、金融风险是怎样产生的？如何从理论上解释金融风险？2 、怎样理解长期资本管理公司破产是一个由制度性缺陷、市场风险和流动性风险所造成的经典案例？3 、在例 4-1 中，当欧洲国家相关企业提出中国绍兴纺织企业向他们购买纺织设备，将终止使用美元支付的惯例，转为以欧元计价结算时，能否估计出１年之内因汇率波动产生的最大损失，若能是多少？4 、在例 4-1 中，能否找到一种套期保值方法，来减少思考题 3 估计出１年之内因汇率波动产生的最大损失。

上海交通大学工程硕士算法总结a.知识点:第一章: 复杂度(简单),欧几米德算法(中等),模运算(中等),算数运算第二章: 大师定理(简单,重要!)第三章: DFS算法(简单,重要!),寻找子部件(中等),图的反转第四章: BFS,Dijsktra算法(中等,重要!),bellman-ford(难)第五章: 最小生成树(中等)第六章: DAG(中等),背包问题(中等,重要!),独立集第七章: 线性规划的解(简单),最大流(中等),对偶(中等,重要!)b.作业(知识点):第一章: 0.1(复杂度),0.2(复杂度),1.2(扩展欧几里得),1.31(复杂度) [难],1.35(逆/费马小定理/模运算)[难]第二章: 2.5(大师定理),2.19(大师定理),2.22(大师定理)[难]第三章: 3.7(DFS),3.11(DFS),3.28(SAT)[难],第四章: 4.11(Dijkstra),4.12(Dijkstra),4.16(2分堆/d堆)第五章: 5.20(完美匹配),5.23(最小生成树),5.26(约束问题)第六章: 6.17(背包),6.21(最小覆盖/独立集),6.22(背包)第七章: 7.6(简单规划),7.7(线性规划的解),7.21(临界边)[难],7.23(归约问题)[难]红/黄/绿指重要度依次降低.注意:作业并不涵盖所有考点.c.考试注意点:1.在图计算中运用了标记,删边,反向等技巧,结合最短路径和DFS等进行计算.2.贪心算法出题比较无规则.3.大师定理,线性规划中的对偶,须理解背熟.几乎是必考.知识点详细:1.复杂度-->如无特殊要求复杂度一般只要表示为大O即可.*常数系数忽略*阶乘>指数>多项式>log*复杂度Θ的证明需要证明O和Ω2.加减的复杂度(n),乘除的复杂度(n2,或者n x m),(乘法可以更加简化)模运算(平方)*满足结合律,分配律,交换率*指数运算可记3.欧几里得算法*辗转相除法求最大公因数(练习)*扩展欧几里得算法求逆元4.逆元*ax=1(modN)成立,则x是关于a模N的一个乘法逆元*gcd(a,N)=1(互质)时,存在x和y,使得ax+Ny=1,x就是要找的逆元(因为y是可以约掉的),利用扩展欧几里得来计算此逆元.(可能性大)*练习一遍很重要5.大师定理*把一个问题分解为相同的子问题,通过子问题的规模,数量和合并子问题的代价可以计算整个问题的复杂度,牢记公式.*大师定理排列问题比较容易出到(归并排序).6.无向图*explore(重要)不断寻找相邻的结点,会变例可到达部分,不连通部分不访问.*DFS(重要)对图中每个点做前序和后序标记,然后运行explore,不连通部分也访问.前序和后序有些特殊性质.7.有向图*DFS(O(|V|+|E|))有向图一般运行DFS,运行后可以得到一颗树.有树边,前向边,回边,横跨边*强连通部件特殊的在这里寻找强连通部件运用了explore算法,整个算法复杂度同DFS.*输出图中所有强连通部件的步骤DFS后标记postpost最大点在源点强连通部件中(post排序)把图反转(每条边反转)explore找到一个强连通部件,也就是源点强连通部件删除后可以寻找下一个(post最大的)8.最短距离*BFS,用于距离的算法,对象无权图无权图的最短路径(O(|V|+|E|))*Dijkstra(采用数组是O(|V|2)),BFS扩展算法,适合边正的有权图迪杰斯特拉算法:预期不会有更长边,因而不适用于负边存在的情况有权图的最短路径*Bellman-ford,BFS扩展算法带负权值的最短路径,对比两点间全部可能的路径,可以使用与负边情况,但是注意,负环存在的情况下没有最短距离.*二分堆上述算法的数据结构9.最小生成树*分割性质两个点集间权重最轻的边总是可以加到最小生成树上(如果能连到树上) *Kruskal算法不断寻找最小边然后合并的过程,包括makeset(x)find(x)union(x,y)10.DAG最短路径*使用一维表计算DAG最短路径11.DAG最长递增子序列*一维表(考试要求的是表达式)12.背包问题*多副本背包一维表,每增加一点价值v为一格,每个格有n种可能性,复杂度nv.习题对应硬币问题.*单副本背包二维表,价值和物品建立一个二维表可得.习题中有对应.13.树的独立集*同样只是个表达式,只能用于树而不能用于图.14.网络流*最大流问题,利用剩余图来解决,即explore剩余图起点既可得到一个最大流的所有边.最小分割最大流定理15.二部图匹配*增加源点和汇点并求最大流(归约),注意二部图的源点汇点是人工增加的,说明加点是一种技巧.16.对偶*写出对偶表达式对偶的写法:题目(喜闻乐见):必考题目,第一梯队:1.大师定理相关(遇题请先上公式)2.有向图相关(涉及回边/横向边/树边),注意和这有关系的是DFS的标记.3.动态规划(20分) --> 估计是二维表,形式可能类似于编辑距离的计算,和最长公共子序列/串相关.(注意必须写出表达式,应该要比背包问题复杂)4.线性规划之对偶据说只有3题的剩余自由空间了,第二梯队:可能性比较大的有:第一章==> 复杂度/扩展欧几里得求逆元选一或者一起第四章==> Dijkstra/二分堆/bellman-ford(可以会一起出)第五章==> Kruskal+Disjoint set/割原理/集合覆盖(注意:老师很喜欢集合覆盖类问题,但是没有讲要考,比较危险)剩下提到要考的但是和已知章节冲突,第三梯队:背包/单纯形法/网络流/二部图匹配等关于证明题:if and only if形式要证明if和only if,只要给出形式就给分.也就是说,对于证明题, =/<=>/Θ这几种形式,必须正反证明.如果没有时间,先证明简单的一种.可能会出现复杂度运算的题目.。

计算方法A 上机大作业1. 共轭梯度法求解线性方程组算法原理：由定理3.4.1可知系数矩阵A 是对称正定矩阵的线性方程组Ax=b 的解与求解二次函数1()2TT f x x Ax b x =- 极小点具有等价性，所以可以利用共轭梯度法求解1()2TT f x x Ax b x =-的极小点来达到求解Ax=b 的目的。

共轭梯度法在形式上具有迭代法的特征，在给定初始值情况下，根据迭代公式：(1)()()k k k k x x d α+=+产生的迭代序列(1)(2)(3)x x x ，，，... 在无舍入误差假定下，最多经过n 次迭代，就可求得()f x 的最小值，也就是方程Ax=b 的解。

首先导出最佳步长k α的计算式。

假设迭代点()k x 和搜索方向()k d 已经给定，便可以通过()()()()k k f x d φαα=+的极小化()()min ()()k k f x d φαα=+来求得，根据多元复合函数的求导法则得:()()()'()()k k T k f x d d φαα=∇+令'()0φα=，得到:()()()()k T k k k T k r d d Adα= ，其中()()k k r b Ax =-然后确定搜索方向()k d 。

给定初始向量(0)x 后，由于负梯度方向是函数下降最快的方向，故第一次迭代取搜索方向(0)(0)(0)(0)()dr f x b Ax ==-∇=- 。

令(1)(0)00x x d α=+其中(0)(0)0(0)(0)T T r d d Adα=。

第二次迭代时，从(1)x 出发的搜索方向不再取(1)r ，而是选取(1)(1)(0)0dr d β=+，使得(1)d 与(0)d 是关于矩阵A 的共轭向量，由此可求得参数0β:(1)(0)0(0)(0)T T r Ad d Adβ=-然后从(1)x 出发，沿(1)d 进行搜索得到(2)(1)(1)1x x d α=+设已经求出(1)()()k k k k x x d α+=+，计算(1)(1)k k r b Ax ++=-。

吉林大学网络教育学院2019-2020学年第一学期期末考试《计算方法》大作业答案学生姓名专业层次年级学号学习中心成绩年月日作业完成要求：大作业要求学生手写，提供手写文档的清晰扫描图片，并将图片添加到word文档内，最终wod文档上传平台，不允许学生提交其他格式文件（如JPG，RAR等非word文档格式），如有雷同、抄袭成绩按不及格处理。

一、解线性方程（每小题8分，共80分）1、用矩阵的LU分解算法求解线性方程组X1+2X2+3X3= 02X1+2X2+8X3= -4-3X1-10X2-2X3= -11答：2、用矩阵的Doolittle分解算法求解线性方程组X1+2X2+3X3= 12X1– X2+9X3= 0-3X1+ 4X2+9X3= 1答：3、用矩阵的Doolittle分解算法求解线性方程组2X1+X2+X3= 46X1+4X2+5X3=154X1+3X2+6X3= 13答：4、用高斯消去法求解线性方程组2X1- X2+3X3= 24X1+2X2+5X3= 4-3X1+4X2-3X3= -3答：5、用无回代过程消元法求解线性方程组2X1- X2+3X3= 24X1+2X2+5X3= 4-3X1+4X2-3X3= -3答：6、用主元素消元法求解线性方程组2X1- X2+3X3= 24X1+2X2+5X3= 4-3X1+4X2-3X3= -3答：7、用高斯消去法求解线性方程组1231231232344272266x x x x x x x x x -+=++=-++=答：8、利用Doolittle 分解法解方程组Ax=b ，即解方程组12341231521917334319174262113x x x x -⎡⎤⎡⎤⎡⎤⎢⎥⎢⎥⎢⎥-⎢⎥⎢⎥⎢⎥=⎢⎥⎢⎥⎢⎥--⎢⎥⎢⎥⎢⎥--⎣⎦⎣⎦⎣⎦ 答：9、利用Doolittle 分解法解方程组Ax=b ，即解方程组123421111443306776081011112x x x x ⎡⎤⎡⎤⎡⎤⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥=⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥-⎣⎦⎣⎦⎣⎦ 答：10、用高斯消元法解方程组1237811351341231x x x -⎡⎤⎡⎤⎡⎤⎢⎥⎢⎥⎢⎥-=-⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦⎣⎦⎣⎦答案：二、计算（每小题10分，共20分）1、已知节点x1,x2及节点处函数值f(x1)，f(x2)，构造线性插值多项式p1(x). 答：2、设f(xi)=i(i=0,1,2),构造二次式p2(x)，使满足： p2(xi)=f(xi)(i=0,1,2)答：。

计算方法大作业一、计算方法学后感计算方法学习主要内容为函数逼近论，数值微分，数值积分，误差分析等。

常用方法有迭代法、差分法、插值法、有限元素法等，印象深刻的是三次样条差值函数，大概是做了作业的缘故。

计算也称数值分析，数学不好的人也没有学会太多东西，不过还是有些小体会的：1>主要矛盾和次要矛盾之间的关系。

现实问题中有很多约束条件，需要我们有侧重的保留摒弃，辨析主要矛盾和次要矛盾，从而提出合理假设；2>尺有所短寸有所长。

没有完美无缺的算法，虽然我们看到有不断地改进优化算法，但这些往往都是以牺牲某些优点为代价的。

比如提高精度，往往会导致格式复杂，产生较大运算量；3>原则不能变。

算法也是要讲原则的，比如要谈算法的优劣性前提是要保证算法的可靠性（相容、收敛、稳定等）。

门外汉的感受也就这水准了。

至于计算方法思想在实际生活中的应用，我想来想去就这些了。

（1）天气预报：天气会受各种因素的影响，稍微一些因素发生改变就会产生很大的变化，所以天气预报其实是一件比较困难的工作。

古代人们用占卜或者经验总结等方式来预计天气状况，这倒更像是统计学。

而有了计算机，我们就可以通过数值模拟来预报天气。

具体过程就是：首先根据大气运动列出数学物理方程（偏微分），其次对空间进行网格划分，然后通过观测数据给出初值条件，最后通过数值方法求解这些偏微分方程得到网格点处的数值解。

这也是为什么主持人总是说大概在...地区，大致在...时段，可能有...量级的降水...因为时空是连续的，而网格划分不可能无限密，所得的数值解也存在误差。

（2）等离子体：对等离子体现有的理论描述中，磁流体力学、符拉索夫方程、福克-普朗克方程等都是微分方程，包含很多参量，如果要求出解析解，物理模型往往需要过分简化以至于无法精确和全面的包罗各种效应，所以需要数值计算，这也是等离子体物理学研究中很重要的一个方面。

比如最简单的单粒子模型，它的牛顿洛伦兹方程是这样：。

数值计算第一次大作业实验目的 以Hilbert 矩阵为例，研究处理病态问题可能遇到的困难。

内容 Hilbert 矩阵的定义是,()11/21/31/1/21/31/41/(1)1/31/41/51/(2)1/1/(1)1/(2)1/(21)n i j H h nn n n n n n =⎡⎤⎢⎥+⎢⎥⎢⎥=+⎢⎥⎢⎥⎢⎥++-⎣⎦它是一个对称正定矩阵，而且()n cond H 随着n 的增加迅速增加，其逆矩阵1,()n i j H α-=，这里,2(1)(1)!(1)!(1)[(1)!(1)!]()!()!i j i jn i n j i j i j n i n j α+-+-+-=+----- 1) 画出ln(())~n cond H n 之间的曲线（可以用任何的一种范数）。

你能猜出ln(())~n cond H n 之间有何种关系吗？提出你的猜想并想法验证。

用行范数for n=1:50 for i=1:n for j=1:nA(i,j)=1/(i+j-1);B(i,j)=factorial(n+i-1)*factorial(n+j-1)/((i+j-1)*(factorial(i-1)*factorial(j-1))^2*factorial(n-i )*factorial(n-j));end endresult1=0; for j=1:nresult1=result1+A(1,j); endresult1=log(result1); result2=0; for i=1:n for j=1:nresult2=B(i,j)+result2; endresult(i)=log(result2); endm=max(result);x(n)=result1+m; end plot([1:50],x)对于更大的n 值，由于Hilbert 逆矩阵中的元素过大，溢出，故在此取50以内的n 。

图1 ln(())~n cond H n 关系曲线图猜想ln(())~n cond H n 之间存在线性关系 验证：设ln(()n cond H an b ∞=+ 在以上程序基础上，再添加>>;>> y=x'; >> l=1:40; >> k=l';>> p=polyfit(k,y,1) %一次多项式拟合 p =3.5446 -3.0931% P=polyfit(k,y,2) %二次多项式拟合 p =-0.0008 3.5778 -3.3253 % P=polyfit(k,y,3) %三次多项式拟合0.0000 -0.0033 3.6198 -3.4777% P=polyfit(k,y,4) %四次多项式拟合-0.0000 0.0002 -0.0082 3.6654 -3.5815 % P=polyfit(k,y,5) %五次多项式拟合 p =0.0000 -0.0000 0.0007 -0.0156 3.7107 -3.6542 从上式可以看出，高次项系数相对于一次项和常数项系数要小很多， 所以取ln(() 3.5446 3.0931n cond H n ∞=-2）设D 是n H 的对角线元素开方构成的矩阵。

交通大学成人本科计算机图形学期末大作业As a person, we must have independent thoughts and personality.《计算机图形学》期末大作业学号：姓名：李燕军学习中心：校本部注：将本作业的word文件、最后一题作品的.fla文件和.swf文件一起压缩成一个文件提交一、术语解释（15题×2分= 30分，1-10题英文缩略词要求写出的中文和英文全称，以课程教材范围内为准；11-15题写出概念解释）1、UI：用户界面2、IBR：基于图像的绘制3、VR：Virtual Reality 虚拟现实4、LOD：5、GKS：6、PHIGS：程序员层次交互式7、RSD：光栅扫描显示器8、CAM：Computer Aided Manufacture 计算机辅助制造9、OpenGL：是独立于视窗操作系统或其它操作系统的，亦是网络透明的。

帮助程序员实现在 PC、工作站、超级计算机等硬件设备上的高性能、极具冲击力的高视觉表现力图形处理软件的开发。

10、UCS：user coordinate system用户坐标系描述物体几何模型的坐标系。

有时也称为局域坐标系（local coordinate system LCS）。

用户坐标系也是实数域坐标系。

11、灭点：与平行投影相比透视投影的特点是所有投影线都从空间一点（称为视点或投影中心）投射，离视点近的物体投影大，离视点远的物体投影小，小到极点消失，称为灭点。

12、裁剪：在二维观察中，需要在观察坐标系下根据窗口大小对世界坐标系中的二维图形进行裁剪（clipping），只将位于窗口内的图形变换到视区输出。

13、投影：答：投影就是从投影中心发出射线经过三维物体上的每一点后与投影面相交所形成的交点集合。

14、消隐：真实感图形绘制过程中，由于投影变换失去了深度信息，往往导致图形的二义性。

要消除这类二义性，就必须在绘制时消除被遮挡的不可见的线或面，习惯上称之为消除隐藏线和隐藏面，或简称为消隐，经过消隐得到的投影图称为物体的真实图形15、走样：二、简答题（2题×5分=10分）1、简要说明计算机图形学、图像处理、模式识别三者之间的区别和联系。

P50 – 1%%牛顿插值多项式function [ c, d] = newpoly( x,y )%这里x为n个节点的横坐标所组成的向量，y为纵坐标所组成的向量。

%c为所求的牛顿插值多项式的系数构成的向量。

n=length(x);d=zeros(n, n);d(: , 1)=y';for j=2 : nfor k=j : nd(k, j)=(d(k, j-1) - d(k-1, j-1)) / (x(k)-x(k-j+1));endendc =d(n, n);for k=(n-1) : - 1 : 1c =conv(c, poly(x(k)));m=length(c);c(m)=c(m)+d(k, k);end>> x = 0.2:0.2:1 ;>> y =[ 0.98,0.92,0.81,0.64,0.38] ;>> c= newpoly(x, y )c = -0.5208 0.8333 -1.1042 0.1917 0.9800%%三次样条插值x=[0.2,0.4,0.6,0.8,1.0];y=[0.98,0.92,0.81,0.64,0.38];x0 = [0.2,0.28,1.0,1.08];pp=csape(x,y,'variational');%% 三次样条函数表达式disp(pp.coefs);-1.3393 -0.0000 -0.2464 0.98000.4464 -0.8036 -0.4071 0.9200-1.6964 -0.5357 -0.6750 0.81002.5893 -1.5536 -1.0929 0.6400绘制曲线图x2 = 0:0.01:1.2 ;y11 = polyval(c,x2) ;y22 = ppval(pp,x2);x0 = [0.2,0.28,1.0,1.08];y110 = polyval(c,x0);y220 = ppval(pp,x0);plot(x2,y11,'r',x0,y110,'^',x2,y22,'g',x0,y220,'h')legend('牛顿插值','牛顿插值样点','三次样条插值','三次样条插值样点')P50 -3(1)x = [0,1,4,9,16,25,36,49,64] ;y = 0:8 ;x1 = 0:0.1:64 ;x2 = 0:0.01:1 ;f = lagrange(x,y)%% 得到多项式函数表达式L(x)= - 3.28063e-10*x^8 + 6.71268e-8*x^7 - 0.00000542921*x^6 + 0.000222972*x^5 - 0.00498071*x^4 + 0.0604294*x^3 - 0.38141*x^2 + 1.32574*xy1 = lagrange(x,y,x1) ;y2 = lagrange(x,y,x2) ;(2)(2) x = [0,1,4,9,16,25,36,49,64] ;y = 0:8 ;x1 = 0:0.1:64 ;x2 = 0:0.01:1 ;%% 得到三次样条差值函数表达式pp=csape(x,y,'not-a-knot');disp(pp.coefs);0.0266 -0.2998 1.2732 00.0266 -0.2199 0.7534 1.0000-0.0021 0.0197 0.1529 2.00000.0005 -0.0112 0.1955 3.0000-0.0000 -0.0001 0.1160 4.00000.0000 -0.0014 0.1026 5.00000.0000 -0.0005 0.0825 6.00000.0000 -0.0004 0.0717 7.0000y11 = ppval(pp,x1) ;y22 = ppval(pp,x2) ;绘制图形（1）在[0,64]显然随着次数越高，多项式插值出现误差很大（2）[0,1]在[0,1]区间上三次样条插值和多项式插值基本一致P137-1insucomplex_4_1.m 文件clear ;clc ;%h为步长,可分别令h=1,0.1,0.01,0.001h = [1,0.1,0.01,0.001]for i = 1:4h(i) ,x=0:h(i):1;y=sqrt(x).*log(x+eps);%复化梯形公式T=trapz(x,y);T=vpa(T,7),f=inline('sqrt(x).*log(x)',x);%复化辛普生公式S=quadl(f,0,1);S=vpa(S,7),end>> t = -log(h) ;>> plot(t,T,'rs',t,S,'r*')>> lengend('复合梯形公式','复合梯形公式')。

大作业3 对于常微方程数值解问题 检查各种数值算法的长期行为 观察步长对于收敛效果的影响给定方程组()()()()(0)0(0)x t ay t y t bx t x y b'=⎧⎪'=-⎪⎨=⎪⎪=⎩ 1．证明方程组的解是xOy 平面上的一个椭圆；2．利用①改进的欧拉折线法，②4阶标准龙格－库塔法，选几个不同的步长h ，计算上述方程组的轨道，看看哪种方法和步长能够保持椭圆轨道不变。

（计算的时间步要足够多――至少10000步） 1. 证明：因为'()(),x t ay t ='()(),y t bx t =- 所以''()'()(),x t ay t abx t ==- 所给方程的特征方程为2,r ab =-所以12,,r r ==是一对共轭复根，所以所求通解为12())),x t C t C t =+ 因为(0)0,x =所以2()),x t C t =同理可得：4()))y t b t C t =+ 因为'()(),x t ay t =即24'()))),x t C t ab t aC t ==+所以240C C ==，所以()),x t t =()) y t b t=所以2221 x yab b+=方程组的解是xOy平面上的一个椭圆。

2.设a=4,b=1;（1）用改进的欧拉折线法计算，程序及结果如下：>>a=4;b=1;n=10000;>> H=6*pi;>> [x,y]=eulermend(a,b,H,n);>>plot(x,y)图1 h=6*pi/10000时的椭圆轨道图像>> H=7*pi;>> [x,y]=eulermend(a,b,H,n);>> plot(x,y,'g--')图2 h=7*pi/10000时的椭圆轨道图像>> H=10*pi;>> [x,y]=eulermend(a,b,H,n);>> plot(x,y,'r--')图3 h=10*pi/10000时的椭圆轨道图像>> H=20*pi;>> [x,y]=eulermend(a,b,H,n);>> plot(x,y)图4 h=20*pi/10000时的椭圆轨道图像>> H=30*pi;>> [x,y]=eulermend(a,b,H,n);>> figure(2),plot(x,y,'r-')图5 h=30*pi/10000时的椭圆轨道图像从上述结果可知，改进的欧拉折线法对于不同步长能够保持椭圆轨道不变。