1      subroutine argos_cafe_uhop(lpbc,lpbcs,lseq,isl,isgm,isga,
2     + isq3,xs,nsaloc,chg,lda,rda,uda)
3c
4      implicit none
5c
6#include "argos_cafe_common.fh"
7#include "bitops.fh"
8#include "util.fh"
9#include "mafdecls.fh"
10c
11      logical lpbc,lpbcs
12      real*8 xs(msa,3),chg(mqt,mqp,mset)
13      integer nsaloc
14c
15      integer isl(msa,mis2)
16      integer lda(16,*),isgm(msa),isga(msa),isq3(msa),lseq(mseq)
17      real*8 rda(11,*),uda(4,*)
18c
19      real*8 q,rwx(3),rwy(3),rw,rwi,uhop(4),facthp,cutoff,uh(4)
20      integer isa,jsa,iss,jss,ndxi
21      integer ilda,idon,idong,idons,idonp,iacc,iaccg,iaccs,iaccp
22      integer idown,iaown,idfr,idto,iafr,iato
23      integer jsprev,i
24      logical lskeep
25c
26      integer ndon,nacc
27c
28      facthp=1.0d0
29      cutoff=1.0d3
30      jsprev=0
31      lskeep=.false.
32c
33c      write(*,5555) stime
34c 5555 format('TIME ',f12.6)
35      do 1 ilda=1,nlda
36c
37      uda(1,ilda)=0.0d0
38      uda(2,ilda)=0.0d0
39      uda(3,ilda)=0.0d0
40      uda(4,ilda)=0.0d0
41c
42      uhop(1)=0.0d0
43      uhop(2)=0.0d0
44      uhop(3)=0.0d0
45      uhop(4)=0.0d0
46c
47      idon=lda(2,ilda)
48      idong=lda(1,ilda)
49      idons=lda(3,ilda)
50      idonp=lda(4,ilda)
51      iacc=lda(7,ilda)
52      iaccg=lda(6,ilda)
53      iaccs=lda(8,ilda)
54      iaccp=lda(9,ilda)
55c
56      idown=lda(11,ilda)
57      idfr=lda(12,ilda)
58      idto=lda(13,ilda)
59      iaown=lda(14,ilda)
60      iafr=lda(15,ilda)
61      iato=lda(16,ilda)
62c
63      ndon=idto-idfr+1
64      nacc=iato-iafr+1
65c
66c      write(*,'(i5,a,5i5)') me,' LDA d ',ilda,nlda,idown,idfr,idto
67c
68c     get donor and acceptor atom coordinates
69c
70      call argos_space_getda(idown,idfr,idto,isl,isga,isq3,xs,nsaloc)
71c      write(*,'(i5,a,5i5)') me,' LDA a ',ilda,nlda,iaown,iafr,iato
72      call argos_space_getda(iaown,iafr,iato,isl,isga,isq3,xs,
73     + nsaloc+ndon)
74c
75      do 2001 isa=1,ndon
76c      write(*,3001) isa,isga(nsaloc+isa),isq3(nsaloc+isa),
77c     + (xs(nsaloc+isa,jsa),jsa=1,3)
78c 3001 format(' donor',3i5,3f12.6)
79      isgm(nsaloc+isa)=idons
80 2001 continue
81      do 2002 isa=1,nacc
82c      write(*,3002) isa,isga(nsaloc+ndon+isa),isq3(nsaloc+ndon+isa),
83c     + (xs(nsaloc+ndon+isa,jsa),jsa=1,3)
84c 3002 format(' acceptor',3i5,3f12.6)
85      isgm(nsaloc+ndon+isa)=iaccs
86 2002 continue
87c
88c      write(*,'(a,8i5)') 'DA ',idon,idong,idons,idonp,
89c     + iacc,iaccg,iaccs,iaccp
90c
91c     loop over the atoms in the donor and acceptor residues
92c
93cx      do 2 isa=1,nsaloc
94      do 2 isa=nsaloc+1,nsaloc+ndon+nacc+1
95c
96      if(isgm(isa).eq.idons.or.isgm(isa).eq.iaccs) then
97c
98c     loop over all atoms not in either donor or acceptor residues
99c
100      do 3 jsa=1,nsaloc
101c
102      iss=isgm(isa)
103      jss=isgm(jsa)
104c
105c      write(*,'(4i5)') isa,jsa,iss,jss
106c
107      if(jsprev.ne.jss) then
108      if(lskeep) then
109      uda(1,ilda)=uda(1,ilda)+facthp*uhop(1)
110      uda(2,ilda)=uda(2,ilda)+facthp*uhop(2)
111      uda(3,ilda)=uda(3,ilda)+facthp*uhop(3)
112      uda(4,ilda)=uda(4,ilda)+facthp*uhop(4)
113      endif
114      uhop(1)=0.0d0
115      uhop(2)=0.0d0
116      uhop(3)=0.0d0
117      uhop(4)=0.0d0
118      jsprev=jss
119      lskeep=.false.
120      endif
121c
122      if(isgm(jsa).ne.idons.and.isgm(jsa).ne.iaccs) then
123c
124c     find the periodic boundary offset from the midpoint
125c
126c      write(*,3323) ilda,isa,xs(isa,1),xs(isa,2),xs(isa,3),
127c     + jsa,xs(jsa,1),xs(jsa,2),xs(jsa,3)
128c 3323 format(2i5,3f12.6,/,5x,i5,3f12.6)
129c
130      if(isgm(isa).eq.idons) then
131      rwy(1)=xs(jsa,1)-rda(6,ilda)
132      rwy(2)=xs(jsa,2)-rda(7,ilda)
133      rwy(3)=xs(jsa,3)-rda(8,ilda)
134      else
135      rwy(1)=xs(jsa,1)-rda(9,ilda)
136      rwy(2)=xs(jsa,2)-rda(10,ilda)
137      rwy(3)=xs(jsa,3)-rda(11,ilda)
138      endif
139c
140c      write(*,3335) rwy
141      if(lpbc.or.lpbcs) call argos_cafe_pbc(0,rwy,1,rwx,1,0,1,1)
142c
143c      write(*,3333) (xs(jsa,i),i=1,3)
144c 3333 format('xs(jsa,1:3)',3f12.6)
145c      write(*,3334) (rda(i,ilda),i=6,11)
146c 3334 format(3x,'rda(6:8)',3f12.6,/,'  rda(9:11)',3f12.6)
147c      write(*,3335) rwx
148c 3335 format(8x,'rwx',3f12.6)
149c
150c     evaluate the distance between atoms isa and jsa
151c
152      if(isga(isa).eq.iaccg) then
153      rwx(1)=-rwx(1)+xs(jsa,1)-rda(1,ilda)
154      rwx(2)=-rwx(2)+xs(jsa,2)-rda(2,ilda)
155      rwx(3)=-rwx(3)+xs(jsa,3)-rda(3,ilda)
156c      write(*,3337) (rda(i,ilda),i=1,3)
157c 3337 format(3x,'rda(1:3)',3f12.6)
158      else
159      rwx(1)=-rwx(1)+xs(jsa,1)-xs(isa,1)
160      rwx(2)=-rwx(2)+xs(jsa,2)-xs(isa,2)
161      rwx(3)=-rwx(3)+xs(jsa,3)-xs(isa,3)
162c      write(*,3336) (xs(isa,i),i=1,3)
163c 3336 format('xs(isa,1:3)',3f12.6)
164      endif
165      rw=sqrt(rwx(1)**2+rwx(2)**2+rwx(3)**2)
166      rwi=1.0d0/rw
167c
168      if(rw.le.cutoff) lskeep=.true.
169c
170c      write(*,'(3i5,f12.6)') ilda,isa,jsa,rw
171c
172c     donor energies
173c
174      uh(1)=0.0d0
175      uh(2)=0.0d0
176      uh(3)=0.0d0
177      uh(4)=0.0d0
178c
179      if(iss.eq.idons) then
180c
181c     before hopping
182c
183      q=chg(isq3(isa),3,lseq(iss))*chg(isq3(jsa),3,lseq(jss))
184      uh(1)=q*rwi
185      uhop(1)=uhop(1)+q*rwi
186c      write(*,6665) stime,1,isga(isa),chg(isq3(isa),3,lseq(iss)),
187c     + isga(jsa),chg(isq3(jsa),3,lseq(jss)),uh(1)
188c 6665 format(f12.6,'q',2i5,f12.6,i5,f12.6,e12.5)
189c
190c     after hopping - skip if donor site
191c
192      if(isga(isa).ne.idong) then
193      ndxi=isa
194      if(isga(isa).eq.idong) ndxi=isa+lda(5,ilda)
195      if(isga(isa).eq.idong+lda(5,ilda)) ndxi=isa-lda(5,ilda)
196      q=chg(isq3(ndxi),3,idonp)*chg(isq3(jsa),3,lseq(jss))
197      uh(3)=q*rwi
198      uhop(3)=uhop(3)+q*rwi
199c      write(*,6665) stime,2,isga(isa),chg(isq3(ndxi),3,idonp),
200c     + isga(jsa),chg(isq3(jsa),3,lseq(jss)),uh(3)
201      endif
202      endif
203c
204c     acceptor energies
205c
206      if(iss.eq.iaccs) then
207c
208c     before hopping - skip if acceptor site
209c
210      if(isga(isa).ne.iaccg) then
211      q=chg(isq3(isa),3,lseq(iss))*chg(isq3(jsa),3,lseq(jss))
212      uh(2)=q*rwi
213      uhop(2)=uhop(2)+q*rwi
214c      write(*,6665) stime,3,isga(isa),chg(isq3(isa),3,lseq(iss)),
215c     + isga(jsa),chg(isq3(jsa),3,lseq(jss)),uh(2)
216      endif
217c
218c     after hopping
219c
220      ndxi=isa
221      if(isga(isa).eq.iaccg) ndxi=isa+lda(10,ilda)
222      if(isga(isa).eq.iaccg+lda(10,ilda)) ndxi=isa-lda(10,ilda)
223      q=chg(isq3(ndxi),3,iaccp)*chg(isq3(jsa),3,lseq(isgm(jsa)))
224      uh(4)=q*rwi
225      uhop(4)=uhop(4)+q*rwi
226c      write(*,6665) stime,4,isga(isa),chg(isq3(ndxi),3,iaccp),
227c     + isga(jsa),chg(isq3(jsa),3,lseq(isgm(jsa))),uh(4)
228      endif
229c
230c      write(*,6666) isga(isa),iss,isga(jsa),jss,lskeep,rw,uh
231c 6666 format(4i5,4x,l1,f12.6,4e12.5)
232c
233      endif
234c
235    3 continue
236c
237      endif
238c
239    2 continue
240c
241      if(lskeep) then
242      uda(1,ilda)=uda(1,ilda)+facthp*uhop(1)
243      uda(2,ilda)=uda(2,ilda)+facthp*uhop(2)
244      uda(3,ilda)=uda(3,ilda)+facthp*uhop(3)
245      uda(4,ilda)=uda(4,ilda)+facthp*uhop(4)
246      endif
247c
248c      write(*,'(a,8f12.3)') 'uhop ',(uda(ndxi,ilda),ndxi=1,4)
249c
250c
251    1 continue
252c
253      return
254      end
255c $Id$
256