简单力场下粒子边界条件设置 - Fortran Coder

[Fortran] 纯文本查看 复制代码
        program main

        INTEGER*4 MTIME, MXPARM,NP, IDO, MC, istep

        PARAMETER (MXPARM=50,MTIME=5000000000, NP=50)

        REAL*8 Y(NP),YPRIME(NP)

        REAL*8 PARAM(MXPARM)

        real*8 t_stop, t_start,del_t

        REAL*8 TOL,T1,T2

        COMMON /blk/a

        external DIVPRK,DSET,FCN

        t_start=0d0

        t_stop=7d0

        del_t=0.1d0

        TOL=1.0D-3

ccccc初值:y(1)=x0=-12,y(2)=v0=3

        y(1)=-12

        y(2)=3

      CALL DSET(MXPARM, 0.0D0, PARAM, 1)

        PARAM(10)=1.0D0

        PARAM(4)=MTIME+1

        istep=0

        IDO=1

        T1=T_START

        mt=0

110        CONTINUE

c

        istep=istep+1

        T2=T1+DEL_T

        CALL DIVPRK(IDO,NP,FCN,T1,T2,TOL,PARAM,Y)

        write(*,*)t2,y(2),y(1)

              IF(T2<=T_STOP) GOTO 110

        print *,'end'

        stop

        end

      subroutine FCN (NP, T2, Y, YPRIME)

      INTEGER*4 NP

        REAL*8 Y(NP),YPRIME(NP)

      COMMON /blk/a

        a=2

      YPRIME(1) = Y(2)

      YPRIME(2) = a

      return

      end subroutine