我利用IVPAG做SHOOTINT METHOD解熱傳方程式
1.運算出來最後的值是199.9556度.希望可能在更接近兩百
請問該如何修正可以讓我的答案更好?
2.我發現如果今天我要換一個方程式.將SUBROUTINE FCN中的
H,YPRMINE(2)
改為
h=5.0E-8.
YPRMINE(2)=-h*(20-Y(1))**4
這時將會在第二輪停止並顯示
FATAL ERROR3 FROM IVPAG MAXIMUM NUMBER OF STEPS ALLOWED ,500, HAS BEEN....
我嘗試去修改主程式一開始宣告的 NPARAM 值.但是代入0,4,32(不確定..
就是改動STEP的那項)也都沒有效果.請問該修改哪個部份才能順利求值?
P.S改為4次厚.Y(1)的最後一項也要接近兩百_
先感謝各位
另外附上程式碼
PROGRAM NEWWORK
USE IMSL
IMPLICIT NONE
INTEGER N, NPARAM
INTEGER :: K=1
PARAMETER (N=2, NPARAM=50)
INTEGER IDO, NOUT
REAL A(1,1), PARAM(NPARAM), T, TEND,TOL, Y(N), BC, G(3), S(2),
+ IEND
DATA (G(K), K= 1,2) /10,11/
DATA (S(K),K=1,2) /0,0/
EXTERNAL FCN, FCNJ,SHOOT
CALL SSET (NPARAM, 0.0, PARAM, 1)
IDO = 1
T = 0.0
Y(1) = 40
Y(2) = G(K)
TOL = 1.0D-6
BC = 200
CALL UMACH (2,NOUT)
WRITE (NOUT,99998)
DO WHILE (ABS(BC-Y(1))>0.1)
DO WHILE (K<3)
IEND = 0
IDO = 1
10 CONTINUE
IEND = IEND + 0.5
TEND = IEND
CALL IVPAG(IDO, N, FCN,FCNJ, A, T, TEND, TOL, PARAM, Y)
IF (IEND .LE. 10) THEN
WRITE (NOUT,99999) T, Y
IF(IEND .EQ. 10) IDO = 3
GO TO 10
END IF
S(K)=Y(2)
K=K+1
END DO
CALL SHOOT(G, S, K)
G(1)=G(2)
S(1)=S(2)
G(2)=G(3)
K=2
99998 FORMAT (11X, 'L', 14X, 'T1', 11X, 'T2')
99999 FORMAT (3F15.5)
END DO
END
SUBROUTINE FCN (N, X, Y, YPRIME)
INTEGER N
REAL X, Y(N),YPRIME(N)
PARAMETER h=1.0E-2
YPRIME(1) = Y(2)
YPRIME(2) = -h*(20-Y(1))
RETURN
END
SUBROUTINE FCNJ (N, X, Y, DYPDY)
INTEGER N
REAL X, Y(N), DYPDY(N,*)
RETURN
END
SUBROUTINE SHOOT(G, S, k)
INTEGER K
REAL G(K), S(K)
BCC=200
G(3)=G(2)+(BCC-S(2))/(S(2)-S(1))*(G(2)-G(1))
RETURN
END