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