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