这是一个IVPAG的范例,我照着打.
PROGRAM HOMEWORK
IMPLICIT NONE
INTEGER N, NPARAM
PARAMETER (N=3, NPARAM=50)
C
C 宣告区域变量
C
INTEGER IDO, IEND, NOUT
REAL A(1,1), PARAM(NPARAM), T, TEND, TOL, Y(N)
C
C 宣告使用的副程式
C
EXTERNAL IVPAG, SSET, UMACH
C
C 宣告函数
C
EXTERNAL FCN, FCNJ
C
C 初始设定
C
CALL SSET (NPARAM, 0.0, PARAM, 1)
IDO = 1
T = 0.0
Y(1) = 0.0
Y(2) = 1.0
Y(3) = 1.0
TOL = 1.0E-6
C
C 写入字段名称
C
CALL UMACH (2, NOUT)
WRITE (NOUT,99998)
C
C 对常微分方程积分
C
IEND=0
10 CONTINUE
IEND = IEND + 1
TEND = IEND
C
CALL IVPAG (IDO, N, FCN, FCNJ, A, T, TEND, TOL, PARAM, Y)
IF (IEND .LE. 10) THEN
WRITE (NOUT,99999) T, Y
C
C 循环结束条件
C
IF (IEND .EQ. 10) IDO = 3
C
GO TO 10
C
END IF
C
C 制作输出表格
C
99998 FORMAT (11X, 'T', 11X, 'Y(1)', 11X, 'Y(2)', 11X, 'Y(3)')
99999 FORMAT (4F15.5)
END
C
SUBROUTINE FCN (N, X, Y, YPRIME)
C
C 宣告参数
C
INTEGER N
REAL X, Y(N), YPRIME(N)
C
YPRIME(1) = Y(2)*Y(3)
YPRIME(2) = -Y(1)*Y(3)
YPRIME(3) = -0.51*Y(1)*Y(2)
RETURN
END
C
SUBROUTINE FCNJ (N, X, Y, DYPDY)
C 宣告参数
INTEGER N
REAL X, Y(N), DYPDY(N,*)
C
RETURN
END
在compile都没事.开始build的时候