24小时热门版块排行榜

返回列表

星星阿雅

新虫 (初入文坛)

应助: 0 (幼儿园)
金币: 3.5
帖子: 3
在线: 5.4小时
虫号: 3215265
注册: 2014-05-18

[求助] 做毕业论文，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