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