diff --git a/.travis.yml b/.travis.yml index d250ed9..a582fbd 100644 --- a/.travis.yml +++ b/.travis.yml @@ -1,5 +1,10 @@ language: ruby rvm: - - 1.9.3 - - 2.0.0 - - jruby + - 2.2 + - 2.3 + - 2.4 + - ruby + - ruby-head + - jruby-9.1.9 +install: gem install minitest +script: rake diff --git a/Gemfile b/Gemfile deleted file mode 100644 index 75db47d..0000000 --- a/Gemfile +++ /dev/null @@ -1,6 +0,0 @@ -source "https://rubygems.org" -gem "rake" - -group :test do - gem "test-unit" -end diff --git a/Gemfile.lock b/Gemfile.lock deleted file mode 100644 index 9b744b9..0000000 --- a/Gemfile.lock +++ /dev/null @@ -1,12 +0,0 @@ -GEM - remote: https://rubygems.org/ - specs: - rake (10.1.0) - test-unit (2.0.0.0) - -PLATFORMS - ruby - -DEPENDENCIES - rake - test-unit diff --git a/Rakefile b/Rakefile index cb70dcb..9014442 100644 --- a/Rakefile +++ b/Rakefile @@ -1,7 +1,19 @@ require 'rake/testtask' -task :default => :test - Rake::TestTask.new do |t| - t.test_files = FileList['test/**/*_test.rb'] + t.pattern = ['test/*.rb'] + t.warning = true +end + +Rake::TestTask.new(bench: :loadavg) do |t| + t.pattern = ['test/bench/*.rb'] + t.warning = true + t.description = "Run benchmarks" end + +desc "Show current system load" +task :loadavg do + puts "/proc/loadavg %s" % (File.read("/proc/loadavg") rescue "Unavailable") +end + +task :default => :test diff --git a/lib/simplex.rb b/lib/simplex.rb index 98c893f..573b678 100644 --- a/lib/simplex.rb +++ b/lib/simplex.rb @@ -1,192 +1,164 @@ -require 'matrix' - -class Vector - public :[]= -end - class Simplex DEFAULT_MAX_PIVOTS = 10_000 - class UnboundedProblem < StandardError - end + class Error < RuntimeError; end + class UnboundedProblem < Error; end + class SanityCheck < Error; end + class TooManyPivots < Error; end attr_accessor :max_pivots + # c - coefficients of objective function; size: num_vars + # a - inequality lhs coefficients; 2dim size: num_inequalities, num_vars + # b - inequality rhs constants size: num_inequalities def initialize(c, a, b) - @pivot_count = 0 + num_vars = c.size + num_inequalities = b.size + raise(ArgumentError, "a doesn't match b") unless a.size == num_inequalities + raise(ArgumentError, "a doesn't match c") unless a.first.size == num_vars + @max_pivots = DEFAULT_MAX_PIVOTS - # Problem dimensions - @num_non_slack_vars = a.first.length - @num_constraints = b.length + # Problem dimensions; these never change + @num_non_slack_vars = num_vars + @num_constraints = num_inequalities @num_vars = @num_non_slack_vars + @num_constraints # Set up initial matrix A and vectors b, c - @c = Vector[*c.map {|c1| -1*c1 } + [0]*@num_constraints] - @a = a.map {|a1| Vector[*(a1.clone + [0]*@num_constraints)]} - @b = Vector[*b.clone] - - unless @a.all? {|a| a.size == @c.size } and @b.size == @a.length - raise ArgumentError, "Input arrays have mismatched dimensions" - end - - 0.upto(@num_constraints - 1) {|i| @a[i][@num_non_slack_vars + i] = 1 } + @c = c.map { |flt| -1 * flt } + Array.new(@num_constraints, 0) + @a = a.map.with_index { |ary, i| + if ary.size != @num_non_slack_vars + raise ArgumentError, "a is inconsistent" + end + ary + Array.new(@num_constraints) { |ci| ci == i ? 1 : 0 } + } + @b = b # set initial solution: all non-slack variables = 0 - @x = Vector[*([0]*@num_vars)] @basic_vars = (@num_non_slack_vars...@num_vars).to_a - update_solution - end - - def solution - solve - current_solution - end - - def current_solution - @x.to_a[0...@num_non_slack_vars] + self.update_solution end + # does not modify vector / matrix def update_solution - 0.upto(@num_vars - 1) {|i| @x[i] = 0 } + @x = Array.new(@num_vars, 0) + + @basic_vars.each { |basic_var| + idx = nil + @num_constraints.times { |i| + if @a[i][basic_var] == 1 + idx =i + break + end + } + raise(SanityCheck, "no idx for basic_var #{basic_var} in a") unless idx + @x[basic_var] = @b[idx] + } + end - @basic_vars.each do |basic_var| - row_with_1 = row_indices.detect do |row_ix| - @a[row_ix][basic_var] == 1 - end - @x[basic_var] = @b[row_with_1] - end + def solution + self.solve + self.current_solution end def solve - while can_improve? - @pivot_count += 1 - raise "Too many pivots" if @pivot_count > max_pivots - pivot + count = 0 + while self.can_improve? + count += 1 + raise(TooManyPivots, count.to_s) unless count < @max_pivots + self.pivot end end - def can_improve? - !!entering_variable + def current_solution + @x[0...@num_non_slack_vars] end - def variables - (0...@c.size).to_a + def can_improve? + !self.entering_variable.nil? end + # idx of @c's minimum negative value + # nil when no improvement is possible + # def entering_variable - variables.select { |var| @c[var] < 0 }. - min_by { |var| @c[var] } + (0...@c.size).select { |i| @c[i] < 0 }.min_by { |i| @c[i] } end def pivot - pivot_column = entering_variable - pivot_row = pivot_row(pivot_column) - raise UnboundedProblem unless pivot_row - leaving_var = basic_variable_in_row(pivot_row) - replace_basic_variable(leaving_var => pivot_column) + pivot_column = self.entering_variable or return nil + pivot_row = self.pivot_row(pivot_column) or raise UnboundedProblem + leaving_var = nil + @a[pivot_row].each_with_index { |a, i| + if a == 1 and @basic_vars.include?(i) + leaving_var = i + break + end + } + raise(SanityCheck, "no leaving_var") if leaving_var.nil? + + @basic_vars.delete(leaving_var) + @basic_vars.push(pivot_column) + @basic_vars.sort! pivot_ratio = Rational(1, @a[pivot_row][pivot_column]) # update pivot row - @a[pivot_row] *= pivot_ratio - @b[pivot_row] = pivot_ratio * @b[pivot_row] + @a[pivot_row] = @a[pivot_row].map { |val| val * pivot_ratio } + @b[pivot_row] = @b[pivot_row] * pivot_ratio # update objective - @c -= @c[pivot_column] * @a[pivot_row] + # @c -= @c[pivot_column] * @a[pivot_row] + @c = @c.map.with_index { |val, i| + val - @c[pivot_column] * @a[pivot_row][i] + } # update A and B - (row_indices - [pivot_row]).each do |row_ix| - r = @a[row_ix][pivot_column] - @a[row_ix] -= r * @a[pivot_row] - @b[row_ix] -= r * @b[pivot_row] - end + @num_constraints.times { |i| + next if i == pivot_row + r = @a[i][pivot_column] + @a[i] = @a[i].map.with_index { |val, j| val - r * @a[pivot_row][j] } + @b[i] = @b[i] - r * @b[pivot_row] + } - update_solution - end - - def replace_basic_variable(hash) - from, to = hash.keys.first, hash.values.first - @basic_vars.delete(from) - @basic_vars << to - @basic_vars.sort! + self.update_solution end def pivot_row(column_ix) - row_ix_a_and_b = row_indices.map { |row_ix| - [row_ix, @a[row_ix][column_ix], @b[row_ix]] - }.reject { |_, a, b| - a == 0 - }.reject { |_, a, b| - (b < 0) ^ (a < 0) # negative sign check - } - row_ix, _, _ = *last_min_by(row_ix_a_and_b) { |_, a, b| - Rational(b, a) + min_ratio = nil + idx = nil + @num_constraints.times { |i| + a, b = @a[i][column_ix], @b[i] + next if a == 0 or (b < 0) ^ (a < 0) + ratio = Rational(b, a) + idx, min_ratio = i, ratio if min_ratio.nil? or ratio <= min_ratio } - row_ix - end - - def basic_variable_in_row(pivot_row) - column_indices.detect do |column_ix| - @a[pivot_row][column_ix] == 1 and @basic_vars.include?(column_ix) - end - end - - def row_indices - (0...@a.length).to_a - end - - def column_indices - (0...@a.first.size).to_a + idx end def formatted_tableau - if can_improve? - pivot_column = entering_variable - pivot_row = pivot_row(pivot_column) + if self.can_improve? + pivot_column = self.entering_variable + pivot_row = self.pivot_row(pivot_column) else pivot_row = nil end - num_cols = @c.size + 1 - c = formatted_values(@c.to_a) - b = formatted_values(@b.to_a) - a = @a.to_a.map {|ar| formatted_values(ar.to_a) } + c = @c.to_a.map { |flt| "%2.3f" % flt } + b = @b.to_a.map { |flt| "%2.3f" % flt } + a = @a.to_a.map { |vec| vec.to_a.map { |flt| "%2.3f" % flt } } if pivot_row a[pivot_row][pivot_column] = "*" + a[pivot_row][pivot_column] end max = (c + b + a + ["1234567"]).flatten.map(&:size).max result = [] - result << c.map {|c| c.rjust(max, " ") } + result << c.map { |str| str.rjust(max, " ") } a.zip(b) do |arow, brow| - result << (arow + [brow]).map {|a| a.rjust(max, " ") } + result << (arow + [brow]).map { |val| val.rjust(max, " ") } result.last.insert(arow.length, "|") end - lines = result.map {|b| b.join(" ") } + lines = result.map { |ary| ary.join(" ") } max_line_length = lines.map(&:length).max lines.insert(1, "-"*max_line_length) lines.join("\n") end - - def formatted_values(array) - array.map {|c| "%2.3f" % c } - end - - # like Enumerable#min_by except if multiple values are minimum - # it returns the last - def last_min_by(array) - best_element, best_value = nil, nil - array.each do |element| - value = yield element - if !best_element || value <= best_value - best_element, best_value = element, value - end - end - best_element - end - - def assert(boolean) - raise unless boolean - end - end - diff --git a/metrics/bench b/metrics/bench new file mode 100644 index 0000000..32918eb --- /dev/null +++ b/metrics/bench @@ -0,0 +1,5 @@ +/proc/loadavg 0.01 0.01 0.00 1/77 1173 +Warming up -------------------------------------- + Simplex Array 84.000 i/100ms +Calculating ------------------------------------- + Simplex Array 849.092 (± 1.3%) i/s - 2.604k in 3.067337s diff --git a/test/bench/simplex.rb b/test/bench/simplex.rb new file mode 100644 index 0000000..10728e0 --- /dev/null +++ b/test/bench/simplex.rb @@ -0,0 +1,98 @@ +require 'simplex' +require 'benchmark/ips' + +Benchmark.ips do |b| + b.config time: 3, warmup: 0.5 + + b.report("Simplex Array") { + Simplex.new([1, 1], + [[2, 1], + [1, 2]], + [4, 3]).solution + + Simplex.new([3, 4], + [[1, 1], + [2, 1]], + [4, 5]).solution + + Simplex.new([2, -1], + [[1, 2], + [3, 2],], + [6, 12]).solution + + Simplex.new([60, 90, 300], + [[1, 1, 1], + [1, 3, 0], + [2, 0, 1]], + [600, 600, 900]).solution + + Simplex.new([70, 210, 140], + [[1, 1, 1], + [5, 4, 4], + [40, 20, 30]], + [100, 480, 3200]).solution + + Simplex.new([2, -1, 2], + [[2, 1, 0], + [1, 2, -2], + [0, 1, 2]], + [10, 20, 5]).solution + + Simplex.new([11, 16, 15], + [[1, 2, Rational(3, 2)], + [Rational(2, 3), Rational(2, 3), 1], + [Rational(1, 2), Rational(1, 3), Rational(1, 2)]], + [12_000, 4_600, 2_400]).solution + + Simplex.new([5, 4, 3], + [[2, 3, 1], + [4, 1, 2], + [3, 4, 2]], + [5, 11, 8]).solution + + Simplex.new([3, 2, -4], + [[1, 4, 0], + [2, 4,-2], + [1, 1,-2]], + [5, 6, 2]).solution + + Simplex.new([2, -1, 8], + [[2, -4, 6], + [-1, 3, 4], + [0, 0, 2]], + [3, 2, 1]).solution + + Simplex.new([100_000, 40_000, 18_000], + [[20, 6, 3], + [0, 1, 0], + [-1,-1, 1], + [-9, 1, 1]], + [182, 10, 0, 0]).solution + + Simplex.new([1, 2, 1, 2], + [[1, 0, 1, 0], + [0, 1, 0, 1], + [1, 1, 0, 0], + [0, 0, 1, 1]], + [1, 4, 2, 2]).solution + + Simplex.new([10, -57, -9, -24], + [[0.5, -5.5, -2.5, 9], + [0.5, -1.5, -0.5, 1], + [ 1, 0, 0, 0]], + [0, 0, 1]).solution + + Simplex.new([25, 20], + [[20, 12], + [1, 1]], + [1800, 8*15]).solution + } + +#b.report("Simplex Array") { +#} + +#b.report("Simplex Matrix") { +#} + +# b.compare! +end diff --git a/test/simplex_test.rb b/test/simplex.rb similarity index 92% rename from test/simplex_test.rb rename to test/simplex.rb index 6b406e8..125e67c 100644 --- a/test/simplex_test.rb +++ b/test/simplex.rb @@ -1,8 +1,7 @@ -require 'test/unit' -$:.push(File.expand_path("../../lib", __FILE__)) +require 'minitest/autorun' require 'simplex' -class SimplexTest < Test::Unit::TestCase +class SimplexTest < Minitest::Test def test_2x2 result = Simplex.new( [1, 1], @@ -180,14 +179,14 @@ def test_cycle2 ], [0, 0] ) - assert_raise Simplex::UnboundedProblem do + assert_raises Simplex::UnboundedProblem do simplex.solution end end def test_error_mismatched_dimensions - assert_raise ArgumentError do - result = Simplex.new( + assert_raises ArgumentError do + Simplex.new( [10, -57, -9], [ [0.5, -5.5, -2.5, 9], @@ -198,8 +197,8 @@ def test_error_mismatched_dimensions ) end - assert_raise ArgumentError do - result = Simplex.new( + assert_raises ArgumentError do + Simplex.new( [10, -57, -9, 2], [ [0.5, -5.5, 9, 4], @@ -210,8 +209,8 @@ def test_error_mismatched_dimensions ) end - assert_raise ArgumentError do - result = Simplex.new( + assert_raises ArgumentError do + Simplex.new( [10, -57, -9, 2], [ [0.5, -5.5, 9, 4], @@ -263,16 +262,16 @@ def test_cup_factory # [1, -2] # ) # while simplex.can_improve? - # puts + # puts # puts simplex.formatted_tableau # simplex.pivot # end # p :done - # puts + # puts # puts simplex.formatted_tableau #end - + def test_unbounded simplex = Simplex.new( [1, 1, 1], @@ -282,7 +281,7 @@ def test_unbounded ], [5, 7] ) - assert_raise Simplex::UnboundedProblem do + assert_raises Simplex::UnboundedProblem do simplex.solution end end