《北航数值分析大作业3》.doc

一、算法设计方案 1.使用牛顿迭代法，对原题中给出的，，()的11*21组分别求出原题中方程组的一组解，于是得到一组和对应的。 2.对于已求出的，使用分片二次代数插值法对原题中关于的数表进行插值得到。于是产生了z=f(x,y)的11*21个数值解。 3.从k=1开始逐渐增大k的值，并使用最小二乘法曲面拟合法对z=f(x,y)进行拟合，得到每次的。当时结束计算，输出拟合结果。 4.计算的值并输出结果，以观察逼近的效果。其中。 二、算法实现方案 1、求： （1）Newton法解非线性方程组 ， 其中，t, u, v ,w为待求的未知量，x, y为代入的已知量。 设，给定精度水平和最大迭代次数M，则解该线性方程组的迭代格式为： 迭代终止条件为，若时仍未达到迭代精度，则迭代计算失败。 其中，雅可比矩阵 ， ， （2）分片双二次插值： 根据题目给出的表格， （其中，） 对于给定的,如果满足 则选择为插值节点，相应的插值多项式为 其中， 如果,则在式（2）中取i=1或i=4; 如果,则在式（2）中取u=1或u=4。 在区域上，将(i=0,1,…,10;0.5+0.2j,j=0,1,…20)代入到非线性方程组（1）中，用Newton法解出，再由分片双二次插值得，则有（i=0,1,…,10;j=0,1,…,20），即求出了。 2、求： 乘积型最小二乘拟合曲面： （1）求系数矩阵C： 其中， 计算中涉及到对矩阵求逆，接着在后面将会具体说明列主元的高斯消去法求矩阵的逆的方法。 （2）确定最小的k值，拟合曲面： 设，给定精度水平和最大迭代次数N，则确定最小k值的迭代格式为： 迭代终止条件为，若时仍未达到迭代精度，则迭代计算失败。 待确定满足精度条件的最小k值后，就可以进行曲面拟合计算了。 3、关于列主元的高斯消去法求矩阵的逆： 设非奇异矩阵，且 ，则 ， 对B和I列分块，有 即 其中， 应用列主元的高斯消去法线性解方程组 ， 则即为A的逆。 注：若A不可逆，则此算法失效。 三、源程序 #include "stdafx.h" #include<iostream> #include <stdio.h> #include <cmath> const int M= 500;//迭代最大次数 const double E=1.0e-12;//确定牛顿迭代精度水平 const double E1=1.0e-7;//确定拟合精度水平 const int kmax=9;//k的最大值 const double matrix[6][6]={ {-0.5, -0.34, 0.14, 0.94, 2.06, 3.5}, {-0.42, -0.5, -0.26, 0.3, 1.18, 2.38}, {-0.18, -0.5, -0.5, -0.18, 0.46, 1.42}, {0.22, -0.34, -0.58, -0.5, -0.1, 0.62}, {0.78, -0.02, -0.5, -0.66, -0.5, -0.02}, {1.5, 0.46, -0.26, -0.66, -0.74, -0.5} }; const double mat_t[6]={0, 0.2, 0.4, 0.6, 0.8, 1.0}; const double mat_u[6]={0, 0.4, 0.8, 1.2, 1.6, 2.0}; double U[11][21];//对f(x,y)的值进行存储 double C[kmax+1][kmax+1];//拟合系数矩阵 ////////////两数绝对值取大///////////// double max(double x,double y){ return fabs(x)>fabs(y)?fabs(x):fabs(y); } /////////高斯消元法解线性方程组/////// void linear_solution(double f[4],double ff[4][4], double (&delta)[4]) { int i, j, k, ik; double tmp, mik; for(k=0; k<3; k++) { ik=k; for(i=k;i<4;i++) if(fabs(ff[ik][k])<fabs(ff[i][k])) ik=i; for(j=k;j<4;j++) { tmp=ff[k][j]; ff[k][j]=ff[ik][j]; ff[ik][j]=tmp; } tmp=f[k]; f[k]=f[ik]; f[ik]=tmp; for(i=k+1;i<4;i++) { mik