1 subroutine read_control
8 include 'COMMON.IOUNITS'
10 include 'COMMON.SBRIDGE'
11 include 'COMMON.CONTROL'
12 include 'COMMON.CLUSTER'
13 include 'COMMON.CHAIN'
14 include 'COMMON.HEADER'
15 include 'COMMON.FFIELD'
17 include 'COMMON.INTERACT'
18 character*320 controlcard,ucase
22 integer i,i1,i2,it1,it2
24 read (INP,'(a80)') titel
25 call card_concat(controlcard)
27 call readi(controlcard,'NRES',nres,0)
28 call readi(controlcard,'RESCALE',rescale_mode,2)
29 call reada(controlcard,'DISTCHAINMAX',distchainmax,50.0d0)
30 write (iout,*) "DISTCHAINMAX",distchainmax
31 call readi(controlcard,'PDBOUT',outpdb,0)
32 call readi(controlcard,'MOL2OUT',outmol2,0)
33 refstr=(index(controlcard,'REFSTR').gt.0)
34 write (iout,*) "REFSTR",refstr
35 pdbref=(index(controlcard,'PDBREF').gt.0)
36 iscode=index(controlcard,'ONE_LETTER')
37 tree=(index(controlcard,'MAKE_TREE').gt.0)
38 with_dihed_constr = index(controlcard,"WITH_DIHED_CONSTR").gt.0
39 call readi(controlcard,'CONSTR_DIST',constr_dist,0)
40 write (iout,*) "with_dihed_constr ",with_dihed_constr,
41 & " CONSTR_DIST",constr_dist
43 min_var=(index(controlcard,'MINVAR').gt.0)
44 plot_tree=(index(controlcard,'PLOT_TREE').gt.0)
45 punch_dist=(index(controlcard,'PUNCH_DIST').gt.0)
46 call readi(controlcard,'NCUT',ncut,1)
47 call readi(controlcard,'SYM',symetr,1)
48 write (iout,*) 'sym', symetr
49 call readi(controlcard,'NSTART',nstart,0)
50 call readi(controlcard,'NEND',nend,0)
51 call reada(controlcard,'ECUT',ecut,10.0d0)
52 call reada(controlcard,'PROB',prob_limit,0.99d0)
53 write (iout,*) "Probability limit",prob_limit
54 lgrp=(index(controlcard,'LGRP').gt.0)
55 caonly=(index(controlcard,'CA_ONLY').gt.0)
56 print_dist=(index(controlcard,'PRINT_DIST').gt.0)
57 call multreada(controlcard,'CUTOFF',rcutoff,ncut,-1.0d0)
58 call readi(controlcard,'IOPT',iopt,2)
59 lside = index(controlcard,"SIDE").gt.0
60 efree = index(controlcard,"EFREE").gt.0
61 call readi(controlcard,'NTEMP',nT,1)
62 write (iout,*) "nT",nT
63 call multreada(controlcard,'TEMPER',beta_h,nT,300.0d0)
64 write (iout,*) "nT",nT
65 write (iout,*) 'beta_h',(beta_h(i),i=1,nT)
67 beta_h(i)=1.0d0/(1.987D-3*beta_h(i))
69 write (iout,*) 'beta_h',(beta_h(i),i=1,nT)
70 lprint_cart=index(controlcard,"PRINT_CART") .gt.0
71 lprint_int=index(controlcard,"PRINT_INT") .gt.0
75 c--------------------------------------------------------------------------
78 C Read molecular data.
82 include 'COMMON.IOUNITS'
85 include 'COMMON.INTERACT'
86 include 'COMMON.LOCAL'
87 include 'COMMON.NAMES'
88 include 'COMMON.CHAIN'
89 include 'COMMON.FFIELD'
90 include 'COMMON.SBRIDGE'
91 include 'COMMON.HEADER'
92 include 'COMMON.CONTROL'
93 include 'COMMON.CONTACTS'
94 include 'COMMON.TIME1'
95 include 'COMMON.TORCNSTR'
99 character*4 sequence(maxres)
100 character*800 weightcard
102 double precision x(maxvar)
103 integer itype_pdb(maxres)
105 integer i,j,kkk,i1,i2,it1,it2
109 C Read weights of the subsequent energy terms.
110 call card_concat(weightcard)
111 call reada(weightcard,'WSC',wsc,1.0d0)
112 call reada(weightcard,'WLONG',wsc,wsc)
113 call reada(weightcard,'WSCP',wscp,1.0d0)
114 call reada(weightcard,'WELEC',welec,1.0D0)
115 call reada(weightcard,'WVDWPP',wvdwpp,welec)
116 call reada(weightcard,'WEL_LOC',wel_loc,1.0D0)
117 call reada(weightcard,'WCORR4',wcorr4,0.0D0)
118 call reada(weightcard,'WCORR5',wcorr5,0.0D0)
119 call reada(weightcard,'WCORR6',wcorr6,0.0D0)
120 call reada(weightcard,'WTURN3',wturn3,1.0D0)
121 call reada(weightcard,'WTURN4',wturn4,1.0D0)
122 call reada(weightcard,'WTURN6',wturn6,1.0D0)
123 call reada(weightcard,'WSTRAIN',wstrain,1.0D0)
124 call reada(weightcard,'WSCCOR',wsccor,1.0D0)
125 call reada(weightcard,'WBOND',wbond,1.0D0)
126 call reada(weightcard,'WTOR',wtor,1.0D0)
127 call reada(weightcard,'WTORD',wtor_d,1.0D0)
128 call reada(weightcard,'WANG',wang,1.0D0)
129 call reada(weightcard,'WSCLOC',wscloc,1.0D0)
130 call reada(weightcard,'SCAL14',scal14,0.4D0)
131 call reada(weightcard,'SCALSCP',scalscp,1.0d0)
132 call reada(weightcard,'CUTOFF',cutoff_corr,7.0d0)
133 call reada(weightcard,'DELT_CORR',delt_corr,0.5d0)
134 if (index(weightcard,'SOFT').gt.0) ipot=6
135 call reada(weightcard,"D0CM",d0cm,3.78d0)
136 call reada(weightcard,"AKCM",akcm,15.1d0)
137 call reada(weightcard,"AKTH",akth,11.0d0)
138 call reada(weightcard,"AKCT",akct,12.0d0)
139 call reada(weightcard,"V1SS",v1ss,-1.08d0)
140 call reada(weightcard,"V2SS",v2ss,7.61d0)
141 call reada(weightcard,"V3SS",v3ss,13.7d0)
142 call reada(weightcard,"EBR",ebr,-5.50D0)
143 call reada(weightcard,"ATRISS",atriss,0.301D0)
144 call reada(weightcard,"BTRISS",btriss,0.021D0)
145 call reada(weightcard,"CTRISS",ctriss,1.001D0)
146 call reada(weightcard,"DTRISS",dtriss,1.001D0)
147 write (iout,*) "ATRISS=", atriss
148 write (iout,*) "BTRISS=", btriss
149 write (iout,*) "CTRISS=", ctriss
150 write (iout,*) "DTRISS=", dtriss
151 dyn_ss=(index(weightcard,'DYN_SS').gt.0)
153 dyn_ss_mask(i)=.false.
157 dyn_ssbond_ij(i,j)=1.0d300
160 call reada(weightcard,"HT",Ht,0.0D0)
162 ss_depth=ebr/wsc-0.25*eps(1,1)
163 Ht=Ht/wsc-0.25*eps(1,1)
164 akcm=akcm*wstrain/wsc
165 akth=akth*wstrain/wsc
166 akct=akct*wstrain/wsc
167 v1ss=v1ss*wstrain/wsc
168 v2ss=v2ss*wstrain/wsc
169 v3ss=v3ss*wstrain/wsc
171 ss_depth=ebr/wstrain-0.25*eps(1,1)*wsc/wstrain
173 write (iout,'(/a)') "Disulfide bridge parameters:"
174 write (iout,'(a,f10.2)') 'S-S bridge energy: ',ebr
175 write (iout,'(2(a,f10.2))') 'd0cm:',d0cm,' akcm:',akcm
176 write (iout,'(2(a,f10.2))') 'akth:',akth,' akct:',akct
177 write (iout,'(3(a,f10.2))') 'v1ss:',v1ss,' v2ss:',v2ss,
180 C 12/1/95 Added weight for the multi-body term WCORR
181 call reada(weightcard,'WCORRH',wcorr,1.0D0)
182 if (wcorr4.gt.0.0d0) wcorr=wcorr4
202 write (iout,10) wsc,wscp,welec,wvdwpp,wbond,wang,wscloc,wtor,
203 & wtor_d,wstrain,wel_loc,wcorr,wcorr5,wcorr6,wturn3,
204 & wturn4,wturn6,wsccor
205 10 format (/'Energy-term weights (unscaled):'//
206 & 'WSCC= ',f10.6,' (SC-SC)'/
207 & 'WSCP= ',f10.6,' (SC-p)'/
208 & 'WELEC= ',f10.6,' (p-p electr)'/
209 & 'WVDWPP= ',f10.6,' (p-p VDW)'/
210 & 'WBOND= ',f10.6,' (stretching)'/
211 & 'WANG= ',f10.6,' (bending)'/
212 & 'WSCLOC= ',f10.6,' (SC local)'/
213 & 'WTOR= ',f10.6,' (torsional)'/
214 & 'WTORD= ',f10.6,' (double torsional)'/
215 & 'WSTRAIN=',f10.6,' (SS bridges & dist. cnstr.)'/
216 & 'WEL_LOC=',f10.6,' (multi-body 3-rd order)'/
217 & 'WCORR4= ',f10.6,' (multi-body 4th order)'/
218 & 'WCORR5= ',f10.6,' (multi-body 5th order)'/
219 & 'WCORR6= ',f10.6,' (multi-body 6th order)'/
220 & 'WTURN3= ',f10.6,' (turns, 3rd order)'/
221 & 'WTURN4= ',f10.6,' (turns, 4th order)'/
222 & 'WTURN6= ',f10.6,' (turns, 6th order)'/
223 & 'WSCCOR= ',f10.6,' (SC-backbone torsinal correalations)')
225 if (wcorr4.gt.0.0d0) then
226 write (iout,'(/2a/)') 'Local-electrostatic type correlation ',
227 & 'between contact pairs of peptide groups'
228 write (iout,'(2(a,f5.3/))')
229 & 'Cutoff on 4-6th order correlation terms: ',cutoff_corr,
230 & 'Range of quenching the correlation terms:',2*delt_corr
231 else if (wcorr.gt.0.0d0) then
232 write (iout,'(/2a/)') 'Hydrogen-bonding correlation ',
233 & 'between contact pairs of peptide groups'
235 write (iout,'(a,f8.3)')
236 & 'Scaling factor of 1,4 SC-p interactions:',scal14
237 write (iout,'(a,f8.3)')
238 & 'General scaling factor of SC-p interactions:',scalscp
239 r0_corr=cutoff_corr-delt_corr
241 aad(i,1)=scalscp*aad(i,1)
242 aad(i,2)=scalscp*aad(i,2)
243 bad(i,1)=scalscp*bad(i,1)
244 bad(i,2)=scalscp*bad(i,2)
248 print *,'indpdb=',indpdb,' pdbref=',pdbref
250 C Read sequence if not taken from the pdb file.
251 if (iscode.gt.0) then
252 read (inp,'(80a1)') (sequence(i)(1:1),i=1,nres)
254 read (inp,'(20(1x,a3))') (sequence(i),i=1,nres)
256 C Convert sequence to numeric code
258 itype(i)=rescode(i,sequence(i),iscode)
261 print '(20i4)',(itype(i),i=1,nres)
265 if (itype(i).eq.ntyp1 .or. itype(i+1).eq.ntyp1) then
267 if (itype(i).eq.ntyp1) then
271 else if (iabs(itype(i+1)).ne.20) then
273 else if (iabs(itype(i)).ne.20) then
280 write (iout,*) "ITEL"
282 write (iout,*) i,itype(i),itel(i)
285 print *,'Call Read_Bridge.'
287 C this fragment reads diheadral constrains
288 if (with_dihed_constr) then
290 read (inp,*) ndih_constr
291 if (ndih_constr.gt.0) then
293 write (iout,*) 'FTORS',ftors
294 C ftors is the force constant for torsional quartic constrains
295 read (inp,*) (idih_constr(i),phi0(i),drange(i),i=1,ndih_constr)
297 & 'There are',ndih_constr,' constraints on phi angles.'
299 write (iout,'(i5,2f8.3)') idih_constr(i),phi0(i),drange(i)
302 phi0(i)=deg2rad*phi0(i)
303 drange(i)=deg2rad*drange(i)
305 endif ! endif ndif_constr.gt.0
306 endif ! with_dihed_constr
309 print *,'NNT=',NNT,' NCT=',NCT
310 if (itype(1).eq.ntyp1) nnt=2
311 if (itype(nres).eq.ntyp1) nct=nct-1
312 if (nstart.lt.nnt) nstart=nnt
313 if (nend.gt.nct .or. nend.eq.0) nend=nct
314 write (iout,*) "nstart",nstart," nend",nend
317 c read(inp,'(a)') pdbfile
318 c write (iout,'(2a)') 'PDB data will be read from file ',pdbfile
319 c open(ipdbin,file=pdbfile,status='old',err=33)
321 c 33 write (iout,'(a)') 'Error opening PDB file.'
324 c print *,'Begin reading pdb data'
326 c print *,'Finished reading pdb data'
327 c write (iout,'(a,i3,a,i3)')'nsup=',nsup,' nstart_sup=',nstart_sup
329 c itype_pdb(i)=itype(i)
332 c write (iout,'(a,i3)') 'nsup=',nsup
334 c if (nsup.le.(nct-nnt+1)) then
335 c do i=0,nct-nnt+1-nsup
336 c if (seq_comp(itype(nnt+i),itype_pdb(nstart_sup),nsup)) then
342 c & 'Error - sequences to be superposed do not match.'
345 c do i=0,nsup-(nct-nnt+1)
346 c if (seq_comp(itype(nnt),itype_pdb(nstart_sup+i),nct-nnt+1))
348 c nstart_sup=nstart_sup+i
354 c & 'Error - sequences to be superposed do not match.'
357 c write (iout,*) 'nsup=',nsup,' nstart_sup=',nstart_sup,
358 c & ' nstart_seq=',nstart_seq
362 write (iout,*) "molread: REFSTR",refstr
364 if (.not.pdbref) then
365 call read_angles(inp,*38)
367 38 write (iout,'(a)') 'Error reading reference structure.'
369 call mp_stopall(Error_Msg)
371 stop 'Error reading reference structure'
384 call contact(.true.,ncont_ref,icont_ref)
387 C write (iout,'(/a,i3,a)')
388 C & 'The chain contains',ns,' disulfide-bridging cysteines.'
389 write (iout,'(20i4)') (iss(i),i=1,ns)
391 write(iout,*)"Running with dynamic disulfide-bond formation"
393 write (iout,'(/a/)') 'Pre-formed links are:'
399 write (iout,'(2a,i3,3a,i3,a,3f10.3)')
400 & restyp(it1),'(',i1,') -- ',restyp(it2),'(',i2,')',dhpb(i),
406 if (ns.gt.0.and.dyn_ss) then
410 forcon(i-nss)=forcon(i)
417 dyn_ss_mask(iss(i))=.true.
420 c Read distance restraints
421 if (constr_dist.gt.0) then
422 call read_dist_constr
427 c-----------------------------------------------------------------------------
428 logical function seq_comp(itypea,itypeb,length)
430 integer length,itypea(length),itypeb(length)
433 if (itypea(i).ne.itypeb(i)) then
441 c-----------------------------------------------------------------------------
442 subroutine read_bridge
443 C Read information about disulfide bridges.
446 include 'COMMON.IOUNITS'
449 include 'COMMON.INTERACT'
450 include 'COMMON.LOCAL'
451 include 'COMMON.NAMES'
452 include 'COMMON.CHAIN'
453 include 'COMMON.FFIELD'
454 include 'COMMON.SBRIDGE'
455 include 'COMMON.HEADER'
456 include 'COMMON.CONTROL'
457 include 'COMMON.TIME1'
459 include 'COMMON.INFO'
462 C Read bridging residues.
463 read (inp,*) ns,(iss(i),i=1,ns)
465 C Check whether the specified bridging residues are cystines.
467 if (itype(iss(i)).ne.1) then
468 write (iout,'(2a,i3,a)')
469 & 'Do you REALLY think that the residue ',
470 & restyp(itype(iss(i))),i,
471 & ' can form a disulfide bridge?!!!'
472 write (*,'(2a,i3,a)')
473 & 'Do you REALLY think that the residue ',
474 & restyp(itype(iss(i))),i,
475 & ' can form a disulfide bridge?!!!'
477 call mp_stopall(error_msg)
483 C Read preformed bridges.
485 read (inp,*) nss,(ihpb(i),jhpb(i),i=1,nss)
488 C Check if the residues involved in bridges are in the specified list of
492 if (ihpb(i).eq.ihpb(j).or.ihpb(i).eq.jhpb(j)
493 & .or.jhpb(i).eq.ihpb(j).or.jhpb(i).eq.jhpb(j)) then
494 write (iout,'(a,i3,a)') 'Disulfide pair',i,
495 & ' contains residues present in other pairs.'
496 write (*,'(a,i3,a)') 'Disulfide pair',i,
497 & ' contains residues present in other pairs.'
499 call mp_stopall(error_msg)
506 if (ihpb(i).eq.iss(j)) goto 10
508 write (iout,'(a,i3,a)') 'Pair',i,' contains unknown cystine.'
511 if (jhpb(i).eq.iss(j)) goto 20
513 write (iout,'(a,i3,a)') 'Pair',i,' contains unknown cystine.'
526 c----------------------------------------------------------------------------
527 subroutine read_angles(kanal,*)
532 include 'COMMON.CHAIN'
533 include 'COMMON.IOUNITS'
535 read (kanal,*,err=10,end=10) (theta(i),i=3,nres)
536 read (kanal,*,err=10,end=10) (phi(i),i=4,nres)
537 read (kanal,*,err=10,end=10) (alph(i),i=2,nres-1)
538 read (kanal,*,err=10,end=10) (omeg(i),i=2,nres-1)
540 theta(i)=deg2rad*theta(i)
541 phi(i)=deg2rad*phi(i)
542 alph(i)=deg2rad*alph(i)
543 omeg(i)=deg2rad*omeg(i)
548 c----------------------------------------------------------------------------
549 subroutine reada(rekord,lancuch,wartosc,default)
551 character*(*) rekord,lancuch
552 double precision wartosc,default
555 iread=index(rekord,lancuch)
560 iread=iread+ilen(lancuch)+1
561 read (rekord(iread:),*) wartosc
564 c----------------------------------------------------------------------------
565 subroutine multreada(rekord,lancuch,tablica,dim,default)
568 double precision tablica(dim),default
569 character*(*) rekord,lancuch
575 iread=index(rekord,lancuch)
576 if (iread.eq.0) return
577 iread=iread+ilen(lancuch)+1
578 read (rekord(iread:),*,end=10,err=10) (tablica(i),i=1,dim)
581 c----------------------------------------------------------------------------
582 subroutine readi(rekord,lancuch,wartosc,default)
584 character*(*) rekord,lancuch
585 integer wartosc,default
588 iread=index(rekord,lancuch)
593 iread=iread+ilen(lancuch)+1
594 read (rekord(iread:),*) wartosc
597 C----------------------------------------------------------------------
598 subroutine multreadi(rekord,lancuch,tablica,dim,default)
601 integer tablica(dim),default
602 character*(*) rekord,lancuch
609 iread=index(rekord,lancuch(:ilen(lancuch))//"=")
610 if (iread.eq.0) return
611 iread=iread+ilen(lancuch)+1
612 read (rekord(iread:),*,end=10,err=10) (tablica(i),i=1,dim)
616 c----------------------------------------------------------------------------
617 subroutine card_concat(card)
619 include 'COMMON.IOUNITS'
621 character*80 karta,ucase
623 read (inp,'(a)') karta
626 do while (karta(80:80).eq.'&')
627 card=card(:ilen(card)+1)//karta(:79)
628 read (inp,'(a)') karta
631 card=card(:ilen(card)+1)//karta
634 c----------------------------------------------------------------------------
643 include 'COMMON.IOUNITS'
644 include 'COMMON.CONTROL'
645 integer lenpre,lenpot,ilen
647 character*16 cformat,cprint
649 integer lenint,lenout
650 call getenv('INPUT',prefix)
651 call getenv('OUTPUT',prefout)
652 call getenv('INTIN',prefintin)
653 call getenv('COORD',cformat)
654 call getenv('PRINTCOOR',cprint)
655 call getenv('SCRATCHDIR',scratchdir)
658 if (index(ucase(cformat),'CX').gt.0) then
665 lenint=ilen(prefintin)
666 C Get the names and open the input files
667 open (inp,file=prefix(:ilen(prefix))//'.inp',status='old')
669 write (liczba,'(bz,i3.3)') me
670 outname=prefout(:lenout)//'_clust.out_'//liczba
672 outname=prefout(:lenout)//'_clust.out'
675 intinname=prefintin(:lenint)//'.bx'
676 else if (from_cx) then
677 intinname=prefintin(:lenint)//'.cx'
679 intinname=prefintin(:lenint)//'.int'
681 rmsname=prefintin(:lenint)//'.rms'
682 open (jplot,file=prefout(:ilen(prefout))//'.tex',
684 open (jrms,file=rmsname,status='unknown')
685 open(iout,file=outname,status='unknown')
686 C Get parameter filenames and open the parameter files.
687 call getenv('BONDPAR',bondname)
688 open (ibond,file=bondname,status='old')
689 call getenv('THETPAR',thetname)
690 open (ithep,file=thetname,status='old')
691 call getenv('ROTPAR',rotname)
692 open (irotam,file=rotname,status='old')
693 call getenv('TORPAR',torname)
694 open (itorp,file=torname,status='old')
695 call getenv('TORDPAR',tordname)
696 open (itordp,file=tordname,status='old')
697 call getenv('FOURIER',fouriername)
698 open (ifourier,file=fouriername,status='old')
699 call getenv('ELEPAR',elename)
700 open (ielep,file=elename,status='old')
701 call getenv('SIDEPAR',sidename)
702 open (isidep,file=sidename,status='old')
703 call getenv('SIDEP',sidepname)
704 open (isidep1,file=sidepname,status="old")
705 call getenv('SCCORPAR',sccorname)
706 open (isccor,file=sccorname,status="old")
709 C 8/9/01 In the newest version SCp interaction constants are read from a file
710 C Use -DOLDSCP to use hard-coded constants instead.
712 call getenv('SCPPAR',scpname)
713 open (iscpp,file=scpname,status='old')
717 subroutine read_dist_constr
718 implicit real*8 (a-h,o-z)
723 include 'COMMON.CONTROL'
724 include 'COMMON.CHAIN'
725 include 'COMMON.IOUNITS'
726 include 'COMMON.SBRIDGE'
727 integer ifrag_(2,100),ipair_(2,100)
728 double precision wfrag_(100),wpair_(100)
729 character*500 controlcard
730 logical lprn /.true./
731 write (iout,*) "Calling read_dist_constr"
732 C write (iout,*) "nres",nres," nstart_sup",nstart_sup," nsup",nsup
734 write(iout,*) "TU sie wywalam?"
735 call card_concat(controlcard)
736 write (iout,*) controlcard
738 call readi(controlcard,"NFRAG",nfrag_,0)
739 call readi(controlcard,"NPAIR",npair_,0)
740 call readi(controlcard,"NDIST",ndist_,0)
741 call reada(controlcard,'DIST_CUT',dist_cut,5.0d0)
742 call multreadi(controlcard,"IFRAG",ifrag_(1,1),2*nfrag_,0)
743 call multreadi(controlcard,"IPAIR",ipair_(1,1),2*npair_,0)
744 call multreada(controlcard,"WFRAG",wfrag_(1),nfrag_,0.0d0)
745 call multreada(controlcard,"WPAIR",wpair_(1),npair_,0.0d0)
746 write (iout,*) "NFRAG",nfrag_," NPAIR",npair_," NDIST",ndist_
747 write (iout,*) "IFRAG"
749 write (iout,*) i,ifrag_(1,i),ifrag_(2,i),wfrag_(i)
751 write (iout,*) "IPAIR"
753 write (iout,*) i,ipair_(1,i),ipair_(2,i),wpair_(i)
757 if (ifrag_(1,i).lt.nstart_sup) ifrag_(1,i)=nstart_sup
758 if (ifrag_(2,i).gt.nstart_sup+nsup-1)
759 & ifrag_(2,i)=nstart_sup+nsup-1
760 c write (iout,*) i,ifrag_(1,i),ifrag_(2,i),wfrag_(i)
762 if (wfrag_(i).gt.0.0d0) then
763 do j=ifrag_(1,i),ifrag_(2,i)-1
765 write (iout,*) "j",j," k",k
767 if (constr_dist.eq.1) then
772 forcon(nhpb)=wfrag_(i)
773 else if (constr_dist.eq.2) then
774 if (ddjk.le.dist_cut) then
779 forcon(nhpb)=wfrag_(i)
786 forcon(nhpb)=wfrag_(i)*dexp(-0.5d0*(ddjk/dist_cut)**2)
789 & write (iout,'(a,3i5,f8.2,1pe12.2)') "+dist.constr ",
790 & nhpb,ihpb(nhpb),jhpb(nhpb),dhpb(nhpb),forcon(nhpb)
796 if (wpair_(i).gt.0.0d0) then
804 do j=ifrag_(1,ii),ifrag_(2,ii)
805 do k=ifrag_(1,jj),ifrag_(2,jj)
809 forcon(nhpb)=wpair_(i)
811 write (iout,'(a,3i5,f8.2,f10.1)') "+dist.constr ",
812 & nhpb,ihpb(nhpb),jhpb(nhpb),dhpb(nhpb),forcon(nhpb)
818 if (constr_dist.eq.11) then
819 read (inp,*) ihpb(nhpb+1),jhpb(nhpb+1),dhpb(i),dhpb1(i),
820 & ibecarb(i),forcon(nhpb+1),fordepth(nhpb+1)
821 fordepth(nhpb+1)=fordepth(nhpb+1)/forcon(nhpb+1)
822 C write (iout,'(a,3i5,f8.2,f10.1)') "+dist.constr ",
823 C & nhpb,ihpb(nhpb),jhpb(nhpb),dhpb(nhpb),forcon(nhpb)
825 read (inp,*) ihpb(nhpb+1),jhpb(nhpb+1),forcon(nhpb+1)
827 if (forcon(nhpb+1).gt.0.0d0) then
829 if (ibecarb(i).gt.0) then
833 if (dhpb(nhpb).eq.0.0d0)
834 & dhpb(nhpb)=dist(ihpb(nhpb),jhpb(nhpb))
835 C dhpb(nhpb)=dist(ihpb(nhpb),jhpb(nhpb))
836 write (iout,'(a,3i5,f8.2,f10.1)') "+dist.constr ",
837 & nhpb,ihpb(nhpb),jhpb(nhpb),dhpb(nhpb),forcon(nhpb)