X-Git-Url: http://mmka.chem.univ.gda.pl/gitweb/?a=blobdiff_plain;f=source%2Fcluster%2Fwham%2Fsrc%2Fenergy_p_new.F;fp=source%2Fcluster%2Fwham%2Fsrc%2Fenergy_p_new.F;h=bf9c56307ea89c66406b2c68a86261bcd5cc2799;hb=0a11a2c4ccee14ed99ae44f2565b270ba8d4bbb6;hp=1b69eb358ecee26b9d81ee9a66a11d0772cb2f57;hpb=5eb407964903815242c59de10960f42761139e10;p=unres.git diff --git a/source/cluster/wham/src/energy_p_new.F b/source/cluster/wham/src/energy_p_new.F index 1b69eb3..bf9c563 100644 --- a/source/cluster/wham/src/energy_p_new.F +++ b/source/cluster/wham/src/energy_p_new.F @@ -336,7 +336,7 @@ C implicit real*8 (a-h,o-z) include 'DIMENSIONS' include 'sizesclu.dat' - include "DIMENSIONS.COMPAR" +c include "DIMENSIONS.COMPAR" parameter (accur=1.0d-10) include 'COMMON.GEO' include 'COMMON.VAR' @@ -497,7 +497,7 @@ C implicit real*8 (a-h,o-z) include 'DIMENSIONS' include 'sizesclu.dat' - include "DIMENSIONS.COMPAR" +c include "DIMENSIONS.COMPAR" include 'COMMON.GEO' include 'COMMON.VAR' include 'COMMON.LOCAL' @@ -586,7 +586,7 @@ C implicit real*8 (a-h,o-z) include 'DIMENSIONS' include 'sizesclu.dat' - include "DIMENSIONS.COMPAR" +c include "DIMENSIONS.COMPAR" include 'COMMON.GEO' include 'COMMON.VAR' include 'COMMON.LOCAL' @@ -713,7 +713,7 @@ C implicit real*8 (a-h,o-z) include 'DIMENSIONS' include 'sizesclu.dat' - include "DIMENSIONS.COMPAR" +c include "DIMENSIONS.COMPAR" include 'COMMON.GEO' include 'COMMON.VAR' include 'COMMON.LOCAL' @@ -844,7 +844,7 @@ C implicit real*8 (a-h,o-z) include 'DIMENSIONS' include 'sizesclu.dat' - include "DIMENSIONS.COMPAR" +c include "DIMENSIONS.COMPAR" include 'COMMON.GEO' include 'COMMON.VAR' include 'COMMON.LOCAL' @@ -2794,16 +2794,16 @@ C Evaluate bridge-strain energy and its gradient in virtual-bond and SC vectors. C implicit real*8 (a-h,o-z) include 'DIMENSIONS' - include 'sizesclu.dat' include 'COMMON.SBRIDGE' include 'COMMON.CHAIN' include 'COMMON.DERIV' include 'COMMON.VAR' include 'COMMON.INTERACT' + include 'COMMON.IOUNITS' dimension ggg(3) ehpb=0.0D0 -cd print *,'edis: nhpb=',nhpb,' fbr=',fbr -cd print *,'link_start=',link_start,' link_end=',link_end +cd write(iout,*)'edis: nhpb=',nhpb,' fbr=',fbr +cd write(iout,*)'link_start=',link_start,' link_end=',link_end if (link_end.eq.0) return do i=link_start,link_end C If ihpb(i) and jhpb(i) > NRES, this is a SC-SC distance, otherwise a @@ -2818,43 +2818,85 @@ C iii and jjj point to the residues for which the distance is assigned. iii=ii jjj=jj endif +c write (iout,*) "i",i," ii",ii," iii",iii," jj",jj," jjj",jjj, +c & dhpb(i),dhpb1(i),forcon(i) C 24/11/03 AL: SS bridges handled separately because of introducing a specific C distance and angle dependent SS bond potential. if (ii.gt.nres .and. itype(iii).eq.1 .and. itype(jjj).eq.1) then call ssbond_ene(iii,jjj,eij) ehpb=ehpb+2*eij +cd write (iout,*) "eij",eij + else if (ii.gt.nres .and. jj.gt.nres) then +c Restraints from contact prediction + dd=dist(ii,jj) + if (dhpb1(i).gt.0.0d0) then + ehpb=ehpb+2*forcon(i)*gnmr1(dd,dhpb(i),dhpb1(i)) + fac=forcon(i)*gnmr1prim(dd,dhpb(i),dhpb1(i))/dd +c write (iout,*) "beta nmr", +c & dd,2*forcon(i)*gnmr1(dd,dhpb(i),dhpb1(i)) + else + dd=dist(ii,jj) + rdis=dd-dhpb(i) +C Get the force constant corresponding to this distance. + waga=forcon(i) +C Calculate the contribution to energy. + ehpb=ehpb+waga*rdis*rdis +c write (iout,*) "beta reg",dd,waga*rdis*rdis +C +C Evaluate gradient. +C + fac=waga*rdis/dd + endif + do j=1,3 + ggg(j)=fac*(c(j,jj)-c(j,ii)) + enddo + do j=1,3 + ghpbx(j,iii)=ghpbx(j,iii)-ggg(j) + ghpbx(j,jjj)=ghpbx(j,jjj)+ggg(j) + enddo + do k=1,3 + ghpbc(k,jjj)=ghpbc(k,jjj)+ggg(k) + ghpbc(k,iii)=ghpbc(k,iii)-ggg(k) + enddo else C Calculate the distance between the two points and its difference from the C target distance. - dd=dist(ii,jj) - rdis=dd-dhpb(i) + dd=dist(ii,jj) + if (dhpb1(i).gt.0.0d0) then + ehpb=ehpb+2*forcon(i)*gnmr1(dd,dhpb(i),dhpb1(i)) + fac=forcon(i)*gnmr1prim(dd,dhpb(i),dhpb1(i))/dd +c write (iout,*) "alph nmr", +c & dd,2*forcon(i)*gnmr1(dd,dhpb(i),dhpb1(i)) + else + rdis=dd-dhpb(i) C Get the force constant corresponding to this distance. - waga=forcon(i) + waga=forcon(i) C Calculate the contribution to energy. - ehpb=ehpb+waga*rdis*rdis + ehpb=ehpb+waga*rdis*rdis +c write (iout,*) "alpha reg",dd,waga*rdis*rdis C C Evaluate gradient. C - fac=waga*rdis/dd + fac=waga*rdis/dd + endif cd print *,'i=',i,' ii=',ii,' jj=',jj,' dhpb=',dhpb(i),' dd=',dd, cd & ' waga=',waga,' fac=',fac - do j=1,3 - ggg(j)=fac*(c(j,jj)-c(j,ii)) - enddo + do j=1,3 + ggg(j)=fac*(c(j,jj)-c(j,ii)) + enddo cd print '(i3,3(1pe14.5))',i,(ggg(j),j=1,3) C If this is a SC-SC distance, we need to calculate the contributions to the C Cartesian gradient in the SC vectors (ghpbx). - if (iii.lt.ii) then + if (iii.lt.ii) then do j=1,3 ghpbx(j,iii)=ghpbx(j,iii)-ggg(j) ghpbx(j,jjj)=ghpbx(j,jjj)+ggg(j) enddo - endif - do j=iii,jjj-1 + endif do k=1,3 - ghpbc(k,j)=ghpbc(k,j)+ggg(k) + ghpbc(k,jjj)=ghpbc(k,jjj)+ggg(k) + ghpbc(k,iii)=ghpbc(k,iii)-ggg(k) enddo - enddo endif enddo ehpb=0.5D0*ehpb @@ -4207,7 +4249,7 @@ c write (iout,*) 'i=',i,' gloc=',gloc(i-3,icg) do i=1,ndih_constr itori=idih_constr(i) phii=phi(itori) - difi=phii-phi0(i) + difi=pinorm(phii-phi0(i)) if (difi.gt.drange(i)) then difi=difi-drange(i) edihcnstr=edihcnstr+0.25d0*ftors*difi**4 @@ -4217,10 +4259,10 @@ c write (iout,*) 'i=',i,' gloc=',gloc(i-3,icg) edihcnstr=edihcnstr+0.25d0*ftors*difi**4 gloc(itori-3,icg)=gloc(itori-3,icg)+ftors*difi**3 endif -! write (iout,'(2i5,2f8.3,2e14.5)') i,itori,rad2deg*phii, -! & rad2deg*difi,0.25d0*ftors*difi**4,gloc(itori-3,icg) +c write (iout,'(2i5,2f8.3,2e14.5)') i,itori,rad2deg*phii, +c & rad2deg*difi,0.25d0*ftors*difi**4,gloc(itori-3,icg) enddo -! write (iout,*) 'edihcnstr',edihcnstr + write (iout,*) 'edihcnstr',edihcnstr return end c------------------------------------------------------------------------------ @@ -4289,24 +4331,26 @@ c write (iout,*) 'i=',i,' gloc=',gloc(i-3,icg) enddo ! 6/20/98 - dihedral angle constraints edihcnstr=0.0d0 +c write (iout,*) "Dihedral angle restraint energy" do i=1,ndih_constr - print *,"i",i itori=idih_constr(i) phii=phi(itori) - difi=phii-phi0(i) + difi=pinorm(phii-phi0(i)) +c write (iout,'(2i5,4f8.3,2e14.5)') i,itori,rad2deg*phii, +c & rad2deg*difi,rad2deg*phi0(i),rad2deg*drange(i) if (difi.gt.drange(i)) then difi=difi-drange(i) edihcnstr=edihcnstr+0.25d0*ftors*difi**4 gloc(itori-3,icg)=gloc(itori-3,icg)+ftors*difi**3 +c write (iout,*) 0.25d0*ftors*difi**4 else if (difi.lt.-drange(i)) then difi=difi+drange(i) edihcnstr=edihcnstr+0.25d0*ftors*difi**4 gloc(itori-3,icg)=gloc(itori-3,icg)+ftors*difi**3 +c write (iout,*) 0.25d0*ftors*difi**4 endif -! write (iout,'(2i5,2f8.3,2e14.5)') i,itori,rad2deg*phii, -! & rad2deg*difi,0.25d0*ftors*difi**4,gloc(itori-3,icg) enddo -! write (iout,*) 'edihcnstr',edihcnstr +c write (iout,*) 'edihcnstr',edihcnstr return end c---------------------------------------------------------------------------- @@ -6195,18 +6239,18 @@ c-------------------------------------------------------------------------- logical lprn common /kutas/ lprn CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC -C -C Parallel Antiparallel -C -C o o -C /l\ /j\ -C / \ / \ -C /| o | | o |\ -C \ j|/k\| / \ |/k\|l / -C \ / \ / \ / \ / -C o o o o -C i i -C +C C +C Parallel Antiparallel C +C C +C o o C +C /l\ /j\ C +C / \ / \ C +C /| o | | o |\ C +C \ j|/k\| / \ |/k\|l / C +C \ / \ / \ / \ / C +C o o o o C +C i i C +C C CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC itk=itortyp(itype(k)) s1= scalar2(AEAb1(1,2,imat),CUgb2(1,i)) @@ -6302,18 +6346,18 @@ c---------------------------------------------------------------------------- logical lprn common /kutas/ lprn CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC -C -C Parallel Antiparallel -C -C o o -C \ /l\ /j\ / -C \ / \ / \ / -C o| o | | o |o -C \ j|/k\| \ |/k\|l -C \ / \ \ / \ -C o o -C i i -C +C C +C Parallel Antiparallel C +C C +C o o C +C \ /l\ /j\ / C +C \ / \ / \ / C +C o| o | | o |o C +C \ j|/k\| \ |/k\|l C +C \ / \ \ / \ C +C o o C +C i i C +C C CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC cd write (2,*) 'eello6_graph2: i,',i,' j',j,' k',k,' l',l C AL 7/4/01 s1 would occur in the sixth-order moment, @@ -6486,18 +6530,18 @@ c---------------------------------------------------------------------------- double precision vv(2),pizda(2,2),auxmat(2,2),auxvec(2) logical swap CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC -C -C Parallel Antiparallel -C -C o o -C /l\ / \ /j\ -C / \ / \ / \ -C /| o |o o| o |\ -C j|/k\| / |/k\|l / -C / \ / / \ / -C / o / o -C i i -C +C C +C Parallel Antiparallel C +C C +C o o C +C /l\ / \ /j\ C +C / \ / \ / \ C +C /| o |o o| o |\ C +C j|/k\| / |/k\|l / C +C / \ / / \ / C +C / o / o C +C i i C +C C CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C C 4/7/01 AL Component s1 was removed, because it pertains to the respective @@ -6604,18 +6648,18 @@ c---------------------------------------------------------------------------- & auxvec1(2),auxmat1(2,2) logical swap CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC -C -C Parallel Antiparallel -C -C o o -C /l\ / \ /j\ -C / \ / \ / \ -C /| o |o o| o |\ -C \ j|/k\| \ |/k\|l -C \ / \ \ / \ -C o \ o \ -C i i -C +C C +C Parallel Antiparallel C +C C +C o o C +C /l\ / \ /j\ C +C / \ / \ / \ C +C /| o |o o| o |\ C +C \ j|/k\| \ |/k\|l C +C \ / \ \ / \ C +C o \ o \ C +C i i C +C C CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C C 4/7/01 AL Component s1 was removed, because it pertains to the respective