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