--- /dev/null
+ double precision function qwolynes(ilevel,jfrag)
+ implicit none
+ include 'DIMENSIONS'
+ include 'DIMENSIONS.ZSCOPT'
+ include 'DIMENSIONS.COMPAR'
+ include 'COMMON.IOUNITS'
+ include 'COMMON.COMPAR'
+ include 'COMMON.CHAIN'
+ include 'COMMON.INTERACT'
+ integer ilevel,jfrag
+ integer i,j,jl,k,l,il,kl,nl,np,ip,kp
+ integer nsep /3/
+ double precision dist
+ double precision qq,qqij,qqijCM,dij,d0ij,dijCM,d0ijCM
+ logical lprn /.false./
+ double precision sigm,x
+ sigm(x)=0.25d0*x
+c write (iout,*) "QWolyes: " jfrag",jfrag,
+c & " ilevel",ilevel
+ qq = 0.0d0
+ if (ilevel.eq.0) then
+ if (lprn) write (iout,*) "Q computed for whole molecule"
+ nl=0
+ do il=nnt+nsep,nct
+ do jl=nnt,il-nsep
+ dij=0.0d0
+ dijCM=0.0d0
+ d0ij=0.0d0
+ d0ijCM=0.0d0
+ qqij=0.0d0
+ qqijCM=0.0d0
+ nl=nl+1
+ d0ij=dsqrt((cref(1,jl)-cref(1,il))**2+
+ & (cref(2,jl)-cref(2,il))**2+
+ & (cref(3,jl)-cref(3,il))**2)
+ dij=dist(il,jl)
+ qqij = dexp(-0.5d0*((dij-d0ij)/(sigm(d0ij)))**2)
+ if (itype(il).ne.10 .or. itype(jl).ne.10) then
+ nl=nl+1
+ d0ijCM=dsqrt(
+ & (cref(1,jl+nres)-cref(1,il+nres))**2+
+ & (cref(2,jl+nres)-cref(2,il+nres))**2+
+ & (cref(3,jl+nres)-cref(3,il+nres))**2)
+ dijCM=dist(il+nres,jl+nres)
+ qqijCM = dexp(-0.5d0*((dijCM-d0ijCM)/(sigm(d0ijCM)))**2)
+ endif
+ qq = qq+qqij+qqijCM
+ if (lprn) then
+ write (iout,*) "il",il," jl",jl,
+ & " itype",itype(il),itype(jl)
+ write (iout,*)"d0ij",d0ij," dij",dij," d0ijCM",d0ijCM,
+ & " dijCM",dijCM," qqij",qqij," qqijCM",qqijCM
+ endif
+ enddo
+ enddo
+ qq = qq/nl
+ if (lprn) write (iout,*) "nl",nl," qq",qq
+ else if (ilevel.eq.1) then
+ if (lprn) write (iout,*) "Level",ilevel," fragment",jfrag
+ nl=0
+c write (iout,*) "nlist_frag",nlist_frag(jfrag)
+ do i=2,nlist_frag(jfrag)
+ do j=1,i-1
+ il=list_frag(i,jfrag)
+ jl=list_frag(j,jfrag)
+ if (iabs(il-jl).gt.nsep) then
+ dij=0.0d0
+ dijCM=0.0d0
+ d0ij=0.0d0
+ d0ijCM=0.0d0
+ qqij=0.0d0
+ qqijCM=0.0d0
+ nl=nl+1
+ d0ij=dsqrt((cref(1,jl)-cref(1,il))**2+
+ & (cref(2,jl)-cref(2,il))**2+
+ & (cref(3,jl)-cref(3,il))**2)
+ dij=dist(il,jl)
+ qqij = dexp(-0.5d0*((dij-d0ij)/(sigm(d0ij)))**2)
+ if (itype(il).ne.10 .or. itype(jl).ne.10) then
+ nl=nl+1
+ d0ijCM=dsqrt(
+ & (cref(1,jl+nres)-cref(1,il+nres))**2+
+ & (cref(2,jl+nres)-cref(2,il+nres))**2+
+ & (cref(3,jl+nres)-cref(3,il+nres))**2)
+ dijCM=dist(il+nres,jl+nres)
+ qqijCM = dexp(-0.5d0*((dijCM-d0ijCM)/(sigm(d0ijCM)))**2)
+ endif
+ qq = qq+qqij+qqijCM
+ if (lprn) then
+ write (iout,*) "i",i," j",j," il",il," jl",jl,
+ & " itype",itype(il),itype(jl)
+ write (iout,*)"d0ij",d0ij," dij",dij," d0ijCM",d0ijCM,
+ & " dijCM",dijCM," qqij",qqij," qqijCM",qqijCM
+ endif
+ endif
+ enddo
+ enddo
+ qq = qq/nl
+ if (lprn) write (iout,*) "nl",nl," qq",qq
+ else if (ilevel.eq.2) then
+ np=npiece(jfrag,ilevel)
+ nl=0
+ do i=2,np
+ ip=ipiece(i,jfrag,ilevel)
+ do j=1,nlist_frag(ip)
+ il=list_frag(j,ip)
+ do k=1,i-1
+ kp=ipiece(k,jfrag,ilevel)
+ do l=1,nlist_frag(kp)
+ kl=list_frag(l,kp)
+ if (iabs(kl-il).gt.nsep) then
+ nl=nl+1
+ dij=0.0d0
+ dijCM=0.0d0
+ d0ij=0.0d0
+ d0ijCM=0.0d0
+ qqij=0.0d0
+ qqijCM=0.0d0
+ d0ij=dsqrt((cref(1,kl)-cref(1,il))**2+
+ & (cref(2,kl)-cref(2,il))**2+
+ & (cref(3,kl)-cref(3,il))**2)
+ dij=dist(il,kl)
+ qqij = dexp(-0.5d0*((dij-d0ij)/(sigm(d0ij)))**2)
+ if (itype(il).ne.10 .or. itype(kl).ne.10) then
+ nl=nl+1
+ d0ijCM=dsqrt(
+ & (cref(1,kl+nres)-cref(1,il+nres))**2+
+ & (cref(2,kl+nres)-cref(2,il+nres))**2+
+ & (cref(3,kl+nres)-cref(3,il+nres))**2)
+ dijCM=dist(il+nres,kl+nres)
+ qqijCM = dexp(-0.5d0*((dijCM-d0ijCM)/
+ & (sigm(d0ijCM)))**2)
+ endif
+ qq = qq+qqij+qqijCM
+ if (lprn) then
+ write (iout,*) "i",i," j",j," k",k," l",l," il",il,
+ & " kl",kl," itype",itype(il),itype(kl)
+ write (iout,*) " d0ij",d0ij," dij",dij," d0ijCM",
+ & d0ijCM," dijCM",dijCM," qqij",qqij," qqijCM",qqijCM
+ endif
+ endif
+ enddo ! l
+ enddo ! k
+ enddo ! j
+ enddo ! i
+ qq = qq/nl
+ if (lprn) write (iout,*) "nl",nl," qq",qq
+ else
+ write (iout,*)"Error: Q can be computed only for level 1 and 2."
+ endif
+ qwolynes=1.0d0-qq
+ return
+ end
+c-------------------------------------------------------------------------------
+ subroutine fragment_list
+ implicit none
+ include 'DIMENSIONS'
+ include 'DIMENSIONS.ZSCOPT'
+ include 'DIMENSIONS.COMPAR'
+ include 'COMMON.IOUNITS'
+ include 'COMMON.COMPAR'
+ logical lprn /.true./
+ integer i,ilevel,j,k,jfrag
+ do jfrag=1,nfrag(1)
+ nlist_frag(jfrag)=0
+ do i=1,npiece(jfrag,1)
+ if (lprn) write (iout,*) "jfrag=",jfrag,
+ & "i=",i," fragment",ifrag(1,i,jfrag),
+ & ifrag(2,i,jfrag)
+ do j=ifrag(1,i,jfrag),ifrag(2,i,jfrag)
+ do k=1,nlist_frag(jfrag)
+ if (list_frag(k,jfrag).eq.j) goto 10
+ enddo
+ nlist_frag(jfrag)=nlist_frag(jfrag)+1
+ list_frag(nlist_frag(jfrag),jfrag)=j
+ enddo
+ 10 continue
+ enddo
+ enddo
+ write (iout,*) "Fragment list"
+ do j=1,nfrag(1)
+ write (iout,*)"Fragment",j," list",(list_frag(k,j),
+ & k=1,nlist_frag(j))
+ enddo
+ return
+ end
projects / unres.git / blobdiff