北航数值分析B大作业第三题
北航数值分析大作业 第一题 幂法与反幂法
数 值 分 析(B ) 大 作 业(一)姓名: 学号: 电话:1、算法设计:①求1λ、501λ和s λ的值:s λ:s λ表示矩阵的按模最小特征值,为求得s λ直接对待求矩阵A 应用反幂法即可。
1λ、501λ:若矩阵A 的特征值满足关系 1n λλ<<且1n λλ≠,要求1λ、及501λ时,可按如下方法求解: a . 对矩阵A 用幂法,求得按模最大的特征值1m λ。
b . 按平移量1m λ对矩阵A 进行原点平移得矩阵1m BA I λ=+,对矩阵B 用反幂法求得B 的按模最小特征值2m λ。
c . 321m m m λλλ=-则:113min(,)m m λλλ=,13max(,)n m m λλλ=即为所求。
②求和A 的与数5011140k k λλμλ-=+最接近的特征值ik λ(k=0,1,…39):求矩阵A 的特征值中与P 最接近的特征值的大小,采用原点平移的方法:先求矩阵 B=A-PI 对应的按模最小特征值k β,则k β+P 即为矩阵A 与P 最接近的特征值。
在本次计算实习中则是先求平移矩阵k B A I μ=-,对该矩阵应用反幂法求得s λ,则与k μ最接近的A 的特征值为:s P λ+重复以上过程39次即可求得ik λ(k=0,1,…39)的值。
③求A 的(谱范数)条件数2cond()A 和行列式det A :在(1)中用反幂法求矩阵A 的按模最小特征值时,要用到Doolittle 分解方法,在Doolittle 分解完成后得到的两个矩阵分别为L 和U ,则A 的行列式可由U 阵求出,即:det(A)=det(U)。
求得det(A)不为0,因此A 为非奇异的实对称矩阵,则: max 2()scond A λλ=,max λ和s λ分别为模最大特征值与模最小特征值。
2、程序源代码:#include "Stdio.h"#include "Conio.h"#include "math.h"//****************************************************************************// // 在存储带状矩阵时,下面的几个量在程序中反复用到,为方便编程故把它们定义成宏.// // M :转换后的矩阵的行数,M=R+S+1。
数值计算B大作业---精品模板
课程设计课程名称:数值计算B设计题目:数值计算B大作业学号:姓名:完成时间:题目一:多项式插值某气象观测站在8:00(AM )开始每隔10分钟对天气作如下观测,用三次多项式插值函数(Newton )逼近如下曲线,插值节点数据如上表,并求出9点30分该地区的温度(x=10)。
二、数学原理假设有n+1个不同的节点及函数在节点上的值(x 0,y 0),……(x n ,y n ),插值多项式有如下形式:)())(()()()(n 10n 102010n x -x )(x -x x -x x P x x x x x x -⋯⋯-+⋯⋯+-++=αααα (1) 其中系数i α(i=0,1,2……n )为特定系数,可由插值样条i i n y x P =)((i=0,1,2……n )确定。
根据均差的定义,把x 看成[a,b ]上的一点,可得f(x)= f (0x )+f [10x x ,](0x -x )f [x , 0x ]= f[10x x ,]+f [x,10x x ,] (1x -x )……f[x , 0x ,…x 1-n ]= f [x, 0x ,…x n ]+ f [x , 0x ,…x n ](x —x n )综合以上式子,把后一式代入前一式,可得到:f (x )= f [0x ]+f[10x x ,](0x -x )+ f[210x x x ,,](0x -x )(1x -x )+ …+ f[x, 0x ,…x n ](0x -x )…(x —x 1-n )+ f [x , 0x ,…x n ,x ])(x 1n +ω= N n (x )+)(x n R 其中N n (x )= f[0x ]+f [10x x ,](0x -x )+ f [210x x x ,,](0x -x )(1x -x )+…+ f [x , 0x ,…x n ](0x -x )…(x —x 1-n ) (2))(x n R = f (x)— N n (x )= f [x, 0x ,…x n ,x ])(x 1n +ω (3) )(x 1n +ω=(0x -x )…(x —x n )Newton 插值的系数i α(i=0,1,2……n )可以用差商表示。
数值分析大作业三、四、五、六、七
大作业 三1. 给定初值0x 及容许误差,编制牛顿法解方程f (x )=0的通用程序.解:Matlab 程序如下:函数m 文件:function Fu=fu(x) Fu=x^3/3-x; end函数m 文件:function Fu=dfu(x) Fu=x^2-1; end用Newton 法求根的通用程序 clear;x0=input('请输入初值x0:'); ep=input('请输入容许误差:'); flag=1;while flag==1x1=x0-fu(x0)/dfu(x0); if abs(x1-x0)<ep flag=0; end x0=x1; endfprintf('方程的一个近似解为:%f\n',x0); 寻找最大δ值的程序: cleareps=input('请输入搜索精度:'); ep=input('请输入容许误差:'); flag=1; k=0; x0=0; while flag==1 sigma=k*eps; x0=sigma; k=k+1; m=0; flag1=1;while flag1==1 && m<=10^3 x1=x0-fu(x0)/dfu(x0);if abs(x1-x0)<epm=m+1; x0=x1; endif flag1==1||abs(x0)>=ep flag=0; end endfprintf('最大的sigma 值为:%f\n',sigma);2.求下列方程的非零根5130.6651()ln 05130.665114000.0918x x f x x +⎛⎫=-= ⎪-⨯⎝⎭解:Matlab 程序为:(1)主程序 clear clc format long x0=765; N=100;errorlim=10^(-5); x=x0-f(x0)/subs(df(),x0); n=1; while n<Nx=x0-f(x0)/subs(df(),x0); if abs(x-x0)>errorlim n=n+1; else break ; end x0=x; enddisp(['迭代次数: n=',num2str(n)])disp(['所求非零根: 正根x1=',num2str(x),' 负根x2=',num2str(-x)]) (2)子函数 非线性函数f function y=f(x)y=log((513+*x)/*x))-x/(1400*; end(3)子函数 非线性函数的一阶导数df function y=df() syms x1y=log((513+*x1)/*x1))-x1/(1400*; y=diff(y); end运行结果如下:所求非零根: 正根x1= 负根x2=大作业 四试编写MATLAB 函数实现Newton 插值,要求能输出插值多项式. 对函数21()14f x x =+在区间[-5,5]上实现10次多项式插值.分析:(1)输出插值多项式。
北航数值分析大作业三
一、题目:关于x, y, t, u, v, w 的下列方程组0.5cos 2.670.5sin 1.070.5cos 3.740.5sin 0.79t u v w x t u v w y t u v w x t u v w y +++-=⎧⎪+++-=⎪⎨+++-=⎪⎪+++-=⎩1、试用数值方法求出f(x, y)在区域 {(,)|00.8,0.5 1.5}D x y x y =≤≤≤≤上的一个近似表达式,0(,)kr s rsr s p x y cx y ==∑要求(,)p x y 一最小的k 值达到以下的精度10202700((,)(,))10i j i j i j f x y p x y σ-===-≤∑∑其中,0.08,0.50.05i j x i y j ==+。
2、计算****(,),(,)i j i j f x y p x y (i = 1, 2, …,8;j = 1, 2,…,5)的值,以观察(,)p x y 逼近(,)f x y 的效果,其中,*i x =0.1i , *j y =0.5+0.2j 。
说明:1、用迭代方法求解非线性方程组时,要求近似解向量()k x 满足()(1)()12||||/||||10k k k x x x --∞∞-≤2、作二元插值时,要使用分片二次代数插值。
3、要由程序自动确定最小的k 值。
4、打印以下内容:●算法的设计方案。
●全部源程序(要求注明主程序和每个子程序的功能)。
●数表:,,i j x y (,)i j f x y (i = 0,1,2,…,10;j = 0,1,2,…,20)。
●选择过程的,k σ值。
●达到精度要求时的,k σ值以及(,)p x y 中的系数rs c (r = 0,1,…,k;s = 0,1,…,k )。
●数表:**,,i j x y ****(,),(,)i j i j f x y p x y (i = 1, 2, ...,8;j = 1, 2, (5)。
北航研究生数值分析试题
∗⎞ ⎟的 A1 ⎠
矩阵。
三、(12 分)试用高斯列主元素法求解线性方程组
⎡ 1 3 −2 −4 ⎤ ⎡ x1 ⎤ ⎡3 ⎤ ⎢ 2 6 −7 −10 ⎥ ⎢ x ⎥ ⎢ −2 ⎥ ⎢ ⎥⎢ 2⎥ = ⎢ ⎥ ⎢ −1 −1 5 9 ⎥ ⎢ x3 ⎥ ⎢14 ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥ ⎢ x4 ⎦ ⎥ ⎣ −6 ⎦ ⎣ −3 −5 0 15 ⎦ ⎣ 四、(12 分)利用矩阵 A 的三角分解 A = LU 求解下列方程组 ⎛ 1 2 1 ⎞ ⎛ x1 ⎞ ⎛ 0 ⎞ ⎜ ⎟⎜ ⎟ ⎜ ⎟ ⎜ 2 2 3 ⎟ ⎜ x2 ⎟ = ⎜ 3 ⎟ ⎜ −1 −3 0 ⎟ ⎜ x ⎟ ⎜ 2 ⎟ ⎝ ⎠⎝ 3 ⎠ ⎝ ⎠
第一章
1、近似数 x = 0.231 关于真值 x = 0.229 有( (1)1;(2)2;(3)3;(4)4。
∗
绪论
一、选择题(四个选项中仅有一项符合题目要求,每小题 3 分,共计 15 分) )位有效数字。
2、取 3 ≈ 1.732 计算 x = ( 3 − 1) ,下列方法中哪种最好?(
4
)
Ax
∞和
A ∞ 的值分别为(
)
3
(1) 8 , 8 ;
(2) 8 , 7 ;
(3) 8 , 6 ;
(4) 7 , 7 。
5 、若解线性代数方程组的 Gauss 部分选主元方法第二步得到的系数矩阵的第三列向量为
(2
6 3 2 −5 4 2 ) ,则第三步主行是(
T
) (4) 第 6 行。
(1) 第 2 行;
1 − cos x , sin x
x ≠ 0且 x << 1 ;
(2)
1 1− x , − 1+ 2x 1+ x
北航数值分析报告大作业第三题(fortran)
“数值分析“计算实习大作业第三题——SY1415215孔维鹏一、计算说明1、将x i=0.08i,y j=0.5+0.05j分别代入方程组(A.3)得到关于t,u,v,w的的方程组,调用离散牛顿迭代子函数求出与x i,y j对应的t i,u j。
2、调用分片二次代数插值子函数在点(t i,u j)处插值得到z(x i,y j)=f(x i,y j),得到数表(x i,y j,f(x i,y j))。
3、对于k=1,2,3,4⋯,分别调用最小二乘拟合子函数计算系数矩阵c rs及误差σ,直到满足精度,即求得最小的k值及系数矩阵c rs。
4、将x i∗=0.1i,y j∗=0.5+0.2j分别代入方程组(A.3)得到关于t∗,u∗,v∗,w∗的的方程组,调用离散牛顿迭代子函数求出与x i∗,y j∗对应的t i∗,u j∗,调用分片二次代数插值子函数在点(t i∗,u j∗)处插值得到z∗(x i∗,y j∗)=f(x i∗,y j∗);调用步骤3中求得的系数矩阵c rs求得p(x i∗,y j∗),打印数表(x i∗,y j∗,f(x i∗,y j∗),p(x i∗,y j∗))。
二、源程序(FORTRAN)PROGRAM SY1415215DIMENSION X(11),Y(21),T(6),U(6),Z(6,6),UX(11,21),TY(11,21),FXY(11,21),C(6,6) DIMENSION X1(8),Y1(5),FXY1(8,5),PXY1(8,5),UX1(8,5),TY1(8,5)REAL(8) X,Y,T,U,Z,FXY,UX,TY,C,E,X1,Y1,FXY1,PXY1,UX1,TY1OPEN (1,FILE='第三题计算结果.TXT')DO I=1,11X(I)=0.08*(I-1)ENDDODO I=1,21Y(I)=0.5+0.05*(I-1)ENDDO!*****求解非线性方程组,得到z=f(t,u)的函数*******DO I=1,11DO J=1,21CALL DISNEWTON_NONLINEAR(X(I),Y(J),UX(I,J),TY(I,J)) ENDDOENDDO!*************分片二次插值得到z=f(x,y)***********DO I=1,11DO J=1,21CALL INTERPOLATION(UX(I,J),TY(I,J),FXY(I,J))ENDDOENDDOWRITE (1,"('数表(x,y,f(x,y)):')")WRITE (1,"(3X,'X',7X,'Y',10X,'F(X,Y)')")DO I=1,11DO J=1,21WRITE(1,'(1X,F5.2,2X,F5.3,2X,E20.13)') X(I),Y(J),FXY(I,J) ENDDOWRITE (1,"('')")ENDDO!***********最小二乘拟合得到P(x,y)**************N=11M=21WRITE (1,'(" ","K和σ分别为:")')DO K=1,20CALL LSFITTING(X,Y,FXY,C,N,M,K,K,E)WRITE (1,'(I3,2X,E20.13)') K-1,EIF(E<10E-7) EXITENDDOWRITE(1,'(" ")')WRITE(1,'("系数矩阵Crs(按行)为:")')DO I=1,KDO J=1,KWRITE (1,'(E20.13,2X,\)') C(I,J)ENDDOWRITE (1,"('')")WRITE (*,"('')")ENDDODO I=1,8X1(I)=0.1*IENDDODO J=1,5Y1(J)=0.5+0.2*JENDDODO I=1,8DO J=1,5CALL DISNEWTON_NONLINEAR(X1(I),Y1(J),UX1(I,J),TY1(I,J))ENDDOENDDODO I=1,8DO J=1,5CALL INTERPOLATION(UX1(I,J),TY1(I,J),FXY1(I,J))ENDDOENDDOPXY1=0DO I=1,8DO J=1,5DO II=1,KDO JJ=1,KPXY1(I,J)=PXY1(I,J)+C(II,JJ)*(X1(I)**(II-1))*(Y1(J)**(JJ-1)) ENDDOENDDOENDDOENDDOWRITE(1,'(" ")')WRITE(1,'("数表(x,y,f(x,y),p(x,y)):")')WRITE(1,"(2X,'X',6X,'Y',12X,'F(X,Y)',14X,'P(X,Y)')")DO I=1,8DO J=1,5WRITE(1,'(F5.3,2X,F5.3,2X,E20.13,2X,E20.13)') X1(I),Y1(J),FXY1(I,J),PXY1(I,J) ENDDOWRITE (1,"('')")ENDDOCLOSE (1)END!***********用离散牛顿法求解非线性方程组****************SUBROUTINE DISNEWTON_NONLINEAR(X1,Y1,U,T)PARAMETER (N=4)REAL EPS !EPS为迭代精度,M为最大迭代次数DIMENSION X(N),H(N),Y(N),JA(N,N),E(N),XK(N)REAL(8) JA,X,H,Y,E,XK,U,T,V,W,X1,Y1,E1,E2F1(T,U,V,W)=0.5*COS(T)+U+V+W-X1-2.67F2(T,U,V,W)=T+0.5*SIN(U)+V+W-Y1-1.07F3(T,U,V,W)=0.5*T+U+COS(V)+W-X1-3.74F4(T,U,V,W)=T+0.5*U+V+SIN(W)-Y1-0.79EPS=10E-12M=100X=1.0DO K=1,MH=1!计算Y=F(x)Y(1)=F1(X(1),X(2),X(3),X(4))Y(2)=F2(X(1),X(2),X(3),X(4))Y(3)=F3(X(1),X(2),X(3),X(4))Y(4)=F4(X(1),X(2),X(3),X(4))!计算JA(N,N)E=0.0DO I=1,NDO J=1,NDO JJ=1,NIF(JJ==J) THENE(JJ)=X(JJ)+H(JJ)ELSEE(JJ)=X(JJ)ENDIFENDDOIF(I==1) THENJA(I,J)=(F1(E(1),E(2),E(3),E(4))-F1(X(1),X(2),X(3),X(4)))/H(J) ELSEIF(I==2) THENJA(I,J)=(F2(E(1),E(2),E(3),E(4))-F2(X(1),X(2),X(3),X(4)))/H(J) ELSEIF(I==3) THENJA(I,J)=(F3(E(1),E(2),E(3),E(4))-F3(X(1),X(2),X(3),X(4)))/H(J) ELSEIF(I==4) THENJA(I,J)=(F4(E(1),E(2),E(3),E(4))-F4(X(1),X(2),X(3),X(4)))/H(J) ENDIFENDDOENDDO!求解线性方程组CALL GAUSS(JA,XK,-Y,N)!判断精度CALL NORM(XK,N,E1)CALL NORM(X,N,E2)IF(E1/E2<=EPS) THENT=X(1)U=X(2)EXITELSEDO I=1,NX(I)=X(I)+XK(I)ENDDOENDIFENDDORETURNEND!**********列主元高斯消去法求解线性方程组********* SUBROUTINE GAUSS(A,X,B,N)DIMENSION A(N,N),B(N),X(N),T(N,N),TB(N)REAL M(N,N)REAL(8) A,B,X,T!消元过程DO K=1,N-1TA=A(K,K)TL=KDO L=K+1,NIF ((A(L,K)>TA).OR.(A(L,K)==TA)) THENTA=A(L,K)TL=LDO J=K,NT(K,J)=A(K,J)A(K,J)=A(TL,J)A(TL,J)=T(K,J)ENDDOTB(K)=B(K)B(K)=B(TL)B(TL)=TB(K)ENDIFENDDODO I=K+1,NM(I,K)=A(I,K)/A(K,K)A(I,K)=0DO J=K+1,NA(I,J)=A(I,J)-M(I,K)*A(K,J)ENDDOB(I)=B(I)-M(I,K)*B(K)ENDDOENDDO!回代过程X(N)=B(N)/A(N,N)DO K=N-1,1,-1S=0.0DO J=K+1,NS=S+A(K,J)*X(J)ENDDOX(K)=(B(K)-S)/A(K,K)ENDDORETURNEND!***********求向量的无穷数************ SUBROUTINE NORM(X,N,A)DIMENSION X(N)REAL(8) X,AA=ABS(X(1))DO I=2,NIF(ABS(X(I))>ABS(X(I-1))) THENA=ABS(X(I))ENDIFENDDORETURNEND!**************分片二次代数插值************** SUBROUTINE INTERPOLATION(U,V,W) PARAMETER (N=6,M=6)DIMENSION X(N),Y(M),Z(M,N),LK(3),LR(3)REAL(8) X,Y,Z,H,TREAL(8) U,V,W,LK,LR !U,V分别为插值点处的坐标,W为插值结果INTEGER R!**********************数据赋值********************** DATA Y/0.0,0.2,0.4,0.6,0.8,1.0/DATA X/0.0,0.4,0.8,1.2,1.6,2.0/DATA Z/-0.5,-0.42,-0.18,0.22,0.78,1.5,&&-0.34,-0.5,-0.5,-0.34,-0.02,0.46,&&0.14,-0.26,-0.5,-0.58,-0.5,-0.26,&&0.94,0.3,-0.18,-0.5,-0.66,-0.66,&&2.06,1.18,0.46,-0.1,-0.5,-0.74,&&3.5,2.38,1.42,0.62,-0.02,-0.5/H=0.4T=0.2!******************计算K,R*************************IF(U<=X(2)+H/2) THENK=2ELSEIF(U>X(N-1)-H/2) THENK=N-1ELSEDO I=3,N-2IF((U>X(I)-H/2).AND.(U<=X(I)+H/2)) THENK=IENDIFENDDOENDIFIF(V<=Y(2)+T/2) THENR=2ELSEIF(V>Y(M-1)-T/2) THENR=M-1ELSEDO J=3,M-2IF((V>Y(J)-T/2).AND.(V<=Y(J)+T/2)) THENR=JENDIFENDDOENDIFI=KJ=RLK(1)=(U-X(I))*(U-X(I+1))/(X(I-1)-X(I))/(X(I-1)-X(I+1))LK(2)=(U-X(I-1))*(U-X(I+1))/(X(I)-X(I-1))/(X(I)-X(I+1)) LK(3)=(U-X(I))*(U-X(I-1))/(X(I+1)-X(I))/(X(I+1)-X(I-1)) LR(1)=(V-Y(J))*(V-Y(J+1))/(Y(J-1)-Y(J))/(Y(J-1)-Y(J+1)) LR(2)=(V-Y(J-1))*(V-Y(J+1))/(Y(J)-Y(J-1))/(Y(J)-Y(J+1)) LR(3)=(V-Y(J))*(V-Y(J-1))/(Y(J+1)-Y(J))/(Y(J+1)-Y(J-1))W=0DO K=1,3DO R=1,3W=W+LK(K)*LR(R)*Z(J+R-2,I+K-2)ENDDOENDDORETURNEND!*******************最小二乘拟合子函数************** SUBROUTINE LSFITTING(X,Y,Z,A,N,M,P,Q,DT1)INTEGER P,QDIMENSION X(N),Y(M),Z(N,M),A(P,Q)DIMENSION APX(20),APY(20),BX(20),BY(20),U(20,20),V(20,M) DIMENSION T(20),T1(20),T2(20)REAL(8) X,Y,Z,A,DT1DO I=1,PDO J=1,QA(I,J)=0.0ENDDOENDDOIF(P>N) P=NIF(P>20) P=20IF(Q>M) Q=MIF(Q>20) Q=20XX=0YY=0D1=NAPX(1)=0.0DO I=1,NAPX(1)=APX(1)+X(I)ENDDOAPX(1)=APX(1)/D1DO J=1,MV(1,J)=0.0DO I=1,NV(1,J)=V(1,J)+Z(I,J)ENDDOV(1,J)=V(1,J)/D1ENDDOIF(P>1) THEND2=0.0APX(2)=0.0DO I=1,NG=X(I)-APX(1)D2=D2+G*GAPX(2)=APX(2)+(X(I)-XX)*G*G ENDDOAPX(2)=APX(2)/D2BX(2)=D2/D1DO J=1,MV(2,J)=0.0DO I=1,NG=X(I)-APX(1)V(2,J)=V(2,J)+Z(I,J)*GENDDOV(2,J)=V(2,J)/D2ENDDOD1=D2ENDIFDO K=3,PD2=0.0APX(K)=0.0DO J=1,MV(K,J)=0.0ENDDODO I=1,NG1=1.0G2=X(I)-APX(1)DO J=3,KG=(X(I)-APX(J-1))*G2-BX(J-1)*G1 G1=G2G2=GENDDOD2=D2+G*GAPX(K)=APX(K)+X(I)*G*GDO J=1,MV(K,J)=V(K,J)+Z(I,J)*GENDDOENDDODO J=1,MV(K,J)=V(K,J)/D2ENDDOAPX(K)=APX(K)/D2BX(K)=D2/D1D1=D2ENDDOD1=MAPY(1)=0.0DO I=1,MAPY(1)=APY(1)+Y(I)ENDDOAPY(1)=APY(1)/D1DO J=1,PU(J,1)=0.0DO I=1,MU(J,1)=U(J,1)+V(J,I) ENDDOU(J,1)=U(J,1)/D1ENDDOIF(Q>1)THEND2=0.0APY(2)=0.0DO I=1,MG=Y(I)-APY(1)D2=D2+G*GAPY(2)=APY(2)+(Y(I))*G*G ENDDOAPY(2)=APY(2)/D2BY(2)=D2/D1DO J=1,PU(J,2)=0.0DO I=1,MG=Y(I)-APY(1)U(J,2)=U(J,2)+V(J,I)*GENDDOU(J,2)=U(J,2)/D2ENDDOD1=D2ENDIFDO K=3,QD2=0.0APY(K)=0.0DO J=1,PU(J,K)=0.0ENDDODO I=1,MG1=1.0G2=Y(I)-APY(1)DO J=3,KG=(Y(I)-APY(J-1))*G2-BY(J-1)*G1 G1=G2G2=GENDDOD2=D2+G*GAPY(K)=APY(K)+Y(I)*G*GDO J=1,PU(J,K)=U(J,K)+V(J,I)*GENDDOENDDODO J=1,PU(J,K)=U(J,K)/D2ENDDOAPY(K)=APY(K)/D2BY(K)=D2/D1D1=D2ENDDOV(1,1)=1.0V(2,1)=-APY(1)V(2,2)=1.0DO I=1,PDO J=1,QA(I,J)=0.0ENDDOENDDODO I=3,QV(I,I)=V(I-1,I-1)V(I,I-1)=-APY(I-1)*V(I-1,I-1)+V(I-1,I-2)IF(I>=4) THENDO K=I-2,2,-1V(I,K)=-APY(I-1)*V(I-1,K)+V(I-1,K-1)-BY(I-1)*V(I-2,K) ENDDOENDIFV(I,1)=-APY(I-1)*V(I-1,1)-BY(I-1)*V(I-2,1)ENDDODO I=1,PIF(I==1) THENT(1)=1.0T1(1)=1.0ELSEIF(I==2) THENT(1)=-APX(1)T(2)=1.0T2(1)=T(1)T2(2)=T(2)ELSET(I)=T2(I-1)T(I-1)=-APX(I-1)*T2(I-1)+T2(I-2)IF(I>=4) THENDO K=I-2,2,-1T(K)=-APX(I-1)*T2(K)+T2(K-1)-BX(I-1)*T1(K)ENDDOENDIFT(1)=-APX(I-1)*T2(1)-BX(I-1)*T1(1)T2(I)=T(I)DO K=I-1,1,-1T1(K)=T2(K)T2(K)=T(K)ENDDOENDIFDO J=1,QDO K=I,1,-1DO L=J,1,-1A(K,L)=A(K,L)+U(I,J)*T(K)*V(J,L)ENDDOENDDOENDDOENDDODT1=0.0DO I=1,NX1=X(I)DO J=1,MY1=Y(J)X2=1.0DD=0.0DO K=1,PG=A(K,Q)DO KK=Q-1,1,-1G=G*Y1+A(K,KK)ENDDOG=G*X2DD=DD+GX2=X2*X1ENDDODT=DD-Z(I,J)DT1=DT1+DT*DTENDDOENDDORETURNEND三、计算结果数表(x,y,f(x,y)):X Y UX TY F(X,Y) 0.00 0.500 1.345 0.243 0.17E+000.00 0.550 1.322 0.269 0.66E+000.00 0.600 1.299 0.295 0.35E+000.00 0.650 1.277 0.322 0.94E+000.00 0.700 1.255 0.350 0.30E-020.00 0.750 1.235 0.377 -0.87E-010.00 0.800 1.215 0.406 -0.58E+000.00 0.850 1.196 0.434 -0.72E+000.00 0.900 1.177 0.463 -0.54E+000.00 0.950 1.159 0.492 -0.86E+000.00 1.050 1.125 0.550 -0.74E+00 0.00 1.100 1.109 0.580 -0.06E+00 0.00 1.150 1.093 0.609 -0.00E+00 0.00 1.200 1.079 0.639 -0.18E+00 0.00 1.250 1.064 0.669 -0.52E+00 0.00 1.300 1.050 0.699 -0.19E+00 0.00 1.350 1.037 0.729 -0.48E+00 0.00 1.400 1.024 0.759 -0.68E+00 0.00 1.450 1.011 0.790 -0.52E+00 0.00 1.500 1.000 0.820 -0.29E+000.08 0.500 1.415 0.228 0.67E+00 0.08 0.550 1.391 0.253 0.08E+00 0.08 0.600 1.368 0.279 0.02E+00 0.08 0.650 1.346 0.306 0.47E+00 0.08 0.700 1.325 0.333 0.57E+00 0.08 0.750 1.304 0.360 0.48E-01 0.08 0.800 1.284 0.388 -0.73E-01 0.08 0.850 1.265 0.416 -0.16E+00 0.08 0.900 1.246 0.444 -0.29E+00 0.08 0.950 1.229 0.473 -0.36E+00 0.08 1.000 1.211 0.502 -0.08E+00 0.08 1.050 1.194 0.531 -0.29E+00 0.08 1.100 1.178 0.560 -0.78E+00 0.08 1.150 1.163 0.589 -0.93E+00 0.08 1.200 1.148 0.619 -0.44E+00 0.08 1.250 1.133 0.649 -0.92E+00 0.08 1.300 1.119 0.679 -0.71E+000.08 1.400 1.093 0.739 -0.37E+00 0.08 1.450 1.080 0.769 -0.83E+00 0.08 1.500 1.068 0.799 -0.92E+000.16 0.500 1.483 0.214 0.31E+00 0.16 0.550 1.460 0.239 0.64E+00 0.16 0.600 1.437 0.264 0.91E+00 0.16 0.650 1.414 0.290 0.06E+00 0.16 0.700 1.393 0.316 0.70E+00 0.16 0.750 1.372 0.343 0.59E+00 0.16 0.800 1.352 0.370 0.12E+00 0.16 0.850 1.333 0.398 0.77E-02 0.16 0.900 1.315 0.426 -0.83E-01 0.16 0.950 1.297 0.454 -0.58E+00 0.16 1.000 1.279 0.483 -0.20E+00 0.16 1.050 1.262 0.512 -0.11E+00 0.16 1.100 1.246 0.541 -0.74E+00 0.16 1.150 1.231 0.570 -0.09E+00 0.16 1.200 1.216 0.600 -0.59E+00 0.16 1.250 1.201 0.629 -0.66E+00 0.16 1.300 1.187 0.659 -0.71E+00 0.16 1.350 1.174 0.689 -0.32E+00 0.16 1.400 1.161 0.718 -0.56E+00 0.16 1.450 1.148 0.748 -0.31E+00 0.16 1.500 1.136 0.778 -0.75E+000.24 0.500 1.551 0.201 0.66E+01 0.24 0.550 1.527 0.225 0.03E+000.24 0.650 1.482 0.275 0.64E+00 0.24 0.700 1.460 0.301 0.47E+00 0.24 0.750 1.439 0.327 0.34E+00 0.24 0.800 1.419 0.354 0.24E+00 0.24 0.850 1.400 0.381 0.69E+00 0.24 0.900 1.381 0.409 0.04E-01 0.24 0.950 1.363 0.437 -0.42E-01 0.24 1.000 1.346 0.465 -0.06E+00 0.24 1.050 1.329 0.494 -0.59E+00 0.24 1.100 1.313 0.523 -0.83E+00 0.24 1.150 1.297 0.552 -0.15E+00 0.24 1.200 1.282 0.581 -0.19E+00 0.24 1.250 1.267 0.610 -0.84E+00 0.24 1.300 1.253 0.640 -0.66E+00 0.24 1.350 1.240 0.669 -0.30E+00 0.24 1.400 1.227 0.699 -0.86E+00 0.24 1.450 1.214 0.729 -0.84E+00 0.24 1.500 1.202 0.759 -0.77E+000.32 0.500 1.617 0.188 0.28E+01 0.32 0.550 1.593 0.212 0.49E+01 0.32 0.600 1.570 0.236 0.68E+00 0.32 0.650 1.547 0.261 0.75E+00 0.32 0.700 1.526 0.286 0.60E+00 0.32 0.750 1.505 0.312 0.77E+00 0.32 0.800 1.485 0.339 0.05E+00 0.32 0.850 1.466 0.365 0.99E+00 0.32 0.900 1.447 0.393 0.27E+000.32 1.000 1.411 0.448 -0.01E-02 0.32 1.050 1.395 0.477 -0.41E-01 0.32 1.100 1.378 0.505 -0.18E+00 0.32 1.150 1.363 0.534 -0.25E+00 0.32 1.200 1.347 0.563 -0.29E+00 0.32 1.250 1.333 0.592 -0.90E+00 0.32 1.300 1.319 0.621 -0.00E+00 0.32 1.350 1.305 0.650 -0.40E+00 0.32 1.400 1.292 0.680 -0.54E+00 0.32 1.450 1.279 0.710 -0.79E+00 0.32 1.500 1.267 0.739 -0.91E+000.40 0.500 1.681 0.177 0.91E+01 0.40 0.550 1.658 0.199 0.00E+01 0.40 0.600 1.634 0.223 0.83E+01 0.40 0.650 1.612 0.247 0.02E+01 0.40 0.700 1.591 0.272 0.94E+00 0.40 0.750 1.570 0.298 0.49E+00 0.40 0.800 1.550 0.324 0.94E+00 0.40 0.850 1.530 0.350 0.40E+00 0.40 0.900 1.512 0.377 0.33E+00 0.40 0.950 1.493 0.405 0.99E+00 0.40 1.000 1.476 0.432 0.68E+00 0.40 1.050 1.459 0.460 0.08E-01 0.40 1.100 1.443 0.488 -0.84E-01 0.40 1.150 1.427 0.517 -0.98E+00 0.40 1.200 1.412 0.545 -0.27E+00 0.40 1.250 1.397 0.574 -0.06E+000.40 1.350 1.369 0.632 -0.66E+00 0.40 1.400 1.356 0.662 -0.37E+00 0.40 1.450 1.343 0.691 -0.43E+00 0.40 1.500 1.331 0.721 -0.12E+000.48 0.500 1.745 0.166 0.69E+01 0.48 0.550 1.721 0.188 0.02E+01 0.48 0.600 1.698 0.211 0.74E+01 0.48 0.650 1.676 0.235 0.40E+01 0.48 0.700 1.654 0.259 0.23E+01 0.48 0.750 1.634 0.284 0.56E+00 0.48 0.800 1.613 0.310 0.28E+00 0.48 0.850 1.594 0.336 0.49E+00 0.48 0.900 1.575 0.363 0.31E+00 0.48 0.950 1.557 0.390 0.66E+00 0.48 1.000 1.539 0.417 0.30E+00 0.48 1.050 1.522 0.444 0.34E+00 0.48 1.100 1.506 0.472 0.07E-01 0.48 1.150 1.490 0.500 -0.62E-01 0.48 1.200 1.475 0.529 -0.45E+00 0.48 1.250 1.460 0.557 -0.86E+00 0.48 1.300 1.446 0.586 -0.39E+00 0.48 1.350 1.432 0.615 -0.22E+00 0.48 1.400 1.419 0.644 -0.67E+00 0.48 1.450 1.406 0.674 -0.55E+00 0.48 1.500 1.394 0.703 -0.14E+000.56 0.500 1.808 0.156 0.48E+010.56 0.600 1.761 0.200 0.10E+01 0.56 0.650 1.739 0.223 0.68E+01 0.56 0.700 1.717 0.247 0.94E+01 0.56 0.750 1.696 0.272 0.33E+01 0.56 0.800 1.676 0.297 0.11E+00 0.56 0.850 1.657 0.323 0.63E+00 0.56 0.900 1.638 0.349 0.97E+00 0.56 0.950 1.620 0.375 0.52E+00 0.56 1.000 1.602 0.402 0.56E+00 0.56 1.050 1.585 0.429 0.47E+00 0.56 1.100 1.568 0.457 0.20E+00 0.56 1.150 1.552 0.485 0.13E+00 0.56 1.200 1.537 0.513 0.09E-01 0.56 1.250 1.522 0.541 -0.47E-01 0.56 1.300 1.508 0.570 -0.99E+00 0.56 1.350 1.494 0.599 -0.82E+00 0.56 1.400 1.481 0.627 -0.26E+00 0.56 1.450 1.468 0.657 -0.71E+00 0.56 1.500 1.455 0.686 -0.98E+000.64 0.500 1.870 0.147 0.74E+01 0.64 0.550 1.846 0.168 0.10E+01 0.64 0.600 1.823 0.190 0.54E+01 0.64 0.650 1.801 0.213 0.42E+01 0.64 0.700 1.779 0.236 0.56E+01 0.64 0.750 1.758 0.260 0.03E+01 0.64 0.800 1.738 0.285 0.42E+01 0.64 0.850 1.718 0.310 0.41E+010.64 0.950 1.681 0.362 0.36E+00 0.64 1.000 1.664 0.388 0.18E+00 0.64 1.050 1.646 0.415 0.28E+00 0.64 1.100 1.630 0.443 0.07E+00 0.64 1.150 1.614 0.470 0.66E+00 0.64 1.200 1.598 0.498 0.09E+00 0.64 1.250 1.584 0.526 0.50E-01 0.64 1.300 1.569 0.554 -0.88E-01 0.64 1.350 1.555 0.583 -0.76E+00 0.64 1.400 1.542 0.611 -0.66E+00 0.64 1.450 1.529 0.640 -0.33E+00 0.64 1.500 1.516 0.669 -0.56E+000.72 0.500 1.931 0.139 0.94E+01 0.72 0.550 1.907 0.159 0.84E+01 0.72 0.600 1.884 0.181 0.36E+01 0.72 0.650 1.862 0.203 0.40E+01 0.72 0.700 1.840 0.226 0.47E+01 0.72 0.750 1.819 0.249 0.56E+01 0.72 0.800 1.799 0.273 0.19E+01 0.72 0.850 1.779 0.298 0.37E+01 0.72 0.900 1.760 0.323 0.86E+01 0.72 0.950 1.742 0.349 0.76E+00 0.72 1.000 1.724 0.375 0.24E+00 0.72 1.050 1.707 0.402 0.55E+00 0.72 1.100 1.691 0.429 0.97E+00 0.72 1.150 1.675 0.456 0.27E+00 0.72 1.200 1.659 0.484 0.31E+000.72 1.250 1.644 0.511 0.51E+00 0.72 1.300 1.630 0.539 0.49E+00 0.72 1.350 1.616 0.568 0.72E-02 0.72 1.400 1.602 0.596 -0.69E-01 0.72 1.450 1.589 0.625 -0.67E+00 0.72 1.500 1.576 0.653 -0.20E+000.80 0.500 1.992 0.131 0.31E+01 0.80 0.550 1.968 0.151 0.44E+01 0.80 0.600 1.945 0.172 0.41E+01 0.80 0.650 1.922 0.193 0.45E+01 0.80 0.700 1.900 0.216 0.00E+01 0.80 0.750 1.879 0.239 0.10E+01 0.80 0.800 1.859 0.263 0.16E+01 0.80 0.850 1.840 0.287 0.52E+01 0.80 0.900 1.821 0.312 0.02E+01 0.80 0.950 1.802 0.337 0.38E+01 0.80 1.000 1.784 0.363 0.89E+01 0.80 1.050 1.767 0.389 0.28E+00 0.80 1.100 1.751 0.416 0.09E+00 0.80 1.150 1.734 0.443 0.23E+00 0.80 1.200 1.719 0.470 0.93E+00 0.80 1.250 1.704 0.498 0.15E+00 0.80 1.300 1.689 0.525 0.86E+00 0.80 1.350 1.675 0.553 0.64E+00 0.80 1.400 1.662 0.582 0.74E-01 0.80 1.450 1.649 0.610 -0.37E-01 0.80 1.500 1.636 0.638 -0.81E+00K和σ分别为:0 0.93E+031 0.61E+012 0.92E-023 0.53E-034 0.16E-055 0.77E-07系数矩阵Crs(按行)为:0.00E+01 -0.83E+01 0.56E+00 0.97E+00 -0.03E+00 0.70E-010.91E+01 -0.99E+00 -0.96E+01 0.17E+01 -0.66E+00 0.10E-010.77E+00 0.42E+01 -0.10E+00 -0.81E+00 0.81E+00 -0.62E-01-0.25E+00 -0.21E+00 0.97E+00 -0.18E+00 0.49E+00 -0.63E-010.34E+00 -0.56E+00 0.69E-01 0.51E+00 -0.77E-01 0.27E-01-0.94E-01 0.94E+00 -0.58E+00 0.69E-01 -0.50E-01 0.53E-02数表(x,y,f(x,y),p(x,y)):X Y F(X,Y) P(X,Y)0.100 0.700 0.58E+00 0.05E+000.100 1.100 -0.66E+00 -0.26E+00 0.100 1.300 -0.68E+00 -0.31E+00 0.100 1.500 -0.52E+00 -0.49E+000.200 0.700 0.54E+00 0.19E+00 0.200 0.900 -0.63E-01 -0.65E-01 0.200 1.100 -0.90E+00 -0.90E+00 0.200 1.300 -0.84E+00 -0.90E+00 0.200 1.500 -0.03E+00 -0.04E+000.300 0.700 0.82E+00 0.09E+00 0.300 0.900 0.48E+00 0.11E+00 0.300 1.100 -0.63E+00 -0.88E+00 0.300 1.300 -0.72E+00 -0.96E+00 0.300 1.500 -0.34E+00 -0.84E+000.400 0.700 0.79E+00 0.89E+00 0.400 0.900 0.56E+00 0.63E+00 0.400 1.100 -0.83E-01 -0.04E-01 0.400 1.300 -0.72E+00 -0.71E+00 0.400 1.500 -0.85E+00 -0.07E+000.500 0.700 0.56E+01 0.92E+01 0.500 0.900 0.51E+00 0.23E+00 0.500 1.100 0.59E+00 0.27E+00 0.500 1.300 -0.53E+00 -0.11E+00 0.500 1.500 -0.67E+00 -0.33E+000.600 0.900 0.14E+00 0.75E+00 0.600 1.100 0.19E+00 0.32E+00 0.600 1.300 -0.70E-01 -0.82E-01 0.600 1.500 -0.08E+00 -0.75E+000.700 0.700 0.89E+01 0.29E+01 0.700 0.900 0.91E+01 0.11E+01 0.700 1.100 0.60E+00 0.97E+00 0.700 1.300 0.22E-01 0.06E-01 0.700 1.500 -0.53E+00 -0.80E+000.800 0.700 0.09E+01 0.06E+01 0.800 0.900 0.32E+01 0.50E+01 0.800 1.100 0.03E+00 0.79E+00 0.800 1.300 0.25E+00 0.50E+00 0.800 1.500 -0.14E+00 -0.28E+00。
BUAA数值分析大作业三
北京航空航天大学2020届研究生《数值分析》实验作业第九题院系:xx学院学号:姓名:2020年11月Q9:方程组A.4一、 算法设计方案(一)总体思路1.题目要求∑∑===k i kj s r rsy x cy x p 00),(对f(x, y) 进行拟合,可选用乘积型最小二乘拟合。
),(i i y x 与),(i i y x f 的数表由方程组与表A-1得到。
2.),(**j i y x f 与1使用相同方法求得,),(**j i y x p 由计算得出的p(x,y)直接带入),(**j i y x 求得。
1. ),(i i y x 与),(i i y x f 的数表的获得对区域D ={ (x,y)|1≤x ≤1.24,1.0≤y ≤1.16}上的f (x , y )值可通过xi=1+0.008i ,yj=1+0.008j ,得到),(i i y x 共31×21组。
将每组带入A4方程组,即可获得五个二元函数组,通过简单牛顿迭代法求解这五个二元数组可获得z1~z5有关x,y 的表达式。
再将),(i i y x 分别带入z1~z5表达式即可获得f(x,y)值。
2.乘积型最小二乘曲面拟合2.1使用乘积型最小二乘拟合,根据k 值不用,有基函数矩阵如下:⎪⎪⎪⎭⎫ ⎝⎛=k i i k x x x x B 0000 , ⎪⎪⎪⎭⎫ ⎝⎛=k j jk y y y y G 0000数表矩阵如下:⎪⎪⎪⎭⎫⎝⎛=),(),(),(),(0000j i i j y x f y x f y x f y x f U记C=[rs c ],则系数rs c 的表达式矩阵为:11-)(-=G G UG B B B C T TT )(通过求解如下线性方程,即可得到系数矩阵C 。
UG B G G C B B T T T =)()(2.2计算),(),,(****j i j i y x p y x f (i =1,2,…,31 ; j =1,2,…,21) 的值),(**j i y x f 的计算与),(j i y x f 相同。
北航数值分析大作业3(学硕)
《数值分析》作业三院系:机械学院学号:SY1307145姓名:龙安林2013年11 月24 日1. 算法设计1) 开始;2) 计算数组[][]0.08,0.050.5,0,1,2,,10;0,1,2,,20x i i y j j i j ==+=⋯=⋯(); 3) 将点[][],0,1,2,,10;0,1,2,,20x i y j i j =⋯=⋯(),()带入非线性方程组: 0.5cos 2.670.5sin 1.070.5cos 3.740.5sin 0.79t u v w x t u v w y t u v w x t u v w y +++-=⎧⎪+++-=⎪⎨+++-=⎪⎪+++-=⎩ 得出相应的点,t u (); 4) 选择拉格朗日插值法,将,t u ()作为中间变量,在题目所给出的二维数表中进行二次代数插值,得到[][],)(z f x i y j =;5) 输出数表:[][][][]()()0,1,2,,10;0,1,2,,20,,,x i y j f x i y j i j =⋯=⋯; 6) 令k=0;7) 以()()(),,,0,1,r r r s x x y y r s ϕψ===…,k 为拟合基函数,将上述数表作为拟合条件,对于给定的k 值,得到矩阵B 、G 、U ;8) 令-1-1(),()T T T A B B B U C AG G G ==,用选主元的LU 分解法分别计算矩阵A 和C 的各列,最后得到系数矩阵C ;9) 以公式:()()()00,k ki j rs r i s j s r p x y C x y ϕψ===∑∑计算每个点的拟合值;10) 利用公式:()()()2102000,,i j i j i j f x y p x y σ===-∑∑计算拟合误差,当σ≤10-7时,循环结束,否则k=k+1,转(6);11) 令[][]()**0.10.50.2 1,2,81,2,5x i i y j j i j ==+=⋯=⋯;,;,;12) 计算()()()******,,,,,i j i j i jf x y p x y delta x y ,输出数表,观察逼近效果; 13) 结束。
北航数值分析计算实习第一题编程
i − t + s +1,t t − k + s +1, k t = max(1,i − r ,k − s )
∑c
c
) / cs +1, k
[i = k + 1, k + 2,⋯ , min( k + r , n); k < n]
(2) 求解 Ly = b,Ux = y (数组 b 先是存放原方程右端向量,后来存放中间向量 y)
0 b a2
b c
c b a3 b c
⋯ ⋯ ⋯ ⋯ ⋯
c b a499 b c
c b a500 b 0
c ⎤ b ⎥ ⎥ a501 ⎥ ⎥ 0 ⎥ 0 ⎥ ⎦
在数组 C 中检索矩阵 A 的带内元素 aij 的方法是: A 的带内元素 aij =C 中的元素 ci − j + s +1, j
2
数值分析计算实习题目一
i −1
bi := bi −
பைடு நூலகம்
i − t + s +1,t t t = max(1,i − r )
∑c
b
(i = 2,3,⋯ , n)
xn := bn / cs +1, n
min( i + s )
xi := (bi −
t = i +1
∑c
i −t + s +1,t t
x ) / cs +1,i
(i = n − 1, n − 2,⋯ ,1)
3、Doolittle 分解求解 n 元带状线性方程组(doolittle()函数)
按照上述对带状矩阵 A 的存储方法和元素 aij 的检索方法,并且把三角分解的结果 ukj 和 lik 分 别存放在 akj 和 aik 原先的存储单元内,那么用 Doolittle 分解法求解 n 元带状线性方程组的算法 可重新表述如下(其中“:=”表示赋值) : (1) 作分解 A = LU 。 对于 k=1,2, ……,n 执行
北航研究生数值分析大作业三
数值分析—计算实习作业三学院:17系专业:精密仪器及机械姓名:张大军学号:DY1417114一、程序设计方案程序设计方案流程图如图1所示。
(注:由本人独立完成,并且有几处算法很巧妙,但同时也有许多不足,可以优化和模块化,由于时间原因只实现了调试通过)图1.程序设计方案流程图二、程序源代码#include <iostream.h>#include <iomanip.h>#include <math.h>#include<stdio.h>#include <conio.h>#define M 10000#define N 4#define E 1.0e-12int zuixiaci;static double c[9][9];static double bijin[8][5];int main(){double X[N]={0,0,0,1};double T[11][21],U[11][21],xianshi[11][21];double diertX[N];double F[N];double f[N][N];double Max1=0,Max2=0;int k,i,j,t,tt=0,yao=0;void qiuF(double * X,double *F,int i,int j);void qiuF2(double *X,double *F,int i,int j);void qiuf(double * X,double (*f)[N]);void qiudiertX(double (*a)[N],double*b,double*X); double gouzaohs(double t,double u); void solve_C(double (*U)[21]); void pp(double (*U)[21],int k);for(i=0;i<11;i++)for(j=0;j<21;j++){for(k=0;k<M;k++){qiuF(X,F,i,j);qiuf(X,f);qiudiertX(f,F,diertX);for(t=0;t<N;t++){X[t]=X[t]+diertX[t];}Max1=0,Max2=0;for(t=0;t<N;t++){if(Max1<fabs(X[t]))Max1=fabs(X[t]);if(Max2<fabs(diertX[t]))Max2=fabs(diertX[t]);}if((Max2/Max1)<=E){k=M;yao=1;T[i][j]=X[0];U[i][j]=X[1];xianshi[i][j]=gouzaohs(X[0],X[1]);cout<<setiosflags(ios::scientific)<<setprecision(12);cout<<setprecision(2)<<"("<<setw(5)<<0.08*i<<","<<setw(5)<< (0.5+0.05*j)<<",";cout<<setprecision(12)<<setw(21)<<xianshi[i][j]<<") ";if(tt==3){tt=0;cout<<'\n';cout<<'\n';}else{tt++;}}}if(yao==0)cout<<"迭代不成功"<<endl; yao=0;}cout<<endl;solve_C(xianshi);pp(xianshi,zuixiaci);tt=0;for(i=1;i<9;i++)for(j=1;j<6;j++){for(k=0;k<M;k++){qiuF2(X,F,i,j);qiuf(X,f);qiudiertX(f,F,diertX);for(t=0;t<N;t++){X[t]=X[t]+diertX[t];}Max1=0,Max2=0;for(t=0;t<N;t++){if(Max1<fabs(X[t]))Max1=fabs(X[t]);if(Max2<fabs(diertX[t]))Max2=fabs(diertX[t]);}if((Max2/Max1)<=E){k=M;yao=1;xianshi[i-1][j-1]=gouzaohs(X[0],X[1]);cout<<setiosflags(ios::scientific)<<setprecision(12);cout<<setprecision(2)<<"("<<setw(5)<<0.1*i<<","<<setw(5)<<( 0.5+0.2*j)<<",";cout<<setprecision(12)<<setw(21)<<xianshi[i-1][j-1]<<","<<set w(21)<<bijin[i-1][j-1]<<") ";if(tt==2){tt=0;cout<<'\n';}else{tt++;}}}if(yao==0)cout<<"迭代不成功"<<endl;yao=0;}cout<<endl;return 1;}void qiuF(double *X,double *F,int i,int j){F[0]=-(0.5*cos(X[0])+X[1]+X[2]+X[3]-0.08*i-2.67);F[1]=-(X[0]+0.5*sin(X[1])+X[2]+X[3]-(0.5+0.05*j)-1.07);F[2]=-(0.5*X[0]+X[1]+cos(X[2])+X[3]-0.08*i-3.74);F[3]=-(X[0]+0.5*X[1]+X[2]+sin(X[3])-(0.5+0.05*j)-0.79); }void qiuF2(double *X,double *F,int i,int j){F[0]=-(0.5*cos(X[0])+X[1]+X[2]+X[3]-0.1*i-2.67);F[1]=-(X[0]+0.5*sin(X[1])+X[2]+X[3]-(0.5+0.2*j)-1.07);F[2]=-(0.5*X[0]+X[1]+cos(X[2])+X[3]-0.1*i-3.74);F[3]=-(X[0]+0.5*X[1]+X[2]+sin(X[3])-(0.5+0.2*j)-0.79); }void qiuf(double *X,double (*f)[N]){f[0][0]=-0.5*sin(X[0]);f[0][1]=1;f[0][2]=1;f[0][3]=1;f[1][0]=1;f[1][1]=0.5*cos(X[1]);f[1][2]=1;f[1][3]=1;f[2][0]=0.5;f[2][1]=1;f[2][2]=-sin(X[2]);f[2][3]=1;f[3][0]=1;f[3][1]=0.5;f[3][2]=1;f[3][3]=cos(X[3]);}//求解关于变化X的线性方程组void qiudiertX(double (*a)[N],double*b,double*X) {double H[N][N]={0},l[N]={0};double B;double sum;int i,j,m,k,z;for(k=0;k<N-1;k++){for(j=k;j<N;j++){l[j]=a[k][j];}z=k;for(m=k;m<N;m++){if(fabs(a[z][k])<fabs(a[m][k]))z=m;}for(j=k;j<N;j++){a[k][j]=a[z][j];a[z][j]=l[j];}B=b[k];b[k]=b[z];b[z]=B;for(i=k+1;i<N;i++){H[i][k]=a[i][k]/a[k][k];for(j=k+1;j<N;j++)a[i][j]=a[i][j]-H[i][k]*a[k][j];b[i]=b[i]-H[i][k]*b[k];}}if(a[N-1][N-1]==0){cout<<"算法失效,停止计算"<<endl; }else{X[N-1]=b[N-1]/a[N-1][N-1];for(k=N-2;k>=0;k--){sum=0;for(j=k+1;j<N;j++){sum=sum+a[k][j]*X[j];}X[k]=(b[k]-sum)/a[k][k];}}}//作二元差值,使用分片二次代数插值double gouzaohs(double t,double u){double T[6]={0,0.2,0.4,0.6,0.8,1},U[6]={0,0.4,0.8,1.2,1.6,2};double Z[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};double g=0,sum=0,sum1=1,sum2=1;int i=0,j=0,k=0,r=0,kk=0,rr=0;for(i=1;(i<6)&&(T[i]-0.1<t);i++){}for(j=1;(j<6)&&(U[j]-0.2<u);j++){}if(i==1)i=2;if(i==6)i=5;if(j==1)j=2;if(j==6)j=5;sum=0;for(k=i-2;k<i+1;k++)for(r=j-2;r<j+1;r++){sum1=1;sum2=1;for(kk=i-2;kk<i+1;kk++){if(k!=kk){sum1=sum1*(t-T[kk])/(T[k]-T[kk]);}}for(rr=j-2;rr<j+1;rr++){if(r!=rr){sum2=sum2*(u-U[rr])/(U[r]-U[rr]);}}sum=sum+sum1*sum2*Z[k][r];}g=sum;return g;}//求r*s阶矩阵A与s*t阶矩阵B的乘积矩阵Cvoid Multi(double *a, double *b, double *c, int la, int lb, int lc, int r, int s, int t){int i, j, k;for (i=0; i<r; i++)for (j=0; j<t; j++){*(c+i*lc+j)=0;for (k=0; k<s; k++)*(c+i*lc+j)+=*(a+i*la+k)*(*(b+k*lb+j));}}//求n阶方阵A的逆矩阵Bdouble Inverse(double *a, double *b, int la, int lb, int n){int i, j, k;double temp;for(i=0; i<n; i++)for(j=0; j<n; j++)if (i==j)*(b+i*lb+j)=1;else*(b+i*lb+j)=0;for (k=0; k<n; k++){j=k;for (i=k+1; i<n; i++)if (fabs(*(a+i*la+k))>fabs(*(a+j*la+k))) j=i;if (j!=k)for (i=0; i<n; i++){temp=*(a+j*la+i);*(a+j*la+i)=*(a+k*la+i);*(a+k*la+i)=temp;temp=*(b+j*lb+i);*(b+j*lb+i)=*(b+k*lb+i);*(b+k*lb+i)=temp;}if (*(a+k*la+k)==0)return 0;if ((temp=*(a+k*la+k))!=1)for (i=0; i<n; i++){*(a+k*la+i)/=temp;*(b+k*lb+i)/=temp;}for (i=0; i<n; i++)if ((*(a+i*la+k)!=0) && (i!=k)){temp=*(a+i*la+k);for (j=0; j<n; j++){*(a+i*la+j)-=temp*(*(a+k*la+j));*(b+i*lb+j)-=temp*(*(b+k*lb+j));}}}return 0;}void solve_C(double (*U)[21]){int i,j,r,s,k;double t1[21][21], t2[21][21], t3[21][21],d[9][9],e[9][9];double B[11][9], B_T[9][11], G[21][9], G_T[9][21],P[11][21];double temp, FangCha;for(i=0;i<9;i++){for(j=0;j<11;j++){B[j][i]=pow(0.08*j,i);B_T[i][j]=pow(0.08*j,i);}for(j=0;j<21;j++){G[j][i]=pow(0.5+0.05*j,i);G_T[i][j]=pow(0.5+0.05*j,i);}}for (k=0; k<9; k++){FangCha=0;Multi(B_T[0], B[0], t1[0], 11, 9, 21, k+1, 11, k+1);Inverse(t1[0], c[0], 21, 9, k+1);Multi(e[0], c[0], d[0], 9, 9, 9, k+1, k+1, k+1);Multi(c[0], B_T[0], t1[0], 9, 11, 21, k+1, k+1, 11);Multi(t1[0], U[0], t2[0], 21, 21, 21, k+1, 11, 21);Multi(G_T[0], G[0], t1[0], 21, 9, 21, k+1, 21, k+1);Inverse(t1[0], c[0], 21, 9, k+1);Multi(G[0], c[0], t3[0], 9, 9, 21, 21, k+1, k+1);Multi(t2[0], t3[0], c[0], 21, 21, 9, k+1, 21, k+1);for(i=0;i<11;i++)for(j=0;j<21;j++){temp=0;for(r=0;r<k+1;r++)for(s=0;s<k+1;s++)temp+=c[r][s]*B[i][r]*G[j][s];P[i][j]=temp;FangCha+=(U[i][j]-temp)*(U[i][j]-temp);}cout<<"k="<<setw(5)<<k<<";"<<setw(5)<<"Sigma="<<FangCha<<" ;\n"<<'\n';if(FangCha<=1.0e-7){zuixiaci=k;cout<<"达到精度要求时: k="<<setw(5)<<k<<";"<<setw(5)<<"Sigma="<<FangCha<<";\n";cout<<" 系数c(r,s)如下:\n";for(i=0;i<k+1;i++){for(j=0;j<k+1;j++){cout<<"C("<<i<<","<<j<<")="<<setw(21)<<c[i][j]<<"; ";}cout<<endl<<'\n';}cout<<endl;return;}}cout<<"经过8次拟合没有达到所需精度;"<<endl;//最高可拟合10次return;}void pp(double (*U)[21],int k){int i,j,r,s;double B[8][9],G[5][9],temp;for(i=0;i<k+1;i++){for(j=0;j<8;j++){B[j][i]=pow(0.1*(j+1),i);}for(j=0;j<5;j++){G[j][i]=pow(0.5+0.2*(j+1),i);}}for(i=0;i<8;i++)for(j=0;j<5;j++){temp=0;for(r=0;r<k+1;r++)for(s=0;s<k+1;s++)temp+=c[r][s]*B[i][r]*G[j][s];bijin[i][j]=temp;}}三、程序运行结果显示程序运行结果显示如图2。
北航数分大作业三
一、算法的设计方案1、对于已给出的非线性方程组,其解集可采用牛顿迭代法进行求解。
在每次迭代过程中,将x ,y 的值固定,如此便可得到一组关于t ,u ,v ,w 的解。
因此可以建立一组(x ,y )和(t ,u )一一对应的关系。
2、采用分片二次插值对题目中所给出的z ,t ,u 二维数表进行处理。
于是在 0≤t ≤1, 0≤u ≤2 的矩形区域就建立了 z 与(t,u)的一一对应关系。
其中选择(m ,n )满足,2322m i m h h t t t m -<≤+≤≤,,2322n j n u u u n ττ-<≤+≤≤。
3、对i x i *08.0= 10,...,2,1,0=i ,j y j 05.05.0+= 20,...,2,1,0=j 。
分别使用前两步算法,可得到一组2111)(],[⨯=j i y x f j i z 的数表。
4、采用最小二乘拟合,设∑∑===k r k s s r rs y x cy x p 00)(,m=10 n=20,M=N=K 。
插值基函数10,...,1,0,)(==i x x r i i r φ k r ,...,1,0=,ks j y y sj j s ,...,1,020,...,1,0,)(===ψ。
、)1()1()1()1(][][+⨯++⨯+==k n sj k m r i y G x B U 即为上面所求的Z[11][21]。
为避免计算过程中出现矩阵求逆,将U B B B A T T 1)(-=改为U B A B B T T =)(,再利用高斯消去法以)(B B T作为系数矩阵,U B T 的每一列作为非线性部分,分别解出A 的每一列。
在将1)(-=G G AG A T 改为AG G G C T =)(,然后利用高斯消元法以)(G G T 作为系数矩阵,AG 的每一行作为非线性部分,分解出C 的每一行。
如此便得到了最小二乘拟合的系数矩阵C 。
北航数值分析B自测题1-3章
D.2
10.若f (x)在区间[a, b]内单调有限,二分k 次后区间记为[ak , bk ],且每次 取xk+1 = ak 近似代替精确解x∗ ,则最小的绝对误差限为( ) A.|xk+1 − x∗ | ≤ C.|xk+1 − x∗ | ≤ 二、 填空题 1.设π 的近似数π ∗ 有4位有效数字,则其相对误差限为 √ 2. x∗ 的相对误差约是x∗ 的相对误差的 倍。
(2)a=6,b=2; (3)a=2,b=3; (4)a=-1,b=2.
3.设矩阵A ∈ Rn∗n , Q ∈ Rn∗n ,且QT Q=E,则下列关系式不成立的是( ) (1) A
2
= AQ ; (2) QA
2 F
= A
F
; (3) Qx
2
= x , 其中x ∈ Rn ;
2
(4)cond∞ ( ) = cond∞ ( Q).
8
3 −1 4 1 4.设矩阵A=−1 2 −2 ,x=−1 ,则 Ax 2 −3 −2 1 ( ) (1)8,8; (2)8,7; (3)8,6; (4)7,7.
A. 7 3
B.− 7 3
7 C. 6
D.− 7 6
8.追赶法适用于求解( )线性方程组. A.上三角 B.下三角 C.对角 D.三对角 二.填空 题:
2 1 0 1.设A=1 2 a,为使A 可分解为A=LLT ,其中L是 对角元素为正的下 0 a 2 三角矩阵 范围是 ,则a的取值 .
1 (s∗ ) s∗
≈
27 110×80
= 0.31% = max 1 ≤ j ≤ n{1, 2, 4} = 4;
北航数值分析大作业3
数值分析第三次作业1. 设计方案对Fredholm 积分方程,用迭代法进行求解:()'(())u x A u x =,其中11(())()(,)()A u x g x K x y u y dy -=-⋅⎰对于公式中的积分部分用数值积分方法。
复化梯形积分法,取2601个节点,取迭代次数上限为50次。
实际计算迭代次数为18次,最后算得误差为r= 0.97E-10。
复化Simpson 积分法,取迭代次数上限为50次,取2*41+1,即83个节点时能满足精度要求。
实际计算迭代次数为17次,最后的误差为 r= 0.97E-10。
Guass 积分法选择的Gauss —Legendre 法,取迭代次数上限为50次,直接选择8个节点,满足精度要求。
实际计算迭代次数为24次,最后算得误差为r= 0.87E-10。
2. 全部源程序 module integral implicit none contains!//////////复化梯形 subroutine trapezoid(m) implicit none integer :: i,j,k,mreal*8 :: x(m+1),u(m+1) real*8 :: sum,sum1,g,r,h real*8 :: e=1.0e-10h=2./m do i=1,m+1x(i)=-1.+(i-1)*h end dou=0.02 do k=1,50 do i=1,m+1 sum1=0.g=dexp(x(i)*4.)+(dexp(x(i)+4.)-dexp(-4.-x(i)))/(x(i)+4.) do j=2,m sum1=sum1+dexp(x(i)*x(j))*u(j) end do sum=h/2.*(dexp(x(i)*-1.)*u(1)+dexp(x(i)*1.)*u(m+1)+2*sum1) u(i)=g-sum end dor=h/2.*((dexp(x(1)*4)-u(1))**2+(dexp(x(m+1)*4)-u(m+1))**2) do i=2,mr=r+h*(dexp(x(i)*4)-u(i))**2end doif(dabs(r)<=e) exitend dowrite(*,*) kopen(1,file="trapezoid.txt")do i=1,m+1write(1,'(3(f18.12))') x(i),u(i),dexp(x(i)*4.)end dowrite(1,'(4x,a2,e9.2)') "r=",rclose(1)returnend subroutine trapezoid!///////////复化simpsonsubroutine simpson(m)implicit noneinteger :: i,j,k,mreal*8 :: x(2*m+1),u(2*m+1)real*8 :: sum,sum1,sum2,g,r,hreal*8 :: e=1.0e-10h=2./(2.*m)do i=1,2*m+1x(i)=-1.+(i-1)*hend dou=0.02do k=1,50do i=1,2*m+1sum1=0.sum2=0.g=dexp(x(i)*4.)+(dexp(x(i)+4.)-dexp(-4.-x(i)))/(x(i)+4.)do j=1,msum1=sum1+dexp(x(i)*x(2*j))*u(2*j)end dodo j=1,m-1sum2=sum2+dexp(x(i)*x(2*j+1))*u(2*j+1)sum=h/3.*(dexp(x(i)*-1.)*u(1)+dexp(x(i)*1.)*u(2*m+1)+4*sum1+2*sum2) u(i)=g-sumend dor=h/3.*((dexp(x(1)*4)-u(1))**2+(dexp(x(2*m+1)*4)-u(2*m+1))**2)do i=1,mr=r+4.*h/3.*(dexp(x(2*i)*4)-u(2*i))**2end dodo i=1,m-1r=r+2.*h/3.*(dexp(x(2*i+1)*4)-u(2*i+1))**2end doif(dabs(r)<=e) exitend dowrite(*,*) kopen(2,file="simpson.txt")do i=1,2*m+1write(2,'(3(f18.12))') x(i),u(i),dexp(x(i)*4.)end dowrite(2,'(4x,a2,e9.2)') "r=",rclose(2)returnend subroutine simpson!///////////Gauss_Legendre法subroutine Gaussimplicit noneinteger,parameter :: m=8integer :: i,j,kreal*8 :: x(m),u(m),a(m)real*8 :: sum,g,rreal*8 :: e=1.0e-10data x /-0.9602898565,-0.7966664774,-0.5255324099,-0.1834346425,&0.1834346425,0.5255324099,0.7966664774,0.9602898565/data a /0.1012285363,0.2223810345,0.3137066459,0.3626837834,&0.3626837834,0.3137066459,0.2223810345,0.1012285363/u=0.02do k=1,50do i=1,mg=dexp(x(i)*4.)+(dexp(x(i)+4.)-dexp(-4.-x(i)))/(x(i)+4.)do j=1,msum=sum+dexp(x(i)*x(j))*u(j)*a(j)end dou(i)=g-sumend dor=0.do i=1,mr=r+a(i)*(dexp(x(i)*4)-u(i))**2end doif(dabs(r)<=e) exitend dowrite(*,*) kopen(3,file="Gauss.txt")do i=1,mwrite(3,'(3(f18.12))') x(i),u(i),dexp(x(i)*4.)end dowrite(3,'(4x,a2,e9.2)') "r=",rclose(3)returnend subroutine Gaussend module!//////////主程序program mainuse integralimplicit noneinteger :: code1=2600integer :: code2=41call trapezoid(code1)call simpson(code2)call Gaussend program3.各种积分方法的节点和数值解(由于数据太多,在打印时用了较计算时少的有效数字)复化Simpson法4.各方法所得曲线(由于所取节点太多,且精度高,所以图中很难看出各曲线的区别。
北航数值分析大作业题目三
1、 算法的设计方案: (一)、总体方案设计:
(1)解非线性方程组。将给定的当作已知量代入题目给定的非线性方程组,求得与相对应 的数组t[i][j],u[i][j]。 (2)分片二次代数插值。通过分片二次代数插值运算,得到与数组t[11][21],u[11][21]] 对应的数组z[11][21],得到二元函数z=。 (3)曲面拟合。利用x[i],y[j],z[11][21]建立二维函数表,再根据精度的要求选择适当k 值,并得到曲面拟合的系数矩阵C[r][s]。 (4)观察和的逼近效果。观察逼近效果只需要重复上面(1)和(2)的过程,得到与新的 插值节点对应的,再与对应的比较即可,这里求解可以直接使用(3)中的C[r][s]和k。
{ temp=0; for(l=k+1;l<=3;l++) {temp=temp+dF[k][l]*dx[l]/dF[k][k];} dx[k]=-F[k]/dF[k][k]-temp; } temp=0; for(l=0;l<=3;l++) /*求解矩阵范数,用无穷范数*/ { if(temp<fabs(dx[l])) temp=fabs(dx[l]); } fx=temp; temp=0; for(l=0;l<=3;l++) { if(temp<fabs(X[l])) temp=fabs(X[l]); } fX=temp; if(fabs(fx/fX)<Epsilon1) /*判断是否成立*/ { t[i][j]=X[0]; u[i][j]=X[1]; goto loop4;} else { for(l=0;l<=3;l++) {X[l]=X[l]+dx[l];} n=n+1; goto loop3;} } loop3:{if(n<M) /*判断是否超出规定迭代次数*/ goto loop1; else printf("迭代不成功\n"); goto loop4; } loop4:{continue;} } } } void fpeccz(double t[11][21],double u[11][21])/*分片二次代数插值子程序*/ { int s[11][21],r[11][21]; int i,j,i1,j1,m; double z0[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},
北航数值分析B第一次上机作业算法作业
数值分析第一次上机作业 一、算法方案设计 (1)存储矩阵A (参考课本25页):矩阵A 501×501为大型带状矩阵,上半带宽S=2,下半带宽R=2,参照课本,可将其用循环语句存储在一个5行501列的二维数组A[5][501]中,使得矩阵A 的第j 列存放数组A 的第j 列带状元素,并使得矩阵的主对角元素存放在数组的第三行。
检索矩阵的元素时只要按照公式:a ij =数组a[i-j+3][j]即可。
(2)求λ1,λ501,λs :第一步,用幂法求出矩阵A 按模最大特征值,再对其判断正负,如果是正数,则该特征值为λ501,如果是负数,则该特征值为λ1。
幂法实现过程具体为:第二步,用反幂法求矩阵A 按模最小特征向量λS 。
班级:ZY16131 学号:ZY16131姓名:第三步,由于λ1≤λ2≤…..≤λ501,可采用原点平移法对矩阵A平移λ1,得到矩阵(-λ1I+A),记为矩阵B,再用用幂法求出B的模最大特征值λB,则λ501=λB+λ1。
(3)距离μk最近的特征值:还是用原点平移法,将矩阵A平移μk个单位,再用反幂法求出平移后矩阵模最小特征值ηk,矩阵A与μk最接近的特征值为:λi k = ηk + μk =ηk +(4)A的谱范数条件数与detA:a、求矩阵A的行列式值:在用反幂法时需要对矩阵A进行Doolittle三角分解,A=LU,根据det(A)=det (LU)=det(L)*det(U),其中det(L)=1,det(A)即为U矩阵的对角线元素的乘积。
b、求A的(谱范数)条件数cond(A)2:由于A是非奇异实对称阵,从而cond(A)2 =∣λ1∣/∣λs∣。
二、运行结果分析(1)取初始向量为 u_0[501]={1,1…1},计算结果截图如下:(2)讨论初始迭代向量取值对计算结果的影响在编程实现算法的过程中遇到了很多问题,多次尝试才得以解决。
例如,对各个函数的定义后需要对矩阵A进行初始化;为简化程序,方便改变变量值,应该尽可能地将某些多级变量写成函数的形式,只需要对初级变量赋值即可。
- 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。
数值分析B大作业第三题学生姓名学生学号SY0905xxxx专业名称飞行器设计院(系)名称航空科学与工程学院2012年12月数值分析大作业第三题一、编程语言选择及方案说明:本文选择Visual C++语言进行编程,针对复化梯形积分法、复化simpson 积分法和高斯法分别编写了三个函数trapezoid()、simpson()以及Guass()。
解方程组的方法为高斯列主元法,即本程序中的Guass_solve()。
最后在同一坐标系中使用Matlab 将保存在“result.txt ”文件中的本程序计算结果绘图。
二、算法设计方案:1. 化原方程为方程组(k=0,1,…...)for i=0,……,k将原积分方程化为:-0()()kx xji j i i j i u x e e u g x λ=+=∑(2.1)ix =-1+2i/k ,则由上式得到k+1个线性方程,联立可得:(I )U k k A b +=,其中: U=(u 0 u 1……u k )T, b k [i]=g(i x ) , A k [i][ j]=x x i jj e λ。
2. 确定方程组A k /b k三种不同方法,将有不同的λi ,χi ,进而方程组中A k ,b k 也不尽相同。
对于复化梯形积分法,取k+1个等距节点,则x i =-1+2i/k, ,i=0,…….,k;λ0=λk =1/k, 其余λi =2/k 。
对于复化simpson 积分法,取2k+1个等距节点,则x i =-1+i/k, i=0,1,……, 2k ; λ0=λ2k =1/(3k) , λ2i-1=4/(3k), i=1,…….,k; λ2i =2/(3k),i=1,…….,k-1。
对于高斯方法,选择ρ(x)≡1,即Guass_legendre 求积公式,对应的x i 和λi 参见教材p170表6-6,进而求解方程组求得u i 。
(仅可用于节点数k<10)求解方程组采用Guass 列主元消去法,由于(I )A +主对角元素严格大于同列其他元素,故其本身即是列主元,无需进行变换,Guass 消去法即是Guass 列主元法。
3. 判断误差for i =0,1,……,k ,u(x i )=e (4*xi)2=(())ki i i u x u σ=-∑,While 1-10e σ≤ 结束,输出计算结果并绘图。
三、源程序:#include<iostream> #include<fstream.h> #include<iomanip.h> #include<stdlib.h> #include<math.h> #include<time.h> #define Err 1.0e-10 #define MAX 3000double u[MAX]={0},x[MAX]={0},coe[MAX]={0};//积分数值解及插值系数向量 double A[MAX][MAX]={0},B[MAX]={0}; //方程组 double err=0; //误差 int trapezoid(int); //复化梯形法 int simpson(int); //复化simpson 法 int Guass(int); //Guass_legendre 法 void Guass_solve(int); //Guass 列主元法解方程组void print(int); //按要求打印结果//Gauss_legender法求积系数及节点表double X[9][9]={0,0,0,0,0,0,0,0,0,-0.5773502692,0.5773502692,0,0,0,0,0,0,0,-0.7745966692,0,0.7745966692,0,0,0,0,0,0,-0.8611363116,-0.3399810436,0.3399810436,0.8611363116,0,0,0,0,0,-0.9061798459,-5384693101,0,5384693101,0.9061798459,0,0,0,0,-0.9324695142,-0.6612093865,-0.2386191861,0.2386191861,0.6612093865,0.932469542,0,0,0,-0.9491079123,-0.7415311856,-0.4058451514,0,0.4058451514,0.7415311856,0.949107912,0,0,-0.9602898565,-0.7966664774,-0.5255324099,-0.1834346425,0.1834346425,0.5255324099,0.7966664774,0.9602898565,0,-0.9681602395,-0.8360311073,-0.6133714327,-0.3242534234,0,0.3242534 234,0.6133714327,0.8360311073,0.9681602395};double C[9][9]={2,0,0,0,0,0,0,0,0,1,1,0,0,0,0,0,0,0,0.5555555556,0.8888888889,0.5555555556,0,0,0,0,0,0,0.3478548451,0.6521451549,0.6521451549,0.3478548451,0,0,0,0,0,0.2369268851,0.4786286705,0.5688888889,0.4786286705,0.2369268851,0,0,0,0,0.1713244924,0.3607615730,0.4679139346,0.4679139346,0.3607615730,0.1713244924,0,0,0,0.1294849662,0.2797053915,0.3818300505,0.4179591837,0.3818300505,0.2797053915,0.1294849662,0, 0,0.1012285363,0.2223810345,0.3137066459,0.3626837834,0.3626837834,0.3137066459,0.2223810345, 0.1012285363,0,0.0812743884,0.1806481607,0.2606106964,0.3123470770,0.3302393550,0.3123470770,0.2606106964,0.1806481607,0.0812743884};using namespace std;int main()//主函数{ int flag=0,k;/*==============复化梯形法===============*/for(k=1;k<MAX;k++){ flag=trapezoid(k);if(flag>0){::cout<<endl<<"用复化梯形法计算,节点个数k="<<flag+1<<"时到达要求精度,结果见'result.txt'!";print(flag);break;}}/*===========复化simpson法==============*/flag=0;for(k=1;k<MAX;k++){ flag=simpson(k);if(flag){ ::cout<<endl<<"复化simpson法节点个数k="<<2*flag+1<<"时到达要求精度,结果见'result.txt'!";print(flag);break;}}/*==============Guass法================*/flag=0;for(k=1;k<MAX;k++){ flag=Guass(k);if(flag>0){ ::cout<<endl<<"用Guass法计算,节点个数k="<<flag<<"时达到要求精度,结果见'result.txt'!";print(flag-1);break;}}return 0;}//========复化梯形解法========================int trapezoid(int k){ double h=2.0/k;err=0;for(int i=0;i<k+1;i++) {coe[i]=h;x[i]=-1+h*i;} coe[0]=coe[k]=h/2;Guass_solve(k);for(i=0;i<k+1;i++)err+=coe[i]*pow((exp(4.0*x[i])-u[i]),2);if(err<=Err) return k; else return 0;}//========复化simpson法=================int simpson(int k){ double h=1.0/k; err=0;for(int i=0;i<2*k+1;i++) x[i]=-1+h*i;for(i=0;i<k+1;i++) {coe[2*i]=2.0*h/3;coe[2*i+1]=4.0*h/3;}coe[0]=coe[2*k]=h/3;Guass_solve(2*k);for(i=0;i<2*k+1;i++) err+=coe[i]*pow((exp(4.0*x[i])-u[i]),2);if(err<=Err) return k; else return 0;}//=====================================int Guass(int k){ for(int i=0;i<k;i++) {x[i]=X[k-1][i];coe[i]=C[k-1][i];}Guass_solve(k-1);//Jacobi_solve(k-1);for(i=0,err=0;i<k;i++) err+=coe[i]*pow((exp(4*x[i])-u[i]),2);if(err<=Err) return k; else return 0;}//=========================解方程组void Guass_solve(int n)//高斯列主元法解方程{for(int i=0;i<n+1;i++){ for(int j=0;j<n+1;j++) A[i][j]= coe[j]*exp(x[i]*x[j]);A[i][i]+=1; B[i]=exp(4.0*x[i])+(exp(x[i]+4)-exp(-x[i]-4))/(x[i]+4);}//第一步:消元 //主对角元素远远大于其他,即主对角已经是主元,不必再选for(int k=0;k<n;k++){ for(int i=k+1;i<n+1;i++){ double lik=A[i][k]/A[k][k];for(int h=k;h<n+1;h++) A[i][h]-=lik*A[k][h]; B[i]-=lik*B[k];}}//第二步:回代========================================double sum=0;u[n]=B[n]/A[n][n];for(i=n-1;i>=0;i--){ sum=0;for(int h=i+1;h<n+1;h++) sum+=A[i][h]*u[h];u[i]=(B[i]-sum)/A[i][i];}}//=========================打印结果void print(int k){::ofstream out;out.open("result.txt",ios::app);for(int i=0;i<k+1;i++)::cout<<endl<<setiosflags(ios::scientific)<<setprecision(12)<<"x"<<i<<"="<<x[i]<<",u"<<i<<"="<<u[i];out<<"\n===============================================\n";::cout<<"\n===============================================\n";out.close();}四、运行结果1.复化梯形法结果:x0=-1.000000000000e+000,u0=1.831678781434e-002 x1=-9.992310649750e-001,u1=1.837320828901e-002 x2=-9.984621299500e-001,u2=1.842980256461e-002 x3=-9.976931949250e-001,u3=1.848657117658e-002x4=-9.969242599000e-001,u4=1.854351466189e-002 x5=-9.961553248750e-001,u5=1.860063355930e-002 x6=-9.953863898501e-001,u6=1.865792840912e-002 x7=-9.946174548251e-001,u7=1.871539975341e-002 x8=-9.938485198001e-001,u8=1.877304813585e-002 x9=-9.930795847751e-001,u9=1.883087410179e-002 x10=-9.923106497501e-001,u10=1.888887819828e-002 x11=-9.915417147251e-001,u11=1.894706097407e-002 x12=-9.907727797001e-001,u12=1.900542297958e-002 x13=-9.900038446751e-001,u13=1.906396476690e-002 x14=-9.892349096501e-001,u14=1.912268688990e-002 x15=-9.884659746251e-001,u15=1.918158990407e-002 x16=-9.876970396002e-001,u16=1.924067436663e-002 x17=-9.869281045752e-001,u17=1.929994083654e-002 x18=-9.861591695502e-001,u18=1.935938987450e-002 x19=-9.853*********e-001,u19=1.941902204291e-002 x20=-9.846212995002e-001,u20=1.947883790589e-002 x21=-9.838523644752e-001,u21=1.953883802930e-002 x22=-9.830834294502e-001,u22=1.959902298075e-002 x23=-9.823144944252e-001,u23=1.965939332964e-002 x24=-9.815455594002e-001,u24=1.971994964706e-002 x25=-9.807766243752e-001,u25=1.978069250590e-002 x26=-9.800076893502e-001,u26=1.984162248077e-002 x27=-9.792387543253e-001,u27=1.990274014812e-002 x28=-9.784698193003e-001,u28=1.996404608610e-002 x29=-9.777008842753e-001,u29=2.002554087470e-002 x30=-9.769319492503e-001,u30=2.008722509574e-002 x31=-9.761630142253e-001,u31=2.014909933262e-002 x32=-9.753940792003e-001,u32=2.021*********e-002 x33=-9.746251441753e-001,u33=2.027*********e-002 x34=-9.738562091503e-001,u34=2.0335********e-002x35=-9.730872741253e-001,u35=2.0398********e-002x36=-9.723183391003e-001,u36=2.046134130623e-002 x37=-9.715494040754e-001,u37=2.052436799418e-002 x38=-9.707804690504e-001,u38=2.058758883352e-002 x39=-9.700115340254e-001,u39=2.065100442234e-002x40=-9.692425990004e-001,u40=2.0714********e-002 x41=-9.684736639754e-001,u41=2.0778********e-002 x42=-9.677047289504e-001,u42=2.084242569416e-002 x43=-9.669357939254e-001,u43=2.0906********e-002 x44=-9.661668589004e-001,u44=2.0971********e-002 x45=-9.653979238754e-001,u45=2.103562141979e-002 x46=-9.646289888504e-001,u46=2.110041715675e-002x47=-9.638600538255e-001,u47=2.116541249465e-002x48=-9.630911188005e-001,u48=2.123060804837e-002x49=-9.623221837755e-001,u49=2.129600443468e-002 x50=-9.615532487505e-001,u50=2.136160227223e-002 x51=-9.607843137255e-001,u51=2.142740218159e-002 x52=-9.600153787005e-001,u52=2.149340478526e-002 x53=-9.592464436755e-001,u53=2.155961070760e-002 x54=-9.584775086505e-001,u54=2.162602057500e-002 x55=-9.577085736255e-001,u55=2.169263501564e-002 x56=-9.569396386005e-001,u56=2.175945465976e-002 x57=-9.561707035755e-001,u57=2.182648013949e-002 x58=-9.554017685506e-001,u58=2.189371208887e-002 x59=-9.546328335256e-001,u59=2.196115114397e-002x60=-9.538638985006e-001,u60=2.202879794272e-002 x61=-9.530949634756e-001,u61=2.209665312512e-002 x62=-9.523260284506e-001,u62=2.216471733311e-002x63=-9.515570934256e-001,u63=2.223299121053e-002 x64=-9.507881584006e-001,u64=2.230147540333e-002 x65=-9.500192233756e-001,u65=2.237017055935e-002 x66=-9.492502883506e-001,u66=2.243907732850e-002 x67=-9.484813533256e-001,u67=2.250819636262e-002 x68=-9.477124183007e-001,u68=2.257752831560e-002 x69=-9.469434832757e-001,u69=2.264707384334e-002 x70=-9.461745482507e-001,u70=2.271683360372e-002 x71=-9.454056132257e-001,u71=2.278680825675e-002 x72=-9.446366782007e-001,u72=2.285699846439e-002 x73=-9.438677431757e-001,u73=2.292740489060e-002 x74=-9.430988081507e-001,u74=2.299802820153e-002 x75=-9.423298731257e-001,u75=2.306886906519e-002 x76=-9.415609381007e-001,u76=2.313992815181e-002 x77=-9.407920030757e-001,u77=2.321120613364e-002x78=-9.400230680507e-001,u78=2.328270368494e-002 x79=-9.392541330258e-001,u79=2.335442148215e-002 x80=-9.384851980008e-001,u80=2.342636020370e-002 x81=-9.377162629758e-001,u81=2.349852053011e-002x82=-9.369473279508e-001,u82=2.357090314412e-002 x83=-9.361783929258e-001,u83=2.364350873045e-002 x84=-9.354094579008e-001,u84=2.371633797594e-002 x85=-9.346405228758e-001,u85=2.378939156958e-002 x86=-9.338715878508e-001,u86=2.386267020250e-002 x87=-9.331026528258e-001,u87=2.393617456795e-002 x88=-9.323337178008e-001,u88=2.400990536121e-002 x89=-9.315647827759e-001,u89=2.408386327986e-002 x90=-9.307958477509e-001,u90=2.415804902355e-002 x91=-9.300269127259e-001,u91=2.423246329410e-002 x92=-9.292579777009e-001,u92=2.430710679546e-002 x93=-9.284890426759e-001,u93=2.438198023381e-002 x94=-9.277201076509e-001,u94=2.445708431742e-002 x95=-9.269511726259e-001,u95=2.453241975686e-002x96=-9.261822376009e-001,u96=2.460798726475e-002 x97=-9.254133025759e-001,u97=2.468378755603e-002 x98=-9.246443675509e-001,u98=2.475982134778e-002 x99=-9.238754325260e-001,u99=2.483608935929e-002 x100=-9.231064975010e-001,u100=2.491259231205e-002 x101=-9.223375624760e-001,u101=2.498933092983e-002 x102=-9.215686274510e-001,u102=2.506630593855e-002 x103=-9.207996924260e-001,u103=2.514351806649e-002 x104=-9.200307574010e-001,u104=2.522096804402e-002 x105=-9.192618223760e-001,u105=2.529865660385e-002 x106=-9.184928873510e-001,u106=2.537658448097e-002 x107=-9.177239523260e-001,u107=2.545475241257e-002 x108=-9.169550173010e-001,u108=2.553316113812e-002 x109=-9.161860822760e-001,u109=2.561181139942e-002 x110=-9.154171472511e-001,u110=2.569070394050e-002 x111=-9.146482122261e-001,u111=2.576983950770e-002 x112=-9.138792772011e-001,u112=2.584921884966e-002 x113=-9.131103421761e-001,u113=2.592884271733e-002 x114=-9.123414071511e-001,u114=2.600871186400e-002 x115=-9.115724721261e-001,u115=2.608882704522e-002 x116=-9.108035371011e-001,u116=2.616918901888e-002 x117=-9.100346020761e-001,u117=2.624979854527e-002 x118=-9.092656670511e-001,u118=2.633065638694e-002 x119=-9.084967320261e-001,u119=2.641176330883e-002 x120=-9.0772********e-001,u120=2.649312007825e-002 x121=-9.069588619762e-001,u121=2.657472746486e-002 x122=-9.061899269512e-001,u122=2.665658624063e-002 x123=-9.054209919262e-001,u123=2.673869718001e-002 x124=-9.046520569012e-001,u124=2.682106105979e-002 x125=-9.038831218762e-001,u125=2.690367865912e-002 x126=-9.0311********e-001,u126=2.698655075963e-002 x127=-9.023*********e-001,u127=2.706967814525e-002 x128=-9.015763168012e-001,u128=2.715306160243e-002 x129=-9.008073817762e-001,u129=2.723670192001e-002 x130=-9.000384467512e-001,u130=2.732059988919e-002 x131=-8.992695117263e-001,u131=2.740475630372e-002 x132=-8.985005767013e-001,u132=2.748917195968e-002 x133=-8.977316416763e-001,u133=2.757384765573e-002 x134=-8.969627066513e-001,u134=2.765878419286e-002 x135=-8.961937716263e-001,u135=2.774398237466e-002 x136=-8.954248366013e-001,u136=2.782944300708e-002 x137=-8.946559015763e-001,u137=2.791516689855e-002 x138=-8.938869665513e-001,u138=2.800115486016e-002 x139=-8.931180315263e-001,u139=2.808740770527e-002 x140=-8.923490965013e-001,u140=2.817392624990e-002 x141=-8.915801614764e-001,u141=2.826071131254e-002 x142=-8.908112264514e-001,u142=2.834776371420e-002 x143=-8.900422914264e-001,u143=2.843508427834e-002 x144=-8.892733564014e-001,u144=2.852*********e-002 x145=-8.885044213764e-001,u145=2.861053320113e-002 x146=-8.877354863514e-001,u146=2.869866321956e-002 x147=-8.869665513264e-001,u147=2.878706472009e-002 x148=-8.861976163014e-001,u148=2.887573853907e-002 x149=-8.854286812764e-001,u149=2.896468551536e-002 x150=-8.846597462514e-001,u150=2.905390649044e-002 x151=-8.838908112265e-001,u151=2.914340230832e-002 x152=-8.831218762015e-001,u152=2.923317381566e-002 x153=-8.823529411765e-001,u153=2.932322186177e-002 x154=-8.815840061515e-001,u154=2.941354729844e-002 x155=-8.808150711265e-001,u155=2.950415098023e-002 x156=-8.800461361015e-001,u156=2.959503376427e-002 x157=-8.792772010765e-001,u157=2.968619651034e-002 x158=-8.785082660515e-001,u158=2.977764008081e-002 x159=-8.777393310265e-001,u159=2.986936534083e-002 x160=-8.769703960015e-001,u160=2.996137315805e-002 x161=-8.762014609765e-001,u161=3.005366440297e-002 x162=-8.754325259516e-001,u162=3.014623994869e-002 x163=-8.746635909266e-001,u163=3.023*********e-002 x164=-8.738946559016e-001,u164=3.033224744819e-002 x165=-8.731257208766e-001,u165=3.042568116173e-002 x166=-8.723567858516e-001,u166=3.0519********e-002 x167=-8.715878508266e-001,u167=3.061341293576e-002 x168=-8.708189158016e-001,u168=3.070771277232e-002 x169=-8.700499807766e-001,u169=3.080230309709e-002 x170=-8.692810457516e-001,u170=3.089718480491e-002 x171=-8.685121107266e-001,u171=3.0992********e-002 x172=-8.677431757017e-001,u172=3.108782596298e-002 x173=-8.669742406767e-001,u173=3.118358721671e-002 x174=-8.662053056517e-001,u174=3.127964346054e-002 x175=-8.654363706267e-001,u175=3.137599560321e-002 x176=-8.646674356017e-001,u176=3.147264455625e-002 x177=-8.638985005767e-001,u177=3.156959123389e-002 x178=-8.631295655517e-001,u178=3.166683655335e-002 x179=-8.623606305267e-001,u179=3.176438143459e-002 x180=-8.615916955017e-001,u180=3.186222680037e-002 x181=-8.608227604767e-001,u181=3.196037357637e-002 x182=-8.600538254517e-001,u182=3.205882269105e-002 x183=-8.592848904268e-001,u183=3.215757507577e-002 x184=-8.585159554018e-001,u184=3.225663166479e-002 x185=-8.577470203768e-001,u185=3.235599339517e-002 x186=-8.569780853518e-001,u186=3.245566120688e-002 x187=-8.562091503268e-001,u187=3.255563604285e-002 x188=-8.554402153018e-001,u188=3.265591884884e-002 x189=-8.546712802768e-001,u189=3.275651057355e-002 x190=-8.539023452518e-001,u190=3.285741216858e-002 x191=-8.531334102268e-001,u191=3.295862458853e-002 x192=-8.523644752018e-001,u192=3.306014879085e-002 x193=-8.515955401769e-001,u193=3.316198573601e-002 x194=-8.508266051519e-001,u194=3.326413638740e-002 x195=-8.500576701269e-001,u195=3.336660171141e-002 x196=-8.492887351019e-001,u196=3.346938267734e-002 x197=-8.485198000769e-001,u197=3.357248025757e-002 x198=-8.477508650519e-001,u198=3.367589542744e-002 x199=-8.469819300269e-001,u199=3.377962916523e-002 x200=-8.462129950019e-001,u200=3.388368245231e-002 x201=-8.454440599769e-001,u201=3.398805627306e-002 x202=-8.446751249519e-001,u202=3.409275161488e-002 x203=-8.439061899270e-001,u203=3.419776946816e-002 x204=-8.431372549020e-001,u204=3.430311082649e-002 x205=-8.423683198770e-001,u205=3.440877668636e-002 x206=-8.415993848520e-001,u206=3.451476804737e-002 x207=-8.408304498270e-001,u207=3.462108591230e-002 x208=-8.400615148020e-001,u208=3.472773128689e-002 x209=-8.392925797770e-001,u209=3.483470518004e-002 x210=-8.385236447520e-001,u210=3.494200860379e-002 x211=-8.377547097270e-001,u211=3.504964257318e-002 x212=-8.369857747020e-001,u212=3.515760810651e-002 x213=-8.362168396770e-001,u213=3.526590622515e-002 x214=-8.354479046521e-001,u214=3.537453795362e-002 x215=-8.346789696271e-001,u215=3.548350431960e-002 x216=-8.339100346021e-001,u216=3.559280635395e-002 x217=-8.331410995771e-001,u217=3.570244509065e-002 x218=-8.323721645521e-001,u218=3.581242156696e-002 x219=-8.316032295271e-001,u219=3.592273682328e-002 x220=-8.308342945021e-001,u220=3.603339190320e-002 x221=-8.300653594771e-001,u221=3.614438785354e-002 x222=-8.292964244521e-001,u222=3.625572572438e-002 x223=-8.285274894271e-001,u223=3.636740656895e-002 x224=-8.277585544022e-001,u224=3.647943144384e-002 x225=-8.269896193772e-001,u225=3.659180140877e-002 x226=-8.262206843522e-001,u226=3.670451752684e-002 x227=-8.254517493272e-001,u227=3.681758086436e-002 x228=-8.246828143022e-001,u228=3.693099249092e-002 x229=-8.239138792772e-001,u229=3.704475347944e-002 x230=-8.231449442522e-001,u230=3.715886490608e-002 x231=-8.223760092272e-001,u231=3.727332785042e-002 x232=-8.216070742022e-001,u232=3.738814339530e-002 x233=-8.208381391772e-001,u233=3.750331262687e-002x234=-8.200692041522e-001,u234=3.761883663469e-002 x235=-8.193002691273e-001,u235=3.773471651157e-002 x236=-8.185313341023e-001,u236=3.785095335389e-002 x237=-8.177623990773e-001,u237=3.796754826118e-002 x238=-8.169934640523e-001,u238=3.808450233650e-002 x239=-8.162245290273e-001,u239=3.820181668626e-002 x240=-8.154555940023e-001,u240=3.831949242025e-002 x241=-8.146866589773e-001,u241=3.843753065172e-002 x242=-8.139177239523e-001,u242=3.855593249740e-002 x243=-8.131487889273e-001,u243=3.867469907736e-002 x244=-8.123798539023e-001,u244=3.879383151512e-002 x245=-8.116109188774e-001,u245=3.891333093773e-002 x246=-8.108419838524e-001,u246=3.903319847570e-002 x247=-8.100730488274e-001,u247=3.915343526300e-002 x248=-8.0930********e-001,u248=3.927404243711e-002 x249=-8.0853********e-001,u249=3.939502113896e-002 x250=-8.0776********e-001,u250=3.951637251307e-002 x251=-8.069973087274e-001,u251=3.963809770745e-002 x252=-8.062283737024e-001,u252=3.976019787368e-002 x253=-8.054594386774e-001,u253=3.988267416682e-002 x254=-8.046905036524e-001,u254=4.000552774549e-002 x255=-8.0392********e-001,u255=4.012875977197e-002 x256=-8.0315********e-001,u256=4.025*********e-002 x257=-8.023*********e-001,u257=4.0376********e-002 x258=-8.016147635525e-001,u258=4.050073821426e-002 x259=-8.008458285275e-001,u259=4.062549572595e-002 x260=-8.000768935025e-001,u260=4.075063755048e-002 x261=-7.993079584775e-001,u261=4.0876********e-002 x262=-7.985390234525e-001,u262=4.100207887722e-002 x263=-7.977700884275e-001,u263=4.112838075812e-002 x264=-7.970011534025e-001,u264=4.125507170927e-002 x265=-7.962322183775e-001,u265=4.138215292922e-002 x266=-7.954632833526e-001,u266=4.150962562017e-002 x267=-7.946943483276e-001,u267=4.163749098801e-002 x268=-7.939254133026e-001,u268=4.176575024247e-002 x269=-7.931564782776e-001,u269=4.189440459681e-002 x270=-7.923875432526e-001,u270=4.202345526820e-002 x271=-7.916186082276e-001,u271=4.215290347748e-002 x272=-7.908496732026e-001,u272=4.228275044923e-002 x273=-7.900807381776e-001,u273=4.241299741187e-002 x274=-7.893118031526e-001,u274=4.254364559756e-002 x275=-7.885428681276e-001,u275=4.267469624224e-002 x276=-7.877739331027e-001,u276=4.280615058570e-002 x277=-7.870049980777e-001,u277=4.293800987152e-002 x278=-7.862360630527e-001,u278=4.307027534713e-002 x279=-7.854671280277e-001,u279=4.320294826380e-002 x280=-7.846981930027e-001,u280=4.333602987664e-002 x281=-7.839292579777e-001,u281=4.346952144461e-002 x282=-7.831603229527e-001,u282=4.360342423062e-002 x283=-7.823913879277e-001,u283=4.373773950137e-002 x284=-7.816224529027e-001,u284=4.387246852758e-002 x285=-7.808535178777e-001,u285=4.400761258377e-002 x286=-7.800845828527e-001,u286=4.414317294844e-002 x287=-7.793156478278e-001,u287=4.427915090404e-002 x288=-7.785467128028e-001,u288=4.441554773693e-002 x289=-7.777777777778e-001,u289=4.455236473748e-002 x290=-7.770088427528e-001,u290=4.468960320002e-002 x291=-7.762399077278e-001,u291=4.482726442280e-002 x292=-7.754709727028e-001,u292=4.496534970824e-002 x293=-7.747020376778e-001,u293=4.510386036255e-002 x294=-7.739331026528e-001,u294=4.524279769613e-002 x295=-7.731641676278e-001,u295=4.538216302335e-002 x296=-7.723952326028e-001,u296=4.552195766261e-002 x297=-7.716262975779e-001,u297=4.566218293645e-002 x298=-7.708573625529e-001,u298=4.580284017142e-002 x299=-7.700884275279e-001,u299=4.594393069819e-002 x300=-7.693194925029e-001,u300=4.608545585144e-002 x301=-7.685505574779e-001,u301=4.622741697009e-002 x302=-7.677816224529e-001,u302=4.636981539711e-002 x303=-7.670126874279e-001,u303=4.651265247962e-002 x304=-7.662437524029e-001,u304=4.665592956895e-002 x305=-7.654748173779e-001,u305=4.679964802045e-002 x306=-7.647058823529e-001,u306=4.694380919380e-002 x307=-7.639369473280e-001,u307=4.708841445274e-002 x308=-7.631680123030e-001,u308=4.723346516530e-002 x309=-7.623990772780e-001,u309=4.737896270368e-002 x310=-7.616301422530e-001,u310=4.752490844438e-002 x311=-7.608612072280e-001,u311=4.767130376798e-002 x312=-7.600922722030e-001,u312=4.781815005953e-002 x313=-7.593233371780e-001,u313=4.796544870815e-002 x314=-7.585544021530e-001,u314=4.811320110735e-002 x315=-7.577854671280e-001,u315=4.826140865489e-002 x316=-7.570165321030e-001,u316=4.841007275288e-002 x317=-7.562475970780e-001,u317=4.855919480765e-002 x318=-7.554786620531e-001,u318=4.870877622999e-002 x319=-7.547097270281e-001,u319=4.885881843500e-002 x320=-7.539407920031e-001,u320=4.900932284208e-002 x321=-7.531718569781e-001,u321=4.916029087496e-002 x322=-7.524029219531e-001,u322=4.931172396200e-002 x323=-7.516339869281e-001,u323=4.946362353567e-002 x324=-7.508650519031e-001,u324=4.961599103299e-002 x325=-7.500961168781e-001,u325=4.976882789542e-002 x326=-7.493271818531e-001,u326=4.992213556887e-002 x327=-7.485582468281e-001,u327=5.007591550361e-002 x328=-7.477893118032e-001,u328=5.023*********e-002 x329=-7.470203767782e-001,u329=5.038489798068e-002 x330=-7.462514417532e-001,u330=5.054010344600e-002 x331=-7.454825067282e-001,u331=5.069578701884e-002 x332=-7.447135717032e-001,u332=5.0851********e-002 x333=-7.439446366782e-001,u333=5.100859438246e-002 x334=-7.431757016532e-001,u334=5.116572113253e-002 x335=-7.424067666282e-001,u335=5.132333190854e-002 x336=-7.416378316032e-001,u336=5.148142820151e-002 x337=-7.408688965782e-001,u337=5.164001150705e-002 x338=-7.400999615532e-001,u338=5.179908332544e-002 x339=-7.393310265283e-001,u339=5.195864516151e-002 x340=-7.385620915033e-001,u340=5.211869852473e-002 x341=-7.377931564783e-001,u341=5.227924492931e-002 x342=-7.370242214533e-001,u342=5.244028589402e-002 x343=-7.362552864283e-001,u343=5.260182294231e-002 x344=-7.354863514033e-001,u344=5.276385760240e-002 x345=-7.347174163783e-001,u345=5.292639140722e-002 x346=-7.339484813533e-001,u346=5.308942589429e-002 x347=-7.331795463283e-001,u347=5.325296260603e-002 x348=-7.324106113033e-001,u348=5.341700308946e-002 x349=-7.316416762784e-001,u349=5.358154889652e-002 x350=-7.308727412534e-001,u350=5.374660158380e-002 x351=-7.301038062284e-001,u351=5.391216271280e-002 x352=-7.293348712034e-001,u352=5.407823384966e-002 x353=-7.285659361784e-001,u353=5.424481656557e-002 x354=-7.277970011534e-001,u354=5.441191243636e-002 x355=-7.270280661284e-001,u355=5.457952304285e-002 x356=-7.262591311034e-001,u356=5.474764997066e-002 x357=-7.254901960784e-001,u357=5.491629481027e-002 x358=-7.247212610534e-001,u358=5.508545915715e-002 x359=-7.239523260285e-001,u359=5.525514461164e-002 x360=-7.231833910035e-001,u360=5.542535277899e-002 x361=-7.224144559785e-001,u361=5.559608526940e-002 x362=-7.216455209535e-001,u362=5.576734369804e-002 x363=-7.208765859285e-001,u363=5.593912968509e-002 x364=-7.201076509035e-001,u364=5.611144485561e-002 x365=-7.193387158785e-001,u365=5.628429083985e-002 x366=-7.185697808535e-001,u366=5.645766927291e-002 x367=-7.178008458285e-001,u367=5.663158179498e-002 x368=-7.170319108035e-001,u368=5.680603005136e-002 x369=-7.162629757785e-001,u369=5.698101569233e-002 x370=-7.154940407536e-001,u370=5.715654037330e-002 x371=-7.147251057286e-001,u371=5.733260575482e-002 x372=-7.139561707036e-001,u372=5.750921350246e-002 x373=-7.131872356786e-001,u373=5.768636528700e-002 x374=-7.124183006536e-001,u374=5.786406278436e-002 x375=-7.116493656286e-001,u375=5.804230767553e-002 x376=-7.108804306036e-001,u376=5.822110164679e-002 x377=-7.101114955786e-001,u377=5.840044638962e-002 x378=-7.0934********e-001,u378=5.858034360058e-002 x379=-7.0857********e-001,u379=5.876079498158e-002 x380=-7.078046905037e-001,u380=5.894180223974e-002 x381=-7.070357554787e-001,u381=5.912336708739e-002 x382=-7.062668204537e-001,u382=5.930549124227e-002 x383=-7.054978854287e-001,u383=5.948817642727e-002 x384=-7.0472********e-001,u384=5.967142437058e-002 x385=-7.0396********e-001,u385=5.985523680583e-002 x386=-7.0319********e-001,u386=6.003961547192e-002 x387=-7.024*********e-001,u387=6.022*********e-002 x388=-7.016532103037e-001,u388=6.0410********e-002 x389=-7.008842752787e-001,u389=6.0596********e-002 x390=-7.001153402537e-001,u390=6.078282741069e-002 x391=-6.993464052288e-001,u391=6.097006350266e-002 x392=-6.985774702038e-001,u392=6.115787637203e-002 x393=-6.978085351788e-001,u393=6.134626779549e-002 x394=-6.970396001538e-001,u394=6.153523955530e-002 x395=-6.962706651288e-001,u395=6.172479343917e-002 x396=-6.955017301038e-001,u396=6.191493124037e-002 x397=-6.947327950788e-001,u397=6.210565475756e-002 x398=-6.939638600538e-001,u398=6.229696579513e-002 x399=-6.931949250288e-001,u399=6.248886616284e-002 x400=-6.924259900038e-001,u400=6.268135767612e-002 x401=-6.916570549789e-001,u401=6.287444215608e-002 x402=-6.908881199539e-001,u402=6.306812142928e-002 x403=-6.901191849289e-001,u403=6.326239732793e-002 x404=-6.893502499039e-001,u404=6.345727169002e-002 x405=-6.885813148789e-001,u405=6.365274635902e-002 x406=-6.878123798539e-001,u406=6.384882318422e-002 x407=-6.870434448289e-001,u407=6.404550402057e-002 x408=-6.862745098039e-001,u408=6.424279072866e-002 x409=-6.855055747789e-001,u409=6.444068517492e-002 x410=-6.847366397539e-001,u410=6.463918923139e-002 x411=-6.839677047290e-001,u411=6.483830477611e-002 x412=-6.831987697040e-001,u412=6.503803369265e-002 x413=-6.824298346790e-001,u413=6.523837787052e-002 x414=-6.816608996540e-001,u414=6.543933920506e-002 x415=-6.808919646290e-001,u415=6.564091959736e-002 x416=-6.801230296040e-001,u416=6.584312095443e-002 x417=-6.793540945790e-001,u417=6.604594518918e-002 x418=-6.785851595540e-001,u418=6.624939422036e-002 x419=-6.778162245290e-001,u419=6.645346997257e-002 x420=-6.770472895040e-001,u420=6.665817437656e-002 x421=-6.762783544790e-001,u421=6.686350936876e-002 x422=-6.755094194541e-001,u422=6.706947689176e-002 x423=-6.747404844291e-001,u423=6.727607889400e-002 x424=-6.739715494041e-001,u424=6.748331733005e-002 x425=-6.732026143791e-001,u425=6.769119416038e-002 x426=-6.724336793541e-001,u426=6.789971135160e-002 x427=-6.716647443291e-001,u427=6.810887087629e-002 x428=-6.708958093041e-001,u428=6.831867471320e-002 x429=-6.701268742791e-001,u429=6.852*********e-002 x430=-6.693579392541e-001,u430=6.874022326887e-002 x431=-6.685890042291e-001,u431=6.895197197555e-002x432=-6.678200692042e-001,u432=6.916437297038e-002 x433=-6.670511341792e-001,u433=6.937742826271e-002 x434=-6.662821991542e-001,u434=6.959113986807e-002 x435=-6.655132641292e-001,u435=6.980550980823e-002 x436=-6.647443291042e-001,u436=7.002054011121e-002 x437=-6.639753940792e-001,u437=7.023*********e-002 x438=-6.632064590542e-001,u438=7.0452********e-002 x439=-6.624375240292e-001,u439=7.066961357065e-002 x440=-6.616685890042e-001,u440=7.0887********e-002 x441=-6.608996539792e-001,u441=7.110566848611e-002 x442=-6.601307189542e-001,u442=7.132470390486e-002 x443=-6.593617839293e-001,u443=7.154441405834e-002 x444=-6.585928489043e-001,u444=7.176480102508e-002 x445=-6.578239138793e-001,u445=7.198586689000e-002 x446=-6.570549788543e-001,u446=7.220761374438e-002 x447=-6.562860438293e-001,u447=7.243004368605e-002 x448=-6.555171088043e-001,u448=7.265315881925e-002 x449=-6.547481737793e-001,u449=7.287696125465e-002 x450=-6.539792387543e-001,u450=7.310145310953e-002 x451=-6.532103037293e-001,u451=7.332663650764e-002 x452=-6.524413687043e-001,u452=7.355251357923e-002 x453=-6.516724336794e-001,u453=7.377908646114e-002 x454=-6.509034986544e-001,u454=7.400635729686e-002 x455=-6.501345636294e-001,u455=7.423432823636e-002 x456=-6.493656286044e-001,u456=7.446300143637e-002 x457=-6.485966935794e-001,u457=7.469237906013e-002 x458=-6.478277585544e-001,u458=7.492246327762e-002 x459=-6.470588235294e-001,u459=7.515325626552e-002 x460=-6.462898885044e-001,u460=7.538476020720e-002 x461=-6.455209534794e-001,u461=7.561697729268e-002 x462=-6.447520184544e-001,u462=7.584990971882e-002 x463=-6.439830834295e-001,u463=7.608355968923e-002 x464=-6.432141484045e-001,u464=7.631792941431e-002 x465=-6.424452133795e-001,u465=7.655302111118e-002 x466=-6.416762783545e-001,u466=7.678883700394e-002 x467=-6.409073433295e-001,u467=7.702537932339e-002 x468=-6.401384083045e-001,u468=7.726265030733e-002 x469=-6.393694732795e-001,u469=7.750065220040e-002 x470=-6.386005382545e-001,u470=7.773938725416e-002 x471=-6.378316032295e-001,u471=7.797885772707e-002 x472=-6.370626682045e-001,u472=7.821906588459e-002 x473=-6.362937331795e-001,u473=7.846001399913e-002 x474=-6.355247981546e-001,u474=7.870170435018e-002 x475=-6.347558631296e-001,u475=7.894413922411e-002 x476=-6.339869281046e-001,u476=7.918732091441e-002 x477=-6.332179930796e-001,u477=7.943125172171e-002 x478=-6.324490580546e-001,u478=7.967593395361e-002 x479=-6.316801230296e-001,u479=7.992136992485e-002 x480=-6.309111880046e-001,u480=8.016756195729e-002 x481=-6.301422529796e-001,u481=8.0414********e-002 x482=-6.293733179546e-001,u482=8.0662********e-002 x483=-6.286043829296e-001,u483=8.0910********e-002 x484=-6.278354479047e-001,u484=8.115993738802e-002 x485=-6.270665128797e-001,u485=8.140994480601e-002 x486=-6.262975778547e-001,u486=8.166072236752e-002 x487=-6.255286428297e-001,u487=8.191227244492e-002 x488=-6.247597078047e-001,u488=8.216459741798e-002 x489=-6.239907727797e-001,u489=8.241769967372e-002 x490=-6.232218377547e-001,u490=8.267158160656e-002 x491=-6.224529027297e-001,u491=8.292624561831e-002 x492=-6.216839677047e-001,u492=8.318169411815e-002 x493=-6.209150326797e-001,u493=8.343792952263e-002 x494=-6.201460976547e-001,u494=8.369495425582e-002 x495=-6.193771626298e-001,u495=8.395277074926e-002 x496=-6.186082276048e-001,u496=8.421138144193e-002 x497=-6.178392925798e-001,u497=8.447078878035e-002 x498=-6.170703575548e-001,u498=8.473099521863e-002 x499=-6.163014225298e-001,u499=8.499200321829e-002 x500=-6.155324875048e-001,u500=8.525381524857e-002 x501=-6.147635524798e-001,u501=8.551643378625e-002 x502=-6.139946174548e-001,u502=8.577986131584e-002 x503=-6.132256824298e-001,u503=8.604410032935e-002 x504=-6.124567474048e-001,u504=8.630915332656e-002 x505=-6.116878123799e-001,u505=8.657502281494e-002 x506=-6.109188773549e-001,u506=8.684171130969e-002 x507=-6.101499423299e-001,u507=8.710922133376e-002 x508=-6.0938********e-001,u508=8.737755541778e-002 x509=-6.086120722799e-001,u509=8.764671610032e-002 x510=-6.078431372549e-001,u510=8.791670592772e-002 x511=-6.070742022299e-001,u511=8.818752745413e-002 x512=-6.063052672049e-001,u512=8.845918324151e-002 x513=-6.0553********e-001,u513=8.873167585989e-002 x514=-6.0476********e-001,u514=8.900500788704e-002 x515=-6.039984621300e-001,u515=8.927918190882e-002 x516=-6.032295271050e-001,u516=8.955420051891e-002 x517=-6.024*********e-001,u517=8.983006631909e-002 x518=-6.016916570550e-001,u518=9.010*********e-002 x519=-6.009227220300e-001,u519=9.038434993677e-002 x520=-6.001537870050e-001,u520=9.0662********e-002 x521=-5.993848519800e-001,u521=9.094205373643e-002 x522=-5.986159169550e-001,u522=9.122219479450e-002 x523=-5.978469819300e-001,u523=9.150319882226e-002 x524=-5.970780469050e-001,u524=9.178506847806e-002 x525=-5.963091118800e-001,u525=9.206780642850e-002 x526=-5.955401768551e-001,u526=9.235141534832e-002 x527=-5.947712418301e-001,u527=9.263589792049e-002 x528=-5.940023068051e-001,u528=9.292125683639e-002 x529=-5.932333717801e-001,u529=9.320749479545e-002 x530=-5.924644367551e-001,u530=9.349461450563e-002 x531=-5.916955017301e-001,u531=9.378261868311e-002 x532=-5.909265667051e-001,u532=9.407151005249e-002 x533=-5.901576316801e-001,u533=9.436129134680e-002 x534=-5.893886966551e-001,u534=9.465196530737e-002 x535=-5.886197616301e-001,u535=9.494353468405e-002 x536=-5.878508266052e-001,u536=9.523600223518e-002 x537=-5.870818915802e-001,u537=9.552937072759e-002 x538=-5.863129565552e-001,u538=9.582364293657e-002 x539=-5.855440215302e-001,u539=9.611882164603e-002 x540=-5.847750865052e-001,u540=9.641490964840e-002 x541=-5.840061514802e-001,u541=9.671190974477e-002 x542=-5.832372164552e-001,u542=9.700982474482e-002 x543=-5.824682814302e-001,u543=9.730865746690e-002 x544=-5.816993464052e-001,u544=9.760841073804e-002 x545=-5.809304113802e-001,u545=9.790908739396e-002 x546=-5.801614763552e-001,u546=9.821069027916e-002 x547=-5.793925413303e-001,u547=9.851322224685e-002 x548=-5.786236063053e-001,u548=9.881668615905e-002 x549=-5.778546712803e-001,u549=9.912108488661e-002 x550=-5.770857362553e-001,u550=9.942642130923e-002 x551=-5.763168012303e-001,u551=9.973269831541e-002 x552=-5.755478662053e-001,u552=1.000399188027e-001 x553=-5.747789311803e-001,u553=1.003480856773e-001 x554=-5.740099961553e-001,u554=1.006572018547e-001 x555=-5.732410611303e-001,u555=1.009672702593e-001 x556=-5.724721261053e-001,u556=1.012782938241e-001 x557=-5.717031910804e-001,u557=1.015902754917e-001 x558=-5.709342560554e-001,u558=1.019032182135e-001 x559=-5.701653210304e-001,u559=1.022*********e-001 x560=-5.693963860054e-001,u560=1.025*********e-001 x561=-5.686274509804e-001,u561=1.028*********e-001 x562=-5.678585159554e-001,u562=1.0316********e-001 x563=-5.670895809304e-001,u563=1.034824515723e-001 x564=-5.663206459054e-001,u564=1.038012231100e-001 x565=-5.655517108804e-001,u565=1.0412********e-001 x566=-5.647827758554e-001,u566=1.044417151209e-001 x567=-5.640138408304e-001,u567=1.0476********e-001 x568=-5.632449058055e-001,u568=1.050861592589e-001 x569=-5.624759707805e-001,u569=1.054098709904e-001 x570=-5.617070357555e-001,u570=1.0573********e-001 x571=-5.609381007305e-001,u571=1.060602890906e-001 x572=-5.601691657055e-001,u572=1.063870016123e-001 x573=-5.594002306805e-001,u573=1.067147205664e-001 x574=-5.586312956555e-001,u574=1.070434490530e-001 x575=-5.578623606305e-001,u575=1.0737********e-001 x576=-5.570934256055e-001,u576=1.0770********e-001 x577=-5.563244905805e-001,u577=1.080357228549e-001 x578=-5.555555555556e-001,u578=1.0836********e-001 x579=-5.547866205306e-001,u579=1.0870********e-001 x580=-5.540176855056e-001,u580=1.0903********e-001 x581=-5.532487504806e-001,u581=1.0937********e-001 x582=-5.524798154556e-001,u582=1.097099953205e-001 x583=-5.517108804306e-001,u583=1.100479506995e-001 x584=-5.509419454056e-001,u584=1.103869471442e-001 x585=-5.501730103806e-001,u585=1.107269878615e-001 x586=-5.494040753556e-001,u586=1.110680760685e-001 x587=-5.486351403306e-001,u587=1.114102149918e-001 x588=-5.478662053057e-001,u588=1.117534078682e-001 x589=-5.470972702807e-001,u589=1.120976579443e-001 x590=-5.463283352557e-001,u590=1.124429684769e-001 x591=-5.455594002307e-001,u591=1.127893427325e-001 x592=-5.447904652057e-001,u592=1.131367839882e-001 x593=-5.440215301807e-001,u593=1.134852955306e-001 x594=-5.432525951557e-001,u594=1.138348806568e-001 x595=-5.424836601307e-001,u595=1.141855426740e-001 x596=-5.417147251057e-001,u596=1.145372848994e-001 x597=-5.409457900807e-001,u597=1.148901106608e-001 x598=-5.401768550557e-001,u598=1.152440232958e-001 x599=-5.394079200308e-001,u599=1.155990261526e-001 x600=-5.386389850058e-001,u600=1.159551225895e-001 x601=-5.378700499808e-001,u601=1.163123159755e-001 x602=-5.371011149558e-001,u602=1.166706096894e-001 x603=-5.363321799308e-001,u603=1.170300071210e-001 x604=-5.355632449058e-001,u604=1.173905116702e-001 x605=-5.347943098808e-001,u605=1.177521267473e-001 x606=-5.340253748558e-001,u606=1.181148557736e-001 x607=-5.332564398308e-001,u607=1.184787021803e-001 x608=-5.324875048058e-001,u608=1.188436694096e-001 x609=-5.317185697809e-001,u609=1.192097609141e-001 x610=-5.309496347559e-001,u610=1.195769801572e-001 x611=-5.301806997309e-001,u611=1.199453306128e-001 x612=-5.294117647059e-001,u612=1.203148157657e-001 x613=-5.286428296809e-001,u613=1.206854391112e-001 x614=-5.278738946559e-001,u614=1.210572041555e-001 x615=-5.271049596309e-001,u615=1.214301144156e-001 x616=-5.263360246059e-001,u616=1.218041734194e-001 x617=-5.255670895809e-001,u617=1.221793847054e-001 x618=-5.247981545559e-001,u618=1.225557518233e-001 x619=-5.240292195309e-001,u619=1.229332783336e-001 x620=-5.232602845060e-001,u620=1.233119678078e-001 x621=-5.224913494810e-001,u621=1.236918238284e-001 x622=-5.217224144560e-001,u622=1.240728499889e-001 x623=-5.209534794310e-001,u623=1.244550498938e-001 x624=-5.201845444060e-001,u624=1.248384271590e-001 x625=-5.194156093810e-001,u625=1.252229854112e-001 x626=-5.186466743560e-001,u626=1.256087282885e-001 x627=-5.178777393310e-001,u627=1.259956594400e-001 x628=-5.171088043060e-001,u628=1.263837825262e-001 x629=-5.163398692810e-001,u629=1.267731012189e-001。