当前位置：文档之家 › 数值分析第二次作业

数值分析第二次作业

  • 格式：doc
  • 大小：301.50 KB
  • 文档页数：36

北航2007级研究生
《数值分析》计算实习作业
（第二题）
院系：航空科学与工程学院
学号： SY0705118
姓名：刘德博
日期： 2007年12月09日
计算实习说明书
目的：训练运用计算机进行科学与工程计算的能力。

要求：
1．独立进行算法设计、程序设计和上机运算，并得出正确的结果。

2．编制程序时全部采用双精度，要求按题目的要求设计输出，并执行打印。

3．只能根据题目给出的信息并且只允许一次计算得出全部结果。

题目：
第二题 关于x ，y ，t ，u ，v ，w 的下列方程组
0.5cos 2.67
0.5sin 1.070.5cos 3.740.5sin 0.79
t u v w x t u v w y t u v w x t u v w y +++-=⎧⎪+++-=⎪
⎨
+++-=⎪⎪+++-=⎩ 以及关于z ，t ，u 的下列二维数表
确定了一个二元函数(,)z f x y =。

1.试用数值方法求出(,)f x y 在区域{(,)|00.8,0.5 1.5}D x y x y =≤≤≤≤上的一个近似表达式：
,0
(,)k
r s rs
r s p x y c
x y ==
∑
要求(,)p x y 一最小的k 值达到以下精度：
1020
2700
((,)(,))10i j i j i j f x y p x y σ-===-≤∑∑
其中0.08,0.50.05i j x i y j ==+。

2.计算****
(,),(,)(1
,2,,8;1,2,,5)i j i j f x y p x y i j == 的值，以观察(,)p x y 逼近(,)f x y 的效果，其中**0.1,0.50.2i j x i y j ==+。

说明：
1．用迭代方法求解非线性方程组时，要求近似解向量()k x 满足：
()(1)()12||||/||||10k k k x x x --∞∞-≤
2．作二元插值时，要使用分片二次代数插值。

3．要求由程序自动确定最小的k 值。

4．打印以下内容： （1）算法的设计方案。

（2）全部源程序（要求注明主程序和每个子程序的功能）。

（3）数表：,,(,)(0,1,,10;0,1,,20)i j i j x y f x y i j == 。

（4）选择过程的,k σ值。

（5）达到精度要求时的,k σ值以及(,)p x y 中的系数(0,1,,;rs c r k s ==
0,1,,)k
（6）数表：******
,,(,),(,)(1
,2,,8;1,2,,5)i j i j i j x y f x y p x y i j == 。

5．采用f 型输出**,,,i j i j x y x y 的准确值，其余实型采用e 型输出并至少显示12位有效数字。

一、程序算法的设计
根据题目要求，并结合题目的特点，算法设计方案如下表所示：
二、全部源程序
本人编程所用的软件是Fortran，所编写的程序如下：
三、计算结果
以下为“计算结果.txt”文件中的内容，即为计算结果：
四、心得体会
本次作业从难度上来说要比第一次的大，首先要求解非线性方程，其次就是分片二次插值，再次就是曲面的拟合，而且拟合要自动确定最小的k值。

本题主要是分析好做题的思路，由于本人提前已经开始编写解非线性方程的程序，由于没有分析好题意，以为将t，u代入方程去求解v，w，x，y，所以采用简单迭代法编程，虽然也能解出值来，但是却不能往下做了。

后来才弄明白是要将x，y视为已知，然后确定t，u，v，w。

在用插值法求z。

看来在做一件事的时候，首先要弄明白怎么去做再去下手。

在本次编程过程中遇到很多问题。

比如说精度问题，对结果有明显的影响，这是在程序的后期阶段，由于程序总是不能收敛时才发现的。

再如，在求拟合系数时是解方程还是求矩阵的逆，本人结合矩阵理论的知识，同时翻阅了关于矩阵求拟的资料，觉得对于方阵通过行变换来求逆的根本原理与高斯消元相同，而且比较简练，所以最终采用这种方法。

再如，Fortran里的数组下标默认均是从1开始，而不是从0开始，所以在求拟合多项式的时候，就用到了加减1的情况，使程序可读性有所降低。

在C语言中，有这样的规定就是尽量不用或少用goto语句，在C里也很容易将goto语句用恰当的循环语句来替换，但是在Fortran中，有些不容易替换，在程序编写好以后，有两处用到goto语句，本人都用do while语句替换了，但是在循环前也都添加了程序执行一次的充分条件。

这也算是程序方面新的认识吧。

这次编程觉得收获很大，以后虽然没有像这样的作业了，但是本人觉得数值分析很有意义，也将会一直坚持下去，学会更多的数值分析方法及其程序的实现。

合集下载

数值分析第二次作业及答案

数值分析第二次作业及答案

数值分析第二次作业及答案1. 用矩阵的直接三角分解法（LU 分解）解方程组1234102050101312431701037x x x x ⎡⎤⎡⎤⎡⎤⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥=⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦⎣⎦⎣⎦。

212223243132333441424344212223243132333441424344110201020101011124310103110201101,112110101l uu u l l u u l l l u luu u l l u u l l l u ⎡⎤⎡⎤⎡⎤⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥=⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦⎣⎦⎣⎦⎡⎤⎡⎤⎡⎤⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⇒=⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎢⎥⎣⎦⎣⎦⎣⎦解：设111122223333444410201012125115102053101310136212117216420101724y y x x y y x x y y x x y y x x ⎡⎤⎢⎥⎢⎥=⎢⎥⎢⎥⎥⎣⎦==⎡⎤⎧⎡⎤⎧⎡⎤⎡⎤⎡⎤⎡⎤⎢⎥⎪⎢⎥⎪⎢⎥⎢⎥⎢⎥⎢⎥==⎪⎪⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥=⇒=⇒⎨⎨⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥==⎪⎪⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎪⎪==⎢⎥⎢⎥⎣⎦⎣⎦⎣⎦⎣⎦⎣⎦⎩⎣⎦⎩由2.矩阵第一行乘以一数，成为211A λλ⎡⎤=⎢⎥⎣⎦,证明当23λ=±时，()cond A ∞有最小值。

112131213 21121223A AA A λλλλλλλλ--∞∞⎧⎛⎫≥- ⎪⎪⎛⎫⎪=⇒==⇒=+ ⎪⎨⎪ ⎪⎝⎭⎪<- ⎪⎪⎩⎝⎭解：123632() min ()7.22343cond A AA cond A λλλλλ-∞∞∞∞⎧+≥⎪⎪∴====⎨⎪+<⎪⎩故当时， 3. 设有方程组AX b =，其中121011221,,302223A b ⎡⎤⎢⎥-⎡⎤⎢⎥⎢⎥⎢⎥==⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦⎢⎥-⎢⎥⎣⎦已知它有解12130X ⎡⎤⎢⎥⎢⎥⎢⎥=-⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦。

数值分析第二次大作业——编写ln函数实验报告

数值分析第二次大作业——编写ln函数实验报告


为一个整体，将 ln(a)进行 Taylor 级数分解，进而按式计算可以求得 ln(x)的值。由于我们选 取级数的前 n 项和近似 ln(a)，则有：
Rn (a
1)

(1)n
(a 1)n1 n 1
…；
|
Rn
(a
1)
||
(a 1)n1 n 1
|
。
2
自 92 乔晖 2009011414
精度为 20 位，即 h4 1020 。
int * GetNumber() {return number;}
//获得大数中的整数数组
bool GetSgn() const {return sgn;}
//获得大数的符号
bool IsZero();
//判断大数是否为 0
CBigNumber Reverse() {sgn = !sgn;return *this;} //取相反数
数值分析——编写 ln 函数
实验报告
自 92 乔晖 2009011414
自 92 乔晖 2009011414
数值分析——编写 ln 函数实验报告
目录
一、 实验任务 ................................................... 2
二、 编译环境 ................................................... 2
1u
x 1
35
5
自 92 乔晖 2009011414
数值分析——编写 ln 函数实验报告
其前 n 项和，近似可得 ln(x)的值。
基于以上原理，可以设计迭代算法，累加 2u2i1 得到前 n 项和。每次计算后，判断计 2i 1
数值分析第二次大作业

数值分析第二次大作业

《数值分析》计算实习报告第二题院系：机械工程及自动化学院学号：姓名：2017年11月一、题目要求试求矩阵A =[a ij ]10×10的全部特征值，并对其中的每一个实特征值求相应的特征向量，已知a ij ={sin (0.5i +0.2j ) i ≠j1.52cos (i +1.2j ) i =j(i,j =1,2, (10)说明：1.用带双步位移的QR 方法求矩阵特征值，要求迭代的精度水平为ε=10−12。

2.打印以下内容： （1）全部源程序；（2）矩阵A 经过拟上三角化后所得的矩阵A (n−1)； （3）对矩阵A (n−1)实行QR 方法迭代结束后所得的矩阵；（4）矩阵A 的全部特征值λi =(R i ,I i ) (i =1,2,⋯,10)，其中R i =Re(λi )，I i = Im(λi ) 。

若λi 是实数，则令I i =0； （5）A 的相应于实特征值的特征向量。

3.采用e 型数输出实型数，并且至少显示12位有效数字。

二、算法设计思路和方案1. 将矩阵A 拟上三角化得到矩阵A (n−1)为了减少计算量，一般先利用Householder 矩阵对矩阵A 作相似变换，把A 化为拟上三角矩阵A (n−1)，然后用QR 方法计算A (n−1)的全部特征值，而A (n−1)的特征值就是A 的特征值。

具体算法如下：记(1)A A =，()r A 的第r 列至第n 列的元素为(r)(1,2,,;,1,,)ij a i n j r r n ==+。

对于1,2,,2r n =-执行（1）若()(2,3,,)r ira i r r n =++全为零，则令(1)()r r A A +=，转(5)；否则转（2）。

（2）计算r d =()()1,sgn r r r r r c a d +=- (若()1,0r r r a +=，则取r r c d =)2()1,r r r r r r h c c a +=-（3）令()()()()1,2,0,,0,,,,T r r r n r r rr r rnru ac aaR ++=-∈。

数值分析大作业

数值分析大作业

第二次计算实验：SVD及其应用梁杰存2014310739航博1431.方法求矩阵A奇异值分解一个途径是求解A T A的特征值，但因为舍入误差容易丢掉小奇异值。

因此通常先将矩阵上双对角化，即构造正交阵Q和W,使得Q T AW=B(upper−bidiagnal)。

这一过程可以通过逐次Household变换或逐次Given’s变换完成，还有一种基于待定系数法思想的Lanczos算法。

由于Linpack中SVD算法需要输入上双对角矩阵，本文采用Lanczos 算法实现上双对角化。

1.1.隐式零移位QR法（implicit zero-shift QR）与传统移位QR迭代算法不同，隐式零移位QR算法不进行移位，并且第一步构造右乘Given’s变换矩阵GR(1,2)将上双对角矩阵B的(1,2)位置上的元素12b消零，而不是传统方法中引入一个非零元素21b。

但这一步可能会使原来为零的b12变为非零。

第二步左乘Given’s阵GL(1,2)使得12b为0，但可能会使为零b13变为非零。

与上述步骤类似，将b13变为0后可能会使b23非零。

如下图所示，重复上述步骤最终将恢复为上双对角矩阵，即完成一步隐式零移位QR迭代。

反复迭代，矩阵B将趋近与对角阵阵，对角元即特征值。

图1隐式QR迭代1.2.分而治之(Divide-and-conquer)分而治之算法将上双对角阵B分成有两个互相独立对角块矩阵与另一矩阵之和，即：B=B100B2+b m vvT=Q1Σ1Q1T00Q2Σ1Q2T+b m vv T =Q100Q2(Σ100Σ1+b m uuT)Q1T00Q2T所以矩阵B的特征值与矩阵D+ρu u T的特征值相同，其中D=Σ100Σ1为对角阵，又：det D+ρu u T−λI=det⁡((D−λI)(I+ρD−λ−1u u T))由于D−λI非奇异，则det I+ρD−λ−1u u T=1+ρu T D−λ−1u=1+ρu i2d i−λ=0ni=1在每个d i与d i+1之间分布着一个特征值，可用牛顿法快速找到该特征值。