projects
/
unres.git
/ blobdiff
commit
grep
author
committer
pickaxe
?
search:
re
summary
|
shortlog
|
log
|
commit
|
commitdiff
|
tree
raw
|
inline
| side by side
pdbread & boxshift
[unres.git]
/
source
/
unres
/
src-HCD-5D
/
unres.F
diff --git
a/source/unres/src-HCD-5D/unres.F
b/source/unres/src-HCD-5D/unres.F
index
5da3f8e
..
08caade
100644
(file)
--- a/
source/unres/src-HCD-5D/unres.F
+++ b/
source/unres/src-HCD-5D/unres.F
@@
-6,7
+6,7
@@
C Program to carry out conformational search of proteins in an united-residue C
C approximation. C
C C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C approximation. C
C C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
- implicit real*8 (a-h,o-z)
+ implicit none
include 'DIMENSIONS'
include 'DIMENSIONS'
@@
-20,7
+20,7
@@
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
include 'COMMON.GEO'
include 'COMMON.HEADER'
include 'COMMON.CONTROL'
include 'COMMON.GEO'
include 'COMMON.HEADER'
include 'COMMON.CONTROL'
- include 'COMMON.CONTACTS'
+c include 'COMMON.CONTACTS'
include 'COMMON.CHAIN'
include 'COMMON.VAR'
include 'COMMON.IOUNITS'
include 'COMMON.CHAIN'
include 'COMMON.VAR'
include 'COMMON.IOUNITS'
@@
-46,7
+46,9
@@
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
& 'Mesoscopic molecular dynamics (MD) ',
& 'Not used 13',
& 'Replica exchange molecular dynamics (REMD)'/
& 'Mesoscopic molecular dynamics (MD) ',
& 'Not used 13',
& 'Replica exchange molecular dynamics (REMD)'/
+ integer ilen
external ilen
external ilen
+ integer ierr
c call memmon_print_usage()
c call memmon_print_usage()
@@
-71,8
+73,8
@@
c write (iout,*) "After readrtns"
call flush(iout)
C
if (modecalc.eq.-2) then
call flush(iout)
C
if (modecalc.eq.-2) then
- call test
- stop
+c call test
+c stop
else if (modecalc.eq.-1) then
write(iout,*) "call check_sc_map next"
call check_bond
else if (modecalc.eq.-1) then
write(iout,*) "call check_sc_map next"
call check_bond
@@
-134,6
+136,7
@@
c call memmon_print_usage()
end
c--------------------------------------------------------------------------
subroutine exec_MD
end
c--------------------------------------------------------------------------
subroutine exec_MD
+ implicit none
include 'DIMENSIONS'
#ifdef MPI
include "mpif.h"
include 'DIMENSIONS'
#ifdef MPI
include "mpif.h"
@@
-156,6
+159,7
@@
c endif
c---------------------------------------------------------------------------
#ifdef MPI
subroutine exec_MREMD
c---------------------------------------------------------------------------
#ifdef MPI
subroutine exec_MREMD
+ implicit none
include 'DIMENSIONS'
#ifdef MPI
include "mpif.h"
include 'DIMENSIONS'
#ifdef MPI
include "mpif.h"
@@
-164,6
+168,7
@@
c---------------------------------------------------------------------------
include 'COMMON.CONTROL'
include 'COMMON.IOUNITS'
include 'COMMON.REMD'
include 'COMMON.CONTROL'
include 'COMMON.IOUNITS'
include 'COMMON.REMD'
+ integer i
if (me.eq.king .or. .not. out1file)
& write (iout,*) "Calling chainbuild"
call chainbuild
if (me.eq.king .or. .not. out1file)
& write (iout,*) "Calling chainbuild"
call chainbuild
@@
-182,7
+187,7
@@
c---------------------------------------------------------------------------
#endif
c---------------------------------------------------------------------------
subroutine exec_eeval_or_minim
#endif
c---------------------------------------------------------------------------
subroutine exec_eeval_or_minim
- implicit real*8 (a-h,o-z)
+ implicit none
include 'DIMENSIONS'
#ifdef MPI
include 'mpif.h'
include 'DIMENSIONS'
#ifdef MPI
include 'mpif.h'
@@
-194,7
+199,7
@@
c---------------------------------------------------------------------------
include 'COMMON.GEO'
include 'COMMON.HEADER'
include 'COMMON.CONTROL'
include 'COMMON.GEO'
include 'COMMON.HEADER'
include 'COMMON.CONTROL'
- include 'COMMON.CONTACTS'
+c include 'COMMON.CONTACTS'
include 'COMMON.CHAIN'
include 'COMMON.VAR'
include 'COMMON.IOUNITS'
include 'COMMON.CHAIN'
include 'COMMON.VAR'
include 'COMMON.IOUNITS'
@@
-202,10
+207,22
@@
c---------------------------------------------------------------------------
include 'COMMON.REMD'
include 'COMMON.MD'
include 'COMMON.SBRIDGE'
include 'COMMON.REMD'
include 'COMMON.MD'
include 'COMMON.SBRIDGE'
+ integer i,icall,iretcode,nfun
common /srutu/ icall
common /srutu/ icall
- double precision energy(0:n_ene)
+ integer nharp,iharp(4,maxres/3)
+ integer nft_sc
+ logical fail
+ double precision energy(0:n_ene),etot,etota
double precision energy_long(0:n_ene),energy_short(0:n_ene)
double precision energy_long(0:n_ene),energy_short(0:n_ene)
+ double precision rms,frac,frac_nn,co
double precision varia(maxvar)
double precision varia(maxvar)
+ double precision time00,time1,time_ene,evals
+#ifdef LBFGS
+ character*9 status
+ integer niter
+ common /lbfgstat/ status,niter,nfun
+#endif
+ integer ilen
if (indpdb.eq.0) call chainbuild
if (indpdb.ne.0) then
dc(1,0)=c(1,1)
if (indpdb.eq.0) call chainbuild
if (indpdb.ne.0) then
dc(1,0)=c(1,1)
@@
-247,20
+264,28
@@
c call flush(iout)
time_ene=tcpu()-time00
#endif
write (iout,*) "Time for energy evaluation",time_ene
time_ene=tcpu()-time00
#endif
write (iout,*) "Time for energy evaluation",time_ene
- print *,"after etotal"
+c print *,"after etotal"
etota = energy(0)
etot =etota
call enerprint(energy(0))
call hairpin(.true.,nharp,iharp)
etota = energy(0)
etot =etota
call enerprint(energy(0))
call hairpin(.true.,nharp,iharp)
- print *,'after hairpin'
+c print *,'after hairpin'
call secondary2(.true.)
call secondary2(.true.)
- print *,'after secondary'
+c print *,'after secondary'
if (minim) then
crc overlap test
if (minim) then
crc overlap test
+ if (indpdb.ne.0 .and. .not.dccart) then
+ call bond_regular
+ call chainbuild_extconf
+ call etotal(energy(0))
+ write (iout,*) "After bond regularization"
+ call enerprint(energy(0))
+ endif
+
if (overlapsc) then
if (overlapsc) then
- print *, 'Calling OVERLAP_SC'
+ write (iout,*) 'Calling OVERLAP_SC'
call overlap_sc(fail)
call overlap_sc(fail)
- print *,"After overlap_sc"
+ write (iout,*) "After overlap_sc"
endif
if (searchsc) then
endif
if (searchsc) then
@@
-278,12
+303,8
@@
crc overlap test
#endif
call minim_dc(etot,iretcode,nfun)
else
#endif
call minim_dc(etot,iretcode,nfun)
else
- if (indpdb.ne.0) then
- call bond_regular
- call chainbuild_extconf
- endif
call geom_to_var(nvar,varia)
call geom_to_var(nvar,varia)
- print *,'Calling MINIMIZE.'
+c print *,'Calling MINIMIZE.'
#ifdef MPI
time1=MPI_WTIME()
#else
#ifdef MPI
time1=MPI_WTIME()
#else
@@
-291,7
+312,11
@@
crc overlap test
#endif
call minimize(etot,varia,iretcode,nfun)
endif
#endif
call minimize(etot,varia,iretcode,nfun)
endif
+#ifdef LBFGS
+ print *,'LBFGS return code is',status,' eval ',nfun
+#else
print *,'SUMSL return code is',iretcode,' eval ',nfun
print *,'SUMSL return code is',iretcode,' eval ',nfun
+#endif
#ifdef MPI
evals=nfun/(MPI_WTIME()-time1)
#else
#ifdef MPI
evals=nfun/(MPI_WTIME()-time1)
#else
@@
-308,11
+333,22
@@
crc overlap test
call enerprint(energy(0))
call intout
call enerprint(energy(0))
call intout
- call briefout(0,etot)
+ if (out_int) call briefout(0,etot)
+ if (out_cart) then
+ cartname=prefix(:ilen(prefix))//'.x'
+ potE=etot
+ call cartoutx(0.0d0)
+ endif
if (refstr) call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
if (refstr) call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
+#ifdef LBFGS
+ write (iout,'(a,a9)') 'LBFGS return code:',status
+ write (iout,'(a,i20)') '# of energy evaluations:',nfun+1
+ write (iout,'(a,f16.3)')'# of energy evaluations/sec:',evals
+#else
write (iout,'(a,i3)') 'SUMSL return code:',iretcode
write (iout,'(a,i20)') '# of energy evaluations:',nfun+1
write (iout,'(a,f16.3)')'# of energy evaluations/sec:',evals
write (iout,'(a,i3)') 'SUMSL return code:',iretcode
write (iout,'(a,i20)') '# of energy evaluations:',nfun+1
write (iout,'(a,f16.3)')'# of energy evaluations/sec:',evals
+#endif
else
print *,'refstr=',refstr
if (refstr) call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
else
print *,'refstr=',refstr
if (refstr) call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
@@
-324,7
+360,7
@@
crc overlap test
end
c---------------------------------------------------------------------------
subroutine exec_regularize
end
c---------------------------------------------------------------------------
subroutine exec_regularize
- implicit real*8 (a-h,o-z)
+ implicit none
include 'DIMENSIONS'
#ifdef MPI
include 'mpif.h'
include 'DIMENSIONS'
#ifdef MPI
include 'mpif.h'
@@
-336,7
+372,7
@@
c---------------------------------------------------------------------------
include 'COMMON.GEO'
include 'COMMON.HEADER'
include 'COMMON.CONTROL'
include 'COMMON.GEO'
include 'COMMON.HEADER'
include 'COMMON.CONTROL'
- include 'COMMON.CONTACTS'
+c include 'COMMON.CONTACTS'
include 'COMMON.CHAIN'
include 'COMMON.VAR'
include 'COMMON.IOUNITS'
include 'COMMON.CHAIN'
include 'COMMON.VAR'
include 'COMMON.IOUNITS'
@@
-345,6
+381,13
@@
c---------------------------------------------------------------------------
include 'COMMON.MD'
include 'COMMON.SBRIDGE'
double precision energy(0:n_ene)
include 'COMMON.MD'
include 'COMMON.SBRIDGE'
double precision energy(0:n_ene)
+ double precision etot,rms,frac,frac_nn,co
+ integer iretcode
+#ifdef LBFGS
+ character*9 status
+ integer niter,nfun
+ common /lbfgstat/ status,niter,nfun
+#endif
call gen_dist_constr
call sc_conf
call gen_dist_constr
call sc_conf
@@
-359,11
+402,16
@@
c---------------------------------------------------------------------------
if (outpdb) call pdbout(etot,titel(:50),ipdb)
if (outmol2) call mol2out(etot,titel(:32))
if (refstr) call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
if (outpdb) call pdbout(etot,titel(:50),ipdb)
if (outmol2) call mol2out(etot,titel(:32))
if (refstr) call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
+#ifdef LBFGS
+ write (iout,'(a,a9)') 'LBFGS return code:',status
+#else
write (iout,'(a,i3)') 'SUMSL return code:',iretcode
write (iout,'(a,i3)') 'SUMSL return code:',iretcode
+#endif
return
end
c---------------------------------------------------------------------------
subroutine exec_thread
return
end
c---------------------------------------------------------------------------
subroutine exec_thread
+ implicit none
include 'DIMENSIONS'
#ifdef MP
include "mpif.h"
include 'DIMENSIONS'
#ifdef MP
include "mpif.h"
@@
-374,9
+422,10
@@
c---------------------------------------------------------------------------
end
c---------------------------------------------------------------------------
subroutine exec_MC
end
c---------------------------------------------------------------------------
subroutine exec_MC
- implicit real*8 (a-h,o-z)
+ implicit none
include 'DIMENSIONS'
character*10 nodeinfo
include 'DIMENSIONS'
character*10 nodeinfo
+ integer ipar
double precision varia(maxvar)
#ifdef MPI
include "mpif.h"
double precision varia(maxvar)
#ifdef MPI
include "mpif.h"
@@
-405,11
+454,12
@@
c---------------------------------------------------------------------------
end
c---------------------------------------------------------------------------
subroutine exec_mult_eeval_or_minim
end
c---------------------------------------------------------------------------
subroutine exec_mult_eeval_or_minim
- implicit real*8 (a-h,o-z)
+ implicit none
include 'DIMENSIONS'
#ifdef MPI
include 'mpif.h'
include 'DIMENSIONS'
#ifdef MPI
include 'mpif.h'
- dimension muster(mpi_status_size)
+ integer muster(mpi_status_size)
+ integer ierr,ierror
#endif
include 'COMMON.SETUP'
include 'COMMON.TIME1'
#endif
include 'COMMON.SETUP'
include 'COMMON.TIME1'
@@
-418,7
+468,7
@@
c---------------------------------------------------------------------------
include 'COMMON.GEO'
include 'COMMON.HEADER'
include 'COMMON.CONTROL'
include 'COMMON.GEO'
include 'COMMON.HEADER'
include 'COMMON.CONTROL'
- include 'COMMON.CONTACTS'
+c include 'COMMON.CONTACTS'
include 'COMMON.CHAIN'
include 'COMMON.VAR'
include 'COMMON.IOUNITS'
include 'COMMON.CHAIN'
include 'COMMON.VAR'
include 'COMMON.IOUNITS'
@@
-427,8
+477,11
@@
c---------------------------------------------------------------------------
include 'COMMON.MD'
include 'COMMON.SBRIDGE'
double precision varia(maxvar)
include 'COMMON.MD'
include 'COMMON.SBRIDGE'
double precision varia(maxvar)
- dimension ind(6)
- double precision energy(0:max_ene)
+ integer i,j,iconf,ind(6)
+ integer n,it,man,nf_mcmf,nmin,imm,mm,nft
+ double precision energy(0:max_ene),ene,etot,ene0
+ double precision rms,frac,frac_nn,co
+ double precision time
logical eof
eof=.false.
#ifdef MPI
logical eof
eof=.false.
#ifdef MPI
@@
-702,7
+755,7
@@
cjlee end
end
c---------------------------------------------------------------------------
subroutine exec_checkgrad
end
c---------------------------------------------------------------------------
subroutine exec_checkgrad
- implicit real*8 (a-h,o-z)
+ implicit none
include 'DIMENSIONS'
#ifdef MPI
include 'mpif.h'
include 'DIMENSIONS'
#ifdef MPI
include 'mpif.h'
@@
-714,14
+767,16
@@
c---------------------------------------------------------------------------
include 'COMMON.GEO'
include 'COMMON.HEADER'
include 'COMMON.CONTROL'
include 'COMMON.GEO'
include 'COMMON.HEADER'
include 'COMMON.CONTROL'
- include 'COMMON.CONTACTS'
+c include 'COMMON.CONTACTS'
include 'COMMON.CHAIN'
include 'COMMON.VAR'
include 'COMMON.IOUNITS'
include 'COMMON.FFIELD'
include 'COMMON.REMD'
include 'COMMON.MD'
include 'COMMON.CHAIN'
include 'COMMON.VAR'
include 'COMMON.IOUNITS'
include 'COMMON.FFIELD'
include 'COMMON.REMD'
include 'COMMON.MD'
+ include 'COMMON.QRESTR'
include 'COMMON.SBRIDGE'
include 'COMMON.SBRIDGE'
+ integer icall
common /srutu/ icall
double precision energy(0:max_ene)
c print *,"A TU?"
common /srutu/ icall
double precision energy(0:max_ene)
c print *,"A TU?"
@@
-798,7
+853,8
@@
c enddo
goto (10,20,30) icheckgrad
10 call check_ecartint
return
goto (10,20,30) icheckgrad
10 call check_ecartint
return
- 20 call check_cartgrad
+ 20 write (iout,*)
+ & "Checking the gradient of Cartesian coordinates disabled."
return
30 call check_eint
return
return
30 call check_eint
return
@@
-812,6
+868,7
@@
C Energy maps
end
c---------------------------------------------------------------------------
subroutine exec_CSA
end
c---------------------------------------------------------------------------
subroutine exec_CSA
+ implicit none
#ifdef MPI
include "mpif.h"
#endif
#ifdef MPI
include "mpif.h"
#endif
@@
-828,16
+885,18
@@
C This method works only with parallel machines!
end
c---------------------------------------------------------------------------
subroutine exec_softreg
end
c---------------------------------------------------------------------------
subroutine exec_softreg
+ implicit none
include 'DIMENSIONS'
include 'COMMON.IOUNITS'
include 'COMMON.CONTROL'
include 'DIMENSIONS'
include 'COMMON.IOUNITS'
include 'COMMON.CONTROL'
- double precision energy(0:max_ene)
+ double precision energy(0:max_ene),etot
+ double precision rms,frac,frac_nn,co
call chainbuild
call etotal(energy(0))
call enerprint(energy(0))
if (.not.lsecondary) then
write(iout,*) 'Calling secondary structure recognition'
call chainbuild
call etotal(energy(0))
call enerprint(energy(0))
if (.not.lsecondary) then
write(iout,*) 'Calling secondary structure recognition'
- call secondary2(debug)
+ call secondary2(.true.)
else
write(iout,*) 'Using secondary structure supplied in pdb'
endif
else
write(iout,*) 'Using secondary structure supplied in pdb'
endif