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 C 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),ftors(i),
298 & 'There are',ndih_constr,' constraints on phi angles.'
300 write (iout,'(i5,3f8.3)') idih_constr(i),phi0(i),drange(i),
304 phi0(i)=deg2rad*phi0(i)
305 drange(i)=deg2rad*drange(i)
307 endif ! endif ndif_constr.gt.0
308 endif ! with_dihed_constr
311 print *,'NNT=',NNT,' NCT=',NCT
312 if (itype(1).eq.ntyp1) nnt=2
313 if (itype(nres).eq.ntyp1) nct=nct-1
314 if (nstart.lt.nnt) nstart=nnt
315 if (nend.gt.nct .or. nend.eq.0) nend=nct
316 write (iout,*) "nstart",nstart," nend",nend
319 c read(inp,'(a)') pdbfile
320 c write (iout,'(2a)') 'PDB data will be read from file ',pdbfile
321 c open(ipdbin,file=pdbfile,status='old',err=33)
323 c 33 write (iout,'(a)') 'Error opening PDB file.'
326 c print *,'Begin reading pdb data'
328 c print *,'Finished reading pdb data'
329 c write (iout,'(a,i3,a,i3)')'nsup=',nsup,' nstart_sup=',nstart_sup
331 c itype_pdb(i)=itype(i)
334 c write (iout,'(a,i3)') 'nsup=',nsup
336 c if (nsup.le.(nct-nnt+1)) then
337 c do i=0,nct-nnt+1-nsup
338 c if (seq_comp(itype(nnt+i),itype_pdb(nstart_sup),nsup)) then
344 c & 'Error - sequences to be superposed do not match.'
347 c do i=0,nsup-(nct-nnt+1)
348 c if (seq_comp(itype(nnt),itype_pdb(nstart_sup+i),nct-nnt+1))
350 c nstart_sup=nstart_sup+i
356 c & 'Error - sequences to be superposed do not match.'
359 c write (iout,*) 'nsup=',nsup,' nstart_sup=',nstart_sup,
360 c & ' nstart_seq=',nstart_seq
364 write (iout,*) "molread: REFSTR",refstr
366 if (.not.pdbref) then
367 call read_angles(inp,*38)
369 38 write (iout,'(a)') 'Error reading reference structure.'
371 call mp_stopall(Error_Msg)
373 stop 'Error reading reference structure'
386 call contact(.true.,ncont_ref,icont_ref)
389 C write (iout,'(/a,i3,a)')
390 C & 'The chain contains',ns,' disulfide-bridging cysteines.'
391 write (iout,'(20i4)') (iss(i),i=1,ns)
393 write(iout,*)"Running with dynamic disulfide-bond formation"
395 write (iout,'(/a/)') 'Pre-formed links are:'
401 write (iout,'(2a,i3,3a,i3,a,3f10.3)')
402 & restyp(it1),'(',i1,') -- ',restyp(it2),'(',i2,')',dhpb(i),
408 if (ns.gt.0.and.dyn_ss) then
412 forcon(i-nss)=forcon(i)
419 dyn_ss_mask(iss(i))=.true.
422 c Read distance restraints
423 if (constr_dist.gt.0) then
424 call read_dist_constr
429 c-----------------------------------------------------------------------------
430 logical function seq_comp(itypea,itypeb,length)
432 integer length,itypea(length),itypeb(length)
435 if (itypea(i).ne.itypeb(i)) then
443 c-----------------------------------------------------------------------------
444 subroutine read_bridge
445 C Read information about disulfide bridges.
448 include 'COMMON.IOUNITS'
451 include 'COMMON.INTERACT'
452 include 'COMMON.LOCAL'
453 include 'COMMON.NAMES'
454 include 'COMMON.CHAIN'
455 include 'COMMON.FFIELD'
456 include 'COMMON.SBRIDGE'
457 include 'COMMON.HEADER'
458 include 'COMMON.CONTROL'
459 include 'COMMON.TIME1'
461 include 'COMMON.INFO'
464 C Read bridging residues.
465 read (inp,*) ns,(iss(i),i=1,ns)
467 C Check whether the specified bridging residues are cystines.
469 if (itype(iss(i)).ne.1) then
470 write (iout,'(2a,i3,a)')
471 & 'Do you REALLY think that the residue ',
472 & restyp(itype(iss(i))),i,
473 & ' can form a disulfide bridge?!!!'
474 write (*,'(2a,i3,a)')
475 & 'Do you REALLY think that the residue ',
476 & restyp(itype(iss(i))),i,
477 & ' can form a disulfide bridge?!!!'
479 call mp_stopall(error_msg)
485 C Read preformed bridges.
487 read (inp,*) nss,(ihpb(i),jhpb(i),i=1,nss)
490 C Check if the residues involved in bridges are in the specified list of
494 if (ihpb(i).eq.ihpb(j).or.ihpb(i).eq.jhpb(j)
495 & .or.jhpb(i).eq.ihpb(j).or.jhpb(i).eq.jhpb(j)) then
496 write (iout,'(a,i3,a)') 'Disulfide pair',i,
497 & ' contains residues present in other pairs.'
498 write (*,'(a,i3,a)') 'Disulfide pair',i,
499 & ' contains residues present in other pairs.'
501 call mp_stopall(error_msg)
508 if (ihpb(i).eq.iss(j)) goto 10
510 write (iout,'(a,i3,a)') 'Pair',i,' contains unknown cystine.'
513 if (jhpb(i).eq.iss(j)) goto 20
515 write (iout,'(a,i3,a)') 'Pair',i,' contains unknown cystine.'
528 c----------------------------------------------------------------------------
529 subroutine read_angles(kanal,*)
534 include 'COMMON.CHAIN'
535 include 'COMMON.IOUNITS'
537 read (kanal,*,err=10,end=10) (theta(i),i=3,nres)
538 read (kanal,*,err=10,end=10) (phi(i),i=4,nres)
539 read (kanal,*,err=10,end=10) (alph(i),i=2,nres-1)
540 read (kanal,*,err=10,end=10) (omeg(i),i=2,nres-1)
542 theta(i)=deg2rad*theta(i)
543 phi(i)=deg2rad*phi(i)
544 alph(i)=deg2rad*alph(i)
545 omeg(i)=deg2rad*omeg(i)
550 c----------------------------------------------------------------------------
551 subroutine reada(rekord,lancuch,wartosc,default)
553 character*(*) rekord,lancuch
554 double precision wartosc,default
557 iread=index(rekord,lancuch)
562 iread=iread+ilen(lancuch)+1
563 read (rekord(iread:),*) wartosc
566 c----------------------------------------------------------------------------
567 subroutine multreada(rekord,lancuch,tablica,dim,default)
570 double precision tablica(dim),default
571 character*(*) rekord,lancuch
577 iread=index(rekord,lancuch)
578 if (iread.eq.0) return
579 iread=iread+ilen(lancuch)+1
580 read (rekord(iread:),*,end=10,err=10) (tablica(i),i=1,dim)
583 c----------------------------------------------------------------------------
584 subroutine readi(rekord,lancuch,wartosc,default)
586 character*(*) rekord,lancuch
587 integer wartosc,default
590 iread=index(rekord,lancuch)
595 iread=iread+ilen(lancuch)+1
596 read (rekord(iread:),*) wartosc
599 C----------------------------------------------------------------------
600 subroutine multreadi(rekord,lancuch,tablica,dim,default)
603 integer tablica(dim),default
604 character*(*) rekord,lancuch
611 iread=index(rekord,lancuch(:ilen(lancuch))//"=")
612 if (iread.eq.0) return
613 iread=iread+ilen(lancuch)+1
614 read (rekord(iread:),*,end=10,err=10) (tablica(i),i=1,dim)
618 c----------------------------------------------------------------------------
619 subroutine card_concat(card)
621 include 'COMMON.IOUNITS'
623 character*80 karta,ucase
625 read (inp,'(a)') karta
628 do while (karta(80:80).eq.'&')
629 card=card(:ilen(card)+1)//karta(:79)
630 read (inp,'(a)') karta
633 card=card(:ilen(card)+1)//karta
636 c----------------------------------------------------------------------------
645 include 'COMMON.IOUNITS'
646 include 'COMMON.CONTROL'
647 integer lenpre,lenpot,ilen
649 character*16 cformat,cprint
651 integer lenint,lenout
652 call getenv('INPUT',prefix)
653 call getenv('OUTPUT',prefout)
654 call getenv('INTIN',prefintin)
655 call getenv('COORD',cformat)
656 call getenv('PRINTCOOR',cprint)
657 call getenv('SCRATCHDIR',scratchdir)
660 if (index(ucase(cformat),'CX').gt.0) then
667 lenint=ilen(prefintin)
668 C Get the names and open the input files
669 open (inp,file=prefix(:ilen(prefix))//'.inp',status='old')
671 write (liczba,'(bz,i3.3)') me
672 outname=prefout(:lenout)//'_clust.out_'//liczba
674 outname=prefout(:lenout)//'_clust.out'
677 intinname=prefintin(:lenint)//'.bx'
678 else if (from_cx) then
679 intinname=prefintin(:lenint)//'.cx'
681 intinname=prefintin(:lenint)//'.int'
683 rmsname=prefintin(:lenint)//'.rms'
684 open (jplot,file=prefout(:ilen(prefout))//'.tex',
686 open (jrms,file=rmsname,status='unknown')
687 open(iout,file=outname,status='unknown')
688 C Get parameter filenames and open the parameter files.
689 call getenv('BONDPAR',bondname)
690 open (ibond,file=bondname,status='old')
691 call getenv('THETPAR',thetname)
692 open (ithep,file=thetname,status='old')
693 call getenv('ROTPAR',rotname)
694 open (irotam,file=rotname,status='old')
695 call getenv('TORPAR',torname)
696 open (itorp,file=torname,status='old')
697 call getenv('TORDPAR',tordname)
698 open (itordp,file=tordname,status='old')
699 call getenv('FOURIER',fouriername)
700 open (ifourier,file=fouriername,status='old')
701 call getenv('ELEPAR',elename)
702 open (ielep,file=elename,status='old')
703 call getenv('SIDEPAR',sidename)
704 open (isidep,file=sidename,status='old')
705 call getenv('SIDEP',sidepname)
706 open (isidep1,file=sidepname,status="old")
707 call getenv('SCCORPAR',sccorname)
708 open (isccor,file=sccorname,status="old")
711 C 8/9/01 In the newest version SCp interaction constants are read from a file
712 C Use -DOLDSCP to use hard-coded constants instead.
714 call getenv('SCPPAR',scpname)
715 open (iscpp,file=scpname,status='old')
719 subroutine read_dist_constr
720 implicit real*8 (a-h,o-z)
725 include 'COMMON.CONTROL'
726 include 'COMMON.CHAIN'
727 include 'COMMON.IOUNITS'
728 include 'COMMON.SBRIDGE'
729 integer ifrag_(2,100),ipair_(2,100)
730 double precision wfrag_(100),wpair_(100)
731 character*500 controlcard
732 logical lprn /.true./
733 write (iout,*) "Calling read_dist_constr"
734 C write (iout,*) "nres",nres," nstart_sup",nstart_sup," nsup",nsup
736 write(iout,*) "TU sie wywalam?"
737 call card_concat(controlcard)
738 write (iout,*) controlcard
740 call readi(controlcard,"NFRAG",nfrag_,0)
741 call readi(controlcard,"NPAIR",npair_,0)
742 call readi(controlcard,"NDIST",ndist_,0)
743 call reada(controlcard,'DIST_CUT',dist_cut,5.0d0)
744 call multreadi(controlcard,"IFRAG",ifrag_(1,1),2*nfrag_,0)
745 call multreadi(controlcard,"IPAIR",ipair_(1,1),2*npair_,0)
746 call multreada(controlcard,"WFRAG",wfrag_(1),nfrag_,0.0d0)
747 call multreada(controlcard,"WPAIR",wpair_(1),npair_,0.0d0)
748 write (iout,*) "NFRAG",nfrag_," NPAIR",npair_," NDIST",ndist_
749 write (iout,*) "IFRAG"
751 write (iout,*) i,ifrag_(1,i),ifrag_(2,i),wfrag_(i)
753 write (iout,*) "IPAIR"
755 write (iout,*) i,ipair_(1,i),ipair_(2,i),wpair_(i)
759 if (ifrag_(1,i).lt.nstart_sup) ifrag_(1,i)=nstart_sup
760 if (ifrag_(2,i).gt.nstart_sup+nsup-1)
761 & ifrag_(2,i)=nstart_sup+nsup-1
762 c write (iout,*) i,ifrag_(1,i),ifrag_(2,i),wfrag_(i)
764 if (wfrag_(i).gt.0.0d0) then
765 do j=ifrag_(1,i),ifrag_(2,i)-1
767 write (iout,*) "j",j," k",k
769 if (constr_dist.eq.1) then
774 forcon(nhpb)=wfrag_(i)
775 else if (constr_dist.eq.2) then
776 if (ddjk.le.dist_cut) then
781 forcon(nhpb)=wfrag_(i)
788 forcon(nhpb)=wfrag_(i)*dexp(-0.5d0*(ddjk/dist_cut)**2)
791 & write (iout,'(a,3i5,f8.2,1pe12.2)') "+dist.constr ",
792 & nhpb,ihpb(nhpb),jhpb(nhpb),dhpb(nhpb),forcon(nhpb)
798 if (wpair_(i).gt.0.0d0) then
806 do j=ifrag_(1,ii),ifrag_(2,ii)
807 do k=ifrag_(1,jj),ifrag_(2,jj)
811 forcon(nhpb)=wpair_(i)
813 write (iout,'(a,3i5,f8.2,f10.1)') "+dist.constr ",
814 & nhpb,ihpb(nhpb),jhpb(nhpb),dhpb(nhpb),forcon(nhpb)
820 if (constr_dist.eq.11) then
821 read (inp,*) ihpb(nhpb+1),jhpb(nhpb+1),dhpb(i),dhpb1(i),
822 & ibecarb(i),forcon(nhpb+1),fordepth(nhpb+1)
823 fordepth(nhpb+1)=fordepth(nhpb+1)/forcon(nhpb+1)
824 C write (iout,'(a,3i5,f8.2,f10.1)') "+dist.constr ",
825 C & nhpb,ihpb(nhpb),jhpb(nhpb),dhpb(nhpb),forcon(nhpb)
827 read (inp,*) ihpb(nhpb+1),jhpb(nhpb+1),forcon(nhpb+1)
829 if (forcon(nhpb+1).gt.0.0d0) then
831 if (ibecarb(i).gt.0) then
835 if (dhpb(nhpb).eq.0.0d0)
836 & dhpb(nhpb)=dist(ihpb(nhpb),jhpb(nhpb))
837 C dhpb(nhpb)=dist(ihpb(nhpb),jhpb(nhpb))
838 write (iout,'(a,3i5,f8.2,f10.1)') "+dist.constr ",
839 & nhpb,ihpb(nhpb),jhpb(nhpb),dhpb(nhpb),forcon(nhpb)