source/unres/src_MD-M/initialize_p.F

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