1   
  2   
  3   
  4   
  5   
  6   
  7   
  8   
  9   
 10   
 11   
 12   
 13   
 14   
 15   
 16   
 17   
 18   
 19   
 20   
 21   
 22   
 23   
 24  """A line search algorithm based on interpolation. 
 25   
 26  This file is part of the U{minfx optimisation library<https://sourceforge.net/projects/minfx>}. 
 27  """ 
 28   
 29   
 30  from copy import deepcopy 
 31  from numpy import dot, sqrt 
 32   
 33   
 34  from minfx.line_search.interpolate import cubic_ext, quadratic_fafbga 
 35   
 36   
 37  quadratic = quadratic_fafbga 
 38   
 39   
 41      """A line search algorithm based on interpolation. 
 42   
 43      Pages 56-57, from 'Numerical Optimization' by Jorge Nocedal and Stephen J. Wright, 1999, 2nd ed. 
 44   
 45      Requires the gradient function. 
 46   
 47   
 48      @param func:            The function to minimise. 
 49      @type func:             func 
 50      @param args:            The tuple of arguments to supply to the functions func. 
 51      @type args:             tuple 
 52      @param x:               The parameter vector at minimisation step k. 
 53      @type x:                numpy array 
 54      @param f:               The function value at the point x. 
 55      @type f:                float 
 56      @param g:               The function gradient vector at the point x. 
 57      @type g:                numpy array 
 58      @param p:               The descent direction. 
 59      @type p:                numpy array 
 60      @keyword a_init:        Initial step length. 
 61      @type a_init:           flaot 
 62      @keyword mu:            Constant determining the slope for the sufficient decrease condition (0 < mu < 1). 
 63      @type mu:               float 
 64      @keyword print_flag:    The higher the value, the greater the amount of info printed out. 
 65      @type print_flag:       int 
 66      """ 
 67   
 68       
 69      i = 1 
 70      f_count = 0 
 71      a0 = {} 
 72      a0['a'] = 0.0 
 73      a0['phi'] = f 
 74      a0['phi_prime'] = dot(g, p) 
 75   
 76       
 77      a = {} 
 78      a['a'] = a_init 
 79      a['phi'] = func(*(x + a['a']*p,)+args) 
 80      f_count = f_count + 1 
 81   
 82      if print_flag: 
 83          print("\n<Line search initial values>") 
 84          print_data("Pre (a0)", i-1, a0) 
 85   
 86       
 87      if a0['phi_prime'] >= 0.0: 
 88          raise NameError("The gradient at point 0 of this line search is positive, ie p is not a descent direction and the line search will not work.") 
 89   
 90       
 91      if a['phi'] <= a0['phi'] + mu * a['a'] * a0['phi_prime']: 
 92          return a['a'], f_count 
 93   
 94       
 95      a_last = deepcopy(a) 
 96   
 97       
 98      a_new = - 0.5 * a0['phi_prime'] * a['a']**2 / (a['phi'] - a0['phi'] - a0['phi_prime']*a['a']) 
 99      a['a'] = a_new 
100      a['phi'] = func(*(x + a['a']*p,)+args) 
101      f_count = f_count + 1 
102   
103       
104      if a['phi'] <= a0['phi'] + mu * a['a'] * a0['phi_prime']: 
105          return a['a'], f_count 
106   
107      while True: 
108          if print_flag: 
109              print("<Line search iteration i = " + repr(i) + " >") 
110              print_data("Initial (a)", i, a) 
111              print_data("Initial (a_last)", i, a_last) 
112   
113          beta1 = 1.0 / (a_last['a']**2 * a['a']**2 * (a['a'] - a_last['a'])) 
114          beta2 = a['phi'] - a0['phi'] - a0['phi_prime'] * a['a'] 
115          beta3 = a_last['phi'] - a0['phi'] - a0['phi_prime'] * a_last['a'] 
116   
117          fact_a = beta1 * (a_last['a']**2 * beta2 - a['a']**2 * beta3) 
118          fact_b = beta1 * (-a_last['a']**3 * beta2 + a['a']**3 * beta3) 
119   
120          a_new = {} 
121          a_new['a'] = (-fact_b + sqrt(fact_b**2 - 3.0 * fact_a * a0['phi_prime'])) / (3.0 * fact_a) 
122          a_new['phi'] = func(*(x + a_new['a']*p,)+args) 
123          f_count = f_count + 1 
124   
125           
126          if a_new['phi'] <= a0['phi'] + mu * a_new['a'] * a0['phi_prime']: 
127              return a_new['a'], f_count 
128   
129           
130          if a['a'] - a_new['a'] > 0.5 * a['a'] or 1.0 - a_new['a']/a['a'] < 0.9: 
131              a_new['a'] = 0.5 * a['a'] 
132              a_new['phi'] = func(*(x + a_new['a']*p,)+args) 
133              f_count = f_count + 1 
134   
135           
136          a_last = deepcopy(a) 
137          a = deepcopy(a_new) 
 138   
139   
141      """Temp func for debugging.""" 
142   
143      print(text + " data printout:") 
144      print("   Iteration:      " + repr(k)) 
145      print("   a:              " + repr(a['a'])) 
146      print("   phi:            " + repr(a['phi'])) 
147      print("   phi_prime:      " + repr(a['phi_prime'])) 
 148