[問題] 關於fortran中imsl的一些問題
大家好，小弟有一些imsl使用上的問題
在程式裡，小弟的目的是想先用程式裡寫好的方程式算出所要的初值
再把算好的初值放到xguess裡去給程式收歛
但是debug完後，程式說
Error: This symbol must be a defined parameter or an argument of an inquiry fun
ction that evaluates to a compile-time constant.
Error: This is not a valid initialization expression.
結果只能事先先算好，再把算好的值(數字)放到xguess裡讓程式收歛
這樣程式才會跑出結果
也就是說程式只收歛常數(數字)而已
後來又試過以下方法
1.試過先用副程式算好，再call到主程式裡
2.用READ去讀檔裡的數字
但結果也是出現debug那兩行
不知道有沒有什麼方式可以在同一程式裡算好
再自動代入xguess裡去算
要不然要一組一組輸入很沒效率(因為小弟有上千組數據等著收)
下面是我的程式碼，先謝過大家了
PROGRAM MAIN
USE IMSL
IMPLICIT NONE
INTEGER, PARAMETER :: ITMAX=100
INTEGER, PARAMETER :: NROOT=1
REAL, PARAMETER :: EPS=1.E-6
REAL, PARAMETER :: ERRABS=1.E-6
REAL, PARAMETER :: ERRREL=1.E-6
REAL, PARAMETER :: ETA=1.E-6
INTEGER INFO(NROOT)
REAL, EXTERNAL :: FS,FR
REAL :: XS(NROOT)
REAL :: XR(NROOT)
real :: VR=1./10.
REAL :: XGUESS1(NROOT)=(/vr/)
REAL :: XGUESS2(NROOT)=(/0.09/)
REAL Z,ZSIM,ZREF
CALL ZREAL (FS,ERRABS,ERRREL,EPS,ETA,NROOT,ITMAX,XGUESS1,XS,INFO)
CALL ZREAL (FR,ERRABS,ERRREL,EPS,ETA,NROOT,ITMAX,XGUESS2,XR,INFO)
WRITE(*,"('VR(0)=',F20.18)") XS
WRITE(*,"('VR(r)=',F20.18)") XR
STOP
END
REAL FUNCTION FS(XS)
IMPLICIT NONE
REAL:: XS
REAL ,PARAMETER :: BE1=0.65392
REAL ,PARAMETER :: GA1=0.060167
REAL GB,GC,GD,TR,PR,TEMP,T,P
INTEGER, PARAMETER :: H=4
REAL :: W1(H),Y1(H),Z1(H),Z
REAL, PARAMETER :: TC= 190.56
REAL, PARAMETER :: PC= 45.99
T=100.
P=20.
W1(1)=0.1181193
W1(2)=0.265728
W1(3)=0.15479
W1(4)=0.030323
Y1(1)=0.0236744
Y1(2)=0.0186984
Y1(3)=0.0
Y1(4)=0.042724
Z1(1)=0.0000155488
Z1(2)=0.0000623689
TR=T/TC
PR=P/PC
GB=W1(1)-W1(2)/TR-W1(3)/TR**2-W1(4)/TR**3
GC=Y1(1)-Y1(2)/TR+Y1(3)/TR**3
GD=Z1(1)+Z1(2)/TR
FS=TR/PR*(1/XS+GB/XS**2+GC/XS**3+GD/XS**6+Y1(4)/TR**3/XS**3*(BE1+GA1/XS**2)*E
XP(-GA1/XS**2))-1
RETURN
END FUNCTION
REAL FUNCTION FR(XR)
IMPLICIT NONE
REAL:: XR
REAL ,PARAMETER :: BE2=1.226
REAL ,PARAMETER :: GA2=0.03754
REAL GB,GC,GD,TR,PR,TEMP,T,P
INTEGER, PARAMETER :: H=4
REAL :: W2(H),Y2(H),Z2(H),Z
REAL, PARAMETER :: TC= 190.56
REAL, PARAMETER :: PC= 45.99
T=100.
P=20.
W2(1)=0.2026579
W2(2)=0.331511
W2(3)=0.027655
W2(4)=0.203488
Y2(1)=0.0313385
Y2(2)=0.0503618
Y2(3)=0.016901
Y2(4)=0.041577
Z2(1)=0.000048736
Z2(2)=0.00000740336
TR=T/TC
PR=P/PC
GB=W2(1)-W2(2)/TR-W2(3)/TR**2-W2(4)/TR**3
GC=Y2(1)-Y2(2)/TR+Y2(3)/TR**3
GD=Z2(1)+Z2(2)/TR
FR=TR/PR*(1/XR+GB/XR**2+GC/XR**3+GD/XR**6+Y2(4)/TR**3/XR**3*(BE2+GA2/XR**2)*E
XP(-GA2/XR**2))-1
RETURN
END FUNCTION
--
※ 發信站: 批踢踢實業坊(ptt.cc)
◆ From: 140.124.50.204
※ 編輯: m710387 來自: 140.124.50.204 (09/19 14:32)