X-Git-Url: http://mmka.chem.univ.gda.pl/gitweb/?a=blobdiff_plain;f=source%2Fwham%2Fsrc-M%2Fenergy_p_new.F;h=810275a816ae31a92c1335c24e4a4c7573371593;hb=a1bbedcfb50a59f5ae082fe9016af63119b205cc;hp=a4f5fb43ebd3b362cd6d95869b25872fc025c3c6;hpb=0d0c28ca14abcf99ccb778cdba15a9646736571d;p=unres.git diff --git a/source/wham/src-M/energy_p_new.F b/source/wham/src-M/energy_p_new.F index a4f5fb4..810275a 100644 --- a/source/wham/src-M/energy_p_new.F +++ b/source/wham/src-M/energy_p_new.F @@ -44,11 +44,13 @@ C Gay-Berne potential (shifted LJ, angular dependence). goto 106 C Gay-Berne-Vorobjev potential (shifted LJ, angular dependence). 105 call egbv(evdw,evdw_t) +C write(iout,*) 'po elektostatyce' C C Calculate electrostatic (H-bonding) energy of the main chain. C - 106 call eelec(ees,evdw1,eel_loc,eello_turn3,eello_turn4) -C + 106 call eelec(ees,evdw1,eel_loc,eello_turn3,eello_turn4) +C write(iout,*) 'po eelec' + C Calculate excluded-volume interaction energy between peptide groups C and side chains. C @@ -56,8 +58,9 @@ C c c Calculate the bond-stretching energy c + call ebond(estr) -c write (iout,*) "estr",estr +C write (iout,*) "estr",estr C C Calculate the disulfide-bridge and other energy and the contributions C from other distance constraints. @@ -68,12 +71,12 @@ C C Calculate the virtual-bond-angle energy. C call ebend(ebe) -cd print *,'Bend energy finished.' +C print *,'Bend energy finished.' C C Calculate the SC local energy. C call esc(escloc) -cd print *,'SCLOC energy finished.' +C print *,'SCLOC energy finished.' C C Calculate the virtual-bond torsional energy. C @@ -1898,13 +1901,14 @@ cd enddo do i=1,nres num_cont_hb(i)=0 enddo -cd print '(a)','Enter EELEC' -cd write (iout,*) 'iatel_s=',iatel_s,' iatel_e=',iatel_e +C print '(a)','Enter EELEC' +C write (iout,*) 'iatel_s=',iatel_s,' iatel_e=',iatel_e do i=1,nres gel_loc_loc(i)=0.0d0 gcorr_loc(i)=0.0d0 enddo do i=iatel_s,iatel_e + if (i.eq.1) cycle if (itype(i).eq.ntyp1 .or. itype(i+1).eq.ntyp1 & .or. itype(i+2).eq.ntyp1 & .or. itype(i-1).eq.ntyp1 @@ -1926,8 +1930,9 @@ cd write (iout,*) 'iatel_s=',iatel_s,' iatel_e=',iatel_e zmedi=mod(zmedi,boxzsize) if (zmedi.lt.0) zmedi=zmedi+boxzsize num_conti=0 -c write (iout,*) 'i',i,' ielstart',ielstart(i),' ielend',ielend(i) +C write (iout,*) 'i',i,' ielstart',ielstart(i),' ielend',ielend(i) do j=ielstart(i),ielend(i) + if (j.eq.1) cycle if (itype(j).eq.ntyp1 .or. itype(j+1).eq.ntyp1 & .or.itype(j+2).eq.ntyp1 & .or.itype(j-1).eq.ntyp1 @@ -2027,7 +2032,7 @@ c write (iout,'(a6,2i5,0pf7.3,2i5,2e11.3)') c &'evdw1',i,j,evdwij c &,iteli,itelj,aaa,evdw1 -c write (iout,'(a6,2i5,0pf7.3)') 'ees',i,j,eesij + write (iout,'(a6,2i5,0pf7.3)') 'ees',i,j,eesij c write(iout,'(2(2i3,2x),7(1pd12.4)/2(3(1pd12.4),5x)/)') c & iteli,i,itelj,j,aaa,bbb,ael6i,ael3i, c & 1.0D0/dsqrt(rrmij),evdwij,eesij, @@ -3272,7 +3277,8 @@ c write (iout,*) i,vbld(i),vbldp0,diff,AKP*diff*diff gradb(j,i-1)=AKP*diff*dc(j,i-1)/vbld(i) enddo endif - + write (iout,'(a7,i5,4f7.3)') + & "estr bb",i,vbld(i),vbldp0,diff,AKP*diff*diff enddo estr=0.5d0*AKP*estr+estr1 c @@ -3357,6 +3363,7 @@ c write (*,'(a,i2)') 'EBEND ICG=',icg c write (iout,*) ithet_start,ithet_end do i=ithet_start,ithet_end C if (itype(i-1).eq.ntyp1) cycle + if (i.le.2) cycle if ((itype(i-1).eq.ntyp1).or.itype(i-2).eq.ntyp1 & .or.itype(i).eq.ntyp1) cycle C Zero the energy function and its derivative at 0 or pi. @@ -3374,6 +3381,10 @@ C Zero the energy function and its derivative at 0 or pi. ichir21=isign(1,itype(i)) ichir22=isign(1,itype(i)) endif + if (i.eq.3) then + y(1)=0.0D0 + y(2)=0.0D0 + else if (i.gt.3 .and. itype(i-3).ne.ntyp1) then #ifdef OSF @@ -3390,6 +3401,7 @@ C Zero the energy function and its derivative at 0 or pi. y(1)=0.0D0 y(2)=0.0D0 endif + endif if (i.lt.nres .and. itype(i+1).ne.ntyp1) then #ifdef OSF phii1=phi(i+1) @@ -3457,6 +3469,8 @@ C Derivatives of the "mean" values in gamma1 and gamma2. & E_theta,E_tc) endif etheta=etheta+ethetai +c write (iout,'(a6,i5,0pf7.3,f7.3,i5)') +c & 'ebend',i,ethetai,theta(i),itype(i) c write (iout,'(2i3,3f8.3,f10.5)') i,it,rad2deg*theta(i), c & rad2deg*phii,rad2deg*phii1,ethetai if (i.gt.3) gloc(i-3,icg)=gloc(i-3,icg)+wang*E_tc*dthetg1 @@ -3605,7 +3619,9 @@ C etheta=0.0D0 c write (iout,*) "ithetyp",(ithetyp(i),i=1,ntyp1) do i=ithet_start,ithet_end +C if (i.eq.2) cycle C if (itype(i-1).eq.ntyp1) cycle + if (i.le.2) cycle if ((itype(i-1).eq.ntyp1).or.itype(i-2).eq.ntyp1 & .or.itype(i).eq.ntyp1) cycle if (iabs(itype(i+1)).eq.20) iblock=2 @@ -3619,6 +3635,14 @@ C if (itype(i-1).eq.ntyp1) cycle coskt(k)=dcos(k*theti2) sinkt(k)=dsin(k*theti2) enddo + if (i.eq.3) then + phii=0.0d0 + ityp1=nthetyp+1 + do k=1,nsingle + cosph1(k)=0.0d0 + sinph1(k)=0.0d0 + enddo + else if (i.gt.3 .and. itype(i-3).ne.ntyp1) then #ifdef OSF phii=phi(i) @@ -3639,6 +3663,7 @@ C if (itype(i-1).eq.ntyp1) cycle sinph1(k)=0.0d0 enddo endif + endif if (i.lt.nres .and. itype(i+1).ne.ntyp1) then #ifdef OSF phii1=phi(i+1) @@ -3799,14 +3824,14 @@ C ALPHA and OMEGA. common /sccalc/ time11,time12,time112,theti,it,nlobit delta=0.02d0*pi escloc=0.0D0 -c write (iout,'(a)') 'ESC' +C write (iout,*) 'ESC' do i=loc_start,loc_end it=itype(i) if (it.eq.ntyp1) cycle if (it.eq.10) goto 1 nlobit=nlob(iabs(it)) c print *,'i=',i,' it=',it,' nlobit=',nlobit -c write (iout,*) 'i=',i,' ssa=',ssa,' ssad=',ssad +C write (iout,*) 'i=',i,' ssa=',ssa,' ssad=',ssad theti=theta(i+1)-pipol x(1)=dtan(theti) x(2)=alph(i) @@ -3842,8 +3867,8 @@ c write (iout,*) "i",i," x",x(1),x(2),x(3) dersc(k)=ss*dersc(k)+(1.0d0-ss)*dersc0(k) enddo dersc(2)=dersc(2)+ssd*(escloci-esclocbi) -c write (iout,*) 'i=',i,x(2)*rad2deg,escloci0,escloci, -c & esclocbi,ss,ssd + write (iout,*) 'i=',i,x(2)*rad2deg,escloci0,escloci, + & esclocbi,ss,ssd escloci=ss*escloci+(1.0d0-ss)*esclocbi c escloci=esclocbi c write (iout,*) escloci @@ -3877,15 +3902,17 @@ c write (iout,*) escloci enddo dersc(2)=dersc(2)+ssd*(escloci-esclocbi) c write (iout,*) 'i=',i,x(2)*rad2deg,escloci0,escloci, -c & esclocbi,ss,ssd +c & esclocbi,ss,ssd escloci=ss*escloci+(1.0d0-ss)*esclocbi -c write (iout,*) escloci +C write (iout,*) 'i=',i, escloci else call enesc(x,escloci,dersc,ddummy,.false.) endif escloc=escloc+escloci -c write (iout,*) 'i=',i,' escloci=',escloci,' dersc=',dersc +C write (iout,*) 'i=',i,' escloci=',escloci,' dersc=',dersc + write (iout,'(a6,i5,0pf7.3)') + & 'escloc',i,escloci gloc(nphi+i-1,icg)=gloc(nphi+i-1,icg)+ & wscloc*dersc(1) @@ -4575,6 +4602,7 @@ C Set lprn=.true. for debugging c lprn=.true. etors=0.0D0 do i=iphi_start,iphi_end + if (i.le.2) cycle if (itype(i-2).eq.ntyp1.or. itype(i-1).eq.ntyp1 & .or. itype(i).eq.ntyp1 .or. itype(i-3).eq.ntyp1) cycle C if (itype(i-2).eq.ntyp1 .or. itype(i-1).eq.ntyp1 @@ -4678,6 +4706,7 @@ C Set lprn=.true. for debugging c lprn=.true. etors_d=0.0D0 do i=iphi_start,iphi_end-1 + if (i.le.3) cycle C if (itype(i-2).eq.ntyp1.or. itype(i-1).eq.ntyp1 C & .or. itype(i).eq.ntyp1 .or. itype(i+1).eq.ntyp1) cycle if ((itype(i-2).eq.ntyp1).or.itype(i-3).eq.ntyp1.or.