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