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