source/unres/src_MD-M/unres.F

   1 CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
   2 C                                                                              C
   3 C                                U N R E S                                     C
   4 C                                                                              C
   5 C Program to carry out conformational search of proteins in an united-residue  C
   6 C approximation.                                                               C
   7 C                                                                              C
   8 CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
   9       implicit real*8 (a-h,o-z)
  10       include 'DIMENSIONS'
  11
  12
  13 #ifdef MPI
  14       include 'mpif.h'
  15       include 'COMMON.SETUP'
  16 #endif
  17       include 'COMMON.TIME1'
  18       include 'COMMON.INTERACT'
  19       include 'COMMON.NAMES'
  20       include 'COMMON.GEO'
  21       include 'COMMON.HEADER'
  22       include 'COMMON.CONTROL'
  23       include 'COMMON.CONTACTS'
  24       include 'COMMON.CHAIN'
  25       include 'COMMON.VAR'
  26       include 'COMMON.IOUNITS'
  27       include 'COMMON.FFIELD'
  28       include 'COMMON.REMD'
  29       include 'COMMON.MD'
  30       include 'COMMON.SBRIDGE'
  31       double precision hrtime,mintime,sectime
  32       character*64 text_mode_calc(-2:14) /'test',
  33      & 'SC rotamer distribution',
  34      & 'Energy evaluation or minimization',
  35      & 'Regularization of PDB structure',
  36      & 'Threading of a sequence on PDB structures',
  37      & 'Monte Carlo (with minimization) ',
  38      & 'Energy minimization of multiple conformations',
  39      & 'Checking energy gradient',
  40      & 'Entropic sampling Monte Carlo (with minimization)',
  41      & 'Energy map',
  42      & 'CSA calculations',
  43      & 'Not used 9',
  44      & 'Not used 10',
  45      & 'Soft regularization of PDB structure',
  46      & 'Mesoscopic molecular dynamics (MD) ',
  47      & 'Not used 13',
  48      & 'Replica exchange molecular dynamics (REMD)'/
  49       external ilen
  50
  51 c      call memmon_print_usage()
  52
  53       call init_task
  54       if (me.eq.king)
  55      & write(iout,*)'### LAST MODIFIED  11/22/15 10:36PM by czarek'
  56       if (me.eq.king) call cinfo
  57 C Read force field parameters and job setup data
  58       call readrtns
  59       call flush(iout)
  60 C
  61       if (me.eq.king .or. .not. out1file) then
  62        write (iout,'(2a/)')
  63      & text_mode_calc(modecalc)(:ilen(text_mode_calc(modecalc))),
  64      & ' calculation.'
  65        if (minim) write (iout,'(a)')
  66      &  'Conformations will be energy-minimized.'
  67        write (iout,'(80(1h*)/)')
  68       endif
  69       call flush(iout)
  70 C
  71       if (modecalc.eq.-2) then
  72         call test
  73         stop
  74       else if (modecalc.eq.-1) then
  75         write(iout,*) "call check_sc_map next"
  76         call check_bond
  77         stop
  78       endif
  79 #ifdef MPI
  80       if (fg_rank.gt.0) then
  81 C Fine-grain slaves just do energy and gradient components.
  82         call ergastulum ! slave workhouse in Latin
  83       else
  84 #endif
  85       if (modecalc.eq.0) then
  86         call exec_eeval_or_minim
  87       else if (modecalc.eq.1) then
  88         call exec_regularize
  89       else if (modecalc.eq.2) then
  90         call exec_thread
  91       else if (modecalc.eq.3 .or. modecalc .eq.6) then
  92         call exec_MC
  93       else if (modecalc.eq.4) then
  94         call exec_mult_eeval_or_minim
  95       else if (modecalc.eq.5) then
  96         call exec_checkgrad
  97       else if (ModeCalc.eq.7) then
  98         call exec_map
  99       else if (ModeCalc.eq.8) then
 100         call exec_CSA
 101       else if (modecalc.eq.11) then
 102         call exec_softreg
 103       else if (modecalc.eq.12) then
 104         call exec_MD
 105       else if (modecalc.eq.14) then
 106 #ifdef MPI
 107         call exec_MREMD
 108 #else
 109         write (iout,*) "Need a parallel version to run MREMD."
 110         stop
 111 #endif
 112       else
 113         write (iout,'(a)') 'This calculation type is not supported',
 114      &   ModeCalc
 115       endif
 116 #ifdef MPI
 117       endif
 118 C Finish task.
 119       if (fg_rank.eq.0) call finish_task
 120 c      call memmon_print_usage()
 121 #ifdef TIMING
 122        call print_detailed_timing
 123 #endif
 124       call MPI_Finalize(ierr)
 125       stop 'Bye Bye...'
 126 #else
 127       call dajczas(tcpu(),hrtime,mintime,sectime)
 128       stop '********** Program terminated normally.'
 129 #endif
 130       end
 131 c--------------------------------------------------------------------------
 132       subroutine exec_MD
 133       include 'DIMENSIONS'
 134 #ifdef MPI
 135       include "mpif.h"
 136 #endif
 137       include 'COMMON.SETUP'
 138       include 'COMMON.CONTROL'
 139       include 'COMMON.IOUNITS'
 140       if (me.eq.king .or. .not. out1file)
 141      &   write (iout,*) "Calling chainbuild"
 142       call chainbuild
 143       call MD
 144       return
 145       end
 146 c---------------------------------------------------------------------------
 147 #ifdef MPI
 148       subroutine exec_MREMD
 149       include 'DIMENSIONS'
 150 #ifdef MPI
 151       include "mpif.h"
 152 #endif
 153       include 'COMMON.SETUP'
 154       include 'COMMON.CONTROL'
 155       include 'COMMON.IOUNITS'
 156       include 'COMMON.REMD'
 157       if (me.eq.king .or. .not. out1file)
 158      &   write (iout,*) "Calling chainbuild"
 159       call chainbuild
 160       if (me.eq.king .or. .not. out1file)
 161      &   write (iout,*) "Calling REMD"
 162       if (remd_mlist) then
 163         call MREMD
 164       else
 165         do i=1,nrep
 166           remd_m(i)=1
 167         enddo
 168         call MREMD
 169       endif
 170       return
 171       end
 172 #endif
 173 c---------------------------------------------------------------------------
 174       subroutine exec_eeval_or_minim
 175       implicit real*8 (a-h,o-z)
 176       include 'DIMENSIONS'
 177 #ifdef MPI
 178       include 'mpif.h'
 179 #endif
 180       include 'COMMON.SETUP'
 181       include 'COMMON.TIME1'
 182       include 'COMMON.INTERACT'
 183       include 'COMMON.NAMES'
 184       include 'COMMON.GEO'
 185       include 'COMMON.HEADER'
 186       include 'COMMON.CONTROL'
 187       include 'COMMON.CONTACTS'
 188       include 'COMMON.CHAIN'
 189       include 'COMMON.VAR'
 190       include 'COMMON.IOUNITS'
 191       include 'COMMON.FFIELD'
 192       include 'COMMON.REMD'
 193       include 'COMMON.MD'
 194       include 'COMMON.SBRIDGE'
 195       common /srutu/ icall
 196       double precision energy(0:n_ene)
 197       double precision energy_long(0:n_ene),energy_short(0:n_ene)
 198       double precision varia(maxvar)
 199       if (indpdb.eq.0)     call chainbuild
 200       print *,'dc',c(1,1)
 201       if (indpdb.ne.0) then
 202       dc(1,0)=c(1,1)
 203       dc(2,0)=c(2,1)
 204       dc(3,0)=c(3,1)
 205       endif
 206 #ifdef MPI
 207       time00=MPI_Wtime()
 208 #else
 209       time00=tcpu()
 210 #endif
 211       call chainbuild_cart
 212       print *,'dc',dc(1,0),dc(2,0),dc(3,0)
 213       if (split_ene) then
 214        print *,"Processor",myrank," after chainbuild"
 215        icall=1
 216        call etotal_long(energy_long(0))
 217        write (iout,*) "Printing long range energy"
 218        call enerprint(energy_long(0))
 219        call etotal_short(energy_short(0))
 220        write (iout,*) "Printing short range energy"
 221        call enerprint(energy_short(0))
 222        do i=0,n_ene
 223          energy(i)=energy_long(i)+energy_short(i)
 224          write (iout,*) i,energy_long(i),energy_short(i),energy(i)
 225        enddo
 226        write (iout,*) "Printing long+short range energy"
 227        call enerprint(energy(0))
 228       endif
 229       call etotal(energy(0))
 230 #ifdef MPI
 231       time_ene=MPI_Wtime()-time00
 232 #else
 233       time_ene=tcpu()-time00
 234 #endif
 235       write (iout,*) "Time for energy evaluation",time_ene
 236       print *,"after etotal"
 237       etota = energy(0)
 238       etot =etota
 239       call enerprint(energy(0))
 240       call hairpin(.true.,nharp,iharp)
 241         print *,'after hairpin'
 242       call secondary2(.true.)
 243         print *,'after secondary'
 244       if (minim) then
 245 crc overlap test
 246         if (overlapsc) then
 247           print *, 'Calling OVERLAP_SC'
 248           call overlap_sc(fail)
 249         endif
 250
 251         if (searchsc) then
 252           call sc_move(2,nres-1,10,1d10,nft_sc,etot)
 253           print *,'SC_move',nft_sc,etot
 254           write(iout,*) 'SC_move',nft_sc,etot
 255         endif
 256
 257         if (dccart) then
 258           print *, 'Calling MINIM_DC'
 259 #ifdef MPI
 260           time1=MPI_WTIME()
 261 #else
 262           time1=tcpu()
 263 #endif
 264           call minim_dc(etot,iretcode,nfun)
 265         else
 266           if (indpdb.ne.0) then
 267             call bond_regular
 268             call chainbuild
 269           endif
 270           call geom_to_var(nvar,varia)
 271           print *,'Calling MINIMIZE.'
 272 #ifdef MPI
 273           time1=MPI_WTIME()
 274 #else
 275           time1=tcpu()
 276 #endif
 277           call minimize(etot,varia,iretcode,nfun)
 278         endif
 279         print *,'SUMSL return code is',iretcode,' eval ',nfun
 280 #ifdef MPI
 281         evals=nfun/(MPI_WTIME()-time1)
 282 #else
 283         evals=nfun/(tcpu()-time1)
 284 #endif
 285         print *,'# eval/s',evals
 286         print *,'refstr=',refstr
 287         call hairpin(.false.,nharp,iharp)
 288         print *,'after hairpin'
 289         call secondary2(.true.)
 290         print *,'after secondary'
 291         call etotal(energy(0))
 292         etot = energy(0)
 293         call enerprint(energy(0))
 294
 295         call intout
 296         call briefout(0,etot)
 297         if (refstr) call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
 298           write (iout,'(a,i3)') 'SUMSL return code:',iretcode
 299           write (iout,'(a,i20)') '# of energy evaluations:',nfun+1
 300           write (iout,'(a,f16.3)')'# of energy evaluations/sec:',evals
 301       else
 302         print *,'refstr=',refstr
 303         if (refstr) call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
 304         call briefout(0,etot)
 305       endif
 306       if (outpdb) call pdbout(etot,titel(:50),ipdb)
 307       if (outmol2) call mol2out(etot,titel(:32))
 308       return
 309       end
 310 c---------------------------------------------------------------------------
 311       subroutine exec_regularize
 312       implicit real*8 (a-h,o-z)
 313       include 'DIMENSIONS'
 314 #ifdef MPI
 315       include 'mpif.h'
 316 #endif
 317       include 'COMMON.SETUP'
 318       include 'COMMON.TIME1'
 319       include 'COMMON.INTERACT'
 320       include 'COMMON.NAMES'
 321       include 'COMMON.GEO'
 322       include 'COMMON.HEADER'
 323       include 'COMMON.CONTROL'
 324       include 'COMMON.CONTACTS'
 325       include 'COMMON.CHAIN'
 326       include 'COMMON.VAR'
 327       include 'COMMON.IOUNITS'
 328       include 'COMMON.FFIELD'
 329       include 'COMMON.REMD'
 330       include 'COMMON.MD'
 331       include 'COMMON.SBRIDGE'
 332       double precision energy(0:n_ene)
 333
 334       call gen_dist_constr
 335       call sc_conf
 336       call intout
 337       call regularize(nct-nnt+1,etot,rms,cref(1,nnt,1),iretcode)
 338       call etotal(energy(0))
 339       energy(0)=energy(0)-energy(14)
 340       etot=energy(0)
 341       call enerprint(energy(0))
 342       call intout
 343       call briefout(0,etot)
 344       if (outpdb) call pdbout(etot,titel(:50),ipdb)
 345       if (outmol2) call mol2out(etot,titel(:32))
 346       if (refstr) call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
 347       write (iout,'(a,i3)') 'SUMSL return code:',iretcode
 348       return
 349       end
 350 c---------------------------------------------------------------------------
 351       subroutine exec_thread
 352       include 'DIMENSIONS'
 353 #ifdef MP
 354       include "mpif.h"
 355 #endif
 356       include "COMMON.SETUP"
 357       call thread_seq
 358       return
 359       end
 360 c---------------------------------------------------------------------------
 361       subroutine exec_MC
 362       implicit real*8 (a-h,o-z)
 363       include 'DIMENSIONS'
 364       character*10 nodeinfo
 365       double precision varia(maxvar)
 366 #ifdef MPI
 367       include "mpif.h"
 368 #endif
 369       include "COMMON.SETUP"
 370       include 'COMMON.CONTROL'
 371       call mcm_setup
 372       if (minim) then
 373 #ifdef MPI
 374         if (modecalc.eq.3) then
 375           call do_mcm(ipar)
 376         else
 377           call entmcm
 378         endif
 379 #else
 380         if (modecalc.eq.3) then
 381           call do_mcm(ipar)
 382         else
 383           call entmcm
 384         endif
 385 #endif
 386       else
 387         call monte_carlo
 388       endif
 389       return
 390       end
 391 c---------------------------------------------------------------------------
 392       subroutine exec_mult_eeval_or_minim
 393       implicit real*8 (a-h,o-z)
 394       include 'DIMENSIONS'
 395 #ifdef MPI
 396       include 'mpif.h'
 397       dimension muster(mpi_status_size)
 398 #endif
 399       include 'COMMON.SETUP'
 400       include 'COMMON.TIME1'
 401       include 'COMMON.INTERACT'
 402       include 'COMMON.NAMES'
 403       include 'COMMON.GEO'
 404       include 'COMMON.HEADER'
 405       include 'COMMON.CONTROL'
 406       include 'COMMON.CONTACTS'
 407       include 'COMMON.CHAIN'
 408       include 'COMMON.VAR'
 409       include 'COMMON.IOUNITS'
 410       include 'COMMON.FFIELD'
 411       include 'COMMON.REMD'
 412       include 'COMMON.MD'
 413       include 'COMMON.SBRIDGE'
 414       double precision varia(maxvar)
 415       dimension ind(6)
 416       double precision energy(0:max_ene)
 417       logical eof
 418       eof=.false.
 419 #ifdef MPI
 420       if(me.ne.king) then
 421         call minim_mcmf
 422         return
 423       endif
 424
 425       close (intin)
 426       open(intin,file=intinname,status='old')
 427       write (istat,'(a5,20a12)')"#    ",
 428      &  (wname(print_order(i)),i=1,nprint_ene)
 429       if (refstr) then
 430         write (istat,'(a5,20a12)')"#    ",
 431      &   (ename(print_order(i)),i=1,nprint_ene),
 432      &   "ETOT total","RMSD","nat.contact","nnt.contact"
 433       else
 434         write (istat,'(a5,20a12)')"#    ",
 435      &    (ename(print_order(i)),i=1,nprint_ene),"ETOT total"
 436       endif
 437
 438       if (.not.minim) then
 439         do while (.not. eof)
 440           if (read_cart) then
 441             read (intin,'(e15.10,e15.5)',end=1100,err=1100) time,ene
 442             call read_x(intin,*11)
 443 #ifdef MPI
 444 c Broadcast the order to compute internal coordinates to the slaves.
 445             if (nfgtasks.gt.1)
 446      &        call MPI_Bcast(6,1,MPI_INTEGER,king,FG_COMM,IERROR)
 447 #endif
 448             call int_from_cart1(.false.)
 449           else
 450             read (intin,'(i5)',end=1100,err=1100) iconf
 451             call read_angles(intin,*11)
 452             call geom_to_var(nvar,varia)
 453             call chainbuild
 454           endif
 455           write (iout,'(a,i7)') 'Conformation #',iconf
 456           call etotal(energy(0))
 457           call briefout(iconf,energy(0))
 458           call enerprint(energy(0))
 459           etot=energy(0)
 460           if (refstr) then
 461             call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
 462             write (istat,'(i5,20(f12.3))') iconf,
 463      &      (energy(print_order(i)),i=1,nprint_ene),etot,
 464      &       rms,frac,frac_nn,co
 465 cjlee end
 466           else
 467             write (istat,'(i5,16(f12.3))') iconf,
 468      &     (energy(print_order(i)),i=1,nprint_ene),etot
 469           endif
 470         enddo
 471 1100    continue
 472         goto 1101
 473       endif
 474
 475       mm=0
 476       imm=0
 477       nft=0
 478       ene0=0.0d0
 479       n=0
 480       iconf=0
 481 c      do n=1,nzsc
 482       do while (.not. eof)
 483         mm=mm+1
 484         if (mm.lt.nodes) then
 485           if (read_cart) then
 486             read (intin,'(e15.10,e15.5)',end=11,err=11) time,ene
 487             call read_x(intin,*11)
 488 #ifdef MPI
 489 c Broadcast the order to compute internal coordinates to the slaves.
 490             if (nfgtasks.gt.1)
 491      &        call MPI_Bcast(6,1,MPI_INTEGER,king,FG_COMM,IERROR)
 492 #endif
 493             call int_from_cart1(.false.)
 494           else
 495             read (intin,'(i5)',end=11,err=11) iconf
 496             call read_angles(intin,*11)
 497             call geom_to_var(nvar,varia)
 498             call chainbuild
 499           endif
 500           write (iout,'(a,i7)') 'Conformation #',iconf
 501           n=n+1
 502          imm=imm+1
 503          ind(1)=1
 504          ind(2)=n
 505          ind(3)=0
 506          ind(4)=0
 507          ind(5)=0
 508          ind(6)=0
 509          ene0=0.0d0
 510          call mpi_send(ind,6,mpi_integer,mm,idint,CG_COMM,
 511      *                  ierr)
 512          call mpi_send(varia,nvar,mpi_double_precision,mm,
 513      *                  idreal,CG_COMM,ierr)
 514          call mpi_send(ene0,1,mpi_double_precision,mm,
 515      *                  idreal,CG_COMM,ierr)
 516 c         print *,'task ',n,' sent to worker ',mm,nvar
 517         else
 518          call mpi_recv(ind,6,mpi_integer,mpi_any_source,idint,
 519      *                 CG_COMM,muster,ierr)
 520          man=muster(mpi_source)
 521 c         print *,'receiving result from worker ',man,' (',iii1,iii,')'
 522          call mpi_recv(varia,nvar,mpi_double_precision,
 523      *               man,idreal,CG_COMM,muster,ierr)
 524          call mpi_recv(ene,1,
 525      *               mpi_double_precision,man,idreal,
 526      *               CG_COMM,muster,ierr)
 527          call mpi_recv(ene0,1,
 528      *               mpi_double_precision,man,idreal,
 529      *               CG_COMM,muster,ierr)
 530 c         print *,'result received from worker ',man,' sending now'
 531
 532           call var_to_geom(nvar,varia)
 533           call chainbuild
 534           call etotal(energy(0))
 535           iconf=ind(2)
 536           write (iout,*)
 537           write (iout,*)
 538           write (iout,'(a,2i7)') 'Conformation #',iconf,ind(5)
 539
 540           etot=energy(0)
 541           call enerprint(energy(0))
 542           call briefout(it,etot)
 543 c          if (minim) call briefout(it,etot)
 544           if (refstr) then
 545             call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
 546             write (istat,'(i5,19(f12.3))') iconf,
 547      &     (energy(print_order(i)),i=1,nprint_ene),etot,
 548      &     rms,frac,frac_nn,co
 549           else
 550             write (istat,'(i5,15(f12.3))') iconf,
 551      &     (energy(print_order(i)),i=1,nprint_ene),etot
 552           endif
 553
 554           imm=imm-1
 555           if (read_cart) then
 556             read (intin,'(e15.10,e15.5)',end=1101,err=1101) time,ene
 557             call read_x(intin,*11)
 558 #ifdef MPI
 559 c Broadcast the order to compute internal coordinates to the slaves.
 560             if (nfgtasks.gt.1)
 561      &        call MPI_Bcast(6,1,MPI_INTEGER,king,FG_COMM,IERROR)
 562 #endif
 563             call int_from_cart1(.false.)
 564           else
 565             read (intin,'(i5)',end=1101,err=1101) iconf
 566             call read_angles(intin,*11)
 567             call geom_to_var(nvar,varia)
 568             call chainbuild
 569           endif
 570           n=n+1
 571           imm=imm+1
 572           ind(1)=1
 573           ind(2)=n
 574           ind(3)=0
 575           ind(4)=0
 576           ind(5)=0
 577           ind(6)=0
 578           call mpi_send(ind,6,mpi_integer,man,idint,CG_COMM,
 579      *                  ierr)
 580           call mpi_send(varia,nvar,mpi_double_precision,man,
 581      *                  idreal,CG_COMM,ierr)
 582           call mpi_send(ene0,1,mpi_double_precision,man,
 583      *                  idreal,CG_COMM,ierr)
 584           nf_mcmf=nf_mcmf+ind(4)
 585           nmin=nmin+1
 586         endif
 587       enddo
 588 11    continue
 589       do j=1,imm
 590         call mpi_recv(ind,6,mpi_integer,mpi_any_source,idint,
 591      *               CG_COMM,muster,ierr)
 592         man=muster(mpi_source)
 593         call mpi_recv(varia,nvar,mpi_double_precision,
 594      *               man,idreal,CG_COMM,muster,ierr)
 595         call mpi_recv(ene,1,
 596      *               mpi_double_precision,man,idreal,
 597      *               CG_COMM,muster,ierr)
 598         call mpi_recv(ene0,1,
 599      *               mpi_double_precision,man,idreal,
 600      *               CG_COMM,muster,ierr)
 601
 602         call var_to_geom(nvar,varia)
 603         call chainbuild
 604         call etotal(energy(0))
 605         iconf=ind(2)
 606         write (iout,*)
 607         write (iout,*)
 608         write (iout,'(a,2i7)') 'Conformation #',iconf,ind(5)
 609
 610         etot=energy(0)
 611         call enerprint(energy(0))
 612         call briefout(it,etot)
 613         if (refstr) then
 614           call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
 615           write (istat,'(i5,19(f12.3))') iconf,
 616      &   (energy(print_order(i)),i=1,nprint_ene),etot,
 617      &   rms,frac,frac_nn,co
 618         else
 619           write (istat,'(i5,15(f12.3))') iconf,
 620      &    (energy(print_order(i)),i=1,nprint_ene),etot
 621         endif
 622         nmin=nmin+1
 623       enddo
 624 1101  continue
 625       do i=1, nodes-1
 626          ind(1)=0
 627          ind(2)=0
 628          ind(3)=0
 629          ind(4)=0
 630          ind(5)=0
 631          ind(6)=0
 632          call mpi_send(ind,6,mpi_integer,i,idint,CG_COMM,
 633      *                  ierr)
 634       enddo
 635 #else
 636       close (intin)
 637       open(intin,file=intinname,status='old')
 638       write (istat,'(a5,20a12)')"#    ",
 639      &   (wname(print_order(i)),i=1,nprint_ene)
 640       write (istat,'("#    ",20(1pe12.4))')
 641      &   (weights(print_order(i)),i=1,nprint_ene)
 642       if (refstr) then
 643         write (istat,'(a5,20a12)')"#    ",
 644      &   (ename(print_order(i)),i=1,nprint_ene),
 645      &   "ETOT total","RMSD","nat.contact","nnt.contact"
 646       else
 647         write (istat,'(a5,14a12)')"#    ",
 648      &   (ename(print_order(i)),i=1,nprint_ene),"ETOT total"
 649       endif
 650       do while (.not. eof)
 651           if (read_cart) then
 652             read (intin,'(e15.10,e15.5)',end=11,err=11) time,ene
 653             call read_x(intin,*11)
 654 #ifdef MPI
 655 c Broadcast the order to compute internal coordinates to the slaves.
 656             if (nfgtasks.gt.1)
 657      &        call MPI_Bcast(6,1,MPI_INTEGER,king,FG_COMM,IERROR)
 658 #endif
 659             call int_from_cart1(.false.)
 660           else
 661             read (intin,'(i5)',end=11,err=11) iconf
 662             call read_angles(intin,*11)
 663             call geom_to_var(nvar,varia)
 664             call chainbuild
 665           endif
 666         write (iout,'(a,i7)') 'Conformation #',iconf
 667         if (minim) call minimize(etot,varia,iretcode,nfun)
 668         call etotal(energy(0))
 669
 670         etot=energy(0)
 671         call enerprint(energy(0))
 672         if (minim) call briefout(it,etot)
 673         if (refstr) then
 674           call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
 675           write (istat,'(i5,18(f12.3))') iconf,
 676      &   (energy(print_order(i)),i=1,nprint_ene),
 677      &   etot,rms,frac,frac_nn,co
 678 cjlee end
 679         else
 680           write (istat,'(i5,14(f12.3))') iconf,
 681      &   (energy(print_order(i)),i=1,nprint_ene),etot
 682         endif
 683       enddo
 684    11 continue
 685 #endif
 686       return
 687       end
 688 c---------------------------------------------------------------------------
 689       subroutine exec_checkgrad
 690       implicit real*8 (a-h,o-z)
 691       include 'DIMENSIONS'
 692 #ifdef MPI
 693       include 'mpif.h'
 694 #endif
 695       include 'COMMON.SETUP'
 696       include 'COMMON.TIME1'
 697       include 'COMMON.INTERACT'
 698       include 'COMMON.NAMES'
 699       include 'COMMON.GEO'
 700       include 'COMMON.HEADER'
 701       include 'COMMON.CONTROL'
 702       include 'COMMON.CONTACTS'
 703       include 'COMMON.CHAIN'
 704       include 'COMMON.VAR'
 705       include 'COMMON.IOUNITS'
 706       include 'COMMON.FFIELD'
 707       include 'COMMON.REMD'
 708       include 'COMMON.MD'
 709       include 'COMMON.SBRIDGE'
 710       common /srutu/ icall
 711       double precision energy(0:max_ene)
 712 c      do i=2,nres
 713 c        vbld(i)=vbld(i)+ran_number(-0.1d0,0.1d0)
 714 c        if (itype(i).ne.10)
 715 c     &      vbld(i+nres)=vbld(i+nres)+ran_number(-0.001d0,0.001d0)
 716 c      enddo
 717       if (indpdb.eq.0) call chainbuild
 718 c      do i=0,nres
 719 c        do j=1,3
 720 c          dc(j,i)=dc(j,i)+ran_number(-0.2d0,0.2d0)
 721 c        enddo
 722 c      enddo
 723 c      do i=1,nres-1
 724 c        if (itype(i).ne.10) then
 725 c          do j=1,3
 726 c            dc(j,i+nres)=dc(j,i+nres)+ran_number(-0.2d0,0.2d0)
 727 c          enddo
 728 c        endif
 729 c      enddo
 730 c      do j=1,3
 731 c        dc(j,0)=ran_number(-0.2d0,0.2d0)
 732 c      enddo
 733       if (constr_homology.eq.0) then
 734       usampl=.true.
 735       totT=1.d0
 736       eq_time=0.0d0
 737       call read_fragments
 738       endif
 739       call chainbuild_cart
 740       call cartprint
 741       call intout
 742       icall=1
 743       call etotal(energy(0))
 744       etot = energy(0)
 745       call enerprint(energy(0))
 746       write (iout,*) "Uconst",Uconst," Uconst_back",uconst_back
 747       print *,'icheckgrad=',icheckgrad
 748       goto (10,20,30) icheckgrad
 749   10  call check_ecartint
 750       return
 751   20  call check_cartgrad
 752       return
 753   30  call check_eint
 754       return
 755       end
 756 c---------------------------------------------------------------------------
 757       subroutine exec_map
 758 C Energy maps
 759       call map_read
 760       call map
 761       return
 762       end
 763 c---------------------------------------------------------------------------
 764       subroutine exec_CSA
 765 #ifdef MPI
 766       include "mpif.h"
 767 #endif
 768       include 'DIMENSIONS'
 769       include 'COMMON.IOUNITS'
 770 C Conformational Space Annealling programmed by Jooyoung Lee.
 771 C This method works only with parallel machines!
 772 #ifdef MPI
 773       call together
 774 #else
 775       write (iout,*) "CSA works on parallel machines only"
 776 #endif
 777       return
 778       end
 779 c---------------------------------------------------------------------------
 780       subroutine exec_softreg
 781       include 'DIMENSIONS'
 782       include 'COMMON.IOUNITS'
 783       include 'COMMON.CONTROL'
 784       double precision energy(0:max_ene)
 785       call chainbuild
 786       call etotal(energy(0))
 787       call enerprint(energy(0))
 788       if (.not.lsecondary) then
 789         write(iout,*) 'Calling secondary structure recognition'
 790         call secondary2(debug)
 791       else
 792         write(iout,*) 'Using secondary structure supplied in pdb'
 793       endif
 794
 795       call softreg
 796
 797       call etotal(energy(0))
 798       etot=energy(0)
 799       call enerprint(energy(0))
 800       call intout
 801       call briefout(0,etot)
 802       call secondary2(.true.)
 803       if (refstr) call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
 804       return
 805       end