| 查看: 559 | 回复: 0 | |||
[求助]
做毕业论文，fortran程序执行遇到问题，急求各位大侠回复！！！在线等！！
|
|
程序如下：执行时的代码见图片，感觉open时候就出问题了，本人小白，求大家指导！！！谢谢！！ parameter(m=713,n=81,mnl=81,ks=0) real x(m,n),egvt(m,mnl),ecof(mnl,n),er(mnl,4) open(21,file='F:\890\890\test.txt',STATUS='old') do i=1,m read(21,*) (x(i,j),j=1,n) enddo close(21) call eof(m,n,mnl,x,ks,er,egvt,ecof) open(22,file='F:\890\890\er.txt',STATUS='old') do i=1,mnl write(22,"(4f30.8)"  (er(i,j),j=1,4)enddo close(22) open(23,file='F:\890\890\ecof.txt',STATUS='old') do j=1,n write(23,"(<mnl>f20.6)" (ecof(i,j),i=1,mnl)enddo close(23) open(24,file='F:\890\890\egvt.txt',STATUS='old') do i=1,m write(24,"(<mnl>f20.6)" (egvt(i,j),j=1,mnl)enddo close(24) print*,"-------------------------------" print*,"这是来自李建平老师的EOF分解程序" print*,"http://bbs.06climate.com整理制作" print*,"-------------------------------" end !*************************************************************************** subroutine eof(m,n,mnl,f,ks,er,egvt,ecof) implicit none integer*4 ::m,n,mnl,ks real*4 ::f(m,n),er(mnl,4),egvt(m,mnl),ecof(mnl,n) real*4,allocatable :: cov(:,  ,s(:,,d(call transf(m,n,f,ks) allocate(cov(mnl,mnl)) call crossproduct(m,n,mnl,f,cov) allocate(s(mnl,mnl)) allocate(d(mnl)) call jacobi(mnl,cov,s,d,0.00001) call arrang(mnl,d,s,er) call tcoeff(m,n,mnl,f,s,er,egvt,ecof) return end subroutine transf(m,n,f,ks) implicit none integer*4 :: m,n,ks real*4 :: f(m,n) real*4,allocatable ::fw( ,wn(integer*4 :: i0,i,j real*4 :: af,sf,vf allocate(fw(n)) allocate(wn(m)) i0=0 do i=1,m do j=1,n fw(j)=f(i,j) enddo call meanvar(n,fw,af,sf,vf) if(sf.eq.0.)then i0=i0+1 wn(i0)=i endif enddo if(i0.ne.0)then write(*,*)'************************* fault ******' write(*,*)' the original field has invalid data.' write(*,*)'There are totally ',i0,' invalid data.' write(*,*)'The array wn(m) stores the positions.' write(*,*)'You must pick off those invalid data ---$' write(*,*)'$---and then reinput a new field eofs.' write(*,*)'************************ fault ********' endif if(ks.eq.-1)return if(ks.eq.0)then do i=1,m do j=1,n fw(j)=f(i,j) enddo call meanvar(n,fw,af,sf,vf) do j=1,n f(i,j)=f(i,j)-af enddo enddo return endif if(ks.eq.1)then do i=1,m do j=1,n fw(j)=f(i,j) enddo call meanvar(n,fw,af,sf,vf) if(sf==0.0)then do j=1,n f(i,j)=0.0 enddo else do j=1,n f(i,j)=(f(i,j)-af)/sf enddo endif enddo endif return end subroutine crossproduct(m,n,mnl,f,cov) implicit none integer*4 ::m,n,mnl real*4 :: f(m,n),cov(mnl,mnl) integer*4 :: i,j,is,js if((n-m)<0)then do i=1,mnl do j=i,mnl cov(i,j)=0.0 do is=1,m cov(i,j)=cov(i,j)+f(is,i)*f(is,j) enddo cov(j,i)=cov(i,j) enddo enddo else do i=1,mnl do j=i,mnl cov(i,j)=0.0 do js=1,n cov(i,j)=cov(i,j)+f(i,js)*f(j,js) enddo cov(j,i)=cov(i,j) enddo enddo endif return end subroutine jacobi(m,a,s,d,eps) implicit none integer*4 :: m real*4 :: a(m,m),s(m,m),d(m),eps integer*4 :: i,j,i1,l,iq,iq1,ip real*4 :: g,s1,s2,s3,v1,v2,v3,u,st,ct s=0.0 do i=1,m s(i,i)=1.0 enddo g=0. do i=2,m i1=i-1 do j=1,i1 g=g+2.0*a(i,j)*a(i,j) enddo enddo s1=sqrt(g) s2=eps/float(m)*s1 s3=s1 l=0 50 s3=s3/float(m) 60 do iq=2,m iq1=iq-1 do ip=1,iq1 if(abs(a(ip,iq)).lt.s3) exit l=1 v1=a(ip,ip) v2=a(ip,iq) v3=a(iq,iq) u=0.5*(v1-v3) if(u.eq.0.0) g=1. if(abs(u).GE.1e-10) g=-sign(1.,u)*v2/sqrt(v2*v2+u*u) st=g/sqrt(2.*(1.+sqrt(1.-g*g))) ct=sqrt(1.-st*st) do i=1,m g=a(i,ip)*ct-a(i,iq)*st a(i,iq)=a(i,ip)*st+a(i,iq)*ct a(i,ip)=g g=s(i,ip)*ct-s(i,iq)*st s(i,iq)=s(i,ip)*st+s(i,iq)*ct s(i,ip)=g enddo do i=1,m a(ip,i)=a(i,ip) a(iq,i)=a(i,iq) enddo g=2.*v2*st*ct a(ip,ip)=v1*ct*ct+v3*st*st-g a(iq,iq)=v1*st*st+v3*ct*ct+g a(ip,iq)=(v1-v3)*st*ct+v2*(ct*ct-st*st) a(iq,ip)=a(ip,iq) enddo enddo if((l-1)==0)then l=0 go to 60 else go to 150 endif 150 if(s3.gt.s2) then goto 50 end if do i=1,m d(i)=a(i,i) enddo return end subroutine arrang(mnl,d,s,er) implicit none integer*4 ::mnl real*4 :: d(mnl),s(mnl,mnl),er(mnl,4) integer*4 ::i,mnl1,k1,k2 real*4 ::c,tr tr=0.0 do i=1,mnl tr=tr+d(i) er(i,1)=d(i) enddo mnl1=mnl-1 do k1=mnl1,1,-1 do k2=k1,mnl1 if(er(k2,1).lt.er(k2+1,1)) then c=er(k2+1,1) er(k2+1,1)=er(k2,1) er(k2,1)=c do i=1,mnl c=s(i,k2+1) s(i,k2+1)=s(i,k2) s(i,k2)=c enddo endif enddo enddo er(1,2)=er(1,1) do i=2,mnl er(i,2)=er(i-1,2)+er(i,1) enddo do i=1,mnl er(i,3)=er(i,1)/tr er(i,4)=er(i,2)/tr enddo return end subroutine tcoeff(m,n,mnl,f,s,er,egvt,ecof) implicit none integer*4 :: m,n,mnl real*4::f(m,n),s(mnl,mnl),er(mnl,4),egvt(m,mnl),ecof(mnl,n) real*4,allocatable:: v( integer*4 :: i,j,js,is real*4 ::c allocate(v(mnl)) do j=1,mnl do i=1,m egvt(i,j)=0. enddo do i=1,n ecof(j,i)=0. enddo enddo do j=1,mnl c=0. do i=1,mnl c=c+s(i,j)*s(i,j) enddo c=sqrt(c) do i=1,mnl s(i,j)=s(i,j)/c enddo enddo if(m.le.n) then do js=1,mnl do i=1,m egvt(i,js)=s(i,js) enddo enddo do j=1,n do i=1,m v(i)=f(i,j) enddo do is=1,mnl do i=1,m ecof(is,j)=ecof(is,j)+v(i)*s(i,is) enddo enddo enddo else do i=1,m do j=1,n v(j)=f(i,j) enddo do js=1,mnl do j=1,n egvt(i,js)=egvt(i,js)+v(j)*s(j,js) enddo enddo enddo do js=1,mnl do j=1,n ecof(js,j)=s(j,js)*sqrt(abs(er(js,1))) enddo do i=1,m egvt(i,js)=egvt(i,js)/sqrt(abs(er(js,1))) enddo enddo endif return end subroutine meanvar(n,x,ax,sx,vx) implicit none integer*4 ::n real*4 :: x(n) real*4 :: ax,vx,sx integer*4 :: i ax=0. vx=0. sx=0. do i=1,n ax=ax+x(i) enddo ax=ax/float(n) do i=1,n vx=vx+(x(i)-ax)**2 enddo vx=vx/float(n) sx=sqrt(vx) return end  执行.png |
回复此楼» 猜你喜欢
filecode
已经有13人回复
-
应该是93bebmhtak前后十一个字符比较关键
已经有35人回复
-
长年满屏的广告，版主太不责任了。基金也等的急
已经有3人回复
-
不应该看fileCode
已经有4人回复
-
FileCode能看出啥？
已经有19人回复
-
2019年青年基金涵评意见，大家看看几个A，几个B?
已经有5人回复
-
职称评审，求友友推荐见刊最快的期刊
已经有4人回复
-
奇怪，两个人的filecode固定段从头到尾一模一样
已经有10人回复
-
综述论文作为代表作会不会影响评审专家的印象分？
已经有13人回复
-
是这周出结果还是下周出结果？
已经有3人回复
找到一些相关的精华帖子，希望有用哦~
请教个问题！谢谢
已经有3人回复
- 急：求用Fortran或者其他语言编写一个程序，处理下面的数据！！！！！！！！！！
已经有37人回复
- 使用music_code运行出现错误
已经有4人回复
- matlab 运行ode45出错，不知道什么原因
已经有4人回复
- VS2012和IVF2013中调用已有程序显示“无法启动程序，系统找不到指定文件”怎么回事？
已经有4人回复
- fortran程序在SSH计算中心上运行的问题！
已经有13人回复
- Fortran编程过程中遇到错误，求大侠帮忙看下
已经有5人回复
- CFX求解时程序出错问题（急！！！）
已经有5人回复
- 求大家帮帮忙，有fortran的主程序和几个模块，如何让他们运行，急！！！
已经有5人回复
- 急求fortran运行错误原因，在线等
已经有7人回复
- 请教关于fortran的ex后缀可执行文件，谢谢
已经有3人回复
- 急求可用的fortran编译器
已经有13人回复
- 写了一个fortran90的小程序，编译通不过，请大侠帮忙
已经有59人回复
- 紧急求助，在线等各位达人，看看编辑回复的这句话是什么意思
已经有14人回复
- 【求助】急求一个NAG FORTRAN LIBRARY里的一个Fortran程序，名字是”F02FJF“！谢谢！
已经有10人回复
点击这里搜索更多相关资源科研从小木虫开始，人人为我，我为人人
