added source code
[unres.git] / source / cluster / wham / src-M / contact.f
1       subroutine contact(lprint,ncont,icont)
2       include 'DIMENSIONS'
3       include 'COMMON.IOUNITS'
4       include 'COMMON.CHAIN'
5       include 'COMMON.INTERACT'
6       include 'COMMON.FFIELD'
7       include 'COMMON.NAMES'
8       real*8 facont /1.569D0/  ! facont = (2/(1-sqrt(1-1/4)))**(1/6)
9       integer ncont,icont(2,maxcont)
10       logical lprint
11       ncont=0
12       kkk=3
13 c     print *,'nnt=',nnt,' nct=',nct
14       do i=nnt+kkk,nct
15         iti=itype(i)
16         do j=nnt,i-kkk
17           itj=itype(j)
18           if (ipot.ne.4) then
19 c           rcomp=sigmaii(iti,itj)+1.0D0
20             rcomp=facont*sigmaii(iti,itj)
21           else 
22 c           rcomp=sigma(iti,itj)+1.0D0
23             rcomp=facont*sigma(iti,itj)
24           endif
25 c         rcomp=6.5D0
26 c         print *,'rcomp=',rcomp,' dist=',dist(nres+i,nres+j)
27           if (dist(nres+i,nres+j).lt.rcomp) then
28             ncont=ncont+1
29             icont(1,ncont)=i
30             icont(2,ncont)=j
31           endif
32         enddo
33       enddo
34       if (lprint) then
35         write (iout,'(a)') 'Contact map:'
36         do i=1,ncont
37           i1=icont(1,i)
38           i2=icont(2,i)
39           it1=itype(i1)
40           it2=itype(i2)
41           write (iout,'(i3,2x,a,i4,2x,a,i4)') 
42      &     i,restyp(it1),i1,restyp(it2),i2 
43         enddo
44       endif
45       return
46       end
47 c----------------------------------------------------------------------------
48       double precision function contact_fract(ncont,ncont_ref,
49      &                                     icont,icont_ref)
50       include 'DIMENSIONS'
51       include 'COMMON.IOUNITS'
52       integer ncont,ncont_ref,icont(2,maxcont),icont_ref(2,maxcont)
53       nmatch=0
54 c     print *,'ncont=',ncont,' ncont_ref=',ncont_ref 
55 c     write (iout,'(20i4)') (icont_ref(1,i),i=1,ncont_ref)
56 c     write (iout,'(20i4)') (icont_ref(2,i),i=1,ncont_ref)
57 c     write (iout,'(20i4)') (icont(1,i),i=1,ncont)
58 c     write (iout,'(20i4)') (icont(2,i),i=1,ncont)
59       do i=1,ncont
60         do j=1,ncont_ref
61           if (icont(1,i).eq.icont_ref(1,j) .and. 
62      &        icont(2,i).eq.icont_ref(2,j)) nmatch=nmatch+1
63         enddo
64       enddo
65 c     print *,' nmatch=',nmatch
66 c     contact_fract=dfloat(nmatch)/dfloat(max0(ncont,ncont_ref))
67       contact_fract=dfloat(nmatch)/dfloat(ncont_ref)
68       return
69       end