1      subroutine argos_cafe_lww(lself,lpbc,xwm,idt,iwfr,iwto,jwfr,jwto,
2     + nlocw,mpairs,npairs,lwwjpt,lwwin,lwwj,list,rwx,rw)
3c $Id$
4      implicit none
5c
6#include "argos_cafe_common.fh"
7#include "bitops.fh"
8#include "mafdecls.fh"
9c
10      integer nlocw
11      real*8 rwx(mscr,3),rw(mscr)
12      real*8 xwm(mwm,3)
13      integer lwwj(*),lwwjpt(nlocw,2)
14      integer lwwin(nlocw,2)
15      integer list(mwm,2),idt(mwm)
16      logical lself,lpbc
17c
18      integer iwm,jwm,iwmfr,iwfr,iwmto,iwto,jwmfr,jwfr,jwmto,jwto
19      integer ilist,nlist,mpairs,npairs
20c
21      iwmfr=iwfr
22      if(lself) then
23      lwwin(nlocw,1)=0
24      lwwin(nlocw,2)=0
25      lwwjpt(nlocw,1)=0
26      lwwjpt(nlocw,2)=0
27      iwmto=iwto-1
28      jwmfr=0
29      jwmto=iwto
30      else
31      iwmto=iwto
32      jwmfr=jwfr
33      jwmto=jwto
34      endif
35c
36      npairs=0
37      do 1 iwm=iwmfr,iwmto
38      if(lself) jwmfr=iwm+1
39c
40      do 3 jwm=jwmfr,jwmto
41      rwx(jwm,1)=xwm(iwm,1)-xwm(jwm,1)
42      rwx(jwm,2)=xwm(iwm,2)-xwm(jwm,2)
43      rwx(jwm,3)=xwm(iwm,3)-xwm(jwm,3)
44    3 continue
45      if(lpbc) call argos_cafe_pbc(1,rwx,mscr,rwx,mscr,0,jwmfr,jwmto)
46      do 8 jwm=jwmfr,jwmto
47      rw(jwm)=rwx(jwm,1)**2+rwx(jwm,2)**2+rwx(jwm,3)**2
48    8 continue
49      do 4 jwm=jwmfr,jwmto
50      list(jwm,1)=0
51      list(jwm,2)=0
52      if(rw(jwm).lt.rshrt2.and.
53     + (iand(idt(iwm),mdynam).eq.ldynam.or.
54     + iand(idt(jwm),mdynam).eq.ldynam)) list(jwm,1)=1
55      if(rw(jwm).lt.rrest2 .and.
56     + (iand(idt(iwm),mrestr).eq.lrestr.or.
57     + iand(idt(jwm),mrestr).eq.lrestr)) list(jwm,1)=1
58    4 continue
59      if(npww.eq.2) then
60      do 2 jwm=jwmfr,jwmto
61      if(rw(jwm).lt.rlong2.and.
62     + (iand(idt(iwm),mdynam).eq.ldynam.or.
63     + iand(idt(jwm),mdynam).eq.ldynam)) list(jwm,2)=1
64      if(list(jwm,1).eq.1) list(jwm,2)=0
65    2 continue
66      endif
67c
68c     short range pairlist
69c
70      nlist=0
71      do 9 jwm=jwmfr,jwmto
72      if(list(jwm,1).eq.1) then
73      nlist=nlist+1
74      list(nlist,1)=jwm
75      endif
76    9 continue
77      if(npairs+nlist.gt.mpairs)
78     + call md_abort('Insufficient memory for pairlist',0)
79      do 5 ilist=1,nlist
80      lwwj(npairs+ilist)=list(ilist,1)
81    5 continue
82      lwwjpt(iwm-iwmfr+1,1)=npairs+1
83      npairs=npairs+nlist
84      lsww=lsww+nlist
85      lwwin(iwm-iwmfr+1,1)=nlist
86c
87c     long range pairs
88c
89      nlist=0
90      if(npww.eq.2) then
91      do 12 jwm=jwmfr,jwmto
92      if(list(jwm,2).eq.1) then
93      nlist=nlist+1
94      list(nlist,2)=jwm
95      endif
96   12 continue
97      if(npairs+nlist.gt.mpairs)
98     + call md_abort('Insufficient memory for pairlist',0)
99      do 13 ilist=1,nlist
100      lwwj(npairs+ilist)=list(ilist,2)
101   13 continue
102      lwwjpt(iwm-iwmfr+1,2)=npairs+1
103      npairs=npairs+nlist
104      llww=llww+nlist
105      lwwin(iwm-iwmfr+1,2)=nlist
106      endif
107    1 continue
108c
109      return
110      end
111