计算实习报告3.docVIP

数值分析〔B〕大作业〔三〕

姓名：学号：：

1、算法设计：

①整体思路

这一次的作业和前面两次相比，感觉要难很多，最明显的就是书上没有完整的算法实现流程。再加上插值和逼近这一段在课上也没听太懂什么意思，所以刚开始思考这个程序时感觉很吃力。不过经过后来好好看书和向同学请教，渐渐对插值和分片元二次拟合有点懂了，也大概知道了解决本次问题分别要完成的子问题。

要完本钱次任务，主体上有三局部问题要解决：

〔1〕解非线性方程组。将D内的)当作量代入题目给定的非线性方程组，得到与相对应的t[i][j],u[i][j]。

〔2〕分片双二次插值。采用分片双二次插值，得到与t[11][21],u[11][21]]对应的z[11][21]，也即建立与的对应关系，得到二元函数。

〔3〕曲面拟合。利用x[i],y[j],z[11][21]建立二维函数表，再根据精度的要求选择适当k值，并得到曲面拟合的系数矩阵C[k][k]。

〔4〕观察逼近效果。观察逼近效果只需要重复上面〔1〕和〔2〕的过程，得到与新的插值节点对应的，再与对应的比拟即可。

程序示意图如下所示：

②各局部问题的算法实现

〔1〕解非线性方程组

非线性方程组常用数值解法有简单迭代法和牛顿法。牛顿法收敛快，一般都能到达平方收敛，因此这里选择Newton法解非线性方程组。

牛顿法解方程组的解，可采用如下算法：

1〕在附近选取，给定精度水平和最大迭代次数M。

2〕对于执行

=1\*GB3①计算和。

=2\*GB3②求解关于的线性方程组

=3\*GB3③假设，那么取，并停止计算；否那么转=4\*GB3④。

=4\*GB3④计算。

=5\*GB3⑤假设，那么继续，否那么，输出M次迭代不成功的信息，并停止计算。

注：第③步中的求取，用到了解线性方程组，使用的是列主元素Gauss消去法。

另外，本次作业中解非线性方程组实际求取的是解向量和，方程组中的、是值，为当前节点的值。

〔2〕分片双二次插值

给定数表以及需要插值的节点，进行分片二次插值的算法：

设数表中的点为：，需要插值的节点为。

1)根据选择插值节点：

假设或，插值节点对应取或，

假设或，插值节点对应取或。

假设

那么选择为插值节点。

2〕计算

插值多项式的公式为：

注：在〔1〕中通过解非线性方程组得到的解向量和与插值节点有着一一对应关系，而此题题目给出的函数表为关于的函数关系。因此，为得到关于的函数关系，此题中应先根据给定数表对进行插值，再利用与的一一对应关系得到与的对应关系。

〔3〕曲面拟合

根据插值得到的数表进行曲面拟合的过程：

根据拟合节点和基底函数写出矩阵B和G：

计算。

在这里，为了简化计算和编程、防止矩阵求逆，记：

，

对上面两式进行变形，得到如下两个线性方程组：

，

通过解上述两个线性方程组，那么有：

对于每一个，。

拟合需要到达的精度条件为：

。

其中对应着插值得到的数表中的值。

让k逐步增加，每一次重复执行以上几步，直到

成立。此时的k值就是所需最小的k。

注：曲面拟合过程中用到的数表为前面插值得到的数表，即。

在2)中，求矩阵A和D的过程用的是列主元素Gauss消去法。

〔4〕观察逼近效果。

观察逼近效果主要要做的是，通过新的插值节点建立新的插值数表，同时求出对应的，比拟即可。

2、程序源代码：

#includestdio.h

#includeconio.h

#includemath.h

#defineEpsilon11e-12//解线性方程组时近似解向量的精度

#defineMaxM200//解线性方程组时的最大迭代次数

#defineK10//预估的能到达精度要求的k值，K=k

#definegivenX11//要求的插值节点xi个数

#definegivenY21//要求的插值节点yi个数

#definetestX8//观察逼近效果时，插值节点xi个数

#definetestY5//观察逼近效果时，插值节点yi个数

/***************************************************************/

/*函数功能描述：