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