source/unres/src_CSA_DiL/initialize_p.F

   1       block data
   2       implicit real*8 (a-h,o-z)
   3       include 'DIMENSIONS'
   4       include 'COMMON.MCM'
   5 c      include 'COMMON.MD'
   6       data MovTypID
   7      &  /'pool','chain regrow','multi-bond','phi','theta','side chain',
   8      &   'total'/
   9 c Conversion from poises to molecular unit and the gas constant
  10 c      data cPoise /2.9361d0/, Rb /0.001986d0/
  11       end
  12 c--------------------------------------------------------------------------
  13       subroutine initialize
  14 C
  15 C Define constants and zero out tables.
  16 C
  17       implicit real*8 (a-h,o-z)
  18       include 'DIMENSIONS'
  19 #ifdef MPI
  20       include 'mpif.h'
  21 #endif
  22 #ifndef ISNAN
  23       external proc_proc
  24 #ifdef WINPGI
  25 cMS$ATTRIBUTES C ::  proc_proc
  26 #endif
  27 #endif
  28       include 'COMMON.IOUNITS'
  29       include 'COMMON.CHAIN'
  30       include 'COMMON.INTERACT'
  31       include 'COMMON.GEO'
  32       include 'COMMON.LOCAL'
  33       include 'COMMON.TORSION'
  34       include 'COMMON.FFIELD'
  35       include 'COMMON.SBRIDGE'
  36       include 'COMMON.MCM'
  37       include 'COMMON.MINIM'
  38       include 'COMMON.DERIV'
  39       include 'COMMON.SPLITELE'
  40 c Common blocks from the diagonalization routines
  41       COMMON /IOFILE/ IR,IW,IP,IJK,IPK,IDAF,NAV,IODA(400)
  42       COMMON /MACHSW/ KDIAG,ICORFL,IXDR
  43       logical mask_r
  44 c      real*8 text1 /'initial_i'/
  45
  46       mask_r=.false.
  47 #ifndef ISNAN
  48 c NaNQ initialization
  49       i=-1
  50       arg=100.0d0
  51       rr=dacos(arg)
  52 #ifdef WINPGI
  53       idumm=proc_proc(rr,i)
  54 #else
  55       call proc_proc(rr,i)
  56 #endif
  57 #endif
  58
  59       kdiag=0
  60       icorfl=0
  61       iw=2
  62 C
  63 C The following is just to define auxiliary variables used in angle conversion
  64 C
  65       pi=4.0D0*datan(1.0D0)
  66       dwapi=2.0D0*pi
  67       dwapi3=dwapi/3.0D0
  68       pipol=0.5D0*pi
  69       deg2rad=pi/180.0D0
  70       rad2deg=1.0D0/deg2rad
  71       angmin=10.0D0*deg2rad
  72 C
  73 C Define I/O units.
  74 C
  75       inp=    1
  76       iout=   2
  77       ipdbin= 3
  78       ipdb=   7
  79       icart = 30
  80       imol2=  4
  81       igeom=  8
  82       intin=  9
  83       ithep= 11
  84       ithep_pdb=51
  85       irotam=12
  86       irotam_pdb=52
  87       itorp= 13
  88       itordp= 23
  89       ielep= 14
  90       isidep=15
  91       iscpp=25
  92       icbase=16
  93       ifourier=20
  94       istat= 17
  95       irest1=55
  96       irest2=56
  97       iifrag=57
  98       ientin=18
  99       ientout=19
 100       ibond = 28
 101       isccor = 29
 102 crc for write_rmsbank1
 103       izs1=21
 104 cdr  include secondary structure prediction bias
 105       isecpred=27
 106 C
 107 C CSA I/O units (separated from others especially for Jooyoung)
 108 C
 109       icsa_rbank=30
 110       icsa_seed=31
 111       icsa_history=32
 112       icsa_bank=33
 113       icsa_bank1=34
 114       icsa_alpha=35
 115       icsa_alpha1=36
 116       icsa_bankt=37
 117       icsa_int=39
 118       icsa_bank_reminimized=38
 119       icsa_native_int=41
 120       icsa_in=40
 121 crc for ifc error 118
 122       icsa_pdb=42
 123 C
 124 C Set default weights of the energy terms.
 125 C
 126       wlong=1.0D0
 127       welec=1.0D0
 128       wtor =1.0D0
 129       wang =1.0D0
 130       wscloc=1.0D0
 131       wstrain=1.0D0
 132 C
 133 C Zero out tables.
 134 C
 135       print '(a,$)','Inside initialize'
 136 c      call memmon_print_usage()
 137       do i=1,maxres2
 138         do j=1,3
 139           c(j,i)=0.0D0
 140           dc(j,i)=0.0D0
 141         enddo
 142       enddo
 143       do i=1,maxres
 144         do j=1,3
 145          xloc(j,i)=0.0D0
 146         enddo
 147       enddo
 148       do i=1,ntyp
 149         do j=1,ntyp
 150           aa(i,j)=0.0D0
 151            bb(i,j)=0.0D0
 152            augm(i,j)=0.0D0
 153            sigma(i,j)=0.0D0
 154            r0(i,j)=0.0D0
 155            chi(i,j)=0.0D0
 156         enddo
 157         do j=1,2
 158           bad(i,j)=0.0D0
 159         enddo
 160         chip(i)=0.0D0
 161         alp(i)=0.0D0
 162         sigma0(i)=0.0D0
 163         sigii(i)=0.0D0
 164         rr0(i)=0.0D0
 165         a0thet(i)=0.0D0
 166         do j=1,2
 167          do ichir1=-1,1
 168           do ichir2=-1,1
 169          athet(j,i,ichir1,ichir2)=0.0D0
 170          bthet(j,i,ichir1,ichir2)=0.0D0
 171           enddo
 172          enddo
 173         enddo
 174         do j=0,3
 175            polthet(j,i)=0.0D0
 176         enddo
 177         do j=1,3
 178           gthet(j,i)=0.0D0
 179         enddo
 180          theta0(i)=0.0D0
 181          sig0(i)=0.0D0
 182          sigc0(i)=0.0D0
 183          do j=1,maxlob
 184           bsc(j,i)=0.0D0
 185            do k=1,3
 186             censc(k,j,i)=0.0D0
 187           enddo
 188           do k=1,3
 189            do l=1,3
 190             gaussc(l,k,j,i)=0.0D0
 191             enddo
 192           enddo
 193             nlob(i)=0
 194         enddo
 195       enddo
 196       nlob(ntyp1)=0
 197       dsc(ntyp1)=0.0D0
 198       do i=-maxtor,maxtor
 199         itortyp(i)=0
 200          do iblock=1,2
 201          do j=-maxtor,maxtor
 202            do k=1,maxterm
 203             v1(k,j,i,iblock)=0.0D0
 204             v2(k,j,i,iblock)=0.0D0
 205           enddo
 206         enddo
 207       enddo
 208       enddo
 209       do iblock=1,2
 210        do i=-maxtor,maxtor
 211         do j=-maxtor,maxtor
 212          do k=-maxtor,maxtor
 213           do l=1,maxtermd_1
 214             v1c(1,l,i,j,k,iblock)=0.0D0
 215             v1s(1,l,i,j,k,iblock)=0.0D0
 216             v1c(2,l,i,j,k,iblock)=0.0D0
 217             v1s(2,l,i,j,k,iblock)=0.0D0
 218           enddo !l
 219           do l=1,maxtermd_2
 220            do m=1,maxtermd_2
 221             v2c(m,l,i,j,k,iblock)=0.0D0
 222             v2s(m,l,i,j,k,iblock)=0.0D0
 223            enddo !m
 224           enddo !l
 225         enddo !k
 226        enddo !j
 227       enddo !i
 228       enddo !i
 229       do i=1,maxres
 230         itype(i)=0
 231         itel(i)=0
 232       enddo
 233 C Initialize the bridge arrays
 234       ns=0
 235       nss=0
 236       nhpb=0
 237       do i=1,maxss
 238         iss(i)=0
 239       enddo
 240       do i=1,maxdim
 241        dhpb(i)=0.0D0
 242       enddo
 243       do i=1,maxres
 244         ihpb(i)=0
 245         jhpb(i)=0
 246       enddo
 247 C
 248 C Initialize timing.
 249 C
 250       call set_timers
 251 C
 252 C Initialize variables used in minimization.
 253 C
 254 c     maxfun=5000
 255 c     maxit=2000
 256       maxfun=500
 257       maxit=200
 258       tolf=1.0D-2
 259       rtolf=5.0D-4
 260 C
 261 C Initialize the variables responsible for the mode of gradient storage.
 262 C
 263       nfl=0
 264       icg=1
 265 C
 266 C Initialize constants used to split the energy into long- and short-range
 267 C components
 268 C
 269       r_cut=2.0d0
 270       rlamb=0.3d0
 271 #ifndef SPLITELE
 272       nprint_ene=nprint_ene-1
 273 #endif
 274       return
 275       end
 276 c-------------------------------------------------------------------------
 277       block data nazwy
 278       implicit real*8 (a-h,o-z)
 279       include 'DIMENSIONS'
 280       include 'COMMON.NAMES'
 281       include 'COMMON.FFIELD'
 282       data restyp /
 283      &'DD' ,'DPR','DLY','DAR','DHI','DAS','DGL','DSG','DGN','DSN','DTH',
 284      &'DYY','DAL','DTY','DTR','DVA','DLE','DIL','DPN','MED','DCY','ZER',
 285      &'CYS','MET','PHE','ILE','LEU','VAL','TRP','TYR','ALA','GLY','THR',
 286      &'SER','GLN','ASN','GLU','ASP','HIS','ARG','LYS','PRO','D'/
 287       data onelet /
 288      &'z','p','k','r','h','d','e','n','q','s','t','g',
 289      &'a','y','w','v','l','i','f','m','c','x',
 290      &'C','M','F','I','L','V','W','Y','A','G','T',
 291      &'S','Q','N','E','D','H','R','K','P','X'/
 292       data potname /'LJ','LJK','BP','GB','GBV'/
 293       data ename /
 294      &   "EVDW SC-SC","EVDW2 SC-p","EES p-p","ECORR4 ","ECORR5 ",
 295      &   "ECORR6 ","EELLO ","ETURN3 ","ETURN4 ","ETURN6 ",
 296      &   "EBE bend","ESC SCloc","ETORS ","ETORSD ","EHPB ","EVDWPP ",
 297      &   "ESTR ","EVDW2_14 ","UCONST ", "      ","ESCCOR"," "," ",
 298      &   "DFA DIS","DFA TOR","DFA NEI","DFA BET"/
 299       data wname /
 300      &   "WSC","WSCP","WELEC","WCORR","WCORR5","WCORR6","WEL_LOC",
 301      &   "WTURN3","WTURN4","WTURN6","WANG","WSCLOC","WTOR","WTORD",
 302      &   "WSTRAIN","WVDWPP","WBOND","SCAL14","     ","    ","WSCCOR",
 303      &   " "," ","WDFAD","WDFAT","WDFAN","WDFAB"/
 304       data nprint_ene /24/
 305       data print_order/1,2,3,11,12,13,14,4,5,6,7,8,9,10,19,18,15,17,16,
 306      & 21,24,25,26,27,0,0,0/
 307       end
 308 c---------------------------------------------------------------------------
 309       subroutine init_int_table
 310       implicit real*8 (a-h,o-z)
 311       include 'DIMENSIONS'
 312 #ifdef MPI
 313       include 'mpif.h'
 314       integer blocklengths(15),displs(15)
 315 #endif
 316       include 'COMMON.CONTROL'
 317       include 'COMMON.SETUP'
 318       include 'COMMON.CHAIN'
 319       include 'COMMON.INTERACT'
 320       include 'COMMON.LOCAL'
 321       include 'COMMON.SBRIDGE'
 322       include 'COMMON.TORCNSTR'
 323       include 'COMMON.IOUNITS'
 324       include 'COMMON.DERIV'
 325       include 'COMMON.CONTACTS'
 326       common /przechowalnia/ iturn3_start_all(0:max_fg_procs),
 327      & iturn3_end_all(0:max_fg_procs),iturn4_start_all(0:max_fg_procs),
 328      & iturn4_end_all(0:max_fg_procs),iatel_s_all(0:max_fg_procs),
 329      &iatel_e_all(0:max_fg_procs),ielstart_all(maxres,0:max_fg_procs-1),
 330      & ielend_all(maxres,0:max_fg_procs-1),
 331      & ntask_cont_from_all(0:max_fg_procs-1),
 332      & itask_cont_from_all(0:max_fg_procs-1,0:max_fg_procs-1),
 333      & ntask_cont_to_all(0:max_fg_procs-1),
 334      & itask_cont_to_all(0:max_fg_procs-1,0:max_fg_procs-1)
 335       integer FG_GROUP,CONT_FROM_GROUP,CONT_TO_GROUP
 336       logical scheck,lprint,flag
 337 #ifdef MPI
 338       integer my_sc_int(0:max_fg_Procs-1),my_sc_intt(0:max_fg_Procs),
 339      & my_ele_int(0:max_fg_Procs-1),my_ele_intt(0:max_fg_Procs)
 340 C... Determine the numbers of start and end SC-SC interaction
 341 C... to deal with by current processor.
 342       do i=0,nfgtasks-1
 343         itask_cont_from(i)=fg_rank
 344         itask_cont_to(i)=fg_rank
 345       enddo
 346       lprint=.false.
 347       if (lprint)
 348      &write (iout,*) 'INIT_INT_TABLE nres=',nres,' nnt=',nnt,' nct=',nct
 349       n_sc_int_tot=(nct-nnt+1)*(nct-nnt)/2-nss
 350       call int_bounds(n_sc_int_tot,my_sc_inds,my_sc_inde)
 351       if (lprint)
 352      &  write (iout,*) 'Processor',fg_rank,' CG group',kolor,
 353      &  ' absolute rank',MyRank,
 354      &  ' n_sc_int_tot',n_sc_int_tot,' my_sc_inds=',my_sc_inds,
 355      &  ' my_sc_inde',my_sc_inde
 356       ind_sctint=0
 357       iatsc_s=0
 358       iatsc_e=0
 359 #endif
 360 c      lprint=.false.
 361       do i=1,maxres
 362         nint_gr(i)=0
 363         nscp_gr(i)=0
 364         do j=1,maxint_gr
 365           istart(i,1)=0
 366           iend(i,1)=0
 367           ielstart(i)=0
 368           ielend(i)=0
 369           iscpstart(i,1)=0
 370           iscpend(i,1)=0
 371         enddo
 372       enddo
 373       ind_scint=0
 374       ind_scint_old=0
 375 cd    write (iout,*) 'ns=',ns,' nss=',nss,' ihpb,jhpb',
 376 cd   &   (ihpb(i),jhpb(i),i=1,nss)
 377       do i=nnt,nct-1
 378         scheck=.false.
 379         do ii=1,nss
 380           if (ihpb(ii).eq.i+nres) then
 381             scheck=.true.
 382             jj=jhpb(ii)-nres
 383             goto 10
 384           endif
 385         enddo
 386    10   continue
 387 cd      write (iout,*) 'i=',i,' scheck=',scheck,' jj=',jj
 388         if (scheck) then
 389           if (jj.eq.i+1) then
 390 #ifdef MPI
 391 c            write (iout,*) 'jj=i+1'
 392             call int_partition(ind_scint,my_sc_inds,my_sc_inde,i,
 393      & iatsc_s,iatsc_e,i+2,nct,nint_gr(i),istart(i,1),iend(i,1),*12)
 394 #else
 395             nint_gr(i)=1
 396             istart(i,1)=i+2
 397             iend(i,1)=nct
 398 #endif
 399           else if (jj.eq.nct) then
 400 #ifdef MPI
 401 c            write (iout,*) 'jj=nct'
 402             call int_partition(ind_scint,my_sc_inds,my_sc_inde,i,
 403      &  iatsc_s,iatsc_e,i+1,nct-1,nint_gr(i),istart(i,1),iend(i,1),*12)
 404 #else
 405             nint_gr(i)=1
 406             istart(i,1)=i+1
 407             iend(i,1)=nct-1
 408 #endif
 409           else
 410 #ifdef MPI
 411             call int_partition(ind_scint,my_sc_inds,my_sc_inde,i,
 412      & iatsc_s,iatsc_e,i+1,jj-1,nint_gr(i),istart(i,1),iend(i,1),*12)
 413             ii=nint_gr(i)+1
 414             call int_partition(ind_scint,my_sc_inds,my_sc_inde,i,
 415      & iatsc_s,iatsc_e,jj+1,nct,nint_gr(i),istart(i,ii),iend(i,ii),*12)
 416 #else
 417             nint_gr(i)=2
 418             istart(i,1)=i+1
 419             iend(i,1)=jj-1
 420             istart(i,2)=jj+1
 421             iend(i,2)=nct
 422 #endif
 423           endif
 424         else
 425 #ifdef MPI
 426           call int_partition(ind_scint,my_sc_inds,my_sc_inde,i,
 427      &    iatsc_s,iatsc_e,i+1,nct,nint_gr(i),istart(i,1),iend(i,1),*12)
 428 #else
 429           nint_gr(i)=1
 430           istart(i,1)=i+1
 431           iend(i,1)=nct
 432           ind_scint=ind_scint+nct-i
 433 #endif
 434         endif
 435 #ifdef MPI
 436         ind_scint_old=ind_scint
 437 #endif
 438       enddo
 439    12 continue
 440 #ifndef MPI
 441       iatsc_s=nnt
 442       iatsc_e=nct-1
 443 #endif
 444 #ifdef MPI
 445       if (lprint) write (*,*) 'Processor',fg_rank,' CG Group',kolor,
 446      &   ' absolute rank',myrank,' iatsc_s=',iatsc_s,' iatsc_e=',iatsc_e
 447 #endif
 448       if (lprint) then
 449       write (iout,'(a)') 'Interaction array:'
 450       do i=iatsc_s,iatsc_e
 451         write (iout,'(i3,2(2x,2i3))')
 452      & i,(istart(i,iint),iend(i,iint),iint=1,nint_gr(i))
 453       enddo
 454       endif
 455       ispp=4
 456 #ifdef MPI
 457 C Now partition the electrostatic-interaction array
 458       npept=nct-nnt
 459       nele_int_tot=(npept-ispp)*(npept-ispp+1)/2
 460       call int_bounds(nele_int_tot,my_ele_inds,my_ele_inde)
 461       if (lprint)
 462      & write (*,*) 'Processor',fg_rank,' CG group',kolor,
 463      &  ' absolute rank',MyRank,
 464      &  ' nele_int_tot',nele_int_tot,' my_ele_inds=',my_ele_inds,
 465      &               ' my_ele_inde',my_ele_inde
 466       iatel_s=0
 467       iatel_e=0
 468       ind_eleint=0
 469       ind_eleint_old=0
 470       do i=nnt,nct-3
 471         ijunk=0
 472         call int_partition(ind_eleint,my_ele_inds,my_ele_inde,i,
 473      &    iatel_s,iatel_e,i+ispp,nct-1,ijunk,ielstart(i),ielend(i),*13)
 474       enddo ! i
 475    13 continue
 476       if (iatel_s.eq.0) iatel_s=1
 477       nele_int_tot_vdw=(npept-2)*(npept-2+1)/2
 478 c      write (iout,*) "nele_int_tot_vdw",nele_int_tot_vdw
 479       call int_bounds(nele_int_tot_vdw,my_ele_inds_vdw,my_ele_inde_vdw)
 480 c      write (iout,*) "my_ele_inds_vdw",my_ele_inds_vdw,
 481 c     & " my_ele_inde_vdw",my_ele_inde_vdw
 482       ind_eleint_vdw=0
 483       ind_eleint_vdw_old=0
 484       iatel_s_vdw=0
 485       iatel_e_vdw=0
 486       do i=nnt,nct-3
 487         ijunk=0
 488         call int_partition(ind_eleint_vdw,my_ele_inds_vdw,
 489      &    my_ele_inde_vdw,i,
 490      &    iatel_s_vdw,iatel_e_vdw,i+2,nct-1,ijunk,ielstart_vdw(i),
 491      &    ielend_vdw(i),*15)
 492 c        write (iout,*) i," ielstart_vdw",ielstart_vdw(i),
 493 c     &   " ielend_vdw",ielend_vdw(i)
 494       enddo ! i
 495       if (iatel_s_vdw.eq.0) iatel_s_vdw=1
 496    15 continue
 497 #else
 498       iatel_s=nnt
 499       iatel_e=nct-5
 500       do i=iatel_s,iatel_e
 501         ielstart(i)=i+4
 502         ielend(i)=nct-1
 503       enddo
 504       iatel_s_vdw=nnt
 505       iatel_e_vdw=nct-3
 506       do i=iatel_s_vdw,iatel_e_vdw
 507         ielstart_vdw(i)=i+2
 508         ielend_vdw(i)=nct-1
 509       enddo
 510 #endif
 511       if (lprint) then
 512         write (*,'(a)') 'Processor',fg_rank,' CG group',kolor,
 513      &  ' absolute rank',MyRank
 514         write (iout,*) 'Electrostatic interaction array:'
 515         do i=iatel_s,iatel_e
 516           write (iout,'(i3,2(2x,2i3))') i,ielstart(i),ielend(i)
 517         enddo
 518       endif ! lprint
 519 c     iscp=3
 520       iscp=2
 521 C Partition the SC-p interaction array
 522 #ifdef MPI
 523       nscp_int_tot=(npept-iscp+1)*(npept-iscp+1)
 524       call int_bounds(nscp_int_tot,my_scp_inds,my_scp_inde)
 525       if (lprint) write (iout,*) 'Processor',fg_rank,' CG group',kolor,
 526      &  ' absolute rank',myrank,
 527      &  ' nscp_int_tot',nscp_int_tot,' my_scp_inds=',my_scp_inds,
 528      &               ' my_scp_inde',my_scp_inde
 529       iatscp_s=0
 530       iatscp_e=0
 531       ind_scpint=0
 532       ind_scpint_old=0
 533       do i=nnt,nct-1
 534         if (i.lt.nnt+iscp) then
 535 cd        write (iout,*) 'i.le.nnt+iscp'
 536           call int_partition(ind_scpint,my_scp_inds,my_scp_inde,i,
 537      &      iatscp_s,iatscp_e,i+iscp,nct,nscp_gr(i),iscpstart(i,1),
 538      &      iscpend(i,1),*14)
 539         else if (i.gt.nct-iscp) then
 540 cd        write (iout,*) 'i.gt.nct-iscp'
 541           call int_partition(ind_scpint,my_scp_inds,my_scp_inde,i,
 542      &      iatscp_s,iatscp_e,nnt,i-iscp,nscp_gr(i),iscpstart(i,1),
 543      &      iscpend(i,1),*14)
 544         else
 545           call int_partition(ind_scpint,my_scp_inds,my_scp_inde,i,
 546      &      iatscp_s,iatscp_e,nnt,i-iscp,nscp_gr(i),iscpstart(i,1),
 547      &      iscpend(i,1),*14)
 548           ii=nscp_gr(i)+1
 549           call int_partition(ind_scpint,my_scp_inds,my_scp_inde,i,
 550      &      iatscp_s,iatscp_e,i+iscp,nct,nscp_gr(i),iscpstart(i,ii),
 551      &      iscpend(i,ii),*14)
 552         endif
 553       enddo ! i
 554    14 continue
 555 #else
 556       iatscp_s=nnt
 557       iatscp_e=nct-1
 558       do i=nnt,nct-1
 559         if (i.lt.nnt+iscp) then
 560           nscp_gr(i)=1
 561           iscpstart(i,1)=i+iscp
 562           iscpend(i,1)=nct
 563         elseif (i.gt.nct-iscp) then
 564           nscp_gr(i)=1
 565           iscpstart(i,1)=nnt
 566           iscpend(i,1)=i-iscp
 567         else
 568           nscp_gr(i)=2
 569           iscpstart(i,1)=nnt
 570           iscpend(i,1)=i-iscp
 571           iscpstart(i,2)=i+iscp
 572           iscpend(i,2)=nct
 573         endif
 574       enddo ! i
 575 #endif
 576       if (lprint) then
 577         write (iout,'(a)') 'SC-p interaction array:'
 578         do i=iatscp_s,iatscp_e
 579           write (iout,'(i3,2(2x,2i3))')
 580      &         i,(iscpstart(i,j),iscpend(i,j),j=1,nscp_gr(i))
 581         enddo
 582       endif ! lprint
 583 C Partition local interactions
 584 #ifdef MPI
 585       call int_bounds(nres-2,loc_start,loc_end)
 586       loc_start=loc_start+1
 587       loc_end=loc_end+1
 588       call int_bounds(nres-2,ithet_start,ithet_end)
 589       ithet_start=ithet_start+2
 590       ithet_end=ithet_end+2
 591       call int_bounds(nct-nnt-2,iturn3_start,iturn3_end)
 592       iturn3_start=iturn3_start+nnt
 593       iphi_start=iturn3_start+2
 594       iturn3_end=iturn3_end+nnt
 595       iphi_end=iturn3_end+2
 596       iturn3_start=iturn3_start-1
 597       iturn3_end=iturn3_end-1
 598       call int_bounds(nres-3,iphi1_start,iphi1_end)
 599       iphi1_start=iphi1_start+3
 600       iphi1_end=iphi1_end+3
 601       call int_bounds(nct-nnt-3,iturn4_start,iturn4_end)
 602       iturn4_start=iturn4_start+nnt
 603       iphid_start=iturn4_start+2
 604       iturn4_end=iturn4_end+nnt
 605       iphid_end=iturn4_end+2
 606       iturn4_start=iturn4_start-1
 607       iturn4_end=iturn4_end-1
 608       call int_bounds(nres-2,ibond_start,ibond_end)
 609       ibond_start=ibond_start+1
 610       ibond_end=ibond_end+1
 611       call int_bounds(nct-nnt,ibondp_start,ibondp_end)
 612       ibondp_start=ibondp_start+nnt
 613       ibondp_end=ibondp_end+nnt
 614       call int_bounds1(nres-1,ivec_start,ivec_end)
 615       print *,"Processor",myrank,fg_rank,fg_rank1,
 616      &  " ivec_start",ivec_start," ivec_end",ivec_end
 617       iset_start=loc_start+2
 618       iset_end=loc_end+2
 619       if (ndih_constr.eq.0) then
 620         idihconstr_start=1
 621         idihconstr_end=0
 622       else
 623         call int_bounds(ndih_constr,idihconstr_start,idihconstr_end)
 624       endif
 625       nsumgrad=(nres-nnt)*(nres-nnt+1)/2
 626       nlen=nres-nnt+1
 627       call int_bounds(nsumgrad,ngrad_start,ngrad_end)
 628       igrad_start=((2*nlen+1)
 629      &    -sqrt(float((2*nlen-1)**2-8*(ngrad_start-1))))/2
 630       jgrad_start(igrad_start)=
 631      &    ngrad_start-(2*nlen-igrad_start)*(igrad_start-1)/2
 632      &    +igrad_start
 633       jgrad_end(igrad_start)=nres
 634       igrad_end=((2*nlen+1)
 635      &    -sqrt(float((2*nlen-1)**2-8*(ngrad_end-1))))/2
 636       if (igrad_end.gt.igrad_start) jgrad_start(igrad_end)=igrad_end+1
 637       jgrad_end(igrad_end)=ngrad_end-(2*nlen-igrad_end)*(igrad_end-1)/2
 638      &    +igrad_end
 639       do i=igrad_start+1,igrad_end-1
 640         jgrad_start(i)=i+1
 641         jgrad_end(i)=nres
 642       enddo
 643       if (lprint) then
 644         write (*,*) 'Processor:',fg_rank,' CG group',kolor,
 645      & ' absolute rank',myrank,
 646      & ' loc_start',loc_start,' loc_end',loc_end,
 647      & ' ithet_start',ithet_start,' ithet_end',ithet_end,
 648      & ' iphi_start',iphi_start,' iphi_end',iphi_end,
 649      & ' iphid_start',iphid_start,' iphid_end',iphid_end,
 650      & ' ibond_start',ibond_start,' ibond_end',ibond_end,
 651      & ' ibondp_start',ibondp_start,' ibondp_end',ibondp_end,
 652      & ' iturn3_start',iturn3_start,' iturn3_end',iturn3_end,
 653      & ' iturn4_start',iturn4_start,' iturn4_end',iturn4_end,
 654      & ' ivec_start',ivec_start,' ivec_end',ivec_end,
 655      & ' iset_start',iset_start,' iset_end',iset_end,
 656      & ' idihconstr_start',idihconstr_start,' idihconstr_end',
 657      &   idihconstr_end
 658        write (*,*) 'Processor:',fg_rank,myrank,' igrad_start',
 659      &   igrad_start,' igrad_end',igrad_end,' ngrad_start',ngrad_start,
 660      &   ' ngrad_end',ngrad_end
 661        do i=igrad_start,igrad_end
 662          write(*,*) 'Processor:',fg_rank,myrank,i,
 663      &    jgrad_start(i),jgrad_end(i)
 664        enddo
 665       endif
 666       if (nfgtasks.gt.1) then
 667         call MPI_Allgather(ivec_start,1,MPI_INTEGER,ivec_displ(0),1,
 668      &    MPI_INTEGER,FG_COMM1,IERROR)
 669         iaux=ivec_end-ivec_start+1
 670         call MPI_Allgather(iaux,1,MPI_INTEGER,ivec_count(0),1,
 671      &    MPI_INTEGER,FG_COMM1,IERROR)
 672         call MPI_Allgather(iset_start-2,1,MPI_INTEGER,iset_displ(0),1,
 673      &    MPI_INTEGER,FG_COMM,IERROR)
 674         iaux=iset_end-iset_start+1
 675         call MPI_Allgather(iaux,1,MPI_INTEGER,iset_count(0),1,
 676      &    MPI_INTEGER,FG_COMM,IERROR)
 677         call MPI_Allgather(ibond_start,1,MPI_INTEGER,ibond_displ(0),1,
 678      &    MPI_INTEGER,FG_COMM,IERROR)
 679         iaux=ibond_end-ibond_start+1
 680         call MPI_Allgather(iaux,1,MPI_INTEGER,ibond_count(0),1,
 681      &    MPI_INTEGER,FG_COMM,IERROR)
 682         call MPI_Allgather(ithet_start,1,MPI_INTEGER,ithet_displ(0),1,
 683      &    MPI_INTEGER,FG_COMM,IERROR)
 684         iaux=ithet_end-ithet_start+1
 685         call MPI_Allgather(iaux,1,MPI_INTEGER,ithet_count(0),1,
 686      &    MPI_INTEGER,FG_COMM,IERROR)
 687         call MPI_Allgather(iphi_start,1,MPI_INTEGER,iphi_displ(0),1,
 688      &    MPI_INTEGER,FG_COMM,IERROR)
 689         iaux=iphi_end-iphi_start+1
 690         call MPI_Allgather(iaux,1,MPI_INTEGER,iphi_count(0),1,
 691      &    MPI_INTEGER,FG_COMM,IERROR)
 692         call MPI_Allgather(iphi1_start,1,MPI_INTEGER,iphi1_displ(0),1,
 693      &    MPI_INTEGER,FG_COMM,IERROR)
 694         iaux=iphi1_end-iphi1_start+1
 695         call MPI_Allgather(iaux,1,MPI_INTEGER,iphi1_count(0),1,
 696      &    MPI_INTEGER,FG_COMM,IERROR)
 697         do i=0,maxprocs-1
 698           do j=1,maxres
 699             ielstart_all(j,i)=0
 700             ielend_all(j,i)=0
 701           enddo
 702         enddo
 703         call MPI_Allgather(iturn3_start,1,MPI_INTEGER,
 704      &    iturn3_start_all(0),1,MPI_INTEGER,FG_COMM,IERROR)
 705         call MPI_Allgather(iturn4_start,1,MPI_INTEGER,
 706      &    iturn4_start_all(0),1,MPI_INTEGER,FG_COMM,IERROR)
 707         call MPI_Allgather(iturn3_end,1,MPI_INTEGER,
 708      &    iturn3_end_all(0),1,MPI_INTEGER,FG_COMM,IERROR)
 709         call MPI_Allgather(iturn4_end,1,MPI_INTEGER,
 710      &    iturn4_end_all(0),1,MPI_INTEGER,FG_COMM,IERROR)
 711         call MPI_Allgather(iatel_s,1,MPI_INTEGER,
 712      &    iatel_s_all(0),1,MPI_INTEGER,FG_COMM,IERROR)
 713         call MPI_Allgather(iatel_e,1,MPI_INTEGER,
 714      &    iatel_e_all(0),1,MPI_INTEGER,FG_COMM,IERROR)
 715         call MPI_Allgather(ielstart(1),maxres,MPI_INTEGER,
 716      &    ielstart_all(1,0),maxres,MPI_INTEGER,FG_COMM,IERROR)
 717         call MPI_Allgather(ielend(1),maxres,MPI_INTEGER,
 718      &    ielend_all(1,0),maxres,MPI_INTEGER,FG_COMM,IERROR)
 719         if (lprint) then
 720         write (iout,*) "iatel_s_all",(iatel_s_all(i),i=0,nfgtasks)
 721         write (iout,*) "iatel_e_all",(iatel_e_all(i),i=0,nfgtasks)
 722         write (iout,*) "iturn3_start_all",
 723      &    (iturn3_start_all(i),i=0,nfgtasks-1)
 724         write (iout,*) "iturn3_end_all",
 725      &    (iturn3_end_all(i),i=0,nfgtasks-1)
 726         write (iout,*) "iturn4_start_all",
 727      &    (iturn4_start_all(i),i=0,nfgtasks-1)
 728         write (iout,*) "iturn4_end_all",
 729      &    (iturn4_end_all(i),i=0,nfgtasks-1)
 730         write (iout,*) "The ielstart_all array"
 731         do i=nnt,nct
 732           write (iout,'(20i4)') i,(ielstart_all(i,j),j=0,nfgtasks-1)
 733         enddo
 734         write (iout,*) "The ielend_all array"
 735         do i=nnt,nct
 736           write (iout,'(20i4)') i,(ielend_all(i,j),j=0,nfgtasks-1)
 737         enddo
 738         call flush(iout)
 739         endif
 740         ntask_cont_from=0
 741         ntask_cont_to=0
 742         itask_cont_from(0)=fg_rank
 743         itask_cont_to(0)=fg_rank
 744         flag=.false.
 745         do ii=iturn3_start,iturn3_end
 746           call add_int(ii,ii+2,iturn3_sent(1,ii),
 747      &                 ntask_cont_to,itask_cont_to,flag)
 748         enddo
 749         do ii=iturn4_start,iturn4_end
 750           call add_int(ii,ii+3,iturn4_sent(1,ii),
 751      &                 ntask_cont_to,itask_cont_to,flag)
 752         enddo
 753         do ii=iturn3_start,iturn3_end
 754           call add_int_from(ii,ii+2,ntask_cont_from,itask_cont_from)
 755         enddo
 756         do ii=iturn4_start,iturn4_end
 757           call add_int_from(ii,ii+3,ntask_cont_from,itask_cont_from)
 758         enddo
 759         if (lprint) then
 760         write (iout,*) "After turn3 ntask_cont_from",ntask_cont_from,
 761      &   " ntask_cont_to",ntask_cont_to
 762         write (iout,*) "itask_cont_from",
 763      &    (itask_cont_from(i),i=1,ntask_cont_from)
 764         write (iout,*) "itask_cont_to",
 765      &    (itask_cont_to(i),i=1,ntask_cont_to)
 766         call flush(iout)
 767         endif
 768 c        write (iout,*) "Loop forward"
 769 c        call flush(iout)
 770         do i=iatel_s,iatel_e
 771 c          write (iout,*) "from loop i=",i
 772 c          call flush(iout)
 773           do j=ielstart(i),ielend(i)
 774             call add_int_from(i,j,ntask_cont_from,itask_cont_from)
 775           enddo
 776         enddo
 777 c        write (iout,*) "Loop backward iatel_e-1",iatel_e-1,
 778 c     &     " iatel_e",iatel_e
 779 c        call flush(iout)
 780         nat_sent=0
 781         do i=iatel_s,iatel_e
 782 c          write (iout,*) "i",i," ielstart",ielstart(i),
 783 c     &      " ielend",ielend(i)
 784 c          call flush(iout)
 785           flag=.false.
 786           do j=ielstart(i),ielend(i)
 787             call add_int(i,j,iint_sent(1,j,nat_sent+1),ntask_cont_to,
 788      &                   itask_cont_to,flag)
 789           enddo
 790           if (flag) then
 791             nat_sent=nat_sent+1
 792             iat_sent(nat_sent)=i
 793           endif
 794         enddo
 795         if (lprint) then
 796         write (iout,*)"After longrange ntask_cont_from",ntask_cont_from,
 797      &   " ntask_cont_to",ntask_cont_to
 798         write (iout,*) "itask_cont_from",
 799      &    (itask_cont_from(i),i=1,ntask_cont_from)
 800         write (iout,*) "itask_cont_to",
 801      &    (itask_cont_to(i),i=1,ntask_cont_to)
 802         call flush(iout)
 803         write (iout,*) "iint_sent"
 804         do i=1,nat_sent
 805           ii=iat_sent(i)
 806           write (iout,'(20i4)') ii,(j,(iint_sent(k,j,i),k=1,4),
 807      &      j=ielstart(ii),ielend(ii))
 808         enddo
 809         write (iout,*) "iturn3_sent iturn3_start",iturn3_start,
 810      &    " iturn3_end",iturn3_end
 811         write (iout,'(20i4)') (i,(iturn3_sent(j,i),j=1,4),
 812      &      i=iturn3_start,iturn3_end)
 813         write (iout,*) "iturn4_sent iturn4_start",iturn4_start,
 814      &    " iturn4_end",iturn4_end
 815         write (iout,'(20i4)') (i,(iturn4_sent(j,i),j=1,4),
 816      &      i=iturn4_start,iturn4_end)
 817         call flush(iout)
 818         endif
 819         call MPI_Gather(ntask_cont_from,1,MPI_INTEGER,
 820      &   ntask_cont_from_all,1,MPI_INTEGER,king,FG_COMM,IERR)
 821 c        write (iout,*) "Gather ntask_cont_from ended"
 822 c        call flush(iout)
 823         call MPI_Gather(itask_cont_from(0),max_fg_procs,MPI_INTEGER,
 824      &   itask_cont_from_all(0,0),max_fg_procs,MPI_INTEGER,king,
 825      &   FG_COMM,IERR)
 826 c        write (iout,*) "Gather itask_cont_from ended"
 827 c        call flush(iout)
 828         call MPI_Gather(ntask_cont_to,1,MPI_INTEGER,ntask_cont_to_all,
 829      &   1,MPI_INTEGER,king,FG_COMM,IERR)
 830 c        write (iout,*) "Gather ntask_cont_to ended"
 831 c        call flush(iout)
 832         call MPI_Gather(itask_cont_to,max_fg_procs,MPI_INTEGER,
 833      &   itask_cont_to_all,max_fg_procs,MPI_INTEGER,king,FG_COMM,IERR)
 834 c        write (iout,*) "Gather itask_cont_to ended"
 835 c        call flush(iout)
 836         if (fg_rank.eq.king) then
 837           write (iout,*)"Contact receive task map (proc, #tasks, tasks)"
 838           do i=0,nfgtasks-1
 839             write (iout,'(20i4)') i,ntask_cont_from_all(i),
 840      &       (itask_cont_from_all(j,i),j=1,ntask_cont_from_all(i))
 841           enddo
 842           write (iout,*)
 843           call flush(iout)
 844           write (iout,*) "Contact send task map (proc, #tasks, tasks)"
 845           do i=0,nfgtasks-1
 846             write (iout,'(20i4)') i,ntask_cont_to_all(i),
 847      &       (itask_cont_to_all(j,i),j=1,ntask_cont_to_all(i))
 848           enddo
 849           write (iout,*)
 850           call flush(iout)
 851 C Check if every send will have a matching receive
 852           ncheck_to=0
 853           ncheck_from=0
 854           do i=0,nfgtasks-1
 855             ncheck_to=ncheck_to+ntask_cont_to_all(i)
 856             ncheck_from=ncheck_from+ntask_cont_from_all(i)
 857           enddo
 858           write (iout,*) "Control sums",ncheck_from,ncheck_to
 859           if (ncheck_from.ne.ncheck_to) then
 860             write (iout,*) "Error: #receive differs from #send."
 861             write (iout,*) "Terminating program...!"
 862             call flush(iout)
 863             flag=.false.
 864           else
 865             flag=.true.
 866             do i=0,nfgtasks-1
 867               do j=1,ntask_cont_to_all(i)
 868                 ii=itask_cont_to_all(j,i)
 869                 do k=1,ntask_cont_from_all(ii)
 870                   if (itask_cont_from_all(k,ii).eq.i) then
 871                     if(lprint)write(iout,*)"Matching send/receive",i,ii
 872                     exit
 873                   endif
 874                 enddo
 875                 if (k.eq.ntask_cont_from_all(ii)+1) then
 876                   flag=.false.
 877                   write (iout,*) "Error: send by",j," to",ii,
 878      &            " would have no matching receive"
 879                 endif
 880               enddo
 881             enddo
 882           endif
 883           if (.not.flag) then
 884             write (iout,*) "Unmatched sends; terminating program"
 885             call flush(iout)
 886           endif
 887         endif
 888         call MPI_Bcast(flag,1,MPI_LOGICAL,king,FG_COMM,IERROR)
 889 c        write (iout,*) "flag broadcast ended flag=",flag
 890 c        call flush(iout)
 891         if (.not.flag) then
 892           call MPI_Finalize(IERROR)
 893           stop "Error in INIT_INT_TABLE: unmatched send/receive."
 894         endif
 895         call MPI_Comm_group(FG_COMM,fg_group,IERR)
 896 c        write (iout,*) "MPI_Comm_group ended"
 897 c        call flush(iout)
 898         call MPI_Group_incl(fg_group,ntask_cont_from+1,
 899      &    itask_cont_from(0),CONT_FROM_GROUP,IERR)
 900         call MPI_Group_incl(fg_group,ntask_cont_to+1,itask_cont_to(0),
 901      &    CONT_TO_GROUP,IERR)
 902         do i=1,nat_sent
 903           ii=iat_sent(i)
 904           iaux=4*(ielend(ii)-ielstart(ii)+1)
 905           call MPI_Group_translate_ranks(fg_group,iaux,
 906      &      iint_sent(1,ielstart(ii),i),CONT_TO_GROUP,
 907      &      iint_sent_local(1,ielstart(ii),i),IERR )
 908 c          write (iout,*) "Ranks translated i=",i
 909 c          call flush(iout)
 910         enddo
 911         iaux=4*(iturn3_end-iturn3_start+1)
 912         call MPI_Group_translate_ranks(fg_group,iaux,
 913      &      iturn3_sent(1,iturn3_start),CONT_TO_GROUP,
 914      &      iturn3_sent_local(1,iturn3_start),IERR)
 915         iaux=4*(iturn4_end-iturn4_start+1)
 916         call MPI_Group_translate_ranks(fg_group,iaux,
 917      &      iturn4_sent(1,iturn4_start),CONT_TO_GROUP,
 918      &      iturn4_sent_local(1,iturn4_start),IERR)
 919         if (lprint) then
 920         write (iout,*) "iint_sent_local"
 921         do i=1,nat_sent
 922           ii=iat_sent(i)
 923           write (iout,'(20i4)') ii,(j,(iint_sent_local(k,j,i),k=1,4),
 924      &      j=ielstart(ii),ielend(ii))
 925           call flush(iout)
 926         enddo
 927         write (iout,*) "iturn3_sent_local iturn3_start",iturn3_start,
 928      &    " iturn3_end",iturn3_end
 929         write (iout,'(20i4)') (i,(iturn3_sent_local(j,i),j=1,4),
 930      &      i=iturn3_start,iturn3_end)
 931         write (iout,*) "iturn4_sent_local iturn4_start",iturn4_start,
 932      &    " iturn4_end",iturn4_end
 933         write (iout,'(20i4)') (i,(iturn4_sent_local(j,i),j=1,4),
 934      &      i=iturn4_start,iturn4_end)
 935         call flush(iout)
 936         endif
 937         call MPI_Group_free(fg_group,ierr)
 938         call MPI_Group_free(cont_from_group,ierr)
 939         call MPI_Group_free(cont_to_group,ierr)
 940         call MPI_Type_contiguous(3,MPI_DOUBLE_PRECISION,MPI_UYZ,IERROR)
 941         call MPI_Type_commit(MPI_UYZ,IERROR)
 942         call MPI_Type_contiguous(18,MPI_DOUBLE_PRECISION,MPI_UYZGRAD,
 943      &    IERROR)
 944         call MPI_Type_commit(MPI_UYZGRAD,IERROR)
 945         call MPI_Type_contiguous(2,MPI_DOUBLE_PRECISION,MPI_MU,IERROR)
 946         call MPI_Type_commit(MPI_MU,IERROR)
 947         call MPI_Type_contiguous(4,MPI_DOUBLE_PRECISION,MPI_MAT1,IERROR)
 948         call MPI_Type_commit(MPI_MAT1,IERROR)
 949         call MPI_Type_contiguous(8,MPI_DOUBLE_PRECISION,MPI_MAT2,IERROR)
 950         call MPI_Type_commit(MPI_MAT2,IERROR)
 951         call MPI_Type_contiguous(6,MPI_DOUBLE_PRECISION,MPI_THET,IERROR)
 952         call MPI_Type_commit(MPI_THET,IERROR)
 953         call MPI_Type_contiguous(9,MPI_DOUBLE_PRECISION,MPI_GAM,IERROR)
 954         call MPI_Type_commit(MPI_GAM,IERROR)
 955 #ifndef MATGATHER
 956 c 9/22/08 Derived types to send matrices which appear in correlation terms
 957         do i=0,nfgtasks-1
 958           if (ivec_count(i).eq.ivec_count(0)) then
 959             lentyp(i)=0
 960           else
 961             lentyp(i)=1
 962           endif
 963         enddo
 964         do ind_typ=lentyp(0),lentyp(nfgtasks-1)
 965         if (ind_typ.eq.0) then
 966           ichunk=ivec_count(0)
 967         else
 968           ichunk=ivec_count(1)
 969         endif
 970 c        do i=1,4
 971 c          blocklengths(i)=4
 972 c        enddo
 973 c        displs(1)=0
 974 c        do i=2,4
 975 c          displs(i)=displs(i-1)+blocklengths(i-1)*maxres
 976 c        enddo
 977 c        do i=1,4
 978 c          blocklengths(i)=blocklengths(i)*ichunk
 979 c        enddo
 980 c        write (iout,*) "blocklengths and displs"
 981 c        do i=1,4
 982 c          write (iout,*) i,blocklengths(i),displs(i)
 983 c        enddo
 984 c        call flush(iout)
 985 c        call MPI_Type_indexed(4,blocklengths(1),displs(1),
 986 c     &    MPI_DOUBLE_PRECISION,MPI_ROTAT1(ind_typ),IERROR)
 987 c        call MPI_Type_commit(MPI_ROTAT1(ind_typ),IERROR)
 988 c        write (iout,*) "MPI_ROTAT1",MPI_ROTAT1
 989 c        do i=1,4
 990 c          blocklengths(i)=2
 991 c        enddo
 992 c        displs(1)=0
 993 c        do i=2,4
 994 c          displs(i)=displs(i-1)+blocklengths(i-1)*maxres
 995 c        enddo
 996 c        do i=1,4
 997 c          blocklengths(i)=blocklengths(i)*ichunk
 998 c        enddo
 999 c        write (iout,*) "blocklengths and displs"
1000 c        do i=1,4
1001 c          write (iout,*) i,blocklengths(i),displs(i)
1002 c        enddo
1003 c        call flush(iout)
1004 c        call MPI_Type_indexed(4,blocklengths(1),displs(1),
1005 c     &    MPI_DOUBLE_PRECISION,MPI_ROTAT2(ind_typ),IERROR)
1006 c        call MPI_Type_commit(MPI_ROTAT2(ind_typ),IERROR)
1007 c        write (iout,*) "MPI_ROTAT2",MPI_ROTAT2
1008         do i=1,8
1009           blocklengths(i)=2
1010         enddo
1011         displs(1)=0
1012         do i=2,8
1013           displs(i)=displs(i-1)+blocklengths(i-1)*maxres
1014         enddo
1015         do i=1,15
1016           blocklengths(i)=blocklengths(i)*ichunk
1017         enddo
1018         call MPI_Type_indexed(8,blocklengths,displs,
1019      &    MPI_DOUBLE_PRECISION,MPI_PRECOMP11(ind_typ),IERROR)
1020         call MPI_Type_commit(MPI_PRECOMP11(ind_typ),IERROR)
1021         do i=1,8
1022           blocklengths(i)=4
1023         enddo
1024         displs(1)=0
1025         do i=2,8
1026           displs(i)=displs(i-1)+blocklengths(i-1)*maxres
1027         enddo
1028         do i=1,15
1029           blocklengths(i)=blocklengths(i)*ichunk
1030         enddo
1031         call MPI_Type_indexed(8,blocklengths,displs,
1032      &    MPI_DOUBLE_PRECISION,MPI_PRECOMP12(ind_typ),IERROR)
1033         call MPI_Type_commit(MPI_PRECOMP12(ind_typ),IERROR)
1034         do i=1,6
1035           blocklengths(i)=4
1036         enddo
1037         displs(1)=0
1038         do i=2,6
1039           displs(i)=displs(i-1)+blocklengths(i-1)*maxres
1040         enddo
1041         do i=1,6
1042           blocklengths(i)=blocklengths(i)*ichunk
1043         enddo
1044         call MPI_Type_indexed(6,blocklengths,displs,
1045      &    MPI_DOUBLE_PRECISION,MPI_PRECOMP22(ind_typ),IERROR)
1046         call MPI_Type_commit(MPI_PRECOMP22(ind_typ),IERROR)
1047         do i=1,2
1048           blocklengths(i)=8
1049         enddo
1050         displs(1)=0
1051         do i=2,2
1052           displs(i)=displs(i-1)+blocklengths(i-1)*maxres
1053         enddo
1054         do i=1,2
1055           blocklengths(i)=blocklengths(i)*ichunk
1056         enddo
1057         call MPI_Type_indexed(2,blocklengths,displs,
1058      &    MPI_DOUBLE_PRECISION,MPI_PRECOMP23(ind_typ),IERROR)
1059         call MPI_Type_commit(MPI_PRECOMP23(ind_typ),IERROR)
1060         do i=1,4
1061           blocklengths(i)=1
1062         enddo
1063         displs(1)=0
1064         do i=2,4
1065           displs(i)=displs(i-1)+blocklengths(i-1)*maxres
1066         enddo
1067         do i=1,4
1068           blocklengths(i)=blocklengths(i)*ichunk
1069         enddo
1070         call MPI_Type_indexed(4,blocklengths,displs,
1071      &    MPI_DOUBLE_PRECISION,MPI_ROTAT_OLD(ind_typ),IERROR)
1072         call MPI_Type_commit(MPI_ROTAT_OLD(ind_typ),IERROR)
1073         enddo
1074 #endif
1075       endif
1076       iint_start=ivec_start+1
1077       iint_end=ivec_end+1
1078       do i=0,nfgtasks-1
1079           iint_count(i)=ivec_count(i)
1080           iint_displ(i)=ivec_displ(i)
1081           ivec_displ(i)=ivec_displ(i)-1
1082           iset_displ(i)=iset_displ(i)-1
1083           ithet_displ(i)=ithet_displ(i)-1
1084           iphi_displ(i)=iphi_displ(i)-1
1085           iphi1_displ(i)=iphi1_displ(i)-1
1086           ibond_displ(i)=ibond_displ(i)-1
1087       enddo
1088       if (nfgtasks.gt.1 .and. fg_rank.eq.king
1089      &     .and. (me.eq.0 .or. out1file)) then
1090         write (iout,*) "IVEC_DISPL, IVEC_COUNT, ISET_START, ISET_COUNT"
1091         do i=0,nfgtasks-1
1092           write (iout,*) i,ivec_displ(i),ivec_count(i),iset_displ(i),
1093      &      iset_count(i)
1094         enddo
1095         write (iout,*) "iphi_start",iphi_start," iphi_end",iphi_end,
1096      &    " iphi1_start",iphi1_start," iphi1_end",iphi1_end
1097         write (iout,*)"IPHI_COUNT, IPHI_DISPL, IPHI1_COUNT, IPHI1_DISPL"
1098         do i=0,nfgtasks-1
1099           write (iout,*) i,iphi_count(i),iphi_displ(i),iphi1_count(i),
1100      &      iphi1_displ(i)
1101         enddo
1102         write(iout,'(i10,a,i10,a,i10,a/a,i3,a)') n_sc_int_tot,' SC-SC ',
1103      & nele_int_tot,' electrostatic and ',nscp_int_tot,
1104      & ' SC-p interactions','were distributed among',nfgtasks,
1105      & ' fine-grain processors.'
1106       endif
1107 #else
1108       loc_start=2
1109       loc_end=nres-1
1110       ithet_start=3
1111       ithet_end=nres
1112       iturn3_start=nnt
1113       iturn3_end=nct-3
1114       iturn4_start=nnt
1115       iturn4_end=nct-4
1116       iphi_start=nnt+3
1117       iphi_end=nct
1118       iphi1_start=4
1119       iphi1_end=nres
1120       idihconstr_start=1
1121       idihconstr_end=ndih_constr
1122       iphid_start=iphi_start
1123       iphid_end=iphi_end-1
1124       ibond_start=2
1125       ibond_end=nres-1
1126       ibondp_start=nnt+1
1127       ibondp_end=nct
1128       ivec_start=1
1129       ivec_end=nres-1
1130       iset_start=3
1131       iset_end=nres+1
1132       iint_start=2
1133       iint_end=nres-1
1134 #endif
1135       return
1136       end
1137 #ifdef MPI
1138 c---------------------------------------------------------------------------
1139       subroutine add_int(ii,jj,itask,ntask_cont_to,itask_cont_to,flag)
1140       implicit none
1141       include "DIMENSIONS"
1142       include "COMMON.INTERACT"
1143       include "COMMON.SETUP"
1144       include "COMMON.IOUNITS"
1145       integer ii,jj,itask(4),ntask_cont_to,
1146      & itask_cont_to(0:max_fg_procs-1)
1147       logical flag
1148       integer iturn3_start_all,iturn3_end_all,iturn4_start_all,
1149      & iturn4_end_all,iatel_s_all,iatel_e_all,ielstart_all,ielend_all
1150       common /przechowalnia/ iturn3_start_all(0:max_fg_procs),
1151      & iturn3_end_all(0:max_fg_procs),iturn4_start_all(0:max_fg_procs),
1152      & iturn4_end_all(0:max_fg_procs),iatel_s_all(0:max_fg_procs),
1153      &iatel_e_all(0:max_fg_procs),ielstart_all(maxres,0:max_fg_procs-1),
1154      & ielend_all(maxres,0:max_fg_procs-1)
1155       integer iproc,isent,k,l
1156 c Determines whether to send interaction ii,jj to other processors; a given
1157 c interaction can be sent to at most 2 processors.
1158 c Sets flag=.true. if interaction ii,jj needs to be sent to at least
1159 c one processor, otherwise flag is unchanged from the input value.
1160       isent=0
1161       itask(1)=fg_rank
1162       itask(2)=fg_rank
1163       itask(3)=fg_rank
1164       itask(4)=fg_rank
1165 c      write (iout,*) "ii",ii," jj",jj
1166 c Loop over processors to check if anybody could need interaction ii,jj
1167       do iproc=0,fg_rank-1
1168 c Check if the interaction matches any turn3 at iproc
1169         do k=iturn3_start_all(iproc),iturn3_end_all(iproc)
1170           l=k+2
1171           if (k.eq.ii-1 .and. l.eq.jj-1 .or. k.eq.ii-1 .and. l.eq.jj+1
1172      &   .or. k.eq.ii+1 .and. l.eq.jj+1 .or. k.eq.ii+1 .and. l.eq.jj-1)
1173      &    then
1174 c            write (iout,*) "turn3 to iproc",iproc," ij",ii,jj,"kl",k,l
1175 c            call flush(iout)
1176             flag=.true.
1177             if (iproc.ne.itask(1).and.iproc.ne.itask(2)
1178      &        .and.iproc.ne.itask(3).and.iproc.ne.itask(4)) then
1179               isent=isent+1
1180               itask(isent)=iproc
1181               call add_task(iproc,ntask_cont_to,itask_cont_to)
1182             endif
1183           endif
1184         enddo
1185 C Check if the interaction matches any turn4 at iproc
1186         do k=iturn4_start_all(iproc),iturn4_end_all(iproc)
1187           l=k+3
1188           if (k.eq.ii-1 .and. l.eq.jj-1 .or. k.eq.ii-1 .and. l.eq.jj+1
1189      &   .or. k.eq.ii+1 .and. l.eq.jj+1 .or. k.eq.ii+1 .and. l.eq.jj-1)
1190      &    then
1191 c            write (iout,*) "turn3 to iproc",iproc," ij",ii,jj," kl",k,l
1192 c            call flush(iout)
1193             flag=.true.
1194             if (iproc.ne.itask(1).and.iproc.ne.itask(2)
1195      &        .and.iproc.ne.itask(3).and.iproc.ne.itask(4)) then
1196               isent=isent+1
1197               itask(isent)=iproc
1198               call add_task(iproc,ntask_cont_to,itask_cont_to)
1199             endif
1200           endif
1201         enddo
1202         if (iatel_s_all(iproc).gt.0 .and. iatel_e_all(iproc).gt.0 .and.
1203      &  iatel_s_all(iproc).le.ii-1 .and. iatel_e_all(iproc).ge.ii-1)then
1204           if (ielstart_all(ii-1,iproc).le.jj-1.and.
1205      &        ielend_all(ii-1,iproc).ge.jj-1) then
1206             flag=.true.
1207             if (iproc.ne.itask(1).and.iproc.ne.itask(2)
1208      &        .and.iproc.ne.itask(3).and.iproc.ne.itask(4)) then
1209               isent=isent+1
1210               itask(isent)=iproc
1211               call add_task(iproc,ntask_cont_to,itask_cont_to)
1212             endif
1213           endif
1214           if (ielstart_all(ii-1,iproc).le.jj+1.and.
1215      &        ielend_all(ii-1,iproc).ge.jj+1) then
1216             flag=.true.
1217             if (iproc.ne.itask(1).and.iproc.ne.itask(2)
1218      &        .and.iproc.ne.itask(3).and.iproc.ne.itask(4)) then
1219               isent=isent+1
1220               itask(isent)=iproc
1221               call add_task(iproc,ntask_cont_to,itask_cont_to)
1222             endif
1223           endif
1224         endif
1225       enddo
1226       return
1227       end
1228 c---------------------------------------------------------------------------
1229       subroutine add_int_from(ii,jj,ntask_cont_from,itask_cont_from)
1230       implicit none
1231       include "DIMENSIONS"
1232       include "COMMON.INTERACT"
1233       include "COMMON.SETUP"
1234       include "COMMON.IOUNITS"
1235       integer ii,jj,itask(2),ntask_cont_from,
1236      & itask_cont_from(0:max_fg_procs-1)
1237       logical flag
1238       integer iturn3_start_all,iturn3_end_all,iturn4_start_all,
1239      & iturn4_end_all,iatel_s_all,iatel_e_all,ielstart_all,ielend_all
1240       common /przechowalnia/ iturn3_start_all(0:max_fg_procs),
1241      & iturn3_end_all(0:max_fg_procs),iturn4_start_all(0:max_fg_procs),
1242      & iturn4_end_all(0:max_fg_procs),iatel_s_all(0:max_fg_procs),
1243      &iatel_e_all(0:max_fg_procs),ielstart_all(maxres,0:max_fg_procs-1),
1244      & ielend_all(maxres,0:max_fg_procs-1)
1245       integer iproc,k,l
1246       do iproc=fg_rank+1,nfgtasks-1
1247         do k=iturn3_start_all(iproc),iturn3_end_all(iproc)
1248           l=k+2
1249           if (k.eq.ii+1 .and. l.eq.jj+1 .or. k.eq.ii+1.and.l.eq.jj-1
1250      &   .or. k.eq.ii-1 .and. l.eq.jj-1 .or. k.eq.ii-1 .and. l.eq.jj+1)
1251      &    then
1252 c            write (iout,*)"turn3 from iproc",iproc," ij",ii,jj," kl",k,l
1253             call add_task(iproc,ntask_cont_from,itask_cont_from)
1254           endif
1255         enddo
1256         do k=iturn4_start_all(iproc),iturn4_end_all(iproc)
1257           l=k+3
1258           if (k.eq.ii+1 .and. l.eq.jj+1 .or. k.eq.ii+1.and.l.eq.jj-1
1259      &   .or. k.eq.ii-1 .and. l.eq.jj-1 .or. k.eq.ii-1 .and. l.eq.jj+1)
1260      &    then
1261 c            write (iout,*)"turn4 from iproc",iproc," ij",ii,jj," kl",k,l
1262             call add_task(iproc,ntask_cont_from,itask_cont_from)
1263           endif
1264         enddo
1265         if (iatel_s_all(iproc).gt.0 .and. iatel_e_all(iproc).gt.0) then
1266           if (ii+1.ge.iatel_s_all(iproc).and.ii+1.le.iatel_e_all(iproc))
1267      &    then
1268             if (jj+1.ge.ielstart_all(ii+1,iproc).and.
1269      &          jj+1.le.ielend_all(ii+1,iproc)) then
1270               call add_task(iproc,ntask_cont_from,itask_cont_from)
1271             endif
1272             if (jj-1.ge.ielstart_all(ii+1,iproc).and.
1273      &          jj-1.le.ielend_all(ii+1,iproc)) then
1274               call add_task(iproc,ntask_cont_from,itask_cont_from)
1275             endif
1276           endif
1277           if (ii-1.ge.iatel_s_all(iproc).and.ii-1.le.iatel_e_all(iproc))
1278      &    then
1279             if (jj-1.ge.ielstart_all(ii-1,iproc).and.
1280      &          jj-1.le.ielend_all(ii-1,iproc)) then
1281               call add_task(iproc,ntask_cont_from,itask_cont_from)
1282             endif
1283             if (jj+1.ge.ielstart_all(ii-1,iproc).and.
1284      &          jj+1.le.ielend_all(ii-1,iproc)) then
1285                call add_task(iproc,ntask_cont_from,itask_cont_from)
1286             endif
1287           endif
1288         endif
1289       enddo
1290       return
1291       end
1292 c---------------------------------------------------------------------------
1293       subroutine add_task(iproc,ntask_cont,itask_cont)
1294       implicit none
1295       include "DIMENSIONS"
1296       integer iproc,ntask_cont,itask_cont(0:max_fg_procs-1)
1297       integer ii
1298       do ii=1,ntask_cont
1299         if (itask_cont(ii).eq.iproc) return
1300       enddo
1301       ntask_cont=ntask_cont+1
1302       itask_cont(ntask_cont)=iproc
1303       return
1304       end
1305 c---------------------------------------------------------------------------
1306       subroutine int_bounds(total_ints,lower_bound,upper_bound)
1307       implicit real*8 (a-h,o-z)
1308       include 'DIMENSIONS'
1309       include 'mpif.h'
1310       include 'COMMON.SETUP'
1311       integer total_ints,lower_bound,upper_bound
1312       integer int4proc(0:max_fg_procs),sint4proc(0:max_fg_procs)
1313       nint=total_ints/nfgtasks
1314       do i=1,nfgtasks
1315         int4proc(i-1)=nint
1316       enddo
1317       nexcess=total_ints-nint*nfgtasks
1318       do i=1,nexcess
1319         int4proc(nfgtasks-i)=int4proc(nfgtasks-i)+1
1320       enddo
1321       lower_bound=0
1322       do i=0,fg_rank-1
1323         lower_bound=lower_bound+int4proc(i)
1324       enddo
1325       upper_bound=lower_bound+int4proc(fg_rank)
1326       lower_bound=lower_bound+1
1327       return
1328       end
1329 c---------------------------------------------------------------------------
1330       subroutine int_bounds1(total_ints,lower_bound,upper_bound)
1331       implicit real*8 (a-h,o-z)
1332       include 'DIMENSIONS'
1333       include 'mpif.h'
1334       include 'COMMON.SETUP'
1335       integer total_ints,lower_bound,upper_bound
1336       integer int4proc(0:max_fg_procs),sint4proc(0:max_fg_procs)
1337       nint=total_ints/nfgtasks1
1338       do i=1,nfgtasks1
1339         int4proc(i-1)=nint
1340       enddo
1341       nexcess=total_ints-nint*nfgtasks1
1342       do i=1,nexcess
1343         int4proc(nfgtasks1-i)=int4proc(nfgtasks1-i)+1
1344       enddo
1345       lower_bound=0
1346       do i=0,fg_rank1-1
1347         lower_bound=lower_bound+int4proc(i)
1348       enddo
1349       upper_bound=lower_bound+int4proc(fg_rank1)
1350       lower_bound=lower_bound+1
1351       return
1352       end
1353 c---------------------------------------------------------------------------
1354       subroutine int_partition(int_index,lower_index,upper_index,atom,
1355      & at_start,at_end,first_atom,last_atom,int_gr,jat_start,jat_end,*)
1356       implicit real*8 (a-h,o-z)
1357       include 'DIMENSIONS'
1358       include 'COMMON.IOUNITS'
1359       integer int_index,lower_index,upper_index,atom,at_start,at_end,
1360      & first_atom,last_atom,int_gr,jat_start,jat_end
1361       logical lprn
1362       lprn=.false.
1363       if (lprn) write (iout,*) 'int_index=',int_index
1364       int_index_old=int_index
1365       int_index=int_index+last_atom-first_atom+1
1366       if (lprn)
1367      &   write (iout,*) 'int_index=',int_index,
1368      &               ' int_index_old',int_index_old,
1369      &               ' lower_index=',lower_index,
1370      &               ' upper_index=',upper_index,
1371      &               ' atom=',atom,' first_atom=',first_atom,
1372      &               ' last_atom=',last_atom
1373       if (int_index.ge.lower_index) then
1374         int_gr=int_gr+1
1375         if (at_start.eq.0) then
1376           at_start=atom
1377           jat_start=first_atom-1+lower_index-int_index_old
1378         else
1379           jat_start=first_atom
1380         endif
1381         if (lprn) write (iout,*) 'jat_start',jat_start
1382         if (int_index.ge.upper_index) then
1383           at_end=atom
1384           jat_end=first_atom-1+upper_index-int_index_old
1385           return1
1386         else
1387           jat_end=last_atom
1388         endif
1389         if (lprn) write (iout,*) 'jat_end',jat_end
1390       endif
1391       return
1392       end
1393 #endif
1394 c------------------------------------------------------------------------------
1395       subroutine hpb_partition
1396       implicit real*8 (a-h,o-z)
1397       include 'DIMENSIONS'
1398 #ifdef MPI
1399       include 'mpif.h'
1400 #endif
1401       include 'COMMON.SBRIDGE'
1402       include 'COMMON.IOUNITS'
1403       include 'COMMON.SETUP'
1404       include 'COMMON.CONTROL'
1405 #ifdef MPI
1406       call int_bounds(nhpb,link_start,link_end)
1407       if (.not. out1file)
1408      &  write (iout,*) 'Processor',fg_rank,' CG group',kolor,
1409      &  ' absolute rank',MyRank,
1410      &  ' nhpb',nhpb,' link_start=',link_start,
1411      &  ' link_end',link_end
1412 #else
1413       link_start=1
1414       link_end=nhpb
1415 #endif
1416       return
1417       end