subroutine minimize(etot,x,iretcode,nfun)
- implicit real*8 (a-h,o-z)
+#ifdef LBFGS
+ use minima
+ use inform
+ use output
+ use iounit
+ use scales
+#endif
+ implicit none
include 'DIMENSIONS'
+#ifndef LBFGS
+ integer liv,lv
parameter (liv=60,lv=(77+maxvar*(maxvar+17)/2))
+#endif
*********************************************************************
* OPTIMIZE sets up SUMSL or DFP and provides a simple interface for *
* the calling subprogram. *
include 'COMMON.VAR'
include 'COMMON.GEO'
include 'COMMON.MINIM'
+ integer icall
common /srutu/ icall
- dimension iv(liv)
- double precision minval,x(maxvar),d(maxvar),v(1:lv),xx(maxvar)
- double precision energia(0:n_ene)
+#ifdef LBFGS
+ double precision grdmin
+ external funcgrad
+ external optsave
+#else
+ integer iv(liv)
+ double precision v(1:lv)
+ common /przechowalnia/ v
+ integer idum
+ double precision rdum
+ double precision fdum
external func,gradient,fdum
external func_restr,grad_restr
logical not_done,change,reduce
+#endif
+ double precision x(maxvar),d(maxvar),xx(maxvar)
+ double precision energia(0:n_ene)
+ integer i,nvar_restr,nfun,iretcode
+ double precision etot
c common /przechowalnia/ v
+#ifdef LBFGS
+ maxiter=maxmin
+ coordtype='RIGIDBODY'
+ grdmin=tolf
+ jout=iout
+ jprint=print_min_stat
+ iwrite=0
+ if (.not. allocated(scale)) allocate (scale(nvar))
+c
+c set scaling parameter for function and derivative values;
+c use square root of median eigenvalue of typical Hessian
+c
+ set_scale = .true.
+c nvar = 0
+ do i = 1, nvar
+c if (use(i)) then
+c do j = 1, 3
+c nvar = nvar + 1
+ scale(i) = 12.0d0
+c end do
+c end if
+ end do
+c write (iout,*) "Calling lbfgs"
+ write (iout,*) 'Calling LBFGS, minimization in angles'
+ call var_to_geom(nvar,x)
+ call chainbuild_extconf
+ call etotal(energia(0))
+ call enerprint(energia(0))
+ call lbfgs (nvar,x,etot,grdmin,funcgrad,optsave)
+ deallocate(scale)
+ write (iout,*) "Minimized energy",etot
+#else
icall = 1
NOT_DONE=.TRUE.
do i=nphi+1,nvar
d(i)=1.0D-1
enddo
-cd print *,'Calling SUMSL'
-c call var_to_geom(nvar,x)
-c call chainbuild
-c call etotal(energia(0))
+ write (iout,*) 'Calling SUMSL'
+ call var_to_geom(nvar,x)
+ call chainbuild_extconf
+ call intout
+ call etotal(energia(0))
+ call enerprint(energia(0))
c etot = energia(0)
IF (mask_r) THEN
call x2xx(x,xx,nvar_restr)
c else
c not_done=.false.
c endif
- call chainbuild
+ call chainbuild_extconf
c call etotal(energia(0))
c etot=energia(0)
c call enerprint(energia(0))
c write (*,*) 'Processor',MyID,' leaves MINIMIZE.'
c ENDDO ! NOT_DONE
-
+#endif
return
end
#ifdef MPI
c----------------------------------------------------------------------------
subroutine ergastulum
- implicit real*8 (a-h,o-z)
+ implicit none
include 'DIMENSIONS'
#ifdef MPI
include "mpif.h"
+ double precision time00
+ integer ierr,ierror
#endif
include 'COMMON.SETUP'
include 'COMMON.DERIV'
include 'COMMON.FFIELD'
include 'COMMON.INTERACT'
include 'COMMON.MD'
+#ifdef FIVEDIAG
+ include 'COMMON.LAGRANGE.5diag'
+#else
+ include 'COMMON.LAGRANGE'
+#endif
include 'COMMON.TIME1'
double precision z(maxres6),d_a_tmp(maxres6)
double precision edum(0:n_ene),time_order(0:10)
double precision Gcopy(maxres2,maxres2)
common /przechowalnia/ Gcopy
integer icall /0/
+ integer i,j,iorder
C Workers wait for variables and NF, and NFL from the boss
iorder=0
do while (iorder.ge.0)
call sum_gradient
c write (2,*) "After sum_gradient"
c write (2,*) "dimen",dimen," dimen3",dimen3
-c call flush(2)
+c call flush(2
+#ifndef FIVEDIAG
else if (iorder.eq.4) then
call ginv_mult(z,d_a_tmp)
else if (iorder.eq.5) then
c call flush(2)
c write (iout,*) "My chunk of ginv_block"
c call MATOUT2(my_ng_count,dimen3,maxres2,maxers2,ginv_block)
+#endif
else if (iorder.eq.6) then
call int_from_cart1(.false.)
else if (iorder.eq.7) then
call chainbuild_cart
else if (iorder.eq.8) then
call intcartderiv
+#ifndef FIVEDIAG
else if (iorder.eq.9) then
call fricmat_mult(z,d_a_tmp)
+#endif
else if (iorder.eq.10) then
call setup_fricmat
endif
end
#endif
************************************************************************
+#ifdef LBFGS
+ double precision function funcgrad(x,g)
+ implicit none
+ include 'DIMENSIONS'
+ include 'COMMON.CONTROL'
+ include 'COMMON.CHAIN'
+ include 'COMMON.DERIV'
+ include 'COMMON.VAR'
+ include 'COMMON.INTERACT'
+ include 'COMMON.FFIELD'
+ include 'COMMON.MD'
+ include 'COMMON.QRESTR'
+ include 'COMMON.IOUNITS'
+ include 'COMMON.GEO'
+ double precision energia(0:n_ene)
+ double precision x(nvar),g(nvar)
+ integer i
+c if (jjj.gt.0) then
+c write (iout,*) "in func x"
+c write (iout,'(10f8.3)') (rad2deg*x(i),i=1,n)
+c endif
+ call var_to_geom(nvar,x)
+ call zerograd
+ call chainbuild_extconf
+ call etotal(energia(0))
+ call sum_gradient
+ funcgrad=energia(0)
+ call cart2intgrad(nvar,g)
+C
+C Add the components corresponding to local energy terms.
+C
+c Add the usampl contributions
+ if (usampl) then
+ do i=1,nres-3
+ gloc(i,icg)=gloc(i,icg)+dugamma(i)
+ enddo
+ do i=1,nres-2
+ gloc(nphi+i,icg)=gloc(nphi+i,icg)+dutheta(i)
+ enddo
+ endif
+ do i=1,nvar
+cd write (iout,*) 'i=',i,'g=',g(i),' gloc=',gloc(i,icg)
+ g(i)=g(i)+gloc(i,icg)
+ enddo
+ return
+ end
+#else
subroutine func(n,x,nf,f,uiparm,urparm,ufparm)
implicit real*8 (a-h,o-z)
include 'DIMENSIONS'
return
end
c-------------------------------------------------------
+#endif
subroutine x2xx(x,xx,n)
implicit real*8 (a-h,o-z)
include 'DIMENSIONS'
include 'COMMON.VAR'
include 'COMMON.CHAIN'
include 'COMMON.INTERACT'
+ include 'COMMON.IOUNITS'
double precision xx(maxvar),x(maxvar)
+c write (iout,*) "nvar",nvar
do i=1,nvar
varall(i)=x(i)
enddo
- ig=0
- igall=0
- do i=4,nres
- igall=igall+1
- if (mask_phi(i).eq.1) then
- ig=ig+1
+ ig=0
+ igall=0
+ do i=4,nres
+ igall=igall+1
+ if (mask_phi(i).eq.1) then
+ ig=ig+1
xx(ig)=x(igall)
- endif
- enddo
-
- do i=3,nres
- igall=igall+1
- if (mask_theta(i).eq.1) then
- ig=ig+1
+ endif
+ enddo
+
+ do i=3,nres
+ igall=igall+1
+ if (mask_theta(i).eq.1) then
+ ig=ig+1
xx(ig)=x(igall)
- endif
+ endif
enddo
- do ij=1,2
- do i=2,nres-1
- if (itype(i).ne.10) then
- igall=igall+1
- if (mask_side(i).eq.1) then
- ig=ig+1
+ do ij=1,2
+ do i=2,nres-1
+ if (itype(i).ne.10 .and. itype(i).ne.ntyp1) then
+ igall=igall+1
+ if (mask_side(i).eq.1) then
+ ig=ig+1
xx(ig)=x(igall)
- endif
- endif
- enddo
+c write (iout,*) "ij",ij," i",i," ig",ig," igall",igall
+c write (iout,*) "x",x(igall)," xx",xx(ig)
+ endif
+ endif
+ enddo
enddo
n=ig
include 'COMMON.VAR'
include 'COMMON.CHAIN'
include 'COMMON.INTERACT'
+ include 'COMMON.IOUNITS'
double precision xx(maxvar),x(maxvar)
do i=1,nvar
x(i)=varall(i)
enddo
- ig=0
- igall=0
- do i=4,nres
- igall=igall+1
- if (mask_phi(i).eq.1) then
- ig=ig+1
+ ig=0
+ igall=0
+ do i=4,nres
+ igall=igall+1
+ if (mask_phi(i).eq.1) then
+ ig=ig+1
x(igall)=xx(ig)
- endif
- enddo
-
- do i=3,nres
- igall=igall+1
- if (mask_theta(i).eq.1) then
- ig=ig+1
+ endif
+ enddo
+
+ do i=3,nres
+ igall=igall+1
+ if (mask_theta(i).eq.1) then
+ ig=ig+1
x(igall)=xx(ig)
- endif
+ endif
enddo
- do ij=1,2
- do i=2,nres-1
- if (itype(i).ne.10) then
- igall=igall+1
- if (mask_side(i).eq.1) then
- ig=ig+1
+ do ij=1,2
+ do i=2,nres-1
+ if (itype(i).ne.10 .and. itype(i).ne.ntyp1) then
+ igall=igall+1
+ if (mask_side(i).eq.1) then
+ ig=ig+1
x(igall)=xx(ig)
- endif
- endif
- enddo
+c write (iout,*) "ij",ij," i",i," ig",ig," igall",igall
+c write (iout,*) "x",x(igall)," xx",xx(ig)
+ endif
+ endif
+ enddo
enddo
return
c----------------------------------------------------------
subroutine minim_dc(etot,iretcode,nfun)
+#ifdef LBFGS
+ use minima
+ use inform
+ use output
+ use iounit
+ use scales
+#endif
implicit real*8 (a-h,o-z)
include 'DIMENSIONS'
+#ifndef LBFGS
parameter (liv=60,lv=(77+maxvar*(maxvar+17)/2))
+#endif
#ifdef MPI
include 'mpif.h'
#endif
include 'COMMON.GEO'
include 'COMMON.MINIM'
include 'COMMON.CHAIN'
+ double precision minval,x(maxvar),d(maxvar),xx(maxvar)
+#ifdef LBFGS
+ double precision grdmin
+ double precision funcgrad_dc
+ external funcgrad_dc,optsave
+#else
dimension iv(liv)
- double precision minval,x(maxvar),d(maxvar),v(1:lv),xx(maxvar)
+ double precision v(1:lv)
common /przechowalnia/ v
-
- double precision energia(0:n_ene)
external func_dc,grad_dc,fdum
- logical not_done,change,reduce
+#endif
double precision g(maxvar),f1
-
+ integer nvarx
+ double precision energia(0:n_ene)
+#ifdef LBFGS
+ maxiter=maxmin
+ coordtype='CARTESIAN'
+ grdmin=tolf
+ jout=iout
+ jprint=print_min_stat
+ iwrite=0
+#else
call deflt(2,iv,liv,lv,v)
* 12 means fresh start, dont call deflt
iv(1)=12
do i=1,6*nres
d(i)=1.0D-1
enddo
-
+#endif
k=0
do i=1,nres-1
do j=1,3
enddo
endif
enddo
+ nvarx=k
+ write (iout,*) "Variables set up nvarx",nvarx
+ write (iout,*) "Before energy minimization"
+ call etotal(energia(0))
+ call enerprint(energia(0))
+#ifdef LBFGS
+c
+c From tinker
+c
+c perform dynamic allocation of some global arrays
+c
+ if (.not. allocated(scale)) allocate (scale(nvarx))
+c
+c set scaling parameter for function and derivative values;
+c use square root of median eigenvalue of typical Hessian
+c
+ set_scale = .true.
+c nvar = 0
+ do i = 1, nvarx
+c if (use(i)) then
+c do j = 1, 3
+c nvar = nvar + 1
+ scale(i) = 12.0d0
+c end do
+c end if
+ end do
+c write (iout,*) "minim_dc Calling lbfgs"
+ call lbfgs (nvarx,x,etot,grdmin,funcgrad_dc,optsave)
+ deallocate(scale)
+c write (iout,*) "minim_dc After lbfgs"
+#else
c-----
c write (iout,*) "checkgrad before SUMSL"
c icheckgrad=1
c write (iout,*) "checkgrad after SUMSL"
c call exec_checkgrad
c-----
-
+#endif
k=0
do i=1,nres-1
do j=1,3
cd x(i)=x(i)-1.0D-5
cd print '(i5,2f15.5)',i,g(i),(f1-f)/1.0D-5
cd enddo
-
+#ifndef LBFGS
etot=v(10)
iretcode=iv(1)
nfun=iv(6)
+#endif
return
end
+#ifdef LBFGS
+ double precision function funcgrad_dc(x,g)
+ implicit real*8 (a-h,o-z)
+ include 'DIMENSIONS'
+#ifdef MPI
+ include 'mpif.h'
+#endif
+ include 'COMMON.SETUP'
+ include 'COMMON.CHAIN'
+ include 'COMMON.DERIV'
+ include 'COMMON.VAR'
+ include 'COMMON.INTERACT'
+ include 'COMMON.FFIELD'
+ include 'COMMON.MD'
+ include 'COMMON.IOUNITS'
+ integer k
+ dimension x(maxvar),g(maxvar)
+ double precision energia(0:n_ene)
+ common /gacia/ nf
+c
+ nf=nf+1
+ k=0
+ do i=1,nres-1
+ do j=1,3
+ k=k+1
+ dc(j,i)=x(k)
+ enddo
+ enddo
+ do i=2,nres-1
+ if (ialph(i,1).gt.0) then
+ do j=1,3
+ k=k+1
+ dc(j,i+nres)=x(k)
+ enddo
+ endif
+ enddo
+ call chainbuild_cart
+ call zerograd
+ call etotal(energia(0))
+c write (iout,*) "energia",energia(0)
+ funcgrad_dc=energia(0)
+C
+C Evaluate the derivatives of virtual bond lengths and SC vectors in variables.
+C
+ call cartgrad
+ k=0
+ do i=1,nres-1
+ do j=1,3
+ k=k+1
+ g(k)=gcart(j,i)
+ enddo
+ enddo
+ do i=2,nres-1
+ if (ialph(i,1).gt.0) then
+ do j=1,3
+ k=k+1
+ g(k)=gxcart(j,i)
+ enddo
+ endif
+ enddo
+ return
+ end
+#else
subroutine func_dc(n,x,nf,f,uiparm,urparm,ufparm)
implicit real*8 (a-h,o-z)
include 'DIMENSIONS'
double precision ufparm
external ufparm
integer uiparm(1)
- real*8 urparm(1)
+ real*8 urparm(1)
dimension x(maxvar)
nfl=nf
cbad icg=mod(nf,2)+1
return
end
+#endif