X-Git-Url: http://mmka.chem.univ.gda.pl/gitweb/?a=blobdiff_plain;f=source%2Funres%2Fio_config.F90;h=61cb7b474b80e9ffbc68a6808efb912211efd6dc;hb=94d31d7bcffc5412d9a88a0aef46b68349b60cf3;hp=5907752bf42ef0eee67a8af626e38c822965a862;hpb=8cd7afe8e87106d4cdd4b6361171845086971607;p=unres4.git diff --git a/source/unres/io_config.F90 b/source/unres/io_config.F90 index 5907752..61cb7b4 100644 --- a/source/unres/io_config.F90 +++ b/source/unres/io_config.F90 @@ -867,8 +867,13 @@ enddo enddo endif + + + + if (.not.allocated(ichargecat)) allocate (ichargecat(ntyp_molec(5))) + if (oldion.eq.1) then do i=1,ntyp_molec(5) - read(iion,*) msc(i,5),restok(i,5) + read(iion,*) msc(i,5),restok(i,5),ichargecat(i) print *,msc(i,5),restok(i,5) enddo ip(5)=0.2 @@ -879,6 +884,23 @@ read (iion,*) (catprm(i,k),i=1,ncatprotparm) enddo print *, catprm + endif + allocate(catnuclprm(14,ntyp_molec(2),ntyp_molec(5))) + do i=1,ntyp_molec(5) + do j=1,ntyp_molec(2) + write(iout,*) i,j + read(iionnucl,*) (catnuclprm(k,j,i),k=1,14) + enddo + enddo + write(*,'(3(5x,a6)11(7x,a6))') "w1 ","w2 ","epslj ","pis1 ", & + "sigma0","epsi0 ","chi1 ","chip1 ","sig ","b1 ","b2 ", & + "b3 ","b4 ","chis1 " + do i=1,ntyp_molec(5) + do j=1,ntyp_molec(2) + write(*,'(3(f10.3,x),11(f12.6,x),a3,2a)') (catnuclprm(k,j,i),k=1,14), & + restyp(i,5),"-",restyp(j,2) + enddo + enddo ! read (iion,*) (vcatprm(k),k=1,ncatprotpram) !---------------------------------------------------- allocate(a0thet(-ntyp:ntyp),theta0(-ntyp:ntyp)) @@ -1573,6 +1595,8 @@ itype2loc(-i)=-itype2loc(i) enddo #else + allocate(iloctyp(-nloctyp:nloctyp)) + allocate(itype2loc(-ntyp1:ntyp1)) iloctyp(0)=10 iloctyp(1)=9 iloctyp(2)=20 @@ -1825,7 +1849,7 @@ allocate(ccold(2,2,-nloctyp-1:nloctyp+1)) allocate(ddold(2,2,-nloctyp-1:nloctyp+1)) allocate(eeold(2,2,-nloctyp-1:nloctyp+1)) - + allocate(b(13,-nloctyp-1:nloctyp+1)) if (lprint) & write (iout,*) "Coefficients of the expansion of Eloc(l1,l2)" do i=0,nloctyp-1 @@ -1841,35 +1865,35 @@ b(4,-i)=-b(4,i) b(5,-i)=-b(5,i) endif -c B1(1,i) = b(3) -c B1(2,i) = b(5) -c B1(1,-i) = b(3) -c B1(2,-i) = -b(5) -c b1(1,i)=0.0d0 -c b1(2,i)=0.0d0 -c B1tilde(1,i) = b(3) -c 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 -c B2(1,i) = b(2) -c B2(2,i) = b(4) -c B2(1,-i) =b(2) -c B2(2,-i) =-b(4) -cc B1tilde(1,i) = b(3,i) -cc B1tilde(2,i) =-b(5,i) -C B1tilde(1,-i) =-b(3,i) -C B1tilde(2,-i) =b(5,i) -cc b1tilde(1,i)=0.0d0 -cc b1tilde(2,i)=0.0d0 -cc B2(1,i) = b(2,i) -cc B2(2,i) = b(4,i) -C B2(1,-i) =b(2,i) -C B2(2,-i) =-b(4,i) - -c b2(1,i)=0.0d0 -c b2(2,i)=0.0d0 +!c B1(1,i) = b(3) +!c B1(2,i) = b(5) +!c B1(1,-i) = b(3) +!c B1(2,-i) = -b(5) +!c b1(1,i)=0.0d0 +!c b1(2,i)=0.0d0 +!c B1tilde(1,i) = b(3) +!c! 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 +!c B2(1,i) = b(2) +!c B2(2,i) = b(4) +!c B2(1,-i) =b(2) +!c B2(2,-i) =-b(4) +!cc B1tilde(1,i) = b(3,i) +!cc B1tilde(2,i) =-b(5,i) +!c B1tilde(1,-i) =-b(3,i) +!c B1tilde(2,-i) =b(5,i) +!cc b1tilde(1,i)=0.0d0 +!cc b1tilde(2,i)=0.0d0 +!cc B2(1,i) = b(2,i) +!cc B2(2,i) = b(4,i) +!c B2(1,-i) =b(2,i) +!c B2(2,-i) =-b(4,i) + +!c b2(1,i)=0.0d0 +!c b2(2,i)=0.0d0 CCold(1,1,i)= b(7,i) CCold(2,2,i)=-b(7,i) CCold(2,1,i)= b(9,i) @@ -1878,27 +1902,27 @@ c b2(2,i)=0.0d0 CCold(2,2,-i)=-b(7,i) CCold(2,1,-i)=-b(9,i) CCold(1,2,-i)=-b(9,i) -c CC(1,1,i)=0.0d0 -c CC(2,2,i)=0.0d0 -c CC(2,1,i)=0.0d0 -c CC(1,2,i)=0.0d0 -c Ctilde(1,1,i)= CCold(1,1,i) -c Ctilde(1,2,i)= CCold(1,2,i) -c Ctilde(2,1,i)=-CCold(2,1,i) -c Ctilde(2,2,i)=-CCold(2,2,i) -c CC(1,1,i)=0.0d0 -c CC(2,2,i)=0.0d0 -c CC(2,1,i)=0.0d0 -c CC(1,2,i)=0.0d0 -c Ctilde(1,1,i)= CCold(1,1,i) -c Ctilde(1,2,i)= CCold(1,2,i) -c Ctilde(2,1,i)=-CCold(2,1,i) -c Ctilde(2,2,i)=-CCold(2,2,i) - -c Ctilde(1,1,i)=0.0d0 -c Ctilde(1,2,i)=0.0d0 -c Ctilde(2,1,i)=0.0d0 -c Ctilde(2,2,i)=0.0d0 +!c CC(1,1,i)=0.0d0 +!c CC(2,2,i)=0.0d0 +!c CC(2,1,i)=0.0d0 +!c CC(1,2,i)=0.0d0 +!c Ctilde(1,1,i)= CCold(1,1,i) +!c Ctilde(1,2,i)= CCold(1,2,i) +!c Ctilde(2,1,i)=-CCold(2,1,i) +!c Ctilde(2,2,i)=-CCold(2,2,i) +!c CC(1,1,i)=0.0d0 +!c CC(2,2,i)=0.0d0 +!c CC(2,1,i)=0.0d0 +!c CC(1,2,i)=0.0d0 +!c Ctilde(1,1,i)= CCold(1,1,i) +!c Ctilde(1,2,i)= CCold(1,2,i) +!c Ctilde(2,1,i)=-CCold(2,1,i) +!c Ctilde(2,2,i)=-CCold(2,2,i) + +!c Ctilde(1,1,i)=0.0d0 +!c Ctilde(1,2,i)=0.0d0 +!c Ctilde(2,1,i)=0.0d0 +!c Ctilde(2,2,i)=0.0d0 DDold(1,1,i)= b(6,i) DDold(2,2,i)=-b(6,i) DDold(2,1,i)= b(8,i) @@ -1907,19 +1931,19 @@ c Ctilde(2,2,i)=0.0d0 DDold(2,2,-i)=-b(6,i) DDold(2,1,-i)=-b(8,i) DDold(1,2,-i)=-b(8,i) -c DD(1,1,i)=0.0d0 -c DD(2,2,i)=0.0d0 -c DD(2,1,i)=0.0d0 -c DD(1,2,i)=0.0d0 -c Dtilde(1,1,i)= DD(1,1,i) -c Dtilde(1,2,i)= DD(1,2,i) -c Dtilde(2,1,i)=-DD(2,1,i) -c Dtilde(2,2,i)=-DD(2,2,i) - -c Dtilde(1,1,i)=0.0d0 -c Dtilde(1,2,i)=0.0d0 -c Dtilde(2,1,i)=0.0d0 -c Dtilde(2,2,i)=0.0d0 +!c DD(1,1,i)=0.0d0 +!c DD(2,2,i)=0.0d0 +!c DD(2,1,i)=0.0d0 +!c DD(1,2,i)=0.0d0 +!c Dtilde(1,1,i)= DD(1,1,i) +!c Dtilde(1,2,i)= DD(1,2,i) +!c Dtilde(2,1,i)=-DD(2,1,i) +!c Dtilde(2,2,i)=-DD(2,2,i) + +!c Dtilde(1,1,i)=0.0d0 +!c Dtilde(1,2,i)=0.0d0 +!c Dtilde(2,1,i)=0.0d0 +!c Dtilde(2,2,i)=0.0d0 EEold(1,1,i)= b(10,i)+b(11,i) EEold(2,2,i)=-b(10,i)+b(11,i) EEold(2,1,i)= b(12,i)-b(13,i) @@ -1930,11 +1954,11 @@ c Dtilde(2,2,i)=0.0d0 EEold(1,2,-i)=-b(12,i)-b(13,i) write(iout,*) "TU DOCHODZE" print *,"JESTEM" -c ee(1,1,i)=1.0d0 -c ee(2,2,i)=1.0d0 -c ee(2,1,i)=0.0d0 -c ee(1,2,i)=0.0d0 -c ee(2,1,i)=ee(1,2,i) +!c ee(1,1,i)=1.0d0 +!c ee(2,2,i)=1.0d0 +!c ee(2,1,i)=0.0d0 +!c ee(1,2,i)=0.0d0 +!c ee(2,1,i)=ee(1,2,i) enddo if (lprint) then write (iout,*) @@ -2185,6 +2209,12 @@ c ee(2,1,i)=ee(1,2,i) enddo enddo endif +#ifndef NEWCORR + do i=1,ntyp1 + itype2loc(i)=itortyp(i) + enddo +#endif + ELSE IF (TOR_MODE.eq.1) THEN !C read valence-torsional parameters @@ -2644,23 +2674,45 @@ c ee(2,1,i)=ee(1,2,i) allocate(wstate(4,ntyp,ntyp)) allocate(dhead(2,2,ntyp,ntyp)) allocate(nstate(ntyp,ntyp)) + allocate(debaykap(ntyp,ntyp)) + if (.not.allocated(sigma)) allocate(sigma(0:ntyp1,0:ntyp1)) if (.not.allocated(chi)) allocate(chi(ntyp1,ntyp1)) !(ntyp,ntyp) + do i=1,ntyp do j=1,i ! write (*,*) "Im in ALAB", i, " ", j read(isidep,*) & - eps(i,j),sigma(i,j),chi(i,j),chi(j,i),chipp(i,j),chipp(j,i), & - (alphasur(k,i,j),k=1,4),sigmap1(i,j),sigmap2(i,j), & - chis(i,j),chis(j,i), & - nstate(i,j),(wstate(k,i,j),k=1,4), & - dhead(1,1,i,j),dhead(1,2,i,j),dhead(2,1,i,j),dhead(2,2,i,j),& - dtail(1,i,j),dtail(2,i,j), & - epshead(i,j),sig0head(i,j), & - rborn(i,j),rborn(j,i),(wqdip(k,i,j),k=1,2),wquad(i,j), & - alphapol(i,j),alphapol(j,i), & - (alphiso(k,i,j),k=1,4),sigiso1(i,j),sigiso2(i,j),epsintab(i,j) -! print *,eps(i,j),sigma(i,j),"SIGMAP",i,j,sigmap1(i,j),sigmap2(j,i) + eps(i,j),sigma(i,j),chi(i,j),chi(j,i),chipp(i,j),chipp(j,i), & !6 w tej linii + (alphasur(k,i,j),k=1,4),sigmap1(i,j),sigmap2(i,j), & !6 w tej linii + chis(i,j),chis(j,i), & !2 w tej linii + nstate(i,j),(wstate(k,i,j),k=1,4), & !5 w tej lini - 1 integer pierwszy + dhead(1,1,i,j),dhead(1,2,i,j),dhead(2,1,i,j),dhead(2,2,i,j),& ! 4 w tej linii + dtail(1,i,j),dtail(2,i,j), & ! 2 w tej lini + epshead(i,j),sig0head(i,j), & ! 2 w tej linii + rborn(i,j),rborn(j,i),(wqdip(k,i,j),k=1,2),wquad(i,j), & ! 5 w tej linii + alphapol(i,j),alphapol(j,i), & ! 2 w tej linii + (alphiso(k,i,j),k=1,4),sigiso1(i,j),sigiso2(i,j),epsintab(i,j),debaykap(i,j) !8 w tej linii + IF ((LaTeX).and.(i.gt.24)) then + write (2,'(2a4,1h&,14(f8.2,1h&),23(f8.2,1h&))') restyp(i,1),restyp(j,1), & + eps(i,j),sigma(i,j),chi(i,j),chi(j,i),chipp(i,j),chipp(j,i), & !6 w tej linii + (alphasur(k,i,j),k=1,4),sigmap1(i,j),sigmap2(i,j), & !6 w tej linii + chis(i,j),chis(j,i) !2 w tej linii + endif + print *,eps(i,j),sigma(i,j),"SIGMAP",i,j,sigmap1(i,j),sigmap2(j,i), wqdip(1,i,j) + END DO + END DO + do i=1,ntyp + do j=1,i + IF ((LaTeX).and.(i.gt.24)) then + write (2,'(2a4,1h&,14(f8.2,1h&),23(f8.2,1h&))') restyp(i,1),restyp(j,1), & + dhead(1,1,i,j),dhead(2,1,i,j),& ! 2 w tej linii + dtail(1,i,j),dtail(2,i,j), & ! 2 w tej lini + epshead(i,j),sig0head(i,j), & ! 2 w tej linii + rborn(i,j),rborn(j,i), & ! 3 w tej linii + alphapol(i,j),alphapol(j,i), & ! 2 w tej linii + (alphiso(k,i,j),k=1,4),sigiso1(i,j),sigiso2(i,j),epsintab(i,j),debaykap(i,j) !8 w tej linii + endif END DO END DO DO i = 1, ntyp @@ -2695,10 +2747,13 @@ c ee(2,1,i)=ee(1,2,i) wquad(i,j) = wquad(j,i) epsintab(i,j) = epsintab(j,i) + debaykap(i,j)=debaykap(j,i) ! if (epsintab(i,j).ne.1.0) print *,"WHAT?",i,j,epsintab(i,j) END DO END DO endif + + ! goto 50 !--------------------- GBV potential ----------------------------------- case(5) @@ -3199,6 +3254,158 @@ c ee(2,1,i)=ee(1,2,i) ! v1ss=0.0d0 ! v2ss=0.0d0 ! v3ss=0.0d0 + +! Ions by Aga + + allocate(alphapolcat(ntyp,ntyp),epsheadcat(ntyp,ntyp),sig0headcat(ntyp,ntyp)) + allocate(sigiso1cat(ntyp,ntyp),rborn1cat(ntyp,ntyp),rborn2cat(ntyp,ntyp),sigmap1cat(ntyp,ntyp)) + allocate(sigmap2cat(ntyp,ntyp),sigiso2cat(ntyp,ntyp)) + allocate(chis1cat(ntyp,ntyp),chis2cat(ntyp,ntyp),wquadcat(ntyp,ntyp),chipp1cat(ntyp,ntyp),chipp2cat(ntyp,ntyp)) + allocate(epsintabcat(ntyp,ntyp)) + allocate(dtailcat(2,ntyp,ntyp)) + allocate(alphasurcat(4,ntyp,ntyp),alphisocat(4,ntyp,ntyp)) + allocate(wqdipcat(2,ntyp,ntyp)) + allocate(wstatecat(4,ntyp,ntyp)) + allocate(dheadcat(2,2,ntyp,ntyp)) + allocate(nstatecat(ntyp,ntyp)) + allocate(debaykapcat(ntyp,ntyp)) + + if (.not.allocated(epscat)) allocate (epscat(0:ntyp1,0:ntyp1)) + if (.not.allocated(sigmacat)) allocate(sigmacat(0:ntyp1,0:ntyp1)) +! if (.not.allocated(chicat)) allocate(chicat(ntyp1,ntyp1)) !(ntyp,ntyp) + if (.not.allocated(chi1cat)) allocate(chi1cat(ntyp1,ntyp1)) !(ntyp,ntyp) + if (.not.allocated(chi2cat)) allocate(chi2cat(ntyp1,ntyp1)) !(ntyp,ntyp) + + + if (.not.allocated(ichargecat)) allocate (ichargecat(ntyp_molec(5))) +! i to SC, j to jon, isideocat - nazwa pliku z ktorego czytam parametry + if (oldion.eq.0) then + if (.not.allocated(icharge)) then ! this mean you are oprating in old sc-sc mode + allocate(icharge(1:ntyp1)) + read(iion,*) (icharge(i),i=1,ntyp) + else + read(iion,*) ijunk + endif + + do i=1,ntyp_molec(5) + read(iion,*) msc(i,5),restok(i,5),ichargecat(i) + print *,msc(i,5),restok(i,5) + enddo + ip(5)=0.2 +!DIR$ NOUNROLL + do j=1,ntyp_molec(5) + do i=1,ntyp +! do j=1,ntyp_molec(5) +! write (*,*) "Im in ALAB", i, " ", j + read(iion,*) & + epscat(i,j),sigmacat(i,j), & +! chicat(i,j),chicat(j,i),chippcat(i,j),chippcat(j,i), & + chi1cat(i,j),chi2cat(i,j),chipp1cat(i,j),chipp2cat(i,j), & + + (alphasurcat(k,i,j),k=1,4),sigmap1cat(i,j),sigmap2cat(i,j),& +! chiscat(i,j),chiscat(j,i), & + chis1cat(i,j),chis2cat(i,j), & + + nstatecat(i,j),(wstatecat(k,i,j),k=1,4), & !5 w tej lini - 1 integer pierwszy + dheadcat(1,1,i,j),dheadcat(1,2,i,j),dheadcat(2,1,i,j),dheadcat(2,2,i,j),& + dtailcat(1,i,j),dtailcat(2,i,j), & + epsheadcat(i,j),sig0headcat(i,j), & +!wdipcat = w1 , w2 +! rborncat(i,j),rborncat(j,i),& + rborn1cat(i,j),rborn2cat(i,j),& + (wqdipcat(k,i,j),k=1,2), & + alphapolcat(i,j),alphapolcat(j,i), & + (alphisocat(k,i,j),k=1,4),sigiso1cat(i,j),sigiso2cat(i,j),epsintabcat(i,j),debaykapcat(i,j) +! print *,eps(i,j),sigma(i,j),"SIGMAP",i,j,sigmap1(i,j),sigmap2(j,i) +! if (i.eq.1) then +! write (iout,*) 'i= ', i, ' j= ', j +! write (iout,*) 'epsi0= ', epscat(1,j) +! write (iout,*) 'sigma0= ', sigmacat(1,j) +! write (iout,*) 'chi(i,j)= ', chicat(1,j) +! write (iout,*) 'chi(j,i)= ', chicat(j,1) +! write (iout,*) 'chip(i,j)= ', chippcat(1,j) +! write (iout,*) 'chip(j,i)= ', chippcat(j,1) +! write (iout,*) 'alphasur1= ', alphasurcat(1,1,j) +! write (iout,*) 'alphasur2= ', alphasurcat(2,1,j) +! write (iout,*) 'alphasur3= ', alphasurcat(3,1,j) +! write (iout,*) 'alphasur4= ', alphasurcat(4,1,j) +! write (iout,*) 'sig1= ', sigmap1cat(1,j) +! write (iout,*) 'chis(i,j)= ', chiscat(1,j) +! write (iout,*) 'chis(j,i)= ', chiscat(j,1) +! write (iout,*) 'dhead= ', dheadcat(1,1,1,j) +! write (iout,*) 'a1= ', rborncat(j,1) +! write (iout,*) 'a2= ', rborncat(1,j) +! write (iout,*) 'epsin= ', epsintabcat(1,j), epsintabcat(j,1) +! write (iout,*) 'alphapol1= ', alphapolcat(1,j) +! write (iout,*) 'w1= ', wqdipcat(1,1,j) +! write (iout,*) 'w2= ', wqdipcat(2,1,j) +! endif + +! +! If ((i.eq.1).and.(j.eq.27)) then +! write (iout,*) 'SEP' +! Write (iout,*) 'w1= ', wqdipcat(1,1,27) +! Write (iout,*) 'w2= ', wqdipcat(2,1,27) +! endif + + END DO + END DO + allocate(aa_aq_cat(-ntyp:ntyp,ntyp),bb_aq_cat(-ntyp:ntyp,ntyp)) + do i=1,ntyp + do j=1,ntyp_molec(5) + epsij=epscat(i,j) + rrij=sigmacat(i,j) + rrij=rrij**expon + sigeps=dsign(1.0D0,epsij) + epsij=dabs(epsij) + aa_aq_cat(i,j)=epsij*rrij*rrij + bb_aq_cat(i,j)=-sigeps*epsij*rrij + enddo + enddo + + do i=1,ntyp + do j=1,ntyp_molec(5) + if (i.eq.10) then + write (iout,*) 'i= ', i, ' j= ', j + write (iout,*) 'epsi0= ', epscat(i,j) + write (iout,*) 'sigma0= ', sigmacat(i,j) + write (iout,*) 'chi1= ', chi1cat(i,j) + write (iout,*) 'chi1= ', chi2cat(i,j) + write (iout,*) 'chip1= ', chipp1cat(1,j) + write (iout,*) 'chip2= ', chipp2cat(1,j) + write (iout,*) 'alphasur1= ', alphasurcat(1,1,j) + write (iout,*) 'alphasur2= ', alphasurcat(2,1,j) + write (iout,*) 'alphasur3= ', alphasurcat(3,1,j) + write (iout,*) 'alphasur4= ', alphasurcat(4,1,j) + write (iout,*) 'sig1= ', sigmap1cat(1,j) + write (iout,*) 'sig2= ', sigmap2cat(1,j) + write (iout,*) 'chis1= ', chis1cat(1,j) + write (iout,*) 'chis1= ', chis2cat(1,j) + write (iout,*) 'nstatecat(i,j)= ', nstatecat(1,j) + write (iout,*) 'wstatecat(k,i,j),k=1= ',wstatecat(1,1,j) + write (iout,*) 'dhead= ', dheadcat(1,1,1,j) + write (iout,*) 'dhead2= ', dheadcat(1,2,1,j) + write (iout,*) 'a1= ', rborn1cat(i,j) + write (iout,*) 'a2= ', rborn2cat(i,j) + write (iout,*) 'epsin= ', epsintabcat(1,j), epsintabcat(j,1) + write (iout,*) 'alphapol1= ', alphapolcat(1,j) + write (iout,*) 'alphapol2= ', alphapolcat(j,1) + write (iout,*) 'w1= ', wqdipcat(1,i,j) + write (iout,*) 'w2= ', wqdipcat(2,i,j) + write (iout,*) 'debaykapcat(i,j)= ', debaykapcat(1,j) + endif + + If ((i.eq.1).and.(j.eq.27)) then + write (iout,*) 'SEP' + Write (iout,*) 'w1= ', wqdipcat(1,1,27) + Write (iout,*) 'w2= ', wqdipcat(2,1,27) + endif + + enddo + enddo + + endif + if(me.eq.king) then write (iout,'(/a)') "Disulfide bridge parameters:" @@ -3407,6 +3614,8 @@ c ee(2,1,i)=ee(1,2,i) if (index(card,"FREE ENERGY").gt.0) read(card(35:),*) efree_temp ! Fish out the ATOM cards. ! write(iout,*) 'card',card(1:20) +! print *,"ATU ",card(1:6), CARD(3:6) +! print *,card if (index(card(1:4),'ATOM').gt.0) then read (card(12:16),*) atom ! write (iout,*) "! ",atom," !",ires @@ -3428,9 +3637,9 @@ c ee(2,1,i)=ee(1,2,i) enddo else call sccenter(ires_old,iii,sccor) - endif + endif !unres_pdb iii=0 - endif + endif !ind_pdb ! Start new residue. if (res.eq.'Cl-' .or. res.eq.'Na+') then ires=ires_old @@ -3442,7 +3651,7 @@ c ee(2,1,i)=ee(1,2,i) ishift=ishift-1 itype(1,1)=ntyp1 nres_molec(molecule)=nres_molec(molecule)+1 - endif + endif ! Gly ires=ires-ishift+ishift1 ires_old=ires ! write (iout,*) "ishift",ishift," ires",ires,& @@ -3461,7 +3670,7 @@ c ee(2,1,i)=ee(1,2,i) ishift=ishift-(ires-ishift+ishift1-ires_old-1) ires=ires-ishift+ishift1 ires_old=ires - endif + endif ! Na Cl ! print *,'atom',ires,atom if (res.eq.'ACE' .or. res.eq.'NHE') then itype(ires,1)=10 @@ -3485,11 +3694,11 @@ c ee(2,1,i)=ee(1,2,i) ! print *,"ires=",ires,istype(ires) ! endif - endif - endif + endif ! atom.eq.CA + endif !ACE else ires=ires-ishift+ishift1 - endif + endif !ires_old ! write (iout,*) "ires_old",ires_old," ires",ires if (card(27:27).eq."A" .or. card(27:27).eq."B") then ! ishift1=ishift1+1 @@ -3574,6 +3783,7 @@ c ee(2,1,i)=ee(1,2,i) read(card(31:54),'(3f8.3)') (sccor(j,iii),j=1,3) endif endif +! print *,"IONS",ions,card(1:6) else if ((ions).and.(card(1:6).eq.'HETATM')) then if (firstion.eq.0) then firstion=1 @@ -3583,10 +3793,10 @@ c ee(2,1,i)=ee(1,2,i) enddo else call sccenter(ires,iii,sccor) - endif - endif + endif ! unres_pdb + endif !firstion read (card(12:16),*) atom - print *,"HETATOM", atom +! print *,"HETATOM", atom read (card(18:20),'(a3)') res if ((atom(1:2).eq.'NA').or.(atom(1:2).eq.'CL').or.& (atom(1:2).eq.'CA').or.(atom(1:2).eq.'MG') & @@ -3600,7 +3810,7 @@ c ee(2,1,i)=ee(1,2,i) res=res(2:3)//' ' itype(ires,molecule)=rescode(ires,res,0,molecule) read(card(31:54),'(3f8.3)') (c(j,ires),j=1,3) - endif + endif! NA endif !atom enddo 10 write (iout,'(a,i5)') ' Number of residues found: ',ires @@ -3674,7 +3884,8 @@ c ee(2,1,i)=ee(1,2,i) ! e2(3)=0.0d0 ! endif do j=1,3 - c(j,i)=c(j,i+1)-1.9d0*e2(j) +! c(j,i)=c(j,i+1)-1.9d0*e2(j) + c(j,i)=c(j,i-1)+1.9d0*(-e1(j)+e2(j))/sqrt(2.0d0) enddo ! else !unres_pdb do j=1,3 @@ -3789,7 +4000,8 @@ c ee(2,1,i)=ee(1,2,i) e2(3)=0.0d0 endif do j=1,3 - c(j,1)=c(j,2)-1.9d0*e2(j) +! c(j,1)=c(j,2)-1.9d0*e2(j) + c(j,1)=c(j,2)+1.9d0*(e1(j)-e2(j))/sqrt(2.0d0) enddo else do j=1,3 @@ -3882,6 +4094,18 @@ c ee(2,1,i)=ee(1,2,i) istype(i)=istype_temp(i) enddo enddo + if ((itype(1,1).eq.ntyp1).and.itype(2,5).ne.0) then +! I have only ions now dummy atoms in the system + molnum(1)=5 + itype(1,5)=itype(2,5) + itype(1,1)=0 + do i=2,nres + itype(i,5)=itype(i+1,5) + enddo + itype(nres,5)=0 + nres=nres-1 + nres_molec(1)=nres_molec(1)-1 + endif ! if (itype(1,1).eq.ntyp1) then ! nsup=nsup-1 ! nstart_sup=2 @@ -4206,6 +4430,7 @@ c ee(2,1,i)=ee(1,2,i) with_theta_constr = index(controlcard,"WITH_THETA_CONSTR").gt.0 protein=index(controlcard,"PROTEIN").gt.0 ions=index(controlcard,"IONS").gt.0 + call readi(controlcard,'OLDION',oldion,1) nucleic=index(controlcard,"NUCLEIC").gt.0 write (iout,*) "with_theta_constr ",with_theta_constr AFMlog=(index(controlcard,'AFM')) @@ -4223,6 +4448,8 @@ c ee(2,1,i)=ee(1,2,i) ! elemode = 0 is orignal UNRES electrostatics ! elemode = 1 is "Momo" potentials in progress ! elemode = 2 is in development EVALD + + write (iout,*) TUBEmode,"TUBEMODE" if (TUBEmode.gt.0) then call reada(controlcard,"XTUBE",tubecenter(1),0.0d0) @@ -4237,7 +4464,7 @@ c ee(2,1,i)=ee(1,2,i) endif ! CUTOFFF ON ELECTROSTATICS - call reada(controlcard,"R_CUT_ELE",r_cut_ele,15.0d0) + call reada(controlcard,"R_CUT_ELE",r_cut_ele,25.0d0) call reada(controlcard,"LAMBDA_ELE",rlamb_ele,0.3d0) write(iout,*) "R_CUT_ELE=",r_cut_ele ! Lipidic parameters @@ -4378,6 +4605,7 @@ c ee(2,1,i)=ee(1,2,i) ! enddo buff_shield=1.0d0 endif + itime_mat=0 return end subroutine read_control !----------------------------------------------------------------------------- @@ -5128,6 +5356,8 @@ c ee(2,1,i)=ee(1,2,i) open (itube,file=tubename,status='old') call getenv_loc('IONPAR',ionname) open (iion,file=ionname,status='old') + call getenv_loc('IONPAR_NUCL',ionnuclname) + open (iionnucl,file=ionnuclname,status='old') #else open(1,file=pref_orig(:ilen(pref_orig))//'.inp',status='old',& readonly) @@ -5187,6 +5417,8 @@ c ee(2,1,i)=ee(1,2,i) open (itube,file=tubename,status='old',action='read') call getenv_loc('IONPAR',ionname) open (iion,file=ionname,status='old',action='read') + call getenv_loc('IONPAR_NUCL',ionnuclname) + open (iionnucl,file=ionnuclname,status='old',action='read') #ifndef CRYST_SC call getenv_loc('ROTPARPDB',rotname_pdb)