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