source/cluster/wham/src/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 #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.NAMES"
  22       include "COMMON.TIME1"
  23 C
  24 C The following is just to define auxiliary variables used in angle conversion
  25 C
  26       pi=4.0D0*datan(1.0D0)
  27       dwapi=2.0D0*pi
  28       dwapi3=dwapi/3.0D0
  29       pipol=0.5D0*pi
  30       deg2rad=pi/180.0D0
  31       rad2deg=1.0D0/deg2rad
  32       angmin=10.0D0*deg2rad
  33       Rgas = 1.987D-3
  34 C
  35 C Define I/O units.
  36 C
  37       inp=    1
  38       iout=   2
  39       ipdbin= 3
  40       ipdb=   7
  41       imol2= 18
  42       jplot= 19
  43       jstatin=10
  44       imol2=  4
  45       igeom=  8
  46       intin=  9
  47       ithep= 11
  48       irotam=12
  49       itorp= 13
  50       itordp= 23
  51       ielep= 14
  52       isidep=15
  53       isidep1=22
  54       iscpp=25
  55       icbase=16
  56       ifourier=20
  57       istat= 17
  58       ibond=28
  59       isccor=29
  60       jrms=30
  61 C
  62 C Set default weights of the energy terms.
  63 C
  64       wlong=1.0D0
  65       welec=1.0D0
  66       wtor =1.0D0
  67       wang =1.0D0
  68       wscloc=1.0D0
  69       wstrain=1.0D0
  70 C
  71 C Zero out tables.
  72 C
  73       ndih_constr=0
  74       do i=1,maxres2
  75         do j=1,3
  76           c(j,i)=0.0D0
  77           dc(j,i)=0.0D0
  78         enddo
  79       enddo
  80       do i=1,maxres
  81         do j=1,3
  82           xloc(j,i)=0.0D0
  83         enddo
  84       enddo
  85       do i=1,ntyp
  86         do j=1,ntyp
  87           aa(i,j)=0.0D0
  88           bb(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           athet(j,i)=0.0D0
 105           bthet(j,i)=0.0D0
 106         enddo
 107         do j=0,3
 108           polthet(j,i)=0.0D0
 109         enddo
 110         do j=1,3
 111           gthet(j,i)=0.0D0
 112         enddo
 113         theta0(i)=0.0D0
 114         sig0(i)=0.0D0
 115         sigc0(i)=0.0D0
 116         do j=1,maxlob
 117           bsc(j,i)=0.0D0
 118           do k=1,3
 119             censc(k,j,i)=0.0D0
 120           enddo
 121           do k=1,3
 122             do l=1,3
 123               gaussc(l,k,j,i)=0.0D0
 124             enddo
 125           enddo
 126           nlob(i)=0
 127         enddo
 128       enddo
 129       nlob(ntyp1)=0
 130       dsc(ntyp1)=0.0D0
 131       do i=1,maxtor
 132         itortyp(i)=0
 133         do j=1,maxtor
 134           do k=1,maxterm
 135             v1(k,j,i)=0.0D0
 136             v2(k,j,i)=0.0D0
 137           enddo
 138         enddo
 139       enddo
 140       do i=1,maxres
 141         itype(i)=0
 142         itel(i)=0
 143       enddo
 144 C Initialize the bridge arrays
 145       ns=0
 146       nss=0
 147       nhpb=0
 148       do i=1,maxss
 149         iss(i)=0
 150       enddo
 151       do i=1,maxss
 152         dhpb(i)=0.0D0
 153       enddo
 154       do i=1,maxss
 155         ihpb(i)=0
 156         jhpb(i)=0
 157       enddo
 158       do i=1,maxres
 159         dyn_ss_mask(i)=.false.
 160       enddo
 161 C
 162 C Initialize timing.
 163 C
 164       call set_timers
 165 C
 166 C Initialize variables used in minimization.
 167 C
 168 c     maxfun=5000
 169 c     maxit=2000
 170       maxfun=500
 171       maxit=200
 172       tolf=1.0D-2
 173       rtolf=5.0D-4
 174 C
 175 C Initialize the variables responsible for the mode of gradient storage.
 176 C
 177       nfl=0
 178       icg=1
 179       do i=1,14
 180         do j=1,14
 181           if (print_order(i).eq.j) then
 182             iw(print_order(i))=j
 183             goto 1121
 184           endif
 185         enddo
 186 1121    continue
 187       enddo
 188       calc_grad=.false.
 189 C Set timers and counters for the respective routines
 190       t_func = 0.0d0
 191       t_grad = 0.0d0
 192       t_fhel = 0.0d0
 193       t_fbet = 0.0d0
 194       t_ghel = 0.0d0
 195       t_gbet = 0.0d0
 196       t_viol = 0.0d0
 197       t_gviol = 0.0d0
 198       n_func = 0
 199       n_grad = 0
 200       n_fhel = 0
 201       n_fbet = 0
 202       n_ghel = 0
 203       n_gbet = 0
 204       n_viol = 0
 205       n_gviol = 0
 206       n_map = 0
 207 #ifndef SPLITELE
 208       nprint_ene=nprint_ene-1
 209 #endif
 210       return
 211       end
 212 c-------------------------------------------------------------------------
 213       block data nazwy
 214       implicit real*8 (a-h,o-z)
 215       include 'DIMENSIONS'
 216       include 'sizesclu.dat'
 217       include 'COMMON.NAMES'
 218       include 'COMMON.FFIELD'
 219       data restyp /
 220      &'CYS','MET','PHE','ILE','LEU','VAL','TRP','TYR','ALA','GLY','THR',
 221      &'SER','GLN','ASN','GLU','ASP','HIS','ARG','LYS','PRO','D'/
 222       data onelet /
 223      &'C','M','F','I','L','V','W','Y','A','G','T',
 224      &'S','Q','N','E','D','H','R','K','P','X'/
 225       data potname /'LJ','LJK','BP','GB','GBV'/
 226       data ename /
 227      &   "EVDW SC-SC","EVDW2 SC-p","EES p-p","ECORR4 ","ECORR5 ",
 228      &   "ECORR6 ","EELLO ","ETURN3 ","ETURN4 ","ETURN6 ",
 229      &   "EBE bend","ESC SCloc","ETORS ","ETORSD ","EHPB","EVDWPP",
 230      &   "ESTR","ESCCOR","EVDW2_14",""/
 231       data wname /
 232      &   "WSC","WSCP","WELEC","WCORR","WCORR5","WCORR6","WEL_LOC",
 233      &   "WTURN3","WTURN4","WTURN6","WANG","WSCLOC","WTOR","WTORD",
 234      &   "WHPB","WVDWPP","WBOND","WSCCOR","WSCP14",""/
 235       data nprint_ene /18/
 236       data print_order /1,2,3,17,11,12,13,14,4,5,6,7,8,9,10,16,15,18,19,
 237      &   20/
 238       end
 239 c---------------------------------------------------------------------------
 240       subroutine init_int_table
 241       implicit real*8 (a-h,o-z)
 242       include 'DIMENSIONS'
 243       include 'sizesclu.dat'
 244 #ifdef MPI
 245       include 'mpif.h'
 246 #endif
 247 #ifdef MPL
 248       include 'COMMON.INFO'
 249 #endif
 250       include 'COMMON.CHAIN'
 251       include 'COMMON.INTERACT'
 252       include 'COMMON.LOCAL'
 253       include 'COMMON.SBRIDGE'
 254       include 'COMMON.IOUNITS'
 255       logical scheck,lprint
 256 #ifdef MPL
 257       integer my_sc_int(0:max_fg_Procs-1),my_sc_intt(0:max_fg_Procs),
 258      & my_ele_int(0:max_fg_Procs-1),my_ele_intt(0:max_fg_Procs)
 259 C... Determine the numbers of start and end SC-SC interaction
 260 C... to deal with by current processor.
 261       lprint=.false.
 262       if (lprint)
 263      &write (iout,*) 'INIT_INT_TABLE nres=',nres,' nnt=',nnt,' nct=',nct
 264       n_sc_int_tot=(nct-nnt+1)*(nct-nnt)/2-nss
 265       MyRank=MyID-(MyGroup-1)*fgProcs
 266       call int_bounds(n_sc_int_tot,my_sc_inds,my_sc_inde)
 267       if (lprint)
 268      &  write (iout,*) 'Processor',MyID,' MyRank',MyRank,
 269      &  ' n_sc_int_tot',n_sc_int_tot,' my_sc_inds=',my_sc_inds,
 270      &  ' my_sc_inde',my_sc_inde
 271       ind_sctint=0
 272       iatsc_s=0
 273       iatsc_e=0
 274 #endif
 275       lprint=.false.
 276       do i=1,maxres
 277         nint_gr(i)=0
 278         nscp_gr(i)=0
 279         do j=1,maxint_gr
 280           istart(i,1)=0
 281           iend(i,1)=0
 282           ielstart(i)=0
 283           ielend(i)=0
 284           iscpstart(i,1)=0
 285           iscpend(i,1)=0
 286         enddo
 287       enddo
 288       ind_scint=0
 289       ind_scint_old=0
 290 cd    write (iout,*) 'ns=',ns,' nss=',nss,' ihpb,jhpb',
 291 cd   &   (ihpb(i),jhpb(i),i=1,nss)
 292       do i=nnt,nct-1
 293         scheck=.false.
 294         do ii=1,nss
 295           if (ihpb(ii).eq.i+nres) then
 296             scheck=.true.
 297             if (dyn_ss) go to 10
 298             jj=jhpb(ii)-nres
 299             goto 10
 300           endif
 301         enddo
 302    10   continue
 303 cd      write (iout,*) 'i=',i,' scheck=',scheck,' jj=',jj
 304         if (scheck) then
 305           if (jj.eq.i+1) then
 306 #ifdef MPL
 307             write (iout,*) 'jj=i+1'
 308             call int_partition(ind_scint,my_sc_inds,my_sc_inde,i,
 309      & iatsc_s,iatsc_e,i+2,nct,nint_gr(i),istart(i,1),iend(i,1),*12)
 310 #else
 311             nint_gr(i)=1
 312             istart(i,1)=i+2
 313             iend(i,1)=nct
 314 #endif
 315           else if (jj.eq.nct) then
 316 #ifdef MPL
 317             write (iout,*) 'jj=nct'
 318             call int_partition(ind_scint,my_sc_inds,my_sc_inde,i,
 319      &  iatsc_s,iatsc_e,i+1,nct-1,nint_gr(i),istart(i,1),iend(i,1),*12)
 320 #else
 321             nint_gr(i)=1
 322             istart(i,1)=i+1
 323             iend(i,1)=nct-1
 324 #endif
 325           else
 326 #ifdef MPL
 327             call int_partition(ind_scint,my_sc_inds,my_sc_inde,i,
 328      & iatsc_s,iatsc_e,i+1,jj-1,nint_gr(i),istart(i,1),iend(i,1),*12)
 329             ii=nint_gr(i)+1
 330             call int_partition(ind_scint,my_sc_inds,my_sc_inde,i,
 331      & iatsc_s,iatsc_e,jj+1,nct,nint_gr(i),istart(i,ii),iend(i,ii),*12)
 332 #else
 333             nint_gr(i)=2
 334             istart(i,1)=i+1
 335             iend(i,1)=jj-1
 336             istart(i,2)=jj+1
 337             iend(i,2)=nct
 338 #endif
 339           endif
 340         else
 341 #ifdef MPL
 342           call int_partition(ind_scint,my_sc_inds,my_sc_inde,i,
 343      &    iatsc_s,iatsc_e,i+1,nct,nint_gr(i),istart(i,1),iend(i,1),*12)
 344 #else
 345           nint_gr(i)=1
 346           istart(i,1)=i+1
 347           iend(i,1)=nct
 348           ind_scint=int_scint+nct-i
 349 #endif
 350         endif
 351 #ifdef MPL
 352         ind_scint_old=ind_scint
 353 #endif
 354       enddo
 355    12 continue
 356 #ifndef MPL
 357       iatsc_s=nnt
 358       iatsc_e=nct-1
 359 #endif
 360 #ifdef MPL
 361       if (lprint) then
 362         write (iout,*) 'Processor',MyID,' Group',MyGroup
 363         write (iout,*) 'iatsc_s=',iatsc_s,' iatsc_e=',iatsc_e
 364       endif
 365 #endif
 366       if (lprint) then
 367       write (iout,'(a)') 'Interaction array:'
 368       do i=iatsc_s,iatsc_e
 369         write (iout,'(i3,2(2x,2i3))')
 370      & i,(istart(i,iint),iend(i,iint),iint=1,nint_gr(i))
 371       enddo
 372       endif
 373       ispp=2
 374 #ifdef MPL
 375 C Now partition the electrostatic-interaction array
 376       npept=nct-nnt
 377       nele_int_tot=(npept-ispp)*(npept-ispp+1)/2
 378       call int_bounds(nele_int_tot,my_ele_inds,my_ele_inde)
 379       if (lprint)
 380      & write (iout,*) 'Processor',MyID,' MyRank',MyRank,
 381      &  ' nele_int_tot',nele_int_tot,' my_ele_inds=',my_ele_inds,
 382      &               ' my_ele_inde',my_ele_inde
 383       iatel_s=0
 384       iatel_e=0
 385       ind_eleint=0
 386       ind_eleint_old=0
 387       do i=nnt,nct-3
 388         ijunk=0
 389         call int_partition(ind_eleint,my_ele_inds,my_ele_inde,i,
 390      &    iatel_s,iatel_e,i+ispp,nct-1,ijunk,ielstart(i),ielend(i),*13)
 391       enddo ! i
 392    13 continue
 393 #else
 394       iatel_s=nnt
 395       iatel_e=nct-3
 396       do i=iatel_s,iatel_e
 397         ielstart(i)=i+2
 398         ielend(i)=nct-1
 399       enddo
 400 #endif
 401       if (lprint) then
 402         write (iout,'(a)') 'Electrostatic interaction array:'
 403         do i=iatel_s,iatel_e
 404           write (iout,'(i3,2(2x,2i3))') i,ielstart(i),ielend(i)
 405         enddo
 406       endif ! lprint
 407 c     iscp=3
 408       iscp=2
 409 C Partition the SC-p interaction array
 410 #ifdef MPL
 411       nscp_int_tot=(npept-iscp+1)*(npept-iscp+1)
 412       call int_bounds(nscp_int_tot,my_scp_inds,my_scp_inde)
 413       if (lprint)
 414      & write (iout,*) 'Processor',MyID,' MyRank',MyRank,
 415      &  ' nscp_int_tot',nscp_int_tot,' my_scp_inds=',my_scp_inds,
 416      &               ' my_scp_inde',my_scp_inde
 417       iatscp_s=0
 418       iatscp_e=0
 419       ind_scpint=0
 420       ind_scpint_old=0
 421       do i=nnt,nct-1
 422         if (i.lt.nnt+iscp) then
 423 cd        write (iout,*) 'i.le.nnt+iscp'
 424           call int_partition(ind_scpint,my_scp_inds,my_scp_inde,i,
 425      &      iatscp_s,iatscp_e,i+iscp,nct,nscp_gr(i),iscpstart(i,1),
 426      &      iscpend(i,1),*14)
 427         else if (i.gt.nct-iscp) then
 428 cd        write (iout,*) 'i.gt.nct-iscp'
 429           call int_partition(ind_scpint,my_scp_inds,my_scp_inde,i,
 430      &      iatscp_s,iatscp_e,nnt,i-iscp,nscp_gr(i),iscpstart(i,1),
 431      &      iscpend(i,1),*14)
 432         else
 433           call int_partition(ind_scpint,my_scp_inds,my_scp_inde,i,
 434      &      iatscp_s,iatscp_e,nnt,i-iscp,nscp_gr(i),iscpstart(i,1),
 435      &      iscpend(i,1),*14)
 436           ii=nscp_gr(i)+1
 437           call int_partition(ind_scpint,my_scp_inds,my_scp_inde,i,
 438      &      iatscp_s,iatscp_e,i+iscp,nct,nscp_gr(i),iscpstart(i,ii),
 439      &      iscpend(i,ii),*14)
 440         endif
 441       enddo ! i
 442    14 continue
 443 #else
 444       iatscp_s=nnt
 445       iatscp_e=nct-1
 446       do i=nnt,nct-1
 447         if (i.lt.nnt+iscp) then
 448           nscp_gr(i)=1
 449           iscpstart(i,1)=i+iscp
 450           iscpend(i,1)=nct
 451         elseif (i.gt.nct-iscp) then
 452           nscp_gr(i)=1
 453           iscpstart(i,1)=nnt
 454           iscpend(i,1)=i-iscp
 455         else
 456           nscp_gr(i)=2
 457           iscpstart(i,1)=nnt
 458           iscpend(i,1)=i-iscp
 459           iscpstart(i,2)=i+iscp
 460           iscpend(i,2)=nct
 461         endif
 462       enddo ! i
 463 #endif
 464       if (lprint) then
 465         write (iout,'(a)') 'SC-p interaction array:'
 466         do i=iatscp_s,iatscp_e
 467           write (iout,'(i3,2(2x,2i3))')
 468      &         i,(iscpstart(i,j),iscpend(i,j),j=1,nscp_gr(i))
 469         enddo
 470       endif ! lprint
 471 C Partition local interactions
 472 #ifdef MPL
 473       call int_bounds(nres-2,loc_start,loc_end)
 474       loc_start=loc_start+1
 475       loc_end=loc_end+1
 476       call int_bounds(nres-2,ithet_start,ithet_end)
 477       ithet_start=ithet_start+2
 478       ithet_end=ithet_end+2
 479       call int_bounds(nct-nnt-2,iphi_start,iphi_end)
 480       iphi_start=iphi_start+nnt+2
 481       iphi_end=iphi_end+nnt+2
 482       call int_bounds(nres-3,itau_start,itau_end)
 483       itau_start=itau_start+3
 484       itau_end=itau_end+3
 485       if (lprint) then
 486         write (iout,*) 'Processor:',MyID,
 487      & ' loc_start',loc_start,' loc_end',loc_end,
 488      & ' ithet_start',ithet_start,' ithet_end',ithet_end,
 489      & ' iphi_start',iphi_start,' iphi_end',iphi_end
 490         write (*,*) 'Processor:',MyID,
 491      & ' loc_start',loc_start,' loc_end',loc_end,
 492      & ' ithet_start',ithet_start,' ithet_end',ithet_end,
 493      & ' iphi_start',iphi_start,' iphi_end',iphi_end
 494       endif
 495       if (fgprocs.gt.1 .and. MyID.eq.BossID) then
 496         write(iout,'(i10,a,i10,a,i10,a/a,i3,a)') n_sc_int_tot,' SC-SC ',
 497      & nele_int_tot,' electrostatic and ',nscp_int_tot,
 498      & ' SC-p interactions','were distributed among',fgprocs,
 499      & ' fine-grain processors.'
 500       endif
 501 #else
 502       loc_start=2
 503       loc_end=nres-1
 504       ithet_start=3
 505       ithet_end=nres
 506       iphi_start=nnt+3
 507       iphi_end=nct
 508       itau_start=4
 509       itau_end=nres
 510
 511 #endif
 512       return
 513       end
 514 c---------------------------------------------------------------------------
 515       subroutine int_partition(int_index,lower_index,upper_index,atom,
 516      & at_start,at_end,first_atom,last_atom,int_gr,jat_start,jat_end,*)
 517       implicit real*8 (a-h,o-z)
 518       include 'DIMENSIONS'
 519       include 'COMMON.IOUNITS'
 520       integer int_index,lower_index,upper_index,atom,at_start,at_end,
 521      & first_atom,last_atom,int_gr,jat_start,jat_end
 522       logical lprn
 523       lprn=.false.
 524       if (lprn) write (iout,*) 'int_index=',int_index
 525       int_index_old=int_index
 526       int_index=int_index+last_atom-first_atom+1
 527       if (lprn)
 528      &   write (iout,*) 'int_index=',int_index,
 529      &               ' int_index_old',int_index_old,
 530      &               ' lower_index=',lower_index,
 531      &               ' upper_index=',upper_index,
 532      &               ' atom=',atom,' first_atom=',first_atom,
 533      &               ' last_atom=',last_atom
 534       if (int_index.ge.lower_index) then
 535         int_gr=int_gr+1
 536         if (at_start.eq.0) then
 537           at_start=atom
 538           jat_start=first_atom-1+lower_index-int_index_old
 539         else
 540           jat_start=first_atom
 541         endif
 542         if (lprn) write (iout,*) 'jat_start',jat_start
 543         if (int_index.ge.upper_index) then
 544           at_end=atom
 545           jat_end=first_atom-1+upper_index-int_index_old
 546           return1
 547         else
 548           jat_end=last_atom
 549         endif
 550         if (lprn) write (iout,*) 'jat_end',jat_end
 551       endif
 552       return
 553       end
 554 c------------------------------------------------------------------------------
 555       subroutine hpb_partition
 556       implicit real*8 (a-h,o-z)
 557       include 'DIMENSIONS'
 558       include 'COMMON.SBRIDGE'
 559       include 'COMMON.IOUNITS'
 560 #ifdef MPL
 561       include 'COMMON.INFO'
 562       call int_bounds(nhpb,link_start,link_end)
 563 #else
 564       link_start=1
 565       link_end=nhpb
 566 #endif
 567 cd    write (iout,*) 'Processor',MyID,' MyRank',MyRank,
 568 cd   &  ' nhpb',nhpb,' link_start=',link_start,
 569 cd   &  ' link_end',link_end
 570       return
 571       end