北航数值分析大作业(三)Word文档格式.docx

北航数值分析大作业(三)Word文档格式.docx

《北航数值分析大作业(三)Word文档格式.docx》由会员分享，可在线阅读，更多相关《北航数值分析大作业(三)Word文档格式.docx（25页珍藏版）》请在冰豆网上搜索。

-0.26

0.3

1.18

2.38

-0.18

0.46

1.42

0.6

0.22

-0.58

-0.1

0.62

0.78

-0.02

-0.66

1.0

1.5

-0.74

确定了一个二元函数z=f（x，y）。

1．试用数值方法求出f（x，y）在区域D={（x,y）︱0≤x≤0.8,0.5≤y≤1.5}上的一个近似表达式

要求p（x,y）一最小的k值达到以下的精度

其中xi=0.08i，yj=0.5+0.05j。

2．计算f（xi*,yj*）,p（xi*,yj*）（i=1,2,…,8;

j=1,2,…,5）的值，以观察p（x,y）逼近f（x,y）的效果，其中xi*=0.1i,yj*=0.5+0.2j。

二、算法方案

1．使用C++语言实现，使用牛顿迭代法求解非线性方程组，对，，（）的共计11×

21组分别求出非线性方程组的解，即求出与对应的。

均为11×

21的矩阵。

2．由求出的，使用分片二次代数插值法对题中给出的数表进行插值得到。

即得到的11×

21个数值解。

3．k=0时的多项式拟合必然不符合要求，从k=1开始迭代，使用最小二乘法的曲面拟合法对进行拟合，计算在不符合要求的情况下增大。

当时结束计算，输出结果。

4．由3中得到的系数计算的值，再次使用牛顿迭代法对进行求解得到，再次进行二次插值得到结果，以观察逼近的效果。

其中，，。

三、源程序：

#include<

iostream>

vector>

cmath>

algorithm>

iomanip>

#defineN4//方程组未知个数

#defineM6//z,t,u数表阶数

#defineX_Num11

#defineY_Num21//定义数表大小

#defineEPSL1e-12//定义阶数，精度

#defineEPSL21e-7

usingnamespacestd;

typedefvector<

vector<

double>

>

Mat;

//将二维数组简写为Mat

Equation（Matinput）;

//定义求解非线性方程的函数，同时供Inverse，Zxy函数调用

MatInverse（intrank,Matinput2）;

//定义求解逆矩阵的函数

doubleAccuracy（vector<

X_1,vector<

X_2）;

//定义求解近似向量精度的函数

doubleInterpolation（doubleu_1,doublet_1）;

//定义分片代数二次插值函数

MatCrs（vector<

X,vector<

Y,MatU）;

//最小二乘法求解近似表达式系数

MatZxy（vector<

X1,vector<

Y1）;

//定义非线性方程组，调用Equation，Accuracy和Interpolation完成求解

//所有的output应该调整，是否调整为输出到文件为好

voidoutput（vector<

Final1,vector<

Final2,MatFinal3）;

//定义输出函数,输出矩阵

voidoutput2（MatXi）;

doublevector_u[M]={0,0.4,0.8,1.2,1.6,2};

doublevector_t[M]={0,0.2,0.4,0.6,0.8,1};

doublemat_z[M][M]={

{-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},

};

//定义初始数表z,t,u，此处使用数组，而其它矩阵和向量均使用的为vector以及二维vector

voidmain（）

{

inti,j;

vector<

x,y;

x.resize（X_Num）;

y.resize（Y_Num）;

for（i=0;

i<

X_Num;

i++）

{

x[i]=0.08*i;

}

Y_Num;

y[i]=0.5+0.05*i;

}//定义插值节点

MatZ=Zxy（x,y）;

//求出插值点函数值

output（x,y,Z）;

MatCr=Crs（x,y,Z）;

//求出近似多项式

output2（Cr）;

x.resize（8）;

y.resize（5）;

8;

x[i]=0.1*（i+1）;

for（j=0;

j<

5;

j++）

y[j]=0.5+0.2*（j+1）;

}//新插值节点

MatZ2=Zxy（x,y）;

MatP;

P.resize（8）;

P[i].resize（5）;

intk=Cr.size（）;

//利用上一部的Cr获得P值即可

doubletmp;

intm,n;

for（m=0;

m<

m++）

{

for（n=0;

n<

n++）

{

tmp=0;

for（i=0;

k;

{

for（j=0;

{

tmp=tmp+Cr[i][j]*pow（x[m],i）*pow（y[n],j）;

}

}

P[m][n]=tmp;

}

}//使用近似多项式得到的近似值

for（j=0;

{

cout<

<

setprecision

（2）<

scientific<

x[i]<

"

"

;

y[j]<

插值"

setprecision（12）<

Z2[i][j]<

拟合"

P[i][j]<

endl;

}

system（"

pause"

）;

}

Y1）

inti,j,k;

x,y;

x=X1;

y=Y1;

intX_No=x.size（）;

intY_No=y.size（）;

X_1,X_2;

doublewrong;

X_1.resize（N）;

//过渡用于判断误差

X_2.resize（N）;

//此处调试发现x,y值略有差异

MatA;

//使用牛顿法迭代的带求非线性方程

Matt,u,v,w;

MatZ;

//对应t,u的数表Z

A.resize（N）;

N;

A[i].resize（N+1）;

t.resize（X_No）;

u.resize（X_No）;

v.resize（X_No）;

w.resize（X_No）;

Z.resize（X_No）;

X_No;

t[i].resize（Y_No）;

u[i].resize（Y_No）;

v[i].resize（Y_No）;

w[i].resize（Y_No）;

Z[i].resize（Y_No）;

Y_No;

t[i][j]=u[i][j]=w[i][j]=v[i][j]=1;

}//将待求量赋予初值

A[i][j]=1;

A[2][0]=0.5;

A[3][1]=0.5;

i++）//求解对应x,y的t,u,v,w

j=0;

while（j<

Y_No）

A[0][4]=2.67+x[i]-0.5*cos（t[i][j]）-u[i][j]-v[i][j]-w[i][j];

//此处应做修改

A[1][4]=1.07+y[j]-0.5*sin（u[i][j]）-t[i][j]-v[i][j]-w[i][j];

A[2][4]=3.74+x[i]-cos（v[i][j]）-0.5*t[i][j]-u[i][j]-w[i][j];

A[3][4]=0.79+y[j]-sin（w[i][j]）-0.5*u[i][j]-t[i][j]-v[i][j];

A[0][0]=-0.5*sin（t[i][j]）;

A[1][1]=0.5*cos（u[i][j]）;

A[2][2]=-sin（v[i][j]）;

A[3][3]=cos（w[i][j]）;

vector<

Change=Equation（A）;

//调用求解方程得到第一组增量,此处需要注意Change赋值的问题，得到每组增量后怎么处理

for（k=0;

k<

k++）

X_1[k]=Change[k];

X_2[0]=t[i][j];

X_2[1]=u[i][j];

X_2[2]=v[i][j];

X_2[3]=w[i][j];

wrong=Accuracy（X_1,X_2）;