资源预览内容
第1页 / 共25页
第2页 / 共25页
第3页 / 共25页
第4页 / 共25页
第5页 / 共25页
第6页 / 共25页
第7页 / 共25页
第8页 / 共25页
第9页 / 共25页
第10页 / 共25页
亲,该文档总共25页,到这儿已超出免费预览范围,如果喜欢就下载吧!
资源描述
数值分析第三次大作业姓名:袁二凯学号:SY0917145联系方式:13488854452642101136qq.com2011年 12月 8日一、 算法的设计方案:(一) 、总体方案设计:(1)解非线性方程组。将给定的 当作已知量代入题目给定的非线性方程组,(,)ixy求得与 相对应的数组 tij,uij。(,)ixy(2)分片二次代数插值。通过分片二次代数插值运算,得到与数组 t1121,u1121对应的数组 z1121,得到二元函数 z= 。(,)ifxy(3)曲面拟合。利用 xi,yj,z1121建立二维函数表,再根据精度的要求选择适当 k值,并得到曲面拟合的系数矩阵 Crs。(4)观察和 的逼近效果。观察逼近效果只需要重复上面(1)和(2)的过(,)ipxy程,得到与新的插值节点 对应的 ,再与对应的 比较即可,这里i (,)ifxy(,)ipxy求解 可以直接使用(3)中的 Crs和 k。(,)ixy(二)具体算法设计:(1)解非线性方程组牛顿法解方程组 的解 ,可采用如下算法:()0Fx*1)在 附近选取 ,给定精度水平 和最大迭代次数 M。*x()D02)对于 执行,1kML 计算 和 。()kFx()k 求解关于 的线性方程组()()()()kkxFx 若 ,则取 ,并停止计算;否则转。()()k*() 计算 。(1)()()kkkxx 若 ,则继续,否则,输出 M 次迭代不成功的信息,并停止计算。M(2)分片双二次插值给定已知数表以及需要插值的节点,进行分片二次插值的算法:设已知数表中的点为: ,需要插值的节点为0(,1),ijxihnyjmL。(,)xy1) 根据 选择插值节点 :(,)xy(,)ijxy若 或 ,插值节点对应取 或 ,12h12nh1in若 或 ,插值节点对应取 或 。yyjim若,2,2iijjhxxinyyjm则选择 为插值节点。(,)1,;1,)krxirj2)计算11()(,)()(1,)itktikjtrtjrxliyl j%插值多项式的公式为:1(,)()(,)jikrkrkirjpxylxyf注:本步进行插值运算的是 ,利用 与 的对应关系就可以得到,tu,ij,tu与 的对应关系。z(,)ijxy(3)曲面拟合根据插值得到的数表 进行曲面拟合的过程:,(,)ijijxyf1) 根据拟合节点和基底函数写出矩阵 B 和 G:0101101()()()()kknnxxBLMO0101101()()()()kkmmyyLMO2) 计算 。1TTCBUG在这里,为了简化计算和编程、避免矩阵求逆,记:,1()TA1()TTD对上面两式进行变形,得到如下两个线性方程组:,()TTBAU()TTGD通过解上述两个线性方程组,则有: TCA3) 对于每一个 , 。(,)ijxy*0(,)()krsij sijrpxyxy4) 拟合需要达到的精度条件为:。*270(,)1nmijijijxyu其中 对应着插值得到的数表 中 的值。iju,(,)ijijf(,)ijfxy5) 让 k 逐步增加,每一次重复执行以上几步,直到成立。此时的 k 值就是要求解最小的 k。*270(,)10nmijijijpxyu二、 源程序:#include#include#include #include #include#include #define Epsilon1 1e-12 /*解线性方程组时近似解向量的精度*/#define M 200 /*解线性方程组时的最大迭代次数*/#define N 10 /*求解迭代次数时假设的 k的最大值,用于定义包含 k的存储空间*/void Newton(); /*牛顿法求解非线性方程组子程序*/void fpeccz(); /*分片二次代数插值子程序*/void qmnh(); /*曲面拟合子程序*/void duibi(); /*对比和 p逼近效果的子程序*/double x11,y21,t1121,u1121;/*定义全局变量*/double z1121,C1010;double kz;void Newton(double x11,double y21)/*牛顿法求解非线性方程组子程序*/ double X4,dx4,F4,dF44,temp,m,fx,fX;int i,j,k,l,p,ik,n; for(i=0;i=0;k-)temp=0;for(l=k+1;l0.7) sij=4;else for(m=2;m0.2*m-0.1)&(tij1.4) rij=4; elsefor(m=2;m0.4*m-0.2)&(uij1e-7) for(i=0;i=0;m-)temp=0;for(i=m+1;i=0;m-)temp=0;for(i=m+1;i=k;i+) temp+=Dil*GtG1mi; Dml=(Gtml-temp)/GtG1mm;for(i=0;i=k;i+)for(j=0;j=k;j+)temp=0; for(m=0;m=20;m+) temp+=Aim*Djm; Cij=temp;sigma=0;/*归零,开始计算 sigma值*/for(i=0;i=10;i+)for(j=0;j=20;j+)pij=0;for(i1=0;i1=k;i1+)for(j1=0;j1=k;j1+) pij+=Ci1j1*pow(xi,i1)*pow(yj,j1); sigma+=pow(pij-zij,2);printf(k=%d sigma=%.12en,k,sigma); k=k+1;k-;printf(达到精度要求时的 k和 sigma值:n);printf(k=%d sigma=%.12en,k,sigma); kz=k; /*定义为全局变量,便于 duibi()调用*/void duibi() /*对比和 p逼近效果的子程序*/ int i,j,i1,j1;double p85;for(i=0;i=7;i+) for(j=0;j=4;j+)xi=0.1*(i+1);yj=0.5+0.2*(j+1); /*重新输入节点*/Newton(x,y);fpeccz(t,u); /*求解( x*,y *)*/ for(i=0;i=7;i+) /*求解 p(x *,y *)*/for(j=0;j=4;j+)pij=0;for(i1=0;i1=kz;i1+) for(j1=0;j1=kz;j1+) pij+=Ci1j1*pow(xi,i1)*pow(yj,j1); printf(x%d=%.6f y%d=%.6fn,i+1,xi,j+1,yj);printf(f(x%d,y%d)=%.12en,i+1,j+1,zij);printf(p(x%d,y%d)=%.12en,i+1,j+1,pij);printf(=%.12en,zij-pij); /*数表 xi*,yi*, (xi*,yi*),p(x i*,yi*)*/void main() int i,j;for(i=0;i=10;i+) for(j=0;j=20;j+)xi=0.08*i;yj=0.5+0.05*j; /*输入节点*/Newton(x,y);fpeccz(t,u);for(i=0;i=10;i+) for(j=0;j=20;j+)printf(x%d=%.6f y%d=%.6f z%d%d=%.12e n,i,xi,j,yj,i,j,zij); /*数表:x i,yi, (xi,yi)*/qmnh();for(i=0;i=kz;i+) for(j=0;j=kz;j+) printf(C%d%d=%.12e n,i,j,Cij);/*数表:Crs*/duibi();三、运行结果输出1、数表数表:x i,yi,(xi,yi)x0=0.000000 y0=0.500000 z00=4.465040184807e-001x0=0.000000 y1=0.550000 z01=3.246832629277e-001x0=0.000000 y2=0.600000 z02=2.101596866827e-001x0=0.000000 y3=0.650000 z03=1.030436083160e-001x0=0.000000 y4=0.700000 z04=3.401895562675e-003x0=0.000000 y5=0.750000 z05=-8.873581363800e-002x0=0.000000 y6=0.800000 z06=-1.733716327497e-001x0=0.000000 y7=0.850000 z07=-2.505346114666e-001x0=0.000000 y8=0.900000 z08=-3.202765063876e-001x0=0.000000 y9=0.950000 z09=-3.826680661097e-001x0=0.000000 y10=1.000000 z010=-4.377957667384e-001x0=0.000000 y11=1.050000 z011=-4.857589414438e-001x0=0.000000 y12=1.100000 z012=-5.266672548835e-001x0=0.000000 y13=1.150000 z013=-5.606384797965e-001x0=0.000000 y14=1.200000 z014=-5.877965387677e-001x0=0.000000 y15=1.250000 z015=-6.082697790899e-001x0=0.000000 y16=1.300000 z016=-6.221894528764e-001x0=0.000000 y17=1.350000 z017=-6.296883781856e-001x0=0.000000 y18=1.400000 z018=-6.308997600028e-001x0=0.000000 y19=1.4
网站客服QQ:2055934822
金锄头文库版权所有
经营许可证:蜀ICP备13022795号 | 川公网安备 51140202000112号