X-Git-Url: http://mmka.chem.univ.gda.pl/gitweb/?a=blobdiff_plain;f=source%2Funres%2Fsrc-HCD-5D%2Fenergy_p_new_barrier.F.safe;h=ae8e4490ae447bc1d188f0a2ac4afdf8059dfeb1;hb=5836ecdab5a8b95f079bbf6e07374dee3fce8a26;hp=647987867e99e3d5d2bd20d39d4f0fb101c87788;hpb=a05e8996b6fb955dc21aef7a85db21ee73d0a9d3;p=unres.git diff --git a/source/unres/src-HCD-5D/energy_p_new_barrier.F.safe b/source/unres/src-HCD-5D/energy_p_new_barrier.F.safe index 6479878..ae8e449 100644 --- a/source/unres/src-HCD-5D/energy_p_new_barrier.F.safe +++ b/source/unres/src-HCD-5D/energy_p_new_barrier.F.safe @@ -324,6 +324,7 @@ C else esccor=0.0d0 endif +#ifdef FOURBODY C print *,"PRZED MULIt" c print *,"Processor",myrank," computed Usccorr" C @@ -352,6 +353,7 @@ c write (iout,*) "MULTIBODY_HB ecorr",ecorr,ecorr5,ecorr6,n_corr, c & n_corr1 c call flush(iout) endif +#endif c print *,"Processor",myrank," computed Ucorr" c write (iout,*) "nsaxs",nsaxs," saxs_mode",saxs_mode if (nsaxs.gt.0 .and. saxs_mode.eq.0) then @@ -1314,9 +1316,16 @@ C Bartek write (iout,10) evdw,wsc,evdw2,wscp,ees,welec,evdw1,wvdwpp, & estr,wbond,ebe,wang, & escloc,wscloc,etors,wtor,etors_d,wtor_d,ehpb,wstrain, +#ifdef FOURBODY & ecorr,wcorr, - & ecorr5,wcorr5,ecorr6,wcorr6,eel_loc,wel_loc,eello_turn3,wturn3, - & eello_turn4,wturn4,eello_turn6,wturn6,esccor,wsccor,edihcnstr, + & ecorr5,wcorr5,ecorr6,wcorr6, +#endif + & eel_loc,wel_loc,eello_turn3,wturn3, + & eello_turn4,wturn4, +#ifdef FOURBODY + & eello_turn6,wturn6, +#endif + & esccor,wsccor,edihcnstr, & ethetacnstr,ebr*nss,Uconst,wumb,eliptran,wliptran,Eafmforce, & etube,wtube,esaxs,wsaxs,ehomology_constr, & edfadis,wdfa_dist,edfator,wdfa_tor,edfanei,wdfa_nei, @@ -1334,13 +1343,17 @@ C Bartek & 'ETORSD=',1pE16.6,' WEIGHT=',1pE16.6,' (double torsional)'/ & 'EHBP= ',1pE16.6,' WEIGHT=',1pE16.6, & ' (SS bridges & dist. cnstr.)'/ +#ifdef FOURBODY & 'ECORR4=',1pE16.6,' WEIGHT=',1pE16.6,' (multi-body)'/ & 'ECORR5=',1pE16.6,' WEIGHT=',1pE16.6,' (multi-body)'/ & 'ECORR6=',1pE16.6,' WEIGHT=',1pE16.6,' (multi-body)'/ +#endif & 'EELLO= ',1pE16.6,' WEIGHT=',1pE16.6,' (electrostatic-local)'/ & 'ETURN3=',1pE16.6,' WEIGHT=',1pE16.6,' (turns, 3rd order)'/ & 'ETURN4=',1pE16.6,' WEIGHT=',1pE16.6,' (turns, 4th order)'/ +#ifdef FOURBODY & 'ETURN6=',1pE16.6,' WEIGHT=',1pE16.6,' (turns, 6th order)'/ +#endif & 'ESCCOR=',1pE16.6,' WEIGHT=',1pE16.6,' (backbone-rotamer corr)'/ & 'EDIHC= ',1pE16.6,' (virtual-bond dihedral angle restraints)'/ & 'ETHETC=',1pE16.6,' (virtual-bond angle restraints)'/ @@ -1361,9 +1374,16 @@ C Bartek write (iout,10) evdw,wsc,evdw2,wscp,ees,welec, & estr,wbond,ebe,wang, & escloc,wscloc,etors,wtor,etors_d,wtor_d,ehpb,wstrain, +#ifdef FOURBODY & ecorr,wcorr, - & ecorr5,wcorr5,ecorr6,wcorr6,eel_loc,wel_loc,eello_turn3,wturn3, - & eello_turn4,wturn4,eello_turn6,wturn6,esccor,wsccor,edihcnstr, + & ecorr5,wcorr5,ecorr6,wcorr6, +#endif + & eel_loc,wel_loc,eello_turn3,wturn3, + & eello_turn4,wturn4, +#ifdef FOURBODY + & eello_turn6,wturn6, +#endif + & esccor,wsccor,edihcnstr, & ethetacnstr,ebr*nss,Uconst,wumb,eliptran,wliptran,Eafmforc, & etube,wtube,esaxs,wsaxs,ehomology_constr, & edfadis,wdfa_dist,edfator,wdfa_tor,edfanei,wdfa_nei, @@ -1380,13 +1400,17 @@ C Bartek & 'ETORSD=',1pE16.6,' WEIGHT=',1pE16.6,' (double torsional)'/ & 'EHBP= ',1pE16.6,' WEIGHT=',1pE16.6, & ' (SS bridges & dist. restr.)'/ +#ifdef FOURBODY & 'ECORR4=',1pE16.6,' WEIGHT=',1pE16.6,' (multi-body)'/ & 'ECORR5=',1pE16.6,' WEIGHT=',1pE16.6,' (multi-body)'/ & 'ECORR6=',1pE16.6,' WEIGHT=',1pE16.6,' (multi-body)'/ +#endif & 'EELLO= ',1pE16.6,' WEIGHT=',1pE16.6,' (electrostatic-local)'/ & 'ETURN3=',1pE16.6,' WEIGHT=',1pE16.6,' (turns, 3rd order)'/ & 'ETURN4=',1pE16.6,' WEIGHT=',1pE16.6,' (turns, 4th order)'/ +#ifdef FOURBODY & 'ETURN6=',1pE16.6,' WEIGHT=',1pE16.6,' (turns, 6th order)'/ +#endif & 'ESCCOR=',1pE16.6,' WEIGHT=',1pE16.6,' (backbone-rotamer corr)'/ & 'EDIHC= ',1pE16.6,' (virtual-bond dihedral angle restraints)'/ & 'ETHETC=',1pE16.6,' (virtual-bond angle restraints)'/ @@ -1425,7 +1449,10 @@ C include 'COMMON.SBRIDGE' include 'COMMON.NAMES' include 'COMMON.IOUNITS' +#ifdef FOURBODY include 'COMMON.CONTACTS' + include 'COMMON.CONTMAT' +#endif double precision gg(3) double precision evdw,evdwij integer i,j,k,itypi,itypj,itypi1,num_conti,iint @@ -1491,6 +1518,7 @@ cgrad gvdwc(l,k)=gvdwc(l,k)+gg(l) cgrad enddo cgrad enddo C +#ifdef FOURBODY C 12/1/95, revised on 5/20/97 C C Calculate the contact function. The ith column of the array JCONT will @@ -1546,10 +1574,13 @@ cd write (iout,'(2i3,3f10.5)') cd & i,j,(gacont(kk,num_conti,i),kk=1,3) endif endif +#endif enddo ! j enddo ! iint C Change 12/1/95 +#ifdef FOURBODY num_cont(i)=num_conti +#endif enddo ! i do i=1,nct do j=1,3 @@ -2076,8 +2107,8 @@ c write (iout,*) "j",j," dc_norm", c & dc_norm(1,nres+j),dc_norm(2,nres+j),dc_norm(3,nres+j) rrij=1.0D0/(xj*xj+yj*yj+zj*zj) rij=dsqrt(rrij) - sss=sscale((1.0d0/rij)/sigma(itypi,itypj)) - sssgrad=sscagrad((1.0d0/rij)/sigma(itypi,itypj)) + sss=sscale(1.0d0/rij,r_cut_int) + sssgrad=sscagrad(1.0d0/rij,r_cut_int) c write (iout,'(a7,4f8.3)') c & "ssscale",sss,((1.0d0/rij)/sigma(itypi,itypj)),r_cut,rlamb @@ -2138,7 +2169,7 @@ C Calculate gradient components. fac=rij*fac c print '(2i4,6f8.4)',i,j,sss,sssgrad* c & evdwij,fac,sigma(itypi,itypj),expon - fac=fac+evdwij/sss*sssgrad/sigma(itypi,itypj)*rij + fac=fac+evdwij/sss*sssgrad*rij c fac=0.0d0 C Calculate the radial part of the gradient gg_lipi(3)=eps1*(eps2rt*eps2rt) @@ -2548,7 +2579,7 @@ C include 'COMMON.SBRIDGE' include 'COMMON.NAMES' include 'COMMON.IOUNITS' - include 'COMMON.CONTACTS' +c include 'COMMON.CONTACTS' dimension gg(3) cd print *,'Entering Esoft_sphere nnt=',nnt,' nct=',nct evdw=0.0D0 @@ -2622,7 +2653,7 @@ C include 'COMMON.CHAIN' include 'COMMON.DERIV' include 'COMMON.INTERACT' - include 'COMMON.CONTACTS' +c include 'COMMON.CONTACTS' include 'COMMON.TORSION' include 'COMMON.VECTORS' include 'COMMON.FFIELD' @@ -2703,8 +2734,8 @@ c write (iout,*) 'i',i,' ielstart',ielstart(i),' ielend',ielend(i) zj=zj_safe-zmedi endif rij=xj*xj+yj*yj+zj*zj - sss=sscale(sqrt(rij)) - sssgrad=sscagrad(sqrt(rij)) + sss=sscale(sqrt(rij),r_cut_int) + sssgrad=sscagrad(sqrt(rij),r_cut_int) if (rij.lt.r0ijsq) then evdw1ij=0.25d0*(rij-r0ijsq)**2 fac=rij-r0ijsq @@ -2946,7 +2977,7 @@ C-------------------------------------------------------------------------- include 'COMMON.CHAIN' include 'COMMON.DERIV' include 'COMMON.INTERACT' - include 'COMMON.CONTACTS' + include 'COMMON.CORRMAT' include 'COMMON.TORSION' include 'COMMON.VECTORS' include 'COMMON.FFIELD' @@ -2962,18 +2993,26 @@ c write(iout,*) "itype2loc",itype2loc #else do i=3,nres+1 #endif - if (i.gt. nnt+2 .and. i.lt.nct+2) then + ii=ireschain(i-2) +c write (iout,*) "i",i,i-2," ii",ii + if (ii.eq.0) cycle + innt=chain_border(1,ii) + inct=chain_border(2,ii) +c write (iout,*) "i",i,i-2," ii",ii," innt",innt," inct",inct +c if (i.gt. nnt+2 .and. i.lt.nct+2) then + if (i.gt. innt+2 .and. i.lt.inct+2) then iti = itype2loc(itype(i-2)) else iti=nloctyp 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 + if (i.gt. innt+1 .and. i.lt.inct+1) then iti1 = itype2loc(itype(i-1)) else iti1=nloctyp endif - write(iout,*),"i",i,i-2," iti",iti," iti1",iti1 +c write(iout,*),"i",i,i-2," iti",itype(i-2),iti, +c & " iti1",itype(i-1),iti1 #ifdef NEWCORR cost1=dcos(theta(i-1)) sint1=dsin(theta(i-1)) @@ -3039,7 +3078,8 @@ c b2tilde(2,i-2)=-b2(2,i-2) write (iout,*) 'theta=', theta(i-1) #endif #else - if (i.gt. nnt+2 .and. i.lt.nct+2) then + if (i.gt. innt+2 .and. i.lt.inct+2) then +c if (i.gt. nnt+2 .and. i.lt.nct+2) then iti = itype2loc(itype(i-2)) else iti=nloctyp @@ -3084,12 +3124,14 @@ c write(iout,*) 'b2=',(b2(k,i-2),k=1,2) #endif enddo + mu=0.0d0 #ifdef PARMAT do i=ivec_start+2,ivec_end+2 #else do i=3,nres+1 #endif - if (i .lt. nres+1) then +c if (itype(i-1).eq.ntyp1 .and. itype(i).eq.ntyp1) cycle + if (i .lt. nres+1 .and. itype(i-1).lt.ntyp1) then sin1=dsin(phi(i)) cos1=dcos(phi(i)) sintab(i-2)=sin1 @@ -3126,7 +3168,7 @@ c Ug2(2,1,i-2)=0.0d0 Ug2(2,2,i-2)=0.0d0 endif - if (i .gt. 3 .and. i .lt. nres+1) then + if (i .gt. 3) then obrot_der(1,i-2)=-sin1 obrot_der(2,i-2)= cos1 Ugder(1,1,i-2)= sin1 @@ -3156,7 +3198,8 @@ c Ug2der(2,2,i-2)=0.0d0 endif c if (i.gt. iatel_s+2 .and. i.lt.iatel_e+5) then - if (i.gt. nnt+2 .and. i.lt.nct+2) then +c if (i.gt. nnt+2 .and. i.lt.nct+2) then + if (i.gt.nnt+2 .and.i.lt.nct+2) then iti = itype2loc(itype(i-2)) else iti=nloctyp @@ -3185,6 +3228,7 @@ c & EE(1,2,iti),EE(2,2,i) c write(iout,*) "Macierz EUG", c & eug(1,1,i-2),eug(1,2,i-2),eug(2,1,i-2), c & eug(2,2,i-2) +#ifdef FOURBODY if (wcorr4.gt.0.0d0 .or. wcorr5.gt.0.0d0 .or. wcorr6.gt.0.0d0) & then call matmat2(CC(1,1,i-2),Ug(1,1,i-2),CUg(1,1,i-2)) @@ -3193,6 +3237,7 @@ c & eug(2,2,i-2) call matvec2(Ctilde(1,1,i-1),obrot(1,i-2),Ctobr(1,i-2)) call matvec2(Dtilde(1,1,i-2),obrot2(1,i-2),Dtobr2(1,i-2)) endif +#endif else do k=1,2 Ub2(k,i-2)=0.0d0 @@ -3226,7 +3271,6 @@ c if (i.gt. iatel_s+1 .and. i.lt.iatel_e+4) then c mu(k,i-2)=b1(k,i-1) c mu(k,i-2)=Ub2(k,i-2) enddo -#define MUOUT #ifdef MUOUT write (iout,'(2hmu,i3,3f8.1,12f10.5)') i-2,rad2deg*theta(i-1), & rad2deg*theta(i),rad2deg*phi(i),mu(1,i-2),mu(2,i-2), @@ -3235,10 +3279,10 @@ c mu(k,i-2)=Ub2(k,i-2) & +dsqrt(b1(1,i-1)**2+b1(2,i-1)**2), & ((ee(l,k,i-2),l=1,2),k=1,2) #endif -#undef MUOUT cd write (iout,*) 'mu1',mu1(:,i-2) cd write (iout,*) 'mu2',mu2(:,i-2) cd write (iout,*) 'mu',i-2,mu(:,i-2) +#ifdef FOURBODY if (wcorr4.gt.0.0d0 .or. wcorr5.gt.0.0d0 .or.wcorr6.gt.0.0d0) & then call matmat2(CC(1,1,i-1),Ugder(1,1,i-2),CUgder(1,1,i-2)) @@ -3257,7 +3301,9 @@ C Vectors and matrices dependent on a single virtual-bond dihedral. call matmat2(EUg(1,1,i-2),DD(1,1,i-1),EUgD(1,1,i-2)) call matmat2(EUgder(1,1,i-2),DD(1,1,i-1),EUgDder(1,1,i-2)) endif +#endif enddo +#ifdef FOURBODY C Matrices dependent on two consecutive virtual-bond dihedrals. C The order of matrices is from left to right. if (wcorr4.gt.0.0d0 .or. wcorr5.gt.0.0d0 .or.wcorr6.gt.0.0d0) @@ -3274,6 +3320,7 @@ c do i=max0(ivec_start,2),ivec_end call matmat2(auxmat(1,1),EUg(1,1,i),Ug2DtEUgder(1,1,1,i)) enddo endif +#endif #if defined(MPI) && defined(PARMAT) #ifdef DEBUG c if (fg_rank.eq.0) then @@ -3342,6 +3389,7 @@ c & " ivec_count",(ivec_count(i),i=0,nfgtasks-1) call MPI_Allgatherv(sintab2(ivec_start),ivec_count(fg_rank1), & MPI_DOUBLE_PRECISION,sintab2(1),ivec_count(0),ivec_displ(0), & MPI_DOUBLE_PRECISION,FG_COMM1,IERR) +#ifdef FOURBODY if (wcorr4.gt.0.0d0 .or. wcorr5.gt.0.0d0 .or. wcorr6.gt.0.0d0) & then call MPI_Allgatherv(Ctobr(1,ivec_start),ivec_count(fg_rank1), @@ -3417,6 +3465,7 @@ c & " ivec_count",(ivec_count(i),i=0,nfgtasks-1) & MPI_MAT2,Ug2DtEUgder(1,1,1,1),ivec_count(0),ivec_displ(0), & MPI_MAT2,FG_COMM1,IERR) endif +#endif #else c Passes matrix info through the ring isend=fg_rank1 @@ -3461,6 +3510,7 @@ c call flush(iout) & iprev,6600+irecv,FG_COMM,status,IERR) c write (iout,*) "Gather PRECOMP12" c call flush(iout) +#ifdef FOURBODY if (wcorr4.gt.0.0d0 .or. wcorr5.gt.0.0d0 .or. wcorr6.gt.0.0d0) & then call MPI_SENDRECV(ug2db1t(1,ivec_displ(isend)+1),1, @@ -3480,6 +3530,7 @@ c call flush(iout) & Ug2DtEUgder(1,1,1,ivec_displ(irecv)+1),1, & MPI_PRECOMP23(lenrecv), & iprev,9900+irecv,FG_COMM,status,IERR) +#endif c write (iout,*) "Gather PRECOMP23" c call flush(iout) endif @@ -3555,7 +3606,11 @@ C include 'COMMON.CHAIN' include 'COMMON.DERIV' include 'COMMON.INTERACT' +#ifdef FOURBODY include 'COMMON.CONTACTS' + include 'COMMON.CONTMAT' +#endif + include 'COMMON.CORRMAT' include 'COMMON.TORSION' include 'COMMON.VECTORS' include 'COMMON.FFIELD' @@ -3628,9 +3683,11 @@ cd enddo eello_turn3=0.0d0 eello_turn4=0.0d0 ind=0 +#ifdef FOURBODY do i=1,nres num_cont_hb(i)=0 enddo +#endif cd print '(a)','Enter EELEC' cd write (iout,*) 'iatel_s=',iatel_s,' iatel_e=',iatel_e do i=1,nres @@ -3680,7 +3737,9 @@ c end if num_conti=0 call eelecij(i,i+2,ees,evdw1,eel_loc) if (wturn3.gt.0.0d0) call eturn3(i,eello_turn3) +#ifdef FOURBODY num_cont_hb(i)=num_conti +#endif enddo do i=iturn4_start,iturn4_end if (i.lt.1) cycle @@ -3736,12 +3795,16 @@ c endif zmedi=mod(zmedi,boxzsize) if (zmedi.lt.0) zmedi=zmedi+boxzsize +#ifdef FOURBODY num_conti=num_cont_hb(i) +#endif c write(iout,*) "JESTEM W PETLI" call eelecij(i,i+3,ees,evdw1,eel_loc) if (wturn4.gt.0.0d0 .and. itype(i+2).ne.ntyp1) & call eturn4(i,eello_turn4) +#ifdef FOURBODY num_cont_hb(i)=num_conti +#endif enddo ! i C Loop over all neighbouring boxes C do xshift=-1,1 @@ -3808,7 +3871,9 @@ c go to 166 c endif c write (iout,*) 'i',i,' ielstart',ielstart(i),' ielend',ielend(i) +#ifdef FOURBODY num_conti=num_cont_hb(i) +#endif C I TU KURWA do j=ielstart(i),ielend(i) C do j=16,17 @@ -3824,7 +3889,9 @@ c & .or.itype(j-1).eq.ntyp1 &) cycle call eelecij(i,j,ees,evdw1,eel_loc) enddo ! j +#ifdef FOURBODY num_cont_hb(i)=num_conti +#endif enddo ! i C enddo ! zshift C enddo ! yshift @@ -3855,7 +3922,11 @@ C------------------------------------------------------------------------------- include 'COMMON.CHAIN' include 'COMMON.DERIV' include 'COMMON.INTERACT' +#ifdef FOURBODY include 'COMMON.CONTACTS' + include 'COMMON.CONTMAT' +#endif + include 'COMMON.CORRMAT' include 'COMMON.TORSION' include 'COMMON.VECTORS' include 'COMMON.FFIELD' @@ -3973,8 +4044,8 @@ C yj=yj-ymedi C zj=zj-zmedi rij=xj*xj+yj*yj+zj*zj - sss=sscale(sqrt(rij)) - sssgrad=sscagrad(sqrt(rij)) + sss=sscale(sqrt(rij),r_cut_int) + sssgrad=sscagrad(sqrt(rij),r_cut_int) c if (sss.gt.0.0d0) then rrmij=1.0D0/rij rij=dsqrt(rij) @@ -4645,6 +4716,7 @@ C Remaining derivatives of eello ENDIF C Change 12/26/95 to calculate four-body contributions to H-bonding energy c if (j.gt.i+1 .and. num_conti.le.maxconts) then +#ifdef FOURBODY if (wcorr+wcorr4+wcorr5+wcorr6.gt.0.0d0 & .and. num_conti.le.maxconts) then c write (iout,*) i,j," entered corr" @@ -4835,6 +4907,7 @@ cdiag enddo endif ! num_conti.le.maxconts endif ! fcont.gt.0 endif ! j.gt.i+1 +#endif if (wturn3.gt.0.0d0 .or. wturn4.gt.0.0d0) then do k=1,4 do l=1,3 @@ -4867,7 +4940,7 @@ C Third- and fourth-order contributions from turns include 'COMMON.CHAIN' include 'COMMON.DERIV' include 'COMMON.INTERACT' - include 'COMMON.CONTACTS' + include 'COMMON.CORRMAT' include 'COMMON.TORSION' include 'COMMON.VECTORS' include 'COMMON.FFIELD' @@ -5050,7 +5123,7 @@ C Third- and fourth-order contributions from turns include 'COMMON.CHAIN' include 'COMMON.DERIV' include 'COMMON.INTERACT' - include 'COMMON.CONTACTS' + include 'COMMON.CORRMAT' include 'COMMON.TORSION' include 'COMMON.VECTORS' include 'COMMON.FFIELD' @@ -5621,7 +5694,7 @@ cd write (iout,*) 'iatscp_s=',iatscp_s,' iatscp_e=',iatscp_e C do xshift=-1,1 C do yshift=-1,1 C do zshift=-1,1 - if (energy_dec) write (iout,*) "escp:",r_cut,rlamb + if (energy_dec) write (iout,*) "escp:",r_cut_int,rlamb do i=iatscp_s,iatscp_e if (itype(i).eq.ntyp1 .or. itype(i+1).eq.ntyp1) cycle iteli=itel(i) @@ -5743,11 +5816,11 @@ CHERE IS THE CALCULATION WHICH MIRROR IMAGE IS THE CLOSEST ONE c print *,xj,yj,zj,'polozenie j' rrij=1.0D0/(xj*xj+yj*yj+zj*zj) c print *,rrij - sss=sscale(1.0d0/(dsqrt(rrij))) + sss=sscale(1.0d0/(dsqrt(rrij)),r_cut_int) c print *,r_cut,1.0d0/dsqrt(rrij),sss,'tu patrz' c if (sss.eq.0) print *,'czasem jest OK' if (sss.le.0.0d0) cycle - sssgrad=sscagrad(1.0d0/(dsqrt(rrij))) + sssgrad=sscagrad(1.0d0/(dsqrt(rrij)),r_cut_int) fac=rrij**expon2 e1=fac*fac*aad(itypj,iteli) e2=fac*bad(itypj,iteli) @@ -6174,6 +6247,12 @@ c estr=0.0d0 estr1=0.0d0 do i=ibondp_start,ibondp_end +c 3/4/2020 Adam: removed dummy bond graient if Calpha and SC coords are +c used +#ifdef FIVEDIAG + if (itype(i-1).eq.ntyp1 .or. itype(i).eq.ntyp1) cycle + diff = vbld(i)-vbldp0 +#else if (itype(i-1).eq.ntyp1 .and. itype(i).eq.ntyp1) cycle c estr1=estr1+gnmr1(vbld(i),-1.0d0,distchainmax) c do j=1,3 @@ -6184,15 +6263,16 @@ c if (energy_dec) write(iout,*) c & "estr1",i,gnmr1(vbld(i),-1.0d0,distchainmax) c else C Checking if it involves dummy (NH3+ or COO-) group - if (itype(i-1).eq.ntyp1 .or. itype(i).eq.ntyp1) then + if (itype(i-1).eq.ntyp1 .or. itype(i).eq.ntyp1) then C YES vbldpDUM is the equlibrium length of spring for Dummy atom - diff = vbld(i)-vbldpDUM - if (energy_dec) write(iout,*) "dum_bond",i,diff - else -C NO vbldp0 is the equlibrium lenght of spring for peptide group - diff = vbld(i)-vbldp0 - endif - if (energy_dec) write (iout,'(a7,i5,4f7.3)') + diff = vbld(i)-vbldpDUM + if (energy_dec) write(iout,*) "dum_bond",i,diff + else +C NO vbldp0 is the equlibrium length of spring for peptide group + diff = vbld(i)-vbldp0 + endif +#endif + if (energy_dec) write (iout,'(a7,i5,4f7.3)') & "estr bb",i,vbld(i),vbldp0,diff,AKP*diff*diff estr=estr+diff*diff do j=1,3 @@ -7144,7 +7224,8 @@ c & sumene4, c & dscp1,dscp2,sumene c sumene = enesc(x,xx,yy,zz,cost2tab(i+1),sint2tab(i+1)) escloc = escloc + sumene -c write (2,*) "i",i," escloc",sumene,escloc,it,itype(i) + if (energy_dec) write (2,*) "i",i," itype",itype(i)," it",it, + & " escloc",sumene,escloc,it,itype(i) c & ,zz,xx,yy c#define DEBUG #ifdef DEBUG @@ -8709,6 +8790,7 @@ c write (iout,*) "EBACK_SC_COR",i,v1ij*cosphi+v2ij*sinphi,intertyp return end +#ifdef FOURBODY c---------------------------------------------------------------------------- subroutine multibody(ecorr) C This subroutine calculates multi-body contributions to energy following @@ -8721,6 +8803,8 @@ C contribution equal to sqrt(eps(i,j)*eps(i+1,j+1)) is added. include 'COMMON.DERIV' include 'COMMON.INTERACT' include 'COMMON.CONTACTS' + include 'COMMON.CONTMAT' + include 'COMMON.CORRMAT' double precision gx(3),gx1(3) logical lprn @@ -8775,6 +8859,8 @@ c------------------------------------------------------------------------------ include 'COMMON.DERIV' include 'COMMON.INTERACT' include 'COMMON.CONTACTS' + include 'COMMON.CONTMAT' + include 'COMMON.CORRMAT' include 'COMMON.SHIELD' double precision gx(3),gx1(3) logical lprn @@ -8829,6 +8915,8 @@ C This subroutine calculates multi-body contributions to hydrogen-bonding include 'COMMON.DERIV' include 'COMMON.INTERACT' include 'COMMON.CONTACTS' + include 'COMMON.CONTMAT' + include 'COMMON.CORRMAT' include 'COMMON.CONTROL' include 'COMMON.LOCAL' double precision gx(3),gx1(3),time00 @@ -9122,6 +9210,8 @@ c------------------------------------------------------------------------------ parameter (max_cont=maxconts) parameter (max_dim=26) include "COMMON.CONTACTS" + include 'COMMON.CONTMAT' + include 'COMMON.CORRMAT' double precision zapas(max_dim,maxconts,max_fg_procs), & zapas_recv(max_dim,maxconts,max_fg_procs) common /przechowalnia/ zapas @@ -9193,6 +9283,8 @@ C This subroutine calculates multi-body contributions to hydrogen-bonding include 'COMMON.LOCAL' include 'COMMON.INTERACT' include 'COMMON.CONTACTS' + include 'COMMON.CONTMAT' + include 'COMMON.CORRMAT' include 'COMMON.CHAIN' include 'COMMON.CONTROL' include 'COMMON.SHIELD' @@ -9563,6 +9655,8 @@ c------------------------------------------------------------------------------ parameter (max_cont=maxconts) parameter (max_dim=70) include "COMMON.CONTACTS" + include 'COMMON.CONTMAT' + include 'COMMON.CORRMAT' double precision zapas(max_dim,maxconts,max_fg_procs), & zapas_recv(max_dim,maxconts,max_fg_procs) common /przechowalnia/ zapas @@ -9616,6 +9710,8 @@ c------------------------------------------------------------------------------ include 'COMMON.DERIV' include 'COMMON.INTERACT' include 'COMMON.CONTACTS' + include 'COMMON.CONTMAT' + include 'COMMON.CORRMAT' include 'COMMON.SHIELD' include 'COMMON.CONTROL' double precision gx(3),gx1(3) @@ -9791,6 +9887,8 @@ C--------------------------------------------------------------------------- include 'COMMON.DERIV' include 'COMMON.INTERACT' include 'COMMON.CONTACTS' + include 'COMMON.CONTMAT' + include 'COMMON.CORRMAT' include 'COMMON.TORSION' include 'COMMON.VAR' include 'COMMON.GEO' @@ -9856,6 +9954,8 @@ C include 'COMMON.DERIV' include 'COMMON.INTERACT' include 'COMMON.CONTACTS' + include 'COMMON.CONTMAT' + include 'COMMON.CORRMAT' include 'COMMON.TORSION' include 'COMMON.VAR' include 'COMMON.GEO' @@ -10242,6 +10342,8 @@ C--------------------------------------------------------------------------- include 'COMMON.DERIV' include 'COMMON.INTERACT' include 'COMMON.CONTACTS' + include 'COMMON.CONTMAT' + include 'COMMON.CORRMAT' include 'COMMON.TORSION' include 'COMMON.VAR' include 'COMMON.GEO' @@ -10363,6 +10465,8 @@ C--------------------------------------------------------------------------- include 'COMMON.DERIV' include 'COMMON.INTERACT' include 'COMMON.CONTACTS' + include 'COMMON.CONTMAT' + include 'COMMON.CORRMAT' include 'COMMON.TORSION' include 'COMMON.VAR' include 'COMMON.GEO' @@ -10767,6 +10871,8 @@ c-------------------------------------------------------------------------- include 'COMMON.DERIV' include 'COMMON.INTERACT' include 'COMMON.CONTACTS' + include 'COMMON.CONTMAT' + include 'COMMON.CORRMAT' include 'COMMON.TORSION' include 'COMMON.VAR' include 'COMMON.GEO' @@ -10907,6 +11013,8 @@ c-------------------------------------------------------------------------- include 'COMMON.DERIV' include 'COMMON.INTERACT' include 'COMMON.CONTACTS' + include 'COMMON.CONTMAT' + include 'COMMON.CORRMAT' include 'COMMON.TORSION' include 'COMMON.VAR' include 'COMMON.GEO' @@ -11011,6 +11119,8 @@ c---------------------------------------------------------------------------- include 'COMMON.DERIV' include 'COMMON.INTERACT' include 'COMMON.CONTACTS' + include 'COMMON.CONTMAT' + include 'COMMON.CORRMAT' include 'COMMON.TORSION' include 'COMMON.VAR' include 'COMMON.GEO' @@ -11196,6 +11306,8 @@ c---------------------------------------------------------------------------- include 'COMMON.DERIV' include 'COMMON.INTERACT' include 'COMMON.CONTACTS' + include 'COMMON.CONTMAT' + include 'COMMON.CORRMAT' include 'COMMON.TORSION' include 'COMMON.VAR' include 'COMMON.GEO' @@ -11311,6 +11423,8 @@ c---------------------------------------------------------------------------- include 'COMMON.DERIV' include 'COMMON.INTERACT' include 'COMMON.CONTACTS' + include 'COMMON.CONTMAT' + include 'COMMON.CORRMAT' include 'COMMON.TORSION' include 'COMMON.VAR' include 'COMMON.GEO' @@ -11555,6 +11669,8 @@ c---------------------------------------------------------------------------- include 'COMMON.DERIV' include 'COMMON.INTERACT' include 'COMMON.CONTACTS' + include 'COMMON.CONTMAT' + include 'COMMON.CORRMAT' include 'COMMON.TORSION' include 'COMMON.VAR' include 'COMMON.GEO' @@ -11873,8 +11989,8 @@ cd write (2,*) 'ekont',ekont cd write (2,*) 'eel_turn6',ekont*eel_turn6 return end - C----------------------------------------------------------------------------- +#endif double precision function scalar(u,v) !DIR$ INLINEALWAYS scalar #ifndef OSF