X-Git-Url: http://mmka.chem.univ.gda.pl/gitweb/?a=blobdiff_plain;f=source%2Fwham%2Fsrc-M-NEWCORR%2Fenergy_p_new.F;h=917b6a6738724ee7d4fb3adbb7b3bf93ba515a44;hb=5467f8060e885151d3415a9897913776c337d88a;hp=0dddd7faca0b3ce3bee64dcc75ecde27aaee2d77;hpb=4d3fd4762ad7bfc6b3fdd9915befe7ea8da7f2e0;p=unres.git diff --git a/source/wham/src-M-NEWCORR/energy_p_new.F b/source/wham/src-M-NEWCORR/energy_p_new.F index 0dddd7f..917b6a6 100644 --- a/source/wham/src-M-NEWCORR/energy_p_new.F +++ b/source/wham/src-M-NEWCORR/energy_p_new.F @@ -1718,8 +1718,8 @@ c write (iout,*) 'theta=', theta(i-1) iti1=ntortyp+1 endif cd write (iout,*) '*******i',i,' iti1',iti -cd write (iout,*) 'b1',b1(:,iti) -cd write (iout,*) 'b2',b2(:,iti) +cd write (iout,*) 'b1',b1(:,i-2) +cd write (iout,*) 'b2',b2(:,i-2) cd write (iout,*) 'Ug',Ug(:,:,i-2) c print *,"itilde1 i iti iti1",i,iti,iti1 if (i .gt. iatel_s+2) then @@ -4362,7 +4362,6 @@ C----------------------------------------------------------------------------- include 'COMMON.IOUNITS' include 'COMMON.FFIELD' include 'COMMON.TORCNSTR' - logical lprn C Set lprn=.true. for debugging lprn=.false. c lprn=.true. @@ -4447,10 +4446,11 @@ c------------------------------------------------------------------------------ include 'COMMON.IOUNITS' include 'COMMON.FFIELD' include 'COMMON.TORCNSTR' - logical lprn + logical lprn,energy_dec C Set lprn=.true. for debugging lprn=.false. c lprn=.true. + energy_dec=.false. etors=0.0D0 do i=iphi_start,iphi_end if (itype(i-2).eq.ntyp1 .or. itype(i-1).eq.ntyp1 @@ -4465,6 +4465,7 @@ c lprn=.true. itori1=itortyp(itype(i-1)) phii=phi(i) gloci=0.0D0 + etors_ii=0.0d0 C Regular cosine and sine terms do j=1,nterm(itori,itori1,iblock) v1ij=v1(j,itori,itori1,iblock) @@ -4472,6 +4473,8 @@ C Regular cosine and sine terms cosphi=dcos(j*phii) sinphi=dsin(j*phii) etors=etors+v1ij*cosphi+v2ij*sinphi + if (energy_dec) etors_ii=etors_ii+ + & v1ij*cosphi+v2ij*sinphi gloci=gloci+j*(v2ij*cosphi-v1ij*sinphi) enddo C Lorentz terms @@ -4488,16 +4491,19 @@ C pom=vl2ij*cosphi+vl3ij*sinphi pom1=1.0d0/(pom*pom+1.0d0) etors=etors+vl1ij*pom1 -c if (energy_dec) etors_ii=etors_ii+ -c & vl1ij*pom1 + if (energy_dec) etors_ii=etors_ii+ + & vl1ij*pom1 pom=-pom*pom1*pom1 gloci=gloci+vl1ij*(vl3ij*cosphi-vl2ij*sinphi)*pom enddo C Subtract the constant term etors=etors-v0(itori,itori1,iblock) + if (energy_dec) write (iout,'(a6,i5,0pf7.3)') + & 'etor',i,etors_ii if (lprn) - & write (iout,'(2(a3,2x,i3,2x),2i3,6f8.3/26x,6f8.3/)') + & write (iout,'(2(a3,2x,i3,2x),2i3,f10.2,6f8.3/36x,6f8.3/)') & restyp(itype(i-2)),i-2,restyp(itype(i-1)),i-1,itori,itori1, + & rad2deg*phii, & (v1(j,itori,itori1,1),j=1,6),(v2(j,itori,itori1,1),j=1,6) gloc(i-3,icg)=gloc(i-3,icg)+wtor*fact*gloci c write (iout,*) 'i=',i,' gloc=',gloc(i-3,icg) @@ -4550,8 +4556,8 @@ C 6/23/01 Compute double torsional energy include 'COMMON.TORCNSTR' logical lprn C Set lprn=.true. for debugging - lprn=.false. -c lprn=.true. +c lprn=.false. + lprn=.true. etors_d=0.0D0 do i=iphi_start,iphi_end-1 if (itype(i-2).eq.ntyp1.or. itype(i-1).eq.ntyp1