北航数值分析大作业题目三.docVIP

《数值分析》第三次大作业

算法的设计方案：

（一）、总体方案设计：

（1）解非线性方程组。将给定的当作已知量代入题目给定的非线性方程组，求得与相对应的数组t[i][j],u[i][j]。

（2）分片二次代数插值。通过分片二次代数插值运算，得到与数组t[11][21],u[11][21]]对应的数组z[11][21]，得到二元函数z=。

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

（4）观察和的逼近效果。观察逼近效果只需要重复上面（1）和（2）的过程，得到与新的插值节点对应的，再与对应的比较即可，这里求解可以直接使用（3）中的C[r][s]和k。

（二）具体算法设计：

（1）解非线性方程组

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

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

2）对于执行

=1\*GB3①计算和。

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

=3\*GB3③若，则取，并停止计算；否则转=4\*GB3④。

=4\*GB3④计算。

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

（2）分片双二次插值

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

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

1)根据选择插值节点：

若或，插值节点对应取或，

若或，插值节点对应取或。

若

则选择为插值节点。

2）计算

插值多项式的公式为：

注：本步进行插值运算的是，利用与的对应关系就可以得到与的对应关系。

（3）曲面拟合

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

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

计算。

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

，

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

，

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

对于每一个，。

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

。

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

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

成立。此时的k值就是要求解最小的k。

源程序：

#includestdio.h

#includeiostream

#includestdlib.h

#includemath.h

#includefloat.h

#includeiomanip

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

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

#defineN10/*求解迭代次数时假设的k的最大值，用于定义包含k的存储空间*/

voidNewton();/*牛顿法求解非线性方程组子程序*/

voidfpeccz();/*分片二次代数插值子程序*/

voidqmnh();/*曲面拟合子程序*/

voidduibi();/*对比??和p逼近效果的子程序*/

doublex[11],y[21],t[11][21],u[11][21];/*定义全局变量*/

doublez[11][21],C[10][10];

doublekz;

voidNewton(doublex[11],doubley[21])/*牛顿法求解非线性方程组子程序*/

{

doubleX[4],dx[4],F[4],dF[4][4],temp,m,fx,fX;

inti,j,k,l,p,ik,n;

for(i=0;i=10;i++)

{

for(j=0;j=20;j++)

{

X[0]=1;/*选取迭代初始向量，四个分别代表t，u，v，w*/

X[1]=1;

X[2]=1;

X[3]=1;

n=0;

loop1:{F[0]=0.5*cos(X[0])+X[1]+X[2]+X[3]-x[i]-2.67;

F[1]=X[0]+0.5*sin(X[1])+X[2]+X[3]-y[j]-1.07;

F[2]=0.5*X[0]+X[1]+cos(X[2])+X[3]-x[i]-3.74;

F[3]=X[0]+0.5*X[1]+X[2]+sin(X[3])-y[j]-0