include 'COMMON.SBRIDGE'
include 'COMMON.MD'
include 'COMMON.SETUP'
+ include 'COMMON.CONTROL'
+ include 'COMMON.SHIELD'
character*1 t1,t2,t3
character*1 onelett(4) /"G","A","P","D"/
character*1 toronelet(-2:2) /"p","a","G","A","P"/
vblinv2=vblinv*vblinv
c
c Read the virtual-bond parameters, masses, and moments of inertia
-c and Stokes' radii of the peptide group and side chains
+c and Stokes radii of the peptide group and side chains
c
#ifdef CRYST_BOND
read (ibond,*) vbldp0,vbldpdum,akp,mp,ip,pstok
#else
read (ibond,*) junk,vbldp0,vbldpdum,akp,rjunk,mp,ip,pstok
do i=1,ntyp
+ print *,i
read (ibond,*) nbondterm(i),(vbldsc0(j,i),aksc(j,i),abond0(j,i),
& j=1,nbondterm(i)),msc(i),isc(i),restok(i)
dsc(i) = vbldsc0(1,i)
& vbldsc0(j,i),aksc(j,i),abond0(j,i)
enddo
enddo
+#ifdef AIX
+ call flush_(iout)
+#else
+ call flush(iout)
+#endif
endif
+C reading lipid parameters
+ read(iliptranpar,*,end=120,err=120) pepliptran
+ do i=1,ntyp
+ read(iliptranpar,*,end=120,err=120) liptranene(i)
+ enddo
+ close(iliptranpar)
+ if (lprint) then
+ write (iout,*) "Lipid transfer parameters"
+ write (iout,'(a5,f10.5)') "pept",pepliptran
+ do i=1,ntyp
+ write (iout,'(a5,f10.5)') restyp(i),liptranene(i)
+ enddo
+#ifdef AIX
+ call flush_(iout)
+#else
+ call flush(iout)
+#endif
+ endif
#ifdef CRYST_THETA
C
C Read the parameters of the probability distribution/energy expression
write (iout,'(a3,i4,2x,5(1pe14.5))') restyp(i),i,theta0(i),
& sig0(i),(gthet(j,i),j=1,3)
enddo
+ call flush(iout)
else
write (iout,'(a)')
& 'Parameters of the virtual-bond valence angles:'
write (iout,'(a3,1h&,2x,5(f8.3,1h&))') restyp(i),theta0(i),
& 100*sig0(i),gthet(1,i)*0.1D0,gthet(2,i),gthet(3,i)*10.0D0
enddo
+#ifdef AIX
+ call flush_(iout)
+#else
+ call flush(iout)
+#endif
endif
endif
#else
enddo
enddo
enddo
- call flush(iout)
+#ifdef AIX
+ call flush_(iout)
+#else
+ call flush(iout)
+#endif
endif
- write (2,*) "Start reading THETA_PDB",ithep_pdb
+ write (iout,*) "Start reading THETA_PDB",ithep_pdb
do i=1,ntyp
c write (2,*) 'i=',i
read (ithep_pdb,*,err=111,end=111)
gthet(j,i)=gthet(j,-i)
enddo
enddo
- write (2,*) "End reading THETA_PDB"
+ write (iout,*) "End reading THETA_PDB"
close (ithep_pdb)
#endif
close(ithep)
endif
endif
enddo
+#ifdef AIX
+ call flush_(iout)
+#else
+ call flush(iout)
+#endif
endif
#else
C
C Read the parameters of the probability distribution/energy expression
C of the side chains.
C
- write (2,*) "Start reading ROTAM_PDB"
+ write (iout,*) "Start reading ROTAM_PDB"
do i=1,ntyp
read (irotam_pdb,'(3x,i3,f8.3)',end=112,err=112) nlob(i),dsc(i)
if (i.eq.10) then
endif
enddo
close (irotam_pdb)
- write (2,*) "End reading ROTAM_PDB"
+c write (iout,*) "End reading ROTAM_PDB"
+c#ifdef AIX
+c call flush_(iout)
+c#else
+c call flush(iout)
+c#endif
#endif
close(irotam)
write (iout,'(6x,6f10.5)') (v2(k,i,j),k=1,nterm_old)
enddo
enddo
+#ifdef AIX
+ call flush_(iout)
+#else
+ call flush(iout)
+#endif
endif
#else
C
enddo
enddo
close (itorp)
+c write (iout,*) "End reading torsional parameters"
+c#ifdef AIX
+c call flush_(iout)
+c#else
+c call flush(iout)
+c#endif
if (lprint) then
write (iout,'(/a/)') 'Torsional constants:'
+ do iblock=1,2
+ write (iout,*) "IBLOCK",iblock
do i=1,ntortyp
do j=1,ntortyp
write (iout,*) 'ityp',i,' jtyp',j
enddo
enddo
enddo
+ enddo
+#ifdef AIX
+ call flush_(iout)
+#else
+ call flush(iout)
+#endif
endif
C
enddo
enddo
enddo
+#ifdef AIX
+ call flush_(iout)
+#else
+ call flush(iout)
+#endif
endif
#endif
C Read of Side-chain backbone correlation parameters
CCC
C
read (isccor,*,end=119,err=119) nsccortyp
+c write (iout,*) "Reading sccor parameters",nsccortyp
+c#ifdef AIX
+c call flush_(iout)
+c#else
+c call flush(iout)
+c#endif
#ifdef SCCORPDB
read (isccor,*,end=119,err=119) (isccortyp(i),i=1,ntyp)
do i=-ntyp,-1
v0ijsccor=v0ijsccor+vlor1sccor(k,i,j)/
&(1+vlor3sccor(k,i,j)**2)
enddo
- v0sccor(i,j,iblock)=v0ijsccor
+ v0sccor(l,i,j)=v0ijsccor
enddo
enddo
enddo
close (isccor)
#endif
+c write (iout,*) "sccor parameters read"
+c#ifdef AIX
+c call flush_(iout)
+c#else
+c call flush(iout)
+c#endif
if (lprint) then
- write (iout,'(/a/)') 'Torsional constants:'
+ write (iout,'(/a/)') 'SC-torsional constants:'
+ do l=1,maxinter
+ write (iout,*) "Torsional type",l
do i=1,nsccortyp
do j=1,nsccortyp
write (iout,*) 'ityp',i,' jtyp',j
enddo
enddo
enddo
+ enddo
+#ifdef AIX
+ call flush_(iout)
+#else
+ call flush(iout)
+#endif
endif
C
c B2(2,i) = b(4)
c B2(1,-i) =b(2)
c B2(2,-i) =-b(4)
+ B1tilde(1,i) = b(3,i)
+ B1tilde(2,i) =-b(5,i)
+C B1tilde(1,-i) =-b(3,i)
+C B1tilde(2,-i) =b(5,i)
+ b1tilde(1,i)=0.0d0
+ b1tilde(2,i)=0.0d0
+ B2(1,i) = b(2,i)
+ 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
EEold(2,1,-i)=-b(12,i)+b(13,i)
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
write (iout,'(2f10.5)') EEold(j,1,i),EEold(j,2,i)
enddo
enddo
+#ifdef AIX
+ call flush_(iout)
+#else
+ call flush(iout)
+#endif
endif
c lprint=.false.
c lprint=.false.
enddo
enddo
+#ifdef AIX
+ if (lprint) call flush_(iout)
+#else
+ if (lprint) call flush(iout)
+#endif
C
C Read side-chain interaction parameters.
C
& ', exponents are ',expon,2*expon
goto (10,20,30,30,40) ipot
C----------------------- LJ potential ---------------------------------
- 10 read (isidep,*,end=116,err=116)((eps(i,j),j=i,ntyp),i=1,ntyp),
+ 10 read (isidep,*,end=117,err=117)((eps(i,j),j=i,ntyp),i=1,ntyp),
& (sigma0(i),i=1,ntyp)
if (lprint) then
write (iout,'(/a/)') 'Parameters of the LJ potential:'
endif
goto 50
C----------------------- LJK potential --------------------------------
- 20 read (isidep,*,end=116,err=116)((eps(i,j),j=i,ntyp),i=1,ntyp),
+ 20 read (isidep,*,end=117,err=117)((eps(i,j),j=i,ntyp),i=1,ntyp),
& (sigma0(i),i=1,ntyp),(rr0(i),i=1,ntyp)
if (lprint) then
write (iout,'(/a/)') 'Parameters of the LJK potential:'
goto 50
C---------------------- GB or BP potential -----------------------------
30 do i=1,ntyp
- read (isidep,*,end=116,err=116)(eps(i,j),j=i,ntyp)
+ read (isidep,*,end=117,err=117)(eps(i,j),j=i,ntyp)
enddo
read (isidep,*,end=116,err=116)(sigma0(i),i=1,ntyp)
read (isidep,*,end=116,err=116)(sigii(i),i=1,ntyp)
read (isidep,*,end=116,err=116)(chip(i),i=1,ntyp)
read (isidep,*,end=116,err=116)(alp(i),i=1,ntyp)
+C now we start reading lipid
+ do i=1,ntyp
+ read (isidep,*,end=1161,err=1161)(epslip(i,j),j=i,ntyp)
+C print *,"WARNING!!"
+C do j=1,ntyp
+C epslip(i,j)=epslip(i,j)+0.05d0
+C enddo
+ enddo
+ write(iout,*) epslip(1,1),"OK?"
C For the GB potential convert sigma'**2 into chi'
if (ipot.eq.4) then
do i=1,ntyp
endif
goto 50
C--------------------- GBV potential -----------------------------------
- 40 read (isidep,*,end=116,err=116)((eps(i,j),j=i,ntyp),i=1,ntyp),
+ 40 read (isidep,*,end=117,err=117)((eps(i,j),j=i,ntyp),i=1,ntyp),
& (sigma0(i),i=1,ntyp),(rr0(i),i=1,ntyp),(sigii(i),i=1,ntyp),
& (chip(i),i=1,ntyp),(alp(i),i=1,ntyp)
if (lprint) then
do i=2,ntyp
do j=1,i-1
eps(i,j)=eps(j,i)
+ epslip(i,j)=epslip(j,i)
enddo
enddo
do i=1,ntyp
epsij=eps(i,j)
sigeps=dsign(1.0D0,epsij)
epsij=dabs(epsij)
- aa(i,j)=epsij*rrij*rrij
- bb(i,j)=-sigeps*epsij*rrij
- aa(j,i)=aa(i,j)
- bb(j,i)=bb(i,j)
+ aa_aq(i,j)=epsij*rrij*rrij
+ bb_aq(i,j)=-sigeps*epsij*rrij
+ aa_aq(j,i)=aa_aq(i,j)
+ bb_aq(j,i)=bb_aq(i,j)
+ epsijlip=epslip(i,j)
+ sigeps=dsign(1.0D0,epsijlip)
+ epsijlip=dabs(epsijlip)
+ aa_lip(i,j)=epsijlip*rrij*rrij
+ bb_lip(i,j)=-sigeps*epsijlip*rrij
+ aa_lip(j,i)=aa_lip(i,j)
+ bb_lip(j,i)=bb_lip(i,j)
if (ipot.gt.2) then
sigt1sq=sigma0(i)**2
sigt2sq=sigma0(j)**2
endif
if (lprint) then
write (iout,'(2(a3,2x),3(1pe10.3),5(0pf8.3))')
- & restyp(i),restyp(j),aa(i,j),bb(i,j),augm(i,j),
+ & restyp(i),restyp(j),aa_aq(i,j),bb_aq(i,j),augm(i,j),
& sigma(i,j),r0(i,j),chi(i,j),chi(j,i)
endif
enddo
enddo
+#ifdef AIX
+ if (lprint) call flush_(iout)
+#else
+ if (lprint) call flush(iout)
+#endif
#ifdef OLDSCP
C
C Define the SC-p interaction constants (hard-coded; old style)
write (iout,'(4f8.3,4e12.4)') eps_scp(i,1),rscp(i,1),
& eps_scp(i,2),rscp(i,2),aad(i,1),bad(i,1),aad(i,2),bad(i,2)
enddo
+#ifdef AIX
+ call flush_(iout)
+#else
+ call flush(iout)
+#endif
endif
c lprint=.false.
#endif
C
C Define the constants of the disulfide bridge
C
- ebr=-5.50D0
+C ebr=-12.00D0
c
c Old arbitrary potential - commented out.
c
c energy surface of diethyl disulfide.
c A. Liwo and U. Kozlowska, 11/24/03
c
- D0CM = 3.78d0
- AKCM = 15.1d0
- AKTH = 11.0d0
- AKCT = 12.0d0
- V1SS =-1.08d0
- V2SS = 7.61d0
- V3SS = 13.7d0
+C D0CM = 3.78d0
+C AKCM = 15.1d0
+C AKTH = 11.0d0
+C AKCT = 12.0d0
+C V1SS =-1.08d0
+C V2SS = 7.61d0
+C V3SS = 13.7d0
c akcm=0.0d0
c akth=0.0d0
c akct=0.0d0
c v2ss=0.0d0
c v3ss=0.0d0
- if(me.eq.king) then
- write (iout,'(/a)') "Disulfide bridge parameters:"
- write (iout,'(a,f10.2)') 'S-S bridge energy: ',ebr
- write (iout,'(2(a,f10.2))') 'd0cm:',d0cm,' akcm:',akcm
- write (iout,'(2(a,f10.2))') 'akth:',akth,' akct:',akct
- write (iout,'(3(a,f10.2))') 'v1ss:',v1ss,' v2ss:',v2ss,
- & ' v3ss:',v3ss
- endif
+C if(me.eq.king) then
+C write (iout,'(/a)') "Disulfide bridge parameters:"
+C write (iout,'(a,f10.2)') 'S-S bridge energy: ',ebr
+C write (iout,'(2(a,f10.2))') 'd0cm:',d0cm,' akcm:',akcm
+C write (iout,'(2(a,f10.2))') 'akth:',akth,' akct:',akct
+C write (iout,'(3(a,f10.2))') 'v1ss:',v1ss,' v2ss:',v2ss,
+C & ' v3ss:',v3ss
+C endif
+C set the variables used for shielding effect
+C write (iout,*) "SHIELD MODE",shield_mode
+C if (shield_mode.gt.0) then
+C VSolvSphere the volume of solving sphere
+C print *,pi,"pi"
+C rpp(1,1) is the energy r0 for peptide group contact and will be used for it
+C there will be no distinction between proline peptide group and normal peptide
+C group in case of shielding parameters
+C VSolvSphere=4.0/3.0*pi*rpp(1,1)**3
+C VSolvSphere_div=VSolvSphere-4.0/3.0*pi*(rpp(1,1)/2.0)**3
+C write (iout,*) VSolvSphere,VSolvSphere_div
+C long axis of side chain
+C do i=1,ntyp
+C long_r_sidechain(i)=vbldsc0(1,i)
+C short_r_sidechain(i)=sigma0(i)
+C enddo
+C lets set the buffor value
+C buff_shield=1.0d0
+C endif
return
111 write (iout,*) "Error reading bending energy parameters."
+#ifdef AIX
+ call flush_(iout)
+#else
+ call flush(iout)
+#endif
goto 999
112 write (iout,*) "Error reading rotamer energy parameters."
+#ifdef AIX
+ call flush_(iout)
+#else
+ call flush(iout)
+#endif
goto 999
113 write (iout,*) "Error reading torsional energy parameters."
+#ifdef AIX
+ call flush_(iout)
+#else
+ call flush(iout)
+#endif
goto 999
114 write (iout,*) "Error reading double torsional energy parameters."
+#ifdef AIX
+ call flush_(iout)
+#else
+ call flush(iout)
+#endif
goto 999
115 write (iout,*)
& "Error reading cumulant (multibody energy) parameters."
+#ifdef AIX
+ call flush_(iout)
+#else
+ call flush(iout)
+#endif
goto 999
116 write (iout,*) "Error reading electrostatic energy parameters."
+#ifdef AIX
+ call flush_(iout)
+#else
+ call flush(iout)
+#endif
+ goto 999
+ 1161 write (iout,*) "Error reading electrostatic energy parameters.Lip"
+#ifdef AIX
+ call flush_(iout)
+#else
+ call flush(iout)
+#endif
goto 999
117 write (iout,*) "Error reading side chain interaction parameters."
+#ifdef AIX
+ call flush_(iout)
+#else
+ call flush(iout)
+#endif
goto 999
118 write (iout,*) "Error reading SCp interaction parameters."
+#ifdef AIX
+ call flush_(iout)
+#else
+ call flush(iout)
+#endif
goto 999
119 write (iout,*) "Error reading SCCOR parameters"
+#ifdef AIX
+ call flush_(iout)
+#else
+ call flush(iout)
+#endif
+ goto 999
+ 120 write (iout,*) "Error reading lipid parameters"
+#ifdef AIX
+ call flush_(iout)
+#else
+ call flush(iout)
+#endif
999 continue
#ifdef MPI
call MPI_Finalize(Ierror)
#else
call getenv(var,val)
#endif
-
+C set the variables used for shielding effect
+C if (shield_mode.gt.0) then
+C VSolvSphere the volume of solving sphere
+C print *,pi,"pi"
+C rpp(1,1) is the energy r0 for peptide group contact and will be used for it
+C there will be no distinction between proline peptide group and normal peptide
+C group in case of shielding parameters
+C VSolvSphere=4.0/3.0*pi*rpp(1,1)**3
+C VSolvSphere_div=VSolvSphere-4.0/3.0*pi*(rpp(1,1)/2.0)**3
+C long axis of side chain
+C do i=1,ntyp
+C long_r_sidechain(i)=vbldsc0(1,i)
+C short_r_sidechain(i)=sigma0(i)
+C enddo
+C lets set the buffor value
+C buff_shield=1.0d0
+C endif
return
end
projects / unres.git / blobdiff