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