      PROGRAM H3SS2
C...
C...  H3 - HEAT CONDUCTION IN A CIRCULAR FIN, STEADY STATE
C...       SOLUTION, VERSION 2
C...
C...  DOUBLE PRECISION CODING IS USED
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
C...
C...  COMMON AREA WITH PROBLEM PARAMETERS
      PARAMETER(NR=6)
      COMMON/P/  ALPHA,   BETA,     UA,     U0,  U(NR),
     +              DR,     R0,     R1,  R(NR),      B,  
     +           RHOCP
C...
C...  OPEN AN OUTPUT FILE
      NO=6
      OPEN(NO,FILE='OUTPUT', STATUS='UNKNOWN')
C...
C...  PROBLEM PARAMETERS
C...
C...     FIN INSIDE RADIUS
         R0=0.5D0
C...
C...     FIN OUTSIDE RADIUS
         R1=2.0D0
C...
C...     FIN THICKNESS
         B=0.2D0
C...
C...     HEAT TRANSFER COEFFICIENT
         H=0.01D0
C...
C...     RHO*CP
         RHOCP=2.0D0
C...
C...     THERMAL DIFFUSIVITY
         ALPHA=0.01D0
C...
C...     AMBIENT TEMPERATURE
         UA=25.0D0
C...
C...     INSIDE TEMPERATURE
         U0=150.0D0
C...
C...  RECIPROCAL TIME FOR HEAT TRANSFER
      BETA=(2.0D0*H)/(B*RHOCP)
C...
C...  INCREMENT FOR RADIAL GRID
      DR=(R1-R0)/DFLOAT(NR-1)
C...
C...  RADIAL GRID
      DO 1 I=1,NR
         R(I)=R0+DFLOAT(I-1)*DR
1     CONTINUE
C...
C...  STEADY STATE SOLUTION WITH DIRICHLET BOUNDARY CONDITION
C...
C...     CONSTANT X0
         X0=DSQRT(BETA/ALPHA)*R0
         X1=DSQRT(BETA/ALPHA)*R1
C...
C...     EVALUATE THE SOLUTION AT A SERIES OF RADIAL POSITIONS
         DO 2 I=1,NR
         X=DSQRT(BETA/ALPHA)*R(I)
C...
C...     ANALYTICAL SOLUTION
         U(I)=DSQRT(X0/X)*DSINH(X1-X)/DSINH(X1-X0)
C...
C...     CONVERSION FROM DIMENSIONLESS TO DIMENSIONAL TEMPERATURE
         U(I)=UA+(U0-UA)*U(I)
2        CONTINUE
C...
C...  PRINT THE SOLUTION
      WRITE(NO,3)(R(I),I=1,NR),(U(I),I=1,NR)
3     FORMAT(/,' Dirichlet boundary conditions',/,
     +         '      r ',6F9.2,/,'   u(r) ',6F9.3,/)
      STOP
      END
      DOUBLE PRECISION FUNCTION DFLOAT(I)
C...
C...  FUNCTION DFLOAT CONVERTS A SINGLE PRECISION INTEGER INTO A
C...  DOUBLE PRECISION FLOATING POINT.  THIS FUNCTION IS PROVIDED
C...  IN CASE THE USER'S FORTRAN COMPILER DOES NOT INCLUDE DFLOAT.
C...
      DFLOAT=DBLE(FLOAT(I))
      RETURN
      END

