北航数值分析大作业(三).docx
《北航数值分析大作业(三).docx》由会员分享，可在线阅读，更多相关《北航数值分析大作业(三).docx（25页珍藏版）》请在冰豆网上搜索。
《数值分析A》
计算实习题目三
姓名:
学号:
学院:
2014年12月
一、题目
关于x，y，t，u，v，w的下列方程组
0.5cost+u+v+w-x=2.67
t+0.5sinu+v+w-y=1.07
0.5t+u+cosv+w-x=3.74
t+0.5u+v+sinw-y=0.79
以及关于z，t，u的下列二维数表
zu
t
0
0.4
0.8
1.2
1.6
2
0
-0.5
-0.34
0.14
0.94
2.06
3.5
0.2
-0.42
-0.5
-0.26
0.3
1.18
2.38
0.4
-0.18
-0.5
-0.5
-0.18
0.46
1.42
0.6
0.22
-0.34
-0.58
-0.5
-0.1
0.62
0.8
0.78
-0.02
-0.5
-0.66
-0.5
-0.02
1.0
1.5
0.46
-0.26
-0.66
-0.74
-0.5
确定了一个二元函数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
#include
#include
#include
#include
#defineN4//方程组未知个数
#defineM6//z,t,u数表阶数
#defineX_Num11
#defineY_Num21//定义数表大小
#defineEPSL1e-12//定义阶数，精度
#defineEPSL21e-7
usingnamespacestd;
typedefvector>Mat;//将二维数组简写为Mat
vectorEquation（Matinput）;//定义求解非线性方程的函数，同时供Inverse，Zxy函数调用
MatInverse（intrank,Matinput2）;//定义求解逆矩阵的函数
doubleAccuracy（vectorX_1,vectorX_2）;//定义求解近似向量精度的函数
doubleInterpolation（doubleu_1,doublet_1）;//定义分片代数二次插值函数
MatCrs（vectorX,vectorY,MatU）;//最小二乘法求解近似表达式系数
MatZxy（vectorX1,vectorY1）;//定义非线性方程组，调用Equation，Accuracy和Interpolation完成求解
//所有的output应该调整，是否调整为输出到文件为好
voidoutput（vectorFinal1,vectorFinal2,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;
vectorx,y;
x.resize（X_Num）;
y.resize（Y_Num）;
for（i=0;i {
x[i]=0.08*i;
}
for（i=0;i {
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）;
for（i=0;i<8;i++）
{
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）;
for（i=0;i<8;i++）
{
P[i].resize（5）;
}
intk=Cr.size（）;//利用上一部的Cr获得P值即可
doubletmp;
intm,n;
for（m=0;m<8;m++）
{
for（n=0;n<5;n++）
{
tmp=0;
for（i=0;i {
for（j=0;j {
tmp=tmp+Cr[i][j]*pow（x[m],i）*pow（y[n],j）;
}
}
P[m][n]=tmp;
}
}//使用近似多项式得到的近似值
for（i=0;i<8;i++）
{
for（j=0;j<5;j++）
{
cout<（2）< cout<（2）< cout<<"插值"< cout<<"拟合"< }
}
system（"pause"）;
}
MatZxy（vectorX1,vectorY1）
{
inti,j,k;
vectorx,y;
x=X1;
y=Y1;
intX_No=x.size（）;
intY_No=y.size（）;
vectorX_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）;
for（i=0;i {
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）;
for（i=0;i {
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）;
}
for（i=0;i {
for（j=0;j {
t[i][j]=u[i][j]=w[i][j]=v[i][j]=1;
}
}//将待求量赋予初值
for（i=0;i {
for（j=0;j {
A[i][j]=1;
}
}
A[2][0]=0.5;
A[3][1]=0.5;
for（i=0;i {
j=0;
while（j {
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]）;
vectorChange=Equation（A）;//调用求解方程得到第一组增量,此处需要注意Change赋值的问题，得到每组增量后怎么处理
for（k=0;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];
wro
升级会员 