Parent

Minimization::Brent

Direct port of Brent algorithm found on GSL. See Unidimensional for methods.

Usage

 min=Minimization::Brent.new(-1000,20000  , proc {|x| (x+1)**2}
 min.expected=1.5  # Expected value
 min.iterate
 min.x_minimum
 min.f_minimum
 min.log

Constants

GSL_SQRT_DBL_EPSILON
(Not documented)

Public Class Methods

new(lower,upper, proc) click to toggle source

(Not documented)

# File lib/minimization.rb, line 215
    def initialize(lower,upper, proc)
      super

      @do_bracketing=true

      # Init

      golden = 0.3819660;      #golden = (3 - sqrt(5))/2

      v = @lower + golden * (@upper - @lower);
      w = v;

      @x_minimum = v ;
      @f_minimum = f(v) ;
      @x_lower=@lower
      @x_upper=@upper
      @f_lower = f(@lower) ;
      @f_upper = f(@lower) ;

      @v = v;
      @w = w;

      @d = 0;
      @e = 0;
      @f_v=f(v)
      @f_w=@f_v
    end

Public Instance Methods

bracketing() click to toggle source

(Not documented)

# File lib/minimization.rb, line 249
    def bracketing
      eval_max=10
      f_left = @f_lower;
      f_right = @f_upper;
      x_left = @x_lower;
      x_right= @x_upper;
      golden = 0.3819660;      # golden = (3 - sqrt(5))/2 */
      nb_eval=0

      if (f_right >= f_left)
        x_center = (x_right - x_left) * golden + x_left;
        nb_eval+=1;
        f_center=f(x_center)
      else
        x_center = x_right ;
        f_center = f_right ;
        x_right = (x_center - x_left).quo(golden) + x_left;
        nb_eval+=1;
        f_right=f(x_right);
      end


      begin
        @log << ["B#{nb_eval}", x_left, x_right, f_left, f_right, (x_left-x_right).abs, (f_left-f_right).abs]
        if (f_center < f_left )
          if (f_center < f_right)
            @x_lower = x_left;
            @x_upper = x_right;
            @x_minimum = x_center;
            @f_lower = f_left;
            @f_upper = f_right;
            @f_minimum = f_center;
            return true;
          elsif (f_center > f_right)
            x_left = x_center;
            f_left = f_center;
            x_center = x_right;
            f_center = f_right;
            x_right = (x_center - x_left).quo(golden) + x_left;
            nb_eval+=1;
            f_right=f(x_right);
          else # f_center == f_right */
            x_right = x_center;
            f_right = f_center;
            x_center = (x_right - x_left).quo(golden) + x_left;
            nb_eval+=1;
            f_center=f(x_center);
          end
        else # f_center >= f_left */
          x_right = x_center;
          f_right = f_center;
          x_center = (x_right - x_left) * golden + x_left;
          nb_eval+=1;
          f_center=f(x_center);
        end
      end while ((nb_eval < eval_max) and
      ((x_right - x_left) > GSL_SQRT_DBL_EPSILON * ( (x_right + x_left) * 0.5 ) + GSL_SQRT_DBL_EPSILON))
      @x_lower = x_left;
      @x_upper = x_right;
      @x_minimum = x_center;
      @f_lower = f_left;
      @f_upper = f_right;
      @f_minimum = f_center;
      return false;

    end
brent_iterate() click to toggle source

Generate one iteration.

# File lib/minimization.rb, line 334
    def brent_iterate
      x_left = @x_lower;
      x_right = @x_upper;

      z = @x_minimum;
      d = @e;
      e = @d;
      v = @v;
      w = @w;
      f_v = @f_v;
      f_w = @f_w;
      f_z = @f_minimum;

      golden = 0.3819660;      # golden = (3 - sqrt(5))/2 */

      w_lower = (z - x_left)
      w_upper = (x_right - z)

      tolerance =  GSL_SQRT_DBL_EPSILON * z.abs

      midpoint = 0.5 * (x_left + x_right)
      _p,q,r=0,0,0
      if (e.abs > tolerance)

        # fit parabola */

        r = (z - w) * (f_z - f_v);
        q = (z - v) * (f_z - f_w);
        _p = (z - v) * q - (z - w) * r;
        q = 2 * (q - r);

        if (q > 0)
          _p = -_p
        else
          q = -q;
        end
        r = e;
        e = d;
      end

      if (_p.abs < (0.5 * q * r).abs and _p < q * w_lower and _p < q * w_upper)
        t2 = 2 * tolerance ;

        d = _p.quo(q);
        u = z + d;

        if ((u - x_left) < t2 or (x_right - u) < t2)
          d = (z < midpoint) ? tolerance : -tolerance ;
        end
      else

        e = (z < midpoint) ? x_right - z : -(z - x_left) ;
        d = golden * e;
      end

      if ( d.abs >= tolerance)
        u = z + d;
      else
        u = z + ((d > 0) ? tolerance : -tolerance) ;
      end

      @e = e;
      @d = d;

      f_u=f(u)

      if (f_u <= f_z)
        if (u < z)
          @x_upper = z;
          @f_upper = f_z;
        else
          @x_lower = z;
          @f_lower = f_z;
        end
        @v = w;
        @f_v = f_w;
        @w = z;
        @f_w = f_z;
        @x_minimum = u;
        @f_minimum = f_u;
        return true;
      else
        if (u < z)
          @x_lower = u;
          @f_lower = f_u;
          return true;
        else
          @x_upper = u;
          @f_upper = f_u;
          return true;
        end

        if (f_u <= f_w or w == z)
          @v = w;
          @f_v = f_w;
          @w = u;
          @f_w = f_u;
          return true;
        elsif f_u <= f_v or v == z or v == w
          @v = u;
          @f_v = f_u;
          return true;
        end

      end
      return false

    end
expected=(v) click to toggle source

(Not documented)

# File lib/minimization.rb, line 243
    def expected=(v)
      @x_minimum=v
      @f_minimum=f(v)
      @do_bracketing=false
    end
iterate() click to toggle source

Start the minimization process If you want to control manually the process, use brent_iterate

# File lib/minimization.rb, line 317
    def iterate
      k=0
      bracketing if @do_bracketing
      while k<@max_iteration and (@x_lower-@x_upper).abs>@epsilon
        k+=1
        result=brent_iterate
        raise FailedIteration,"Error on iteration" if !result
        begin 
          @log << [k, @x_lower, @x_upper, @f_lower, @f_upper, (@x_lower-@x_upper).abs, (@f_lower-@f_upper).abs]
        rescue =>@e
          @log << [k, @e.to_s,nil,nil,nil,nil,nil]
        end
      end
      @iterations=k
      return true
    end

Disabled; run with --debug to generate this.

[Validate]

Generated with the Darkfish Rdoc Generator 1.1.6.