C
C *** TRANSMISSION BOUNDARY ***
C
      SUBROUTINE TMB(OM,ML,NL,DNS,VS,THK,CC,NMODE
     *              ,MMD,ZZ,MD4,ND4,DV,DP,DK,WK,IWK)
C
C ** INPUT DATA **
C OM     : OMEGA
C ML     : MAIN NL
C NL     : NUMBER OF LAYERS
C DNS(I) : DENSITY OF THE Ith LAYER 
C VS (I) : VS OF THE Ith LAYER
C THK(I) : THICKNESS OF THE Ith LAYER
C CC(J)  : PHASE VELOCITY OF THE Jth LOVE MODE
C NMODE  : OUTPUT HIGHEST MODE
C MMD    : MAX MODE NUMBER IN MAIN DIMENSION
C ZZ(I)  : DEPTH OF Ith NODE FROM THE FREE SURFACE
C MD4    : MAIN ND4
C ND4    : NUMBER OF TOTAL NODE
C DV(I,J): MATRIX DV
C DP(I,J): MATRIX DP
C WK(I)  : WORK MATRIX
C IWK(I) : WORK MATRIX
C ** OUTPUT DATA **
C DK(I,J): MATRIX DK (STIFFNESS MATRIX DUE TO TRANSMISSION BOUNDARY)
C
      IMPLICIT REAL*8(A-H,O-Z)
      DIMENSION CC(0:MMD),DNS(ML),VS(ML),THK(ML),ZZ(MD4)
     *         ,WK(MD4),IWK(MD4) 
      COMPLEX*16 DV(MD4,MD4),DP(MD4,MD4),DK(MD4,MD4)
C
C *** PHASE VELOCITY ***
C
      NDC=100
      CALL PHAS(OM,ML,NL,DNS,VS,THK,NDC,MMD,CC,NMODE)
C
C *** DV AND DP ***
C
      CALL DVDP(OM,ML,NL,DNS,VS,THK,CC,NMODE
     *         ,MMD,ZZ,MD4,ND4,DV,DP)
C
C *** DK ***
C
      ESP=1.D-10
      CALL $ZIF1M(DV,ND4,MD4,ESP,RDET,IDET,WK,IWK,IER)
      IF(IER.NE.0)THEN
        WRITE(*,*)'FROM TMB; IER.NE.0, THUS STOP'
        WRITE(*,*)'IER=',IER
        STOP
      END IF
C
      DO 100 I=1,ND4
        DO 110 J=1,ND4
          DK(I,J)=(0.D0,0.D0)
  110   CONTINUE
  100 CONTINUE
C
      DO 200 I=1,ND4
        DO 210 J=1,ND4
          DO 220 K=1,ND4
            DK(I,J)=DK(I,J)-DV(I,K)*DP(K,J)
  220     CONTINUE
  210   CONTINUE
  200 CONTINUE
C
      RETURN
      END
