C     ***************************************************
C     **  1-D DIFFUSION EQUATION (DU/DT+D/DX(DU/DX)=0) **
C     **                                               **
C     **  IMPLICIT METHOD                              **
C     **  EULER BACKWARD METHOD                        **
C     **                                               **
C     **                                               **
C     **  1999/08/23  VERSION 1.0                      **
C     **                                               **
C     **  FLUID RESEARCHING INSTITUTE                  **
C     ***************************************************
C
C     ******************
C     ** MAIN ROUTINE **
C     ******************
C
C
         CALL CONST
C
         CALL INIT
C
         CALL SOLVE
C
C
         STOP
         END
C
C
C     *************************************************
C     ** SOLVE DIFFUSION EQ. BY TIME MARCHING METHOD **
C     *************************************************
C
      SUBROUTINE SOLVE
C
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
      PARAMETER (MX=256)
      COMMON /DIVIDE/ NXP1,NX,NXM1
      COMMON /TMCND/ DT,NSTEP
      COMMON /RESULT/ X(MX),U(MX)
C
C     ***** EXPLICIT TIME MARCHING LOOP *****
C
      DO 500 ISTEP=1,NSTEP
C
C
      WRITE(*,*)'           '
      WRITE(*,600) NSTEP,ISTEP
  600 FORMAT('  STEP=',2I8)
C
C
          CALL BOUND
C
          CALL INTEG
C
  500 CONTINUE
C
          CALL OUTPUT
C
C
      RETURN
      END
C
C     *****************************
C     ** SUBROUTINE CONSTANT SET **
C     *****************************
C
      SUBROUTINE CONST
C
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
      PARAMETER (MX=256)
      COMMON /DIVIDE/ NXP1,NX,NXM1
      COMMON /TMCND/ DT,NSTEP
      COMMON /LENGTH/ DX,DXI
C
C     ***** NUMBER OF DIVIDE *****
C
         NX=11
         DX=1.D0/DBLE(NX-1)
         DXI=1.D0/DX
C
         NXM1=NX-1
         NXP1=NX+1
C
         DT=2.4D0*DX**2
C
C     ***** TIME ITERATION *****
C
         NSTEP=2
C
C
      RETURN
      END
C
C     *****************************
C     ** SUBROUTINE INITIAL DATA **
C     *****************************
C
      SUBROUTINE INIT
C
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
      PARAMETER (MX=256)
      COMMON /DIVIDE/ NXP1,NX,NXM1
      COMMON /TMCND/ DT,NSTEP
      COMMON /LENGTH/ DX,DXI
      COMMON /RESULT/ X(MX),U(MX)
C
        PI=4.D0*DATAN(1.D0)
C
C     *****  INITIAL DISTRIBUTION  *****
C
         X(1)=0.D0
C
      DO 100 I=2,NX
C
      IF(I.LT.6) THEN
       U(I)=2.D0*X(I)
      ELSE
       U(I)=2.D0-2.D0*X(I)
      END IF
C
  100 CONTINUE
C
C
      RETURN
      END
C
C     ***********************
C     ** SUBROUTINE OUTPUT **
C     ***********************
C
      SUBROUTINE OUTPUT
C
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
      PARAMETER (MX=256)
      COMMON /DIVIDE/ NXP1,NX,NXM1
      COMMON /LENGTH/ DX,DXI
      COMMON /RESULT/ X(MX),U(MX)
C
C
      WRITE(2,600)(X(I),U(I),I=1,NX)
  600 FORMAT(2F12.7)
C
C
      RETURN
      END
C
C     ***********************************
C     ** SUBROUTINE BOUNDARY CONDITION **
C     ***********************************
C
      SUBROUTINE BOUND
C
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
      PARAMETER (MX=256)
      COMMON /DIVIDE/ NXP1,NX,NXM1
      COMMON /RESULT/ X(MX),U(MX)
C
C
        U(NX)=0.D0
        U( 1)=0.D0
C
C
      RETURN
      END
C
C     *******************************
C     ** SOLVED DIFFUSION EQUATION **
C     *******************************
C
      SUBROUTINE INTEG
C
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
      PARAMETER (MX=256)
      COMMON /DIVIDE/ NXP1,NX,NXM1
      COMMON /RESULT/ X(MX),U(MX)
      COMMON /TMCND/ DT,NSTEP
      COMMON /LENGTH/ DX,DXI
      DIMENSION UO(MX)
C
C
      DO 100 I=1,NX
C
  100    UO(I) = U(I)
C
        COEF=1.D0+2.D0*DT*DXI*DXI
C
C     ***** COMPUTE EQUATION *****
C
      DO 300 LOOP=1,100
C
         ERR=0.D0
C
      DO 120 I=2,NXM1
C
        DIFF=(U(I+1)+U(I-1))*DXI*DXI
C
        DU=(UO(I)+DT*DIFF)/COEF-U(I)
C
        ERR=ERR+DU*DU
C
        U(I)=U(I)+1.D0*DU
C
  120 CONTINUE
C
        CALL BOUND
C
        ERR=DSQRT(ERR/DBLE(NX-2))
        IF(ERR.LT.1.D-5) GOTO 99
C
  300 CONTINUE
C
   99 CONTINUE
C
      WRITE(*,*) LOOP,ERR
C
C
      RETURN
      END