From c6588e178e83ce42d404f65fecf3dc0aa71764c5 Mon Sep 17 00:00:00 2001 From: Adam Sieradzan Date: Tue, 29 Sep 2015 19:48:56 +0200 Subject: [PATCH] ESSENTIAL CHANGE - BUG FIX in ENERGY to have old Correlation --- source/unres/src_MD-M/energy_p_new_barrier.F | 14 +++++++++++++- source/wham/src-M/energy_p_new.F | 20 +++++++++++++++++--- 2 files changed, 30 insertions(+), 4 deletions(-) diff --git a/source/unres/src_MD-M/energy_p_new_barrier.F b/source/unres/src_MD-M/energy_p_new_barrier.F index 3192819..3e1eb6e 100644 --- a/source/unres/src_MD-M/energy_p_new_barrier.F +++ b/source/unres/src_MD-M/energy_p_new_barrier.F @@ -2719,6 +2719,17 @@ c write(iout,*) 'b1=',b1(1,i-2) c write (iout,*) 'theta=', theta(i-1) enddo #else + if (i.gt. nnt+2 .and. i.lt.nct+2) then + iti = itortyp(itype(i-2)) + else + iti=ntortyp+1 + endif +c if (i.gt. iatel_s+1 .and. i.lt.iatel_e+4) then + if (i.gt. nnt+1 .and. i.lt.nct+1) then + iti1 = itortyp(itype(i-1)) + else + iti1=ntortyp+1 + endif b1(1,i-2)=b(3,iti) b1(2,i-2)=b(5,iti) b2(1,i-2)=b(2,iti) @@ -2873,6 +2884,7 @@ c if (i.gt. iatel_s+1 .and. i.lt.iatel_e+4) then do k=1,2 mu(k,i-2)=Ub2(k,i-2)+b1(k,i-1) enddo +C write (iout,*) 'mumu',i,b1(1,i-1),Ub2(1,i-2) c write (iout,*) 'mu ',mu(:,i-2),i-2 cd write (iout,*) 'mu1',mu1(:,i-2) cd write (iout,*) 'mu2',mu2(:,i-2) @@ -3987,7 +3999,7 @@ C Contribution to the local-electrostatic energy coming from the i-j pair & +a33*muij(4) c write (iout,*) 'i',i,' j',j,itype(i),itype(j), c & ' eel_loc_ij',eel_loc_ij -c write(iout,*) 'muije=',muij(1),muij(2),muij(3),muij(4) +C write(iout,*) 'muije=',i,j,muij(1),muij(2),muij(3),muij(4) C Calculate patrial derivative for theta angle #ifdef NEWCORR geel_loc_ij=a22*gmuij1(1) diff --git a/source/wham/src-M/energy_p_new.F b/source/wham/src-M/energy_p_new.F index c228c20..b8e07c5 100644 --- a/source/wham/src-M/energy_p_new.F +++ b/source/wham/src-M/energy_p_new.F @@ -913,7 +913,7 @@ C lipbufthick is thickenes of lipid buffore & +aa_aq(itypi,itypj)*(2.0d0-sslipi-sslipj)/2.0d0 bb=bb_lip(itypi,itypj)*(sslipi+sslipj)/2.0d0 & +bb_aq(itypi,itypj)*(2.0d0-sslipi-sslipj)/2.0d0 - write(iout,*),aa,aa_lip(itypi,itypj),aa_aq(itypi,itypj) +C write(iout,*),aa,aa_lip(itypi,itypj),aa_aq(itypi,itypj) C checking the distance dist_init=(xj-xi)**2+(yj-yi)**2+(zj-zi)**2 xj_safe=xj @@ -1843,6 +1843,8 @@ c print *,"itilde3 i iti iti1",i,iti,iti1 do k=1,2 mu(k,i-2)=Ub2(k,i-2)+b1(k,iti1) enddo +C write (iout,*) 'mumu',i,b1(1,iti),Ub2(1,i-2) + C Vectors and matrices dependent on a single virtual-bond dihedral. call matvec2(DD(1,1,iti),b1tilde(1,iti1),auxvec(1)) call matvec2(Ug2(1,1,i-2),auxvec(1),Ug2Db1t(1,i-2)) @@ -2464,8 +2466,9 @@ C Contribution to the local-electrostatic energy coming from the i-j pair eel_loc_ij=a22*muij(1)+a23*muij(2)+a32*muij(3) & +a33*muij(4) c write (iout,*) 'i',i,' j',j,' eel_loc_ij',eel_loc_ij -c write (iout,'(a6,2i5,0pf7.3)') -c & 'eelloc',i,j,eel_loc_ij +C write (iout,'(a6,2i5,0pf7.3)') +C & 'eelloc',i,j,eel_loc_ij +C write(iout,*) 'muije=',i,j,muij(1),muij(2),muij(3),muij(4) c write (iout,*) a22,muij(1),a23,muij(2),a32,muij(3) eel_loc=eel_loc+eel_loc_ij C Partial derivatives in virtual-bond dihedral angles gamma @@ -2803,6 +2806,16 @@ C Cartesian derivatives enddo endif else if (j.eq.i+3 .and. itype(i+2).ne.ntyp1) then + if (itype(i).eq.ntyp1 .or. itype(i+1).eq.ntyp1 +C changes suggested by Ana to avoid out of bounds + & .or.((i+5).gt.nres) + & .or.((i-1).le.0) +C end of changes suggested by Ana + & .or. itype(i+3).eq.ntyp1 + & .or. itype(i+4).eq.ntyp1 + & .or. itype(i+5).eq.ntyp1 + & .or. itype(i).eq.ntyp1 + & .or. itype(i-1).eq.ntyp1) goto 178 CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C C Fourth-order contributions @@ -2942,6 +2955,7 @@ C Remaining derivatives of this turn contribution gcorr4_turn(l,j1)=gcorr4_turn(l,j1)-(s1+s2+s3) enddo endif + 178 continue endif return end -- 1.7.9.5