c
c @(#)bndtphys.f	1.3
c
      subroutine bndtphys
c
c
c
      include'mladv.inc'
      include'trc.inc'

      common/parm/kz,m,mm1,mm2,mp1,msq,tim,tmax,tstart,hy,rea,dt,
     +   xbasin,yy(mx0),icalc,kzp1,kzm1,iftop,ifbot,iftvv,ifrel,theta,
     +   ifsbl,ifeva

      if (itphys.eq.0) return

      if ((itphys .eq. 2) .or. (itphys.eq. 3)) then

         icurtrace = 1

         do 5 ibnd=1,mbnd

            call findkml(bnddpth(ibnd),zml,l,kml)

            if (itphys .eq. 2) then
               krem = kml
              else
               krem = keuph
             endif

            do 10 k=1,krem

               ip=ipofb(ibnd)

               xn = bndtrace(ibnd,k,1) + trcqgofk(ip,k)

               if (xn .lt. 0.0) then
                   xn = 0.0
                 endif

               if (k .eq. kml) then
                   z1 = zml(kml) + dzml(1) / 2.0
                   dz1 = bnddpth(ibnd) + z1
                   xn = xn * dz1 / dzml(1)
                 endif

               if (k .eq. kml+1) then
                  z1 = -bnddpth(ibnd) - zml(kml+1)
                  dz1 = z1 + dzml(1) / 2.0
                  xn = xn * dz1 / dzml(1)
                endif

               bndtrace(ibnd,l+1,1) = bndtrace(ibnd,l+1,1) + xn
     
               bndtrace(ibnd,k,1) = -trcqgofk(ip,k)

   10        continue

    5     continue

        else
          write(*,*)'Invalid biological boundary condition model'
          stop
       endif

      return
      end
