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