*DECK PRO2P4
      SUBROUTINE INITAL
C...
C...  SUBROUTINE INITAL IS CALLED BY MAIN PROGRAM PRO1P3 TO DEFINE THE
C...  INITIAL CONDITION IN THE NUMOL SOLUTION OF EQUATIONS (2.11) TO
C...  (2.14)
C...
      COMMON/T/   TIME,    NSTOP,    NORUN
     1      /Y/  T(51)
     2      /F/ TT(51)
     3      /R/  R(51),   TR(51),  TRR(51),
     4              R0,       T0,       DR
     5      /I/      N
C...
C...  SURFACE TEMPERATURE
      T0=25.
C...
C...  RADIUS
      R0=1.0
C...
C...  GRID SPACING
      N=51
      DR=R0/FLOAT(N-1)
C...
C...  INITIAL CONDITION (2.12)
      DO 1 I=1,N
         R(I)=FLOAT(I-1)*DR
         T(I)=0.
1     CONTINUE
      RETURN
      END
      SUBROUTINE DERV
C...
C...  SUBROUTINE DERV IS CALLED BY INTEGRATOR RKF45 TO DEFINE THE
C...  ODES IN THE NUMOL SOLUTION OF EQUATIONS (2.11) TO (2.14)
C...
      COMMON/T/   TIME,    NSTOP,    NORUN
     1      /Y/  T(51)
     2      /F/ TT(51)
     3      /R/  R(51),   TR(51),  TRR(51),
     4              R0,       T0,       DR
     5      /I/      N
C...
C...  BOUNDARY CONDITION (2.14)
       T(N)=T0
      TT(N)=0.
C...
C...  DERIVATIVE TR
      CALL DSS002(0.,R0,N,T,TR)
C...
C...  BOUNDARY CONDITION (2.13)
      TR(1)=0.
C...
C...  DERIVATIVE TRR
      CALL DSS002(0.,R0,N,TR,TRR)
C...
C...  EQUATION (2.11)
      DO 2 I=1,N-1
C...
C...     R = 0 (EQUATIONS (2.11) AND (2.15))
         IF(I.EQ.1)THEN
            TT(I)=2.0*TRR(I)
C...
C...     R NE 0 (EQUATION (2.11))
         ELSE IF(I.NE.1)THEN
            TT(I)=TRR(I)+1.0/R(I)*TR(I)
         END IF
2     CONTINUE
      RETURN
      END
      SUBROUTINE PRINT(NI,NO)
C...
C...  SUBROUTINE PRINT IS CALLED BY MAIN PROGRAM PRO1P3 TO PRINT THE
C...  NUMOL SOLUTION OF EQUATIONS (2.11) TO (2.14)
C...
      COMMON/T/   TIME,    NSTOP,    NORUN
     1      /Y/  T(51)
     2      /F/ TT(51)
     3      /R/  R(51),   TR(51),  TRR(51),
     4              R0,       T0,       DR
     5      /I/      N
C...
C...  PRINT THE NUMOL SOLUTION TO EQUATIONS (2.11) TO (2.14)
      WRITE(NO,2)TIME,(T(I),I=1,N,10)
2     FORMAT(' TIME = ',F6.2,/,15X,'  R=0',5X,'R=0.2',5X,'R=0.4',5X,
     1                             'R=0.6',5X,'R=0.8',5X,'  R=1',/,
     2                         '    T(R,T)',6F10.3,/)
      RETURN
      END
