source/unres/src_MD/initialize_p.F

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