source/maxlik/src-Fmatch/zscorez.F

   1       program ENE_OPT
   2 c Optimize the UNRES energy function by minimization of a quartic target
   3 c function or by the VMC method.
   4       implicit none
   5 #ifndef ISNAN
   6       external proc_proc
   7 #endif
   8 #ifdef WINPGI
   9 cMS$ATTRIBUTES C ::  proc_proc
  10 #endif
  11       include "DIMENSIONS.ZSCOPT"
  12 #ifdef MPI
  13       include "mpif.h"
  14       integer IERROR,ERRCODE,kolor,key
  15       include "COMMON.MPI"
  16 #endif
  17       include "COMMON.IOUNITS"
  18       include "COMMON.OPTIM"
  19       integer nvarr,iparm
  20       double precision rr,x(max_paropt)
  21       integer idumm
  22       integer i
  23 c      print *,"Starting..."
  24 #ifdef MPI
  25 c      print *,"Initializing MPI..."
  26       call MPI_Init( IERROR )
  27       call MPI_Comm_rank( MPI_COMM_WORLD, me, IERROR )
  28       call MPI_Comm_size( MPI_COMM_WORLD, nprocs, IERROR )
  29 c      print *,"Finished initializing MPI..."
  30       Master = 0
  31       Master1 = 1
  32 c      print *,"Me",me," Master",master," Ierror",ierror
  33       if (ierror.gt.0) then
  34         write(iout,*) "SEVERE ERROR - Can't initialize MPI."
  35         call mpi_finalize(ierror)
  36         stop
  37       endif
  38 #endif
  39 #ifndef ISNAN
  40 c NaNQ initialization
  41       i=-1
  42       rr=dacos(100.0d0)
  43 #ifdef WINPGI
  44       idumm=proc_proc(rr,i)
  45 #else
  46       call proc_proc(rr,i)
  47 #endif
  48 #endif
  49       call initialize
  50 c      print *,"calling openunits"
  51       call openunits
  52 c      print *,"openunits called"
  53       call read_general_data(*10)
  54       write (iout,'(80(1h-)/10x,
  55      & "Maximum likelihood & force matching optimization",
  56      & " of UNRES energy function",
  57      & " v. 01/12/20"/80(1h-))')
  58       call flush(iout)
  59       call cinfo
  60 c      call promienie(*10)
  61       write (iout,*) "Finished READ_GENERAL_DATA"
  62       call flush(iout)
  63       do iparm=1,nparmset
  64         call parmread(iparm,*10)
  65       enddo
  66       call read_pmf_data(*10)
  67       write (iout,*) "Finished parmread"
  68       call flush(iout)
  69       call read_optim_parm(*10)
  70       call print_general_data(*10)
  71       call read_protein_data(*10)
  72       write (iout,*) "Finished READ_PROTEIN_DATA"
  73       call flush(iout)
  74       call read_database(*10)
  75       write (iout,*) "Finished READ_DATABASE"
  76       call flush(iout)
  77 #ifdef MPI
  78 c      write (iout,*) Me,' calling PROC_GROUPS'
  79       call proc_groups
  80 c      write (iout,*) Me,' calling WORK_PARTITION_MAP'
  81 c      call work_partition_map(nvarr)
  82 #endif
  83 #ifdef NEWCORR
  84       call pmfread(*10)
  85       write (iout,*) "Exit pmfread"
  86       call flush(iout)
  87 #endif
  88 c      write (iout,*) "calling proc_data"
  89 c      call flush(iout)
  90       call proc_data(nvarr,x,*10)
  91 c      write (iout,*) "exit proc_data"
  92 c      call flush(iout)
  93       call read_thermal
  94 #ifdef MPI
  95       if (me.eq.Master) then
  96 #endif
  97       call maxlikopt(nvarr,x)
  98 #ifdef MPI
  99       call jebadelko(nvarr)
 100       else
 101         call jebadelko(nvarr)
 102       endif
 103       call bilans
 104       call MPI_Finalize( IERROR )
 105 #else
 106       call bilans
 107 #endif
 108       stop
 109    10 write (iout,*) "Error termination of the program"
 110       call MPI_Finalize( IERROR )
 111       stop
 112       end
 113 c------------------------------------------------------------------------------
 114       subroutine bilans
 115       implicit none
 116       integer ntasks
 117       parameter (ntasks=11)
 118       include "DIMENSIONS.ZSCOPT"
 119 #ifdef MPI
 120       include "mpif.h"
 121       include "COMMON.MPI"
 122       integer IERROR
 123       double precision ttask_all(ntasks)
 124       integer nctask_all(ntasks)
 125 #endif
 126       include "COMMON.IOUNITS"
 127       integer i
 128       double precision ttask
 129       integer nctask
 130       common /timing/ ttask(ntasks),nctask(ntasks)
 131       character*16 tname(ntasks) /"function","gradient",9*''/
 132
 133       write (iout,*)
 134       write (iout,'(80(1h-))')
 135 #ifdef MPI
 136       write (iout,*) "Routine call info from the processor ",me," ..."
 137 #else
 138       write (iout,*) "Routine call info ..."
 139 #endif
 140       write (iout,*)
 141       write (iout,'(65(1h-))')
 142       write (iout,100) "task","   # calls","     total time",
 143      &  "  ave.time/call"
 144       write (iout,'(65(1h-))')
 145       do i=1,ntasks
 146         write (iout,200) tname(i),nctask(i),ttask(i),
 147      &  ttask(i)/(nctask(i)+1.0d-10)
 148       enddo
 149       write (iout,*)
 150 #ifdef MPI
 151       call MPI_Reduce(ttask(1),ttask_all(1),ntasks,
 152      &  MPI_DOUBLE_PRECISION, MPI_SUM, Master, WHAM_COMM,IERROR)
 153       call MPI_Reduce(nctask(1),nctask_all(1),ntasks,
 154      &  MPI_INTEGER, MPI_SUM, Master, WHAM_COMM,IERROR)
 155       if (Me.eq.Master) then
 156       write (iout,'(80(1h-))')
 157       write (iout,*) "Routine call info from all processors ..."
 158       write (iout,*)
 159       write (iout,'(65(1h-))')
 160       write (iout,100) "task","   # calls","     total time",
 161      &  "  ave.time/call"
 162       write (iout,'(65(1h-))')
 163       do i=1,ntasks
 164         write (iout,200) tname(i),nctask_all(i),ttask_all(i),
 165      &    ttask_all(i)/(nctask_all(i)+1.0d-10)
 166       enddo
 167       write (iout,*)
 168       endif
 169 #endif
 170       return
 171   100 format(a,t21,a,t31,a,t46,a)
 172   200 format(a,t21,i10,f15.2,f15.8)
 173       end
 174 c------------------------------------------------------------------------------
 175 #ifdef MPI
 176       subroutine proc_groups
 177 C Split the processors into the Master and Workers group, if needed.
 178       implicit none
 179       include "DIMENSIONS"
 180       include "DIMENSIONS.ZSCOPT"
 181       include "mpif.h"
 182       include "COMMON.IOUNITS"
 183       include "COMMON.MPI"
 184       include "COMMON.VMCPAR"
 185       integer n,chunk,iprot,i,j,ii,remainder
 186       integer kolor,key,ierror,errcode
 187       logical lprint
 188       lprint=.true.
 189 C No secondary structure constraints.
 190       Comm1 = WHAM_COMM
 191       Me1 = Me
 192       Nprocs1 = Nprocs
 193       return
 194       end
 195 c-------------------------------------------------------------------------------
 196       subroutine work_partition(lprint)
 197 c Split the conformations between processors
 198       implicit none
 199       include "DIMENSIONS"
 200       include "DIMENSIONS.ZSCOPT"
 201       include "mpif.h"
 202       include "COMMON.CLASSES"
 203       include "COMMON.IOUNITS"
 204       include "COMMON.MPI"
 205       include "COMMON.VMCPAR"
 206       integer n,chunk,iprot,i,j,ii,remainder
 207       integer kolor,key,ierror,errcode
 208       logical lprint
 209 C
 210 C Divide conformations between processors; for each proteins the first and
 211 C the last conformation to handle by ith processor is stored in
 212 C indstart(i,iprot) and indend(i,iprot), respectively.
 213 C
 214 C First try to assign equal number of conformations to each processor.
 215 C
 216       do iprot=1,nprot
 217         n=ntot_work(iprot)
 218         write (iout,*) "Protein",iprot," n=",n
 219         indstart(0,iprot)=1
 220         chunk = N/nprocs1
 221         scount(0,iprot) = chunk
 222 c        print *,"i",0," indstart",indstart(0,iprot)," scount",
 223 c     &     scount(0,iprot)
 224         do i=1,nprocs1-1
 225           indstart(i,iprot)=chunk+indstart(i-1,iprot)
 226           scount(i,iprot)=scount(i-1,iprot)
 227 c          print *,"i",i," indstart",indstart(i,iprot)," scount",
 228 c     &     scount(i,iprot)
 229         enddo
 230 C
 231 C Determine how many conformations remained yet unassigned.
 232 C
 233         remainder=N-(indstart(nprocs1-1,iprot)
 234      &    +scount(nprocs1-1,iprot)-1)
 235 c        print *,"remainder",remainder
 236 C
 237 C Assign the remainder conformations to consecutive processors, starting
 238 C from the lowest rank; this continues until the list is exhausted.
 239 C
 240         if (remainder .gt. 0) then
 241           do i=1,remainder
 242             scount(i-1,iprot) = scount(i-1,iprot) + 1
 243             indstart(i,iprot) = indstart(i,iprot) + i
 244           enddo
 245           do i=remainder+1,nprocs1-1
 246             indstart(i,iprot) = indstart(i,iprot) + remainder
 247           enddo
 248         endif
 249
 250         indstart(nprocs1,iprot)=N+1
 251         scount(nprocs1,iprot)=0
 252
 253         do i=0,NProcs1
 254           indend(i,iprot)=indstart(i,iprot)+scount(i,iprot)-1
 255           idispl(i,iprot)=indstart(i,iprot)-1
 256         enddo
 257
 258         N=0
 259         do i=0,Nprocs1-1
 260           N=N+indend(i,iprot)-indstart(i,iprot)+1
 261         enddo
 262
 263 c        print *,"N",n," NTOT",ntot_work(iprot)
 264         if (N.ne.ntot_work(iprot)) then
 265           write (iout,*) "!!! Checksum error on processor",me
 266           call flush(iout)
 267           call MPI_Abort( WHAM_COMM, Ierror, Errcode )
 268         endif
 269
 270         ii=0
 271         do i=1,ntot_work(iprot)
 272           if (i.ge.indstart(me1,iprot) .and. i.le.indend(me1,iprot))
 273      &    then
 274             ii=ii+1
 275             i2ii(i,iprot)=ii
 276           else
 277             i2ii(i,iprot)=0
 278           endif
 279 c          write (iout,*) "i",i," iprot",iprot," i2ii",i2ii(i,iprot)
 280         enddo
 281
 282       enddo ! iprot
 283       if (lprint) then
 284         write (iout,*) "Partition of work between processors"
 285         do iprot=1,nprot
 286           write (iout,*) "Protein",iprot
 287 #ifdef DEBUG
 288           write (iout,*) "The i2ii array"
 289           do j=1,ntot_work(iprot)
 290             write (iout,*) j,i2ii(j,iprot)
 291           enddo
 292 #endif
 293           do i=0,nprocs1-1
 294             write (iout,'(a,i5,a,i7,a,i7,a,i7)')
 295      &        "Processor",i," indstart",indstart(i,iprot),
 296      &        " indend",indend(i,iprot)," count",scount(i,iprot)
 297           enddo
 298         enddo
 299       endif
 300       return
 301       end
 302 c------------------------------------------------------------------------------
 303       subroutine jebadelko(nvarr)
 304       implicit none
 305       include "DIMENSIONS"
 306       include "DIMENSIONS.ZSCOPT"
 307       include "mpif.h"
 308       include "COMMON.IOUNITS"
 309       include "COMMON.MPI"
 310       include "COMMON.VMCPAR"
 311       include "COMMON.CLASSES"
 312       include "COMMON.OPTIM"
 313       include "COMMON.WEIGHTS"
 314       include "COMMON.WEIGHTDER"
 315       include "COMMON.ENERGIES"
 316       include "COMMON.TIME1"
 317       include "COMMON.PROTNAME"
 318       include "COMMON.PROTFILES"
 319       include "COMMON.TORSION"
 320       include "COMMON.COMPAR"
 321       integer What, TAG, IERROR, status(MPI_STATUS_SIZE), istop, iprot,
 322      &  nvarr, nf, errcode, ider
 323       integer count
 324       double precision x(max_paropt), g(max_paropt), viol
 325       integer uiparm,i,j
 326       integer ibatch,ib
 327       external fdum
 328       double precision rdum,rdif,chisquare_force
 329       double precision tcpu,t1,t1w,t1_ini,t1w_ini
 330       logical lprn,op,lder
 331       integer iun
 332       lprn=.false.
 333       ider=1
 334       write(iout,*) "Processor",me,me1," called JEBADELKO"
 335       call flush(iout)
 336       if (me.eq.Master) then
 337         istop=0
 338         call func1(nvarr,istop,x)
 339       else
 340 #ifdef DEBUG
 341       write (iout,*) "ELOWEST at slave  starting JEBADELKO"
 342       do iprot=1,nprot
 343         do ibatch=1,natlike(iprot)+2
 344           do ib=1,nbeta(ibatch,iprot)
 345             write (iout,*) "iprot",iprot," ibatch",ibatch," ib",ib,
 346      &         " elowest",elowest(ib,ibatch,iprot)
 347           enddo
 348         enddo
 349       enddo
 350 #endif
 351         istop=1
 352         t1w_ini = MPI_WTIME()
 353         t1_ini = tcpu()
 354         do while (istop.ne.0)
 355 #ifdef DEBUG
 356       write (iout,*) "ELOWEST at slave  calling FUNC1 from JBADELKO"
 357       do iprot=1,nprot
 358         do ibatch=1,natlike(iprot)+2
 359           do ib=1,nbeta(ibatch,iprot)
 360             write (iout,*) "iprot",iprot," ibatch",ibatch," ib",ib,
 361      &         " elowest",elowest(ib,ibatch,iprot)
 362           enddo
 363         enddo
 364       enddo
 365 #endif
 366           call func1(nvarr,istop,x)
 367 c          write (iout,*) "slave: after func1"
 368 c          call flush(iout)
 369 #ifdef NEWCORR
 370           if (istop.eq.1 .and. mod_fourier(nloctyp).gt.0) then
 371             rdum = rdif(nvarr,x,g,ider)
 372 c            write (iout,*) "slave: after rdif"
 373 c            call flush(iout)
 374           endif
 375 #endif
 376         if (istop.eq.1 .or. istop.eq.2 .or. istop.eq.3) then
 377 c          write (iout,*) "slave: calling chisquare_force"
 378           rdum = chisquare_force(nvarr,x,g,ider)
 379 c          write (iout,*) "slave: after chisquare_force"
 380 c          call flush(iout)
 381         endif
 382
 383         enddo
 384         t1w = mpi_wtime() - t1w_ini
 385         t1 = tcpu() - t1_ini
 386         write (iout,*)
 387         write (iout,*) "CPU time",t1," wall clock time",t1w
 388         call flush(iout)
 389       endif
 390 #ifdef DEBUG
 391       write (iout,*)
 392       write (iout,*) "Energies of processor",me
 393       do j=1,nprot
 394         write (iout,*)
 395         write (iout,*) "Protein ",protname(iprot)
 396         do i=1,ntot_work(j)
 397           if (i.ge.indstart(me1,j).and.i.le.indend(me1,j)) then
 398             write(iout,*)i,e_total(i,j),rmstb(i,j),iscore(i,0,j)
 399           endif
 400         enddo
 401       enddo
 402 #endif
 403       do iprot=1,nprot
 404         if (maxlik(iprot)) then
 405         write (iout,*) "Deleting scratchfile",bprotfiles(iprot)
 406         inquire (file=bprotfiles(iprot),number=iun,opened=op)
 407         write (iout,*) "unit",iun," icbase",icbase
 408         if (.not. op) then
 409           open (icbase,file=bprotfiles(iprot),status="old",
 410      &    form="unformatted",access="direct",recl=lenrec(iprot))
 411           close(icbase,status="delete")
 412         else
 413           close(iun,status="delete")
 414         endif
 415         call flush(iout)
 416         if (.not.mod_other_params) then
 417         write (iout,*) "Deleting scratchfile",benefiles(iprot)
 418         inquire (file=benefiles(iprot),number=iun,opened=op)
 419         write (iout,*) "unit",iun," ientout",icbase
 420         if (.not.op) then
 421           open (ientout,file=benefiles(iprot),status="old",
 422      &    form="unformatted",access="direct",recl=lenrec_ene(iprot))
 423           close(ientout,status="delete")
 424         else
 425           close (iun,status="delete")
 426         endif
 427         endif
 428         call flush(iout)
 429         endif
 430         if (fmatch(iprot)) then
 431         write (iout,*) "Deleting scratchfile",bprotfiles_MD(iprot)
 432         inquire (file=bprotfiles_MD(iprot),number=iun,opened=op)
 433         write (iout,*) "unit",iun," icbase",icbase
 434         if (.not. op) then
 435           open (icbase,file=bprotfiles_MD(iprot),status="old",
 436      &    form="unformatted",access="direct",recl=lenrec_MD(iprot))
 437           close(icbase,status="delete")
 438         else
 439           close(iun,status="delete")
 440         endif
 441         call flush(iout)
 442         if (.not.mod_other_params) then
 443         write (iout,*) "Deleting scratchfile",bforcefiles(iprot)
 444         inquire (file=bforcefiles(iprot),number=iun,opened=op)
 445         write (iout,*) "unit",iun," ientout",icbase
 446         if (.not.op) then
 447           open (ientout,file=bforcefiles(iprot),status="old",
 448      &    form="unformatted",access="direct",recl=lenrec_forces(iprot))
 449           close(ientout,status="delete")
 450         else
 451           close (iun,status="delete")
 452         endif
 453         endif
 454         call flush(iout)
 455         endif
 456       enddo
 457       write (iout,*) "Processor",me,"leaves JEBADELKO"
 458       call flush(iout)
 459       return
 460       end
 461 #endif