3 c ###################################################
4 c ## COPYRIGHT (C) 1990 by Jay William Ponder ##
5 c ## All Rights Reserved ##
6 c ###################################################
8 c ##################################################################
10 c ## program minimize -- low storage BFGS Cartesian optimizer ##
12 c ##################################################################
15 c "minimize" performs energy minimization in Cartesian coordinate
16 c space using a low storage BFGS nonlinear optimization
32 real*8 minimum,minimiz1
33 real*8 grdmin,gnorm,grms
35 real*8, allocatable :: xx(:)
36 real*8, allocatable :: derivs(:,:)
37 logical exist,analytic
45 c set up the structure and mechanics calculation
51 c perform the setup functions needed for optimization
55 c use either analytical or numerical gradients
60 c get termination criterion as RMS gradient per atom
63 call nextarg (string,exist)
64 if (exist) read (string,*,err=10,end=10) grdmin
66 if (grdmin .le. 0.0d0) then
68 20 format (/,' Enter RMS Gradient per Atom Criterion',
70 read (input,30) grdmin
73 if (grdmin .le. 0.0d0) grdmin = 0.01d0
75 c write out a copy of coordinates for later update
78 minfile = filename(1:leng)//'.xyz'
79 call version (minfile,'new')
80 open (unit=imin,file=minfile,status='new')
85 c perform dynamic allocation of some global arrays
87 if (.not. allocated(scale)) allocate (scale(3*n))
89 c set scaling parameter for function and derivative values;
90 c use square root of median eigenvalue of typical Hessian
103 c perform dynamic allocation of some local arrays
106 allocate (derivs(3,n))
108 c convert atomic coordinates to optimization parameters
114 xx(nvar) = x(i) * scale(nvar)
116 xx(nvar) = y(i) * scale(nvar)
118 xx(nvar) = z(i) * scale(nvar)
122 c make the call to the optimization routine
124 call lbfgs (nvar,xx,minimum,grdmin,minimiz1,optsave)
126 c convert optimization parameters to atomic coordinates
132 x(i) = xx(nvar) / scale(nvar)
134 y(i) = xx(nvar) / scale(nvar)
136 z(i) = xx(nvar) / scale(nvar)
140 c compute the final function and RMS gradient values
143 call gradient (minimum,derivs)
146 call numgrad (energy,derivs,eps)
148 if (use_rattle) call shakef (derivs)
153 gnorm = gnorm + derivs(j,i)**2
158 grms = gnorm / sqrt(dble(nvar/3))
160 c perform deallocation of some local arrays
165 c write out the final function and gradient values
167 if (digits .ge. 8) then
168 if (grms .gt. 1.0d-8) then
169 write (iout,40) minimum,grms,gnorm
170 40 format (/,' Final Function Value :',2x,f20.8,
171 & /,' Final RMS Gradient :',4x,f20.8,
172 & /,' Final Gradient Norm :',3x,f20.8)
174 write (iout,50) minimum,grms,gnorm
175 50 format (/,' Final Function Value :',2x,f20.8,
176 & /,' Final RMS Gradient :',4x,d20.8,
177 & /,' Final Gradient Norm :',3x,d20.8)
179 else if (digits .ge. 6) then
180 if (grms .gt. 1.0d-6) then
181 write (iout,60) minimum,grms,gnorm
182 60 format (/,' Final Function Value :',2x,f18.6,
183 & /,' Final RMS Gradient :',4x,f18.6,
184 & /,' Final Gradient Norm :',3x,f18.6)
186 write (iout,70) minimum,grms,gnorm
187 70 format (/,' Final Function Value :',2x,f18.6,
188 & /,' Final RMS Gradient :',4x,d18.6,
189 & /,' Final Gradient Norm :',3x,d18.6)
192 if (grms .gt. 1.0d-4) then
193 write (iout,80) minimum,grms,gnorm
194 80 format (/,' Final Function Value :',2x,f16.4,
195 & /,' Final RMS Gradient :',4x,f16.4,
196 & /,' Final Gradient Norm :',3x,f16.4)
198 write (iout,90) minimum,grms,gnorm
199 90 format (/,' Final Function Value :',2x,f16.4,
200 & /,' Final RMS Gradient :',4x,d16.4,
201 & /,' Final Gradient Norm :',3x,d16.4)
205 c write the final coordinates into a file
207 if (use_bounds) call bounds
209 open (unit=imin,file=minfile,status='old')
214 c perform any final tasks before program exit
220 c ###############################################################
222 c ## function minimiz1 -- energy and gradient for minimize ##
224 c ###############################################################
227 c "minimiz1" is a service routine that computes the energy and
228 c gradient for a low storage BFGS optimization in Cartesian
232 function minimiz1 (xx,g)
243 real*8, allocatable :: derivs(:,:)
248 c use either analytical or numerical gradients
253 c convert optimization parameters to atomic coordinates
259 x(i) = xx(nvar) / scale(nvar)
261 y(i) = xx(nvar) / scale(nvar)
263 z(i) = xx(nvar) / scale(nvar)
267 c adjust atomic coordinates to satisfy distance constraints
269 if (use_rattle) call shake (x,y,z)
271 c perform dynamic allocation of some local arrays
273 allocate (derivs(3,n))
275 c compute and store the energy and gradient
278 call gradient (e,derivs)
281 call numgrad (energy,derivs,eps)
285 c adjust gradient to remove components along constraints
287 if (use_rattle) call shakef (derivs)
289 c convert coordinates and gradient to optimization parameters
295 xx(nvar) = x(i) * scale(nvar)
296 g(nvar) = derivs(1,i) / scale(nvar)
298 xx(nvar) = y(i) * scale(nvar)
299 g(nvar) = derivs(2,i) / scale(nvar)
301 xx(nvar) = z(i) * scale(nvar)
302 g(nvar) = derivs(3,i) / scale(nvar)
306 c perform deallocation of some local arrays