在fortran中子例程的主程序中将输出写入相同的txt文件中

Question

我已经写了一个包含6个子例程的代码，首先我在主程序中写了输出，我得到了正确的文本文件，但现在我还想在同一个文本文件中添加一个变量的另一个子例程的输出。这是我的代码示例

 ! code for the particle tracing in magnetic field    
      implicit none
      include "mydata.inc"
      External ad,bsstep
      real*8 x(ndim,nv),b(4)

      real*8 pwE(ndim,3),pwB(ndim,3),ep

      real*8 ke,gamma,kew,epx,hmin,rr,epx0
      real*8 eps(3),eta0(5)
      real*8 m,d_eta
      real*8 t,tw,dt,rt
      integer i,j,k,nok,nbad,filen
      character*80 fname1
 !     real,parameter :: pi

      t=0.d0
      tw=0.d0
      m=me
      x=0.d0
      pwE=0.d0
      pwB=0.d0

 !     

      !=====loop================!
      k=1
 120  call ini_v(x,ke,pwB)
      x(2,3)=eta0(k)
      call dispersion(pwE,pwB,eps,ep)

      write(fname1,1005) k

      open(1,file=fname1)

      write(1,*)'time,eta,phi,pitch angle,ke,deta'

      !====normalization========!
      ke=ke/(m*c**2)
      gamma=ke+1
      x(1,2)=cos(alpha)*dsqrt(ke**2+2*ke)
      x(2,2)=sin(alpha)*dsqrt(ke**2+2*ke)
      dt=0.001

      epx=1.d-6       ! required accuracy
    epx0=epx
    hmin=0.d0
      !================================

      rt=0
      kew=ke*(m*c**2)/1.0d6/e !ke in the units of electron volts

 100  write(*,*)rt
      if(rt.gt.15)then
          goto 110
      endif
      if(rt.eq.0)then
      write(1,1000)t/(2*pi),
     &(x(i,3),i=1,2),alpha*180.d0/pi,kew,d_eta

      else
      write(1,1000)t/(2*pi),
     &(x(i,3),i=1,2),alpha*180.d0/pi,kew,d_eta
      endif

      !=========adaptivestepsize method==============!

      call odeint(x,t,t+dt,epx,dt,hmin,nok,nbad,ad,bsstep,
     &pwE,pwB,gamma,ep)

      x(1:2,3)=mod(x(1:2,3),2*pi) !this is the range of angles
      alpha=dacos(x(1,2)/dsqrt(x(1,2)**2+x(2,2)**2))
      call go(KE,x(1:2,2),gamma)
      KEw=KE*(m*c**2)/1.0d6/e
      t=t+dt
      rt=rt+dt
      goto 100

 110  write(*,*) eta0(k)

      if(k.le.size(eta0)) then
          k=k+1
          t=0
          goto 120

      endif 

 1000 format(6(E10.4,2x),4(f15.6,2x))
 1005  format("output",i2.2,".txt")     


      stop
      end

我想要从中获得值的子例程是这样的

=

! subroutine to calculate all the derivatives for this code
      subroutine ad(step,kin,kout,pwE,pwB,gamma,ep)
      include "mydata.inc"
      real*8 step
      real*8 kin(ndim,nv),kout(ndim,nv)
      real*8 pwE(ndim,3),pwB(ndim,3),ep
      real*8 gamma
      real*8 tm(3)







      kout(1,1)=kin(2,2)*cos(kin(2,3)-kin(1,3))/gamma
      kout(2,1)=kin(2,2)*sin(kin(2,3)-kin(1,3))/gamma
      kout(3,1)=kin(1,2)/gamma



      kout(1,2)=-pwE(3,1)*sin(kin(1,3))
     &+(pwB(1,1)*kin(2,2)/(2.*gamma))
     &*((1+ep)*sin(kin(2,3))+(1-ep)*sin(kin(2,3)-2*kin(1,3)))

      tm(1)=-pwB(1,1)*kin(1,2)/2./gamma
      tm(2)=pwB(1,1)*emicwave/wavek/cos(si)
      tm(3)=tan(si)*pwE(3,1)

      kout(2,2)=sin(kin(2,3))*(tm(1)*(1+ep)+1./2.*((1+ep)*tm(2)-tm(3)))
     &+sin(kin(2,3)-2*kin(1,3))
     &*(tm(1)*(1-ep)+1./2.*((1-ep)*tm(2)+tm(3))) 



      kout(1,3)=EMICwave-kin(1,2)*wavek*cos(si)/gamma
     &-kin(2,2)*wavek*sin(si)*cos(kin(2,3)-kin(1,3))/gamma

!----------------------d(eta)/dt --------------------------!
!     kin(2,3)=eta, kout(2,3)=d(eta)/dt
!   time derivatives of phase(phi)
!----------------------------------------------------------!

      tm(1)=-2*tm(1)
      kout(2,3)=EMICwave-kin(1,2)*wavek*cos(si)/gamma
     &-kin(2,2)*wavek*sin(si)*cos(kin(2,3)-kin(1,3))/gamma
     &+1./gamma-cos(kin(2,3))/(2.*kin(2,2))
     &*((1+ep)*tm(1)-(1+ep)*tm(2)+tm(3))
     &-cos(kin(2,3)-2*kin(1,3))/(2.*kin(2,2))
     &*((1-ep)*tm(1)-(1-ep)*tm(2)-tm(3))
     &-pwB(3,1)/gamma*cos(kin(1,3))



    return
    end

Answer 1

首先要做的是:我会远离小于10的文件单元--您永远不会知道您的编译器是否恰好使用这个特定的数字来表示标准输出或类似的东西(大多数使用单元0表示标准错误，使用5和6表示标准输入和标准输出，但不能保证)。

最安全的方法是在打开时使用newunit参数：

 integer :: unt

 open(newunit=unt, file=filename)

 write(unt, *) "Some text"

 close(unt)

如果在调用子例程时文件是打开的，那么理论上该子例程应该只需要该文件单元来将输出写入该文件(假设子例程不是pure)。

如果您愿意，可以将文件i/o传递给一个模块，然后在主程序和子例程中对其执行use操作。