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