source/cluster/wham/src-M/initialize_p.F

   1       subroutine initialize
   2 C
   3 C Define constants and zero out tables.
   4 C
   5       implicit real*8 (a-h,o-z)
   6       include 'DIMENSIONS'
   7       include 'sizesclu.dat'
   8       include 'COMMON.IOUNITS'
   9       include 'COMMON.CHAIN'
  10       include 'COMMON.INTERACT'
  11       include 'COMMON.GEO'
  12       include 'COMMON.LOCAL'
  13       include 'COMMON.TORSION'
  14       include 'COMMON.FFIELD'
  15       include 'COMMON.SBRIDGE'
  16       include 'COMMON.MINIM'
  17       include 'COMMON.DERIV'
  18       include "COMMON.NAMES"
  19       include "COMMON.TIME1"
  20 C
  21 C The following is just to define auxiliary variables used in angle conversion
  22 C
  23       pi=4.0D0*datan(1.0D0)
  24       dwapi=2.0D0*pi
  25       dwapi3=dwapi/3.0D0
  26       pipol=0.5D0*pi
  27       deg2rad=pi/180.0D0
  28       rad2deg=1.0D0/deg2rad
  29       angmin=10.0D0*deg2rad
  30       Rgas = 1.987D-3
  31 C
  32 C Define I/O units.
  33 C
  34       inp=    1
  35       iout=   2
  36       ipdbin= 3
  37       ipdb=   7
  38       imol2= 18
  39       jplot= 19
  40       jstatin=10
  41       imol2=  4
  42       igeom=  8
  43       intin=  9
  44       ithep= 11
  45       irotam=12
  46       itorp= 13
  47       itordp= 23
  48       ielep= 14
  49       isidep=15
  50       isidep1=22
  51       iscpp=25
  52       icbase=16
  53       ifourier=20
  54       istat= 17
  55       ibond=28
  56       isccor=29
  57       jrms=30
  58       iliptran=60
  59 C
  60 C Set default weights of the energy terms.
  61 C
  62       wlong=1.0D0
  63       welec=1.0D0
  64       wtor =1.0D0
  65       wang =1.0D0
  66       wscloc=1.0D0
  67       wstrain=1.0D0
  68 C
  69 C Zero out tables.
  70 C
  71       ndih_constr=0
  72       do i=1,maxres2
  73         do j=1,3
  74           c(j,i)=0.0D0
  75           dc(j,i)=0.0D0
  76         enddo
  77       enddo
  78       do i=1,maxres
  79         do j=1,3
  80           xloc(j,i)=0.0D0
  81         enddo
  82       enddo
  83       do i=1,ntyp
  84         do j=1,ntyp
  85           aa_aq(i,j)=0.0D0
  86           bb_aq(i,j)=0.0D0
  87           aa_lip(i,j)=0.0D0
  88           bb_lip(i,j)=0.0D0
  89           augm(i,j)=0.0D0
  90           sigma(i,j)=0.0D0
  91           r0(i,j)=0.0D0
  92           chi(i,j)=0.0D0
  93         enddo
  94         do j=1,2
  95           bad(i,j)=0.0D0
  96         enddo
  97         chip(i)=0.0D0
  98         alp(i)=0.0D0
  99         sigma0(i)=0.0D0
 100         sigii(i)=0.0D0
 101         rr0(i)=0.0D0
 102         a0thet(i)=0.0D0
 103         do j=1,2
 104          do ichir1=-1,1
 105           do ichir2=-1,1
 106           athet(j,i,ichir1,ichir2)=0.0D0
 107           bthet(j,i,ichir1,ichir2)=0.0D0
 108           enddo
 109          enddo
 110         enddo
 111         do j=0,3
 112           polthet(j,i)=0.0D0
 113         enddo
 114         do j=1,3
 115           gthet(j,i)=0.0D0
 116         enddo
 117         theta0(i)=0.0D0
 118         sig0(i)=0.0D0
 119         sigc0(i)=0.0D0
 120         do j=1,maxlob
 121           bsc(j,i)=0.0D0
 122           do k=1,3
 123             censc(k,j,i)=0.0D0
 124           enddo
 125           do k=1,3
 126             do l=1,3
 127               gaussc(l,k,j,i)=0.0D0
 128             enddo
 129           enddo
 130           nlob(i)=0
 131         enddo
 132       enddo
 133       nlob(ntyp1)=0
 134       dsc(ntyp1)=0.0D0
 135       do i=-maxtor,maxtor
 136         itortyp(i)=0
 137        do iblock=1,2
 138         do j=-maxtor,maxtor
 139           do k=1,maxterm
 140             v1(k,j,i,iblock)=0.0D0
 141             v2(k,j,i,iblock)=0.0D0
 142            enddo
 143          enddo
 144        enddo
 145       enddo
 146       do iblock=1,2
 147        do i=-maxtor,maxtor
 148         do j=-maxtor,maxtor
 149          do k=-maxtor,maxtor
 150           do l=1,maxtermd_1
 151             v1c(1,l,i,j,k,iblock)=0.0D0
 152             v1s(1,l,i,j,k,iblock)=0.0D0
 153             v1c(2,l,i,j,k,iblock)=0.0D0
 154             v1s(2,l,i,j,k,iblock)=0.0D0
 155           enddo !l
 156           do l=1,maxtermd_2
 157            do m=1,maxtermd_2
 158             v2c(m,l,i,j,k,iblock)=0.0D0
 159             v2s(m,l,i,j,k,iblock)=0.0D0
 160            enddo !m
 161           enddo !l
 162         enddo !k
 163        enddo !j
 164       enddo !i
 165       enddo !iblock
 166       do i=1,maxres
 167         itype(i)=0
 168         itel(i)=0
 169       enddo
 170 C Initialize the bridge arrays
 171       ns=0
 172       nss=0
 173       nhpb=0
 174       do i=1,maxss
 175         iss(i)=0
 176       enddo
 177       do i=1,maxss
 178         dhpb(i)=0.0D0
 179       enddo
 180       do i=1,maxss
 181         ihpb(i)=0
 182         jhpb(i)=0
 183       enddo
 184 C
 185 C Initialize timing.
 186 C
 187       call set_timers
 188 C
 189 C Initialize variables used in minimization.
 190 C
 191 c     maxfun=5000
 192 c     maxit=2000
 193       maxfun=500
 194       maxit=200
 195       tolf=1.0D-2
 196       rtolf=5.0D-4
 197 C
 198 C Initialize the variables responsible for the mode of gradient storage.
 199 C
 200       nfl=0
 201       icg=1
 202       do i=1,14
 203         do j=1,14
 204           if (print_order(i).eq.j) then
 205             iw(print_order(i))=j
 206             goto 1121
 207           endif
 208         enddo
 209 1121    continue
 210       enddo
 211       calc_grad=.false.
 212 C Set timers and counters for the respective routines
 213       t_func = 0.0d0
 214       t_grad = 0.0d0
 215       t_fhel = 0.0d0
 216       t_fbet = 0.0d0
 217       t_ghel = 0.0d0
 218       t_gbet = 0.0d0
 219       t_viol = 0.0d0
 220       t_gviol = 0.0d0
 221       n_func = 0
 222       n_grad = 0
 223       n_fhel = 0
 224       n_fbet = 0
 225       n_ghel = 0
 226       n_gbet = 0
 227       n_viol = 0
 228       n_gviol = 0
 229       n_map = 0
 230 #ifndef SPLITELE
 231       nprint_ene=nprint_ene-1
 232 #endif
 233       return
 234       end
 235 c-------------------------------------------------------------------------
 236       block data nazwy
 237       implicit real*8 (a-h,o-z)
 238       include 'DIMENSIONS'
 239       include 'sizesclu.dat'
 240       include 'COMMON.NAMES'
 241       include 'COMMON.FFIELD'
 242       data restyp /
 243      &'DD','DAU','DAI','DDB','DSM','DPR','DLY','DAR','DHI','DAS','DGL',
 244      & 'DSG','DGN','DSN','DTH',
 245      &'DYY','DAL','DTY','DTR','DVA','DLE','DIL','DPN','MED','DCY','ZER',
 246      &'CYS','MET','PHE','ILE','LEU','VAL','TRP','TYR','ALA','GLY','THR',
 247      &'SER','GLN','ASN','GLU','ASP','HIS','ARG','LYS','PRO','SME','DBZ',
 248      &'AIB','ABU','D'/
 249       data onelet /
 250      &'z','z','z','z','z','p','k','r','h','d','e','n','q','s','t','g',
 251      &'a','y','w','v','l','i','f','m','c','x',
 252      &'C','M','F','I','L','V','W','Y','A','G','T',
 253      &'S','Q','N','E','D','H','R','K','P','z','z','z','z','X'/
 254       data potname /'LJ','LJK','BP','GB','GBV'/
 255       data ename /
 256      &   "EVDW SC-SC","EVDW2 SC-p","EES p-p","ECORR4 ","ECORR5 ",
 257      &   "ECORR6 ","EELLO ","ETURN3 ","ETURN4 ","ETURN6 ",
 258      &   "EBE bend","ESC SCloc","ETORS ","ETORSD ","EHPB","EVDWPP",
 259      &   "EVDW2_14","ESTR","ESCCOR","EDIHC","EVDW_T","ELIPTRAN",
 260      &   "EAFM","ETHETC","EMPTY"/
 261       data wname /
 262      &   "WSC","WSCP","WELEC","WCORR","WCORR5","WCORR6","WEL_LOC",
 263      &   "WTURN3","WTURN4","WTURN6","WANG","WSCLOC","WTOR","WTORD",
 264      &   "WHPB","WVDWPP","WSCP14","WBOND","WSCCOR","WDIHC","WSC",
 265      &   "WLIPTRAN","WAFM","WTHETC","WSHIELD"/
 266       data nprint_ene /22/
 267       data print_order /1,2,3,18,11,12,13,14,4,5,6,7,8,9,10,19,
 268      &  16,15,17,20,21,24,22,23,1/
 269       end
 270 c---------------------------------------------------------------------------
 271       subroutine init_int_table
 272       implicit real*8 (a-h,o-z)
 273       include 'DIMENSIONS'
 274       include 'COMMON.CHAIN'
 275       include 'COMMON.INTERACT'
 276       include 'COMMON.LOCAL'
 277       include 'COMMON.SBRIDGE'
 278       include 'COMMON.IOUNITS'
 279       include "COMMON.TORCNSTR"
 280       logical scheck,lprint
 281       lprint=.false.
 282       do i=1,maxres
 283         nint_gr(i)=0
 284         nscp_gr(i)=0
 285         do j=1,maxint_gr
 286           istart(i,1)=0
 287           iend(i,1)=0
 288           ielstart(i)=0
 289           ielend(i)=0
 290           iscpstart(i,1)=0
 291           iscpend(i,1)=0
 292         enddo
 293       enddo
 294       ind_scint=0
 295       ind_scint_old=0
 296 cd    write (iout,*) 'ns=',ns,' nss=',nss,' ihpb,jhpb',
 297 cd   &   (ihpb(i),jhpb(i),i=1,nss)
 298       do i=nnt,nct-1
 299         scheck=.false.
 300         if (dyn_ss) goto 10
 301         do ii=1,nss
 302           if (ihpb(ii).eq.i+nres) then
 303             scheck=.true.
 304             jj=jhpb(ii)-nres
 305             goto 10
 306           endif
 307         enddo
 308    10   continue
 309 cd      write (iout,*) 'i=',i,' scheck=',scheck,' jj=',jj
 310         if (scheck) then
 311           if (jj.eq.i+1) then
 312             nint_gr(i)=1
 313             istart(i,1)=i+2
 314             iend(i,1)=nct
 315           else if (jj.eq.nct) then
 316             nint_gr(i)=1
 317             istart(i,1)=i+1
 318             iend(i,1)=nct-1
 319           else
 320             nint_gr(i)=2
 321             istart(i,1)=i+1
 322             iend(i,1)=jj-1
 323             istart(i,2)=jj+1
 324             iend(i,2)=nct
 325           endif
 326         else
 327           nint_gr(i)=1
 328           istart(i,1)=i+1
 329           iend(i,1)=nct
 330           ind_scint=int_scint+nct-i
 331         endif
 332       enddo
 333    12 continue
 334       iatsc_s=nnt
 335       iatsc_e=nct-1
 336       if (lprint) then
 337       write (iout,'(a)') 'Interaction array:'
 338       do i=iatsc_s,iatsc_e
 339         write (iout,'(i3,2(2x,2i3))')
 340      & i,(istart(i,iint),iend(i,iint),iint=1,nint_gr(i))
 341       enddo
 342       endif
 343       ispp=2
 344       iatel_s=nnt
 345       iatel_e=nct-3
 346       do i=iatel_s,iatel_e
 347         ielstart(i)=i+4
 348         ielend(i)=nct-1
 349       enddo
 350       if (lprint) then
 351         write (iout,'(a)') 'Electrostatic interaction array:'
 352         do i=iatel_s,iatel_e
 353           write (iout,'(i3,2(2x,2i3))') i,ielstart(i),ielend(i)
 354         enddo
 355       endif ! lprint
 356 c     iscp=3
 357       iscp=2
 358 C Partition the SC-p interaction array
 359       iatscp_s=nnt
 360       iatscp_e=nct-1
 361       do i=nnt,nct-1
 362         if (i.lt.nnt+iscp) then
 363           nscp_gr(i)=1
 364           iscpstart(i,1)=i+iscp
 365           iscpend(i,1)=nct
 366         elseif (i.gt.nct-iscp) then
 367           nscp_gr(i)=1
 368           iscpstart(i,1)=nnt
 369           iscpend(i,1)=i-iscp
 370         else
 371           nscp_gr(i)=2
 372           iscpstart(i,1)=nnt
 373           iscpend(i,1)=i-iscp
 374           iscpstart(i,2)=i+iscp
 375           iscpend(i,2)=nct
 376         endif
 377       enddo ! i
 378       if (lprint) then
 379         write (iout,'(a)') 'SC-p interaction array:'
 380         do i=iatscp_s,iatscp_e
 381           write (iout,'(i3,2(2x,2i3))')
 382      &         i,(iscpstart(i,j),iscpend(i,j),j=1,nscp_gr(i))
 383         enddo
 384       endif ! lprint
 385 C Partition local interactions
 386       loc_start=2
 387       loc_end=nres-1
 388       ithet_start=3
 389       ithet_end=nres
 390       iturn3_start=nnt
 391       iturn3_end=nct-3
 392       iturn4_start=nnt
 393       iturn4_end=nct-4
 394       iphi_start=nnt+3
 395       iphi_end=nct
 396       idihconstr_start=1
 397       idihconstr_end=ndih_constr
 398       ithetaconstr_start=1
 399       ithetaconstr_end=ntheta_constr
 400       itau_start=4
 401       itau_end=nres
 402       return
 403       end
 404 c---------------------------------------------------------------------------
 405       subroutine int_partition(int_index,lower_index,upper_index,atom,
 406      & at_start,at_end,first_atom,last_atom,int_gr,jat_start,jat_end,*)
 407       implicit real*8 (a-h,o-z)
 408       include 'DIMENSIONS'
 409       include 'COMMON.IOUNITS'
 410       integer int_index,lower_index,upper_index,atom,at_start,at_end,
 411      & first_atom,last_atom,int_gr,jat_start,jat_end
 412       logical lprn
 413       lprn=.false.
 414       if (lprn) write (iout,*) 'int_index=',int_index
 415       int_index_old=int_index
 416       int_index=int_index+last_atom-first_atom+1
 417       if (lprn)
 418      &   write (iout,*) 'int_index=',int_index,
 419      &               ' int_index_old',int_index_old,
 420      &               ' lower_index=',lower_index,
 421      &               ' upper_index=',upper_index,
 422      &               ' atom=',atom,' first_atom=',first_atom,
 423      &               ' last_atom=',last_atom
 424       if (int_index.ge.lower_index) then
 425         int_gr=int_gr+1
 426         if (at_start.eq.0) then
 427           at_start=atom
 428           jat_start=first_atom-1+lower_index-int_index_old
 429         else
 430           jat_start=first_atom
 431         endif
 432         if (lprn) write (iout,*) 'jat_start',jat_start
 433         if (int_index.ge.upper_index) then
 434           at_end=atom
 435           jat_end=first_atom-1+upper_index-int_index_old
 436           return1
 437         else
 438           jat_end=last_atom
 439         endif
 440         if (lprn) write (iout,*) 'jat_end',jat_end
 441       endif
 442       return
 443       end
 444 c------------------------------------------------------------------------------
 445       subroutine hpb_partition
 446       implicit real*8 (a-h,o-z)
 447       include 'DIMENSIONS'
 448       include 'COMMON.SBRIDGE'
 449       include 'COMMON.IOUNITS'
 450       link_start=1
 451       link_end=nhpb
 452 cd    write (iout,*) 'Processor',MyID,' MyRank',MyRank,
 453 cd   &  ' nhpb',nhpb,' link_start=',link_start,
 454 cd   &  ' link_end',link_end
 455       return
 456       end