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