X-Git-Url: http://mmka.chem.univ.gda.pl/gitweb/?a=blobdiff_plain;f=source%2Funres%2Fsrc_MD%2Fmd-diff%2Fmd%2Finit_md.f;fp=source%2Funres%2Fsrc_MD%2Fmd-diff%2Fmd%2Finit_md.f;h=0000000000000000000000000000000000000000;hb=0a11a2c4ccee14ed99ae44f2565b270ba8d4bbb6;hp=a7e2037ac96a60c78038972f47254656020e5baa;hpb=5eb407964903815242c59de10960f42761139e10;p=unres.git diff --git a/source/unres/src_MD/md-diff/md/init_md.f b/source/unres/src_MD/md-diff/md/init_md.f deleted file mode 100644 index a7e2037..0000000 --- a/source/unres/src_MD/md-diff/md/init_md.f +++ /dev/null @@ -1,189 +0,0 @@ -c--------------------------------------------------------- - subroutine init_MD -c Set up the initial conditions of a MD simulation - implicit real*8 (a-h,o-z) - include 'DIMENSIONS' -#ifdef MP - include 'mpif.h' - include 'COMMON.INFO' - include 'COMMON.SETUP' - character*4 liczba - character*16 form -#endif - include 'COMMON.CONTROL' - include 'COMMON.VAR' - include 'COMMON.MD' - include 'COMMON.LANGEVIN' - include 'COMMON.CHAIN' - include 'COMMON.DERIV' - include 'COMMON.GEO' - include 'COMMON.LOCAL' - include 'COMMON.INTERACT' - include 'COMMON.IOUNITS' - include 'COMMON.NAMES' - real*8 energia(0:n_ene),energia_long(0:n_ene), - & energia_short(0:n_ene),vcm(3),incr(3) - double precision cm(3),L(3),xv,sigv,lowb,highb - character*256 qstr - integer ilen - external ilen - character*50 tytul - common /gucio/ cm - d_time0=d_time - write(iout,*) "d_time", d_time -c Compute the standard deviations of stochastic forces for Langevin dynamics -c if the friction coefficients do not depend on surface area - if (lang.gt.0 .and. .not.surfarea) then - do i=nnt,nct-1 - stdforcp(i)=stdfp*dsqrt(gamp) - enddo - do i=nnt,nct - stdforcsc(i)=stdfsc(itype(i))*dsqrt(gamsc(itype(i))) - enddo - endif -c Open the pdb file for snapshotshots -#ifdef MPI - if (nprocs.eq.1) then - npos=3 - else - npos = dlog10(dfloat(nprocs-1))+1 - endif - if (npos.lt.3) npos=3 - write (liczba,'(i1)') npos - form = '(bz,i'//liczba(:ilen(liczba))//'.'//liczba(:ilen(liczba)) - & //')' - write (liczba,form) myrank - if(mdpdb) then - open(ipdb, - & file=prefix(:ilen(prefix))//"_MD"//liczba(:ilen(liczba)) - & //".pdb") - else - cartname=prefix(:ilen(prefix))//"_MD"//liczba(:ilen(liczba)) - & //".cx" - endif -#else - if(mdpdb) then - open(ipdb,file=prefix(:ilen(prefix))//"_MD.pdb") - else - cartname=prefix(:ilen(prefix))//"_MD.cx" - endif -#endif - if (usampl) then - write (qstr,'(256(1h ))') - ipos=1 - do i=1,nfrag - iq = qinfrag(i)*10 - iw = wfrag(i)/100 - if (iw.gt.0) then - write (iout,*) "Frag",qinfrag(i),wfrag(i),iq,iw - write (qstr(ipos:ipos+6),'(2h_f,i1,1h_,i1,1h_,i1)') i,iq,iw - ipos=ipos+7 - endif - enddo - do i=1,npair - iq = qinpair(i)*10 - iw = wpair(i)/100 - if (iw.gt.0) then - write (iout,*) "Pair",i,qinpair(i),wpair(i),iq,iw - write (qstr(ipos:ipos+6),'(2h_p,i1,1h_,i1,1h_,i1)') i,iq,iw - ipos=ipos+7 - endif - enddo - pdbname=pdbname(:ilen(pdbname)-4)//qstr(:ipos-1)//'.pdb' - cartname=cartname(:ilen(cartname)-2)//qstr(:ipos-1)//'.x' - statname=statname(:ilen(statname)-5)//qstr(:ipos-1)//'.stat' - endif - icg=1 - if (rest) then - write(iout,*) "Initial state will be read from file ", - & rest2name(:ilen(rest2name)) - call readrst - else -c Generate initial velocities - write(iout,*) "Initial velocities randomly generated" - call random_vel - totT=0.0d0 - endif -c rest2name = prefix(:ilen(prefix))//'.rst' - write(iout,*) "Initial backbone velocities" - do i=nnt,nct-1 - write(iout,*) (d_t(j,i),j=1,3) - enddo - write(iout,*) "Initial side-chain velocities" - do i=nnt,nct - write(iout,*) (d_t(j,i+nres),j=1,3) - enddo -c Zeroing the total angular momentum of the system - write(iout,*) "Calling the zero-angular - & momentum subroutine" - call inertia_tensor -c Getting the potential energy and forces and velocities and accelerations - call vcm_vel(vcm) -c write (iout,*) "velocity of the center of the mass:" -c write (iout,*) (vcm(j),j=1,3) - do j=1,3 - d_t(j,0)=d_t(j,0)-vcm(j) - enddo -c Removing the velocity of the center of mass - call vcm_vel(vcm) - write (iout,*) "vcm right after adjustment:" - write (iout,*) (vcm(j),j=1,3) - if (.not.rest) then - call chainbuild - endif - call chainbuild_cart - call kinetic(EK) - if (tbf) then - call verlet_bath(EK) - endif - kinetic_T=2.0d0/(dimen*Rb)*EK - call cartprint - call intout - call zerograd - call etotal(energia) - potE=energia(0) - call cartgrad - call lagrangian - call max_accel - if (amax*d_time .gt. dvmax) d_time=d_time*dvmax/amax - write(iout,*) "Potential energy" - write(iout,*) (energia(i),i=0,n_ene) - potE=energia(0)-energia(20) - totE=EK+potE - call statout(itime) - write (iout,*) "Initial:", - & " Kinetic energy",EK," potential energy",potE, - & " total energy",totE," maximum acceleration ", - & amax - if (large) then - write (iout,*) "Initial velocities" - do i=0,nres - write (iout,'(i3,3f10.5,3x,3f10.5)') i,(d_t(j,i),j=1,3), - & (d_t(j,i+nres),j=1,3) - enddo - write (iout,*) "Initial accelerations" - do i=0,nres - write (iout,'(i3,3f10.5,3x,3f10.5)') i,(d_a(j,i),j=1,3), - & (d_a(j,i+nres),j=1,3) - enddo - endif - do i=0,2*nres - do j=1,3 - dc_old(j,i)=dc(j,i) - d_t_old(j,i)=d_t(j,i) - d_a_old(j,i)=d_a(j,i) - enddo -c write (iout,*) "dc_old",i,(dc_old(j,i),j=1,3) - enddo - if (RESPA) then - call zerograd - call etotal_long(energia_long) - call cartgrad - call lagrangian -c call etotal_short(energia_short) -c write (iout,*) "energia_long",energia_long(0), -c & " energia_short",energia_short(0), -c & " total",energia_long(0)+energia_short(0) - endif - return - end