source/unres/src_MD/initialize_p.F

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