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*64 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*64 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,bprotfile_temp
 419       character*3 liczba
 420       character*2 liczba1
 421       integer i,j,ii,jj(maxslice),k,kk(maxslice),l,
 422      &  ll(maxslice),mm(maxslice),if
 423       integer nrec,nlines,iscor,iunit,islice
 424       double precision energ
 425       integer ilen,iroof
 426       external ilen,iroof
 427       double precision rmsdev,energia(0:max_ene),efree,eini,temp
 428       double precision prop(maxQ)
 429       integer ntot_all(maxslice,0:maxprocs-1)
 430       integer iparm,ib,iib,ir,nprop,nthr,npars
 431       double precision etot,time
 432       integer ixdrf,iret
 433       logical lerr,linit
 434
 435       lenrec1=12*(nres+nct-nnt+1)+4*(2*nss+2)+24
 436       lenrec2=12*(nres+nct-nnt+1)+4*(2*nss+2)+24+8*nQ
 437       lenrec=lenrec2+8
 438       write (iout,*) "lenrec",lenrec," lenrec1",lenrec1,
 439      &  " lenrec2",lenrec2
 440
 441       do i=1,nQ
 442         prop(i)=0.0d0
 443       enddo
 444       do islice=1,nslice
 445         ll(islice)=0
 446         mm(islice)=0
 447       enddo
 448       write (iout,*) "nparmset",nparmset
 449       if (hamil_rep) then
 450         npars=1
 451       else
 452         npars=nparmset
 453       endif
 454       do iparm=1,npars
 455
 456       if (replica(iparm)) then
 457         nthr = 1
 458       else
 459         nthr = nT_h(iparm)
 460       endif
 461
 462       do ib=1,nthr
 463       do iR=1,nRR(ib,iparm)
 464
 465       write (iout,*) "iparm",iparm," ib",ib," iR",iR," nQ",nQ
 466       do islice=1,nslice
 467         jj(islice)=0
 468         kk(islice)=0
 469       enddo
 470
 471       IF (NFILE_BIN(iR,ib,iparm).GT.0) THEN
 472 c Read conformations from binary DA files (one per batch) and write them to
 473 c a binary DA scratchfile.
 474         write (liczba,'(bz,i3.3)') me
 475         do if=1,nfile_bin(iR,ib,iparm)
 476           nazwa=protfiles(if,1,iR,ib,iparm)
 477      &     (:ilen(protfiles(if,1,iR,ib,iparm)))//".bx"
 478           open (ientin,file=nazwa,status="old",form="unformatted",
 479      &     access="direct",recl=lenrec2,err=1111)
 480           ii=0
 481           do islice=1,nslice
 482             call opentmp(islice,ientout,bprotfile_temp)
 483             call bxread(nazwa,ii,jj(islice),kk(islice),ll(islice),
 484      &        mm(islice),iR,ib,iparm)
 485             close(ientout)
 486           enddo
 487           close(ientin)
 488         enddo
 489       ENDIF ! NFILE_BIN>0
 490 c
 491       IF (NFILE_ASC(iR,ib,iparm).GT.0) THEN
 492 c Read conformations from multiple ASCII int files and write them to a binary
 493 c DA scratchfile.
 494         do if=1,nfile_asc(iR,ib,iparm)
 495           nazwa=protfiles(if,2,iR,ib,iparm)
 496      &     (:ilen(protfiles(if,2,iR,ib,iparm)))//".x"
 497           open(unit=ientin,file=nazwa,status='old',err=1111)
 498           write(iout,*) "reading ",nazwa(:ilen(nazwa))
 499           ii=0
 500           call xread(nazwa,ii,jj,kk,ll,mm,iR,ib,iparm)
 501         enddo ! if
 502       ENDIF
 503       IF (NFILE_CX(iR,ib,iparm).gt.0) THEN
 504 c Read conformations from cx files and write them to a binary
 505 c DA scratchfile.
 506         do if=1,nfile_cx(iR,ib,iparm)
 507           nazwa=protfiles(if,2,iR,ib,iparm)
 508      &     (:ilen(protfiles(if,2,iR,ib,iparm)))//".cx"
 509           write(iout,*) "reading ",nazwa(:ilen(nazwa))
 510           ii=0
 511           print *,"Calling cxread"
 512           call cxread(nazwa,ii,jj,kk,ll,mm,iR,ib,iparm,
 513      &       *1111)
 514           close(ientout)
 515           write (iout,*) "exit cxread"
 516           call flush(iout)
 517         enddo
 518       ENDIF
 519
 520       do islice=1,nslice
 521         stot(islice)=stot(islice)+jj(islice)
 522       enddo
 523
 524       enddo
 525       enddo
 526       write (iout,*) "IPARM",iparm
 527       enddo
 528
 529       if (nslice.eq.1) then
 530 #ifdef MPI
 531         write (liczba,'(bz,i3.3)') me
 532         bprotfile_temp = scratchdir(:ilen(scratchdir))//"/"//
 533      &    prefix(:ilen(prefix))//liczba//".xbin.tmp"
 534 #else
 535         bprotfile_temp = scratchdir(:ilen(scratchdir))//
 536      &     "/"//prefix(:ilen(prefix))//".xbin.tmp"
 537 #endif
 538         write(iout,*) mm(1)," conformations read",ll(1),
 539      &    " conformations written to ",
 540      &    bprotfile_temp(:ilen(bprotfile_temp))
 541       else
 542         do islice=1,nslice
 543           write (liczba1,'(bz,i2.2)') islice
 544 #ifdef MPI
 545           write (liczba,'(bz,i3.3)') me
 546           bprotfile_temp = scratchdir(:ilen(scratchdir))//"/"//
 547      &      prefix(:ilen(prefix))//liczba//".xbin.tmp"//liczba1
 548 #else
 549           bprotfile_temp = scratchdir(:ilen(scratchdir))//
 550      &       "/"//prefix(:ilen(prefix))//".xbin.tmp"//liczba1
 551 #endif
 552           write(iout,*) mm(islice)," conformations read",ll(islice),
 553      &    " conformations written to ",
 554      &    bprotfile_temp(:ilen(bprotfile_temp))
 555         enddo
 556       endif
 557
 558 #ifdef MPI
 559 c Check if everyone has the same number of conformations
 560       call MPI_Allgather(stot(1),maxslice,MPI_INTEGER,
 561      &  ntot_all(1,0),maxslice,MPI_INTEGER,MPI_Comm_World,IERROR)
 562       lerr=.false.
 563       do i=0,nprocs-1
 564         if (i.ne.me) then
 565           do islice=1,nslice
 566           if (stot(islice).ne.ntot_all(islice,i)) then
 567             write (iout,*) "Number of conformations at processor",i,
 568      &       " differs from that at processor",me,
 569      &       stot(islice),ntot_all(islice,i)," slice",islice
 570             lerr = .true.
 571           endif
 572           enddo
 573         endif
 574       enddo
 575       if (lerr) then
 576         write (iout,*)
 577         write (iout,*) "Numbers of conformations read by processors"
 578         write (iout,*)
 579         do i=0,nprocs-1
 580           write (iout,'(8i10)') i,(ntot_all(islice,i),islice=1,nslice)
 581         enddo
 582         write (iout,*) "Calculation terminated."
 583         call flush(iout)
 584         return1
 585       endif
 586       do islice=1,nslice
 587         ntot(islice)=stot(islice)
 588       enddo
 589       return
 590 #endif
 591  1111 write(iout,*) "Error opening coordinate file ",nazwa(:ilen(nazwa))
 592       call flush(iout)
 593       return1
 594       end
 595 c------------------------------------------------------------------------------
 596       subroutine card_concat(card,to_upper)
 597       implicit none
 598       include 'DIMENSIONS.ZSCOPT'
 599       include "COMMON.IOUNITS"
 600       character*(*) card
 601       character*80 karta,ucase
 602       logical to_upper
 603       integer ilen
 604       external ilen
 605       read (inp,'(a)') karta
 606       if (to_upper) karta=ucase(karta)
 607       card=' '
 608       do while (karta(80:80).eq.'&')
 609         card=card(:ilen(card)+1)//karta(:79)
 610         read (inp,'(a)') karta
 611         if (to_upper) karta=ucase(karta)
 612       enddo
 613       card=card(:ilen(card)+1)//karta
 614       return
 615       end
 616 c------------------------------------------------------------------------------
 617       subroutine readi(rekord,lancuch,wartosc,default)
 618       implicit none
 619       character*(*) rekord,lancuch
 620       integer wartosc,default
 621       integer ilen,iread
 622       external ilen
 623       iread=index(rekord,lancuch(:ilen(lancuch))//"=")
 624       if (iread.eq.0) then
 625         wartosc=default
 626         return
 627       endif
 628       iread=iread+ilen(lancuch)+1
 629       read (rekord(iread:),*) wartosc
 630       return
 631       end
 632 c----------------------------------------------------------------------------
 633       subroutine reada(rekord,lancuch,wartosc,default)
 634       implicit none
 635       character*(*) rekord,lancuch
 636       character*80 aux
 637       double precision wartosc,default
 638       integer ilen,iread
 639       external ilen
 640       iread=index(rekord,lancuch(:ilen(lancuch))//"=")
 641       if (iread.eq.0) then
 642         wartosc=default
 643         return
 644       endif
 645       iread=iread+ilen(lancuch)+1
 646       read (rekord(iread:),*) wartosc
 647       return
 648       end
 649 c----------------------------------------------------------------------------
 650       subroutine multreadi(rekord,lancuch,tablica,dim,default)
 651       implicit none
 652       integer dim,i
 653       integer tablica(dim),default
 654       character*(*) rekord,lancuch
 655       character*80 aux
 656       integer ilen,iread
 657       external ilen
 658       do i=1,dim
 659         tablica(i)=default
 660       enddo
 661       iread=index(rekord,lancuch(:ilen(lancuch))//"=")
 662       if (iread.eq.0) return
 663       iread=iread+ilen(lancuch)+1
 664       read (rekord(iread:),*,end=10,err=10) (tablica(i),i=1,dim)
 665    10 return
 666       end
 667 c----------------------------------------------------------------------------
 668       subroutine multreada(rekord,lancuch,tablica,dim,default)
 669       implicit none
 670       integer dim,i
 671       double precision tablica(dim),default
 672       character*(*) rekord,lancuch
 673       character*80 aux
 674       integer ilen,iread
 675       external ilen
 676       do i=1,dim
 677         tablica(i)=default
 678       enddo
 679       iread=index(rekord,lancuch(:ilen(lancuch))//"=")
 680       if (iread.eq.0) return
 681       iread=iread+ilen(lancuch)+1
 682       read (rekord(iread:),*,end=10,err=10) (tablica(i),i=1,dim)
 683    10 return
 684       end
 685 c----------------------------------------------------------------------------
 686       subroutine reads(rekord,lancuch,wartosc,default)
 687       implicit none
 688       character*(*) rekord,lancuch,wartosc,default
 689       character*80 aux
 690       integer ilen,lenlan,lenrec,iread,ireade
 691       external ilen
 692       logical iblnk
 693       external iblnk
 694       lenlan=ilen(lancuch)
 695       lenrec=ilen(rekord)
 696       iread=index(rekord,lancuch(:lenlan)//"=")
 697 c      print *,"rekord",rekord," lancuch",lancuch
 698 c      print *,"iread",iread," lenlan",lenlan," lenrec",lenrec
 699       if (iread.eq.0) then
 700         wartosc=default
 701         return
 702       endif
 703       iread=iread+lenlan+1
 704 c      print *,"iread",iread
 705 c      print *,"|",rekord(iread:iread),"|",iblnk(rekord(iread:iread))
 706       do while (iread.le.lenrec .and. iblnk(rekord(iread:iread)))
 707         iread=iread+1
 708 c      print *,"|",rekord(iread:iread),"|",iblnk(rekord(iread:iread))
 709       enddo
 710 c      print *,"iread",iread
 711       if (iread.gt.lenrec) then
 712          wartosc=default
 713         return
 714       endif
 715       ireade=iread+1
 716 c      print *,"ireade",ireade
 717       do while (ireade.lt.lenrec .and.
 718      &   .not.iblnk(rekord(ireade:ireade)))
 719         ireade=ireade+1
 720       enddo
 721       wartosc=rekord(iread:ireade)
 722       return
 723       end
 724 c----------------------------------------------------------------------------
 725       subroutine multreads(rekord,lancuch,tablica,dim,default)
 726       implicit none
 727       integer dim,i
 728       character*(*) rekord,lancuch,tablica(dim),default
 729       character*80 aux
 730       integer ilen,lenlan,lenrec,iread,ireade
 731       external ilen
 732       logical iblnk
 733       external iblnk
 734       do i=1,dim
 735         tablica(i)=default
 736       enddo
 737       lenlan=ilen(lancuch)
 738       lenrec=ilen(rekord)
 739       iread=index(rekord,lancuch(:lenlan)//"=")
 740 c      print *,"rekord",rekord," lancuch",lancuch
 741 c      print *,"iread",iread," lenlan",lenlan," lenrec",lenrec
 742       if (iread.eq.0) return
 743       iread=iread+lenlan+1
 744       do i=1,dim
 745 c      print *,"iread",iread
 746 c      print *,"|",rekord(iread:iread),"|",iblnk(rekord(iread:iread))
 747       do while (iread.le.lenrec .and. iblnk(rekord(iread:iread)))
 748         iread=iread+1
 749 c      print *,"|",rekord(iread:iread),"|",iblnk(rekord(iread:iread))
 750       enddo
 751 c      print *,"iread",iread
 752       if (iread.gt.lenrec) return
 753       ireade=iread+1
 754 c      print *,"ireade",ireade
 755       do while (ireade.lt.lenrec .and.
 756      &   .not.iblnk(rekord(ireade:ireade)))
 757         ireade=ireade+1
 758       enddo
 759       tablica(i)=rekord(iread:ireade)
 760       iread=ireade+1
 761       enddo
 762       end
 763 c----------------------------------------------------------------------------
 764       subroutine split_string(rekord,tablica,dim,nsub)
 765       implicit none
 766       integer dim,nsub,i,ii,ll,kk
 767       character*(*) tablica(dim)
 768       character*(*) rekord
 769       integer ilen
 770       external ilen
 771       do i=1,dim
 772         tablica(i)=" "
 773       enddo
 774       ii=1
 775       ll = ilen(rekord)
 776       nsub=0
 777       do i=1,dim
 778 C Find the start of term name
 779         kk = 0
 780         do while (ii.le.ll .and. rekord(ii:ii).eq." ")
 781           ii = ii+1
 782         enddo
 783 C Parse the name into TABLICA(i) until blank found
 784         do while (ii.le.ll .and. rekord(ii:ii).ne." ")
 785           kk = kk+1
 786           tablica(i)(kk:kk)=rekord(ii:ii)
 787           ii = ii+1
 788         enddo
 789         if (kk.gt.0) nsub=nsub+1
 790         if (ii.gt.ll) return
 791       enddo
 792       return
 793       end
 794 c--------------------------------------------------------------------------------
 795       integer function iroof(n,m)
 796       ii = n/m
 797       if (ii*m .lt. n) ii=ii+1
 798       iroof = ii
 799       return
 800       end