五点光滑法的fortran程序

1、精选优质文档-倾情为你奉上五点光滑法的fortran程序 program main use dflibparameter(i_m=40)integer ireal(kind=8)px(-1:i_m),py(-1:i_m)real(kind=8)xc(1:i_m-2),yc(1:i_m-2),xd(1:i_m-2),yd(1:i_m-2)real(kind=8)xa(1:i_m-2),ya(1:i_m-2),xb(1:i_m-2),yb(1:i_m-2)real(kind=8)K1(1:i_m-2),K2(1:i_m-2),K3(1:i_m-2),K4(1:i_m-2)real(kind=8)K

2、5(1:i_m-1)real(kind=8)c1(1:i_m-2),c2(1:i_m-2),c3(1:i_m-2),c4(1:i_m-2)REALAA,BB,DD,EE,diffrealaaa(1:i_m-2),bbb(1:i_m-2),ccc(1:i_m-2),ddd(1:i_m-2)realdiff2(1:i_m-2),right(1:i_m-2),left(1:i_m-2)integer : result!取回画图函数运行状态integer(kind=4) : color!设定颜色用real(kind=8), parameter : X_Start=-10.0! x轴最小范围real(k

3、ind=8), parameter : X_End=20.0! x轴最大范围 real(kind=8), parameter : Y_Top=20.0! y轴最大范围 real(kind=8), parameter : Y_Bottom=-10.0! y轴最小范围type(wxycoord): wt!返回上一次的虚拟坐标位置integer,parameter:cpoints=38 !传入坐标点的数目type(wxycoord) : ppoints(cpoints) !用数组来存储多边形顶点的坐标open(unit=1,file="ppoints.txt")d

4、o i=-1,i_mread(1,*)px(i),py(i)enddoclose(1) open(unit=3,file="diff1.txt")!*求第一点第二点和倒数第一点第二点的导数*!*(需要采用补点的方法,补点采用增量相等的原则来补)*PX(0)=PX(3)-3*(PX(2)-PX(1) PY(0)=PY(3)-3*(PY(2)-PY(1)PX(-1)=PX(2)-3*(PX(1)-PX(0) PY(-1)=PY(2)-3*(PY(1)-PY(0)PX(i_m-1)=PX(i_m-4)-3*(PX(i_m-3)-PX(i_m-2) PY(i_m-1

5、)=PY(i_m-4)-3*(PY(i_m-3)-PY(i_m-2)PX(i_m)=PX(i_m-3)-3*(PX(i_m-2)-PX(i_m-1)PY(i_m)=PY(i_m-3)-3*(PY(i_m-2)-PY(i_m-1)open(unit=2,file="pp.txt")do i=-1,i_mwrite(2,*)px(i),py(i)enddoclose(2) !*求Pi的导数*do 10 i=1,i_m-2K1(i)=(py(i-1)-py(i-2)/(px(i-1)-px(i-2)C1(i)=0.5*(py(i-1)+py(i-2)-K1(i)*

6、(px(i-1)+px(i-2)K2(i)=(py(i)-py(i-1)/(px(i)-px(i-1)C2(i)=0.5*(py(i)+py(i-1)-K2(i)*(px(i)+px(i-1)K3(i)=(py(i+1)-py(i)/(px(i+1)-px(i)C3(i)=0.5*(py(i+1)+py(i)-K3(i)*(px(i+1)+px(i)K4(i)=(py(i+2)-py(i+1)/(px(i+2)-px(i+1)C4(i)=0.5*(py(i+2)+py(i+1)-K4(i)*(px(i+2)+px(i+1)XC(i)=-(C3(i)-C1(i)/(K3(i)-K1(i)YC(i

7、)=K1(i)*XC(i)+C1(i)XD(i)=-(C4(i)-C2(i)/(K4(i)-K2(i)YD(i)=K2(i)*XD(i)+C2(i)write(*,*)XC(i),YC(i),XD(i),YD(i)diff=1000XA(i)=-10YA(i)=K1(i)*XA(i)+C1(i)do 11 while(diff.gt.1e-3.and.XA(i).lt.1000)XB(i)=(px(i)*(py(i)-YA(i)+(C4(i)-py(i)*(px(i)-XA(i)/(py(i)-YA(i)-K4(i)*px(i)+K4(i)*XA(i)YB(i)=K4(i)*XB(i)+C4(

8、i)aa=(py(i-1)-YC(i)*2+(px(i-1)-XC(i)*2bb=(py(i+1)-YD(i)*2+(px(i+1)-XD(i)*2dd=(YC(i)-YA(i)*2+(XC(i)-XA(i)*2ee=(YD(i)-XB(i)*2+(XD(i)-YB(i)*2diff=ABS(aa/bb-dd/ee)XA(i)=XA(i)+0.0001YA(i)=K1(i)*XA(i)+C1(i)11continue K5(i)=(py(i)-YA(i)/(px(i)-XA(i)write(3,800)i,xa(i),diff800 format(1x,i2,4f)10continueK5(i

9、 _m-1)=K5(i_m-2)-0.02open(unit=5,file="k5.txt")do 14 i=1,i_m-1write(5,*)i,k5(i)14continueclose(5)!*求Pi和Pi-1两点之间的曲线函数关系式*open(unit=4,file="diff2.txt")do 12 i=1,i_m-2aaa(i)=-0.6diff2=1000do 13 while(diff2(i).gt.1e-3.and.aaa(i).lt.0.3)bbb(i)=(K5(i+1)-K5(i)-3*aaa(i)*(

10、px(i+1)*2-px(i)*2)/(px(i+1)-px(i)/2ccc(i)=(K5(i+1)+K5(i)-3*aaa(i)*(px(i+1)*2+px(i)*2)-2*bbb(i)*(px(i+1)-px(i)*0.5right(i)=aaa(i)*(px(i+1)*px(i)*3-px(i)*px(i+1)*3)+bbb(i)*(px(i+1)*px(i)*2-px(i)*px(i+1)*2)left(i)=px(i+1)*py(i)-px(i)*py(i+1)ddd(i)=(left(i)-right(i)/(px(i+1)-px(i)diff2(i)=abs(aaa(i)*px(

11、i)*3+bbb(i)*px(i)*2+ccc(i)*px(i)+ddd(i)-py(i)aaa(i)=aaa(i)+1e-613continue write(4,700)i,aaa(i),bbb(i),ccc(i),ddd(i),diff2(i)700 format(1x,i2,5f)12continue!*画图并测试五点光滑法的效果*result=SETBKCOLORRGB(#FFFFFF) !SETBKCOLORRGB()函数调用返回值 或者直接用call语句都可以call CLEARSCREEN($GCLEARSCREEN)result=SetWindow(.true. , X_Sta

12、rt, Y_Top, X_End, Y_Bottom ) ! 设定虚拟坐标范围，Y轴向上为正。 result=SetColorRGB(#)! 利用color来设定颜色call MoveTo_w(X_Start,0.0_8,wt)! result=LineTo_w(X_End,0.0_8)! 画X轴 call MoveTo_w(0.0_8,Y_Top,wt)! result=LineTo_w(0.0_8,Y_Bottom)! 画Y轴 !*画点连接线*result=setcolorRGB(#00ff00) !画笔颜色为绿色do 19 i=1,i_m-3call moveto_w(px(i),py(i),wt) result=lineto_w(px(i+1),py(i+1) 19continue!*画拟合光滑曲线*result=setcolorRGB(#0000ff) !画笔颜色为红色do 29 i=1,