3 C Read the parameters of the probability distributions of the virtual-bond
4 C valence angles and the side chains and energy parameters.
6 C Important! Energy-term weights ARE NOT read here; they are read from the
7 C main input file instead, because NO defaults have yet been set for these
10 implicit real*8 (a-h,o-z)
16 include 'COMMON.IOUNITS'
17 include 'COMMON.CHAIN'
18 include 'COMMON.INTERACT'
20 include 'COMMON.LOCAL'
21 include 'COMMON.TORSION'
22 include 'COMMON.SCCOR'
23 include 'COMMON.SCROT'
24 include 'COMMON.FFIELD'
25 include 'COMMON.NAMES'
26 include 'COMMON.SBRIDGE'
28 include 'COMMON.SETUP'
30 character*1 onelett(4) /"G","A","P","D"/
32 dimension blower(3,3,maxlob)
34 character*3 lancuch,ucase
36 C For printing parameters after they are read set the following in the UNRES
39 C setenv PRINT_PARM YES
41 C To print parameters in LaTeX format rather than as ASCII tables:
45 call getenv_loc("PRINT_PARM",lancuch)
46 lprint = (ucase(lancuch).eq."YES" .or. ucase(lancuch).eq."Y")
47 call getenv_loc("LATEX",lancuch)
48 LaTeX = (ucase(lancuch).eq."YES" .or. ucase(lancuch).eq."Y")
50 dwa16=2.0d0**(1.0d0/6.0d0)
52 C Assign virtual-bond length
57 c Read the virtual-bond parameters, masses, and moments of inertia
58 c and Stokes' radii of the peptide group and side chains
61 read (ibond,*) vbldp0,akp,mp,ip,pstok
64 read (ibond,*) vbldsc0(1,i),aksc(1,i),msc(i),isc(i),restok(i)
69 dsc_inv(i)=1.0D0/dsc(i)
73 read (ibond,*) junk,vbldp0,akp,rjunk,mp,ip,pstok
75 read (ibond,*) nbondterm(i),(vbldsc0(j,i),aksc(j,i),abond0(j,i),
76 & j=1,nbondterm(i)),msc(i),isc(i),restok(i)
81 dsc_inv(i)=1.0D0/dsc(i)
86 write(iout,'(/a/)')"Dynamic constants of the interaction sites:"
87 write (iout,'(a10,a3,6a10)') 'Type','N','VBL','K','A0','mass',
89 write(iout,'(a10,i3,6f10.5)') "p",1,vbldp0,akp,0.0d0,mp,ip,pstok
91 write (iout,'(a10,i3,6f10.5)') restyp(i),nbondterm(i),
92 & vbldsc0(1,i),aksc(1,i),abond0(1,i),msc(i),isc(i),restok(i)
94 write (iout,'(13x,3f10.5)')
95 & vbldsc0(j,i),aksc(j,i),abond0(j,i)
101 C Read the parameters of the probability distribution/energy expression
102 C of the virtual-bond valence angles theta
105 read (ithep,*,err=111,end=111) a0thet(i),(athet(j,i),j=1,2),
107 read (ithep,*,err=111,end=111) (polthet(j,i),j=0,3)
108 read (ithep,*,err=111,end=111) (gthet(j,i),j=1,3)
109 read (ithep,*,err=111,end=111) theta0(i),sig0(i),sigc0(i)
116 & 'Parameters of the virtual-bond valence angles:'
117 write (iout,'(/a/9x,5a/79(1h-))') 'Fourier coefficients:',
118 & ' ATHETA0 ',' A1 ',' A2 ',
121 write(iout,'(a3,i4,2x,5(1pe14.5))') restyp(i),i,
122 & a0thet(i),(athet(j,i),j=1,2),(bthet(j,i),j=1,2)
124 write (iout,'(/a/9x,5a/79(1h-))')
125 & 'Parameters of the expression for sigma(theta_c):',
126 & ' ALPH0 ',' ALPH1 ',' ALPH2 ',
127 & ' ALPH3 ',' SIGMA0C '
129 write (iout,'(a3,i4,2x,5(1pe14.5))') restyp(i),i,
130 & (polthet(j,i),j=0,3),sigc0(i)
132 write (iout,'(/a/9x,5a/79(1h-))')
133 & 'Parameters of the second gaussian:',
134 & ' THETA0 ',' SIGMA0 ',' G1 ',
137 write (iout,'(a3,i4,2x,5(1pe14.5))') restyp(i),i,theta0(i),
138 & sig0(i),(gthet(j,i),j=1,3)
142 & 'Parameters of the virtual-bond valence angles:'
143 write (iout,'(/a/9x,5a/79(1h-))')
144 & 'Coefficients of expansion',
145 & ' theta0 ',' a1*10^2 ',' a2*10^2 ',
146 & ' b1*10^1 ',' b2*10^1 '
148 write(iout,'(a3,1h&,2x,5(f8.3,1h&))') restyp(i),
149 & a0thet(i),(100*athet(j,i),j=1,2),(10*bthet(j,i),j=1,2)
151 write (iout,'(/a/9x,5a/79(1h-))')
152 & 'Parameters of the expression for sigma(theta_c):',
153 & ' alpha0 ',' alph1 ',' alph2 ',
154 & ' alhp3 ',' sigma0c '
156 write (iout,'(a3,1h&,2x,5(1pe12.3,1h&))') restyp(i),
157 & (polthet(j,i),j=0,3),sigc0(i)
159 write (iout,'(/a/9x,5a/79(1h-))')
160 & 'Parameters of the second gaussian:',
161 & ' theta0 ',' sigma0*10^2 ',' G1*10^-1',
164 write (iout,'(a3,1h&,2x,5(f8.3,1h&))') restyp(i),theta0(i),
165 & 100*sig0(i),gthet(1,i)*0.1D0,gthet(2,i),gthet(3,i)*10.0D0
171 C Read the parameters of Utheta determined from ab initio surfaces
172 C Kozlowska et al., J. Phys.: Condens. Matter 19 (2007) 285203
174 read (ithep,*,err=111,end=111) nthetyp,ntheterm,ntheterm2,
175 & ntheterm3,nsingle,ndouble
176 nntheterm=max0(ntheterm,ntheterm2,ntheterm3)
177 read (ithep,*,err=111,end=111) (ithetyp(i),i=1,ntyp1)
183 aathet(l,i,j,k)=0.0d0
187 bbthet(m,l,i,j,k)=0.0d0
188 ccthet(m,l,i,j,k)=0.0d0
189 ddthet(m,l,i,j,k)=0.0d0
190 eethet(m,l,i,j,k)=0.0d0
196 ffthet(mm,m,l,i,j,k)=0.0d0
197 ggthet(mm,m,l,i,j,k)=0.0d0
207 read (ithep,'(3a)',end=111,err=111) res1,res2,res3
208 read (ithep,*,end=111,err=111) aa0thet(i,j,k)
209 read (ithep,*,end=111,err=111)(aathet(l,i,j,k),l=1,ntheterm)
210 read (ithep,*,end=111,err=111)
211 & ((bbthet(lll,ll,i,j,k),lll=1,nsingle),
212 & (ccthet(lll,ll,i,j,k),lll=1,nsingle),
213 & (ddthet(lll,ll,i,j,k),lll=1,nsingle),
214 & (eethet(lll,ll,i,j,k),lll=1,nsingle),ll=1,ntheterm2)
215 read (ithep,*,end=111,err=111)
216 & (((ffthet(llll,lll,ll,i,j,k),ffthet(lll,llll,ll,i,j,k),
217 & ggthet(llll,lll,ll,i,j,k),ggthet(lll,llll,ll,i,j,k),
218 & llll=1,lll-1),lll=2,ndouble),ll=1,ntheterm3)
223 C For dummy ends assign glycine-type coefficients of theta-only terms; the
224 C coefficients of theta-and-gamma-dependent terms are zero.
229 aathet(l,i,j,nthetyp+1)=aathet(l,i,j,1)
230 aathet(l,nthetyp+1,i,j)=aathet(l,1,i,j)
232 aa0thet(i,j,nthetyp+1)=aa0thet(i,j,1)
233 aa0thet(nthetyp+1,i,j)=aa0thet(1,i,j)
236 aathet(l,nthetyp+1,i,nthetyp+1)=aathet(l,1,i,1)
238 aa0thet(nthetyp+1,i,nthetyp+1)=aa0thet(1,i,1)
241 C Control printout of the coefficients of virtual-bond-angle potentials
244 write (iout,'(//a)') 'Parameter of virtual-bond-angle potential'
248 write (iout,'(//4a)')
249 & 'Type ',onelett(i),onelett(j),onelett(k)
250 write (iout,'(//a,10x,a)') " l","a[l]"
251 write (iout,'(i2,1pe15.5)') 0,aa0thet(i,j,k)
252 write (iout,'(i2,1pe15.5)')
253 & (l,aathet(l,i,j,k),l=1,ntheterm)
255 write (iout,'(//2h m,4(9x,a,3h[m,,i1,1h]))')
256 & "b",l,"c",l,"d",l,"e",l
258 write (iout,'(i2,4(1pe15.5))') m,
259 & bbthet(m,l,i,j,k),ccthet(m,l,i,j,k),
260 & ddthet(m,l,i,j,k),eethet(m,l,i,j,k)
264 write (iout,'(//3hm,n,4(6x,a,5h[m,n,,i1,1h]))')
265 & "f+",l,"f-",l,"g+",l,"g-",l
268 write (iout,'(i1,1x,i1,4(1pe15.5))') n,m,
269 & ffthet(n,m,l,i,j,k),ffthet(m,n,l,i,j,k),
270 & ggthet(n,m,l,i,j,k),ggthet(m,n,l,i,j,k)
279 write (2,*) "Start reading THETA_PDB"
281 read (ithep_pdb,*,err=111,end=111) a0thet(i),(athet(j,i),j=1,2),
283 read (ithep_pdb,*,err=111,end=111) (polthet(j,i),j=0,3)
284 read (ithep_pdb,*,err=111,end=111) (gthet(j,i),j=1,3)
285 read (ithep_pdb,*,err=111,end=111) theta0(i),sig0(i),sigc0(i)
288 write (2,*) "End reading THETA_PDB"
294 C Read the parameters of the probability distribution/energy expression
295 C of the side chains.
298 read (irotam,'(3x,i3,f8.3)',end=112,err=112) nlob(i),dsc(i)
302 dsc_inv(i)=1.0D0/dsc(i)
313 read(irotam,*,end=112,err=112)(censc(k,1,i),k=1,3),
314 & ((blower(k,l,1),l=1,k),k=1,3)
316 read (irotam,*,end=112,err=112) bsc(j,i)
317 read (irotam,*,end=112,err=112) (censc(k,j,i),k=1,3),
318 & ((blower(k,l,j),l=1,k),k=1,3)
325 akl=akl+blower(k,m,j)*blower(l,m,j)
336 write (iout,'(/a)') 'Parameters of side-chain local geometry'
341 write (iout,'(/3a,i2,a,f8.3)') 'Residue type: ',restyp(i),
342 & ' # of gaussian lobes:',nlobi,' dsc:',dsc(i)
343 write (iout,'(1h&,a,3(2h&&,f8.3,2h&&))')
344 & 'log h',(bsc(j,i),j=1,nlobi)
345 write (iout,'(1h&,a,3(1h&,f8.3,1h&,f8.3,1h&,f8.3,1h&))')
346 & 'x',((censc(k,j,i),k=1,3),j=1,nlobi)
348 write (iout,'(2h& ,5(2x,1h&,3(f7.3,1h&)))')
349 & ((gaussc(k,l,j,i),l=1,3),j=1,nlobi)
352 write (iout,'(/a,8x,i1,4(25x,i1))') 'Lobe:',(j,j=1,nlobi)
353 write (iout,'(a,f10.4,4(16x,f10.4))')
354 & 'Center ',(bsc(j,i),j=1,nlobi)
355 write (iout,'(5(2x,3f8.4))') ((censc(k,j,i),k=1,3),
364 C Read scrot parameters for potentials determined from all-atom AM1 calculations
365 C added by Urszula Kozlowska 07/11/2007
368 read (irotam,*,end=112,err=112)
370 read (irotam,*,end=112,err=112)
373 read(irotam,*,end=112,err=112) sc_parmin(j,i)
378 C Read the parameters of the probability distribution/energy expression
379 C of the side chains.
382 read (irotam_pdb,'(3x,i3,f8.3)',end=112,err=112) nlob(i),dsc(i)
386 dsc_inv(i)=1.0D0/dsc(i)
397 read(irotam_pdb,*,end=112,err=112)(censc(k,1,i),k=1,3),
398 & ((blower(k,l,1),l=1,k),k=1,3)
400 read (irotam_pdb,*,end=112,err=112) bsc(j,i)
401 read (irotam_pdb,*,end=112,err=112) (censc(k,j,i),k=1,3),
402 & ((blower(k,l,j),l=1,k),k=1,3)
409 akl=akl+blower(k,m,j)*blower(l,m,j)
424 C Read torsional parameters in old format
426 read (itorp,*,end=113,err=113) ntortyp,nterm_old
427 if (lprint)write (iout,*) 'ntortyp,nterm',ntortyp,nterm_old
428 read (itorp,*,end=113,err=113) (itortyp(i),i=1,ntyp)
433 read (itorp,*,end=113,err=113) kk,v1(k,j,i),v2(k,j,i)
439 write (iout,'(/a/)') 'Torsional constants:'
442 write (iout,'(2i3,6f10.5)') i,j,(v1(k,i,j),k=1,nterm_old)
443 write (iout,'(6x,6f10.5)') (v2(k,i,j),k=1,nterm_old)
449 C Read torsional parameters
451 read (itorp,*,end=113,err=113) ntortyp
452 read (itorp,*,end=113,err=113) (itortyp(i),i=1,ntyp)
453 c write (iout,*) 'ntortyp',ntortyp
456 read (itorp,*,end=113,err=113) nterm(i,j),nlor(i,j)
460 read (itorp,*,end=113,err=113) kk,v1(k,i,j),v2(k,i,j)
461 v0ij=v0ij+si*v1(k,i,j)
465 read (itorp,*,end=113,err=113) kk,vlor1(k,i,j),
466 & vlor2(k,i,j),vlor3(k,i,j)
467 v0ij=v0ij+vlor1(k,i,j)/(1+vlor3(k,i,j)**2)
474 write (iout,'(/a/)') 'Torsional constants:'
477 write (iout,*) 'ityp',i,' jtyp',j
478 write (iout,*) 'Fourier constants'
480 write (iout,'(2(1pe15.5))') v1(k,i,j),v2(k,i,j)
482 write (iout,*) 'Lorenz constants'
484 write (iout,'(3(1pe15.5))')
485 & vlor1(k,i,j),vlor2(k,i,j),vlor3(k,i,j)
491 C 6/23/01 Read parameters for double torsionals
496 read (itordp,'(3a1)',end=114,err=114) t1,t2,t3
497 if (t1.ne.onelett(i) .or. t2.ne.onelett(j)
498 & .or. t3.ne.onelett(k)) then
499 write (iout,*) "Error in double torsional parameter file",
502 call MPI_Finalize(Ierror)
504 stop "Error in double torsional parameter file"
506 read (itordp,*,end=114,err=114) ntermd_1(i,j,k),
508 read (itordp,*,end=114,err=114) (v1c(1,l,i,j,k),l=1,
510 read (itordp,*,end=114,err=114) (v1s(1,l,i,j,k),l=1,
512 read (itordp,*,end=114,err=114) (v1c(2,l,i,j,k),l=1,
514 read (itordp,*,end=114,err=114) (v1s(2,l,i,j,k),l=1,
516 read (itordp,*,end=114,err=114) ((v2c(l,m,i,j,k),
517 & v2c(m,l,i,j,k),v2s(l,m,i,j,k),v2s(m,l,i,j,k),
518 & m=1,l-1),l=1,ntermd_2(i,j,k))
524 write (iout,*) 'Constants for double torsionals'
528 write (iout,*) 'ityp',i,' jtyp',j,' ktyp',k,
529 & ' nsingle',ntermd_1(i,j,k),' ndouble',ntermd_2(i,j,k)
531 write (iout,*) 'Single angles:'
532 do l=1,ntermd_1(i,j,k)
533 write (iout,'(i5,2f10.5,5x,2f10.5)') l,
534 & v1c(1,l,i,j,k),v1s(1,l,i,j,k),
535 & v1c(2,l,i,j,k),v1s(2,l,i,j,k)
538 write (iout,*) 'Pairs of angles:'
539 write (iout,'(3x,20i10)') (l,l=1,ntermd_2(i,j,k))
540 do l=1,ntermd_2(i,j,k)
541 write (iout,'(i5,20f10.5)')
542 & l,(v2c(l,m,i,j,k),m=1,ntermd_2(i,j,k))
545 write (iout,'(3x,20i10)') (l,l=1,ntermd_2(i,j,k))
546 do l=1,ntermd_2(i,j,k)
547 write (iout,'(i5,20f10.5)')
548 & l,(v2s(l,m,i,j,k),m=1,ntermd_2(i,j,k))
556 C Read of Side-chain backbone correlation parameters
557 C Modified 11 May 2012 by Adasko
560 read (isccor,*,end=113,err=113) nsccortyp
561 read (isccor,*,end=113,err=113) (isccortyp(i),i=1,ntyp)
562 c write (iout,*) 'ntortyp',ntortyp
564 cc maxinter is maximum interaction sites
568 read (isccor,*,end=113,err=113) nterm_sccor(i,j),nlor_sccor(i,j)
571 do k=1,nterm_sccor(i,j)
572 read (isccor,*,end=113,err=113) kk,v1sccor(k,i,j),v2sccor(k,i,j)
573 v0ijsccor=v0ijsccor+si*v1sccor(k,i,j)
576 do k=1,nlor_sccor(i,j)
577 read (isccor,*,end=113,err=113) kk,vlor1sccor(k,i,j),
578 & vlor2sccor(k,i,j),vlor3sccor(k,i,j)
579 v0ijsccor=v0ijsccor+vlor1sccor(k,i,j)/
580 &(1+vlor3sccor(k,i,j)**2)
582 v0sccor(i,j)=v0ijsccor
589 write (iout,'(/a/)') 'Torsional constants:'
592 write (iout,*) 'ityp',i,' jtyp',j
593 write (iout,*) 'Fourier constants'
594 do k=1,nterm_sccor(i,j)
595 write (iout,'(2(1pe15.5))') v1sccor(k,i,j),v2sccor(k,i,j)
597 write (iout,*) 'Lorenz constants'
598 do k=1,nlor_sccor(i,j)
599 write (iout,'(3(1pe15.5))')
600 & vlor1sccor(k,i,j),vlor2sccor(k,i,j),vlor3sccor(k,i,j)
607 C 9/18/99 (AL) Read coefficients of the Fourier expansion of the local
608 C interaction energy of the Gly, Ala, and Pro prototypes.
612 write (iout,*) "Coefficients of the cumulants"
614 read (ifourier,*) nloctyp
616 read (ifourier,*,end=115,err=115)
617 read (ifourier,*,end=115,err=115) (b(ii),ii=1,13)
619 write (iout,*) 'Type',i
620 write (iout,'(a,i2,a,f10.5)') ('b(',ii,')=',b(ii),ii=1,13)
646 c Ctilde(1,1,i)=0.0d0
647 c Ctilde(1,2,i)=0.0d0
648 c Ctilde(2,1,i)=0.0d0
649 c Ctilde(2,2,i)=0.0d0
662 c Dtilde(1,1,i)=0.0d0
663 c Dtilde(1,2,i)=0.0d0
664 c Dtilde(2,1,i)=0.0d0
665 c Dtilde(2,2,i)=0.0d0
666 EE(1,1,i)= b(10)+b(11)
667 EE(2,2,i)=-b(10)+b(11)
668 EE(2,1,i)= b(12)-b(13)
669 EE(1,2,i)= b(12)+b(13)
674 c ee(2,1,i)=ee(1,2,i)
678 write (iout,*) 'Type',i
680 write(iout,*) B1(1,i),B1(2,i)
682 write(iout,*) B2(1,i),B2(2,i)
685 write (iout,'(2f10.5)') CC(j,1,i),CC(j,2,i)
689 write (iout,'(2f10.5)') DD(j,1,i),DD(j,2,i)
693 write (iout,'(2f10.5)') EE(j,1,i),EE(j,2,i)
698 C Read electrostatic-interaction parameters
702 write (iout,'(/a)') 'Electrostatic interaction constants:'
703 write (iout,'(1x,a,1x,a,10x,a,11x,a,11x,a,11x,a)')
704 & 'IT','JT','APP','BPP','AEL6','AEL3'
706 read (ielep,*,end=116,err=116) ((epp(i,j),j=1,2),i=1,2)
707 read (ielep,*,end=116,err=116) ((rpp(i,j),j=1,2),i=1,2)
708 read (ielep,*,end=116,err=116) ((elpp6(i,j),j=1,2),i=1,2)
709 read (ielep,*,end=116,err=116) ((elpp3(i,j),j=1,2),i=1,2)
714 app (i,j)=epp(i,j)*rri*rri
715 bpp (i,j)=-2.0D0*epp(i,j)*rri
716 ael6(i,j)=elpp6(i,j)*4.2D0**6
717 ael3(i,j)=elpp3(i,j)*4.2D0**3
718 if (lprint) write(iout,'(2i3,4(1pe15.4))')i,j,app(i,j),bpp(i,j),
719 & ael6(i,j),ael3(i,j)
723 C Read side-chain interaction parameters.
725 read (isidep,*,end=117,err=117) ipot,expon
726 if (ipot.lt.1 .or. ipot.gt.5) then
727 write (iout,'(2a)') 'Error while reading SC interaction',
728 & 'potential file - unknown potential type.'
730 call MPI_Finalize(Ierror)
736 & write(iout,'(/3a,2i3)') 'Potential is ',potname(ipot),
737 & ', exponents are ',expon,2*expon
738 goto (10,20,30,30,40) ipot
739 C----------------------- LJ potential ---------------------------------
740 10 read (isidep,*,end=116,err=116)((eps(i,j),j=i,ntyp),i=1,ntyp),
741 & (sigma0(i),i=1,ntyp)
743 write (iout,'(/a/)') 'Parameters of the LJ potential:'
744 write (iout,'(a/)') 'The epsilon array:'
745 call printmat(ntyp,ntyp,ntyp,iout,restyp,eps)
746 write (iout,'(/a)') 'One-body parameters:'
747 write (iout,'(a,4x,a)') 'residue','sigma'
748 write (iout,'(a3,6x,f10.5)') (restyp(i),sigma0(i),i=1,ntyp)
751 C----------------------- LJK potential --------------------------------
752 20 read (isidep,*,end=116,err=116)((eps(i,j),j=i,ntyp),i=1,ntyp),
753 & (sigma0(i),i=1,ntyp),(rr0(i),i=1,ntyp)
755 write (iout,'(/a/)') 'Parameters of the LJK potential:'
756 write (iout,'(a/)') 'The epsilon array:'
757 call printmat(ntyp,ntyp,ntyp,iout,restyp,eps)
758 write (iout,'(/a)') 'One-body parameters:'
759 write (iout,'(a,4x,2a)') 'residue',' sigma ',' r0 '
760 write (iout,'(a3,6x,2f10.5)') (restyp(i),sigma0(i),rr0(i),
764 C---------------------- GB or BP potential -----------------------------
765 30 read (isidep,*,end=116,err=116)((eps(i,j),j=i,ntyp),i=1,ntyp),
766 & (sigma0(i),i=1,ntyp),(sigii(i),i=1,ntyp),(chip(i),i=1,ntyp),
768 C For the GB potential convert sigma'**2 into chi'
771 chip(i)=(chip(i)-1.0D0)/(chip(i)+1.0D0)
775 write (iout,'(/a/)') 'Parameters of the BP potential:'
776 write (iout,'(a/)') 'The epsilon array:'
777 call printmat(ntyp,ntyp,ntyp,iout,restyp,eps)
778 write (iout,'(/a)') 'One-body parameters:'
779 write (iout,'(a,4x,4a)') 'residue',' sigma ','s||/s_|_^2',
781 write (iout,'(a3,6x,4f10.5)') (restyp(i),sigma0(i),sigii(i),
782 & chip(i),alp(i),i=1,ntyp)
785 C--------------------- GBV potential -----------------------------------
786 40 read (isidep,*,end=116,err=116)((eps(i,j),j=i,ntyp),i=1,ntyp),
787 & (sigma0(i),i=1,ntyp),(rr0(i),i=1,ntyp),(sigii(i),i=1,ntyp),
788 & (chip(i),i=1,ntyp),(alp(i),i=1,ntyp)
790 write (iout,'(/a/)') 'Parameters of the GBV potential:'
791 write (iout,'(a/)') 'The epsilon array:'
792 call printmat(ntyp,ntyp,ntyp,iout,restyp,eps)
793 write (iout,'(/a)') 'One-body parameters:'
794 write (iout,'(a,4x,5a)') 'residue',' sigma ',' r0 ',
795 & 's||/s_|_^2',' chip ',' alph '
796 write (iout,'(a3,6x,5f10.5)') (restyp(i),sigma0(i),rr0(i),
797 & sigii(i),chip(i),alp(i),i=1,ntyp)
801 C-----------------------------------------------------------------------
802 C Calculate the "working" parameters of SC interactions.
810 sigma(i,j)=dsqrt(sigma0(i)**2+sigma0(j)**2)
811 sigma(j,i)=sigma(i,j)
812 rs0(i,j)=dwa16*sigma(i,j)
816 if (lprint) write (iout,'(/a/10x,7a/72(1h-))')
817 & 'Working parameters of the SC interactions:',
818 & ' a ',' b ',' augm ',' sigma ',' r0 ',
823 if (ipot.eq.1 .or. ipot.eq.3 .or. ipot.eq.4) then
832 sigeps=dsign(1.0D0,epsij)
834 aa(i,j)=epsij*rrij*rrij
835 bb(i,j)=-sigeps*epsij*rrij
843 ratsig1=sigt2sq/sigt1sq
844 ratsig2=1.0D0/ratsig1
845 chi(i,j)=(sigii1-1.0D0)/(sigii1+ratsig1)
846 if (j.gt.i) chi(j,i)=(sigii2-1.0D0)/(sigii2+ratsig2)
847 rsum_max=dsqrt(sigii1*sigt1sq+sigii2*sigt2sq)
851 c if (ipot.eq.1 .or. ipot.eq.3 .or. ipot.eq.4) then
852 sigmaii(i,j)=rsum_max
853 sigmaii(j,i)=rsum_max
855 c sigmaii(i,j)=r0(i,j)
856 c sigmaii(j,i)=r0(i,j)
858 cd write (iout,*) i,j,r0(i,j),sigma(i,j),rsum_max
859 if ((ipot.eq.2 .or. ipot.eq.5) .and. r0(i,j).gt.rsum_max) then
860 r_augm=sigma(i,j)*(rrij-sigma(i,j))/rrij
861 augm(i,j)=epsij*r_augm**(2*expon)
862 c augm(i,j)=0.5D0**(2*expon)*aa(i,j)
869 write (iout,'(2(a3,2x),3(1pe10.3),5(0pf8.3))')
870 & restyp(i),restyp(j),aa(i,j),bb(i,j),augm(i,j),
871 & sigma(i,j),r0(i,j),chi(i,j),chi(j,i)
877 C Define the SC-p interaction constants (hard-coded; old style)
880 C "Soft" SC-p repulsion (causes helices to be too flat, but facilitates
882 c aad(i,1)=0.3D0*4.0D0**12
883 C Following line for constants currently implemented
884 C "Hard" SC-p repulsion (gives correct turn spacing in helices)
885 aad(i,1)=1.5D0*4.0D0**12
886 c aad(i,1)=0.17D0*5.6D0**12
888 C "Soft" SC-p repulsion
890 C Following line for constants currently implemented
891 c aad(i,1)=0.3D0*4.0D0**6
892 C "Hard" SC-p repulsion
893 bad(i,1)=3.0D0*4.0D0**6
894 c bad(i,1)=-2.0D0*0.17D0*5.6D0**6
903 C 8/9/01 Read the SC-p interaction constants from file
906 read (iscpp,*,end=118,err=118) (eps_scp(i,j),rscp(i,j),j=1,2)
909 aad(i,1)=dabs(eps_scp(i,1))*rscp(i,1)**12
910 aad(i,2)=dabs(eps_scp(i,2))*rscp(i,2)**12
911 bad(i,1)=-2*eps_scp(i,1)*rscp(i,1)**6
912 bad(i,2)=-2*eps_scp(i,2)*rscp(i,2)**6
916 write (iout,*) "Parameters of SC-p interactions:"
918 write (iout,'(4f8.3,4e12.4)') eps_scp(i,1),rscp(i,1),
919 & eps_scp(i,2),rscp(i,2),aad(i,1),bad(i,1),aad(i,2),bad(i,2)
924 C Define the constants of the disulfide bridge
928 c Old arbitrary potential - commented out.
933 c Constants of the disulfide-bond potential determined based on the RHF/6-31G**
934 c energy surface of diethyl disulfide.
935 c A. Liwo and U. Kozlowska, 11/24/03
952 write (iout,'(/a)') "Disulfide bridge parameters:"
953 write (iout,'(a,f10.2)') 'S-S bridge energy: ',ebr
954 write (iout,'(2(a,f10.2))') 'd0cm:',d0cm,' akcm:',akcm
955 write (iout,'(2(a,f10.2))') 'akth:',akth,' akct:',akct
956 write (iout,'(3(a,f10.2))') 'v1ss:',v1ss,' v2ss:',v2ss,
960 111 write (iout,*) "Error reading bending energy parameters."
962 112 write (iout,*) "Error reading rotamer energy parameters."
964 113 write (iout,*) "Error reading torsional energy parameters."
966 114 write (iout,*) "Error reading double torsional energy parameters."
969 & "Error reading cumulant (multibody energy) parameters."
971 116 write (iout,*) "Error reading electrostatic energy parameters."
973 117 write (iout,*) "Error reading side chain interaction parameters."
975 118 write (iout,*) "Error reading SCp interaction parameters."
977 119 write (iout,*) "Error reading SCCOR parameters"
980 call MPI_Finalize(Ierror)
987 subroutine getenv_loc(var, val)
988 character(*) var, val
994 open (196,file='env',status='old',readonly,shared)
996 c write(*,*)'looking for ',var
997 10 read(196,*,err=11,end=11)line
998 iread=index(line,var)
999 c write(*,*)iread,' ',var,' ',line
1000 if (iread.eq.0) go to 10
1001 c write(*,*)'---> ',line
1007 iread=iread+ilen(var)+1
1008 read (line(iread:),*,err=12,end=12) val
1009 c write(*,*)'OK: ',var,' = ',val
1015 #elif (defined CRAY)
1016 integer lennam,lenval,ierror
1018 c getenv using a POSIX call, useful on the T3D
1019 c Sept 1996, comment out error check on advice of H. Pritchard
1022 if(lennam.le.0) stop '--error calling getenv--'
1023 call pxfgetenv(var,lennam,val,lenval,ierror)
1024 c-HP- if(ierror.ne.0) stop '--error returned by pxfgetenv--'
1026 call getenv(var,val)