X-Git-Url: http://mmka.chem.univ.gda.pl/gitweb/?a=blobdiff_plain;f=source%2Funres%2Fsrc_MD%2Fparmread.F;h=bfb4c220b2dd0cfb06b99bf775744d98ad16ce7f;hb=49b4d93d38312bdd6e5e9c6bb845431396c9b0f5;hp=074db7c4bb5d8874f481ff998bdd5a5c94b3ec5e;hpb=a1174757c39d2d9dc4579febaeb2e7d3d84602a3;p=unres.git diff --git a/source/unres/src_MD/parmread.F b/source/unres/src_MD/parmread.F index 074db7c..bfb4c22 100644 --- a/source/unres/src_MD/parmread.F +++ b/source/unres/src_MD/parmread.F @@ -106,9 +106,9 @@ C read (ithep,*,err=111,end=111) a0thet(i),(athet(j,i,1,1),j=1,2), & (bthet(j,i,1,1),j=1,2) read (ithep,*,err=111,end=111) (polthet(j,i),j=0,3) - read (ithep,*,err=111,end=111) (gthet(j,i),j=1,3) - read (ithep,*,err=111,end=111) theta0(i),sig0(i),sigc0(i) - sigc0(i)=sigc0(i)**2 + read (ithep,*,err=111,end=111) (gthet(j,i),j=1,3) + read (ithep,*,err=111,end=111) theta0(i),sig0(i),sigc0(i) + sigc0(i)=sigc0(i)**2 enddo do i=1,ntyp athet(1,i,1,-1)=athet(1,i,1,1) @@ -144,6 +144,7 @@ C gthet(j,i)=gthet(j,-i) enddo enddo + close (ithep) if (lprint) then if (.not.LaTeX) then @@ -209,6 +210,7 @@ C Kozlowska et al., J. Phys.: Condens. Matter 19 (2007) 285203 C read (ithep,*,err=111,end=111) nthetyp,ntheterm,ntheterm2, & ntheterm3,nsingle,ndouble +C print *, "tu" nntheterm=max0(ntheterm,ntheterm2,ntheterm3) read (ithep,*,err=111,end=111) (ithetyp(i),i=1,ntyp1) do i=-ntyp1,-1 @@ -248,6 +250,7 @@ c VAR:ntethtyp is type of theta potentials type currently 0=glycine c VAR:1=non-glicyne non-proline 2=proline c VAR:negative values for D-aminoacid do i=0,nthetyp +C print *,i do j=-nthetyp,nthetyp do k=-nthetyp,nthetyp read (ithep,'(6a)',end=111,err=111) res1 @@ -460,18 +463,10 @@ C bsc(1,i)=0.0D0 read(irotam,*,end=112,err=112)(censc(k,1,i),k=1,3), & ((blower(k,l,1),l=1,k),k=1,3) - censc(1,1,-i)=censc(1,1,i) - censc(2,1,-i)=censc(2,1,i) - censc(3,1,-i)=-censc(3,1,i) - do j=2,nlob(i) read (irotam,*,end=112,err=112) bsc(j,i) read (irotam,*,end=112,err=112) (censc(k,j,i),k=1,3), & ((blower(k,l,j),l=1,k),k=1,3) - censc(1,j,-i)=censc(1,j,i) - censc(2,j,-i)=censc(2,j,i) - censc(3,j,-i)=-censc(3,j,i) -C BSC is amplitude of Gaussian enddo do j=1,nlob(i) do k=1,3 @@ -482,14 +477,6 @@ C BSC is amplitude of Gaussian enddo gaussc(k,l,j,i)=akl gaussc(l,k,j,i)=akl - if (((k.eq.3).and.(l.ne.3)) - & .or.((l.eq.3).and.(k.ne.3))) then - gaussc(k,l,j,-i)=-akl - gaussc(l,k,j,-i)=-akl - else - gaussc(k,l,j,-i)=akl - gaussc(l,k,j,-i)=akl - endif enddo enddo enddo @@ -617,40 +604,25 @@ C Read torsional parameters C read (itorp,*,end=113,err=113) ntortyp read (itorp,*,end=113,err=113) (itortyp(i),i=1,ntyp) - do iblock=1,2 - do i=-ntyp,-1 - itortyp(i)=-itortyp(-i) - enddo c write (iout,*) 'ntortyp',ntortyp - do i=0,ntortyp-1 - do j=-ntortyp+1,ntortyp-1 - read (itorp,*,end=113,err=113) nterm(i,j,iblock), - & nlor(i,j,iblock) - nterm(-i,-j,iblock)=nterm(i,j,iblock) - nlor(-i,-j,iblock)=nlor(i,j,iblock) + do i=1,ntortyp + do j=1,ntortyp + read (itorp,*,end=113,err=113) nterm(i,j),nlor(i,j) v0ij=0.0d0 si=-1.0d0 - do k=1,nterm(i,j,iblock) - read (itorp,*,end=113,err=113) kk,v1(k,i,j,iblock), - & v2(k,i,j,iblock) - v1(k,-i,-j,iblock)=v1(k,i,j,iblock) - v2(k,-i,-j,iblock)=-v2(k,i,j,iblock) - v0ij=v0ij+si*v1(k,i,j,iblock) + do k=1,nterm(i,j) + read (itorp,*,end=113,err=113) kk,v1(k,i,j),v2(k,i,j) + v0ij=v0ij+si*v1(k,i,j) si=-si -c write(iout,*) i,j,k,iblock,nterm(i,j,iblock) -c write(iout,*) v1(k,-i,-j,iblock),v1(k,i,j,iblock), -c &v2(k,-i,-j,iblock),v2(k,i,j,iblock) enddo - do k=1,nlor(i,j,iblock) + do k=1,nlor(i,j) read (itorp,*,end=113,err=113) kk,vlor1(k,i,j), & vlor2(k,i,j),vlor3(k,i,j) v0ij=v0ij+vlor1(k,i,j)/(1+vlor3(k,i,j)**2) enddo - v0(i,j,iblock)=v0ij - v0(-i,-j,iblock)=v0ij + v0(i,j)=v0ij enddo enddo - enddo close (itorp) if (lprint) then write (iout,'(/a/)') 'Torsional constants:' @@ -658,12 +630,11 @@ c &v2(k,-i,-j,iblock),v2(k,i,j,iblock) do j=1,ntortyp write (iout,*) 'ityp',i,' jtyp',j write (iout,*) 'Fourier constants' - do k=1,nterm(i,j,iblock) - write (iout,'(2(1pe15.5))') v1(k,i,j,iblock), - & v2(k,i,j,iblock) + do k=1,nterm(i,j) + write (iout,'(2(1pe15.5))') v1(k,i,j),v2(k,i,j) enddo write (iout,*) 'Lorenz constants' - do k=1,nlor(i,j,iblock) + do k=1,nlor(i,j) write (iout,'(3(1pe15.5))') & vlor1(k,i,j),vlor2(k,i,j),vlor3(k,i,j) enddo @@ -673,16 +644,15 @@ c &v2(k,-i,-j,iblock),v2(k,i,j,iblock) C C 6/23/01 Read parameters for double torsionals C - do iblock=1,2 - do i=0,ntortyp-1 - do j=-ntortyp+1,ntortyp-1 - do k=-ntortyp+1,ntortyp-1 + do i=1,ntortyp + do j=1,ntortyp + do k=1,ntortyp read (itordp,'(3a1)',end=114,err=114) t1,t2,t3 -c write (iout,*) "OK onelett", -c & i,j,k,t1,t2,t3 + write (iout,*) "OK onelett", + & i,j,k,t1,t2,t3 - if (t1.ne.toronelet(i) .or. t2.ne.toronelet(j) - & .or. t3.ne.toronelet(k)) then + if (t1.ne.onelett(i) .or. t2.ne.onelett(j) + & .or. t3.ne.onelett(k)) then write (iout,*) "Error in double torsional parameter file", & i,j,k,t1,t2,t3 #ifdef MPI @@ -690,44 +660,22 @@ c & i,j,k,t1,t2,t3 #endif stop "Error in double torsional parameter file" endif - read (itordp,*,end=114,err=114) ntermd_1(i,j,k,iblock), - & ntermd_2(i,j,k,iblock) - ntermd_1(-i,-j,-k,iblock)=ntermd_1(i,j,k,iblock) - ntermd_2(-i,-j,-k,iblock)=ntermd_2(i,j,k,iblock) - read (itordp,*,end=114,err=114) (v1c(1,l,i,j,k,iblock),l=1, - & ntermd_1(i,j,k,iblock)) - read (itordp,*,end=114,err=114) (v1s(1,l,i,j,k,iblock),l=1, - & ntermd_1(i,j,k,iblock)) - read (itordp,*,end=114,err=114) (v1c(2,l,i,j,k,iblock),l=1, - & ntermd_1(i,j,k,iblock)) - read (itordp,*,end=114,err=114) (v1s(2,l,i,j,k,iblock),l=1, - & ntermd_1(i,j,k,iblock)) -C Martix of D parameters for one dimesional foureir series - do l=1,ntermd_1(i,j,k,iblock) - v1c(1,l,-i,-j,-k,iblock)=v1c(1,l,i,j,k,iblock) - v1s(1,l,-i,-j,-k,iblock)=-v1s(1,l,i,j,k,iblock) - v1c(2,l,-i,-j,-k,iblock)=v1c(2,l,i,j,k,iblock) - v1s(2,l,-i,-j,-k,iblock)=-v1s(2,l,i,j,k,iblock) -c write(iout,*) "whcodze" , -c & v1s(2,l,-i,-j,-k,iblock),v1s(2,l,i,j,k,iblock) - enddo - read (itordp,*,end=114,err=114) ((v2c(l,m,i,j,k,iblock), - & v2c(m,l,i,j,k,iblock),v2s(l,m,i,j,k,iblock), - & v2s(m,l,i,j,k,iblock), - & m=1,l-1),l=1,ntermd_2(i,j,k,iblock)) -C Martix of D parameters for two dimesional fourier series - do l=1,ntermd_2(i,j,k,iblock) - do m=1,l-1 - v2c(l,m,-i,-j,-k,iblock)=v2c(l,m,i,j,k,iblock) - v2c(m,l,-i,-j,-k,iblock)=v2c(m,l,i,j,k,iblock) - v2s(l,m,-i,-j,-k,iblock)=-v2s(l,m,i,j,k,iblock) - v2s(m,l,-i,-j,-k,iblock)=-v2s(m,l,i,j,k,iblock) - enddo!m - enddo!l - enddo!k - enddo!j - enddo!i - enddo!iblock + read (itordp,*,end=114,err=114) ntermd_1(i,j,k), + & ntermd_2(i,j,k) + read (itordp,*,end=114,err=114) (v1c(1,l,i,j,k),l=1, + & ntermd_1(i,j,k)) + read (itordp,*,end=114,err=114) (v1s(1,l,i,j,k),l=1, + & ntermd_1(i,j,k)) + read (itordp,*,end=114,err=114) (v1c(2,l,i,j,k),l=1, + & ntermd_1(i,j,k)) + read (itordp,*,end=114,err=114) (v1s(2,l,i,j,k),l=1, + & ntermd_1(i,j,k)) + read (itordp,*,end=114,err=114) ((v2c(l,m,i,j,k), + & v2c(m,l,i,j,k),v2s(l,m,i,j,k),v2s(m,l,i,j,k), + & m=1,l-1),l=1,ntermd_2(i,j,k)) + enddo + enddo + enddo if (lprint) then write (iout,*) write (iout,*) 'Constants for double torsionals' @@ -736,29 +684,26 @@ C Martix of D parameters for two dimesional fourier series do j=-ntortyp+1,ntortyp-1 do k=-ntortyp+1,ntortyp-1 write (iout,*) 'ityp',i,' jtyp',j,' ktyp',k, - & ' nsingle',ntermd_1(i,j,k,iblock), - & ' ndouble',ntermd_2(i,j,k,iblock) + & ' nsingle',ntermd_1(i,j,k),' ndouble',ntermd_2(i,j,k) write (iout,*) write (iout,*) 'Single angles:' - do l=1,ntermd_1(i,j,k,iblock) - write (iout,'(i5,2f10.5,5x,2f10.5,5x,2f10.5)') l, - & v1c(1,l,i,j,k,iblock),v1s(1,l,i,j,k,iblock), - & v1c(2,l,i,j,k,iblock),v1s(2,l,i,j,k,iblock), - & v1s(1,l,-i,-j,-k,iblock),v1s(2,l,-i,-j,-k,iblock) + do l=1,ntermd_1(i,j,k) + write (iout,'(i5,2f10.5,5x,2f10.5)') l, + & v1c(1,l,i,j,k),v1s(1,l,i,j,k), + & v1c(2,l,i,j,k),v1s(2,l,i,j,k) enddo write (iout,*) write (iout,*) 'Pairs of angles:' - write (iout,'(3x,20i10)') (l,l=1,ntermd_2(i,j,k,iblock)) - do l=1,ntermd_2(i,j,k,iblock) + write (iout,'(3x,20i10)') (l,l=1,ntermd_2(i,j,k)) + do l=1,ntermd_2(i,j,k) write (iout,'(i5,20f10.5)') - & l,(v2c(l,m,i,j,k,iblock),m=1,ntermd_2(i,j,k,iblock)) + & l,(v2c(l,m,i,j,k),m=1,ntermd_2(i,j,k)) enddo write (iout,*) - write (iout,'(3x,20i10)') (l,l=1,ntermd_2(i,j,k,iblock)) - do l=1,ntermd_2(i,j,k,iblock) + write (iout,'(3x,20i10)') (l,l=1,ntermd_2(i,j,k)) + do l=1,ntermd_2(i,j,k) write (iout,'(i5,20f10.5)') - & l,(v2s(l,m,i,j,k,iblock),m=1,ntermd_2(i,j,k,iblock)), - & (v2s(l,m,-i,-j,-k,iblock),m=1,ntermd_2(i,j,k,iblock)) + & l,(v2s(l,m,i,j,k),m=1,ntermd_2(i,j,k)) enddo write (iout,*) enddo @@ -771,9 +716,9 @@ C Read of Side-chain backbone correlation parameters C Modified 11 May 2012 by Adasko CCC C - read (isccor,*,end=119,err=119) nsccortyp + read (isccor,*,end=1113,err=1113) nsccortyp #ifdef SCCORPDB - read (isccor,*,end=119,err=119) (isccortyp(i),i=1,ntyp) + read (isccor,*,end=1113,err=1113) (isccortyp(i),i=1,ntyp) do i=-ntyp,-1 isccortyp(i)=-isccortyp(-i) enddo @@ -783,8 +728,10 @@ c write (iout,*) 'ntortyp',ntortyp cc maxinter is maximum interaction sites do l=1,maxinter do i=1,nsccortyp - do j=1,nsccortyp - read (isccor,*,end=119,err=119) nterm_sccor(i,j),nlor_sccor(i,j) + do j=1,nsccortyp + read (isccor,*,end=1113,err=1113) nterm_sccor(i,j), + & nlor_sccor(i,j) + print *,i,j,l v0ijsccor=0.0d0 v0ijsccor1=0.0d0 v0ijsccor2=0.0d0 @@ -793,8 +740,8 @@ cc maxinter is maximum interaction sites nterm_sccor(-i,j)=nterm_sccor(i,j) nterm_sccor(-i,-j)=nterm_sccor(i,j) nterm_sccor(i,-j)=nterm_sccor(i,j) - do k=1,nterm_sccor(i,j) - read (isccor,*,end=119,err=119) kk,v1sccor(k,l,i,j) + do k=1,nterm_sccor(i,j) + read (isccor,*,end=119,err=119) kk,v1sccor(k,l,i,j) & ,v2sccor(k,l,i,j) if (j.eq.iscprol) then if (i.eq.isccortyp(10)) then @@ -828,14 +775,16 @@ cc maxinter is maximum interaction sites endif endif endif +C read (isccor,*,end=1113,err=1113) kk,v1sccor(k,l,i,j) +C & ,v2sccor(k,l,i,j) v0ijsccor=v0ijsccor+si*v1sccor(k,l,i,j) v0ijsccor1=v0ijsccor+si*v1sccor(k,l,-i,j) v0ijsccor2=v0ijsccor+si*v1sccor(k,l,i,-j) v0ijsccor3=v0ijsccor+si*v1sccor(k,l,-i,-j) si=-si enddo - do k=1,nlor_sccor(i,j) - read (isccor,*,end=119,err=119) kk,vlor1sccor(k,i,j), + do k=1,nlor_sccor(i,j) + read (isccor,*,end=1113,err=1113) kk,vlor1sccor(k,i,j), & vlor2sccor(k,i,j),vlor3sccor(k,i,j) v0ijsccor=v0ijsccor+vlor1sccor(k,i,j)/ &(1+vlor3sccor(k,i,j)**2) @@ -907,7 +856,7 @@ C write (iout,*) "Coefficients of the cumulants" endif read (ifourier,*) nloctyp - do i=0,nloctyp-1 + do i=1,nloctyp read (ifourier,*,end=115,err=115) read (ifourier,*,end=115,err=115) (b(ii),ii=1,13) if (lprint) then @@ -916,31 +865,22 @@ C endif B1(1,i) = b(3) B1(2,i) = b(5) - B1(1,-i) = b(3) - B1(2,-i) = -b(5) c b1(1,i)=0.0d0 c b1(2,i)=0.0d0 B1tilde(1,i) = b(3) B1tilde(2,i) =-b(5) - B1tilde(1,-i) =-b(3) - B1tilde(2,-i) =b(5) +C B1tilde(1,-i) =-b(3) +C B1tilde(2,-i) =b(5) c b1tilde(1,i)=0.0d0 c b1tilde(2,i)=0.0d0 B2(1,i) = b(2) B2(2,i) = b(4) - B2(1,-i) =b(2) - B2(2,-i) =-b(4) - c b2(1,i)=0.0d0 c b2(2,i)=0.0d0 CC(1,1,i)= b(7) CC(2,2,i)=-b(7) CC(2,1,i)= b(9) CC(1,2,i)= b(9) - CC(1,1,-i)= b(7) - CC(2,2,-i)=-b(7) - CC(2,1,-i)=-b(9) - CC(1,2,-i)=-b(9) c CC(1,1,i)=0.0d0 c CC(2,2,i)=0.0d0 c CC(2,1,i)=0.0d0 @@ -949,11 +889,6 @@ c CC(1,2,i)=0.0d0 Ctilde(1,2,i)=b(9) Ctilde(2,1,i)=-b(9) Ctilde(2,2,i)=b(7) - Ctilde(1,1,-i)=b(7) - Ctilde(1,2,-i)=-b(9) - Ctilde(2,1,-i)=b(9) - Ctilde(2,2,-i)=b(7) - c Ctilde(1,1,i)=0.0d0 c Ctilde(1,2,i)=0.0d0 c Ctilde(2,1,i)=0.0d0 @@ -962,10 +897,6 @@ c Ctilde(2,2,i)=0.0d0 DD(2,2,i)=-b(6) DD(2,1,i)= b(8) DD(1,2,i)= b(8) - DD(1,1,-i)= b(6) - DD(2,2,-i)=-b(6) - DD(2,1,-i)=-b(8) - DD(1,2,-i)=-b(8) c DD(1,1,i)=0.0d0 c DD(2,2,i)=0.0d0 c DD(2,1,i)=0.0d0 @@ -974,11 +905,6 @@ c DD(1,2,i)=0.0d0 Dtilde(1,2,i)=b(8) Dtilde(2,1,i)=-b(8) Dtilde(2,2,i)=b(6) - Dtilde(1,1,-i)=b(6) - Dtilde(1,2,-i)=-b(8) - Dtilde(2,1,-i)=b(8) - Dtilde(2,2,-i)=b(6) - c Dtilde(1,1,i)=0.0d0 c Dtilde(1,2,i)=0.0d0 c Dtilde(2,1,i)=0.0d0 @@ -987,11 +913,6 @@ c Dtilde(2,2,i)=0.0d0 EE(2,2,i)=-b(10)+b(11) EE(2,1,i)= b(12)-b(13) EE(1,2,i)= b(12)+b(13) - EE(1,1,-i)= b(10)+b(11) - EE(2,2,-i)=-b(10)+b(11) - EE(2,1,-i)=-b(12)+b(13) - EE(1,2,-i)=-b(12)-b(13) - c ee(1,1,i)=1.0d0 c ee(2,2,i)=1.0d0 c ee(2,1,i)=0.0d0 @@ -1296,6 +1217,9 @@ c v3ss=0.0d0 goto 999 113 write (iout,*) "Error reading torsional energy parameters." goto 999 + 1113 write (iout,*) + & "Error reading side-chain torsional energy parameters." + goto 999 114 write (iout,*) "Error reading double torsional energy parameters." goto 999 115 write (iout,*)