X-Git-Url: http://mmka.chem.univ.gda.pl/gitweb/?a=blobdiff_plain;f=source%2Funres%2Fsrc_MD-M%2FMD_A-MTS.F;fp=source%2Funres%2Fsrc_MD-M%2FMD_A-MTS.F;h=e0c1cfb3070c751a19bc7645b9f4f20d5eaf4d37;hb=fe9e3a2f1191191edd529c64aba73f2b98723a47;hp=6c6fb143c6ddb18cdde3ba4dac20bbcbc9411f9a;hpb=c89d7d207940ebb57e0042cb35ab65d0f8a4b55e;p=unres.git diff --git a/source/unres/src_MD-M/MD_A-MTS.F b/source/unres/src_MD-M/MD_A-MTS.F index 6c6fb14..e0c1cfb 100644 --- a/source/unres/src_MD-M/MD_A-MTS.F +++ b/source/unres/src_MD-M/MD_A-MTS.F @@ -209,7 +209,7 @@ c Variable time step algorithm. enddo enddo do i=nnt,nct - if (itype(i).ne.10 .and. itype(i).ne.21) then + if (itype(i).ne.10 .and. itype(i).ne.ntyp1) then do j=1,3 ind=ind+1 v_work(ind)=d_t(j,i+nres) @@ -232,6 +232,7 @@ c Variable time step algorithm. if (mod(itime,ntwx).eq.0) then write (tytul,'("time",f8.2)') totT if(mdpdb) then + call check_ecartint call hairpin(.true.,nharp,iharp) call secondary2(.true.) call pdbout(potE,tytul,ipdb) @@ -955,7 +956,7 @@ c forces). enddo enddo do i=nnt,nct - if (itype(i).ne.10 .and. itype(i).ne.21) then + if (itype(i).ne.10 .and. itype(i).ne.ntyp1) then inres=i+nres do j=1,3 d_t(j,inres)=d_t(j,inres)+0.5d0*d_a(j,inres)*d_time @@ -1005,7 +1006,7 @@ c Applying velocity Verlet algorithm - step 1 to coordinates enddo enddo do i=nnt,nct - if (itype(i).ne.10 .and. itype(i).ne.21) then + if (itype(i).ne.10 .and. itype(i).ne.ntyp1) then inres=i+nres do j=1,3 adt=d_a_old(j,inres)*d_time @@ -1049,7 +1050,7 @@ c Step 2 of the velocity Verlet algorithm: update velocities enddo enddo do i=nnt,nct - if (itype(i).ne.10 .and. itype(i).ne.21) then + if (itype(i).ne.10 .and. itype(i).ne.ntyp1) then inres=i+nres do j=1,3 d_t(j,inres)=d_t_new(j,inres)+0.5d0*d_a(j,inres)*d_time @@ -1148,7 +1149,7 @@ c ind=ind+3 enddo do i=nnt,nct - if (itype(i).ne.10 .and. itype(i).ne.21) then + if (itype(i).ne.10 .and. itype(i).ne.ntyp1) then inres=i+nres do j=1,3 adt=(d_a_old(j,inres)+d_af_work(ind+j))*d_time @@ -1213,7 +1214,7 @@ c ind=ind+3 enddo do i=nnt,nct - if (itype(i).ne.10 .and. itype(i).ne.21) then + if (itype(i).ne.10 .and. itype(i).ne.ntyp1) then inres=i+nres do j=1,3 d_t(j,inres)=d_t_new(j,inres)+(0.5d0*(d_a(j,inres) @@ -1260,6 +1261,7 @@ c accel(j)=aux(j)+0.5d0*(d_a(j,i)-d_a_old(j,i)) c if (dabs(accel(j)).gt.amax) amax=dabs(accel(j)) if (dabs(accel(j)).gt.dabs(accel_old(j))) then dacc=dabs(accel(j)-accel_old(j)) +c write (iout,*) i,dacc if (dacc.gt.amax) amax=dacc endif enddo @@ -1278,7 +1280,7 @@ c accel(j)=aux(j) enddo endif do i=nnt,nct - if (itype(i).ne.10 .and. itype(i).ne.21) then + if (itype(i).ne.10 .and. itype(i).ne.ntyp1) then do j=1,3 c accel(j)=accel(j)+d_a(j,i+nres)-d_a_old(j,i+nres) accel_old(j)=accel_old(j)+d_a_old(j,i+nres) @@ -1289,6 +1291,7 @@ c accel(j)=accel(j)+d_a(j,i+nres)-d_a_old(j,i+nres) c if (dabs(accel(j)).gt.amax) amax=dabs(accel(j)) if (dabs(accel(j)).gt.dabs(accel_old(j))) then dacc=dabs(accel(j)-accel_old(j)) +c write (iout,*) "side-chain",i,dacc if (dacc.gt.amax) amax=dacc endif enddo @@ -1331,7 +1334,7 @@ c write (iout,*) "back",i,j,epdriftij enddo endif c Side chains - if (itype(i).ne.10 .and. itype(i).ne.21) then + if (itype(i).ne.10 .and. itype(i).ne.ntyp1) then do j=1,3 epdriftij= & dabs((d_a(j,i+nres)-d_a_old(j,i+nres))*gxcart(j,i)) @@ -1378,7 +1381,7 @@ c write(iout,*) "fact", fact enddo enddo do i=nnt,nct - if (itype(i).ne.10 .and. itype(i).ne.21) then + if (itype(i).ne.10 .and. itype(i).ne.ntyp1) then inres=i+nres do j=1,3 d_t(j,inres)=fact*d_t(j,inres) @@ -1433,7 +1436,8 @@ c if the friction coefficients do not depend on surface area stdforcp(i)=stdfp*dsqrt(gamp) enddo do i=nnt,nct - stdforcsc(i)=stdfsc(itype(i))*dsqrt(gamsc(itype(i))) + stdforcsc(i)=stdfsc(iabs(itype(i))) + & *dsqrt(gamsc(iabs(itype(i)))) enddo endif c Open the pdb file for snapshotshots @@ -1605,7 +1609,7 @@ c Removing the velocity of the center of mass call chainbuild_cart call kinetic(EK) if (tbf) then - call verlet_bath(EK) + call verlet_bath endif kinetic_T=2.0d0/(dimen3*Rb)*EK if(me.eq.king.or..not.out1file)then @@ -1824,7 +1828,7 @@ c Transfer to the d_t vector enddo enddo do i=nnt,nct - if (itype(i).ne.10 .and. itype(i).ne.21) then + if (itype(i).ne.10 .and. itype(i).ne.ntyp1) then do j=1,3 ind=ind+1 d_t(j,i+nres)=d_t_work(ind) @@ -2053,7 +2057,7 @@ c enddo ind=ind+3 enddo do i=nnt,nct - if (itype(i).ne.10 .and. itype(i).ne.21) then + if (itype(i).ne.10 .and. itype(i).ne.ntyp1) then do j=1,3 dc_work(ind+j)=dc_old(j,i+nres) d_t_work(ind+j)=d_t_old(j,i+nres) @@ -2162,7 +2166,7 @@ c ind=ind+3 enddo do i=nnt,nct - if (itype(i).ne.10 .and. itype(i).ne.21) then + if (itype(i).ne.10 .and. itype(i).ne.ntyp1) then inres=i+nres do j=1,3 d_t(j,inres)=d_t_work(ind+j) @@ -2319,7 +2323,7 @@ c enddo ind=ind+3 enddo do i=nnt,nct - if (itype(i).ne.10 .and. itype(i).ne.21) then + if (itype(i).ne.10 .and. itype(i).ne.ntyp1) then do j=1,3 dc_work(ind+j)=dc_old(j,i+nres) d_t_work(ind+j)=d_t_old(j,i+nres) @@ -2368,7 +2372,7 @@ c enddo ind=ind+3 enddo do i=nnt,nct - if (itype(i).ne.10 .and. itype(i).ne.21) then + if (itype(i).ne.10 .and. itype(i).ne.ntyp1) then inres=i+nres do j=1,3 dc(j,inres)=dc_work(ind+j) @@ -2429,7 +2433,7 @@ c ddt2=ddt2+vrand_mat2(i,j)*stochforcvecV(j) ind=ind+3 enddo do i=nnt,nct - if (itype(i).ne.10 .and. itype(i).ne.21) then + if (itype(i).ne.10 .and. itype(i).ne.ntyp1) then inres=i+nres do j=1,3 d_t(j,inres)=d_t_work(ind+j)