不好意思
我最近在改寫一個FORTRAN的程式碼(是關於把汙染物放進河川裡)
我把汙染物的初始濃度 設定為1000000000 (1E9)
然後數據跑出來 很多都是*********星號(代表數據過大)
我覺得 這應該是我給fortran的空間不夠
但我修改之後 還是沒辦法修正結果 還是會跑出*********星號
我想可能是我修改錯誤
想請大家幫忙看一下 我可能錯在哪裡 相關code如下:
PROGRAM DISPERSION
IMPLICIT NONE
REAL*8 :: CC(300), Xstop(10), Xgauge(10), UU(10), EE(10)
REAL*8 :: INI_X, END_X, E, U, DE, DX
REAL*8 :: T, TT, DT
INTEGER :: N, I, II, IT, ITT, ITW, WT, NO_S
CHARACTER(LEN=100) :: filename
CHARACTER(LEN=40) :: ext=".dat"
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! INPUT !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
CHARACTER(LEN=20000) :: casename='test123' ! directory for saveing files
INI_X = 0. ! initial of x
END_X = 1700. ! end of x
DX = 10.
N = INT((END_X-INI_X)/DX)+1 ! number of node
T = 0. ! the same unit as U and E
ITW = 50 ! iteration steps duration to write
TT = 3120. ! end time of T
DT = 0.1
CC(1) = 1E9 ! concentration of first node at time = 0.0
CC(N) = 0. ! concentration of last node at time = 0.0
DE = 0.0000008333333333333333333333
NO_S = 3 ! number of station
!!!!! stations !!!!!
Xstop(1) = 0. ! first station
Xstop(2) = 850. ! second station
Xstop(3) = 1700. ! third station
!!!!! velocity, U !!!!!
UU(1) = 0.49 ! velocity between 1st and 2nd station
UU(2) = 0.49 ! velocity between 2nd and 3rd station
!!!!! diffusion coefficient, E !!!!!
EE(1) = 6.5 ! coefficient between 1st and 2nd station
EE(2) = 6.5 ! coefficient between 2nd and 3rd station
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! END SETTING !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
OPEN(101, FILE='time'//ext, STATUS='UNKNOWN')
DO I = 1, NO_S
II = 1
DO WHILE ((INI_X+float(II-1)*DX-Xstop(I)) .LE. 0.)
II = II+1
ENDDO
Xgauge(I) = II-1
WRITE(filename,'(I3.3)') I ! filename contains the time you write
OPEN(301,file='station'//TRIM(filename)//ext,STATUS='unknown')
ENDDO
DO I=2, N-1
CC(I)=0.
ENDDO
IT = 0
ITT = 0
WT = 0
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! start !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
WRITE(101,102) T
102 FORMAT(F100.95)
WRITE(filename,'(I5.5)') WT ! filename contains the time you write
OPEN(201,file=TRIM(filename)//ext,STATUS='unknown')
WRITE(201,202) (CC(I),I=1,N)
202 FORMAT(2X,10000F8.4)
DO II = 1, NO_S
WRITE(filename,'(I3.3)') II ! filename contains the time you write
OPEN(301,file='station'//TRIM(filename)//ext,STATUS='old')
WRITE(301,302) CC(INT(Xgauge(II)))
302 FORMAT(F8.4)
ENDDO
OPEN(401,file='x'//ext,STATUS='unknown')
DO I = 1, N
WRITE(401,*) (INI_X+(I-1)*DX)
ENDDO
12 T=T+DT
IT=IT+1
ITT = ITT+1
DO WHILE (T .LT. TT+DT)
DO I=2,N-1
IF ((I-1)*DX+INI_X .GE. Xstop(1) .AND. (I-1)*DX+INI_X .LT. Xstop(2))
THEN
U = UU(1)
E = EE(1)
ELSEIF ((I-1)*DX+INI_X .GE. Xstop(2) .AND. (I-1)*DX+INI_X .LT.
Xstop(3)) THEN
U = UU(2)
E = EE(2)
ENDIF
CC(I)=CC(I)+E*DT*(CC(I+1)-2*CC(I)+CC(I-1))/DX/DX-&
U*DT*(CC(I+1)-CC(I-1))/2./DX-&
DE*DT*CC(I)
CC(N+1) = CC(N-1)
ENDDO
IF(IT.LT.ITW) GOTO 12
WRITE(101,102) T
WT = WT+1
WRITE(filename,'(I5.5)') WT ! filename contains the time you write
OPEN(201,file=TRIM(filename)//ext,STATUS='unknown')
WRITE(201,202) (CC(I),I=1,N)
IT = 0
DO II = 1, NO_S
WRITE(filename,'(I3.3)') II ! filename contains the time you write
OPEN(301,file='station'//TRIM(filename)//ext,STATUS='old',ACCESS='APPEND')
WRITE(301,302) CC(INT(Xgauge(II)))
ENDDO
WRITE(*,15) T
15 FORMAT(1X,'output data at time = ',F10.4)
GOTO 12
ENDDO
CLOSE(101); CLOSE(201); CLOSE(301); CLOSE(401);
call system('mkdir ' // casename)
call system('move *.dat ' // casename)
STOP
END