--- /dev/null
+ subroutine initialize
+C
+C Define constants and zero out tables.
+C
+ implicit real*8 (a-h,o-z)
+ include 'DIMENSIONS'
+ include 'DIMENSIONS.ZSCOPT'
+#ifdef MPI
+ include 'mpif.h'
+#endif
+ include 'COMMON.IOUNITS'
+ include 'COMMON.CHAIN'
+ include 'COMMON.INTERACT'
+ include 'COMMON.GEO'
+ include 'COMMON.LOCAL'
+ include 'COMMON.TORSION'
+ include 'COMMON.FFIELD'
+ include 'COMMON.SBRIDGE'
+ include 'COMMON.MINIM'
+ include 'COMMON.DERIV'
+ include "COMMON.WEIGHTS"
+ include "COMMON.NAMES"
+ include "COMMON.TIME1"
+C
+C The following is just to define auxiliary variables used in angle conversion
+C
+ pi=4.0D0*datan(1.0D0)
+ dwapi=2.0D0*pi
+ dwapi3=dwapi/3.0D0
+ pipol=0.5D0*pi
+ deg2rad=pi/180.0D0
+ rad2deg=1.0D0/deg2rad
+ angmin=10.0D0*deg2rad
+C
+C Define I/O units.
+C
+ inp= 1
+ iout= 2
+ ipdbin= 3
+ ipdb= 7
+ imol2= 4
+ igeom= 8
+ intin= 9
+ ithep= 11
+ irotam=12
+ itorp= 13
+ itordp= 23
+ ielep= 14
+ isidep=15
+ isidep1=22
+ iscpp=25
+ icbase=16
+ ifourier=20
+ istat= 17
+ ientin=18
+ ientout=19
+ ibond=28
+ isccor=29
+C
+C WHAM files
+C
+ ihist=30
+ iweight=31
+ izsc=32
+C
+C Set default weights of the energy terms.
+C
+ wlong=1.0D0
+ welec=1.0D0
+ wtor =1.0D0
+ wang =1.0D0
+ wscloc=1.0D0
+ wstrain=1.0D0
+C
+C Zero out tables.
+C
+ ndih_constr=0
+ do i=1,maxres2
+ do j=1,3
+ c(j,i)=0.0D0
+ dc(j,i)=0.0D0
+ enddo
+ enddo
+ do i=1,maxres
+ do j=1,3
+ xloc(j,i)=0.0D0
+ enddo
+ enddo
+ do i=1,ntyp
+ do j=1,ntyp
+ aa(i,j)=0.0D0
+ bb(i,j)=0.0D0
+ augm(i,j)=0.0D0
+ sigma(i,j)=0.0D0
+ r0(i,j)=0.0D0
+ chi(i,j)=0.0D0
+ enddo
+ do j=1,2
+ bad(i,j)=0.0D0
+ enddo
+ chip(i)=0.0D0
+ alp(i)=0.0D0
+ sigma0(i)=0.0D0
+ sigii(i)=0.0D0
+ rr0(i)=0.0D0
+ a0thet(i)=0.0D0
+ do j=1,2
+ athet(j,i)=0.0D0
+ bthet(j,i)=0.0D0
+ enddo
+ do j=0,3
+ polthet(j,i)=0.0D0
+ enddo
+ do j=1,3
+ gthet(j,i)=0.0D0
+ enddo
+ theta0(i)=0.0D0
+ sig0(i)=0.0D0
+ sigc0(i)=0.0D0
+ do j=1,maxlob
+ bsc(j,i)=0.0D0
+ do k=1,3
+ censc(k,j,i)=0.0D0
+ enddo
+ do k=1,3
+ do l=1,3
+ gaussc(l,k,j,i)=0.0D0
+ enddo
+ enddo
+ nlob(i)=0
+ enddo
+ enddo
+ nlob(ntyp1)=0
+ dsc(ntyp1)=0.0D0
+ do i=1,maxtor
+ itortyp(i)=0
+ do j=1,maxtor
+ do k=1,maxterm
+ v1(k,j,i)=0.0D0
+ v2(k,j,i)=0.0D0
+ enddo
+ enddo
+ enddo
+ do i=1,maxres
+ itype(i)=0
+ itel(i)=0
+ enddo
+C Initialize the bridge arrays
+ ns=0
+ nss=0
+ nhpb=0
+ do i=1,maxss
+ iss(i)=0
+ enddo
+ do i=1,maxdim
+ dhpb(i)=0.0D0
+ enddo
+ do i=1,maxres
+ ihpb(i)=0
+ jhpb(i)=0
+ enddo
+C
+C Initialize timing.
+C
+ call set_timers
+C
+C Initialize variables used in minimization.
+C
+c maxfun=5000
+c maxit=2000
+ maxfun=500
+ maxit=200
+ tolf=1.0D-2
+ rtolf=5.0D-4
+C
+C Initialize the variables responsible for the mode of gradient storage.
+C
+ nfl=0
+ icg=1
+ do i=1,14
+ do j=1,14
+ if (print_order(i).eq.j) then
+ iw(print_order(i))=j
+ goto 1121
+ endif
+ enddo
+1121 continue
+ enddo
+ calc_grad=.false.
+C Set timers and counters for the respective routines
+ t_func = 0.0d0
+ t_grad = 0.0d0
+ t_fhel = 0.0d0
+ t_fbet = 0.0d0
+ t_ghel = 0.0d0
+ t_gbet = 0.0d0
+ t_viol = 0.0d0
+ t_gviol = 0.0d0
+ n_func = 0
+ n_grad = 0
+ n_fhel = 0
+ n_fbet = 0
+ n_ghel = 0
+ n_gbet = 0
+ n_viol = 0
+ n_gviol = 0
+ n_map = 0
+#ifndef SPLITELE
+ nprint_ene=nprint_ene-1
+#endif
+ return
+ end
+c-------------------------------------------------------------------------
+ block data nazwy
+ implicit real*8 (a-h,o-z)
+ include 'DIMENSIONS'
+ include 'DIMENSIONS.ZSCOPT'
+ include 'COMMON.NAMES'
+ include 'COMMON.WEIGHTS'
+ include 'COMMON.FFIELD'
+ include 'COMMON.INTERACT'
+ data restyp /
+ &'CYS','MET','PHE','ILE','LEU','VAL','TRP','TYR','ALA','GLY','THR',
+ &'SER','GLN','ASN','GLU','ASP','HIS','ARG','LYS','PRO','D'/
+ data onelet /
+ &'C','M','F','I','L','V','W','Y','A','G','T',
+ &'S','Q','N','E','D','H','R','K','P','X'/
+ data potname /'LJ','LJK','BP','GB','GBV','MM'/
+ data ename /
+ & "EVDW SC-SC","EVDW2 SC-p","EES p-p","ECORR4 ","ECORR5 ",
+ & "ECORR6 ","EELLO ","ETURN3 ","ETURN4 ","ETURN6 ",
+ & "EBE bend","ESC SCloc","ETORS ","ETORSD ","EHPB","EVDWPP",
+ & "EVDW2_14","ESTR","ESCCOR","EDIHC","EVDW_T"/
+ data wname /
+ & "WSC","WSCP","WELEC","WCORR","WCORR5","WCORR6","WEL_LOC",
+ & "WTURN3","WTURN4","WTURN6","WANG","WSCLOC","WTOR","WTORD",
+ & "WHPB","WVDWPP","WSCP14","WBOND","WSCCOR","WDIHC","WSC"/
+ data ww0 /1.0d0,1.0d0,1.0d0,1.0d0,1.0d0,1.0d0,1.0d0,1.0d0,1.0d0,
+ & 1.0d0,1.0d0,1.0d0,1.0d0,1.0d0,1.0d0,1.0d0,0.4d0,1.0d0,1.0d0,
+ & 0.0d0,0.0/
+ data nprint_ene /21/
+ data print_order /1,2,3,18,11,12,13,14,4,5,6,7,8,9,10,19,
+ & 16,15,17,20,21/
+c Dielectric constant of water
+ data eps_out /80.0d0/
+ end
+c---------------------------------------------------------------------------
+ subroutine init_int_table
+ implicit real*8 (a-h,o-z)
+ include 'DIMENSIONS'
+ include 'DIMENSIONS.ZSCOPT'
+#ifdef MPI
+ include 'mpif.h'
+#endif
+#ifdef MP
+ include 'COMMON.INFO'
+#endif
+ include 'COMMON.CHAIN'
+ include 'COMMON.INTERACT'
+ include 'COMMON.LOCAL'
+ include 'COMMON.SBRIDGE'
+ include 'COMMON.IOUNITS'
+ logical scheck,lprint
+#ifdef MPL
+ integer my_sc_int(0:max_fg_Procs-1),my_sc_intt(0:max_fg_Procs),
+ & my_ele_int(0:max_fg_Procs-1),my_ele_intt(0:max_fg_Procs)
+C... Determine the numbers of start and end SC-SC interaction
+C... to deal with by current processor.
+ lprint=.false.
+ if (lprint)
+ &write (iout,*) 'INIT_INT_TABLE nres=',nres,' nnt=',nnt,' nct=',nct
+ n_sc_int_tot=(nct-nnt+1)*(nct-nnt)/2-nss
+ MyRank=MyID-(MyGroup-1)*fgProcs
+ call int_bounds(n_sc_int_tot,my_sc_inds,my_sc_inde)
+ if (lprint)
+ & write (iout,*) 'Processor',MyID,' MyRank',MyRank,
+ & ' n_sc_int_tot',n_sc_int_tot,' my_sc_inds=',my_sc_inds,
+ & ' my_sc_inde',my_sc_inde
+ ind_sctint=0
+ iatsc_s=0
+ iatsc_e=0
+#endif
+ lprint=.false.
+ do i=1,maxres
+ nint_gr(i)=0
+ nscp_gr(i)=0
+ do j=1,maxint_gr
+ istart(i,1)=0
+ iend(i,1)=0
+ ielstart(i)=0
+ ielend(i)=0
+ iscpstart(i,1)=0
+ iscpend(i,1)=0
+ enddo
+ enddo
+ ind_scint=0
+ ind_scint_old=0
+cd write (iout,*) 'ns=',ns,' nss=',nss,' ihpb,jhpb',
+cd & (ihpb(i),jhpb(i),i=1,nss)
+ do i=nnt,nct-1
+ scheck=.false.
+ do ii=1,nss
+ if (ihpb(ii).eq.i+nres) then
+ scheck=.true.
+ jj=jhpb(ii)-nres
+ goto 10
+ endif
+ enddo
+ 10 continue
+cd write (iout,*) 'i=',i,' scheck=',scheck,' jj=',jj
+ if (scheck) then
+ if (jj.eq.i+1) then
+#ifdef MPL
+ write (iout,*) 'jj=i+1'
+ call int_partition(ind_scint,my_sc_inds,my_sc_inde,i,
+ & iatsc_s,iatsc_e,i+2,nct,nint_gr(i),istart(i,1),iend(i,1),*12)
+#else
+ nint_gr(i)=1
+ istart(i,1)=i+2
+ iend(i,1)=nct
+#endif
+ else if (jj.eq.nct) then
+#ifdef MPL
+ write (iout,*) 'jj=nct'
+ call int_partition(ind_scint,my_sc_inds,my_sc_inde,i,
+ & iatsc_s,iatsc_e,i+1,nct-1,nint_gr(i),istart(i,1),iend(i,1),*12)
+#else
+ nint_gr(i)=1
+ istart(i,1)=i+1
+ iend(i,1)=nct-1
+#endif
+ else
+#ifdef MPL
+ call int_partition(ind_scint,my_sc_inds,my_sc_inde,i,
+ & iatsc_s,iatsc_e,i+1,jj-1,nint_gr(i),istart(i,1),iend(i,1),*12)
+ ii=nint_gr(i)+1
+ call int_partition(ind_scint,my_sc_inds,my_sc_inde,i,
+ & iatsc_s,iatsc_e,jj+1,nct,nint_gr(i),istart(i,ii),iend(i,ii),*12)
+#else
+ nint_gr(i)=2
+ istart(i,1)=i+1
+ iend(i,1)=jj-1
+ istart(i,2)=jj+1
+ iend(i,2)=nct
+#endif
+ endif
+ else
+#ifdef MPL
+ call int_partition(ind_scint,my_sc_inds,my_sc_inde,i,
+ & iatsc_s,iatsc_e,i+1,nct,nint_gr(i),istart(i,1),iend(i,1),*12)
+#else
+ nint_gr(i)=1
+ istart(i,1)=i+1
+ iend(i,1)=nct
+ ind_scint=int_scint+nct-i
+#endif
+ endif
+#ifdef MPL
+ ind_scint_old=ind_scint
+#endif
+ enddo
+ 12 continue
+#ifndef MPL
+ iatsc_s=nnt
+ iatsc_e=nct-1
+#endif
+#ifdef MPL
+ if (lprint) then
+ write (iout,*) 'Processor',MyID,' Group',MyGroup
+ write (iout,*) 'iatsc_s=',iatsc_s,' iatsc_e=',iatsc_e
+ endif
+#endif
+ if (lprint) then
+ write (iout,'(a)') 'Interaction array:'
+ do i=iatsc_s,iatsc_e
+ write (iout,'(i3,2(2x,2i3))')
+ & i,(istart(i,iint),iend(i,iint),iint=1,nint_gr(i))
+ enddo
+ endif
+ ispp=2
+#ifdef MPL
+C Now partition the electrostatic-interaction array
+ npept=nct-nnt
+ nele_int_tot=(npept-ispp)*(npept-ispp+1)/2
+ call int_bounds(nele_int_tot,my_ele_inds,my_ele_inde)
+ if (lprint)
+ & write (iout,*) 'Processor',MyID,' MyRank',MyRank,
+ & ' nele_int_tot',nele_int_tot,' my_ele_inds=',my_ele_inds,
+ & ' my_ele_inde',my_ele_inde
+ iatel_s=0
+ iatel_e=0
+ ind_eleint=0
+ ind_eleint_old=0
+ do i=nnt,nct-3
+ ijunk=0
+ call int_partition(ind_eleint,my_ele_inds,my_ele_inde,i,
+ & iatel_s,iatel_e,i+ispp,nct-1,ijunk,ielstart(i),ielend(i),*13)
+ enddo ! i
+ 13 continue
+#else
+ iatel_s=nnt
+ iatel_e=nct-3
+ do i=iatel_s,iatel_e
+ ielstart(i)=i+2
+ ielend(i)=nct-1
+ enddo
+#endif
+ if (lprint) then
+ write (iout,'(a)') 'Electrostatic interaction array:'
+ do i=iatel_s,iatel_e
+ write (iout,'(i3,2(2x,2i3))') i,ielstart(i),ielend(i)
+ enddo
+ endif ! lprint
+c iscp=3
+ iscp=2
+C Partition the SC-p interaction array
+#ifdef MPL
+ nscp_int_tot=(npept-iscp+1)*(npept-iscp+1)
+ call int_bounds(nscp_int_tot,my_scp_inds,my_scp_inde)
+ if (lprint)
+ & write (iout,*) 'Processor',MyID,' MyRank',MyRank,
+ & ' nscp_int_tot',nscp_int_tot,' my_scp_inds=',my_scp_inds,
+ & ' my_scp_inde',my_scp_inde
+ iatscp_s=0
+ iatscp_e=0
+ ind_scpint=0
+ ind_scpint_old=0
+ do i=nnt,nct-1
+ if (i.lt.nnt+iscp) then
+cd write (iout,*) 'i.le.nnt+iscp'
+ call int_partition(ind_scpint,my_scp_inds,my_scp_inde,i,
+ & iatscp_s,iatscp_e,i+iscp,nct,nscp_gr(i),iscpstart(i,1),
+ & iscpend(i,1),*14)
+ else if (i.gt.nct-iscp) then
+cd write (iout,*) 'i.gt.nct-iscp'
+ call int_partition(ind_scpint,my_scp_inds,my_scp_inde,i,
+ & iatscp_s,iatscp_e,nnt,i-iscp,nscp_gr(i),iscpstart(i,1),
+ & iscpend(i,1),*14)
+ else
+ call int_partition(ind_scpint,my_scp_inds,my_scp_inde,i,
+ & iatscp_s,iatscp_e,nnt,i-iscp,nscp_gr(i),iscpstart(i,1),
+ & iscpend(i,1),*14)
+ ii=nscp_gr(i)+1
+ call int_partition(ind_scpint,my_scp_inds,my_scp_inde,i,
+ & iatscp_s,iatscp_e,i+iscp,nct,nscp_gr(i),iscpstart(i,ii),
+ & iscpend(i,ii),*14)
+ endif
+ enddo ! i
+ 14 continue
+#else
+ iatscp_s=nnt
+ iatscp_e=nct-1
+ do i=nnt,nct-1
+ if (i.lt.nnt+iscp) then
+ nscp_gr(i)=1
+ iscpstart(i,1)=i+iscp
+ iscpend(i,1)=nct
+ elseif (i.gt.nct-iscp) then
+ nscp_gr(i)=1
+ iscpstart(i,1)=nnt
+ iscpend(i,1)=i-iscp
+ else
+ nscp_gr(i)=2
+ iscpstart(i,1)=nnt
+ iscpend(i,1)=i-iscp
+ iscpstart(i,2)=i+iscp
+ iscpend(i,2)=nct
+ endif
+ enddo ! i
+#endif
+ if (lprint) then
+ write (iout,'(a)') 'SC-p interaction array:'
+ do i=iatscp_s,iatscp_e
+ write (iout,'(i3,2(2x,2i3))')
+ & i,(iscpstart(i,j),iscpend(i,j),j=1,nscp_gr(i))
+ enddo
+ endif ! lprint
+C Partition local interactions
+#ifdef MPL
+ call int_bounds(nres-2,loc_start,loc_end)
+ loc_start=loc_start+1
+ loc_end=loc_end+1
+ call int_bounds(nres-2,ithet_start,ithet_end)
+ ithet_start=ithet_start+2
+ ithet_end=ithet_end+2
+ call int_bounds(nct-nnt-2,iphi_start,iphi_end)
+ iphi_start=iphi_start+nnt+2
+ iphi_end=iphi_end+nnt+2
+ call int_bounds(nres-3,itau_start,itau_end)
+ itau_start=itau_start+3
+ itau_end=itau_end+3
+ if (lprint) then
+ write (iout,*) 'Processor:',MyID,
+ & ' loc_start',loc_start,' loc_end',loc_end,
+ & ' ithet_start',ithet_start,' ithet_end',ithet_end,
+ & ' iphi_start',iphi_start,' iphi_end',iphi_end
+ write (*,*) 'Processor:',MyID,
+ & ' loc_start',loc_start,' loc_end',loc_end,
+ & ' ithet_start',ithet_start,' ithet_end',ithet_end,
+ & ' iphi_start',iphi_start,' iphi_end',iphi_end
+ endif
+ if (fgprocs.gt.1 .and. MyID.eq.BossID) then
+ write(iout,'(i10,a,i10,a,i10,a/a,i3,a)') n_sc_int_tot,' SC-SC ',
+ & nele_int_tot,' electrostatic and ',nscp_int_tot,
+ & ' SC-p interactions','were distributed among',fgprocs,
+ & ' fine-grain processors.'
+ endif
+#else
+ loc_start=2
+ loc_end=nres-1
+ ithet_start=3
+ ithet_end=nres
+ iphi_start=nnt+3
+ iphi_end=nct
+ itau_start=4
+ itau_end=nres
+#endif
+ return
+ end
+c---------------------------------------------------------------------------
+ subroutine int_partition(int_index,lower_index,upper_index,atom,
+ & at_start,at_end,first_atom,last_atom,int_gr,jat_start,jat_end,*)
+ implicit real*8 (a-h,o-z)
+ include 'DIMENSIONS'
+ include 'COMMON.IOUNITS'
+ integer int_index,lower_index,upper_index,atom,at_start,at_end,
+ & first_atom,last_atom,int_gr,jat_start,jat_end
+ logical lprn
+ lprn=.false.
+ if (lprn) write (iout,*) 'int_index=',int_index
+ int_index_old=int_index
+ int_index=int_index+last_atom-first_atom+1
+ if (lprn)
+ & write (iout,*) 'int_index=',int_index,
+ & ' int_index_old',int_index_old,
+ & ' lower_index=',lower_index,
+ & ' upper_index=',upper_index,
+ & ' atom=',atom,' first_atom=',first_atom,
+ & ' last_atom=',last_atom
+ if (int_index.ge.lower_index) then
+ int_gr=int_gr+1
+ if (at_start.eq.0) then
+ at_start=atom
+ jat_start=first_atom-1+lower_index-int_index_old
+ else
+ jat_start=first_atom
+ endif
+ if (lprn) write (iout,*) 'jat_start',jat_start
+ if (int_index.ge.upper_index) then
+ at_end=atom
+ jat_end=first_atom-1+upper_index-int_index_old
+ return1
+ else
+ jat_end=last_atom
+ endif
+ if (lprn) write (iout,*) 'jat_end',jat_end
+ endif
+ return
+ end
+c------------------------------------------------------------------------------
+ subroutine hpb_partition
+ implicit real*8 (a-h,o-z)
+ include 'DIMENSIONS'
+ include 'COMMON.SBRIDGE'
+ include 'COMMON.IOUNITS'
+#ifdef MPL
+ include 'COMMON.INFO'
+ call int_bounds(nhpb,link_start,link_end)
+#else
+ link_start=1
+ link_end=nhpb
+#endif
+cd write (iout,*) 'Processor',MyID,' MyRank',MyRank,
+cd & ' nhpb',nhpb,' link_start=',link_start,
+cd & ' link_end',link_end
+ return
+ end