PROGRAM SETMBC2
C=======================================================================
C     DATA GENERATING PROGRAM FOR FEM4Q.FOR
C     DOMAIN : TLX BY TLY
C     BOUNDARY CONDITIONS
C        ON FACES OF -X AND -Y: T(X,Y)=0.
C        ON FACE OF +X: DT/DN = 0
C        ON FACE OF +Y: T(X,Y) = SIN ( 0.25*PI*X )
C     EIJI FUKUMORI NOVEMBER 28, 1995 
C=======================================================================
      IMPLICIT REAL*8 ( A-H , O-Z )
      PARAMETER (ND=4,MXE=4000, MXN=5000,MXB=1000,MXG=100,MXGEACH=100)
      PARAMETER ( TLX = 10., TLY = 5., MXX=100 )
C=======================================================================
      DIMENSION NODEX(MXE,ND),XCOORD(MXN),YCOORD(MXN),IBNDT(MXB),
     *          BVT(MXB), ITYPET(MXB),NMOVE(MXG), MOVEND(MXG,MXGEACH)
      DIMENSION XBOTTOM(MXX), XTOP(MXX)
      CHARACTER PROJECT*12,EXFILE*3
      LOGICAL YES
C=======================================================================
      DATA PROJECT / 'MOVINGBC.DAT' /
C=======================================================================
C         NEY: NEMBER OF VERTICAL ELEMENTS (NUMBER OF NODES: NEY+1)
C         NEX: NEMBER OF HORIZONTAL ELEMENTS (NUMBER OF NODES: NEX+1)
C=======================================================================
      MAXIT=700
      ERALLOW=0.0005
      JPOWER = 2
C
      NEX = 30
      NEY = 10
C
      DX = TLX / NEX
      DY = TLY / NEY
      NDX=NEX+1
      NDY=NEY+1
       EXX=1.
       EXY = 0.
       EYY=1.
C
C---------------NDMAXH2 = # OF NODES IN Y DOWN-STREAM HEAD IS ASSIGNED
      H1 = TLY
      NDMAXH2 = 0.25*NEY+1
C=======================================================================
C                    ELEMENT CREATION
      NE = 0
      DO I = 1 , NEX
      DO J = 1 , NEY
      NE = NE + 1
      IF ( NE .GT. MXE ) STOP 'NE > MXE'
      NODEX(NE,1) = NDY*(I-1) + J
      NODEX(NE,2) = NODEX(NE,1) + NDY
      NODEX(NE,3) = NODEX(NE,2) + 1
      NODEX(NE,4) = NODEX(NE,3) - NDY
      END DO
      END DO
C=======================================================================
C                 NODAL COORDINATE CREATION
C------XCOORDINATE VALUES ALONG BOTTOM AND TOP LINES
      TLXTOP = TLX*0.25
      DXTOP = TLXTOP/NEX
      DO I = 1 , NDX
      XBOTTOM(I) = DX*(I-1)
      XTOP(I) = TLX/2-TLXTOP/2+DXTOP*(I-1)
      END DO
C
      NNODE = 0
      DO I = 1 , NDX
      DO J = 1 , NDY
      NNODE = NNODE + 1
      IF ( NNODE .GT. MXN ) STOP 'NNODE > MXN'
      NODE = NDY*(I-1) + J
      XCOORD(NODE) =  F1(XBOTTOM(I),XTOP(I),NDY,J)
      YCOORD(NODE) =  DY*(J-1)
      END DO
      END DO
C=======================================================================
      NBT  = 0
C--------- FACE OF -X
      DO J = 1 , NDY
      NBT = NBT + 1
      IBNDT(NBT) = J
      ITYPET(NBT) = 1
      BVT(NBT) = H1          
      END DO
C--------- FACE OF -Y
C   NO FLOW; Q=0
C--------- FACE OF +X
C== HYDRAULIC HEAD = CONSTANT NODE
      H2 = YCOORD(NDY*NEX+NDMAXH2)
      DO J = 1 , NDMAXH2
      NBT  = NBT  + 1
      IBNDT (NBT) = NDY*NEX + J
      ITYPET(NBT) = 1
      BVT(NBT ) = H2
      END DO
C== SEEPAGE FACE <---------------- ITYPET(NBT) = 3 -----------------
      DO J = NDMAXH2+1, NEY
      NBT  = NBT  + 1
      IBNDT (NBT) = NDY*NEX + J
      ITYPET(NBT) = 3
      BVT(NBT ) = YCOORD (IBNDT(NBT))
      END DO
C--------- FACE OF +Y
C   NO FLOW; Q=0 AND FREE SURFACE
C
C---------- MOVING NODES (GROUPING)
C*******************************************************************
      NGROUP = 0
      DO I = 2 , NDX
      II = I-1
      NGROUP = NGROUP + 1
      NMOVEI = 0
      ISTART = 1
      IF ( I .EQ. NDX ) ISTART = NDMAXH2
      IF ( I .EQ. NDX-1 ) ISTART = NDMAXH2-2
      IF ( I .EQ. NDX-2 ) ISTART = NDMAXH2-3
      IF ( ISTART .LE. 0 ) ISTART = 1
      DO J = ISTART, NDY
      NMOVEI = NMOVEI + 1
      JJ = J - ISTART+1
      MOVEND(II,JJ) = NDY*(I-1)+J
      END DO
      NMOVE(II) = NMOVEI
      END DO
C=======================================================================
      WRITE (*,*) ' NUMBER OF ELEMENTS (NE) = ',NE
      WRITE (*,*) ' NUMBER OF NODAL POINTS (NNODE) = ',NNODE
C=======================================================================
C                      DATA FILE INQUIRY
      EXFILE = 'NEW'
      INQUIRE ( FILE = PROJECT, EXIST = YES )
      IF ( YES ) EXFILE='OLD'
C=======================================================================
C                     MAKING DATA FILES
C---------- 'PROJECT'.JNK
      OPEN ( 1, FILE=PROJECT, STATUS = EXFILE )
C---------- MATERIAL PARAMETER
      WRITE (1,*) EXX, EXY, EYY
C---------- CONVERGENCE PARAMETER
      WRITE (1,*) MAXIT,ERALLOW, JPOWER
C---------- ELEMENT INFORMATION
      WRITE (1,*) NE
      DO IEL = 1 , NE
        WRITE (1,'(10I5)') IEL,(NODEX(IEL,J),J=1,ND)
      END DO
C---------- NODAL COORDINATE INFOMATION
      WRITE (1,*) NNODE
      DO NODE = 1 , NNODE
        WRITE (1,*) NODE,XCOORD(NODE),YCOORD(NODE)
      END DO
C---------- BOUNDARY CONDITIONS
      WRITE (1,*) NBT
      DO I = 1 , NBT
        WRITE (1,*) IBNDT(I), ITYPET(I), BVT(I)
      END DO
C---------- MOVING NODE INFORMATION
      WRITE (1,*) NGROUP
      DO I = 1 , NGROUP
      NNMOVE = NMOVE(I)
      WRITE (1,*) NNMOVE, (MOVEND(I,J), J=1,NNMOVE)
      END DO
      CLOSE (1)
      STOP
      END
C
C
      FUNCTION F1(XBOTTOM,XTOP,NDY,J)
      IMPLICIT REAL*8 ( A-H , O-Z )
      T = DFLOAT(J-1)/DFLOAT(NDY-1)
      F1 = (1.-T)*XBOTTOM + T*XTOP
      RETURN
      END