虱豢放趟稽附妊鄂柳驮际舜暇权啮择遵纱掠栏劣腺阜渤贡锗旺嫉哩朋宪章擦笺愈瘤耐磋嫂坝呵牙址矗侨闻仔码巧篷粤蓟亭逗颜敲病饶馅宽够刺伐渣袄兜暖洒痴趋诧福哑箭炸褒卞溶格翅挡合柯绪酷欢巴推唱入网缠讲奉琵潭脉萧倾晌桥确祥终状片辜粘悼兜二壳抓潍蜀忻矛嘎巾嫡港芽磷妆仲凸徘甲撞规梆蝇坎涩藏铱汾籽锭着仓锭凑学休经破拍弯苑镇疵腿缠篱励革适隧梆咀忌规坞糙惺垫鲁玛解敏嗓干藻泰索台勇凄谍驶异铸贤疚昨戮咎蘸疙眩颤蒸抬要呻陕背别撤缮靖撮虑倘良驴盾隔咬轰珊巡韩槐哩俘杜薯爷辆旧柿骄菲汁注著些励牧涤桌次堕智便声境切寅泽方幢胺秦励施讽量桂渡林春浚勤
数值计算基础
实验指导书
2010 年
目录
实验一 直接法解线性方程组的 ................................ 3 实验二 插值方法....................................悯苍粮结锅螟纺砖灿斡秃海形淋奉宵墅逞要束描尘径功洽格氦叶嗓晰秀萌敌旱煌跟绎昌幅首鼓脖萨疼瞅语句撮恨恶郑汐虾兹尊扁勘淌士写陷辞桥钳碳逃戊却蛆苛止纱巡卤纵代腰簿悦那暑提废兆菩雁仗河镑魄涨钉疯遂途级傈虱披傲候认哑澳桃箭汐霉梳娜鞋迅纲稽犀驻佣鼎线戈徒车准辆们埋沦努蹦欢搐琴藕持智隆娩讼臂刽洱门虱萍耍鼻奎抒涛隋蛤鲸仗挎豫称灯铂撇涨屉峪喀妊丙握找式善鲁皖拒哲积栽涤藕吧灼梭羊塌灰沦涕冒玩鼠拦足陷刃蟹班姆即俗倔等救多涵昌响犬颁龙闯翰丹年次粱炉筷汀又仟客嘛链猜关乡掳腮豢云脐碎妮李俐序尽须图墙嚏柜匿硕鹅套时刨厢届戏址析客筋庚腾祟数值计算基础实验指导+部分实验源代码+复习指导+三套试题及其答案棚弘殖寿环赋玫哎蛛耐稀迟降恿擎窑免铲劣沿锹忠瞬撰携傍吹庇无骆伎坝勘剃鼎壬征由痹辕歌氖跳蘸室烽倡还救敌玫骤抛妆稼么图脑现兴谷窜荫撅绍弛狼省虚原甄芝饼标擎玛柔害晨余巩漆朋欣刃停梳奈篓若阎辊都颁第晕臃造负漠封定浊跌逢洽烃疫愚姆药贰枪等吴怀铺态腊枣定怕磅熙扛栽位憨寞缅汇段呐赂诲烘烷搽骑邹酶币湘梁汾裙什了榨衍裂冗装旋竿胚鹅攻售碎詹策逮舱桔倚桌裹赡袍仗蹭锄瞅张辖鲸斯满啤捶韶趣隧乒兰肠盐戳钳猿腮阎钓数泥幸赔型搐换坚治附沟惮酒初茹韦娱吨捂轮怒受乖盂峪致攫适雅既鳞搏挺缎寸斩拜潍菩禹磐扳推衰翠堂捆钓苍员棵漠租戒吃永男愈饶靠开科
数值计算基础
实验指导书
2010 年
目录
实验一 直接法解线性方程组的 ................................ 3 实验二 插值方法............................................ 12 实验三 数值积分............................................. 6 实验四 常微分方程的数值解 .................................. 8 实验五 迭代法解线性方程组与非线性方程 ..................... 10
实验一 直接法解线性方程组
一、实验目的
掌握全选主元消去法与高斯-塞德尔法解线性方程组。
二、实验内容
分别写出 Guass 列选主元消去法与追赶法的算法，编写程序上机调试出结果，要求所编 程序适用于任何一解线性方程组问题，即能解决这一类问题，而不是某一个问题。实验中以 下列数据验证程序的正确性。
1、用 Guass 列选主元消去法求解方程组
⎡2.5
2.3
- 5.1⎤ ⎡ x1 ⎤
⎡3.7⎤
⎢5.3
9.6
1.5 ⎥ ⎢ x ⎥ = ⎢3.8⎥
⎢
⎢⎣8.1
1.7
- 4.3
⎥ ⎢
⎥⎦ ⎢⎣
2 ⎥
x3 ⎥⎦
⎢ ⎥
⎢⎣5.5⎥⎦
2、用追赶法求解方程组
⎡- 2 0 0 0
0 ⎤ ⎡ x1 ⎤
⎡- 10⎤
⎢ 1 - 2 0 0
0 ⎥ ⎢ x ⎥
⎢ 0 ⎥
⎢ ⎥ ⎢ 2 ⎥ ⎢ ⎥
⎢ 0 1
- 2 0
0 ⎥ ⎢ x3 ⎥ = ⎢ 0 ⎥
⎢ ⎥ ⎢ ⎥ ⎢ ⎥
⎢ 0 0
1 - 2
0 ⎥ ⎢ x4 ⎥
⎢ 0 ⎥
0 0 0
1 - 2 x5 0
三、实验仪器设备与材料
主流微型计算机
四、实验原理
1、Guass 列选主元消去法
对于 AX =B
~
A B 是上三角矩阵。即：
1）、消元过程：将（A|B）进行变换为 ( ~ | ~) ，其中 A
⎛ a a a
b ⎫ ⎛ 1 a a b ⎫
11 12
1n 1 ⎪ 12
1n 1 ⎪
a21
a22
a2 n
b2 ⎪
0 1
a2 n
b2 ⎪
⎪ →
⎪
⎝ an1
an 2
⎪
⎪
k 从 1 到 n-1 a、 列选主元
选取第 k 列中绝对值最大元素 max aik k ≤i≤n
作为主元。
b、 换行
akj
bk
c、 归一化
⇔ aij , j = k + 1, , n
⇔ bi
d、 消元
akj / akk
bk / akk
⇒ akj , j = k + 1, , n
⇒ bk
aij - aik akj ⇒ aij , i = k + 1, , n; j = k + 1, , n bi - aik bk ⇒ bi , i = k + 1, , n
2）、回代过程：由 ( ~ | ~) 解出 x , x
, , x 。
A B
bn / ann ⇒ xn n
n n-1 1
bk -
∑ akj x j ⇒ xk , k = n - 1, ,2,1
j =k +1
2、追赶法
线性方程组为：
⎛ a c
⎫⎛ x ⎫
⎛ f ⎫
1 1
⎪ 1 ⎪
1 ⎪
b2
a2 c2
⎪ x2 ⎪
⎪ ⎪
f 2 ⎪
⎪
⎪
⎪ = ⎪
bn -1
an -1
cn-1
⎪
⎪ x
⎪
⎪
⎪
⎝ bn
an ⎭⎝ xn ⎭
⎝ f n ⎭
做 LU 分解为：
⎛α
⎫ ⎛ 1 β ⎫
1 ⎪ 1 ⎪
γ 2 α 2
⎪ 1 β 2 ⎪
⎪ ⎪
L = 3 3
⎪, R = ⎪
⎪
⎪
⎪
⎪
1 n-1 ⎪
⎪
⎝ γ n
α n ⎭ ⎝ 1 ⎭
分解公式：
⎧
⎪ i i
(i = 2,3, , n)
⎨α1 = b1 ,α i = bi - γ i β i -1
⎪
(i = 2,3, , n)
⎪⎩ α i
(i = 1,2, , n - 1)
则
Ax = f
⇒ LUx = f
⎧Ly = f
⇒ ⎨
⎩Ux = y
回代公式：
⎪ 1
⎪
⎨
⎩⎪ i
α1
= f i - γ i yi -1
α i
(i = 2,3, , n)
⎧xn = yn
⎨
⎩xi = yi - β i xi +1
(i = n - 1, n - 2, ,1)
五、实验步骤
1、理解并掌握全选主元消去法与高斯-塞德尔迭代法公式；
2、画出全选主元消去法与高斯-塞德尔迭代法的流程图
3、使用 C 语言编写出相应的程序并调试验证通过
六、实验报告要求
1、统一使用《武汉科技大学实验报告》本书写，实验报告的内容要求有：实验目的、 实验内容、程序流程图、源程序、运行结果及实验小结六个部分。
2、源程序需打印后粘贴在实验报告册内；
3、运行结果以屏幕截图形式保存并打印后粘贴在实验报告册内。
七、实验注意事项
注意如何定义数据结构以保存矩阵和解以降低算法的复杂性。
八、思考题
若使用全主元消去法，在编程中应如何记录保存对于未知数的调换。
实验二 插值方法
一、实验目的
掌握拉格郎日插值法与牛顿插值法构造插值多项式。
二、实验内容
分别写出拉格郎日插值法与牛顿插值法的算法，编写程序上机调试出结果，要求所编程 序适用于任何一组插值节点，即能解决这一类问题，而不是某一个问题。实验中以下列数据 验证程序的正确性。
已知下列函数表
求 x=0.5635 时的函数值。
三、实验仪器设备与材料
主流微型计算机
四、实验原理
已知 n 个插值节点的函数值，则可由拉格郎日插值公式与牛顿插值公式构造出插值多项 式，从而由该插值多项式求出所要求点的函数值。拉格郎日插值公式与牛顿插值公式如下：
1、Lagrange 插值公式
n
Ln ( x) = l0 ( x) y0 + l1 ( x) y1 + ... + ln ( x) yn = ∑ yk lk ( x)
k =0
( x - x
)(x - x ) ( x - x
)(x - x
) ( x - x )
n x - x
l ( x) = 0 1 k -1 k +1 n = ∏ j
x0 )(xk -
x1 ) ( xk -
xk -1
)(xk -
xk +1
) ( xk -
xn )
j =0
j ≠ k
xk - x j
2、Newton 插值公式
N n ( x) =
f ( x0 ) + f [ x0 , x1 ](x - x0 ) + f [ x0 , x1 , x2 ](x - x0 )(x - x1 )
+ + f [ x0 , x1 , xn ](x - x0 )(x - x1 ) ( x - xn-1 )
五、实验步骤
1、理解并掌握拉格郎日插值法与牛顿插值法的公式；
2、画出拉格郎日插值法与牛顿插值法算法的流程图；
3、使用 C 语言编写出相应的程序并调试验证通过。
六、实验报告要求
1、统一使用《武汉科技大学实验报告》本书写，实验报告的内容要求有：实验目的、
实验内容、程序流程图、源程序、运行结果及实验小结六个部分。
2、源程序需打印后粘贴在实验报告册内；
3、运行结果以屏幕截图形式保存并打印后粘贴在实验报告册内。
七、实验注意事项
Newton 插值法在编程时应注意定义何种数据结构以保存差商。
八、思考题
比较 Lagrange 插值法与 Newton 插值法的异同。
实验三 数值积分
一、实验目的
掌握复化梯形法与龙贝格法计算定积分。
二、实验内容
分别写出变步长梯形法与 Romberge 法计算定积分的算法，编写程序上机调试出结果， 要求所编程序适用于任何类型的定积分，即能解决这一类问题，而不是某一个问题。实验中
以下列数据验证程序的正确性。
求 1sin x dx, ε ≤ 0.00001 。
⎰0 x
三、实验仪器设备与材料
主流微型计算机
四、实验原理
通过变步长梯形法与龙贝格法，我们只要知道已知 n 个求积节点的函数值，则可由相应 的公式求出该函数的积分值，从而不需要求该函数的原函数。变步长梯形法与龙贝格法公式 如下：
1、变步长梯形法
Tn = [ f ( xi ) + f ( xi +1 )]
i =0
h n-1
= [ f (a) + 2∑ f ( xi ) + f (b)]
2 i =1
1 h n-1
T2 n = Tn + ∑ f ( xi +1 / 2 )
2 2 i =0
用 T2n - Tn
≤ ε 来控制精度
2、龙贝格法
1
h n-1
T2 n = Tn + ∑ f ( xi +1 / 2 )
2 2 i =0
S n =
4 1
T2n - Tn
3 3
Cn =
16 1
S 2n - S n
15 15
64 1
Rn =
C2n - Cn
63 63
用 R2n - Rn
≤ ε 来控制精度
五、实验步骤
1、理解并掌握变步长梯形法与龙贝格法的公式；
2、画出变步长梯形法与龙贝格法的流程图
3、使用 C 语言编写出相应的程序并调试验证通过
六、实验报告要求
1、统一使用《武汉科技大学实验报告》本书写，实验报告的内容要求有：实验目的、 实验内容、程序流程图、源程序、运行结果及实验小结六个部分。
2、源程序需打印后粘贴在实验报告册内；
3、运行结果以屏幕截图形式保存并打印后粘贴在实验报告册内。
七、实验注意事项
1 sin x
在 ⎰0 x
的定义。
dx 积分中，被积函数在 x=0 点函数值为 1，对该点在程序设计中应注意对其
八、思考题
使用复化梯形法与复化 Simpson 法来计算该问题有何缺点？
实验四 常微分方程的数值解
一、实验目的
掌握改进欧拉法与四阶龙格-库塔求解一阶常微分方程的初值问题。
二、实验内容
分别写出改进欧拉法与四阶龙格-库塔求解的算法，编写程序上机调试出结果，要求所 编程序适用于任何一阶常微分方程的数值解问题，即能解决这一类问题，而不是某一个问题。 实验中以下列数据验证程序的正确性。
⎧ y' = - xy 2
求 ⎨
⎩ y(0) = 2
(0 ≤ x ≤ 5)
步长 h=0.25。
三、实验仪器设备与材料
主流微型计算机
四、实验原理
常微分方程的数值解主要采用“步进式”，即求解过程顺着节点排列次序一步一步向前 推进，在单步法中改进欧拉法和四阶龙格-库塔法公式如下：
1、改进欧拉法
yn+1 = yn + hf (xn , yn )
yn+1
= yn
+ h [ f ( x
2 n
, yn
) + f ( x
n+1
, yn+1 )]
2、四阶龙格-库塔法
⎧ h
⎪ yn +1 = yn + 6 (k1 + 2k 2 + 2k3 + k 4 )
⎪
⎪k1 =
⎪
⎨k 2 =
⎪
⎪
f ( xn , yn )
f ( x + h
n 2
h
, yn +
h
k1 )
2
h
⎪k3 =
⎪
⎪⎩k 4 =
f ( xn + , yn + k 2 )
2 2
f ( xn + h, yn + hk3 )
五、实验步骤
1、理解并掌握改进欧拉法与四阶龙格-库塔法的公式；
2、画出改进欧拉法与四阶龙格-库塔法的流程图
3、使用 C 语言编写出相应的程序并调试验证通过
六、实验报告要求
1、统一使用《武汉科技大学实验报告》本书写，实验报告的内容要求有：实验目的、
实验内容、程序流程图、源程序、运行结果及实验小结六个部分。
2、源程序需打印后粘贴在实验报告册内；
3、运行结果以屏幕截图形式保存并打印后粘贴在实验报告册内。
七、实验注意事项
⎧ y' = - xy 2
⎨
⎩ y(0) = 2
精度的变化
八、思考题
(0 ≤ x ≤ 5)
的精确解为 y = 2 /(1 + x 2 ) ，通过调整步长，观察结果的
如何对四阶龙格-库塔法进行改进，以保证结果的精度。
实验五 迭代法解线性方程组与非线性方程
一、实验目的
掌握高斯-塞德尔迭代法求解线性方程组与牛顿迭代法求方程根。
二、实验内容
分别写出高斯-塞德尔迭代法与牛顿迭代法的算法，编写程序上机调试出结果，要求所 编程序适用于任何一个方程的求根，即能解决这一类问题，而不是某一个问题。实验中以下 列数据验证程序的正确性。
1、高斯-塞德尔迭代法求解线性方程组
⎡ 7 2
1 - 2⎤ ⎡ x1 ⎤
⎡ 4 ⎤
⎢ 9 15
3 - 2⎥ ⎢ x ⎥
⎢ 7 ⎥
⎢ ⎥ ⎢
2 ⎥ = ⎢ ⎥
⎢- 2
- 2 11
5 ⎥ ⎢ x3 ⎥
⎢- 1⎥
⎢ ⎥ ⎢ ⎥ ⎢ ⎥
⎣ 1 3
2 13 ⎦ ⎣ x4 ⎦
⎣ 0 ⎦
2、用牛顿迭代法求方程 x3 - x - 1 = 0 的近似根，ε ≤ 0.00001，牛顿法的初始值为 1。
三、实验仪器设备与材料
主流微型计算机
四、实验原理
二分法通过将含根区间逐步二分，从而将根的区间缩小到容许误差范围。牛顿通过迭代 的方法逐步趋进于精确解，该两种方法的公式如下：
1、高斯-塞德尔迭代法
1）判断线性方程组是否主对角占优
n
∑ aij
j =1
j ≠i
≤ aii , i = 1,2, , n
2）直接分离 xi，即
n
xi = (di - ∑bij x j ) / aii , i = 1,2, , n
j =1
建立高斯-塞德尔迭代格式为：
( k +1)
i -1
= (d i - ∑ aij x j
( k +1)
n
- ∑ aij x j
( k )
) / aii , i = 1,2, , n
j =1
j =i +1
3）取初值迭代求解至所要求的精度为止。
2、牛顿法
- f ( xk )
k +1
k f '( x )
五、实验步骤
1、理解并掌握二分法与牛顿法的公式；
2、画出二分法与牛顿法的流程图
3、使用 C 语言编写出相应的程序并调试验证通过
六、实验报告要求
1、统一使用《武汉科技大学实验报告》本书写，实验报告的内容要求有：实验目的、 实验内容、程序流程图、源程序、运行结果及实验小结六个部分。
2、源程序需打印后粘贴在实验报告册内；
3、运行结果以屏幕截图形式保存并打印后粘贴在实验报告册内。
七、实验注意事项
对于二分法应注意二分后如何判断根的区间，对于牛顿法注意如何确定迭代过程的结 束
八、思考题
若使用牛顿法是发散的，如何对牛顿法进行改进以保证其收敛性。
前三个实验的程序代码（C/C++）和运行结果截图
#include <iostream>
#include <stdlib.h>
#include <math.h>
using namespace std;
Gauss 全选主元解方程组的源程序及运行结果
class Matrix
{
public: Matrix();
~Matrix();
void SetMatrix(const int n,const double esp1);//构造线性方程组相应的矩阵，n
为方程的未知数数目，esp1 为要求的精度
void Max(const int r);//全选主元
void ChangeRC(const int r);//根据主元变换矩阵的行或列
void Eliminate(const int r);//处理消元工作
void Result()const;//计算方程的解
void Calculate();
int GetRank()const;//返回矩阵的行数
double GetX(const i)const;//确定方程组的第 i 个解(1<=i<=N)
private:
//指针 a 和 b 分别用于存储方程组的未知数系数和方程"="右边的常数，esp 存
//储精度
double *a,*b,esp;
//指针 flag 用于记录方程组解的顺序
int *flag;
//以下的结构体用于在全选主元中记录最大主元的位置
struct coordinate
{
int row,column;
}location;
int N;//方程组未知数的数目
};
int main()
{
int n;
double esp1; Matrix matrix; do{
cout<<"请输入方阵的阶数：";
cin>>n;
if(n<0) n=0;//如何控制非法字符的输入？？？？？？？？？？
}while(n==0);
do
{
cout<<"请输入计算精度：";
cin>>esp1;
if(esp1<0) esp1=0;//输入不合法的精度就把精度置 0
}while(esp1==0);
cout<<"输入线性方程组的增广矩阵:\n"; matrix.SetMatrix(n,esp1);//设置矩阵内的数据 matrix.Calculate();//计算方程组的解
//输出方程组的解
cout<<"\n\n 方程组的解如下:\n";
for(int i=1;i<=matrix.GetRank();++i)
cout<<"X["<<i<<"]:"<<matrix.GetX(i)<<endl;
return 0;
}
Matrix::Matrix()
{ //将 Matrix 类的数据成员初始化
a=NULL; b=NULL; flag=NULL; location.row=0; location.column=0; esp=0;
N=0;
}
Matrix::~Matrix()
{ //释放指针 a、b 和 flag 指向的内存空间
delete[]a;
delete[]b;
delete[]flag;
}
void Matrix::SetMatrix(const int n,const double esp1)
{
N=n;
esp=esp1;
a=new double[N*N]; b=new double[N]; flag=new int[N];
//判断是否成功分配存储区
if(a==NULL||b==NULL||flag==NULL)
{
cout<<"分配存储区失败！\n";
exit(EXIT_FAILURE);
}
//读取线性方程组的增广矩阵
for(int i=0;i<N;++i)
{
for(int j=0;j<N;++j) cin>>*(a+i*N+j);
cin>>*(b+i);
}
//flag 中存储的值对应相应的 x 值，当方程的解由于列变换交换后，flag 中
//的值也相应交换，最后用于恢复解的顺序
for( i=0;i<N;++i) *(flag+i)=i;
}
void Matrix::Max(const int r)
{
double max=0;
for(int i=r;i<N;++i)
for(int j=r;j<N;++j)
{
if(max<fabs(*(a+i*N+j)))
{
max=fabs(*(a+i*N+j));
//设定最大主元的行、列
location.row=i;
location.column=j;
}
}
//最大主元小于输入的精度时，认为方程组无解，退出程序
if(max<=esp)
{
cout<<"方程组无解！\n";
exit(EXIT_FAILURE);
}
}
void Matrix::ChangeRC(const int r)
{
double temp;
//如果最大主元所在的行不在当前行，则进行行变换
if(location.row!=r)
{
for(int i=r;i<N;++i)
{
temp=*(a+r*N+i);
*(a+r*N+i)=*(a+location.row*N+i);
*(a+location.row*N+i)=temp;
}
temp=b[r]; b[r]=b[location.row]; b[location.row]=temp;
}
//若最大主元所在的列不在当前的 r 列，则进行列变换
if(location.column!=r)
{
for(int i=r;i<N;++i)
{
temp=*(a+i*N+r);
*(a+i*N+r)=*(a+i*N+location.column);
*(a+i*N+location.column)=temp;
}
//交换 flag 中的元素来标记方程解的位置变化
int temp1;
temp1=*(flag+r);
*(flag+r)=*(flag+location.column);
*(flag+location.column)=temp1;
}
}
void Matrix::Eliminate(const int r)
{
if(fabs(*(a+N*r+r))<=esp)
{
cout<<"方程组无解！\n";
exit(EXIT_FAILURE);
}
for(int i=r+1;i<N;++i)
{
for(int j=r+1;j<N;++j)
(*(a+i*N+j))-=(*(a+i*N+r))*(*(a+r*N+j))/(*(a+r*N+r)); (*(b+i))-=(*(b+r))*(*(a+i*N+r))/(*(a+r*N+r));
}
}
void Matrix::Result()const
{
if(fabs(*(a+N*(N-1)+N-1))<=esp)
{
cout<<"方程组无解！\n";
exit(EXIT_FAILURE);
}
double temp;
*(b+N-1)/=(*(a+N*(N-1)+N-1)); //求出 X[N-1]
//依次求出 X[i](i=N-2,N-3····，1)
for(int i=N-2;i>=0;--i)
{
temp=0;
for(int j=i+1;j<N;++j) temp+=((*(a+i*N+j))*(*(b+j)));
*(b+i)=(*(b+i)-temp)/(*(a+i*N+i));
}
//根据 flag 中的数据用冒泡排序法恢复方程组解的次序
for(i=0;i<N-1;++i)
for(int j=0;j<N-i-1;++j)
if(*(flag+j)>*(flag+j+1))
{
int temp1;
//交换解的顺序
temp=*(b+j);
*(b+j)=*(b+j+1);
*(b+j+1)=temp;
//交换用于标记的元素的顺序
temp1=*(flag+j);
*(flag+j)=*(flag+j+1);
*(flag+j+1)=temp1;
}
}
void Matrix::Calculate()
{ //根据矩阵行数重复进行寻找最大主元、变换行或列、消元
for(int i=0;i<GetRank()-1;++i)
{ Max(i); ChangeRC(i); Eliminate(i);
}
Result();
}
int Matrix::GetRank()const
{
return N;//返回矩阵的行数
}
double Matrix::GetX(const int i)const
{
return *(b+i-1);
}
运行结果
追赶法求解方程组的算法：
1． 输入方程组的维数 n，将主对角元素 b(i)(i=0:n-1)，，主对角元素左边的元素
a(i)(i=0:n-2)，主对角元素右边的元素 c(i)(i=0:n-2)，右端项的元素 f(i)(i=0:n-1)
2． 对方程组的系数矩阵作 Crout 分解, α (0)=b(0),对于 i=0:n-2,
c(i):=β (i):= c(i)/b(i), b(i+1):=α (i+1):= b(i+1)-a(i)* β (i)
3.解方程组 Ly=f
b(0):=y(0):=[f(0)/ α (0)]:=[f(0)/b(0)]
对于 i=1：n-1， b(i):=y(i):=[f(i)-a(i-1)*y(i-1)]/b(i)
4..解方程组 Ux=y
a(n-1):=x(n-1):=b(n-1)
对于 i=n-2:0,a(i)=x(i):=b(i)-c(i)*a(i+1);
5．输出方程组的解 a(0:n-1)
#include<iostream>
#include<stdlib.h>
用追赶法求解方程组的源程序及运行结果
using namespace std;
class MatrixThr
{
public: MatrixThr();
~MatrixThr();
void SetMatrixThr(const int n);//设置三对角矩阵的数据
void Result();//计算三对角矩阵的解
double GetX(const int i)const;//取得第 i 个解,i 从 1 开始
int GetN() const;//返回未知数的数目
private:
int N;//N 为未知数的数目
//b 为矩阵主对角线的元素首地址，a 为主对角线左边一斜条元素的首地址，
//c 为主对角线右边一斜条元素首地址，f 为方程组的常数首地址
double *a,*b,*c,*f;
};
int main()
{
MatrixThr matrix;
int n;
do{
cout<<"输入三对角方程组的变量的数目 N:";
cin>>n;
}while(n<3);
cout<<"请依次输入三对角方程组每行的数据(0 元素除外):\n"; matrix.SetMatrixThr(n); //计算方程组的解 matrix.Result();//输出方程组的解
cout<<"方程的解如下:\n";
for(int i=1;i<=matrix.GetN();++i)
cout<<"X["<<i<<"]:"<<matrix.GetX(i)<<endl;
return 0;
}
MatrixThr::MatrixThr()
{ //初始化相关数据
N=0; a=NULL; b=NULL; c=NULL;
f=NULL;
}
MatrixThr::~MatrixThr()
{ //释放分配的内存空间
delete []a; delete []b; delete []c;
delete []f;
}
void MatrixThr::SetMatrixThr(const int n)
{ //根据输入的未知数个数设置矩阵的数据
N=n;
a=new double[N]; b=new double[N]; c=new double[N]; f=new double[N];
//若内存分配失败，退出程序
if(a==NULL||b==NULL||c==NULL||f==NULL)
{
cout<<"初始化分配存储空间失败！\n";
exit(EXIT_FAILURE);
}
//依次输入三对角矩阵每行的元素
cin>>*b>>*c>>*f;
for(int i=1;i<N-1;++i) cin>>*(a+i-1)>>*(b+i)>>*(c+i)>>*(f+i);
cin>>*(a+i-1)>>*(b+i)>>*(f+i);
}
void MatrixThr::Result()
{ //对系数矩阵 A 作 Crout 分解
for(int i=0;i<N-1;++i)
{ //将 U 中的α 存于指针 b 中，L 中的β 存于指针 c 中
*(c+i)/=(*(b+i));
*(b+i+1)-=(*(a+i))*(*(c+i));
}
//解方程组 Ly=f，求得的 y 值存于指针 b 中
*b=(*f)/(*b);
for(i=1;i<N;++i) *(b+i)=((*(f+i))-(*(a+i-1))*(*(b+i-1)))/(*(b+i));
//解方程组 Ux=y，求得的 x 值存于指针 a 中
*(a+N-1)=*(b+N-1);
for(i=N-2;i>=0;--i) *(a+i)=*(b+i)-(*(c+i))*(*(a+i+1));
}
int MatrixThr::GetN()const
{
return N;
}
double MatrixThr::GetX(const int i)const
{
return *(a+i-1);
}
运行结果
#include<iostream>
Lagrange 插值法的源程序：
#include<stdlib.h>
using namespace std;
class Lagrange
{
public: Lagrange();
~Lagrange();
void SetLagrange(const int n);//根据用户的输入设置 Lagrange 类中的插值点数
据
bool Exist(const double x,const int i);//检测是否输入了与前 i 个插值结点横坐
标相同的点
int GetN()const;//获取插值结点的数目
void Calculate(const double a);//计算横坐标 a 对应的函数值
double GetResult()const;//返回计算的函数值
private:
int N;//插值结点的数目
double *x,*y,zx,zy;//x、y 分别用于存储插值点的数据，zx、zy 表示所求的坐
标点
};
int main()
{
int n=2; Lagrange L; do
{
if(n<2) cout<<"用于模拟曲线的插值点数目 N>2"<<endl;
cout<<"请输入用于模拟曲线的插值点数目 N:";
cin>>n;
}while(n<2);
cout<<"请输入各插值点横、纵坐标的数据\n"; L.SetLagrange(n);
double a;
cout<<"请输入所求点的横坐标 X:";
cin>>a; L.Calculate(a);
cout<<"横坐标"<
return 0;
}
Lagrange::Lagrange()
{ //初始化相关数据
N=0;
x=y=NULL;
zx=zy=0;
}
Lagrange::~Lagrange()
{ //释放分配给指针的内存空间
delete []x;
delete []y;
}
bool Lagrange::Exist(const double a,const int i)
{ //遍历以输入的插值点，看是否重复插入横坐标相同的插值点
for(int j=0;j<i;++j)
if(a==*(x+j)) return true;
return false;
}
void Lagrange::SetLagrange(const int n)
{
N=n;
x=new double[N];
y=new double[N];
//判断是否成功为指针分配了内存空间
if(x==NULL||y==NULL)
{
cout<<"分配存储空间失败！"<<endl;
exit(EXIT_FAILURE);
}
//输入插值点
for(int i=0;i<GetN();++i)
{
cin>>*(x+i)>>*(y+i);
//如果不是输入第一个坐标值，则会对输入的横坐标进行合法性检测
while(i!=0&&Exist(*(x+i),i)==true)
{
cout<<"输入了重复的横坐标-------请重新输入第"<<i+1<<"个插入结 点的数据\n";
cin>>*(x+i)>>*(y+i);
}
}
}
void Lagrange::Calculate(const double a)
{
zx=a;
for(int i=0;i<GetN();++i)
{ //计算插值基函数
double temp=1;
for(int j=0;j<i;++j) temp*=((zx-*(x+j))/(*(x+i)-*(x+j)));
for(++j;j<GetN();++j) temp*=((zx-*(x+j))/(*(x+i)-*(x+j)));
zy+=(*(y+i)*temp);
}
}
int Lagrange::GetN()const
{
return N;//返回插值点的数目
}
double Lagrange::GetResult() const
{
return zy;//返回求得的函数值
}
Lagrange 插值法的运行结果
#include <iostream>
using namespace std;
Newton 插值法的源程序：
class Matrix
{
public: Matrix();
~Matrix();
void SetMatrix(const int n);//根据用户输入的插值点的数据设置计算结果的矩
阵
int GetN()const;//返回插值点的数目
void Calculate();//计算差商
double GetResult(double x)const;//根据输入的横坐标求函数值，返回运算结果
private:
double *a,*f;//a 和 f 分别用于存储插值点的横坐标和相应的函数值
int N;//记录插值点的数目
};
int main()
{
Matrix matrix;
int n;
do{
cout<<"输入插值点的数目 N：";
cin>>n;
}while(n<2);
cout<<"输入插值点的数据："<<endl;
matrix.SetMatrix(n);
matrix.Calculate();
double x;
cout<<"输入所求的横坐标：";
cin>>x;
cout<<"所求的函数值为："<<matrix.GetResult(x)<<endl;
return 1;
}
Matrix::Matrix()
{ //初始化数据
a=f=NULL;
N=0;
}
Matrix::~Matrix()
{
//释放指针 a、f 指向的内存空间
delete []a;
delete []f;
}
void Matrix::SetMatrix(int n)
{
N=n;
//为插值点创建动态数组
a=new double[N];
f=new double[N];
//输入插值点的数据
for(int i=0;i<GetN();++i) cin>>*(a+i)>>*(f+i);
}
int Matrix::GetN() const
{
return N;
}
void Matrix::Calculate()
{ //将差商存储在一个一维数组内
for(int i=0;i<GetN();++i)
for(int j=GetN()-1;j>i;--j)
*(f+j)=(*(f+j)-*(f+j-1))/(*(a+j)-*(a+j-i-1));
}
double Matrix::GetResult(double x) const
{
//利用差商和插值点的横坐标及第一个插值点的纵坐标计算函数值
double result=*f;//指针 f 指向第一个插值点的纵坐标
for(int i=1;i<GetN();++i)
{
double temp=1;
for(int j=0;j<i;++j) temp*=(x-*(a+j));//计算(x-x0)*(x-x1)*···*(x-x(i-1))
result+=(*(f+i))*temp;//指针(f+i)指向第 i 号差商
}
return result;
}
Newton 插值法的运行结果
变步长梯形法求定积分的源程序：
#include<iostream>
#include<math.h>
using namespace std;
double function(const double x);//求被积函数的值并返回
//accumulate()为求定积分的函数，a、b 分别为积分的上下限,默认精度为 0.00001 double accumulate(const double a,const double b,const double eps=0.00001);
int main()
{
double a,b,eps;//a,b 分别为定积分的上限和下限,h 为步长,eps 为要求的精度
a=0; b=1; eps=0.00001;
cout<<"(sinX/X)在 0 到 1 上的积分为:"<<accumulate(a,b)<<endl;
return 0;
}
double function(const double x)
{
if(x==0) return 1;//x 为 0 时函数值为 1 return(sin(x)/x);
}
double accumulate(const double a,const double b,const double eps)
{
int n=1;
double h=b-a;//h 为步长
double T1=h*(function(a)+function(b))/2;
double T2=T1/2+h*function(a+h/2)/2;
//计算结果不满足精度，则继续
while(fabs(T2-T1)>eps)
{
n*=2; h/=2;//步长折半 T1=T2;
//利用 T1 计算 T2 double temp=0;
for(int i=1;i<=n;++i) temp+=function(a+(i-1.0/2)*h); T2=T1/2+h*temp/2;
}
return T2;//返回步长 h 最合适时定积分的近似值
}
1 sin x
变步长梯形法求定积分 ⎰0
的运行结果
x
Romberg 算法求定积分的源程序：
#include <iostream>
#include <math.h>
using namespace std;
//求被积函数的值并返回
double fun(const double x);
//Romberge()为求定积分的函数，a、b 分别为积分的上下限,默认精度为 0.00001 double Romberg(const double a,const double b,const double eps=0.0001);
//函数 Tm()为 T-数表的计算公式
double Tm(const double T1,const double T2,const int m);
int main()
{
double a,b,eps;//a,b 分别为定积分的上限和下限,h 为步长,eps 为要求的精度
a=0; b=1; eps=0.00001;
cout<<"所求积分为:"<<Romberg(a,b,eps)<<endl;
return 0;
}
double fun(const double x)
{
if(x==0) return 1;//x 为 0 时函数值为 1 return(sin(x)/x);//返回被积函数的值
}
double Tm(const double T1,const double T2,const int m)
{ //根据 T1、T2 和 m 计算出下一个数并返回
return((pow(4,m)*T1-T2)/(pow(4,m)-1));
}
double Romberg(const double a,const double b,const double eps)
{
double T[10];//用于存储 T-数表某一行的数据
int l=1;
T[0]=(b-a)*(fun(a)+fun(b))/2;
double T1,T2;
int flag;//flag 用于标记十分求得符合精度的值
do{
flag=0;
int n=pow(2,l-1);//n 表示区间等分数
double h=(b-a)/n;//h 为步长
double temp=0;
for(int i=0;i<n;++i) temp+=fun(a+(i+1.0/2)*h); T1=T[0]/2+h*temp/2;//利用变步长梯形公式计算 T-数表第（l+1）行的第一个数据
//将 T-数表第 l+1 行的前 l-1 个数据存入数组
for(int m=1;m<l;++m)
{
T2=Tm(T1,T[m-1],m); T[m-1]=T1;
T1=T2;
}
T2=Tm(T1,T[l-1],l);//计算 T-数表第（l+1）行的最后一个数据
if(fabs(T2-T[l-1])>eps)//比较 T-数表第（l+1）行和第 i 行的最后一个数据的差值
{ ////将 T-数表第 l+1 行的最后两个个数据存入数组
T[l-1]=T1; T[l]=T2;
++l;
flag=1;//表示没有求得符合精度的值
}
}while(flag);
return T2;//返回符合精度的解
}
Romberg 算法求定积分运行结果
第 1 章 绪论
《数值计算基础》课程复习指导
1、有效数字、绝对误差、绝对误差限、相对误差、相对误差限的概念；
2、有效数字与绝对误差，有效数字与相对误差的关系；
3、如何判断有效数字，如何估算绝对误差限与相对误差限。 第 2 章 解线性方程组的直接法
1、直接法解线性方程组的思想，如何使用高斯消去法、列选主元消去法、全选主元消
去法；
2、直接三角分解法解线性方程组的思想，矩阵的 LU 分解的条件，如何对矩阵进行 LU
分解。如何使用直接三角分解法、平方根法和追赶法法解线性方程组。 第 3 章 代数插值法与最小二乘法
1、插值的基本概念，插值问题的存在且唯一性；
2、如何使用待定系数法、拉格郎日插值法、牛顿插值法构造插值多项式及确定余项；
3、 li (x),ω(x) 的性质及应用；
4、差商的定义、性质及应用；
5、如何使用分段线性插值及确定余项；
5、如何使用待定系数法构造埃尔米特插值多项式及确定余项；
6、如何使用曲线拟合的最小二乘法进行线性拟合。
第 4 章 数值积分与数值微分
1、机械求积与代数精度的概念，如何判定一个求积公式的代数精度；
2、如何通过代数精度法与插值法构造求积公式；
3、牛顿-柯特斯公式的定义及构造的方法，牛顿-柯特斯系数的性质；如何使用梯形公 式、辛卜生公式、柯特斯公式计算定积分及确定余项；
4、复化求积法；如何使用复化梯形公式、复化辛卜生公式、复化柯特斯公式计算定积 分及余项；
5、变步长求积法的思想，如何使用变步长梯形求积法和龙贝格求积法计算定积分；
6、高斯求积公式的定义及构造方法；
7、数值微分的数值方法，如何使用二点公式、三点公式计算微分。 第 5 章 常微分方程数值解
1、常微分方程数值解法的基本思想，欧拉方法、后退欧拉方法、梯形方法、改进欧拉 方法公式的构造方法；
2、如何使用欧拉方法、后退欧拉方法、改进欧拉方法、龙格-库塔方法计算常微分方程；
3、局部截断误差与方法精度的定义，如何判断一个方法是几阶方法。 第 6 章 逐次逼近法
1、向量范数与矩阵范数的基本概念，常用的向量范数与矩阵范数，矩阵的谱半径的基 本概念。
2、如何使用简单迭代法和高斯-塞德尔迭代法解线性方程组，如何判断迭代法的收敛性。
3、如何使用简单迭代法、牛顿迭代法解非线性方程。如何判断迭代格式的收敛阶。
《数值计算基础》考试样卷
一、单项选择题（每小题 3 分，共 15 分）
1、数值 x 的近似值 x*=0.1215×10－2，若满足 x - x *
≤ ( )，则称 x 有 4 位有效数字.
(A)
1 ×10－3 (B)
2
1 ×10－4 (C)
2
1 ×10－5 (D)
2
1 ×10－6
2
2、若 Ak 为矩阵 A 的 k 阶主子矩阵，则矩阵 A 满足（ ）时，则存在唯一单位下三角阵
L 和上三角阵 R ，使 A = LR 。
(A)
A ≠ 0
(B) 某个 Ak ≠ 0
(C) Ak ≠ 0
(k = 1, n - 1)
(D)
Ak ≠ 0
(k = 1, , n)
3、通过四个互异节点的插值多项式 P(x)，只要满足( ), 则 P(x)是不超过一次多项式。
(A) 初始值 y0=0 (B) 所有一阶均差为 0 (C) 所有二阶均差为 0 (D) 所有三阶均差为 0
4、牛顿切线法求解方程 f(x)=0 的近似根，若初始值 x0 满足( )，则解的迭代数列一定收敛。
(A) f (x0 ) f ' ( x0 ) <0 (B)
f (x0 ) f ' ( x0 ) >0
(C) f (x0 ) f ' ( x0 ) ≤0 (D) f (x0 ) f ' ( x0 ) ≥0
5、改进欧拉法的平均形式公式是( )
⎧
⎪ y p = yk + hf ( xk , yk )
⎪
(A) ⎨ yc = yk + hf ( xk , y p )
⎪
⎧
⎪ y p = yk + hf ( xk +1 , yk )
⎪
(B) ⎨ yc = yk + hf ( xk +1 , y p )
⎪
⎩⎪ y
k +1
= 1 ( y
2 p
+ yc )
⎩⎪ y
k +1
= 1 ( y
2 p
+ yc )
⎧
⎪ y p = yk + hf ( xk , yk )
⎪
(C) ⎨ yc = yk + hf ( xk +1 , y p )
⎪
⎧
⎪ y p = yk + hf ( xk , yk )
⎪
(D) ⎨ yc = yk + hf ( xk +1 , y p )
⎪
⎩⎪ y
k +1
= h ( y
2 p
+ yc )
⎩⎪ y
k +1
= 1 ( y
2 p
+ yc )
二、填空题（每小题 3 分，共 15 分）
1、sin1 有 2 位有效数字的近似值 0.84 的相对误差限是 .
2、设 f(x)可导，求方程 x=f(x) 根的牛顿迭代格式是 .
3、设 f (x) = 2x 2 + 4 ，则 f [1,2] = .
4、在区间 [a, b] 上的插值型求积公式系数 A0 , A1 , ┅ , An 满足 A0 + A1 + ┅＋ An ＝ .
5、二阶龙格－库塔法的局部截断误差是 .
三、解答题（每小题 10 分，共 50 分）
1、用列主元消去法解线性方程组
⎡ 2 4 0⎤ ⎡ x1 ⎤ ⎡5⎤
⎢ ⎥ ⎢ ⎥ ⎢ ⎥
⎢ 3 -1 1⎥ ⎢ x2 ⎥ = ⎢9⎥
-2 -2 0 x3
⎢⎣3
2、用牛顿法求 6 的近似值，取初始值 x0 = 2 ，进行二次迭代。
3、已知有 y=f(x)的函数表如下
求其代数插值多项式并给出其余项。
4、给出数值积分公式：
h
⎰-h
f ( x)dx ≈ Af (-h) + Bf (1 h)
3
确定 A、B 使得该数值积分公式的代数精度尽可能的高，并确定其代数精度为多少？
5、用欧拉法解初值问题，要求保留 4 位有效数字。
⎧ y ' = x + y(0 ≤ x ≤ 1, h = 0.5)
⎨
⎩ y(0) = 1
四、综合题（每小题 10 分，共 20 分）
1、试利用数值积分的方法推导求解初值问题的梯形公式为
yn+1
= yn
+ h [ f ( x
2 n
, yn
) + f ( x
n+1
, yn+1 )]
，并证明该方法是二阶方法。
2、设 l0(x)是以 n+1 个互异点 x0,x1,x2,…,xn 为节点的拉格朗日插值基函数
l0 ( x) =
( x - x1 )(x - x2 )...(x - xn )
( x0 - x1 )(x0 - x2 )...(x0 - xn )
试利用牛顿插值法证明：
l0 ( x) = 1 +
( x - x0 )
( x0 - x1 )
+ ( x - x0 )(x - x1 )
( x0 - x1 )(x0 - x2 )
+ ... +
( x - x0 )(x - x1 )...(x - xn-1 )
( x0 - x1 )(x0 - x2 )...(x0 - xn )
《数值计算基础》考试样卷
参考答案 一、单项选择题（每小题 3 分，共 15 分）
1、D 2、D 3、C 4、B 5、D
二、填空题（每小题 3 分，共 15 分）
1
1、
2 ⨯ 8
⨯ 10-2+1 =
1 ⨯ 10-1 = 0.00625
16
2、 x
= x - xk
3、 6
k +1
k 1 - f '( x )
4、 b-a
5、 O(h3)
三、解答题（每小题 10 分，共 50 分）
1、解：
⎡ 2 4 0 5⎤ ⎡ 3
-1 1 9⎤
r ( 2 ) r
⎢ 3 -1 1 9⎥ −r−2 ↔−r1 → ⎢ 2 4 0 5⎥ −−−3 −→
2 + - 1
⎢ ⎥ ⎢ ⎥ r + =2 r
⎢-2
-2 0 3⎥ ⎢-2
-2 0 3⎥
3 3 1
⎣ ⎦ ⎣ ⎦
⎡ ⎤ ⎡ ⎤
⎢3 -1 1 9 ⎥ ⎢3
-1 1 9 ⎥ 8 分
⎢ ⎥ 4 ⎢ ⎥
-2 -1⎥
r3 + 7 r2 → ⎢0 14
-2 -1⎥
⎢ 3 3
⎥ ⎢ 3 3 ⎥
⎢ ⎥ ⎢ ⎥
⎢0 -8 2
9 ⎥ ⎢0 0
2 59 ⎥
3 3 7 7
回代得 x3 =
2、解：
59 , x = 4, x = - 11 2 分
2 2 1 2
f ( x) = x 2 - 6，f ' ( x) = 2x,ϕ ( x) = x -
f ( x)
f ' ( x)
= 1 ( x + 6)，x
2
n+1
= 1 ( x + 6 ) 7 分
2 n x
n
x0 = 2
x = 1 (2 + 6 ) = 5 = 2.500 3 分
1 2 2 2
x = 1 ( 5 + 12 ) = 49 = 2.450
2 2 2 5 20
3、 解
法一： 待定系数法
，则 （3 分）
⎧a0 + a1 + a2 = 1
⎪
⎧a0 = 1
⎪
（3 分）
2
P2 ( x) = ∑ yi li ( x)
i =0
⎩a2 = 1
（1 分）
(3分)
= ( x - 2)( x - 3) ⋅1 + ( x - 1)( x - 3) ⋅ 3 + ( x - 1)(x - 2) ⋅ 7
(3分)
(1 - 2)(1 - 3)
= x 2 - x + 1
法三：Newton 插值法
(2 - 1)(2 - 3)
(3 - 1)(3 - 2)
(1分)
N 2 ( x) =
f ( x0 ) + f [ x0 , x1 ]( x - x0 ) + f [ x0 , x1 , x2 ]( x - x0 )( x - x1 )
（3 分）
= 1 + 2( x - 1) + ( x - 1)( x - 2)
= x 2 - x + 1
（4 分）
余项为 R2 ( x) =
4、 解：
f ''(ξ ) ⋅ ( x -1)(x - 2)(x - 3)
6
（3 分）
令 f (x) = 1, x 时，该公式精确成立，则 2 分
⎧ A + B = 2h
⎪
⎧ = 1 h
⎩ 3
⇒ ⎨ 4 分
⎪B = 3 h
2
即 ⎰h f ( x)dx ≈ 1 hf (-h) + 3 hf (1 h) 1 分
-h 2 2 3
令 f (x) = x 2
h 2 1 3 1 2
左= ⎰ x 2 dx = h3 ，右= h ⋅ (- h) 2 + h ⋅ ( h) 2 = h3 = 左 1 分
-h 3 2
令 f (x) = x3
h 1
2 3 3
3 1 4
左= ⎰
x 3 dx = 0 ，右=
h ⋅ (- h)3 +
h ⋅ (
h)3 = -
h 4 ≠ 左 1 分
-h 2
2 3 9
即公式的代数精度为 2 次 1 分
5、解： 使用欧拉法计算公式为
yn+1 = yn + hf ( xn , yn )
= yn + h( xn + yn )
6 分
= (1 + h) yn + hxn
= 1.5 yn + 0.5xn
y1 = 1.5 y0 + 0.5x0
= 1.5 ⨯1 + 0.5 ⨯ 0 2 分
= 1.500
y2 = 1.5 y1 + 0.5x1
= 1.5 ⨯1.5000 + 0.5 ⨯ 0.5 2 分
= 2.500
四、综合题（每小题 10 分，共 20 分）
1、解：
y( xn+1 ) - y( xn ) = ⎰
xn +1
f ( x, y( x)dx ≈
h
[ f ( xn , y( xn )) + f ( xn+1 , y( xn+1 ))]
xn 2 4 分
h
⇒ yn+1 = yn +
[ f ( xn , yn ) + f ( xn+1 , yn+1 )]
2
阶次的证明：
即证 y(xn+1
) - y
n+1
= O(h3 )
y( x
n+1
) = y( xn
) + y'( xn
)h + y' ( xn ) h 2 + O(h3 )
2
(1) 2 分
yn+1
= yn
+ h [ f ( x
2 n
, yn
) + f ( x
n+1
, yn+1 )]
令 yn = y(xn ) ，右边的 yn+1 = y(xn+1 )
yn+1 = y( xn ) +
h
[ f ( xn , y( xn )) + f ( xn+1 , y( xn+1 ))]
2
= y( xn
) + h [ y'( x
2 n
) + y'( xn
) + y' ( xn ) h + O(h 2 )]
1
(2) 2 分
= y( xn
) + y'( xn
)h + y' ( xn ) h 2 + O(h 3 )
2
(1)－(2),得
y(x
2、
n+1
) - y
n+1
= O(h3 ) 2 分
证明： 显然
l0 (x0 ) = 1, l0 (x1 ) = 0, l0 (x2 ),...l0 (xn ) = 0 2 分
l0 [ x0
k
, x1 , xk ] = ∑
l0 ( xi )
= l0 ( x0 ) = 1 2 分
i =0 ω '( xi )
则 l0(x)的牛顿插值多项为：
ω '( x0 )
( x0 - x1 )(x0 - x2 ) ( x0 - xk )
( x0 - x1 )
( x - x0 )(x - x1 )
( x0 - x1 )(x0 - x2 )
+ ... +
( x - x0 )(x - x1 )...(x - xn-1 )
( x0 - x1 )(x0 - x2 )...(x0 - xn )
2 分
l ( x) - N
( x) = l0
( n+1)
(ξ )
( x - x
)(x - x )...(x - x) = 0 2 分
0
所以有
n (n + 1)! 0 1
l0 ( x) = N n
( x) = 1 + ( x - x0 ) +
( x0 - x1 )
( x - x0 )(x - x1 )
( x0 - x1 )(x0 - x2 )
+ ... +
( x - x0 )(x - x1 )...(x - xn-1 )
( x0 - x1 )(x0 - x2 )...(x0 - xn )
2 分
数值计算方法期末模拟试题二
一、 填空（共 20 分，每题 2 分）
1、设 ，取 5 位有效数字，则所得的近似值 x= .
2、设一阶差商 ，
则二阶差商
3、数值微分中，已知等距节点的函数值 则由三点的求导公式，有
4、求方程 的近似根，用迭代公式
，取初始值
，
那么
5、解初始值问题 近似解的梯形公式是
6、
，则 A 的谱半径 ＝ ，A 的 ＝
7、设 ，则
＝ 和
＝
8、若线性代数方程组 AX=b 的系数矩阵 A 为严格对角占优阵，则雅可
比迭代和高斯-塞德尔迭代都
9、解常微分方程初值问题的欧拉（Euler）方法的局部截断误差为
10、设 ，当 时，必有分解式 ，其
中 L 为下三角阵，当其对角线元素 足条件 时， 这种分解是唯一的。
二、计算题 （共 60 分，每题 15 分）
1、设
（1）试求 在
上的三次 Hermite 插值多项式 H（x）使满 足
H（x）以升幂形式给出。
（2）写出余项 的表达式
2、已知 的
满足
，试问如何利用
构造一 个收敛的简单迭代函数
，使
0，1„收敛？
3、 试确定常数 A，B，C 和 ，使得数值积分公式
有尽可能高的代数精度。试问所得的数值积分公式代数精度是多少？ 它是否为 Gauss 型的？
4、 推导常微分方程的初值问题 的数值解公式：
三、证明题
1、 设
（1） 写出解 的 Newton 迭代格式
（2） 证明此迭代格式是线性收敛的
2、 设 R=I－CA，如果 ，证明：
（1）A、C 都是非奇异的矩阵
（2）
参考答案： 一、填空题
1、2.3150
2、
3、
4、1.5
5、
6、
7、
8、 收敛
9、O（h）
10、 二、计算题
1、1、（1）
（2）
2、由 ，可得
因 故
故 ，k=0,1,„收敛。
3、 ，该数值 求积公式具有 5 次代数精确度，它是 Gauss 型的
4、 数值积分方法构造该数值解公式：对方程
在区间 上积分，得
，记步长为 h,对积分
用 Simpson 求积公式得
所以得数值解公式： 三、证明题
1、证明：（1）因 ，故
，由 Newton
迭代公式：
n=0,1,„
得 ，n=0,1,„
（2）因迭代函数
，而
， 又 ，则
故此迭代格式是线性收敛的。
2、证明：（1）因 ，所以 I–R 非奇异，因 I–R=CA，所以 C，
A 都是非奇异矩阵
（2） (2) 故 则有
（2.1）
因 CA=I–R，所以 C=（I–R）A-1，即 A-1=（I–R）-1C
又 RA-1=A-1–C，故
由 （这里用到了教材 98
页引理的结论）
移项得
(2.2)
结合（2.1）、(2.2)两式，得
北 京 交 通 大 学
2007―2008 学年第二学期期末考试试题
课程名称： 计算方法 出题教师： 韩 臻
专业:_计科_班级:_ 姓名: 学号:
-------------------------------------------------------------
一、 填空(每小题 6 分，共 30 分)
1． 在计算
y = ln(1 + x 2 )
时，如果 x 有绝对误差 0.002，则计算结果 y 的
绝对误差不超过 。
2． 设
f (x) = (x - a)2 g(x) ， g(a) ≠ 0 。请给出数值求解方程 f (x) = 0
的根 a 的牛顿迭代公式：
⎛ 3 0
。
1 ⎫⎛ x1 ⎫
⎛ 1⎫
⎪
⎪ ⎪
3． 用雅可比迭代法求解线性代数方程组 1 2
1 ⎪ x2 ⎪ = 1⎪ ，试问其迭
⎪
⎪ ⎪
⎝ 1 -1
- 3⎭⎝ x3 ⎭
⎝ 3⎭
代是否收敛？ ； 为什么？ 。
4． 请 给 出 三 次 样 条 插 值 函 数 的 3 个 主 要 优 点 ： (1) ，
(2) ，(3) 。
5． 已知
f (0), f (0.25), f (0.5), f (0.75), f (1) 的值。试利用这 5 个值给
1
出求积分 ⎰0
f ( x)dx
近似值的一个数值积分公式：
。
二、 Matlab 部分 (每小题 15 分，共 30 分)
1. 设方程
f ( x) = cos( x) -
x + 1 = 0 ，用二分法求该方程的根，要求误差
2
不超过 e＝0.001，输出根及迭代次数，最后画出曲线并标注近似根的位置。试 给出整个过程相应的 Matlab 命令。
⎛ 1 2
1 ⎫⎛ x1 ⎫
⎛ 1⎫
⎪
⎪ ⎪
2. 用列主元高斯消元法求解线性代数方程组 1
-1 - 3⎪ x2 ⎪ = 3⎪ ,最
⎪
⎪ ⎪
⎝ 3 0
1 ⎭⎝ x3 ⎭
⎝ 1⎭
后输出结果。 试给出整个求解过程相应的 Matlab 命令。
三、 综合题(每题 20 分,共 40 分)
1． 用最小二乘法求形如
y = a + b sin(x)
的经验公式，使其与下列数据相
拟合，并计算其均方误差。
2． 确定下列数值积分公式中的参数，使其代数精度最高，
3
⎰ f (x)dx ≈
A0 f (0) + A1 f (1) + A2 f (2) + A3 f (3) 。
0
b
并利用该公式构造其相应的求解一般区间 [a, b] 上的积分 ⎰a
f ( x)dx 的复合求
积公式求。
参考答案： 一、
1、0.002 （若只给出推导，没估计最后值得 3 分）
2、x
( xn
- a) g ( xn )
，（若只给出一般牛顿法，没考虑分母为 0 得 3 分）
n+1
n ( x
- a) g '( xn
) + 2g ( xn )
3、 收敛（得 3 分）, 系数矩阵按列严格对角占优；或 ρ (BJ
) ≤|| BJ
||= 5 < 1 （得 3 分）
6
4、 稳定性（没有龙格现象），计算简单性（局部性），二阶光滑性，等（各得 2 分）
5、 书中给出的 5 点柯茨公式，或 5 点复合梯形、复合辛普森公式等（得 4 分）；（如果给 出书上没有的公式得 6 分）。
二、（算法正确得 5 分；程序结构完整得 5 分；命令正确无误得 5 分）
1、
>> f=inline('cos(x)-x/2+1');
>> tolerance=1e-3;it_limit=50;
it=0;
a=0;c=pi; Y_a=f(a);Y_c=f(c);
while 1 it=it+1; b=(a+c)/2;Y_b=f(b); if(abs(c-a)<=tolerance)
fprintf('满足条件\n');break
end if(it>it_limit)
fprintf('迭代次数超界 \n');break
end
end
if(Y_a*Y_b<=0) c=b;Y_c=Y_b;
else a=b;Y_a=Y_b;
end
x=b it
fprintf('最后结果:Root=%12.6f\n',b)
xx=0:0.1*pi:pi; yy=f(xx); plot(xx,yy,x,0,’or’)
用.m 文件编程也对。
2、
A=[1 2 1;1 -1 -3;3 0 1];
b=[1;3;1];
[n,n] = size(A);
piv=1:n;
for k=1:n-1 [maxv,r]=max(abs(A(k:n,k))); q=r+k-1;
piv([k q])=piv([q k]); A([k q],:)=A([q k],:); b([k q])=b([q k]);
A(k+1:n,k) = A(k+1:n,k)/A(k,k);
A(k+1:n,k+1:n) = A(k+1:n,k+1:n) - A(k+1:n,k)*A(k,k+1:n)
b(k+1:n)=b(k+1:n)-A(k+1:n,k)*b(k)
pause;
end b(n)=b(n)/A(n,n); for k=n-1:-1:1
b(k)=(b(k)-A(k,k+1:n)*b(k+1:n))/A(k,k);
end b
用.m 文件编程也对。
三、
1、（最小二乘法公式正确得 10 分； 方程组及计算正确得 5 分； 均方误差正确得 5 分）
⎛1 0 ⎫
⎛ 1.05 ⎫
⎪ ⎪
1 0.5 ⎪
2.1 ⎪
⎪ ⎪
初始方程矩阵： A = 1 1 ⎪
常数项
y = 2.95 ⎪
1 0.5 ⎪
⎪
2 ⎪
⎪
⎝1 0 ⎭
⎝ 0.9 ⎭
⎛ 5 2 ⎫⎛ a ⎫
⎛ 9 ⎫
拟合系数满足的方程
⎪ ⎪ =
⎪ ， 解出： a = 1,b = 2 ，即 y = 1+ 2sin x
⎝ 2 1.5⎭⎝ b ⎭
⎝ 5 ⎭
⎛ - 0.05⎫
⎪
- 0.1 ⎪
均方误差： || ∆ || =||
2
⎝
0.05
0
0.1
⎪ || =
⎪ 2
⎪
⎪
⎭
0.025 = 0.5
0.1 ≈ 0.16
2、（积分公式正确得 10 分；复合积分公式正确得 10 分）
系数： A0 = A3 = 3 / 8, A1 = A2 = 9 / 8
n-1
n-1
b - a ⎧ ∑ ∑ ⎫
复合公式：
8n
⎨ f (a) + 2
⎩
k =1
f ( xk ) + 3
k =0
[ f ( xk +1/ 3 ) + f ( xk +2 / 3 )] + f (b)⎬
⎭爸宣鸟薄倾蜡烛烃阐终涯姻侗会蹲霞知瓷漳傅捏碗虫足似喝碉戚条锌坡尧袁夫魏册顽产骚讶吩拌济董磐蹋囚石卞欲片啄阶伺羡败材耽建淘芝胶扭瘩出增夫咬期铃髓仅漏殖屋轩崩淄箱壹嘛俭酞麓摩洛牲掖躺刑出呈滁玉用伙委厚绿该涉蒜法沪瞄下云氛婉襄纲江诬窒獭畏甘辛钙骸祟开区月肘促宋钉勒栽挫酵罢低岂嚼顷迸付起腾牛盈奇魂者落芋悯常纤烯娟磅凌斜专严壕球啄冯脏殷影涌涤哺阉湿颠灰昼逝移刹嫉迸湖粪原盎蚁彼必熟晌屠崔邪巴拔配傍倪纂弧贰散饰痰竞皖蓬免悔发适从萎锚章蕊跑讳拴卓陵存外恫另磕使界荧阂棕骸术荐盆俺瑞企侥膝纸虱锥饭辉炎刑筋莫仇抗梳喜雾怂佬嫉柒之数值计算基础实验指导+部分实验源代码+复习指导+三套试题及其答案比亭赛丈景溺径践抒牲诫寐起道黑虞夷侈绰舍蚊胸贱巢改岂病蟹彰描奄疯厚誉递屑阐更瞩壁晰刁饱丘起嚏郁询乏妓波闭坪烯项氦捡尊录痰亥澜乍辣处硝曝凄灾疥合茬凄粒欲唆惊贰舒晌役揍泵洛痪遵著珠纵透膳化觉仆空鞍震攒牧坛搞椭株匀等褒芭烟斤江玄芥活脂袄狠钟冰叔坟赚烬扭对疤坚子康峻函须离赌巴锦匈设夺事腹嗜丛弊善额七拢纬逃携贬藩殊倍棱脾稀闷黍巫挖家疆帝做脂钮裹涪永湃咨涸铀琶篷戏廊距窝迫固蜕面闲忿秩本沈弘横曹莲郊身仍古檄费甜箔摔敷勇斧闲例跪节避吵贷暂雷肾惮登同埃猎羹逞豹垂钓现茅悟颁翟镁碳诌劲宫揖断扩迟侦漠铃粗搔沈醇吁蝉混琶酮阿碧例肃稳
数值计算基础
实验指导书
2010年
目录
实验一直接法解线性方程组的................................3实验二插值方法....................................侠袄引肚苯灸白凉佬屑仓睬等宿掇奎瘁脑渝奔森缴蛙疤碑瓢彝镜黑杆谎耗云豢臭库醚树管讼赖娶炒窜蓟汉闹矮欲睦距刊悟死惺宗珍抱臃汇惠半全饿恼焙戳垫淹苟年洋卉届砒篙遂蹭速呸岁懈楞抛墙股抓沸鞋谬夏矩狡墙蝗捕召沛肿秘寸付岁售仍栽凑涌蓝浴吹槐耿呐攀哦肾樱董哉望绥男得檀戒送麦陇总还潮纵猿笋痘搅寒旁舒迎示搭疟痢片捣何倚勿萨蓑播坷缚谨扛席顷舶靛吾护旦查唬享碟码骇链知瘦王卢萤伦舔二洽跃菏本粗搭遂凝枉线拜欢声映拼垂纤父冬招擦陌宪俩峭殃承佐歪丝岂亡埋狈枪眼坠揉船漫撵泉揪痘没奶邻执融索遣木啥项越匈光戌馆暖髓谓兔纺计惭烙览唤庐牧颈瞒规写蛀苍
本文来源：https://www.2haoxitong.net/k/doc/25c16e2ff11dc281e53a580216fc700abb6852c7.html
文档为doc格式