/[escript]/trunk/downunder/py_src/minimizers.py
ViewVC logotype

Diff of /trunk/downunder/py_src/minimizers.py

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 4078 by gross, Thu Nov 15 03:45:24 2012 UTC revision 4079 by gross, Fri Nov 16 07:59:01 2012 UTC
# Line 52  def _zoom(phi, gradphi, phiargs, alpha_l Line 52  def _zoom(phi, gradphi, phiargs, alpha_l
52          alpha=alpha_lo+.5*(alpha_hi-alpha_lo) # should use interpolation...          alpha=alpha_lo+.5*(alpha_hi-alpha_lo) # should use interpolation...
53          args_a=phiargs(alpha)          args_a=phiargs(alpha)
54          phi_a=phi(alpha, *args_a)          phi_a=phi(alpha, *args_a)
55          zoomlogger.debug("iteration %d, alpha=%e, phi(alpha)=%e"%(i,alpha,phi_a))          zoomlogger.debug("Zoom.iteration %d, alpha=%e, phi(alpha)=%e"%(i,alpha,phi_a))
56          if phi_a > phi0+c1*alpha*gphi0 or phi_a >= phi_lo:          if phi_a > phi0+c1*alpha*gphi0 or phi_a >= phi_lo:
57              alpha_hi=alpha              alpha_hi=alpha
58          else:          else:
59              gphi_a=gradphi(alpha, *args_a)              gphi_a=gradphi(alpha, *args_a)
60              zoomlogger.debug("grad(phi(alpha))=%e"%(gphi_a))              zoomlogger.debug("Zoom.grad(phi(alpha))=%e"%(gphi_a))
61              if np.abs(gphi_a) <= -c2*gphi0:              if np.abs(gphi_a) <= -c2*gphi0:
62                  break                  break
63              if gphi_a*(alpha_hi-alpha_lo) >= 0:              if gphi_a*(alpha_hi-alpha_lo) >= 0:
# Line 120  def line_search(f, x, p, gf, fx, alpha_m Line 120  def line_search(f, x, p, gf, fx, alpha_m
120    
121          args_a=phiargs(alpha)          args_a=phiargs(alpha)
122          phi_a=phi(alpha, *args_a)          phi_a=phi(alpha, *args_a)
123          lslogger.debug("iteration %d, alpha=%e, phi(alpha)=%e"%(i,alpha,phi_a))          lslogger.debug("Line Search.iteration %d, alpha=%e, phi(alpha)=%e"%(i,alpha,phi_a))
124          if (phi_a > phi0+c1*alpha*gphi0) or ((phi_a>=old_phi_a) and (i>1)):          if (phi_a > phi0+c1*alpha*gphi0) or ((phi_a>=old_phi_a) and (i>1)):
125              alpha, phi_a, gphi_a = _zoom(phi, gradphi, phiargs, old_alpha, alpha, old_phi_a, phi_a, c1, c2, phi0, gphi0)              alpha, phi_a, gphi_a = _zoom(phi, gradphi, phiargs, old_alpha, alpha, old_phi_a, phi_a, c1, c2, phi0, gphi0)
126              break              break
# Line 271  class MinimizerLBFGS(AbstractMinimizer): Line 271  class MinimizerLBFGS(AbstractMinimizer):
271    
272          while error > self._tol and k < self._imax:          while error > self._tol and k < self._imax:
273              #self.logger.info("\033[1;31miteration %d\033[1;30m, error=%e"%(k,error))              #self.logger.info("\033[1;31miteration %d\033[1;30m, error=%e"%(k,error))
274              self.logger.info("iteration %d, error=%e"%(k,error))              if k >0:
275                    self.logger.info("LBFGS.iteration %d, error=%e"%(k,error))
276                else:
277                    self.logger.info("LBFGS.iteration %d"%k)
278              # determine search direction              # determine search direction
279              p = -self._twoLoop(invH_scale, gf, s_and_y, x, *args)              p = -self._twoLoop(invH_scale, gf, s_and_y, x, *args)
280    
# Line 282  class MinimizerLBFGS(AbstractMinimizer): Line 285  class MinimizerLBFGS(AbstractMinimizer):
285    
286              # determine step length              # determine step length
287              alpha, fx_new, gf_new = line_search(self._f, x, p, gf, fx)              alpha, fx_new, gf_new = line_search(self._f, x, p, gf, fx)
288              self.logger.debug("alpha=%e"%(alpha))              self.logger.debug("LBFGS.alpha=%e"%(alpha))
289              # execute the step              # execute the step
290              x_new = x + alpha*p              x_new = x + alpha*p
291                print gf_new
292              if gf_new is None:              if gf_new is None:
293              args=self._f.getArguments(x_new)              args=self._f.getArguments(x_new)
294                  gf_new=self._f.getGradient(x_new, args)                  gf_new=self._f.getGradient(x_new, args)
295    
296              delta_x=x_new-x              delta_x=x_new-x
297              delta_g=gf_new-gf              delta_g=gf_new-gf
298              s_and_y.append((delta_x,delta_g, self._f.getDualProduct(delta_x, delta_g) ))              rho=self._f.getDualProduct(delta_x, delta_g)
299                            print "rho =",rho
300                if abs(rho)>0 :
301                   s_and_y.append((delta_x,delta_g, rho ))
302                else:
303               raise ZeroDivisionError("LBFGS break down.")
304              # this needs to be reviewed              # this needs to be reviewed
305              if fx_new==0.:              if fx_new==0.:
306                  error=fx                  error=fx
# Line 317  class MinimizerLBFGS(AbstractMinimizer): Line 325  class MinimizerLBFGS(AbstractMinimizer):
325                invH_scale=self._f.getDualProduct(delta_x,delta_g)/denom                invH_scale=self._f.getDualProduct(delta_x,delta_g)/denom
326            else:            else:
327                invH_scale=self._initial_H                invH_scale=self._initial_H
328                self.logger.debug("Break down in H update. Resetting to initial value %s."%self._initial_H)                self.logger.debug("LBFGS.Break down in H update. Resetting to initial value %s."%self._initial_H)
329    
330          if k >= self._imax:          if k >= self._imax:
331              reason=self.MAX_ITERATIONS_REACHED              reason=self.MAX_ITERATIONS_REACHED
332              self.logger.warning("Maximum number of iterations reached!")              self.logger.warning("LBFGS.Maximum number of iterations reached!")
333          else:          else:
334              reason=self.TOLERANCE_REACHED              reason=self.TOLERANCE_REACHED
335              self.logger.warning("Success after %d iterations! Final error=%e"%(k,error))              self.logger.warning("LBFGS.Success after %d iterations! Final error=%e"%(k,error))
336          self._result=x          self._result=x
337          return reason          return reason
338    
# Line 340  class MinimizerLBFGS(AbstractMinimizer): Line 348  class MinimizerLBFGS(AbstractMinimizer):
348              alpha.append(a)              alpha.append(a)
349              q=q-a*y              q=q-a*y
350    
351          if self._f.provides_inverse_Hessian_approximation:              if self._f.provides_inverse_Hessian_approximation:  
352               r= self._f.getInverseHessianApproximation(x, q, *args)               r= self._f.getInverseHessianApproximation(x, q, *args)
353          else:          else:
354           r= invH_scale * q           r= invH_scale * q
Legend:
Removed from v.4078  
changed lines
  Added in v.4079
  ViewVC Help
Powered by ViewVC 1.1.26