From c4eebcc1774bd6bfe98a1f484f697ab334d1782e Mon Sep 17 00:00:00 2001 From: Cezary Czaplewski Date: Tue, 8 Mar 2016 18:03:24 +0900 Subject: [PATCH] corrections to src_MD-M for the same energy as src_MD --- source/unres/src_MD-M/COMMON.IOUNITS | 10 ++-- source/unres/src_MD-M/energy_p_new_barrier.F | 14 +++-- source/unres/src_MD-M/initialize_p.F | 2 + source/unres/src_MD-M/parmread.F | 77 ++++++++++++++++++++++---- source/unres/src_MD-M/readrtns_CSA.F | 45 ++++++++++++++- 5 files changed, 124 insertions(+), 24 deletions(-) diff --git a/source/unres/src_MD-M/COMMON.IOUNITS b/source/unres/src_MD-M/COMMON.IOUNITS index a9ace0b..49b6db3 100644 --- a/source/unres/src_MD-M/COMMON.IOUNITS +++ b/source/unres/src_MD-M/COMMON.IOUNITS @@ -11,11 +11,11 @@ C General I/O units & files integer inp,iout,igeom,intin,ipdb,imol2,ipdbin,ithep,irotam, & itorp,itordp,ifourier,ielep,isidep,iscpp,icbase,istat, & ientin,ientout,izs1,isecpred,ibond,irest2,iifrag,icart, - & irest1,isccor + & irest1,isccor,ithep_pdb,irotam_pdb common /iounits/ inp,iout,igeom,intin,ipdb,imol2,ipdbin,ithep, & irotam,itorp,itordp,ifourier,ielep,isidep,iscpp,icbase, & istat,ientin,ientout,izs1,isecpred,ibond,irest2,iifrag, - & icart,irest1,isccor + & icart,irest1,isccor,ithep_pdb,irotam_pdb character*256 outname,intname,pdbname,mol2name,statname,intinname, & entname,prefix,secpred,rest2name,qname,cartname,tmpdir, & mremd_rst_name,curdir,pref_orig @@ -38,9 +38,11 @@ C CSA I/O units & files & icsa_bank_reminimized,icsa_native_int,icsa_in,icsa_pdb C Parameter files character*256 bondname,thetname,rotname,torname,tordname, - & fouriername,elename,sidename,scpname,sccorname,patname + & fouriername,elename,sidename,scpname,sccorname,patname, + & thetname_pdb,rotname_pdb common /parfiles/ bondname,thetname,rotname,torname,tordname, - & fouriername,elename,sidename,scpname,sccorname,patname + & fouriername,elename,sidename,scpname,sccorname,patname, + & thetname_pdb,rotname_pdb character*3 pot C----------------------------------------------------------------------- C INP - main input file diff --git a/source/unres/src_MD-M/energy_p_new_barrier.F b/source/unres/src_MD-M/energy_p_new_barrier.F index eb1cbd4..9ea982a 100644 --- a/source/unres/src_MD-M/energy_p_new_barrier.F +++ b/source/unres/src_MD-M/energy_p_new_barrier.F @@ -4546,7 +4546,8 @@ C logical lprn /.false./, lprn1 /.false./ etheta=0.0D0 do i=ithet_start,ithet_end - if (itype(i-1).eq.21) cycle + if ((itype(i-1).eq.ntyp1).or.(itype(i-2).eq.ntyp1).or. + &(itype(i).eq.ntyp1)) cycle dethetai=0.0d0 dephii=0.0d0 dephii1=0.0d0 @@ -4556,7 +4557,8 @@ C coskt(k)=dcos(k*theti2) sinkt(k)=dsin(k*theti2) enddo - if (i.gt.3 .and. itype(i-2).ne.21) then +C if (i.gt.3) then + if (i.gt.3 .and. itype(imax0(i-3,1)).ne.ntyp1) then #ifdef OSF phii=phi(i) if (phii.ne.phii) phii=150.0 @@ -4570,13 +4572,13 @@ C enddo else phii=0.0d0 - ityp1=nthetyp+1 + ityp1=ithetyp(itype(i-2)) do k=1,nsingle cosph1(k)=0.0d0 sinph1(k)=0.0d0 enddo endif - if (i.lt.nres .and. itype(i).ne.21) then + if ((i.lt.nres).and. itype(i+1).ne.ntyp1) then #ifdef OSF phii1=phi(i+1) if (phii1.ne.phii1) phii1=150.0 @@ -4591,7 +4593,7 @@ C enddo else phii1=0.0d0 - ityp3=nthetyp+1 + ityp3=ithetyp(itype(i)) do k=1,nsingle cosph2(k)=0.0d0 sinph2(k)=0.0d0 @@ -4701,6 +4703,8 @@ C & i,theta(i)*rad2deg,phii*rad2deg, & phii1*rad2deg,ethetai etheta=etheta+ethetai + if (energy_dec) write (iout,'(a6,i5,0pf7.3)') + & 'ebend',i,ethetai if (i.gt.3) gloc(i-3,icg)=gloc(i-3,icg)+wang*dephii if (i.lt.nres) gloc(i-2,icg)=gloc(i-2,icg)+wang*dephii1 gloc(nphi+i-2,icg)=wang*dethetai diff --git a/source/unres/src_MD-M/initialize_p.F b/source/unres/src_MD-M/initialize_p.F index ae884df..37692d4 100644 --- a/source/unres/src_MD-M/initialize_p.F +++ b/source/unres/src_MD-M/initialize_p.F @@ -80,7 +80,9 @@ C igeom= 8 intin= 9 ithep= 11 + ithep_pdb=51 irotam=12 + irotam_pdb=52 itorp= 13 itordp= 23 ielep= 14 diff --git a/source/unres/src_MD-M/parmread.F b/source/unres/src_MD-M/parmread.F index 64203a2..183da9d 100644 --- a/source/unres/src_MD-M/parmread.F +++ b/source/unres/src_MD-M/parmread.F @@ -276,6 +276,17 @@ C enddo call flush(iout) endif + write (2,*) "Start reading THETA_PDB" + do i=1,ntyp + read (ithep_pdb,*,err=111,end=111) a0thet(i),(athet(j,i),j=1,2), + & (bthet(j,i),j=1,2) + read (ithep_pdb,*,err=111,end=111) (polthet(j,i),j=0,3) + read (ithep_pdb,*,err=111,end=111) (gthet(j,i),j=1,3) + read (ithep_pdb,*,err=111,end=111) theta0(i),sig0(i),sigc0(i) + sigc0(i)=sigc0(i)**2 + enddo + write (2,*) "End reading THETA_PDB" + close (ithep_pdb) #endif close(ithep) #ifdef CRYST_SC @@ -363,6 +374,48 @@ C enddo endif enddo +C +C Read the parameters of the probability distribution/energy expression +C of the side chains. +C + do i=1,ntyp + read (irotam_pdb,'(3x,i3,f8.3)',end=112,err=112) nlob(i),dsc(i) + if (i.eq.10) then + dsc_inv(i)=0.0D0 + else + dsc_inv(i)=1.0D0/dsc(i) + endif + if (i.ne.10) then + do j=1,nlob(i) + do k=1,3 + do l=1,3 + blower(l,k,j)=0.0D0 + enddo + enddo + enddo + bsc(1,i)=0.0D0 + read(irotam_pdb,*,end=112,err=112)(censc(k,1,i),k=1,3), + & ((blower(k,l,1),l=1,k),k=1,3) + do j=2,nlob(i) + read (irotam_pdb,*,end=112,err=112) bsc(j,i) + read (irotam_pdb,*,end=112,err=112) (censc(k,j,i),k=1,3), + & ((blower(k,l,j),l=1,k),k=1,3) + enddo + do j=1,nlob(i) + do k=1,3 + do l=1,k + akl=0.0D0 + do m=1,3 + akl=akl+blower(k,m,j)*blower(l,m,j) + enddo + gaussc(k,l,j,i)=akl + gaussc(l,k,j,i)=akl + enddo + enddo + enddo + endif + enddo + close (irotam_pdb) #endif close(irotam) @@ -590,23 +643,23 @@ cc maxinter is maximum interaction sites c write (iout,*) 'ntortyp',ntortyp maxinter=3 cc maxinter is maximum interaction sites - do l=1,maxinter + do l=1,maxinter do i=1,nsccortyp - do j=1,nsccortyp + do j=1,nsccortyp read (isccor,*,end=113,err=113) & nterm_sccor(i,j),nlor_sccor(i,j) v0ijsccor=0.0d0 si=-1.0d0 - - do k=1,nterm_sccor(i,j) + + do k=1,nterm_sccor(i,j) read (isccor,*,end=113,err=113) kk,v1sccor(k,l,i,j) - & ,v2sccor(k,l,i,j) + & ,v2sccor(k,l,i,j) v0ijsccor=v0ijsccor+si*v1sccor(k,l,i,j) si=-si enddo - do k=1,nlor_sccor(i,j) + do k=1,nlor_sccor(i,j) read (isccor,*,end=113,err=113) kk,vlor1sccor(k,i,j), - & vlor2sccor(k,i,j),vlor3sccor(k,i,j) + & vlor2sccor(k,i,j),vlor3sccor(k,i,j) v0ijsccor=v0ijsccor+vlor1sccor(k,i,j)/ &(1+vlor3sccor(k,i,j)**2) enddo @@ -615,12 +668,12 @@ cc maxinter is maximum interaction sites enddo enddo close (isccor) - + #endif if (lprint) then - write (iout,'(/a/)') 'Torsional constants:' - do i=1,nsccortyp - do j=1,nsccortyp + write (iout,'(/a/)') 'Torsional constants:' + do i=1,nsccortyp + do j=1,nsccortyp write (iout,*) 'ityp',i,' jtyp',j write (iout,*) 'Fourier constants' do k=1,nterm_sccor(i,j) @@ -628,7 +681,7 @@ cc maxinter is maximum interaction sites enddo write (iout,*) 'Lorenz constants' do k=1,nlor_sccor(i,j) - write (iout,'(3(1pe15.5))') + write (iout,'(3(1pe15.5))') & vlor1sccor(k,i,j),vlor2sccor(k,i,j),vlor3sccor(k,i,j) enddo enddo diff --git a/source/unres/src_MD-M/readrtns_CSA.F b/source/unres/src_MD-M/readrtns_CSA.F index 0b08ab7..bd17218 100644 --- a/source/unres/src_MD-M/readrtns_CSA.F +++ b/source/unres/src_MD-M/readrtns_CSA.F @@ -1862,8 +1862,16 @@ C Get parameter filenames and open the parameter files. open (ibond,file=bondname,status='old',readonly,shared) call getenv_loc('THETPAR',thetname) open (ithep,file=thetname,status='old',readonly,shared) +#ifndef CRYST_THETA + call getenv_loc('THETPARPDB',thetname_pdb) + open (ithep_pdb,file=thetname_pdb,status='old',readonly,shared) +#endif call getenv_loc('ROTPAR',rotname) open (irotam,file=rotname,status='old',readonly,shared) +#ifndef CRYST_SC + call getenv_loc('ROTPARPDB',rotname_pdb) + open (irotam_pdb,file=rotname_pdb,status='old',readonly,shared) +#endif call getenv_loc('TORPAR',torname) open (itorp,file=torname,status='old',readonly,shared) call getenv_loc('TORDPAR',tordname) @@ -1888,9 +1896,17 @@ c print *,"Processor",myrank," opened file IBOND" call getenv_loc('THETPAR',thetname) open (ithep,file=thetname,status='old',action='read') c print *,"Processor",myrank," opened file ITHEP" +#ifndef CRYST_THETA + call getenv_loc('THETPARPDB',thetname_pdb) + open (ithep_pdb,file=thetname_pdb,status='old',action='read') +#endif call getenv_loc('ROTPAR',rotname) open (irotam,file=rotname,status='old',action='read') c print *,"Processor",myrank," opened file IROTAM" +#ifndef CRYST_SC + call getenv_loc('ROTPARPDB',rotname_pdb) + open (irotam_pdb,file=rotname_pdb,status='old',action='read') +#endif call getenv_loc('TORPAR',torname) open (itorp,file=torname,status='old',action='read') c print *,"Processor",myrank," opened file ITORP" @@ -1919,8 +1935,16 @@ C Get parameter filenames and open the parameter files. open (ibond,file=bondname,status='old') call getenv_loc('THETPAR',thetname) open (ithep,file=thetname,status='old') +#ifndef CRYST_THETA + call getenv_loc('THETPARPDB',thetname_pdb) + open (ithep_pdb,file=thetname_pdb,status='old') +#endif call getenv_loc('ROTPAR',rotname) open (irotam,file=rotname,status='old') +#ifndef CRYST_SC + call getenv_loc('ROTPARPDB',rotname_pdb) + open (irotam_pdb,file=rotname_pdb,status='old') +#endif call getenv_loc('TORPAR',torname) open (itorp,file=torname,status='old') call getenv_loc('TORDPAR',tordname) @@ -1943,8 +1967,18 @@ C Get parameter filenames and open the parameter files. open (ibond,file=bondname,status='old',readonly) call getenv_loc('THETPAR',thetname) open (ithep,file=thetname,status='old',readonly) +#ifndef CRYST_THETA + call getenv_loc('THETPARPDB',thetname_pdb) + print *,"thetname_pdb ",thetname_pdb + open (ithep_pdb,file=thetname_pdb,status='old',action='read') + print *,ithep_pdb," opened" +#endif call getenv_loc('ROTPAR',rotname) open (irotam,file=rotname,status='old',readonly) +#ifndef CRYST_SC + call getenv_loc('ROTPARPDB',rotname_pdb) + open (irotam_pdb,file=rotname_pdb,status='old',action='read') +#endif call getenv_loc('TORPAR',torname) open (itorp,file=torname,status='old',readonly) call getenv_loc('TORDPAR',tordname) @@ -2388,6 +2422,9 @@ C Initialize random number generator C implicit real*8 (a-h,o-z) include 'DIMENSIONS' +#ifdef AMD64 + integer*8 iseedi8 +#endif #ifdef MPI include 'mpif.h' logical OKRandom, prng_restart @@ -2423,9 +2460,11 @@ C if (fg_rank.eq.0) then seed=seed*(me+1)+1 #ifdef AMD64 - if(me.eq.king) - & write (iout,*) 'MPI: node= ', me, ' iseed= ',iseed - OKRandom = prng_restart(me,iseed) + iseedi8=dint(seed) + if(me.eq.king .or. .not. out1file) + & write (iout,*) 'MPI: node= ', me, ' iseed= ',iseedi8 + write (*,*) 'MPI: node= ', me, ' iseed= ',iseedi8 + OKRandom = prng_restart(me,iseedi8) #else do i=1,4 tmp=65536.0d0**(4-i) -- 1.7.9.5