首页 北航数值分析上机作业第3题

北航数值分析上机作业第3题

举报
开通vip

北航数值分析上机作业第3题北航数值分析上机作业第3题 一 算法设计 1 设计方案(流程图) 在编程学习过程中,可知道本实习题目的训练目的在于对解非线性方程组、函数插值以及曲面拟合的逼近(第2部分介绍)的进一步认识。本实习题目的算法设计方案如下: 开始 i=0,10,j=0,20,对xi,yj赋值 Newton法解非线性方程组,得到tuvw 插值得到z=f(xi,yj),输出结果 曲面拟合,得到拟合系数矩阵C,输出C i=1,8,j=1,5,对x*i,y*j赋值 计算p(x*i,y*j) Newton法解非线性方程组,得到...

北航数值分析上机作业第3题
北航数值 分析 定性数据统计分析pdf销售业绩分析模板建筑结构震害分析销售进度分析表京东商城竞争战略分析 上机作业第3题 一 算法设计 1 设计 方案 气瓶 现场处置方案 .pdf气瓶 现场处置方案 .doc见习基地管理方案.doc关于群访事件的化解方案建筑工地扬尘治理专项方案下载 ( 流程 快递问题件怎么处理流程河南自建厂房流程下载关于规范招聘需求审批流程制作流程表下载邮件下载流程设计 图) 在编程学习过程中,可知道本实习题目的训练目的在于对解非线性方程组、函数插值以及曲面拟合的逼近(第2部分介绍)的进一步认识。本实习题目的算法设计方案如下: 开始 i=0,10,j=0,20,对xi,yj赋值 Newton法解非线性方程组,得到tuvw 插值得到z=f(xi,yj),输出结果 曲面拟合,得到拟合系数矩阵C,输出C i=1,8,j=1,5,对x*i,y*j赋值 计算p(x*i,y*j) Newton法解非线性方程组,得到tuvw 插值得到z=f(x*i,y*j) 输出f(x*i,y*j),p(x*i,y*j) 结束 2 部分算法说明 本部分主要列举了Newton法求解非线性方程组、双二次插值和曲面拟合三种算法,由于在之前的实习题中已对列主元Gauss消去法的描述分析,在此不再赘述。算法说明的文字部分来自教材。 Newton法解方程组 (0)*,(1) 在x附近选取xD,,给定精度水平和最大迭代次数M (2) 对于k=0,1,…,M执行 ()k'()kFx()Fx()? 计算和 '()()()kkk()kFxxFx()(),,,? 求解关于,x的线性方程组 ()()kk*()k,,xx/,? 若,则取xx,,并停止计算;否则转? (1)()()kkk,? 计算 xxx,,, ? 若,则继续;否则,输出M次迭代不成功的信息,停止计算。 kM, 双二次插值 设 xxih,, (0,1,...,)in,i0 yyj,,, (0,1,...,)jm,i0 若(x,y)满足 xhxxh,,,,/2/2 22,,,inii yyy,,,,,,/2/2 22,,,jmjj (,)xy则选择为插值节点,相应的插值多项式为 (1,,1;1,,1)kiiirjjj,,,,,,kr j,1i,1 pxylxlyfxy(,)()()(,), ,,22krkrkirj,,,,11其中 i,1xx,t (1,,1)kiii,,, ()lx,,kxx,ti,,1kttk, j,1yy,t (1,,1)rjjj,,, ()ly,,ryy,tj,,1rttr, xxh,,/2xxh,,/2yy,,,/2如果或者,则取i=1或i=n-1;如果或者1n,11 yy,,,/2,则取j=1或j=m-1。 m,1 曲面拟合 (,)xuy(0,1,...,)im,固定,以为基函数组对数据作最小二乘拟合,得到n+1条,()x,,iijir 拟合曲线 M qxx()(),,,(0,1,...,)jn, ,jrjr,0r T其中是法方程 (,,...,),,,,,1ojjMjj TT(0,1,...,)jn, BBBu,,jj T,,A,,的解,而,。记,Bx,,()uuuu,(,,...,),,rjri1jojjmj,,(1)(1),,,mM(1)(1)Mn,,, ,,,则有 Uu,ij,,(1)(1)mn,,, TT,1ABBBU,() (,())yqx固定x,以为基函数组对数据 作最小二乘拟合得到所求(0,1,...,)jn,,()y,,jjs 的拟合曲面 N*pxyy(,)(),,, ,sjs,s0 T((),(),...,())(),,,,xxxx,其中是法方程 01N TTGGxGqx,()(), T,,qxqxqxqx()((),(),...,()),的解,而,。 Gy,,()01nsj,,(1)(1),,,nN TT,1,,设,则有 ()()GGGy,,sjj,,(1)(1)Nn,,, Mn ,,,,()()xx, ,,srjsjrrj,,00 所以 NM**pxycxy(,)()(),,, ,,rsrs,,00sr 其中 n*c,,,(0,1,...,;0,1,...,)rMsN,, ,rsrjsjj0, **,,c记Cc,,则有系数的 关于同志近三年现实表现材料材料类招标技术评分表图表与交易pdf视力表打印pdf用图表说话 pdf 达式 rsrs,,,,,(1)(1)MN TTT,,11CBBBUGGG,()() 二 源程序及功能说明 本题目程序部分采用Fortran语言编写,分为主程序和多个子程序。主程序通过对变量进行初始赋值,调用子程序实现算法目的的计算,并对部分计算结果进行输出。同时,程序中编写了多个子程序实现不同的功能,并在子程序过程中对部分结果进行了输出。由于个人编程的习惯,在程序中建立了一个module存放部分变量和常量。在本部分中,首先列出了源程序,并在相关的部分进行程序说明(由于个人编程的习惯,本人喜欢把每个程序分别放于一个子文件中,因此以下源程序的顺序是按文件名的排列顺序排列);接着在第二部分分别对主程序和每个子程序功能进行说明。由于在题目中给定的变量角标很多是从0开始,而在Fortran语言中默认数组从1开始,因此在程序中的大部分强制定义从0开始,以避免换算过程造成的失误。 1 源程序 function delta(c,k,z) !求对应的插后值精度 implicit none double precision c(0:10,0:10),z(0:230) double precision x(0:10),y(0:20),tt(0:230) !11*21 double precision delta,pxy,del_sum integer m,n,i,j,k m=20 n=10 do i=0,n x(i)=0.08d0*i end do do j=0,m y(j)=0.5d0+0.05d0*j end do do i=0,n do j=0,m tt(i*21+j)=pxy(c,k,x(i),y(j)) !调用函数 end do end do del_sum=0.0d0 do i=0,230 del_sum=del_sum+(tt(i)-z(i))*(tt(i)-z(i)) end do delta=del_sum end function delta function pxy(c,k,x,y) !拟合近似表达式 implicit none double precision c(0:10,0:10) double precision pxy,fin,x,y integer i,j,k fin=0.0d0 do i=0,k do j=0,k fin=fin+c(i,j)*(x**i)*(y**j) end do end do pxy=fin end function pxy function interpo(t,u) !二元插值函数 P101 implicit none double precision interpo double precision x(0:5),y(0:5),f(0:5,0:5) double precision lk(3),lr(3),sum_interp(3) double precision t,u,h,tao,sum integer n,m,i,j,it,jt data f /-0.5d0,-0.42d0,-0.18d0,0.22d0,0.78d0,1.5d0,& !第1列 -0.34d0,-0.5d0,-0.5d0,-0.34d0,-0.02d0,0.46d0,& !第2列 0.14d0,-0.26d0,-0.5d0,-0.58d0,-0.5d0,-0.26d0,& !第3列 0.94d0,0.3d0,-0.18d0,-0.5d0,-0.66d0,-0.66d0,& !第4列 2.06d0,1.18d0,0.46d0,-0.1d0,-0.5d0,-0.74d0,& !第5列 3.5d0,2.38d0,1.42d0,0.62d0,-0.02d0,-0.5d0/ !第6列 h=0.2d0 !定义步长 tao=0.4d0 n=5 m=5 do i=0,5 !xi,yi x(i)=0.2d0*i y(i)=0.4d0*i end do if (t<=(x(1)+h/2.0)) then !t方向位置判断 it=1 end if do i=2,n-2 if ((t>x(i)-h/2.0).and.(t<=x(i)+h/2.0)) then it=i exit end if end do if (t>x(n-1)-h/2.0) then it=n-1 end if if (u<=y(1)+tao/2.0) then !u方向位置判断 jt=1 end if do j=2,m-2 if ((u>y(j)-tao/2.0).and.(u<=y(j)+tao/2.0)) then jt=j exit end if end do if (u>y(m-1)-tao/2.0) then jt=m-1 end if lk=0.d0 lr=0.d0 lk(1)=(t-x(it))*(t-x(it+1))/(x(it-1)-x(it))/(x(it-1)-x(it+1)) lk(2)=(t-x(it-1))*(t-x(it+1))/(x(it)-x(it-1))/(x(it)-x(it+1)) lk(3)=(t-x(it-1))*(t-x(it))/(x(it+1)-x(it-1))/(x(it+1)-x(it)) lr(1)=(u-y(jt))*(u-y(jt+1))/(y(jt-1)-y(jt))/(y(jt-1)-y(jt+1)) lr(2)=(u-y(jt-1))*(u-y(jt+1))/(y(jt)-y(jt-1))/(y(jt)-y(jt+1)) lr(3)=(u-y(jt-1))*(u-y(jt))/(y(jt+1)-y(jt-1))/(y(jt+1)-y(jt)) sum=0.0d0 sum_interp=0.d0 sum_interp(1)=lk(1)*lr(1)*f(it-1,jt-1)+lk(1)*lr(2)*f(it-1,jt)+lk(1)*lr(3)*f(it-1,jt+1) sum_interp(2)=lk(2)*lr(1)*f(it,jt-1)+lk(2)*lr(2)*f(it,jt)+lk(2)*lr(3)*f(it,jt+1) sum_interp(3)=lk(3)*lr(1)*f(it+1,jt-1)+lk(3)*lr(2)*f(it+1,jt)+lk(3)*lr(3)*f(it+1,jt+1) !共3*3项 do i=1,3 sum=sum+sum_interp(i) end do interpo=sum sum=0.0d0 lk=0.d0 lr=0.d0 end function interpo function max_inf(Delta_val,n) !定义Delta绝对值最大值函数(范数) implicit none integer i,n double precision max_inf double precision Delta_val(n) max_inf=0.d0 do i=1,n if (abs(Delta_val(i))>abs(max_inf)) then max_inf=Delta_val(i) end if end do end function max_inf module param !建立模,存放变量 implicit none double precision equ_x(0:10),equ_y(0:20) !定义方程组的变量 double precision tuvw(4) double precision xx(0:7),yy(0:4) !第2问的变量 integer,parameter::maxits=10000 !定义最大迭代次数 double precision,parameter::epsino=1.0d-12 end module param subroutine inverse(a,k,b) !子程序 求逆 implicit none double precision a(0:10,0:10),b(0:10,0:10),a2(0:11,0:11),e(0:11,0:11),e2(0:11),x(0:11) integer k,k2,k3,n,i,j k2=k n=10 do i=0,n !矩阵赋值,保留原始矩阵 do j=0,n a2(i+1,j+1)=a(i,j) end do end do e=0.d0 do i=0,n !单位对角阵 do j=0,n if (i==j) then e(i+1,j+1)=1.0d0 else e(i+1,j+1)=0.d0 end if end do end do do k3=0,k2 do i=0,n e2(i+1)=e(i+1,k3+1) end do call gauss_cal_inverse(k2,a2,e2,x) !循环调用求解 do i=1,(k2+1) b(i-1,k3)=x(i) end do end do end subroutine inverse subroutine delta2(c,z,tt) !第2问 implicit none double precision c(0:10,0:10),z(0:44),tt(0:44) !由前面的计算,拟合系数矩阵已经存储 double precision xx(0:7),yy(0:4) double precision pxy integer m,n,i,j,k m=8 n=5 k=5 do i=1,m xx(i-1)=0.1d0*i end do do j=1,n yy(j-1)=0.5d0+0.2d0*j end do do i=1,m do j=1,n tt((i-1)*5+j-1)=pxy(c,k,xx(i-1),yy(j-1)) end do end do end subroutine delta2 subroutine equa_cal(x,y,tuvw) !牛顿法求解线性方程组 use param,only:maxits,epsino implicit none double precision x,y,t,u,v,w double precision delta_t,delta_u,delta_v,delta_w double precision A(4,4),Delta(4),B(4),tuvw(4) !A为求导矩阵,Delta(4)为求解的值(排列顺序 double precision max_inf ! t,u,v,w),B(4)为移项后的式,tuvw为所要的值 integer i integer its Delta=0.0d0 t=1.0d0 !选定初始值 u=1.0d0 v=1.0d0 w=1.0d0 tuvw(1)=t !将四个变量放于矩阵中 tuvw(2)=u tuvw(3)=v tuvw(4)=w A(1,1)=-0.5d0*dsin(tuvw(1)) !求导矩阵赋值 A(1,2)=1.0d0 A(1,3)=1.0d0 A(1,4)=1.0d0 A(2,1)=1.0d0 A(2,2)=0.5d0*dcos(tuvw(2)) A(2,3)=1.0d0 A(2,4)=1.0d0 A(3,1)=0.5d0 A(3,2)=1.0d0 A(3,3)=-1.0d0*dsin(tuvw(3)) A(3,4)=1.0d0 A(4,1)=1.0d0 A(4,2)=0.5d0 A(4,3)=1.0d0 A(4,4)=1.0d0*dcos(tuvw(4)) B(1)=-(0.5d0*dcos(tuvw(1))+tuvw(2)+tuvw(3)+tuvw(4)-x-2.67d0) !方程组赋值 B(2)=-(tuvw(1)+0.5d0*dsin(tuvw(2))+tuvw(3)+tuvw(4)-y-1.07d0) B(3)=-(0.5*tuvw(1)+tuvw(2)+dcos(tuvw(3))+tuvw(4)-x-3.74d0) B(4)=-(tuvw(1)+0.5*tuvw(2)+tuvw(3)+dsin(tuvw(4))-y-0.79d0) its=1 call mainrow_cancle(A,B,Delta) !在循环外做一次求解,得到一组Delta值 do while (abs(max_inf(Delta,4))/abs(max_inf(tuvw,4))>epsino.and.its<=maxits) !在规定精度和 its=its+1 ! 迭代次数内迭代 do i=1,4 tuvw(i)=tuvw(i)+Delta(i) end do !----------------------- A=0.d0 B=0.d0 A(1,1)=-0.5d0*dsin(tuvw(1)) !求导矩阵赋值 A(1,2)=1.0d0 A(1,3)=1.0d0 A(1,4)=1.0d0 A(2,1)=1.0d0 A(2,2)=0.5d0*dcos(tuvw(2)) A(2,3)=1.0d0 A(2,4)=1.0d0 A(3,1)=0.5d0 A(3,2)=1.0d0 A(3,3)=-1.0d0*dsin(tuvw(3)) A(3,4)=1.0d0 A(4,1)=1.0d0 A(4,2)=0.5d0 A(4,3)=1.0d0 A(4,4)=1.0d0*dcos(tuvw(4)) B(1)=-(0.5d0*dcos(tuvw(1))+tuvw(2)+tuvw(3)+tuvw(4)-x-2.67d0) !方程组赋值 B(2)=-(tuvw(1)+0.5d0*dsin(tuvw(2))+tuvw(3)+tuvw(4)-y-1.07d0) B(3)=-(0.5*tuvw(1)+tuvw(2)+dcos(tuvw(3))+tuvw(4)-x-3.74d0) B(4)=-(tuvw(1)+0.5*tuvw(2)+tuvw(3)+dsin(tuvw(4))-y-0.79d0) !----------------------- call mainrow_cancle(A,B,Delta) !调用列主元Gauss消去法求解线性方程组 A=0.d0 B=0.d0 if (its==maxits) then print*,'迭代不成功' exit end if end do end subroutine equa_cal subroutine gauss_cal_inverse(kk,a,b,x) !求逆中用到的Gauss implicit none integer kk double precision a(0:11,0:11),b(0:11),x(0:11) !a(0:11,0:11),b(0:11),x(0:11) double precision a2(0:11,0:11),b2(0:11),m(0:11,0:11) double precision mid_max,m1,d1 integer i,j,k,ik,n n=kk+1 do i=1,n do j=1,n !保持原有矩阵型式 a2(i,j)=a(i,j) end do b2(i)=b(i) end do do k=1,n-1 mid_max=a2(k,k) do i=k+1,n if (abs(mid_max)abs(max_mid)) then max_mid=a(ik,k) t=ik end if end do a(t,k)=a(k,k) a(k,k)=max_mid b_mid=b(k) !列选主元 b(k)=b(t) b(t)=b_mid do i=k+1,n max_mid_row(i)=a(k,i) a(k,i)=a(t,i) a(t,i)=max_mid_row(i) end do !----------------- do i=k+1,n m(i)=a(i,k)/a(k,k) do j=k+1,n !j表示列 a(i,j)=a(i,j)-m(i)*a(k,j) end do b(i)=b(i)-m(i)*b(k) end do end do !------return 回代------ x(n)=b(n)/a(n,n) do k=n-1,1,-1 sum=0 do j=k+1,n sum=sum+a(k,j)*x(j) end do x(k)=(b(k)-sum)/a(k,k) end do end subroutine mainrow_cancle subroutine fit(z,c) !曲面拟合 P142 implicit none double precision z(0:230),c(0:10,0:10) double precision x(0:10),y(0:20),u(0:10,0:20), & b(0:10,0:10),b2(0:10,0:10),b3(0:10,0:10), & g(0:20,0:10),g2(0:10,0:10),g3(0:10,0:10), & d2(0:10,0:10),d3(0:10,0:20),d4(0:10,0:10) double precision d1,yi,delta_cal,delta !d1作为中间变量 integer m,n,i,j,k,i2 m=20 n=10 yi=1.0d-7 !精度要求 do i=0,n x(i)=0.08d0*i end do do j=0,m y(j)=0.5d0+0.05d0*j end do do j=0,n !基 b(j,0)=1.0d0 end do do i=1,n do j=0,n b(j,i)=x(j)**i end do end do do j=0,m !基 g(j,0)=1.0d0 end do do i=1,n do j=0,m g(j,i)=y(j)**i end do end do do i=0,n !赋值 do j=0,m u(i,j)=z(i*(m+1)+j) end do end do do i=0,n do j=0,n d1=0.d0 do k=0,n d1=d1+b(k,i)*b(k,j) end do b2(i,j)=d1 end do end do do i=0,n do j=0,n d1=0.d0 do k=0,m d1=d1+g(k,i)*g(k,j) end do g2(i,j)=d1 end do end do open (22,file='k值及误差.txt') !k值及对应误差结果输出 do k=0,n !循环判断 call inverse(b2,k,b3) call inverse(g2,k,g3) do i=0,k do j=0,n d1=0.0d0 do i2=0,k d1=d1+b3(i,i2)*b(j,i2) end do d2(i,j)=d1 end do end do do i=0,k do j=0,m d1=0.d0 do i2=0,n d1=d1+d2(i,i2)*u(i2,j) end do d3(i,j)=d1 end do end do do i=0,k do j=0,k d1=0.0d0 do i2=0,m d1=d1+d3(i,i2)*g(i2,j) end do d4(i,j)=d1 end do end do do i=0,k !矩阵相乘 do j=0,k d1=0.d0 do i2=0,k d1=d1+d4(i,i2)*g3(i2,j) end do c(i,j)=d1 end do end do delta_cal=0.d0 delta_cal=delta(c,k,z) write(22,220) k,delta_cal print*,k,delta_cal open (33,file='Matrix_C.txt') !放于do循环中,将原来的矩阵覆盖 do i=0,k write(33,'(20e20.12)') (c(i,j),j=0,k) !由计算结果,k=5时满足精度要求 end do close(33) if (delta_cal<=yi) then exit end if end do close(22) 220 format (i2,e20.12) end subroutine fit program test3_main !主程序 use param implicit none double precision tu(0:230,0:1),ut(0:1),zt(0:230),c(0:10,0:10) double precision tut(45,0:1),ztt(0:39),tt(45) double precision z,interpo,zp integer i,j,m,n integer mm,nn n=10 m=20 do i=0,10 !x,y的初值确定 equ_x(i)=0.08*i end do do j=0,20 equ_y(j)=0.5d0+0.05*j end do open (11,file='x_y_t_u_z.txt') write(11,100) 'x','y','t','u','z' do i=0,10 !循环求解11*21次,得到i,j对应的t,u,v,w值,并插值得到z do j=0,20 call equa_cal(equ_x(i),equ_y(j),tuvw) tu(i*21+j,0)=tuvw(1) tu(i*21+j,1)=tuvw(2) zt(i*21+j)=interpo(tu(i*21+j,0),tu(i*21+j,1)) write(11,200) equ_x(i),equ_y(j),tuvw(1),tuvw(2),zt(i*21+j) end do end do close(11) open (44,file='f(x,y).txt') !插值结果 write(44,441) 'xi','yi','f(xi,yi)' do i=0,n do j=0,m write(44,440) equ_x(i),equ_y(j),zt(i*21+j) end do end do close(44) call fit(zt,c) !调用子程序,曲面拟合 mm=8 !第2问求解 nn=5 do i=1,mm xx(i-1)=0.1d0*i end do do j=1,nn yy(j-1)=0.5d0+0.2d0*j end do do i=1,mm do j=1,nn call equa_cal(xx(i-1),yy(j-1),tuvw) tut((i-1)*5+j,0)=tuvw(1) tut((i-1)*5+j,1)=tuvw(2) end do end do do i=1,40 zp=interpo(tut(i,0),tut(i,1)) !第2问f(x,y) ztt(i-1)=zp end do call delta2(c,ztt,tt) open (55,file='f2(x,y).txt') write(55,551) 'x*i','y*i','f(xi,yi)','p(xi,yi)' do i=1,mm do j=1,nn write(55,550) xx(i-1),yy(j-1),ztt((i-1)*5+j-1),tt((i-1)*5+j) end do end do close(55) 100 format(2a10,3a18) 200 format(2e12.4,3e20.12) 440 format(2e12.4,e20.12) 441 format(2a10,a18) 550 format(2e20.7,2e20.12) 551 format(2a18,2a18) end program test3_main 2 程序说明 (1) function delta(c,k,z) 本函数用于计算误差的积累,被其它子程序调用实现对精度的判断,以选择达到要求的k值。按照题目给定的要求,定义函数 21020 ,,,fxypxy(,)(,),,,,ijij,,ij00 (2) function pxy(c,k,x,y) 本函数用于计算函数数值拟合结果。函数中以x,y的n次幂为基函数,并用到了曲面拟合得到的拟合系数矩阵C。 (3) function interpo(t,u) xy,本函数用于计算二元函数代数插值。由的值解非线性方程组得到对应的t,u值,ij zfxy,(,)通过插值函数,可得到对应的z值,即。插值多项式中一共含3*3项。在题目的第1、2问中,都调用到此函数。 (4) function max_inf(Delta_val,n) 本函数用于计算行列式中的元素绝对值最大值,即行列式的无穷范数。函数值用于判断迭代结果是否满足题目的精度要求。 (5) module param 由于编程习惯的因素,我习惯在编程开始时建立个module存放在主程序和部分子程序中需要用到的变量和常量。在此题目中,module param定义的量多只是用于主程序中,被主 程序调用。 (6) subroutine inverse(a,k,b) 此子程序用于求矩阵的逆,供子程序subroutine fit(z,c)调用。由于求逆相当于是对n阶矩阵解n次线性方程组,因此在子程序中n次调用了子程序subroutine gauss_cal_inverse(kk,a,b,x)求解方程,并将解得的结果合成矩阵所求的逆。 (7) subroutine delta2(c,z,tt) ****xy,题目第二问中的与第一问的值不同,因此在此子程序中通过对进行xy,xy,ijijij变量赋值,并通过调用函数pxy得到第二问数据的拟合结果。拟合的结果在主程序中以文件形式输出,并与插值结果比较。 (8) subroutine equa_cal(x,y,tuvw) 此子程序的功能在于Newton法解非线性方程组,得到方程组中的t,u,v,w的值。在此子程序中,对所需要的变量行列式赋初始值,并对方程组求导存放于4*4矩阵中,同时对方程组移项变号存放于行列式中。子程序在不满足精度的条件下通过调用子程序subroutine mainrow_cancle(a,b,x)解线性方程组,迭代直至得到满足精度的解,得到的t,u值用于后边的计算。 (9) subroutine gauss_cal_inverse(kk,a,b,x) 此子程序用于求解线性方程组,以使子程序subroutine inverse(a,k,b)中调用此子程序实现得到矩阵的逆的功能。程序中采用列主元Gauss消去法求解方程组。由于此子程序的最终目标在于实现对曲面插值基函数最高阶次的确定和插值系数矩阵的计算,在程序中对线性方程组的维数采用变维数的形式,维数用kk表示,通过对此程序的多次调用实现kk的计算和选择。 (10) subroutine mainrow_cancle(a,b,x) 此子程序用于求解线性方程组,采用列主元Gauss消去法求解。程序在Newton法求解非线性方程组过程中被调用。由于非线性方程组中含4个未知数,所以相别于子程序subroutine gauss_cal_inverse(kk,a,b,x),此程序变量中没有对维数进行定义,而是在程序中直接固定。 (11) subroutine fit(z,c) 此子程序用于计算曲面拟合的系数矩阵。在此程序中,首先对变量进行赋值,同时对矩阵做相乘。在迭代循环中,对矩阵做求逆运算,以及矩阵相乘,同时采用最小二乘 原则 组织架构调整原则组织架构设计原则组织架构设置原则财政预算编制原则问卷调查设计原则 判断 ,是否满足精度要求。程序中对选择过程中的,进行输出,同时在满足精度要求的条件k 下对系数矩阵进行输出。 (12) program test3_main xy,此部分为主程序。程序通过对进行赋值,通过调用子程序实现对非线性方程组的ij zfxy,(,)求解,得到11*21组(t,u,v,w)数据,同时通过插值,得到对应的值,并以文档形 **式对计算结果进行输出。程序通过调用子程序实现对曲面拟合系数矩阵的求解,并对xy,ij **进行赋值,分别通过解非线性方程组和插值的方式得到对应的,通过拟合计算得fxy(,)ij **到,并以文件形式同时输出。 pxy(,)ij 三 结果输出 **1 数表 xyfxy,,(,)ijij xi yi f(xi,yi) 0.0000E+00 0.5000E+00 0.446504018481E+00 0.0000E+00 0.5500E+00 0.324683261166E+00 0.0000E+00 0.6000E+00 0.210159683380E+00 0.0000E+00 0.6500E+00 0.103043603694E+00 0.0000E+00 0.7000E+00 0.340188984706E-02 0.0000E+00 0.7500E+00 -0.887358202231E-01 0.0000E+00 0.8000E+00 -0.173371639982E+00 0.0000E+00 0.8500E+00 -0.250534619127E+00 0.0000E+00 0.9000E+00 -0.320276514262E+00 0.0000E+00 0.9500E+00 -0.382668073988E+00 0.0000E+00 0.1000E+01 -0.437795774417E+00 0.0000E+00 0.1050E+01 -0.485758948724E+00 0.0000E+00 0.1100E+01 -0.526667261575E+00 0.0000E+00 0.1150E+01 -0.560638485713E+00 0.0000E+00 0.1200E+01 -0.587796543731E+00 0.0000E+00 0.1250E+01 -0.608269782929E+00 0.0000E+00 0.1300E+01 -0.622189455424E+00 0.0000E+00 0.1350E+01 -0.629688379283E+00 0.0000E+00 0.1400E+01 -0.630899759496E+00 0.0000E+00 0.1450E+01 -0.625956150287E+00 0.0000E+00 0.1500E+01 -0.614988542456E+00 0.8000E-01 0.5000E+00 0.638015222109E+00 0.8000E-01 0.5500E+00 0.506611749043E+00 0.8000E-01 0.6000E+00 0.382176361692E+00 0.8000E-01 0.6500E+00 0.264863482309E+00 0.8000E-01 0.7000E+00 0.154780190406E+00 0.8000E-01 0.7500E+00 0.519926728023E-01 0.8000E-01 0.8000E+00 -0.434680514793E-01 0.8000E-01 0.8500E+00 -0.131601068428E+00 0.8000E-01 0.9000E+00 -0.212431100097E+00 0.8000E-01 0.9500E+00 -0.286004562782E+00 0.8000E-01 0.1000E+01 -0.352386090434E+00 0.8000E-01 0.1050E+01 -0.411655467506E+00 0.8000E-01 0.1100E+01 -0.463904921838E+00 0.8000E-01 0.1150E+01 -0.509236734168E+00 0.8000E-01 0.1200E+01 -0.547761126397E+00 0.8000E-01 0.1250E+01 -0.579594395565E+00 0.8000E-01 0.1300E+01 -0.604857264739E+00 0.8000E-01 0.1350E+01 -0.623673425643E+00 0.8000E-01 0.1400E+01 -0.636168251031E+00 0.8000E-01 0.1450E+01 -0.642467657463E+00 0.8000E-01 0.1500E+01 -0.642697101484E+00 0.1600E+00 0.5000E+00 0.840081386516E+00 0.1600E+00 0.5500E+00 0.699764154764E+00 0.1600E+00 0.6000E+00 0.566061430001E+00 0.1600E+00 0.6500E+00 0.439171594548E+00 0.1600E+00 0.7000E+00 0.319242123477E+00 0.1600E+00 0.7500E+00 0.206376177056E+00 0.1600E+00 0.8000E+00 0.100638508018E+00 0.1600E+00 0.8500E+00 0.206072387448E-02 0.1600E+00 0.9000E+00 -0.893540410607E-01 0.1600E+00 0.9500E+00 -0.173626985044E+00 0.1600E+00 0.1000E+01 -0.250799972016E+00 0.1600E+00 0.1050E+01 -0.320932284774E+00 0.1600E+00 0.1100E+01 -0.384097749610E+00 0.1600E+00 0.1150E+01 -0.440382189108E+00 0.1600E+00 0.1200E+01 -0.489881164904E+00 0.1600E+00 0.1250E+01 -0.532697976849E+00 0.1600E+00 0.1300E+01 -0.568941889159E+00 0.1600E+00 0.1350E+01 -0.598726557770E+00 0.1600E+00 0.1400E+01 -0.622168636234E+00 0.1600E+00 0.1450E+01 -0.639386540257E+00 0.1600E+00 0.1500E+01 -0.650499353278E+00 0.2400E+00 0.5000E+00 0.105151507734E+01 0.2400E+00 0.5500E+00 0.902927414730E+00 0.2400E+00 0.6000E+00 0.760580249339E+00 0.2400E+00 0.6500E+00 0.624715179428E+00 0.2400E+00 0.7000E+00 0.495519736312E+00 0.2400E+00 0.7500E+00 0.373134022343E+00 0.2400E+00 0.8000E+00 0.257656727924E+00 0.2400E+00 0.8500E+00 0.149150538171E+00 0.2400E+00 0.9000E+00 0.476469654652E-01 0.2400E+00 0.9500E+00 -0.468493443612E-01 0.2400E+00 0.1000E+01 -0.134356781184E+00 0.2400E+00 0.1050E+01 -0.214913365159E+00 0.2400E+00 0.1100E+01 -0.288573720098E+00 0.2400E+00 0.1150E+01 -0.355406383286E+00 0.2400E+00 0.1200E+01 -0.415491413839E+00 0.2400E+00 0.1250E+01 -0.468918265989E+00 0.2400E+00 0.1300E+01 -0.515783897638E+00 0.2400E+00 0.1350E+01 -0.556191088040E+00 0.2400E+00 0.1400E+01 -0.590246941568E+00 0.2400E+00 0.1450E+01 -0.618061557256E+00 0.2400E+00 0.1500E+01 -0.639746846133E+00 0.3200E+00 0.5000E+00 0.127124673150E+01 0.3200E+00 0.5500E+00 0.111500199653E+01 0.3200E+00 0.6000E+00 0.964607704189E+00 0.3200E+00 0.6500E+00 0.820347345257E+00 0.3200E+00 0.7000E+00 0.682447652997E+00 0.3200E+00 0.7500E+00 0.551085182680E+00 0.3200E+00 0.8000E+00 0.426392359478E+00 0.3200E+00 0.8500E+00 0.308462968931E+00 0.3200E+00 0.9000E+00 0.197357102936E+00 0.3200E+00 0.9500E+00 0.931055942708E-01 0.3200E+00 0.1000E+01 -0.428601843970E-02 0.3200E+00 0.1050E+01 -0.948339509092E-01 0.3200E+00 0.1100E+01 -0.178573015179E+00 0.3200E+00 0.1150E+01 -0.255553802874E+00 0.3200E+00 0.1200E+01 -0.325840172789E+00 0.3200E+00 0.1250E+01 -0.389507009624E+00 0.3200E+00 0.1300E+01 -0.446638224310E+00 0.3200E+00 0.1350E+01 -0.497324969757E+00 0.3200E+00 0.1400E+01 -0.541664048802E+00 0.3200E+00 0.1450E+01 -0.579756493853E+00 0.3200E+00 0.1500E+01 -0.611706300008E+00 0.4000E+00 0.5000E+00 0.149832102673E+01 0.4000E+00 0.5500E+00 0.133499860468E+01 0.4000E+00 0.6000E+00 0.117712509493E+01 0.4000E+00 0.6500E+00 0.102502402501E+01 0.4000E+00 0.7000E+00 0.878959992192E+00 0.4000E+00 0.7500E+00 0.739145076979E+00 0.4000E+00 0.8000E+00 0.605744839252E+00 0.4000E+00 0.8500E+00 0.478883828551E+00 0.4000E+00 0.9000E+00 0.358650593313E+00 0.4000E+00 0.9500E+00 0.245102203798E+00 0.4000E+00 0.1000E+01 0.138268318921E+00 0.4000E+00 0.1050E+01 0.381548340039E-01 0.4000E+00 0.1100E+01 -0.552528517584E-01 0.4000E+00 0.1150E+01 -0.141986910389E+00 0.4000E+00 0.1200E+01 -0.222094467461E+00 0.4000E+00 0.1250E+01 -0.295635259426E+00 0.4000E+00 0.1300E+01 -0.362679536888E+00 0.4000E+00 0.1350E+01 -0.423306187853E+00 0.4000E+00 0.1400E+01 -0.477601057837E+00 0.4000E+00 0.1450E+01 -0.525655446285E+00 0.4000E+00 0.1500E+01 -0.567564761030E+00 0.4800E+00 0.5000E+00 0.173189270866E+01 0.4800E+00 0.5500E+00 0.156203454383E+01 0.4800E+00 0.6000E+00 0.139721688339E+01 0.4800E+00 0.6500E+00 0.123780097070E+01 0.4800E+00 0.7000E+00 0.108408749564E+01 0.4800E+00 0.7500E+00 0.936322734519E+00 0.4800E+00 0.8000E+00 0.794704410728E+00 0.4800E+00 0.8500E+00 0.659387159406E+00 0.4800E+00 0.9000E+00 0.530487548152E+00 0.4800E+00 0.9500E+00 0.408088646927E+00 0.4800E+00 0.1000E+01 0.292244163087E+00 0.4800E+00 0.1050E+01 0.182982169313E+00 0.4800E+00 0.1100E+01 0.803084573091E-01 0.4800E+00 0.1150E+01 -0.157904487575E-01 0.4800E+00 0.1200E+01 -0.105344586107E+00 0.4800E+00 0.1250E+01 -0.188398123682E+00 0.4800E+00 0.1300E+01 -0.265007180792E+00 0.4800E+00 0.1350E+01 -0.335237868599E+00 0.4800E+00 0.1400E+01 -0.399164531631E+00 0.4800E+00 0.1450E+01 -0.456868168928E+00 0.4800E+00 0.1500E+01 -0.508435016628E+00 0.5600E+00 0.5000E+00 0.197122174854E+01 0.5600E+00 0.5500E+00 0.179532955996E+01 0.5600E+00 0.6000E+00 0.162406707221E+01 0.5600E+00 0.6500E+00 0.145783054044E+01 0.5600E+00 0.7000E+00 0.129695460646E+01 0.5600E+00 0.7500E+00 0.114171806136E+01 0.5600E+00 0.8000E+00 0.992349488679E+00 0.5600E+00 0.8500E+00 0.849032618359E+00 0.5600E+00 0.9000E+00 0.711911307203E+00 0.5600E+00 0.9500E+00 0.581094113999E+00 0.5600E+00 0.1000E+01 0.456658468676E+00 0.5600E+00 0.1050E+01 0.338654452169E+00 0.5600E+00 0.1100E+01 0.227108212603E+00 0.5600E+00 0.1150E+01 0.122025047040E+00 0.5600E+00 0.1200E+01 0.233921787024E-01 0.5600E+00 0.1250E+01 -0.688187414915E-01 0.5600E+00 0.1300E+01 -0.154649382128E+00 0.5600E+00 0.1350E+01 -0.234152702585E+00 0.5600E+00 0.1400E+01 -0.307391126074E+00 0.5600E+00 0.1450E+01 -0.374434894373E+00 0.5600E+00 0.1500E+01 -0.435360586262E+00 0.6400E+00 0.5000E+00 0.221566781955E+01 0.6400E+00 0.5500E+00 0.203420108776E+01 0.6400E+00 0.6000E+00 0.185695509626E+01 0.6400E+00 0.6500E+00 0.168435811551E+01 0.6400E+00 0.7000E+00 0.151677630268E+01 0.6400E+00 0.7500E+00 0.135451899060E+01 0.6400E+00 0.8000E+00 0.119784403552E+01 0.6400E+00 0.8500E+00 0.104696299790E+01 0.6400E+00 0.9000E+00 0.902046032159E+00 0.6400E+00 0.9500E+00 0.763226426125E+00 0.6400E+00 0.1000E+01 0.630604770752E+00 0.6400E+00 0.1050E+01 0.504252763956E+00 0.6400E+00 0.1100E+01 0.384216665685E+00 0.6400E+00 0.1150E+01 0.270520427776E+00 0.6400E+00 0.1200E+01 0.163168524738E+00 0.6400E+00 0.1250E+01 0.621485118599E-01 0.6400E+00 0.1300E+01 -0.325666640541E-01 0.6400E+00 0.1350E+01 -0.121016577714E+00 0.6400E+00 0.1400E+01 -0.203251440524E+00 0.6400E+00 0.1450E+01 -0.279330398316E+00 0.6400E+00 0.1500E+01 -0.349319993987E+00 0.7200E+00 0.5000E+00 0.246468417214E+01 0.7200E+00 0.5500E+00 0.227805892712E+01 0.7200E+00 0.6000E+00 0.209525119701E+01 0.7200E+00 0.6500E+00 0.191671807283E+01 0.7200E+00 0.7000E+00 0.174285457249E+01 0.7200E+00 0.7500E+00 0.157399837017E+01 0.7200E+00 0.8000E+00 0.141043477742E+01 0.7200E+00 0.8500E+00 0.125240169239E+01 0.7200E+00 0.9000E+00 0.110009435123E+01 0.7200E+00 0.9500E+00 0.953669792929E+00 0.7200E+00 0.1000E+01 0.813250997212E+00 0.7200E+00 0.1050E+01 0.678930685454E+00 0.7200E+00 0.1100E+01 0.550774791740E+00 0.7200E+00 0.1150E+01 0.428825621175E+00 0.7200E+00 0.1200E+01 0.313104717120E+00 0.7200E+00 0.1250E+01 0.203615460798E+00 0.7200E+00 0.1300E+01 0.100345426589E+00 0.7200E+00 0.1350E+01 0.326851531060E-02 0.7200E+00 0.1400E+01 -0.876531138281E-01 0.7200E+00 0.1450E+01 -0.172467293594E+00 0.7200E+00 0.1500E+01 -0.251230264216E+00 0.8000E+00 0.5000E+00 0.271781105246E+01 0.8000E+00 0.5500E+00 0.252639944245E+01 0.8000E+00 0.6000E+00 0.233841132645E+01 0.8000E+00 0.6500E+00 0.215432931549E+01 0.8000E+00 0.7000E+00 0.197457449369E+01 0.8000E+00 0.7500E+00 0.179951051520E+01 0.8000E+00 0.8000E+00 0.162944815596E+01 0.8000E+00 0.8500E+00 0.146464997869E+01 0.8000E+00 0.9000E+00 0.130533490236E+01 0.8000E+00 0.9500E+00 0.115168255603E+01 0.8000E+00 0.1000E+01 0.100383735488E+01 0.8000E+00 0.1050E+01 0.861912272681E+00 0.8000E+00 0.1100E+01 0.725992307271E+00 0.8000E+00 0.1150E+01 0.596137648608E+00 0.8000E+00 0.1200E+01 0.472386566146E+00 0.8000E+00 0.1250E+01 0.354758035484E+00 0.8000E+00 0.1300E+01 0.243254125325E+00 0.8000E+00 0.1350E+01 0.137862165422E+00 0.8000E+00 0.1400E+01 0.385567151676E-01 0.8000E+00 0.1450E+01 -0.546986489666E-01 0.8000E+00 0.1500E+01 -0.141949710438E+00 值 2 选择过程的k,, 0 0.144288072856E+03 1 0.322090910310E+01 2 0.465995962880E-02 3 0.172117546878E-03 4 0.330953468916E-05 5 0.254138772222E-07 , 由计算可知可知,达到精度要求的=5,=0.254138772222E-07。 k 3 pxy(,)中的系数矩阵 0.202123054827E+01 -0.366842687735E+01 0.709248647094E+00 0.848605426028E+00 -0.415897446219E+00 0.674319972750E-01 0.319190893717E+01 -0.741110296047E+00 -0.269712471124E+01 0.163118352392E+01 -0.484720025328E+00 0.606142897886E-01 0.256889784709E+00 0.157991867606E+01 -0.463408198208E+00 -0.813443139195E-01 0.102094203234E+00 -0.210152261425E-01 -0.269260300556E+00 -0.730247711297E+00 0.107614520472E+01 -0.807012951933E+00 0.302872874774E+00 -0.459726490662E-01 0.217459720487E+00 -0.178372284630E+00 -0.724059676286E-01 0.243330626865E+00 -0.141334814136E+00 0.265102556150E-01 -0.559032532792E-01 0.143199192891E+00 -0.136270281044E+00 0.407195391454E-01 0.377507134544E-02 -0.266770853341E-02 ****** 4 数表 xyfxypxy,,(,),(,)ijijij x*i y*i f(xi,yi) p(xi,yi) 0.1000000E+00 0.7000000E+00 0.194720407918E+00 0.194730345568E+00 0.1000000E+00 0.9000000E+00 -0.183037079189E+00 -0.183041851510E+00 0.1000000E+00 0.1100000E+01 -0.445497646915E+00 -0.445500050597E+00 0.1000000E+00 0.1300000E+01 -0.597566707641E+00 -0.597558856792E+00 0.1000000E+00 0.1500000E+01 -0.646459593901E+00 -0.646446096505E+00 0.2000000E+00 0.7000000E+00 0.405979189288E+00 0.405989523192E+00 0.2000000E+00 0.9000000E+00 -0.225159583746E-01 -0.225211325675E-01 0.2000000E+00 0.1100000E+01 -0.338220816040E+00 -0.338224033861E+00 0.2000000E+00 0.1300000E+01 -0.544437831522E+00 -0.544430451918E+00 0.2000000E+00 0.1500000E+01 -0.647361338568E+00 -0.647347996223E+00 0.3000000E+00 0.7000000E+00 0.634777195151E+00 0.634787431284E+00 0.3000000E+00 0.9000000E+00 0.158801168839E+00 0.158796275492E+00 0.3000000E+00 0.1100000E+01 -0.207365694171E+00 -0.207368594549E+00 0.3000000E+00 0.1300000E+01 -0.465357906898E+00 -0.465349925805E+00 0.3000000E+00 0.1500000E+01 -0.620270953075E+00 -0.620257124337E+00 0.4000000E+00 0.7000000E+00 0.878960023174E+00 0.878969838435E+00 0.4000000E+00 0.9000000E+00 0.358650625882E+00 0.358646018554E+00 0.4000000E+00 0.1100000E+01 -0.552528211681E-01 -0.552554541106E-01 0.4000000E+00 0.1300000E+01 -0.362679511503E+00 -0.362671067149E+00 0.4000000E+00 0.1500000E+01 -0.567564743655E+00 -0.567550568365E+00 0.5000000E+00 0.7000000E+00 0.113661091016E+01 0.113662032098E+01 0.5000000E+00 0.9000000E+00 0.574980340948E+00 0.574975813546E+00 0.5000000E+00 0.1100000E+01 0.115992376792E+00 0.115989301133E+00 0.5000000E+00 0.1300000E+01 -0.238568304012E+00 -0.238560425080E+00 0.5000000E+00 0.1500000E+01 -0.491434393656E+00 -0.491420886205E+00 0.6000000E+00 0.7000000E+00 0.140604179891E+01 0.140605064928E+01 0.6000000E+00 0.9000000E+00 0.805941494063E+00 0.805937268290E+00 0.6000000E+00 0.1100000E+01 0.304429221045E+00 0.304425808669E+00 0.6000000E+00 0.1300000E+01 -0.950161300996E-01 -0.950089530931E-01 0.6000000E+00 0.1500000E+01 -0.393902307746E+00 -0.393889823049E+00 0.7000000E+00 0.7000000E+00 0.168578351531E+01 0.168579117406E+01 0.7000000E+00 0.9000000E+00 0.104988115306E+01 0.104987770135E+01 0.7000000E+00 0.1100000E+01 0.508293783940E+00 0.508291018474E+00 0.7000000E+00 0.1300000E+01 0.661487967065E-01 0.661563470006E-01 0.7000000E+00 0.1500000E+01 -0.276834341778E+00 -0.276822028305E+00 0.8000000E+00 0.7000000E+00 0.197457455665E+01 0.197458121264E+01 0.8000000E+00 0.9000000E+00 0.130533496765E+01 0.130533196186E+01 0.8000000E+00 0.1100000E+01 0.725992371111E+00 0.725989281702E+00 0.8000000E+00 0.1300000E+01 0.243254184181E+00 0.243260781141E+00 0.8000000E+00 0.1500000E+01 -0.141949659709E+00 -0.141938772513E+00
本文档为【北航数值分析上机作业第3题】,请使用软件OFFICE或WPS软件打开。作品中的文字与图均可以修改和编辑, 图片更改请在作品中右键图片并更换,文字修改请直接点击文字进行修改,也可以新增和删除文档中的内容。
该文档来自用户分享,如有侵权行为请发邮件ishare@vip.sina.com联系网站客服,我们会及时删除。
[版权声明] 本站所有资料为用户分享产生,若发现您的权利被侵害,请联系客服邮件isharekefu@iask.cn,我们尽快处理。
本作品所展示的图片、画像、字体、音乐的版权可能需版权方额外授权,请谨慎使用。
网站提供的党政主题相关内容(国旗、国徽、党徽..)目的在于配合国家政策宣传,仅限个人学习分享使用,禁止用于任何广告和商用目的。
下载需要: 免费 已有0 人下载
最新资料
资料动态
专题动态
is_196623
暂无简介~
格式:doc
大小:92KB
软件:Word
页数:41
分类:
上传时间:2017-09-30
浏览量:77