source/wham/src-M/readrtns.F

   1       subroutine read_general_data(*)
   2       implicit none
   3       include "DIMENSIONS"
   4       include "DIMENSIONS.ZSCOPT"
   5       include "DIMENSIONS.FREE"
   6       include "COMMON.TORSION"
   7       include "COMMON.INTERACT"
   8       include "COMMON.IOUNITS"
   9       include "COMMON.TIME1"
  10       include "COMMON.PROT"
  11       include "COMMON.PROTFILES"
  12       include "COMMON.CHAIN"
  13       include "COMMON.NAMES"
  14       include "COMMON.FFIELD"
  15       include "COMMON.ENEPS"
  16       include "COMMON.WEIGHTS"
  17       include "COMMON.FREE"
  18       include "COMMON.CONTROL"
  19       include "COMMON.ENERGIES"
  20       include "COMMON.SPLITELE"
  21       include "COMMON.SBRIDGE"
  22       character*800 controlcard
  23       integer i,j,k,ii,n_ene_found
  24       integer ind,itype1,itype2,itypf,itypsc,itypp
  25       integer ilen
  26       external ilen
  27       character*16 ucase
  28       character*16 key
  29       external ucase
  30       double precision pi
  31       call card_concat(controlcard,.true.)
  32       call readi(controlcard,"N_ENE",n_ene,max_ene)
  33       if (n_ene.gt.max_ene) then
  34         write (iout,*) "Error: parameter out of range: N_ENE",n_ene,
  35      &    max_ene
  36         return1
  37       endif
  38       call readi(controlcard,"NPARMSET",nparmset,1)
  39       separate_parset = index(controlcard,"SEPARATE_PARSET").gt.0
  40       call readi(controlcard,"IPARMPRINT",iparmprint,1)
  41       write (iout,*) "PARMPRINT",iparmprint
  42       if (nparmset.gt.max_parm) then
  43         write (iout,*) "Error: parameter out of range: NPARMSET",
  44      &    nparmset, Max_Parm
  45         return1
  46       endif
  47       call readi(controlcard,"MAXIT",maxit,5000)
  48       call reada(controlcard,"FIMIN",fimin,1.0d-3)
  49       call readi(controlcard,"ENSEMBLES",ensembles,0)
  50       hamil_rep=index(controlcard,"HAMIL_REP").gt.0
  51       write (iout,*) "Number of energy parameter sets",nparmset
  52       call multreadi(controlcard,"ISAMPL",isampl,nparmset,1)
  53       write (iout,*) "MaxSlice",MaxSlice
  54       call readi(controlcard,"NSLICE",nslice,1)
  55       call flush(iout)
  56       if (nslice.gt.MaxSlice) then
  57         write (iout,*) "Error: parameter out of range: NSLICE",nslice,
  58      &    MaxSlice
  59         return1
  60       endif
  61       write (iout,*) "Frequency of storing conformations",
  62      & (isampl(i),i=1,nparmset)
  63       write (iout,*) "Maxit",maxit," Fimin",fimin
  64       call readi(controlcard,"NQ",nQ,1)
  65       if (nQ.gt.MaxQ) then
  66         write (iout,*) "Error: parameter out of range: NQ",nq,
  67      &    maxq
  68         return1
  69       endif
  70       indpdb=0
  71       if (index(controlcard,"CLASSIFY").gt.0) indpdb=1
  72       call reada(controlcard,"DELTA",delta,1.0d-2)
  73       call readi(controlcard,"EINICHECK",einicheck,2)
  74       call reada(controlcard,"DELTRMS",deltrms,5.0d-2)
  75       call reada(controlcard,"DELTRGY",deltrgy,5.0d-2)
  76       call readi(controlcard,"RESCALE",rescale_mode,1)
  77       check_conf=index(controlcard,"NO_CHECK_CONF").eq.0
  78       call reada(controlcard,'DISTCHAINMAX',distchainmax,50.0d0)
  79        call reada(controlcard,'BOXX',boxxsize,100.0d0)
  80        call reada(controlcard,'BOXY',boxysize,100.0d0)
  81        call reada(controlcard,'BOXZ',boxzsize,100.0d0)
  82 c Cutoff range for interactions
  83       call reada(controlcard,"R_CUT",r_cut,15.0d0)
  84       call reada(controlcard,"LAMBDA",rlamb,0.3d0)
  85       call reada(controlcard,"LIPTHICK",lipthick,0.0d0)
  86       call reada(controlcard,"LIPAQBUF",lipbufthick,0.0d0)
  87       if (lipthick.gt.0.0d0) then
  88        bordliptop=(boxzsize+lipthick)/2.0
  89        bordlipbot=bordliptop-lipthick
  90 C      endif
  91       if ((bordliptop.gt.boxzsize).or.(bordlipbot.lt.0.0))
  92      & write(iout,*) "WARNING WRONG SIZE OF LIPIDIC PHASE"
  93       buflipbot=bordlipbot+lipbufthick
  94       bufliptop=bordliptop-lipbufthick
  95       if ((lipbufthick*2.0d0).gt.lipthick)
  96      &write(iout,*) "WARNING WRONG SIZE OF LIP AQ BUF"
  97       endif
  98       write(iout,*) "bordliptop=",bordliptop
  99       write(iout,*) "bordlipbot=",bordlipbot
 100       write(iout,*) "bufliptop=",bufliptop
 101       write(iout,*) "buflipbot=",buflipbot
 102       call readi(controlcard,'SYM',symetr,1)
 103       write (iout,*) "DISTCHAINMAX",distchainmax
 104       write (iout,*) "delta",delta
 105       write (iout,*) "einicheck",einicheck
 106       write (iout,*) "rescale_mode",rescale_mode
 107       call flush(iout)
 108       bxfile=index(controlcard,"BXFILE").gt.0
 109       cxfile=index(controlcard,"CXFILE").gt.0
 110       if (nslice .eq. 1 .and. .not.bxfile .and. .not.cxfile)
 111      & bxfile=.true.
 112       histfile=index(controlcard,"HISTFILE").gt.0
 113       histout=index(controlcard,"HISTOUT").gt.0
 114       entfile=index(controlcard,"ENTFILE").gt.0
 115       zscfile=index(controlcard,"ZSCFILE").gt.0
 116       with_dihed_constr = index(controlcard,"WITH_DIHED_CONSTR").gt.0
 117       write (iout,*) "with_dihed_constr ",with_dihed_constr
 118       with_theta_constr = index(controlcard,"WITH_THETA_CONSTR").gt.0
 119       write (iout,*) "with_theta_constr ",with_theta_constr
 120       call readi(controlcard,'CONSTR_DIST',constr_dist,0)
 121       write (iout,*) "with_dihed_constr ",with_dihed_constr,
 122      & " CONSTR_DIST",constr_dist
 123       call readi(controlcard,'CONSTR_HOMOL',constr_homology,0)
 124       write (iout,*) "with_homology_constr ",with_dihed_constr,
 125      & " CONSTR_HOMOLOGY",constr_homology
 126       refstr = index(controlcard,'REFSTR').gt.0
 127       pdbref = index(controlcard,'PDBREF').gt.0
 128       dyn_ss=(index(controlcard,'DYN_SS').gt.0)
 129 C /06/28/2013 Adasko: dyn_ss is keyword allowing to break and create bond
 130 C disulfide bond. Note that in conterary to dynamics this in
 131 C CONTROLCARD. The bond are read in molread_zs.F
 132       call flush(iout)
 133       return
 134       end
 135 c------------------------------------------------------------------------------
 136       subroutine read_efree(*)
 137 C
 138 C Read molecular data
 139 C
 140       implicit none
 141       include 'DIMENSIONS'
 142       include 'DIMENSIONS.ZSCOPT'
 143       include 'DIMENSIONS.COMPAR'
 144       include 'DIMENSIONS.FREE'
 145       include 'COMMON.IOUNITS'
 146       include 'COMMON.TIME1'
 147       include 'COMMON.SBRIDGE'
 148       include 'COMMON.CONTROL'
 149       include 'COMMON.CHAIN'
 150       include 'COMMON.HEADER'
 151       include 'COMMON.GEO'
 152       include 'COMMON.FREE'
 153       character*320 controlcard,ucase
 154       integer iparm,ib,i,j,npars
 155       integer ilen
 156       external ilen
 157
 158       if (hamil_rep) then
 159         npars=1
 160       else
 161         npars=nParmSet
 162       endif
 163
 164       do iparm=1,npars
 165
 166       call card_concat(controlcard,.true.)
 167       call readi(controlcard,'NT',nT_h(iparm),1)
 168       write (iout,*) "iparm",iparm," nt",nT_h(iparm)
 169       call flush(iout)
 170       if (nT_h(iparm).gt.MaxT_h) then
 171         write (iout,*)  "Error: parameter out of range: NT",nT_h(iparm),
 172      &    MaxT_h
 173         return1
 174       endif
 175       replica(iparm)=index(controlcard,"REPLICA").gt.0
 176       umbrella(iparm)=index(controlcard,"UMBRELLA").gt.0
 177       read_iset(iparm)=index(controlcard,"READ_ISET").gt.0
 178       write (iout,*) "nQ",nQ," nT",nT_h(iparm)," replica ",
 179      &  replica(iparm)," umbrella ",umbrella(iparm),
 180      &  " read_iset",read_iset(iparm)
 181       call flush(iout)
 182       do ib=1,nT_h(iparm)
 183         call card_concat(controlcard,.true.)
 184         call readi(controlcard,'NR',nR(ib,iparm),1)
 185         if (umbrella(iparm)) then
 186           nRR(ib,iparm)=1
 187         else
 188           nRR(ib,iparm)=nR(ib,iparm)
 189         endif
 190         if (nR(ib,iparm).gt.MaxR) then
 191           write (iout,*)  "Error: parameter out of range: NR",
 192      &      nR(ib,iparm),MaxR
 193         return1
 194         endif
 195         call reada(controlcard,'TEMP',beta_h(ib,iparm),298.0d0)
 196         beta_h(ib,iparm)=1.0d0/(beta_h(ib,iparm)*1.987D-3)
 197         call multreada(controlcard,'FI',f(1,ib,iparm),nR(ib,iparm),
 198      &    0.0d0)
 199         do i=1,nR(ib,iparm)
 200           call card_concat(controlcard,.true.)
 201           call multreada(controlcard,'KH',KH(1,i,ib,iparm),nQ,
 202      &      100.0d0)
 203           call multreada(controlcard,'Q0',Q0(1,i,ib,iparm),nQ,
 204      &      0.0d0)
 205         enddo
 206       enddo
 207       do ib=1,nT_h(iparm)
 208         write (iout,*) "ib",ib," beta_h",
 209      &    1.0d0/(0.001987*beta_h(ib,iparm))
 210         write (iout,*) "nR",nR(ib,iparm)
 211         write (iout,*) "fi",(f(i,ib,iparm),i=1,nR(ib,iparm))
 212         do i=1,nR(ib,iparm)
 213           write (iout,*) "i",i," Kh",(Kh(j,i,ib,iparm),j=1,nQ),
 214      &      "q0",(q0(j,i,ib,iparm),j=1,nQ)
 215         enddo
 216         call flush(iout)
 217       enddo
 218
 219       write (iout,*) "HOMOL_NSET",homol_nset
 220       enddo
 221
 222       if (hamil_rep) then
 223
 224        do iparm=2,nParmSet
 225           nT_h(iparm)=nT_h(1)
 226          do ib=1,nT_h(iparm)
 227            nR(ib,iparm)=nR(ib,1)
 228            if (umbrella(iparm)) then
 229              nRR(ib,iparm)=1
 230            else
 231              nRR(ib,iparm)=nR(ib,1)
 232            endif
 233            beta_h(ib,iparm)=beta_h(ib,1)
 234            do i=1,nR(ib,iparm)
 235              f(i,ib,iparm)=f(i,ib,1)
 236              do j=1,nQ
 237                KH(j,i,ib,iparm)=KH(j,i,ib,1)
 238                Q0(j,i,ib,iparm)=Q0(j,i,ib,1)
 239              enddo
 240            enddo
 241            replica(iparm)=replica(1)
 242            umbrella(iparm)=umbrella(1)
 243            read_iset(iparm)=read_iset(1)
 244          enddo
 245        enddo
 246
 247       endif
 248
 249       return
 250       end
 251 c-----------------------------------------------------------------------------
 252       subroutine read_protein_data(*)
 253       implicit none
 254       include "DIMENSIONS"
 255       include "DIMENSIONS.ZSCOPT"
 256       include "DIMENSIONS.FREE"
 257 #ifdef MPI
 258       include "mpif.h"
 259       integer IERROR,ERRCODE,STATUS(MPI_STATUS_SIZE)
 260       include "COMMON.MPI"
 261 #endif
 262       include "COMMON.CHAIN"
 263       include "COMMON.IOUNITS"
 264       include "COMMON.PROT"
 265       include "COMMON.PROTFILES"
 266       include "COMMON.NAMES"
 267       include "COMMON.FREE"
 268       include "COMMON.OBCINKA"
 269       character*64 nazwa
 270       character*16000 controlcard
 271       integer i,ii,ib,iR,iparm,ilen,iroof,nthr,npars
 272       external ilen,iroof
 273       if (hamil_rep) then
 274         npars=1
 275       else
 276         npars=nparmset
 277       endif
 278
 279       do iparm=1,npars
 280
 281 C Read names of files with conformation data.
 282       if (replica(iparm)) then
 283         nthr = 1
 284       else
 285         nthr = nT_h(iparm)
 286       endif
 287       do ib=1,nthr
 288       do ii=1,nRR(ib,iparm)
 289       write (iout,*) "Parameter set",iparm," temperature",ib,
 290      & " window",ii
 291       call flush(iout)
 292       call card_concat(controlcard,.true.)
 293       write (iout,*) controlcard(:ilen(controlcard))
 294       call readi(controlcard,"NFILE_BIN",nfile_bin(ii,ib,iparm),0)
 295       call readi(controlcard,"NFILE_ASC",nfile_asc(ii,ib,iparm),0)
 296       call readi(controlcard,"NFILE_CX",nfile_cx(ii,ib,iparm),0)
 297       call readi(controlcard,"REC_START",rec_start(ii,ib,iparm),1)
 298       call readi(controlcard,"REC_END",rec_end(ii,ib,iparm),
 299      & maxstr*isampl(iparm)+rec_start(ii,ib,iparm)-1)
 300       call reada(controlcard,"TIME_START",
 301      &  time_start_collect(ii,ib,iparm),0.0d0)
 302       call reada(controlcard,"TIME_END",time_end_collect(ii,ib,iparm),
 303      &  1.0d10)
 304       write (iout,*) "rec_start",rec_start(ii,ib,iparm),
 305      & " rec_end",rec_end(ii,ib,iparm)
 306       write (iout,*) "time_start",time_start_collect(ii,ib,iparm),
 307      & " time_end",time_end_collect(ii,ib,iparm)
 308       call flush(iout)
 309       if (replica(iparm)) then
 310         call readi(controlcard,"TOTRAJ",totraj(ii,iparm),1)
 311         write (iout,*) "Number of trajectories",totraj(ii,iparm)
 312         call flush(iout)
 313       endif
 314       if (nfile_bin(ii,ib,iparm).lt.2
 315      &    .and. nfile_asc(ii,ib,iparm).eq.0
 316      &    .and. nfile_cx(ii,ib,iparm).eq.0) then
 317         write (iout,*) "Error - no action specified!"
 318         return1
 319       endif
 320       if (nfile_bin(ii,ib,iparm).gt.0) then
 321         call card_concat(controlcard,.false.)
 322         call split_string(controlcard,protfiles(1,1,ii,ib,iparm),
 323      &   maxfile_prot,nfile_bin(ii,ib,iparm))
 324 #ifdef DEBUG
 325         write(iout,*)"nfile_bin",nfile_bin(ii,ib,iparm)
 326         write(iout,*) (protfiles(i,1,ii,ib,iparm),
 327      &    i=1,nfile_bin(ii,ib,iparm))
 328 #endif
 329       endif
 330       if (nfile_asc(ii,ib,iparm).gt.0) then
 331         call card_concat(controlcard,.false.)
 332         call split_string(controlcard,protfiles(1,2,ii,ib,iparm),
 333      &   maxfile_prot,nfile_asc(ii,ib,iparm))
 334 #ifdef DEBUG
 335         write(iout,*) "nfile_asc(ii,ib,iparm)",nfile_asc(ii,ib,iparm)
 336         write(iout,*) (protfiles(i,2,ii,ib,iparm),
 337      &    i=1,nfile_asc(ii,ib,iparm))
 338 #endif
 339       else if (nfile_cx(ii,ib,iparm).gt.0) then
 340         call card_concat(controlcard,.false.)
 341         call split_string(controlcard,protfiles(1,2,ii,ib,iparm),
 342      &   maxfile_prot,nfile_cx(ii,ib,iparm))
 343 #ifdef DEBUG
 344         write(iout,*) "nfile_cx(ii,ib,iparm)",nfile_cx(ii,ib,iparm)
 345         write(iout,*) (protfiles(i,2,ii,ib,iparm),
 346      &    i=1,nfile_cx(ii,ib,iparm))
 347 #endif
 348       endif
 349       call flush(iout)
 350       enddo
 351       enddo
 352
 353       enddo
 354
 355       return
 356       end
 357 c-------------------------------------------------------------------------------
 358       subroutine opentmp(islice,iunit,bprotfile_temp)
 359       implicit none
 360       include "DIMENSIONS"
 361       include "DIMENSIONS.ZSCOPT"
 362       include "DIMENSIONS.FREE"
 363 #ifdef MPI
 364       include "mpif.h"
 365       integer IERROR,ERRCODE,STATUS(MPI_STATUS_SIZE)
 366       include "COMMON.MPI"
 367 #endif
 368       include "COMMON.IOUNITS"
 369       include "COMMON.PROTFILES"
 370       include "COMMON.PROT"
 371       include "COMMON.FREE"
 372       character*64 bprotfile_temp
 373       character*3 liczba,liczba2
 374       character*2 liczba1
 375       integer iunit,islice
 376       integer ilen,iroof
 377       external ilen,iroof
 378       logical lerr
 379
 380       write (liczba1,'(bz,i2.2)') islice
 381       write (liczba,'(bz,i3.3)') me
 382 #ifdef MPI
 383 c      write (iout,*) "separate_parset ",separate_parset,
 384 c     &  " myparm",myparm
 385       if (separate_parset) then
 386       write (liczba2,'(bz,i3.3)') myparm
 387       bprotfile_temp = scratchdir(:ilen(scratchdir))//"/"//
 388      &  prefix(:ilen(prefix))//liczba//"_"//liczba2//".xbin.tmp"//liczba1
 389       open (iunit,file=bprotfile_temp,status="unknown",
 390      &    form="unformatted",access="direct",recl=lenrec)
 391       else
 392       bprotfile_temp = scratchdir(:ilen(scratchdir))//"/"//
 393      &  prefix(:ilen(prefix))//liczba//".xbin.tmp"//liczba1
 394       open (iunit,file=bprotfile_temp,status="unknown",
 395      &    form="unformatted",access="direct",recl=lenrec)
 396       endif
 397 #else
 398       bprotfile_temp = scratchdir(:ilen(scratchdir))//
 399      &  "/"//prefix(:ilen(prefix))//".xbin.tmp"//liczba1
 400       open (iunit,file=bprotfile_temp,status="unknown",
 401      &    form="unformatted",access="direct",recl=lenrec)
 402 #endif
 403 c      write (iout,*) "OpenTmp iunit",iunit," bprotfile_temp",
 404 c     &  bprotfile_temp
 405 c      call flush(iout)
 406       return
 407       end
 408 c-------------------------------------------------------------------------------
 409       subroutine read_database(*)
 410       implicit none
 411       include "DIMENSIONS"
 412       include "DIMENSIONS.ZSCOPT"
 413       include "DIMENSIONS.FREE"
 414 #ifdef MPI
 415       include "mpif.h"
 416       integer IERROR,ERRCODE,STATUS(MPI_STATUS_SIZE)
 417       include "COMMON.MPI"
 418 #endif
 419       include "COMMON.CHAIN"
 420       include "COMMON.IOUNITS"
 421       include "COMMON.PROTFILES"
 422       include "COMMON.NAMES"
 423       include "COMMON.VAR"
 424       include "COMMON.GEO"
 425       include "COMMON.ENEPS"
 426       include "COMMON.PROT"
 427       include "COMMON.INTERACT"
 428       include "COMMON.FREE"
 429       include "COMMON.SBRIDGE"
 430       include "COMMON.OBCINKA"
 431       real*4 csingle(3,maxres2)
 432       character*64 nazwa,bprotfile_temp
 433       character*3 liczba
 434       character*2 liczba1
 435       integer i,j,ii,jj(maxslice),k,kk(maxslice),l,
 436      &  ll(maxslice),mm(maxslice),if
 437       integer nrec,nlines,iscor,iunit,islice
 438       double precision energ
 439       integer ilen,iroof
 440       external ilen,iroof
 441       double precision rmsdev,energia(0:max_ene),efree,eini,temp
 442       double precision prop(maxQ)
 443       integer ntot_all(maxslice,0:maxprocs-1), maxslice_buff
 444       integer iparm,ib,iib,ir,nprop,nthr,npars
 445       double precision etot,time
 446       integer ixdrf,iret
 447       logical lerr,linit
 448
 449       lenrec1=12*(nres+nct-nnt+1)+4*(2*nss+2)+24
 450       lenrec2=12*(nres+nct-nnt+1)+4*(2*nss+2)+24+8*nQ
 451       lenrec=lenrec2+8
 452       write (iout,*) "lenrec",lenrec," lenrec1",lenrec1,
 453      &  " lenrec2",lenrec2
 454
 455       do i=1,nQ
 456         prop(i)=0.0d0
 457       enddo
 458       do islice=1,nslice
 459         ll(islice)=0
 460         mm(islice)=0
 461       enddo
 462       write (iout,*) "nparmset",nparmset
 463       if (hamil_rep) then
 464         npars=1
 465       else
 466         npars=nparmset
 467       endif
 468       do iparm=1,npars
 469
 470       if (replica(iparm)) then
 471         nthr = 1
 472       else
 473         nthr = nT_h(iparm)
 474       endif
 475
 476       do ib=1,nthr
 477       do iR=1,nRR(ib,iparm)
 478
 479       write (iout,*) "iparm",iparm," ib",ib," iR",iR," nQ",nQ
 480       do islice=1,nslice
 481         jj(islice)=0
 482         kk(islice)=0
 483       enddo
 484
 485       IF (NFILE_BIN(iR,ib,iparm).GT.0) THEN
 486 c Read conformations from binary DA files (one per batch) and write them to
 487 c a binary DA scratchfile.
 488         write (liczba,'(bz,i3.3)') me
 489         do if=1,nfile_bin(iR,ib,iparm)
 490           nazwa=protfiles(if,1,iR,ib,iparm)
 491      &     (:ilen(protfiles(if,1,iR,ib,iparm)))//".bx"
 492           open (ientin,file=nazwa,status="old",form="unformatted",
 493      &     access="direct",recl=lenrec2,err=1111)
 494           ii=0
 495           do islice=1,nslice
 496             call opentmp(islice,ientout,bprotfile_temp)
 497             call bxread(nazwa,ii,jj(islice),kk(islice),ll(islice),
 498      &        mm(islice),iR,ib,iparm)
 499             close(ientout)
 500           enddo
 501           close(ientin)
 502         enddo
 503       ENDIF ! NFILE_BIN>0
 504 c
 505       IF (NFILE_ASC(iR,ib,iparm).GT.0) THEN
 506 c Read conformations from multiple ASCII int files and write them to a binary
 507 c DA scratchfile.
 508         do if=1,nfile_asc(iR,ib,iparm)
 509           nazwa=protfiles(if,2,iR,ib,iparm)
 510      &     (:ilen(protfiles(if,2,iR,ib,iparm)))//".x"
 511           open(unit=ientin,file=nazwa,status='old',err=1111)
 512           write(iout,*) "reading ",nazwa(:ilen(nazwa))
 513           ii=0
 514           call xread(nazwa,ii,jj,kk,ll,mm,iR,ib,iparm)
 515         enddo ! if
 516       ENDIF
 517       IF (NFILE_CX(iR,ib,iparm).gt.0) THEN
 518 c Read conformations from cx files and write them to a binary
 519 c DA scratchfile.
 520         do if=1,nfile_cx(iR,ib,iparm)
 521           nazwa=protfiles(if,2,iR,ib,iparm)
 522      &     (:ilen(protfiles(if,2,iR,ib,iparm)))//".cx"
 523           write(iout,*) "reading ",nazwa(:ilen(nazwa))
 524           ii=0
 525           print *,"Calling cxread"
 526           call cxread(nazwa,ii,jj,kk,ll,mm,iR,ib,iparm,
 527      &       *1111)
 528           close(ientout)
 529           write (iout,*) "exit cxread"
 530           call flush(iout)
 531         enddo
 532       ENDIF
 533
 534       do islice=1,nslice
 535         stot(islice)=stot(islice)+jj(islice)
 536       enddo
 537
 538       enddo
 539       enddo
 540       write (iout,*) "IPARM",iparm
 541       enddo
 542
 543       if (nslice.eq.1) then
 544 #ifdef MPI
 545         write (liczba,'(bz,i3.3)') me
 546         bprotfile_temp = scratchdir(:ilen(scratchdir))//"/"//
 547      &    prefix(:ilen(prefix))//liczba//".xbin.tmp"
 548 #else
 549         bprotfile_temp = scratchdir(:ilen(scratchdir))//
 550      &     "/"//prefix(:ilen(prefix))//".xbin.tmp"
 551 #endif
 552         write(iout,*) mm(1)," conformations read",ll(1),
 553      &    " conformations written to ",
 554      &    bprotfile_temp(:ilen(bprotfile_temp))
 555       else
 556         do islice=1,nslice
 557           write (liczba1,'(bz,i2.2)') islice
 558 #ifdef MPI
 559           write (liczba,'(bz,i3.3)') me
 560           bprotfile_temp = scratchdir(:ilen(scratchdir))//"/"//
 561      &      prefix(:ilen(prefix))//liczba//".xbin.tmp"//liczba1
 562 #else
 563           bprotfile_temp = scratchdir(:ilen(scratchdir))//
 564      &       "/"//prefix(:ilen(prefix))//".xbin.tmp"//liczba1
 565 #endif
 566           write(iout,*) mm(islice)," conformations read",ll(islice),
 567      &    " conformations written to ",
 568      &    bprotfile_temp(:ilen(bprotfile_temp))
 569         enddo
 570       endif
 571
 572 #ifdef MPI
 573 c Check if everyone has the same number of conformations
 574
 575 c      call MPI_ALLgather(MPI_IN_PLACE,stot(1),MPI_DATATYPE_NULL,
 576 c     &  ntot_all(1,0),maxslice,MPI_INTEGER,MPI_Comm_World,IERROR)
 577
 578       maxslice_buff=maxslice
 579
 580       call MPI_Allgather(stot(1),maxslice_buff,MPI_INTEGER,
 581      &  ntot_all(1,0),maxslice,MPI_INTEGER,MPI_Comm_World,IERROR)
 582       lerr=.false.
 583       do i=0,nprocs-1
 584         if (i.ne.me) then
 585           do islice=1,nslice
 586           if (stot(islice).ne.ntot_all(islice,i)) then
 587             write (iout,*) "Number of conformations at processor",i,
 588      &       " differs from that at processor",me,
 589      &       stot(islice),ntot_all(islice,i)," slice",islice
 590             lerr = .true.
 591           endif
 592           enddo
 593         endif
 594       enddo
 595       if (lerr) then
 596         write (iout,*)
 597         write (iout,*) "Numbers of conformations read by processors"
 598         write (iout,*)
 599         do i=0,nprocs-1
 600           write (iout,'(8i10)') i,(ntot_all(islice,i),islice=1,nslice)
 601         enddo
 602         write (iout,*) "Calculation terminated."
 603         call flush(iout)
 604         return1
 605       endif
 606       do islice=1,nslice
 607         ntot(islice)=stot(islice)
 608       enddo
 609       return
 610 #endif
 611  1111 write(iout,*) "Error opening coordinate file ",nazwa(:ilen(nazwa))
 612       call flush(iout)
 613       return1
 614       end
 615 c------------------------------------------------------------------------------
 616       subroutine card_concat(card,to_upper)
 617       implicit none
 618       include 'DIMENSIONS.ZSCOPT'
 619       include "COMMON.IOUNITS"
 620       character*(*) card
 621       character*80 karta,ucase
 622       logical to_upper
 623       integer ilen
 624       external ilen
 625       read (inp,'(a)') karta
 626       if (to_upper) karta=ucase(karta)
 627       card=' '
 628       do while (karta(80:80).eq.'&')
 629         card=card(:ilen(card)+1)//karta(:79)
 630         read (inp,'(a)') karta
 631         if (to_upper) karta=ucase(karta)
 632       enddo
 633       card=card(:ilen(card)+1)//karta
 634       return
 635       end
 636 c------------------------------------------------------------------------------
 637       subroutine readi(rekord,lancuch,wartosc,default)
 638       implicit none
 639       character*(*) rekord,lancuch
 640       integer wartosc,default
 641       integer ilen,iread
 642       external ilen
 643       iread=index(rekord,lancuch(:ilen(lancuch))//"=")
 644       if (iread.eq.0) then
 645         wartosc=default
 646         return
 647       endif
 648       iread=iread+ilen(lancuch)+1
 649       read (rekord(iread:),*) wartosc
 650       return
 651       end
 652 c----------------------------------------------------------------------------
 653       subroutine reada(rekord,lancuch,wartosc,default)
 654       implicit none
 655       character*(*) rekord,lancuch
 656       character*80 aux
 657       double precision wartosc,default
 658       integer ilen,iread
 659       external ilen
 660       iread=index(rekord,lancuch(:ilen(lancuch))//"=")
 661       if (iread.eq.0) then
 662         wartosc=default
 663         return
 664       endif
 665       iread=iread+ilen(lancuch)+1
 666       read (rekord(iread:),*) wartosc
 667       return
 668       end
 669 c----------------------------------------------------------------------------
 670       subroutine multreadi(rekord,lancuch,tablica,dim,default)
 671       implicit none
 672       integer dim,i
 673       integer tablica(dim),default
 674       character*(*) rekord,lancuch
 675       character*80 aux
 676       integer ilen,iread
 677       external ilen
 678       do i=1,dim
 679         tablica(i)=default
 680       enddo
 681       iread=index(rekord,lancuch(:ilen(lancuch))//"=")
 682       if (iread.eq.0) return
 683       iread=iread+ilen(lancuch)+1
 684       read (rekord(iread:),*,end=10,err=10) (tablica(i),i=1,dim)
 685    10 return
 686       end
 687 c----------------------------------------------------------------------------
 688       subroutine multreada(rekord,lancuch,tablica,dim,default)
 689       implicit none
 690       integer dim,i
 691       double precision tablica(dim),default
 692       character*(*) rekord,lancuch
 693       character*80 aux
 694       integer ilen,iread
 695       external ilen
 696       do i=1,dim
 697         tablica(i)=default
 698       enddo
 699       iread=index(rekord,lancuch(:ilen(lancuch))//"=")
 700       if (iread.eq.0) return
 701       iread=iread+ilen(lancuch)+1
 702       read (rekord(iread:),*,end=10,err=10) (tablica(i),i=1,dim)
 703    10 return
 704       end
 705 c----------------------------------------------------------------------------
 706       subroutine reads(rekord,lancuch,wartosc,default)
 707       implicit none
 708       character*(*) rekord,lancuch,wartosc,default
 709       character*80 aux
 710       integer ilen,lenlan,lenrec,iread,ireade
 711       external ilen
 712       logical iblnk
 713       external iblnk
 714       lenlan=ilen(lancuch)
 715       lenrec=ilen(rekord)
 716       iread=index(rekord,lancuch(:lenlan)//"=")
 717 c      print *,"rekord",rekord," lancuch",lancuch
 718 c      print *,"iread",iread," lenlan",lenlan," lenrec",lenrec
 719       if (iread.eq.0) then
 720         wartosc=default
 721         return
 722       endif
 723       iread=iread+lenlan+1
 724 c      print *,"iread",iread
 725 c      print *,"|",rekord(iread:iread),"|",iblnk(rekord(iread:iread))
 726       do while (iread.le.lenrec .and. iblnk(rekord(iread:iread)))
 727         iread=iread+1
 728 c      print *,"|",rekord(iread:iread),"|",iblnk(rekord(iread:iread))
 729       enddo
 730 c      print *,"iread",iread
 731       if (iread.gt.lenrec) then
 732          wartosc=default
 733         return
 734       endif
 735       ireade=iread+1
 736 c      print *,"ireade",ireade
 737       do while (ireade.lt.lenrec .and.
 738      &   .not.iblnk(rekord(ireade:ireade)))
 739         ireade=ireade+1
 740       enddo
 741       wartosc=rekord(iread:ireade)
 742       return
 743       end
 744 c----------------------------------------------------------------------------
 745       subroutine multreads(rekord,lancuch,tablica,dim,default)
 746       implicit none
 747       integer dim,i
 748       character*(*) rekord,lancuch,tablica(dim),default
 749       character*80 aux
 750       integer ilen,lenlan,lenrec,iread,ireade
 751       external ilen
 752       logical iblnk
 753       external iblnk
 754       do i=1,dim
 755         tablica(i)=default
 756       enddo
 757       lenlan=ilen(lancuch)
 758       lenrec=ilen(rekord)
 759       iread=index(rekord,lancuch(:lenlan)//"=")
 760 c      print *,"rekord",rekord," lancuch",lancuch
 761 c      print *,"iread",iread," lenlan",lenlan," lenrec",lenrec
 762       if (iread.eq.0) return
 763       iread=iread+lenlan+1
 764       do i=1,dim
 765 c      print *,"iread",iread
 766 c      print *,"|",rekord(iread:iread),"|",iblnk(rekord(iread:iread))
 767       do while (iread.le.lenrec .and. iblnk(rekord(iread:iread)))
 768         iread=iread+1
 769 c      print *,"|",rekord(iread:iread),"|",iblnk(rekord(iread:iread))
 770       enddo
 771 c      print *,"iread",iread
 772       if (iread.gt.lenrec) return
 773       ireade=iread+1
 774 c      print *,"ireade",ireade
 775       do while (ireade.lt.lenrec .and.
 776      &   .not.iblnk(rekord(ireade:ireade)))
 777         ireade=ireade+1
 778       enddo
 779       tablica(i)=rekord(iread:ireade)
 780       iread=ireade+1
 781       enddo
 782       end
 783 c----------------------------------------------------------------------------
 784       subroutine split_string(rekord,tablica,dim,nsub)
 785       implicit none
 786       integer dim,nsub,i,ii,ll,kk
 787       character*(*) tablica(dim)
 788       character*(*) rekord
 789       integer ilen
 790       external ilen
 791       do i=1,dim
 792         tablica(i)=" "
 793       enddo
 794       ii=1
 795       ll = ilen(rekord)
 796       nsub=0
 797       do i=1,dim
 798 C Find the start of term name
 799         kk = 0
 800         do while (ii.le.ll .and. rekord(ii:ii).eq." ")
 801           ii = ii+1
 802         enddo
 803 C Parse the name into TABLICA(i) until blank found
 804         do while (ii.le.ll .and. rekord(ii:ii).ne." ")
 805           kk = kk+1
 806           tablica(i)(kk:kk)=rekord(ii:ii)
 807           ii = ii+1
 808         enddo
 809         if (kk.gt.0) nsub=nsub+1
 810         if (ii.gt.ll) return
 811       enddo
 812       return
 813       end
 814 c--------------------------------------------------------------------------------
 815       integer function iroof(n,m)
 816       ii = n/m
 817       if (ii*m .lt. n) ii=ii+1
 818       iroof = ii
 819       return
 820       end