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