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