      subroutine beckex(dnstya,dnstyb,dnstyt,gmoda,gmodb,dlapa,dlapb,
     &                  xcvta,xcvtb,ex,vxa,vxb)

c  this subroutine written by alain st-amant of the department of
c  pharmaceutical chemistry, university of california, san francisco.
c  all rights reserved.  this is part of the DeFT project.

      implicit real*8(a-h,o-z)

      data tolerance,beta/1.0d-15,0.0042d0/

      data zero,one,two/0.0d0,1.0d0,2.0d0/
      data three,four,six/3.0d0,4.0d0,6.0d0/

      data third,fourthirds/0.33333333333333d0,1.33333333333333d0/

      if(dnstya.gt.tolerance) then
                                    dka=gmoda/dnstya**fourthirds
                              else
                                    dka=zero
                              endif

      if(dnstyb.gt.tolerance) then
                                    dkb=gmodb/dnstyb**fourthirds
                              else
                                    dkb=zero
                              endif

      sinha=log(dka+sqrt(one+dka*dka))
      sinhb=log(dkb+sqrt(one+dkb*dkb))

      denoma=one+six*beta*dka*sinha
      denomb=one+six*beta*dkb*sinhb

      dsinha=one/sqrt(one+dka*dka)
      dsinhb=one/sqrt(one+dkb*dkb)

      eterma=fourthirds*dka*dka*denoma*dnstya**third
      etermb=fourthirds*dkb*dkb*denomb*dnstyb**third

      vterma=one+three*beta*dka*(sinha-dka*dsinha)
      vtermb=one+three*beta*dkb*(sinhb-dkb*dsinhb)

      if(dnstya.gt.tolerance) then
                                    vterma=two*vterma*dlapa
     &                                    /dnstya**fourthirds

                                    pterma=six*beta*(dnstya**third)
     &                                     *((xcvta/dnstya**3)
     &                                       -fourthirds*dka**3)

                                    ex=ex-beta*(dnstya**fourthirds)
     &                                        *dka*dka/denoma/dnstyt
                              else
                                    vterma=zero
                                    pterma=zero
                              endif

      if(dnstyb.gt.tolerance) then
                                    vtermb=two*vtermb*dlapb
     &                                    /dnstyb**fourthirds

                                    ptermb=six*beta*(dnstyb**third)
     &                                     *((xcvtb/dnstyb**3)
     &                                       -fourthirds*dkb**3)

                                    ex=ex-beta*(dnstyb**fourthirds)
     &                                        *dkb*dkb/denomb/dnstyt
                              else
                                    vtermb=zero
                                    ptermb=zero
                              endif

      factora=(three*sinha*(one+two*beta*dka*sinha)
     &        +four*dka*dsinha*(one-three*beta*dka*dka*dsinha))/denoma
     &        +dka*(dsinha**3)
      factorb=(three*sinhb*(one+two*beta*dkb*sinhb)
     &        +four*dkb*dsinhb*(one-three*beta*dkb*dkb*dsinhb))/denomb
     &        +dkb*(dsinhb**3)

      factora=eterma-vterma+pterma*factora
      factorb=etermb-vtermb+ptermb*factorb

      vxa=vxa-beta*factora/(denoma*denoma)
      vxb=vxb-beta*factorb/(denomb*denomb)

      return
      end
