+++ /dev/null
- subroutine readrtns
- implicit real*8 (a-h,o-z)
- include 'DIMENSIONS'
-#ifdef MPI
- include 'mpif.h'
-#endif
- include 'COMMON.SETUP'
- include 'COMMON.CONTROL'
- include 'COMMON.SBRIDGE'
- include 'COMMON.IOUNITS'
- logical file_exist
-C Read force-field parameters except weights
- call parmread
-C Read job setup parameters
- call read_control
-C Read control parameters for energy minimzation if required
- if (minim) call read_minim
-C Read MCM control parameters if required
-c if (modecalc.eq.3 .or. modecalc.eq.6) call mcmread
-C Read MD control parameters if reqjuired
-c if (modecalc.eq.12) call read_MDpar
-C Read MREMD control parameters if required
-c if (modecalc.eq.14) then
-c call read_MDpar
-c call read_REMDpar
-c endif
-C Read MUCA control parameters if required
-c if (lmuca) call read_muca
-C Read CSA control parameters if required (from fort.40 if exists
-C otherwise from general input file)
- if (modecalc.eq.8) then
- inquire (file="fort.40",exist=file_exist)
- if (.not.file_exist) call csaread
- endif
-cfmc if (modecalc.eq.10) call mcmfread
-C Read molecule information, molecule geometry, energy-term weights, and
-C restraints if requested
- call molread
-
-C Print restraint information
-#ifdef MPI
- if (.not. out1file .or. me.eq.king) then
-#endif
- if (nhpb.gt.nss)
- &write (iout,'(a,i5,a)') "The following",nhpb-nss,
- & " distance constraints have been imposed"
- do i=nss+1,nhpb
- write (iout,'(3i6,f10.5)') i-nss,ihpb(i),jhpb(i),forcon(i)
- enddo
-#ifdef MPI
- endif
-#endif
-c print *,"Processor",myrank," leaves READRTNS"
- return
- end
-C-------------------------------------------------------------------------------
- subroutine read_control
-C
-C Read contorl data
-C
- implicit real*8 (a-h,o-z)
- include 'DIMENSIONS'
-#ifdef MP
- include 'mpif.h'
- logical OKRandom, prng_restart
- real*8 r1
-#endif
- include 'COMMON.IOUNITS'
- include 'COMMON.TIME1'
-c include 'COMMON.THREAD'
- include 'COMMON.SBRIDGE'
- include 'COMMON.CONTROL'
- include 'COMMON.MCM'
-c include 'COMMON.MAP'
- include 'COMMON.HEADER'
-c include 'COMMON.CSA'
- include 'COMMON.CHAIN'
-c include 'COMMON.MUCA'
-c include 'COMMON.MD'
- include 'COMMON.FFIELD'
- include 'COMMON.SETUP'
- COMMON /MACHSW/ KDIAG,ICORFL,IXDR
- character*8 diagmeth(0:3) /'Library','EVVRSP','Givens','Jacobi'/
- character*80 ucase
- character*320 controlcard
-
- nglob_csa=0
- eglob_csa=1d99
- nmin_csa=0
- read (INP,'(a)') titel
- call card_concat(controlcard)
-c out1file=index(controlcard,'OUT1FILE').gt.0 .or. fg_rank.gt.0
-c print *,"Processor",me," fg_rank",fg_rank," out1file",out1file
- call reada(controlcard,'SEED',seed,0.0D0)
- call random_init(seed)
-C Set up the time limit (caution! The time must be input in minutes!)
- read_cart=index(controlcard,'READ_CART').gt.0
- call readi(controlcard,'CONSTR_DIST',constr_dist,0)
- call reada(controlcard,'TIMLIM',timlim,960.0D0) ! default 16 hours
- unres_pdb = index(controlcard,'UNRES_PDB') .gt. 0
- call reada(controlcard,'SAFETY',safety,30.0D0) ! default 30 minutes
- call reada(controlcard,'RMSDBC',rmsdbc,3.0D0)
- call reada(controlcard,'RMSDBC1',rmsdbc1,0.5D0)
- call reada(controlcard,'RMSDBC1MAX',rmsdbc1max,1.5D0)
- call reada(controlcard,'RMSDBCM',rmsdbcm,3.0D0)
- call reada(controlcard,'DRMS',drms,0.1D0)
- if(me.eq.king .or. .not. out1file .and. fg_rank.eq.0) then
- write (iout,'(a,f10.1)')'RMSDBC = ',rmsdbc
- write (iout,'(a,f10.1)')'RMSDBC1 = ',rmsdbc1
- write (iout,'(a,f10.1)')'RMSDBC1MAX = ',rmsdbc1max
- write (iout,'(a,f10.1)')'DRMS = ',drms
- write (iout,'(a,f10.1)')'RMSDBCM = ',rmsdbcm
- write (iout,'(a,f10.1)') 'Time limit (min):',timlim
- endif
- call readi(controlcard,'NZ_START',nz_start,0)
- call readi(controlcard,'NZ_END',nz_end,0)
- call readi(controlcard,'IZ_SC',iz_sc,0)
- timlim=60.0D0*timlim
- safety = 60.0d0*safety
- timem=timlim
- modecalc=0
- call reada(controlcard,"T_BATH",t_bath,300.0d0)
- minim=(index(controlcard,'MINIMIZE').gt.0)
- dccart=(index(controlcard,'CART').gt.0)
- overlapsc=(index(controlcard,'OVERLAP').gt.0)
- overlapsc=.not.overlapsc
- searchsc=(index(controlcard,'NOSEARCHSC').gt.0)
- searchsc=.not.searchsc
- sideadd=(index(controlcard,'SIDEADD').gt.0)
- energy_dec=(index(controlcard,'ENERGY_DEC').gt.0)
- outpdb=(index(controlcard,'PDBOUT').gt.0)
- outmol2=(index(controlcard,'MOL2OUT').gt.0)
- pdbref=(index(controlcard,'PDBREF').gt.0)
- refstr=pdbref .or. (index(controlcard,'REFSTR').gt.0)
- indpdb=index(controlcard,'PDBSTART')
- extconf=(index(controlcard,'EXTCONF').gt.0)
- call readi(controlcard,'IPRINT',iprint,0)
- call readi(controlcard,'MAXGEN',maxgen,10000)
- call readi(controlcard,'MAXOVERLAP',maxoverlap,1000)
- call readi(controlcard,"KDIAG",kdiag,0)
- call readi(controlcard,"RESCALE_MODE",rescale_mode,0)
- if(me.eq.king .or. .not. out1file .and. fg_rank.eq.0)
- & write (iout,*) "RESCALE_MODE",rescale_mode
- split_ene=index(controlcard,'SPLIT_ENE').gt.0
- if (index(controlcard,'REGULAR').gt.0.0D0) then
- call reada(controlcard,'WEIDIS',weidis,0.1D0)
- modecalc=1
- refstr=.true.
- endif
- if (index(controlcard,'CHECKGRAD').gt.0) then
- modecalc=5
- if (index(controlcard,'CART').gt.0) then
- icheckgrad=1
- elseif (index(controlcard,'CARINT').gt.0) then
- icheckgrad=2
- else
- icheckgrad=3
- endif
- elseif (index(controlcard,'THREAD').gt.0) then
- modecalc=2
- call readi(controlcard,'THREAD',nthread,0)
- if (nthread.gt.0) then
- call reada(controlcard,'WEIDIS',weidis,0.1D0)
- else
- if (fg_rank.eq.0)
- & write (iout,'(a)')'A number has to follow the THREAD keyword.'
- stop 'Error termination in Read_Control.'
- endif
- else if (index(controlcard,'MCMA').gt.0) then
- modecalc=3
- else if (index(controlcard,'MCEE').gt.0) then
- modecalc=6
- else if (index(controlcard,'MULTCONF').gt.0) then
- modecalc=4
- else if (index(controlcard,'MAP').gt.0) then
- modecalc=7
- call readi(controlcard,'MAP',nmap,0)
- else if (index(controlcard,'CSA').gt.0) then
- modecalc=8
-crc else if (index(controlcard,'ZSCORE').gt.0) then
-crc
-crc ZSCORE is rm from UNRES, modecalc=9 is available
-crc
-crc modecalc=9
-cfcm else if (index(controlcard,'MCMF').gt.0) then
-cfmc modecalc=10
- else if (index(controlcard,'SOFTREG').gt.0) then
- modecalc=11
- else if (index(controlcard,'CHECK_BOND').gt.0) then
- modecalc=-1
- else if (index(controlcard,'TEST').gt.0) then
- modecalc=-2
- else if (index(controlcard,'MD').gt.0) then
- modecalc=12
- else if (index(controlcard,'RE ').gt.0) then
- modecalc=14
- endif
-
- lmuca=index(controlcard,'MUCA').gt.0
- call readi(controlcard,'MUCADYN',mucadyn,0)
- call readi(controlcard,'MUCASMOOTH',muca_smooth,0)
- if (lmuca .and. (me.eq.king .or. .not.out1file ))
- & then
- write (iout,*) 'MUCADYN=',mucadyn
- write (iout,*) 'MUCASMOOTH=',muca_smooth
- endif
-
- iscode=index(controlcard,'ONE_LETTER')
- indphi=index(controlcard,'PHI')
- indback=index(controlcard,'BACK')
- iranconf=index(controlcard,'RAND_CONF')
- i2ndstr=index(controlcard,'USE_SEC_PRED')
- gradout=index(controlcard,'GRADOUT').gt.0
- gnorm_check=index(controlcard,'GNORM_CHECK').gt.0
-
- if(me.eq.king.or..not.out1file)
- & write (iout,'(2a)') diagmeth(kdiag),
- & ' routine used to diagonalize matrices.'
- return
- end
-c------------------------------------------------------------------------------
- subroutine molread
-C
-C Read molecular data.
-C
- implicit real*8 (a-h,o-z)
- include 'DIMENSIONS'
-#ifdef MPI
- include 'mpif.h'
- integer error_msg
-#endif
- include 'COMMON.IOUNITS'
- include 'COMMON.GEO'
- include 'COMMON.VAR'
- include 'COMMON.INTERACT'
- include 'COMMON.LOCAL'
- include 'COMMON.NAMES'
- include 'COMMON.CHAIN'
- include 'COMMON.FFIELD'
- include 'COMMON.SBRIDGE'
- include 'COMMON.HEADER'
- include 'COMMON.CONTROL'
-c include 'COMMON.DBASE'
-c include 'COMMON.THREAD'
- include 'COMMON.CONTACTS'
- include 'COMMON.TORCNSTR'
- include 'COMMON.TIME1'
- include 'COMMON.BOUNDS'
-c include 'COMMON.MD'
-c include 'COMMON.REMD'
- include 'COMMON.SETUP'
- character*4 sequence(maxres)
- integer rescode
- double precision x(maxvar)
- character*256 pdbfile
- character*320 weightcard
- character*80 weightcard_t,ucase
- dimension itype_pdb(maxres)
- common /pizda/ itype_pdb
- logical seq_comp,fail
- double precision energia(0:n_ene)
- integer ilen
- external ilen
-C
-C Body
-C
-C Read weights of the subsequent energy terms.
-
- call card_concat(weightcard)
- call reada(weightcard,'WLONG',wlong,1.0D0)
- call reada(weightcard,'WSC',wsc,wlong)
- call reada(weightcard,'WSCP',wscp,wlong)
- call reada(weightcard,'WELEC',welec,1.0D0)
- call reada(weightcard,'WVDWPP',wvdwpp,welec)
- call reada(weightcard,'WEL_LOC',wel_loc,1.0D0)
- call reada(weightcard,'WCORR4',wcorr4,0.0D0)
- call reada(weightcard,'WCORR5',wcorr5,0.0D0)
- call reada(weightcard,'WCORR6',wcorr6,0.0D0)
- call reada(weightcard,'WTURN3',wturn3,1.0D0)
- call reada(weightcard,'WTURN4',wturn4,1.0D0)
- call reada(weightcard,'WTURN6',wturn6,1.0D0)
- call reada(weightcard,'WSCCOR',wsccor,1.0D0)
- call reada(weightcard,'WSTRAIN',wstrain,1.0D0)
- call reada(weightcard,'WBOND',wbond,1.0D0)
- call reada(weightcard,'WTOR',wtor,1.0D0)
- call reada(weightcard,'WTORD',wtor_d,1.0D0)
- call reada(weightcard,'WANG',wang,1.0D0)
- call reada(weightcard,'WSCLOC',wscloc,1.0D0)
-C Juyong
- call reada(weightcard,'WDFAD',wdfa_dist,0.0d0)
- call reada(weightcard,'WDFAT',wdfa_tor,0.0d0)
- call reada(weightcard,'WDFAN',wdfa_nei,0.0d0)
- call reada(weightcard,'WDFAB',wdfa_beta,0.0d0)
-C
- call reada(weightcard,'SCAL14',scal14,0.4D0)
- call reada(weightcard,'SCALSCP',scalscp,1.0d0)
- call reada(weightcard,'CUTOFF',cutoff_corr,7.0d0)
- call reada(weightcard,'DELT_CORR',delt_corr,0.5d0)
- call reada(weightcard,'TEMP0',temp0,300.0d0)
- if (index(weightcard,'SOFT').gt.0) ipot=6
-C 12/1/95 Added weight for the multi-body term WCORR
- call reada(weightcard,'WCORRH',wcorr,1.0D0)
- if (wcorr4.gt.0.0d0) wcorr=wcorr4
- weights(1)=wsc
- weights(2)=wscp
- weights(3)=welec
- weights(4)=wcorr
- weights(5)=wcorr5
- weights(6)=wcorr6
- weights(7)=wel_loc
- weights(8)=wturn3
- weights(9)=wturn4
- weights(10)=wturn6
- weights(11)=wang
- weights(12)=wscloc
- weights(13)=wtor
- weights(14)=wtor_d
- weights(15)=wstrain
- weights(16)=wvdwpp
- weights(17)=wbond
- weights(18)=scal14
- weights(21)=wsccor
-C JUYONG
- weights(24)=wdfa_dist
- weights(25)=wdfa_tor
- weights(26)=wdfa_nei
- weights(27)=wdfa_beta
-C
-
- if(me.eq.king.or..not.out1file)
- & write (iout,10) wsc,wscp,welec,wvdwpp,wbond,wang,wscloc,wtor,
- & wtor_d,wstrain,wel_loc,wcorr,wcorr5,wcorr6,wsccor,wturn3,
- & wturn4,wturn6,
- & wdfa_dist,wdfa_tor,wdfa_nei,wdfa_beta
-
- 10 format (/'Energy-term weights (unscaled):'//
- & 'WSCC= ',f10.6,' (SC-SC)'/
- & 'WSCP= ',f10.6,' (SC-p)'/
- & 'WELEC= ',f10.6,' (p-p electr)'/
- & 'WVDWPP= ',f10.6,' (p-p VDW)'/
- & 'WBOND= ',f10.6,' (stretching)'/
- & 'WANG= ',f10.6,' (bending)'/
- & 'WSCLOC= ',f10.6,' (SC local)'/
- & 'WTOR= ',f10.6,' (torsional)'/
- & 'WTORD= ',f10.6,' (double torsional)'/
- & 'WSTRAIN=',f10.6,' (SS bridges & dist. cnstr.)'/
- & 'WEL_LOC=',f10.6,' (multi-body 3-rd order)'/
- & 'WCORR4= ',f10.6,' (multi-body 4th order)'/
- & 'WCORR5= ',f10.6,' (multi-body 5th order)'/
- & 'WCORR6= ',f10.6,' (multi-body 6th order)'/
- & 'WSCCOR= ',f10.6,' (back-scloc correlation)'/
- & 'WTURN3= ',f10.6,' (turns, 3rd order)'/
- & 'WTURN4= ',f10.6,' (turns, 4th order)'/
- & 'WTURN6= ',f10.6,' (turns, 6th order)'/
- & 'WDFA_D= ',f10.6,' (DFA, distance)' /
- & 'WDFA_T= ',f10.6,' (DFA, torsional)' /
- & 'WDFA_N= ',f10.6,' (DFA, number of neighbor)' /
- & 'WDFA_B= ',f10.6,' (DFA, beta formation)')
-
- if(me.eq.king.or..not.out1file)then
- if (wcorr4.gt.0.0d0) then
- write (iout,'(/2a/)') 'Local-electrostatic type correlation ',
- & 'between contact pairs of peptide groups'
- write (iout,'(2(a,f5.3/))')
- & 'Cutoff on 4-6th order correlation terms: ',cutoff_corr,
- & 'Range of quenching the correlation terms:',2*delt_corr
- else if (wcorr.gt.0.0d0) then
- write (iout,'(/2a/)') 'Hydrogen-bonding correlation ',
- & 'between contact pairs of peptide groups'
- endif
- write (iout,'(a,f8.3)')
- & 'Scaling factor of 1,4 SC-p interactions:',scal14
- write (iout,'(a,f8.3)')
- & 'General scaling factor of SC-p interactions:',scalscp
- endif
- r0_corr=cutoff_corr-delt_corr
- do i=1,20
- aad(i,1)=scalscp*aad(i,1)
- aad(i,2)=scalscp*aad(i,2)
- bad(i,1)=scalscp*bad(i,1)
- bad(i,2)=scalscp*bad(i,2)
- enddo
-c call rescale_weights(t_bath)
- if(me.eq.king.or..not.out1file)
- & write (iout,22) wsc,wscp,welec,wvdwpp,wbond,wang,wscloc,wtor,
- & wtor_d,wstrain,wel_loc,wcorr,wcorr5,wcorr6,wsccor,wturn3,
- & wturn4,wturn6,
- & wdfa_dist,wdfa_tor,wdfa_nei,wdfa_beta
-
- 22 format (/'Energy-term weights (scaled):'//
- & 'WSCC= ',f10.6,' (SC-SC)'/
- & 'WSCP= ',f10.6,' (SC-p)'/
- & 'WELEC= ',f10.6,' (p-p electr)'/
- & 'WVDWPP= ',f10.6,' (p-p VDW)'/
- & 'WBOND= ',f10.6,' (stretching)'/
- & 'WANG= ',f10.6,' (bending)'/
- & 'WSCLOC= ',f10.6,' (SC local)'/
- & 'WTOR= ',f10.6,' (torsional)'/
- & 'WTORD= ',f10.6,' (double torsional)'/
- & 'WSTRAIN=',f10.6,' (SS bridges & dist. cnstr.)'/
- & 'WEL_LOC=',f10.6,' (multi-body 3-rd order)'/
- & 'WCORR4= ',f10.6,' (multi-body 4th order)'/
- & 'WCORR5= ',f10.6,' (multi-body 5th order)'/
- & 'WCORR6= ',f10.6,' (multi-body 6th order)'/
- & 'WSCCOR= ',f10.6,' (back-scloc correlatkion)'/
- & 'WTURN3= ',f10.6,' (turns, 3rd order)'/
- & 'WTURN4= ',f10.6,' (turns, 4th order)'/
- & 'WTURN6= ',f10.6,' (turns, 6th order)'/
- & 'WDFA_D= ',f10.6,' (DFA, distance)' /
- & 'WDFA_T= ',f10.6,' (DFA, torsional)' /
- & 'WDFA_N= ',f10.6,' (DFA, number of neighbor)' /
- & 'WDFA_B= ',f10.6,' (DFA, beta formation)')
-
- if(me.eq.king.or..not.out1file)
- & write (iout,*) "Reference temperature for weights calculation:",
- & temp0
- call reada(weightcard,"D0CM",d0cm,3.78d0)
- call reada(weightcard,"AKCM",akcm,15.1d0)
- call reada(weightcard,"AKTH",akth,11.0d0)
- call reada(weightcard,"AKCT",akct,12.0d0)
- call reada(weightcard,"V1SS",v1ss,-1.08d0)
- call reada(weightcard,"V2SS",v2ss,7.61d0)
- call reada(weightcard,"V3SS",v3ss,13.7d0)
- call reada(weightcard,"EBR",ebr,-5.50D0)
- if(me.eq.king.or..not.out1file) then
- write (iout,*) "Parameters of the SS-bond potential:"
- write (iout,*) "D0CM",d0cm," AKCM",akcm," AKTH",akth,
- & " AKCT",akct
- write (iout,*) "V1SS",v1ss," V2SS",v2ss," V3SS",v3ss
- write (iout,*) "EBR",ebr
- print *,'indpdb=',indpdb,' pdbref=',pdbref
- endif
- if (indpdb.gt.0 .or. pdbref) then
- read(inp,'(a)') pdbfile
- if(me.eq.king.or..not.out1file)
- & write (iout,'(2a)') 'PDB data will be read from file ',
- & pdbfile(:ilen(pdbfile))
- open(ipdbin,file=pdbfile,status='old',err=33)
- goto 34
- 33 write (iout,'(a)') 'Error opening PDB file.'
- stop
- 34 continue
-c print *,'Begin reading pdb data'
- call readpdb
-c print *,'Finished reading pdb data'
- if(me.eq.king.or..not.out1file)
- & write (iout,'(a,i3,a,i3)')'nsup=',nsup,
- & ' nstart_sup=',nstart_sup
- do i=1,nres
- itype_pdb(i)=itype(i)
- enddo
- close (ipdbin)
- nnt=nstart_sup
- nct=nstart_sup+nsup-1
- call contact(.false.,ncont_ref,icont_ref,co)
-
- if (sideadd) then
-C Following 2 lines for diagnostics; comment out if not needed
- write (iout,*) "Before sideadd"
- call intout
- if(me.eq.king.or..not.out1file)
- & write(iout,*)'Adding sidechains'
- maxsi=1000
- do i=2,nres-1
- iti=itype(i)
- if (iti.ne.10) then
- nsi=0
- fail=.true.
- do while (fail.and.nsi.le.maxsi)
-c call gen_side(iti,theta(i+1),alph(i),omeg(i),fail)
- nsi=nsi+1
- enddo
- if(fail) write(iout,*)'Adding sidechain failed for res ',
- & i,' after ',nsi,' trials'
- endif
- enddo
- endif
-C 10/03/12 Adam: Recalculate coordinates with new side chain positions
- call chainbuild
-C Following 2 lines for diagnostics; comment out if not needed
- write (iout,*) "After sideadd"
- call intout
- endif
-
- if (indpdb.eq.0) then
-C Read sequence if not taken from the pdb file.
- read (inp,*) nres
-c print *,'nres=',nres
- if (iscode.gt.0) then
- read (inp,'(80a1)') (sequence(i)(1:1),i=1,nres)
- else
- read (inp,'(20(1x,a3))') (sequence(i),i=1,nres)
- endif
-C Convert sequence to numeric code
- do i=1,nres
- itype(i)=rescode(i,sequence(i),iscode)
- enddo
-C Assign initial virtual bond lengths
- do i=2,nres
- vbld(i)=vbl
- vbld_inv(i)=vblinv
- enddo
- do i=2,nres-1
- vbld(i+nres)=dsc(iabs(itype(i)))
- vbld_inv(i+nres)=dsc_inv(iabs(itype(i)))
-c write (iout,*) "i",i," itype",itype(i),
-c & " dsc",dsc(itype(i))," vbld",vbld(i),vbld(i+nres)
- enddo
- endif
-c print *,nres
-c print '(20i4)',(itype(i),i=1,nres)
- do i=1,nres
-#ifdef PROCOR
- if (itype(i).eq.21 .or. itype(i+1).eq.21) then
-#else
- if (itype(i).eq.21) then
-#endif
- itel(i)=0
-#ifdef PROCOR
- else if (iabs(itype(i+1)).ne.20) then
-#else
- else if (iabs(itype(i)).ne.20) then
-#endif
- itel(i)=1
- else
- itel(i)=2
- endif
- enddo
- if(me.eq.king.or..not.out1file)then
- write (iout,*) "ITEL"
- do i=1,nres-1
- write (iout,*) i,itype(i),itel(i)
- enddo
- print *,'Call Read_Bridge.'
- endif
- call read_bridge
-C 8/13/98 Set limits to generating the dihedral angles
- do i=1,nres
- phibound(1,i)=-pi
- phibound(2,i)=pi
- enddo
- read (inp,*) ndih_constr
- if (ndih_constr.gt.0) then
- read (inp,*) ftors
- read (inp,*) (idih_constr(i),phi0(i),drange(i),i=1,ndih_constr)
- if(me.eq.king.or..not.out1file)then
- write (iout,*)
- & 'There are',ndih_constr,' constraints on phi angles.'
- do i=1,ndih_constr
- write (iout,'(i5,2f8.3)') idih_constr(i),phi0(i),drange(i)
- enddo
- endif
- do i=1,ndih_constr
- phi0(i)=deg2rad*phi0(i)
- drange(i)=deg2rad*drange(i)
- enddo
- if(me.eq.king.or..not.out1file)
- & write (iout,*) 'FTORS',ftors
- do i=1,ndih_constr
- ii = idih_constr(i)
- phibound(1,ii) = phi0(i)-drange(i)
- phibound(2,ii) = phi0(i)+drange(i)
- enddo
- endif
- nnt=1
-#ifdef MPI
- if (me.eq.king) then
-#endif
- write (iout,'(a)') 'Boundaries in phi angle sampling:'
- do i=1,nres
- write (iout,'(a3,i5,2f10.1)')
- & restyp(itype(i)),i,phibound(1,i)*rad2deg,phibound(2,i)*rad2deg
- enddo
-#ifdef MP
- endif
-#endif
- nct=nres
-cd print *,'NNT=',NNT,' NCT=',NCT
- if (itype(1).eq.21) nnt=2
- if (itype(nres).eq.21) nct=nct-1
-
-C Juyong:READ init_vars
-C Initialize variables!
-C Juyong:READ read_info
-C READ fragment information!!
-C both routines should be in dfa.F file!!
-
- if (.not. (wdfa_dist.eq.0.0 .and. wdfa_tor.eq.0.0 .and.
- & wdfa_nei.eq.0.0 .and. wdfa_beta.eq.0.0)) then
- call init_dfa_vars
- print*, 'init_dfa_vars finished!'
- call read_dfa_info
- print*, 'read_dfa_info finished!'
- endif
-C
-C
-
-
- if (pdbref) then
- if(me.eq.king.or..not.out1file)
- & write (iout,'(a,i3)') 'nsup=',nsup
- nstart_seq=nnt
- if (nsup.le.(nct-nnt+1)) then
- do i=0,nct-nnt+1-nsup
- if (seq_comp(itype(nnt+i),itype_pdb(nstart_sup),nsup)) then
- nstart_seq=nnt+i
- goto 111
- endif
- enddo
- write (iout,'(a)')
- & 'Error - sequences to be superposed do not match.'
- stop
- else
- do i=0,nsup-(nct-nnt+1)
- if (seq_comp(itype(nnt),itype_pdb(nstart_sup+i),nct-nnt+1))
- & then
- nstart_sup=nstart_sup+i
- nsup=nct-nnt+1
- goto 111
- endif
- enddo
- write (iout,'(a)')
- & 'Error - sequences to be superposed do not match.'
- endif
- 111 continue
- if (nsup.eq.0) nsup=nct-nnt
- if (nstart_sup.eq.0) nstart_sup=nnt
- if (nstart_seq.eq.0) nstart_seq=nnt
- if(me.eq.king.or..not.out1file)
- & write (iout,*) 'nsup=',nsup,' nstart_sup=',nstart_sup,
- & ' nstart_seq=',nstart_seq
- endif
-c--- Zscore rms -------
- if (nz_start.eq.0) nz_start=nnt
- if (nz_end.eq.0 .and. nsup.gt.0) then
- nz_end=nnt+nsup-1
- else if (nz_end.eq.0) then
- nz_end=nct
- endif
- if(me.eq.king.or..not.out1file)then
- write (iout,*) 'NZ_START=',nz_start,' NZ_END=',nz_end
- write (iout,*) 'IZ_SC=',iz_sc
- endif
-c----------------------
- call init_int_table
- if (refstr) then
- if (.not.pdbref) then
- call read_angles(inp,*38)
- goto 39
- 38 write (iout,'(a)') 'Error reading reference structure.'
-#ifdef MPI
- call MPI_Finalize(MPI_COMM_WORLD,IERROR)
- stop 'Error reading reference structure'
-#endif
- 39 call chainbuild
- call setup_var
-czscore call geom_to_var(nvar,coord_exp_zs(1,1))
- nstart_sup=nnt
- nstart_seq=nnt
- nsup=nct-nnt+1
- do i=1,2*nres
- do j=1,3
- cref(j,i)=c(j,i)
- enddo
- enddo
- call contact(.true.,ncont_ref,icont_ref,co)
- endif
-c write (iout,*) "constr_dist",constr_dist,nstart_sup,nsup
- call flush(iout)
- if (constr_dist.gt.0) call read_dist_constr
-c write (iout,*) "After read_dist_constr nhpb",nhpb
- call hpb_partition
- if(me.eq.king.or..not.out1file)
- & write (iout,*) 'Contact order:',co
- if (pdbref) then
- if(me.eq.king.or..not.out1file)
- & write (2,*) 'Shifting contacts:',nstart_seq,nstart_sup
- do i=1,ncont_ref
- do j=1,2
- icont_ref(j,i)=icont_ref(j,i)+nstart_seq-nstart_sup
- enddo
- if(me.eq.king.or..not.out1file)
- & write (2,*) i,' ',restyp(itype(icont_ref(1,i))),' ',
- & icont_ref(1,i),' ',
- & restyp(itype(icont_ref(2,i))),' ',icont_ref(2,i)
- enddo
- endif
- endif
- if (indpdb.eq.0 .and. modecalc.ne.2 .and. modecalc.ne.4
- & .and. modecalc.ne.8 .and. modecalc.ne.9 .and.
- & modecalc.ne.10) then
-C If input structure hasn't been supplied from the PDB file read or generate
-C initial geometry.
- if (iranconf.eq.0 .and. .not. extconf) then
- if(me.eq.king.or..not.out1file .and.fg_rank.eq.0)
- & write (iout,'(a)') 'Initial geometry will be read in.'
- if (read_cart) then
- read(inp,'(8f10.5)',end=36,err=36)
- & ((c(l,k),l=1,3),k=1,nres),
- & ((c(l,k+nres),l=1,3),k=nnt,nct)
- call int_from_cart1(.false.)
- do i=1,nres-1
- do j=1,3
- dc(j,i)=c(j,i+1)-c(j,i)
- dc_norm(j,i)=dc_norm(j,i)*vbld_inv(i+1)
- enddo
- enddo
- do i=nnt,nct
- if (itype(i).ne.10) then
- do j=1,3
- dc(j,i+nres)=c(j,i+nres)-c(j,i)
- dc_norm(j,i+nres)=dc_norm(j,i+nres)*vbld_inv(i+nres)
- enddo
- endif
- enddo
- return
- else
- call read_angles(inp,*36)
- endif
- goto 37
- 36 write (iout,'(a)') 'Error reading angle file.'
-#ifdef MPI
- call mpi_finalize( MPI_COMM_WORLD,IERR )
-#endif
- stop 'Error reading angle file.'
- 37 continue
- else if (extconf) then
- if(me.eq.king.or..not.out1file .and. fg_rank.eq.0)
- & write (iout,'(a)') 'Extended chain initial geometry.'
- do i=3,nres
- theta(i)=90d0*deg2rad
- enddo
- do i=4,nres
- phi(i)=180d0*deg2rad
- enddo
- do i=2,nres-1
- alph(i)=110d0*deg2rad
- enddo
- do i=2,nres-1
- omeg(i)=-120d0*deg2rad
- if (itype(i).le.0) omeg(i)=-omeg(i)
- enddo
- else
- if(me.eq.king.or..not.out1file)
- & write (iout,'(a)') 'Random-generated initial geometry.'
-
-
-#ifdef MPI
- if (me.eq.king .or. fg_rank.eq.0 .and. (
- & modecalc.eq.12 .or. modecalc.eq.14) ) then
-#endif
- do itrial=1,100
- itmp=1
-c call gen_rand_conf(itmp,*30)
- goto 40
- 30 write (iout,*) 'Failed to generate random conformation',
- & ', itrial=',itrial
- write (*,*) 'Processor:',me,
- & ' Failed to generate random conformation',
- & ' itrial=',itrial
- call intout
-
-#ifdef AIX
- call flush_(iout)
-#else
- call flush(iout)
-#endif
- enddo
- write (iout,'(a,i3,a)') 'Processor:',me,
- & ' error in generating random conformation.'
- write (*,'(a,i3,a)') 'Processor:',me,
- & ' error in generating random conformation.'
- call flush(iout)
-#ifdef MPI
- call MPI_Abort(mpi_comm_world,error_msg,ierrcode)
- 40 continue
- endif
-#else
- do itrial=1,100
- itmp=1
-c call gen_rand_conf(itmp,*31)
- goto 40
- 31 write (iout,*) 'Failed to generate random conformation',
- & ', itrial=',itrial
- write (*,*) 'Failed to generate random conformation',
- & ', itrial=',itrial
- enddo
- write (iout,'(a,i3,a)') 'Processor:',me,
- & ' error in generating random conformation.'
- write (*,'(a,i3,a)') 'Processor:',me,
- & ' error in generating random conformation.'
- stop
- 40 continue
-#endif
- endif
- elseif (modecalc.eq.4) then
- read (inp,'(a)') intinname
- open (intin,file=intinname,status='old',err=333)
- if (me.eq.king .or. .not.out1file.and.fg_rank.eq.0)
- & write (iout,'(a)') 'intinname',intinname
- write (*,'(a)') 'Processor',myrank,' intinname',intinname
- goto 334
- 333 write (iout,'(2a)') 'Error opening angle file ',intinname
-#ifdef MPI
- call MPI_Finalize(MPI_COMM_WORLD,IERR)
-#endif
- stop 'Error opening angle file.'
- 334 continue
-
- endif
-C Generate distance constraints, if the PDB structure is to be regularized.
- call setup_var
- if (me.eq.king .or. .not. out1file)
- & call intout
- if (ns.gt.0 .and. (me.eq.king .or. .not.out1file) ) then
- write (iout,'(/a,i3,a)')
- & 'The chain contains',ns,' disulfide-bridging cysteines.'
- write (iout,'(20i4)') (iss(i),i=1,ns)
- write (iout,'(/a/)') 'Pre-formed links are:'
- do i=1,nss
- i1=ihpb(i)-nres
- i2=jhpb(i)-nres
- it1=itype(i1)
- it2=itype(i2)
- if (me.eq.king.or..not.out1file)
- & write (iout,'(2a,i3,3a,i3,a,3f10.3)')
- & restyp(it1),'(',i1,') -- ',restyp(it2),'(',i2,')',dhpb(i),
- & ebr,forcon(i)
- enddo
- write (iout,'(a)')
- endif
-c if (i2ndstr.gt.0) call secstrp2dihc
-c call geom_to_var(nvar,x)
-c call etotal(energia(0))
-c call enerprint(energia(0))
-c call briefout(0,etot)
-c stop
-cd write (iout,'(2(a,i3))') 'NNT',NNT,' NCT',NCT
-cd write (iout,'(a)') 'Variable list:'
-cd write (iout,'(i4,f10.5)') (i,rad2deg*x(i),i=1,nvar)
-#ifdef MPI
- if (me.eq.king .or. (fg_rank.eq.0 .and. .not.out1file))
- & write (iout,'(//80(1h*)/20x,a,i4,a/80(1h*)//)')
- & 'Processor',myrank,': end reading molecular data.'
-#endif
- return
- end
-c--------------------------------------------------------------------------
- logical function seq_comp(itypea,itypeb,length)
- implicit none
- integer length,itypea(length),itypeb(length)
- integer i
- do i=1,length
- if (itypea(i).ne.itypeb(i)) then
- seq_comp=.false.
- return
- endif
- enddo
- seq_comp=.true.
- return
- end
-c-----------------------------------------------------------------------------
- subroutine read_bridge
-C Read information about disulfide bridges.
- implicit real*8 (a-h,o-z)
- include 'DIMENSIONS'
-#ifdef MPI
- include 'mpif.h'
-#endif
- include 'COMMON.IOUNITS'
- include 'COMMON.GEO'
- include 'COMMON.VAR'
- include 'COMMON.INTERACT'
- include 'COMMON.LOCAL'
- include 'COMMON.NAMES'
- include 'COMMON.CHAIN'
- include 'COMMON.FFIELD'
- include 'COMMON.SBRIDGE'
- include 'COMMON.HEADER'
- include 'COMMON.CONTROL'
-c include 'COMMON.DBASE'
-c include 'COMMON.THREAD'
- include 'COMMON.TIME1'
- include 'COMMON.SETUP'
-C Read bridging residues.
- read (inp,*) ns,(iss(i),i=1,ns)
- print *,'ns=',ns
- if(me.eq.king.or..not.out1file)
- & write (iout,*) 'ns=',ns,' iss:',(iss(i),i=1,ns)
-C Check whether the specified bridging residues are cystines.
- do i=1,ns
- if (itype(iss(i)).ne.1) then
- if (me.eq.king.or..not.out1file) write (iout,'(2a,i3,a)')
- & 'Do you REALLY think that the residue ',restyp(iss(i)),i,
- & ' can form a disulfide bridge?!!!'
- write (*,'(2a,i3,a)')
- & 'Do you REALLY think that the residue ',restyp(iss(i)),i,
- & ' can form a disulfide bridge?!!!'
-#ifdef MPI
- call MPI_Finalize(MPI_COMM_WORLD,ierror)
- stop
-#endif
- endif
- enddo
-C Read preformed bridges.
- if (ns.gt.0) then
- read (inp,*) nss,(ihpb(i),jhpb(i),i=1,nss)
- write (iout,*) 'nss=',nss,' ihpb,jhpb: ',(ihpb(i),jhpb(i),i=1,nss)
- if (nss.gt.0) then
- nhpb=nss
-C Check if the residues involved in bridges are in the specified list of
-C bridging residues.
- do i=1,nss
- do j=1,i-1
- if (ihpb(i).eq.ihpb(j).or.ihpb(i).eq.jhpb(j)
- & .or.jhpb(i).eq.ihpb(j).or.jhpb(i).eq.jhpb(j)) then
- write (iout,'(a,i3,a)') 'Disulfide pair',i,
- & ' contains residues present in other pairs.'
- write (*,'(a,i3,a)') 'Disulfide pair',i,
- & ' contains residues present in other pairs.'
-#ifdef MPI
- call MPI_Finalize(MPI_COMM_WORLD,ierror)
- stop
-#endif
- endif
- enddo
- do j=1,ns
- if (ihpb(i).eq.iss(j)) goto 10
- enddo
- write (iout,'(a,i3,a)') 'Pair',i,' contains unknown cystine.'
- 10 continue
- do j=1,ns
- if (jhpb(i).eq.iss(j)) goto 20
- enddo
- write (iout,'(a,i3,a)') 'Pair',i,' contains unknown cystine.'
- 20 continue
- dhpb(i)=dbr
- forcon(i)=fbr
- enddo
- do i=1,nss
- ihpb(i)=ihpb(i)+nres
- jhpb(i)=jhpb(i)+nres
- enddo
- endif
- endif
- return
- end
-c----------------------------------------------------------------------------
- subroutine read_x(kanal,*)
- implicit real*8 (a-h,o-z)
- include 'DIMENSIONS'
- include 'COMMON.GEO'
- include 'COMMON.VAR'
- include 'COMMON.CHAIN'
- include 'COMMON.IOUNITS'
- include 'COMMON.CONTROL'
- include 'COMMON.LOCAL'
- include 'COMMON.INTERACT'
-c Read coordinates from input
-c
- read(kanal,'(8f10.5)',end=10,err=10)
- & ((c(l,k),l=1,3),k=1,nres),
- & ((c(l,k+nres),l=1,3),k=nnt,nct)
- do j=1,3
- c(j,nres+1)=c(j,1)
- c(j,2*nres)=c(j,nres)
- enddo
- call int_from_cart1(.false.)
- do i=1,nres-1
- do j=1,3
- dc(j,i)=c(j,i+1)-c(j,i)
- dc_norm(j,i)=dc(j,i)*vbld_inv(i+1)
- enddo
- enddo
- do i=nnt,nct
- if (itype(i).ne.10) then
- do j=1,3
- dc(j,i+nres)=c(j,i+nres)-c(j,i)
- dc_norm(j,i+nres)=dc(j,i+nres)*vbld_inv(i+nres)
- enddo
- endif
- enddo
-
- return
- 10 return1
- end
-c------------------------------------------------------------------------------
- subroutine setup_var
- implicit real*8 (a-h,o-z)
- include 'DIMENSIONS'
- include 'COMMON.IOUNITS'
- include 'COMMON.GEO'
- include 'COMMON.VAR'
- include 'COMMON.INTERACT'
- include 'COMMON.LOCAL'
- include 'COMMON.NAMES'
- include 'COMMON.CHAIN'
- include 'COMMON.FFIELD'
- include 'COMMON.SBRIDGE'
- include 'COMMON.HEADER'
- include 'COMMON.CONTROL'
-c include 'COMMON.DBASE'
-c include 'COMMON.THREAD'
- include 'COMMON.TIME1'
-C Set up variable list.
- ntheta=nres-2
- nphi=nres-3
- nvar=ntheta+nphi
- nside=0
- do i=2,nres-1
- if (itype(i).ne.10) then
- nside=nside+1
- ialph(i,1)=nvar+nside
- ialph(nside,2)=i
- endif
- enddo
- if (indphi.gt.0) then
- nvar=nphi
- else if (indback.gt.0) then
- nvar=nphi+ntheta
- else
- nvar=nvar+2*nside
- endif
-cd write (iout,'(3i4)') (i,ialph(i,1),ialph(i,2),i=2,nres-1)
- return
- end
-c----------------------------------------------------------------------------
- subroutine gen_dist_constr
-C Generate CA distance constraints.
- implicit real*8 (a-h,o-z)
- include 'DIMENSIONS'
- include 'COMMON.IOUNITS'
- include 'COMMON.GEO'
- include 'COMMON.VAR'
- include 'COMMON.INTERACT'
- include 'COMMON.LOCAL'
- include 'COMMON.NAMES'
- include 'COMMON.CHAIN'
- include 'COMMON.FFIELD'
- include 'COMMON.SBRIDGE'
- include 'COMMON.HEADER'
- include 'COMMON.CONTROL'
-c include 'COMMON.DBASE'
-c include 'COMMON.THREAD'
- include 'COMMON.TIME1'
- dimension itype_pdb(maxres)
- common /pizda/ itype_pdb
- character*2 iden
-cd print *,'gen_dist_constr: nnt=',nnt,' nct=',nct
-cd write (2,*) 'gen_dist_constr: nnt=',nnt,' nct=',nct,
-cd & ' nstart_sup',nstart_sup,' nstart_seq',nstart_seq,
-cd & ' nsup',nsup
- do i=nstart_sup,nstart_sup+nsup-1
-cd write (2,*) 'i',i,' seq ',restyp(itype(i+nstart_seq-nstart_sup)),
-cd & ' seq_pdb', restyp(itype_pdb(i))
- do j=i+2,nstart_sup+nsup-1
- nhpb=nhpb+1
- ihpb(nhpb)=i+nstart_seq-nstart_sup
- jhpb(nhpb)=j+nstart_seq-nstart_sup
- forcon(nhpb)=weidis
- dhpb(nhpb)=dist(i,j)
- enddo
- enddo
-cd write (iout,'(a)') 'Distance constraints:'
-cd do i=nss+1,nhpb
-cd ii=ihpb(i)
-cd jj=jhpb(i)
-cd iden='CA'
-cd if (ii.gt.nres) then
-cd iden='SC'
-cd ii=ii-nres
-cd jj=jj-nres
-cd endif
-cd write (iout,'(a,1x,a,i4,3x,a,1x,a,i4,2f10.3)')
-cd & restyp(itype(ii)),iden,ii,restyp(itype(jj)),iden,jj,
-cd & dhpb(i),forcon(i)
-cd enddo
- return
- end
-c----------------------------------------------------------------------------
- subroutine csaread
- implicit real*8 (a-h,o-z)
- include 'DIMENSIONS'
- include 'COMMON.IOUNITS'
- include 'COMMON.GEO'
- include 'COMMON.CSA'
- include 'COMMON.BANK'
- include 'COMMON.CONTROL'
- character*80 ucase
- character*620 mcmcard
- call card_concat(mcmcard)
-
- call readi(mcmcard,'NCONF',nconf,50)
- call readi(mcmcard,'NADD',nadd,0)
- call readi(mcmcard,'JSTART',jstart,1)
- call readi(mcmcard,'JEND',jend,1)
- call readi(mcmcard,'NSTMAX',nstmax,500000)
- call readi(mcmcard,'N0',n0,1)
- call readi(mcmcard,'N1',n1,6)
- call readi(mcmcard,'N2',n2,4)
- call readi(mcmcard,'N3',n3,0)
- call readi(mcmcard,'N4',n4,0)
- call readi(mcmcard,'N5',n5,0)
- call readi(mcmcard,'N6',n6,10)
- call readi(mcmcard,'N7',n7,0)
- call readi(mcmcard,'N8',n8,0)
- call readi(mcmcard,'N9',n9,0)
- call readi(mcmcard,'N14',n14,0)
- call readi(mcmcard,'N15',n15,0)
- call readi(mcmcard,'N16',n16,0)
- call readi(mcmcard,'N17',n17,0)
- call readi(mcmcard,'N18',n18,0)
-
- vdisulf=(index(mcmcard,'DYNSS').gt.0)
-
- call readi(mcmcard,'NDIFF',ndiff,2)
- call reada(mcmcard,'DIFFCUT',diffcut,0.0d0)
- call readi(mcmcard,'IS1',is1,1)
- call readi(mcmcard,'IS2',is2,8)
- call readi(mcmcard,'NRAN0',nran0,4)
- call readi(mcmcard,'NRAN1',nran1,2)
- call readi(mcmcard,'IRR',irr,1)
- call readi(mcmcard,'NSEED',nseed,20)
- call readi(mcmcard,'NTOTAL',ntotal,10000)
- call reada(mcmcard,'CUT1',cut1,2.0d0)
- call reada(mcmcard,'CUT2',cut2,5.0d0)
- call reada(mcmcard,'ESTOP',estop,-300000.0d0)
- call readi(mcmcard,'ICMAX',icmax,1)
- call readi(mcmcard,'NBANKM',nbankm,400)
- call readi(mcmcard,'IUCUT',iucut,2)
- call readi(mcmcard,'IRESTART',irestart,0)
-c!bankt call readi(mcmcard,'NBANKTM',ntbankm,0)
- ntbankm=0
-c!bankt
- call reada(mcmcard,'DELE',dele,20.0d0)
- call reada(mcmcard,'DIFCUT',difcut,720.0d0)
- call readi(mcmcard,'IREF',iref,0)
- call reada(mcmcard,'RMSCUT',rmscut,4.0d0)
- call reada(mcmcard,'PNCCUT',pnccut,0.5d0)
- call readi(mcmcard,'NCONF_IN',nconf_in,0)
- call reada(mcmcard,'RDIH_BIAS',rdih_bias,0.5d0)
- write (iout,*) "NCONF_IN",nconf_in
- tm_score=(index(mcmcard,'TMSCORE').gt.0)
- if (tm_score) write (iout,*) "Using TM_Score instead of DIFF",
- & " for torsional angles"
- return
- end
-
- subroutine read_minim
- implicit real*8 (a-h,o-z)
- include 'DIMENSIONS'
- include 'COMMON.MINIM'
- include 'COMMON.IOUNITS'
- character*80 ucase
- character*320 minimcard
- call card_concat(minimcard)
- call readi(minimcard,'MAXMIN',maxmin,2000)
- call readi(minimcard,'MAXFUN',maxfun,5000)
- call readi(minimcard,'MINMIN',minmin,maxmin)
- call readi(minimcard,'MINFUN',minfun,maxmin)
- call reada(minimcard,'TOLF',tolf,1.0D-2)
- call reada(minimcard,'RTOLF',rtolf,1.0D-4)
- print_min_stat=min0(index(minimcard,'PRINT_MIN_STAT'),1)
- print_min_res=min0(index(minimcard,'PRINT_MIN_RES'),1)
- print_min_ini=min0(index(minimcard,'PRINT_MIN_INI'),1)
- write (iout,'(/80(1h*)/20x,a/80(1h*))')
- & 'Options in energy minimization:'
- write (iout,'(4(a,i5),a,1pe14.5,a,1pe14.5)')
- & 'MaxMin:',MaxMin,' MaxFun:',MaxFun,
- & 'MinMin:',MinMin,' MinFun:',MinFun,
- & ' TolF:',TolF,' RTolF:',RTolF
- return
- end
-c----------------------------------------------------------------------------
- subroutine read_angles(kanal,*)
- implicit real*8 (a-h,o-z)
- include 'DIMENSIONS'
- include 'COMMON.GEO'
- include 'COMMON.VAR'
- include 'COMMON.CHAIN'
- include 'COMMON.IOUNITS'
- include 'COMMON.CONTROL'
-c Read angles from input
-c
- read (kanal,*,err=10,end=10) (theta(i),i=3,nres)
- read (kanal,*,err=10,end=10) (phi(i),i=4,nres)
- read (kanal,*,err=10,end=10) (alph(i),i=2,nres-1)
- read (kanal,*,err=10,end=10) (omeg(i),i=2,nres-1)
-
- do i=1,nres
-c 9/7/01 avoid 180 deg valence angle
- if (theta(i).gt.179.99d0) theta(i)=179.99d0
-c
- theta(i)=deg2rad*theta(i)
- phi(i)=deg2rad*phi(i)
- alph(i)=deg2rad*alph(i)
- omeg(i)=deg2rad*omeg(i)
- enddo
- return
- 10 return1
- end
-c----------------------------------------------------------------------------
- subroutine reada(rekord,lancuch,wartosc,default)
- implicit none
- character*(*) rekord,lancuch
- double precision wartosc,default
- integer ilen,iread
- external ilen
- iread=index(rekord,lancuch)
- if (iread.eq.0) then
- wartosc=default
- return
- endif
- iread=iread+ilen(lancuch)+1
- read (rekord(iread:),*,err=10,end=10) wartosc
- return
- 10 wartosc=default
- return
- end
-c----------------------------------------------------------------------------
- subroutine readi(rekord,lancuch,wartosc,default)
- implicit none
- character*(*) rekord,lancuch
- integer wartosc,default
- integer ilen,iread
- external ilen
- iread=index(rekord,lancuch)
- if (iread.eq.0) then
- wartosc=default
- return
- endif
- iread=iread+ilen(lancuch)+1
- read (rekord(iread:),*,err=10,end=10) wartosc
- return
- 10 wartosc=default
- return
- end
-c----------------------------------------------------------------------------
- subroutine multreadi(rekord,lancuch,tablica,dim,default)
- implicit none
- integer dim,i
- integer tablica(dim),default
- character*(*) rekord,lancuch
- character*80 aux
- integer ilen,iread
- external ilen
- do i=1,dim
- tablica(i)=default
- enddo
- iread=index(rekord,lancuch(:ilen(lancuch))//"=")
- if (iread.eq.0) return
- iread=iread+ilen(lancuch)+1
- read (rekord(iread:),*,end=10,err=10) (tablica(i),i=1,dim)
- 10 return
- end
-c----------------------------------------------------------------------------
- subroutine multreada(rekord,lancuch,tablica,dim,default)
- implicit none
- integer dim,i
- double precision tablica(dim),default
- character*(*) rekord,lancuch
- character*80 aux
- integer ilen,iread
- external ilen
- do i=1,dim
- tablica(i)=default
- enddo
- iread=index(rekord,lancuch(:ilen(lancuch))//"=")
- if (iread.eq.0) return
- iread=iread+ilen(lancuch)+1
- read (rekord(iread:),*,end=10,err=10) (tablica(i),i=1,dim)
- 10 return
- end
-c----------------------------------------------------------------------------
- subroutine openunits
- implicit real*8 (a-h,o-z)
- include 'DIMENSIONS'
-#ifdef MPI
- include 'mpif.h'
- character*16 form,nodename
- integer nodelen
-#endif
- include 'COMMON.SETUP'
- include 'COMMON.IOUNITS'
-c include 'COMMON.MD'
- include 'COMMON.CONTROL'
- integer lenpre,lenpot,ilen,lentmp
- external ilen
- character*3 out1file_text,ucase
- character*3 ll
- external ucase
-c print *,"Processor",myrank,"fg_rank",fg_rank," entered openunits"
- call getenv_loc("PREFIX",prefix)
- pref_orig = prefix
- call getenv_loc("POT",pot)
- call getenv_loc("DIRTMP",tmpdir)
- call getenv_loc("CURDIR",curdir)
- call getenv_loc("OUT1FILE",out1file_text)
-c print *,"Processor",myrank,"fg_rank",fg_rank," did GETENV"
- out1file_text=ucase(out1file_text)
- if (out1file_text(1:1).eq."Y") then
- out1file=.true.
- else
- out1file=fg_rank.gt.0
- endif
- lenpre=ilen(prefix)
- lenpot=ilen(pot)
- lentmp=ilen(tmpdir)
- if (lentmp.gt.0) then
- write (*,'(80(1h!))')
- write (*,'(a,19x,a,19x,a)') "!"," A T T E N T I O N ","!"
- write (*,'(80(1h!))')
- write (*,*)"All output files will be on node /tmp directory."
-#ifdef MPI
- call MPI_GET_PROCESSOR_NAME( nodename, nodelen, IERROR )
- if (me.eq.king) then
- write (*,*) "The master node is ",nodename
- else if (fg_rank.eq.0) then
- write (*,*) "I am the CG slave node ",nodename
- else
- write (*,*) "I am the FG slave node ",nodename
- endif
-#endif
- PREFIX = tmpdir(:lentmp)//'/'//prefix(:lenpre)
- lenpre = lentmp+lenpre+1
- endif
- entname=prefix(:lenpre)//'_'//pot(:lenpot)//'.entr'
-C Get the names and open the input files
-#if defined(WINIFL) || defined(WINPGI)
- open(1,file=pref_orig(:ilen(pref_orig))//
- & '.inp',status='old',readonly,shared)
- open (9,file=prefix(:ilen(prefix))//'.intin',status='unknown')
-C open (18,file=prefix(:ilen(prefix))//'.entin',status='unknown')
-C Get parameter filenames and open the parameter files.
- call getenv_loc('BONDPAR',bondname)
- open (ibond,file=bondname,status='old',readonly,shared)
- call getenv_loc('THETPAR',thetname)
- open (ithep,file=thetname,status='old',readonly,shared)
-#ifndef CRYST_THETA
- call getenv_loc('THETPARPDB',thetname_pdb)
- open (ithep_pdb,file=thetname_pdb,status='old',readonly,shared)
-#endif
- call getenv_loc('ROTPAR',rotname)
- open (irotam,file=rotname,status='old',readonly,shared)
-#ifndef CRYST_SC
- call getenv_loc('ROTPARPDB',rotname_pdb)
- open (irotam_pdb,file=rotname_pdb,status='old',readonly,shared)
-#endif
- call getenv_loc('TORPAR',torname)
- open (itorp,file=torname,status='old',readonly,shared)
- call getenv_loc('TORDPAR',tordname)
- open (itordp,file=tordname,status='old',readonly,shared)
- call getenv_loc('FOURIER',fouriername)
- open (ifourier,file=fouriername,status='old',readonly,shared)
- call getenv_loc('ELEPAR',elename)
- open (ielep,file=elename,status='old',readonly,shared)
- call getenv_loc('SIDEPAR',sidename)
- open (isidep,file=sidename,status='old',readonly,shared)
-#elif (defined CRAY) || (defined AIX)
- open(1,file=pref_orig(:ilen(pref_orig))//'.inp',status='old',
- & action='read')
-c print *,"Processor",myrank," opened file 1"
- open (9,file=prefix(:ilen(prefix))//'.intin',status='unknown')
-c print *,"Processor",myrank," opened file 9"
-C open (18,file=prefix(:ilen(prefix))//'.entin',status='unknown')
-C Get parameter filenames and open the parameter files.
- call getenv_loc('BONDPAR',bondname)
- open (ibond,file=bondname,status='old',action='read')
-c print *,"Processor",myrank," opened file IBOND"
- call getenv_loc('THETPAR',thetname)
- open (ithep,file=thetname,status='old',action='read')
-c print *,"Processor",myrank," opened file ITHEP"
-#ifndef CRYST_THETA
- call getenv_loc('THETPARPDB',thetname_pdb)
- open (ithep_pdb,file=thetname_pdb,status='old',action='read')
-#endif
- call getenv_loc('ROTPAR',rotname)
- open (irotam,file=rotname,status='old',action='read')
-c print *,"Processor",myrank," opened file IROTAM"
-#ifndef CRYST_SC
- call getenv_loc('ROTPARPDB',rotname_pdb)
- open (irotam_pdb,file=rotname_pdb,status='old',action='read')
-#endif
- call getenv_loc('TORPAR',torname)
- open (itorp,file=torname,status='old',action='read')
-c print *,"Processor",myrank," opened file ITORP"
- call getenv_loc('TORDPAR',tordname)
- open (itordp,file=tordname,status='old',action='read')
-c print *,"Processor",myrank," opened file ITORDP"
- call getenv_loc('SCCORPAR',sccorname)
- open (isccor,file=sccorname,status='old',action='read')
-c print *,"Processor",myrank," opened file ISCCOR"
- call getenv_loc('FOURIER',fouriername)
- open (ifourier,file=fouriername,status='old',action='read')
-c print *,"Processor",myrank," opened file IFOURIER"
- call getenv_loc('ELEPAR',elename)
- open (ielep,file=elename,status='old',action='read')
-c print *,"Processor",myrank," opened file IELEP"
- call getenv_loc('SIDEPAR',sidename)
- open (isidep,file=sidename,status='old',action='read')
-c print *,"Processor",myrank," opened file ISIDEP"
-c print *,"Processor",myrank," opened parameter files"
-#elif (defined G77)
- open(1,file=pref_orig(:ilen(pref_orig))//'.inp',status='old')
- open (9,file=prefix(:ilen(prefix))//'.intin',status='unknown')
-C open (18,file=prefix(:ilen(prefix))//'.entin',status='unknown')
-C Get parameter filenames and open the parameter files.
- call getenv_loc('BONDPAR',bondname)
- open (ibond,file=bondname,status='old')
- call getenv_loc('THETPAR',thetname)
- open (ithep,file=thetname,status='old')
-#ifndef CRYST_THETA
- call getenv_loc('THETPARPDB',thetname_pdb)
- open (ithep_pdb,file=thetname_pdb,status='old')
-#endif
- call getenv_loc('ROTPAR',rotname)
- open (irotam,file=rotname,status='old')
-#ifndef CRYST_SC
- call getenv_loc('ROTPARPDB',rotname_pdb)
- open (irotam_pdb,file=rotname_pdb,status='old')
-#endif
- call getenv_loc('TORPAR',torname)
- open (itorp,file=torname,status='old')
- call getenv_loc('TORDPAR',tordname)
- open (itordp,file=tordname,status='old')
- call getenv_loc('SCCORPAR',sccorname)
- open (isccor,file=sccorname,status='old')
- call getenv_loc('FOURIER',fouriername)
- open (ifourier,file=fouriername,status='old')
- call getenv_loc('ELEPAR',elename)
- open (ielep,file=elename,status='old')
- call getenv_loc('SIDEPAR',sidename)
- open (isidep,file=sidename,status='old')
-#else
- open(1,file=pref_orig(:ilen(pref_orig))//'.inp',status='old',
- & readonly)
- open (9,file=prefix(:ilen(prefix))//'.intin',status='unknown')
-C open (18,file=prefix(:ilen(prefix))//'.entin',status='unknown')
-C Get parameter filenames and open the parameter files.
- call getenv_loc('BONDPAR',bondname)
- open (ibond,file=bondname,status='old',readonly)
- call getenv_loc('THETPAR',thetname)
- open (ithep,file=thetname,status='old',readonly)
-#ifndef CRYST_THETA
- call getenv_loc('THETPARPDB',thetname_pdb)
- print *,"thetname_pdb ",thetname_pdb
- open (ithep_pdb,file=thetname_pdb,status='old',readonly)
- print *,ithep_pdb," opened"
-#endif
- call getenv_loc('ROTPAR',rotname)
- open (irotam,file=rotname,status='old',readonly)
-#ifndef CRYST_SC
- call getenv_loc('ROTPARPDB',rotname_pdb)
- open (irotam_pdb,file=rotname_pdb,status='old',readonly)
-#endif
- call getenv_loc('TORPAR',torname)
- open (itorp,file=torname,status='old',readonly)
- call getenv_loc('TORDPAR',tordname)
- open (itordp,file=tordname,status='old',readonly)
- call getenv_loc('SCCORPAR',sccorname)
- open (isccor,file=sccorname,status='old',readonly)
- call getenv_loc('FOURIER',fouriername)
- open (ifourier,file=fouriername,status='old',readonly)
- call getenv_loc('ELEPAR',elename)
- open (ielep,file=elename,status='old',readonly)
- call getenv_loc('SIDEPAR',sidename)
- open (isidep,file=sidename,status='old',readonly)
-#endif
-#ifndef OLDSCP
-C
-C 8/9/01 In the newest version SCp interaction constants are read from a file
-C Use -DOLDSCP to use hard-coded constants instead.
-C
- call getenv_loc('SCPPAR',scpname)
-#if defined(WINIFL) || defined(WINPGI)
- open (iscpp,file=scpname,status='old',readonly,shared)
-#elif (defined CRAY) || (defined AIX)
- open (iscpp,file=scpname,status='old',action='read')
-#elif (defined G77)
- open (iscpp,file=scpname,status='old')
-#else
- open (iscpp,file=scpname,status='old',readonly)
-#endif
-#endif
- call getenv_loc('PATTERN',patname)
-#if defined(WINIFL) || defined(WINPGI)
- open (icbase,file=patname,status='old',readonly,shared)
-#elif (defined CRAY) || (defined AIX)
- open (icbase,file=patname,status='old',action='read')
-#elif (defined G77)
- open (icbase,file=patname,status='old')
-#else
- open (icbase,file=patname,status='old',readonly)
-#endif
-#ifdef MPI
-C Open output file only for CG processes
-c print *,"Processor",myrank," fg_rank",fg_rank
- if (fg_rank.eq.0) then
-
- if (nodes.eq.1) then
- npos=3
- else
- npos = dlog10(dfloat(nodes-1))+1
- endif
- if (npos.lt.3) npos=3
- write (liczba,'(i1)') npos
- form = '(bz,i'//liczba(:ilen(liczba))//'.'//liczba(:ilen(liczba))
- & //')'
- write (liczba,form) me
- outname=prefix(:lenpre)//'.out_'//pot(:lenpot)//
- & liczba(:ilen(liczba))
- intname=prefix(:lenpre)//'_'//pot(:lenpot)//liczba(:ilen(liczba))
- & //'.int'
- pdbname=prefix(:lenpre)//'_'//pot(:lenpot)//liczba(:ilen(liczba))
- & //'.pdb'
- mol2name=prefix(:lenpre)//'_'//pot(:lenpot)//
- & liczba(:ilen(liczba))//'.mol2'
- statname=prefix(:lenpre)//'_'//pot(:lenpot)//
- & liczba(:ilen(liczba))//'.stat'
- if (lentmp.gt.0)
- & call copy_to_tmp(pref_orig(:ilen(pref_orig))//'_'//pot(:lenpot)
- & //liczba(:ilen(liczba))//'.stat')
- rest2name=prefix(:ilen(prefix))//"_"//liczba(:ilen(liczba))
- & //'.rst'
-c if(usampl) then
-c qname=prefix(:lenpre)//'_'//pot(:lenpot)//
-c & liczba(:ilen(liczba))//'.const'
-c endif
-
- endif
-#else
- outname=prefix(:lenpre)//'.out_'//pot(:lenpot)
- intname=prefix(:lenpre)//'_'//pot(:lenpot)//'.int'
- pdbname=prefix(:lenpre)//'_'//pot(:lenpot)//'.pdb'
- mol2name=prefix(:lenpre)//'_'//pot(:lenpot)//'.mol2'
- statname=prefix(:lenpre)//'_'//pot(:lenpot)//'.stat'
- if (lentmp.gt.0)
- & call copy_to_tmp(pref_orig(:ilen(pref_orig))//'_'//pot(:lenpot)//
- & '.stat')
- rest2name=prefix(:ilen(prefix))//'.rst'
-c if(usampl) then
-c qname=prefix(:lenpre)//'_'//pot(:lenpot)//'.const'
-c endif
-#endif
-#if defined(AIX) || defined(PGI)
- if (me.eq.king .or. .not. out1file)
- & open(iout,file=outname,status='unknown')
-#ifdef DEBUG
- if (fg_rank.gt.0) then
- write (liczba,'(i3.3)') myrank/nfgtasks
- write (ll,'(bz,i3.3)') fg_rank
- open(iout,file="debug"//liczba(:ilen(liczba))//"."//ll,
- & status='unknown')
- endif
-#endif
- if(me.eq.king) then
- open(igeom,file=intname,status='unknown',position='append')
- open(ipdb,file=pdbname,status='unknown')
- open(imol2,file=mol2name,status='unknown')
- open(istat,file=statname,status='unknown',position='append')
- else
-c1out open(iout,file=outname,status='unknown')
- endif
-#else
- if (me.eq.king .or. .not.out1file)
- & open(iout,file=outname,status='unknown')
-#ifdef DEBUG
- if (fg_rank.gt.0) then
- write (liczba,'(i3.3)') myrank/nfgtasks
- write (ll,'(bz,i3.3)') fg_rank
- open(iout,file="debug"//liczba(:ilen(liczba))//"."//ll,
- & status='unknown')
- endif
-#endif
- if(me.eq.king) then
- open(igeom,file=intname,status='unknown',access='append')
- open(ipdb,file=pdbname,status='unknown')
- open(imol2,file=mol2name,status='unknown')
- open(istat,file=statname,status='unknown',access='append')
- else
-c1out open(iout,file=outname,status='unknown')
- endif
-#endif
- csa_rbank=prefix(:lenpre)//'.CSA.rbank'
- csa_seed=prefix(:lenpre)//'.CSA.seed'
- csa_history=prefix(:lenpre)//'.CSA.history'
- csa_bank=prefix(:lenpre)//'.CSA.bank'
- csa_bank1=prefix(:lenpre)//'.CSA.bank1'
- csa_alpha=prefix(:lenpre)//'.CSA.alpha'
- csa_alpha1=prefix(:lenpre)//'.CSA.alpha1'
-c!bankt csa_bankt=prefix(:lenpre)//'.CSA.bankt'
- csa_int=prefix(:lenpre)//'.int'
- csa_bank_reminimized=prefix(:lenpre)//'.CSA.bank_reminimized'
- csa_native_int=prefix(:lenpre)//'.CSA.native.int'
- csa_in=prefix(:lenpre)//'.CSA.in'
-c print *,"Processor",myrank,"fg_rank",fg_rank," opened files"
-C Write file names
- if (me.eq.king)then
- write (iout,'(80(1h-))')
- write (iout,'(30x,a)') "FILE ASSIGNMENT"
- write (iout,'(80(1h-))')
- write (iout,*) "Input file : ",
- & pref_orig(:ilen(pref_orig))//'.inp'
- write (iout,*) "Output file : ",
- & outname(:ilen(outname))
- write (iout,*)
- write (iout,*) "Sidechain potential file : ",
- & sidename(:ilen(sidename))
-#ifndef OLDSCP
- write (iout,*) "SCp potential file : ",
- & scpname(:ilen(scpname))
-#endif
- write (iout,*) "Electrostatic potential file : ",
- & elename(:ilen(elename))
- write (iout,*) "Cumulant coefficient file : ",
- & fouriername(:ilen(fouriername))
- write (iout,*) "Torsional parameter file : ",
- & torname(:ilen(torname))
- write (iout,*) "Double torsional parameter file : ",
- & tordname(:ilen(tordname))
- write (iout,*) "SCCOR parameter file : ",
- & sccorname(:ilen(sccorname))
- write (iout,*) "Bond & inertia constant file : ",
- & bondname(:ilen(bondname))
- write (iout,*) "Bending parameter file : ",
- & thetname(:ilen(thetname))
- write (iout,*) "Rotamer parameter file : ",
- & rotname(:ilen(rotname))
- write (iout,*) "Threading database : ",
- & patname(:ilen(patname))
- if (lentmp.ne.0)
- &write (iout,*)" DIRTMP : ",
- & tmpdir(:lentmp)
- write (iout,'(80(1h-))')
- endif
- return
- end
-c----------------------------------------------------------------------------
- subroutine card_concat(card)
- implicit real*8 (a-h,o-z)
- include 'DIMENSIONS'
- include 'COMMON.IOUNITS'
- character*(*) card
- character*80 karta,ucase
- external ilen
- read (inp,'(a)') karta
- karta=ucase(karta)
- card=' '
- do while (karta(80:80).eq.'&')
- card=card(:ilen(card)+1)//karta(:79)
- read (inp,'(a)') karta
- karta=ucase(karta)
- enddo
- card=card(:ilen(card)+1)//karta
- return
- end
-
- subroutine read_dist_constr
- implicit real*8 (a-h,o-z)
- include 'DIMENSIONS'
-#ifdef MPI
- include 'mpif.h'
-#endif
- include 'COMMON.SETUP'
- include 'COMMON.CONTROL'
- include 'COMMON.CHAIN'
- include 'COMMON.IOUNITS'
- include 'COMMON.SBRIDGE'
- integer ifrag_(2,100),ipair_(2,100)
- double precision wfrag_(100),wpair_(100)
- character*500 controlcard
- write (iout,*) "Calling read_dist_constr"
- write (iout,*) "nres",nres," nstart_sup",nstart_sup," nsup",nsup
- call flush(iout)
- call card_concat(controlcard)
- call readi(controlcard,"NFRAG",nfrag_,0)
- call readi(controlcard,"NPAIR",npair_,0)
- call readi(controlcard,"NDIST",ndist_,0)
- call reada(controlcard,'DIST_CUT',dist_cut,5.0d0)
- call multreadi(controlcard,"IFRAG",ifrag_(1,1),2*nfrag_,0)
- call multreadi(controlcard,"IPAIR",ipair_(1,1),2*npair_,0)
- call multreada(controlcard,"WFRAG",wfrag_(1),nfrag_,0.0d0)
- call multreada(controlcard,"WPAIR",wpair_(1),npair_,0.0d0)
-c write (iout,*) "NFRAG",nfrag_," NPAIR",npair_," NDIST",ndist_
-c write (iout,*) "IFRAG"
-c do i=1,nfrag_
-c write (iout,*) i,ifrag_(1,i),ifrag_(2,i),wfrag_(i)
-c enddo
-c write (iout,*) "IPAIR"
-c do i=1,npair_
-c write (iout,*) i,ipair_(1,i),ipair_(2,i),wpair_(i)
-c enddo
- call flush(iout)
- do i=1,nfrag_
- if (ifrag_(1,i).lt.nstart_sup) ifrag_(1,i)=nstart_sup
- if (ifrag_(2,i).gt.nstart_sup+nsup-1)
- & ifrag_(2,i)=nstart_sup+nsup-1
-c write (iout,*) i,ifrag_(1,i),ifrag_(2,i),wfrag_(i)
- call flush(iout)
- if (wfrag_(i).gt.0.0d0) then
- do j=ifrag_(1,i),ifrag_(2,i)-1
- do k=j+1,ifrag_(2,i)
- write (iout,*) "j",j," k",k
- ddjk=dist(j,k)
- if (constr_dist.eq.1) then
- nhpb=nhpb+1
- ihpb(nhpb)=j
- jhpb(nhpb)=k
- dhpb(nhpb)=ddjk
- forcon(nhpb)=wfrag_(i)
- else if (constr_dist.eq.2) then
- if (ddjk.le.dist_cut) then
- nhpb=nhpb+1
- ihpb(nhpb)=j
- jhpb(nhpb)=k
- dhpb(nhpb)=ddjk
- forcon(nhpb)=wfrag_(i)
- endif
- else
- nhpb=nhpb+1
- ihpb(nhpb)=j
- jhpb(nhpb)=k
- dhpb(nhpb)=ddjk
- forcon(nhpb)=wfrag_(i)*dexp(-0.5d0*(ddjk/dist_cut)**2)
- endif
-#ifdef MPI
- if (.not.out1file .or. me.eq.king)
- & write (iout,'(a,3i5,f8.2,1pe12.2)') "+dist.constr ",
- & nhpb,ihpb(nhpb),jhpb(nhpb),dhpb(nhpb),forcon(nhpb)
-#else
- write (iout,'(a,3i5,f8.2,1pe12.2)') "+dist.constr ",
- & nhpb,ihpb(nhpb),jhpb(nhpb),dhpb(nhpb),forcon(nhpb)
-#endif
- enddo
- enddo
- endif
- enddo
- do i=1,npair_
- if (wpair_(i).gt.0.0d0) then
- ii = ipair_(1,i)
- jj = ipair_(2,i)
- if (ii.gt.jj) then
- itemp=ii
- ii=jj
- jj=itemp
- endif
- do j=ifrag_(1,ii),ifrag_(2,ii)
- do k=ifrag_(1,jj),ifrag_(2,jj)
- nhpb=nhpb+1
- ihpb(nhpb)=j
- jhpb(nhpb)=k
- forcon(nhpb)=wpair_(i)
- dhpb(nhpb)=dist(j,k)
-#ifdef MPI
- if (.not.out1file .or. me.eq.king)
- & write (iout,'(a,3i5,f8.2,f10.1)') "+dist.constr ",
- & nhpb,ihpb(nhpb),jhpb(nhpb),dhpb(nhpb),forcon(nhpb)
-#else
- write (iout,'(a,3i5,f8.2,f10.1)') "+dist.constr ",
- & nhpb,ihpb(nhpb),jhpb(nhpb),dhpb(nhpb),forcon(nhpb)
-#endif
- enddo
- enddo
- endif
- enddo
- do i=1,ndist_
- read (inp,*) ihpb(nhpb+1),jhpb(nhpb+1),forcon(nhpb+1)
- if (forcon(nhpb+1).gt.0.0d0) then
- nhpb=nhpb+1
- dhpb(nhpb)=dist(ihpb(nhpb),jhpb(nhpb))
-#ifdef MPI
- if (.not.out1file .or. me.eq.king)
- & write (iout,'(a,3i5,f8.2,f10.1)') "+dist.constr ",
- & nhpb,ihpb(nhpb),jhpb(nhpb),dhpb(nhpb),forcon(nhpb)
-#else
- write (iout,'(a,3i5,f8.2,f10.1)') "+dist.constr ",
- & nhpb,ihpb(nhpb),jhpb(nhpb),dhpb(nhpb),forcon(nhpb)
-#endif
- endif
- enddo
- call flush(iout)
- return
- end
-c-------------------------------------------------------------------------------
-#ifdef WINIFL
- subroutine flush(iu)
- return
- end
-#endif
-#ifdef AIX
- subroutine flush(iu)
- call flush_(iu)
- return
- end
-#endif
-c------------------------------------------------------------------------------
- subroutine copy_to_tmp(source)
- include "DIMENSIONS"
- include "COMMON.IOUNITS"
- character*(*) source
- character* 256 tmpfile
- integer ilen
- external ilen
- logical ex
- tmpfile=curdir(:ilen(curdir))//"/"//source(:ilen(source))
- inquire(file=tmpfile,exist=ex)
- if (ex) then
- write (*,*) "Copying ",tmpfile(:ilen(tmpfile)),
- & " to temporary directory..."
- write (*,*) "/bin/cp "//tmpfile(:ilen(tmpfile))//" "//tmpdir
- call system("/bin/cp "//tmpfile(:ilen(tmpfile))//" "//tmpdir)
- endif
- return
- end
-c------------------------------------------------------------------------------
- subroutine move_from_tmp(source)
- include "DIMENSIONS"
- include "COMMON.IOUNITS"
- character*(*) source
- integer ilen
- external ilen
- write (*,*) "Moving ",source(:ilen(source)),
- & " from temporary directory to working directory"
- write (*,*) "/bin/mv "//source(:ilen(source))//" "//curdir
- call system("/bin/mv "//source(:ilen(source))//" "//curdir)
- return
- end
-c------------------------------------------------------------------------------
- subroutine random_init(seed)
-C
-C Initialize random number generator
-C
- implicit real*8 (a-h,o-z)
- include 'DIMENSIONS'
-#ifdef AMD64
- integer*8 iseedi8
-#endif
-#ifdef MPI
- include 'mpif.h'
- logical OKRandom, prng_restart
- real*8 r1
- integer iseed_array(4)
-#endif
- include 'COMMON.IOUNITS'
- include 'COMMON.TIME1'
-c include 'COMMON.THREAD'
- include 'COMMON.SBRIDGE'
- include 'COMMON.CONTROL'
- include 'COMMON.MCM'
-c include 'COMMON.MAP'
- include 'COMMON.HEADER'
-c include 'COMMON.CSA'
- include 'COMMON.CHAIN'
-c include 'COMMON.MUCA'
-c include 'COMMON.MD'
- include 'COMMON.FFIELD'
- include 'COMMON.SETUP'
- iseed=-dint(dabs(seed))
- if (iseed.eq.0) then
- write (iout,'(/80(1h*)/20x,a/80(1h*))')
- & 'Random seed undefined. The program will stop.'
- write (*,'(/80(1h*)/20x,a/80(1h*))')
- & 'Random seed undefined. The program will stop.'
-#ifdef MPI
- call mpi_finalize(mpi_comm_world,ierr)
-#endif
- stop 'Bad random seed.'
- endif
-#ifdef MPI
- if (fg_rank.eq.0) then
- seed=seed*(me+1)+1
-#ifdef AMD64
- iseedi8=dint(seed)
- if(me.eq.king .or. .not. out1file)
- & write (iout,*) 'MPI: node= ', me, ' iseed= ',iseedi8
- write (*,*) 'MPI: node= ', me, ' iseed= ',iseedi8
- OKRandom = prng_restart(me,iseedi8)
-#else
- do i=1,4
- tmp=65536.0d0**(4-i)
- iseed_array(i) = dint(seed/tmp)
- seed=seed-iseed_array(i)*tmp
- enddo
- if(me.eq.king .or. .not. out1file)
- & write (iout,*) 'MPI: node= ', me, ' iseed(4)= ',
- & (iseed_array(i),i=1,4)
- write (*,*) 'MPI: node= ',me, ' iseed(4)= ',
- & (iseed_array(i),i=1,4)
- OKRandom = prng_restart(me,iseed_array)
-#endif
- if (OKRandom) then
- r1=ran_number(0.0D0,1.0D0)
- if(me.eq.king .or. .not. out1file)
- & write (iout,*) 'ran_num',r1
- if (r1.lt.0.0d0) OKRandom=.false.
- endif
- if (.not.OKRandom) then
- write (iout,*) 'PRNG IS NOT WORKING!!!'
- print *,'PRNG IS NOT WORKING!!!'
- if (me.eq.0) then
- call flush(iout)
- call mpi_abort(mpi_comm_world,error_msg,ierr)
- stop
- else
- write (iout,*) 'too many processors for parallel prng'
- write (*,*) 'too many processors for parallel prng'
- call flush(iout)
- stop
- endif
- endif
- endif
-#else
- call vrndst(iseed)
- write (iout,*) 'ran_num',ran_number(0.0d0,1.0d0)
-#endif
- return
- end
projects / unres.git / blobdiff