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