Algorithm-Bertsekas
view release on metacpan or search on metacpan
lib/Algorithm/Bertsekas.pm view on Meta::CPAN
package Algorithm::Bertsekas;
use strict;
use warnings FATAL => 'all';
use diagnostics;
require Exporter;
our @ISA = qw(Exporter);
our @EXPORT = qw( auction );
our $VERSION = '0.87';
#Variables global to the package
my $maximize_total_benefit;
my $matrix_spaces; # used to print messages on the screen
my $decimals; # the number of digits after the decimal point
my ( $array1_size, $array2_size, $min_size, $max_size, $original_max_size );
my ( $need_transpose, $inicial_price, $iter_count_global, $iter_count_local );
my ( $epsilon_scaling, $max_epsilon_scaling, $max_matrix_value, $target, $output );
my ( %index_correlation, %assignned_object, %assignned_person, %price_object );
my ( %objects_desired_by_this, %locked_list, %seen_person, %seen_assignned_objects );
sub auction { # => default values
my %args = ( matrix_ref => undef, # reference to array: matrix N x M
maximize_total_benefit => 0, # 0: minimize_total_benefit ; 1: maximize_total_benefit
inicial_stepsize => undef, # auction algorithm terminates with a feasible assignment if the problem data are integer and stepsize < 1/min(N,M)
inicial_price => 0,
verbose => 3, # level of verbosity, 0: quiet; 1, 2, 3, 4, 5, 8, 9, 10: debug information.
@_, # argument pair list goes here
);
$max_matrix_value = 0;
$iter_count_global = 0;
$epsilon_scaling = 0;
$need_transpose = 0;
%index_correlation = ();
%assignned_object = ();
%assignned_person = ();
%price_object = ();
%objects_desired_by_this = ();
%locked_list = ();
%seen_person = ();
my @matrix_input = @{$args{matrix_ref}}; # Input: Reference to the input matrix (NxM) = $min_size x $max_size
$array1_size = $#matrix_input + 1;
$array2_size = $#{$matrix_input[0]} + 1;
$min_size = $array1_size < $array2_size ? $array1_size : $array2_size ; # square matrix --> $min_size = $max_size and $array1_size = $array2_size
$max_size = $array1_size < $array2_size ? $array2_size : $array1_size ;
$original_max_size = $max_size;
$target = 'auction-' . $array1_size . 'x' . $array2_size . '-output.txt' ;
if ( $args{verbose} >= 8 ){
print "\n verbose = $args{verbose} ---> print the verbose messages to <$target> file \n";
if ( open ( $output, '>', $target ) ) {
print "\n *** Open <$target> for writing. *** \n";
} else {
$args{verbose} = 7;
warn "\n *** Could not open <$target> for writing: $! *** \n";
}
}
$maximize_total_benefit = $args{maximize_total_benefit};
my $optimal_benefit = 0;
my %assignment_hash; # assignment: a hash representing edges in the mapping, as in the Algorithm::Kuhn::Munkres.
my @output_index; # output_index: an array giving the number of the value assigned, as in the Algorithm::Munkres.
my @matrix;
my @matrix_index;
foreach ( @matrix_input ){ # copy the orginal matrix N x M
push @matrix, [ @$_ ];
}
if ( $max_size <= 1 ){ # matrix_input 1 x 1
$assignment_hash{0} = 0;
$output_index[0] = 0;
$matrix_index[0] = 0;
$optimal_benefit = $matrix_input[0]->[0];
}
$need_transpose = 1 if ( $array1_size > $array2_size ); # will always be chosen N <= M
if ( $need_transpose ){
my $transposed = transpose(\@matrix);
@matrix = @$transposed;
}
get_matrix_info( \@matrix, $args{verbose} );
delete_multiple_columns( \@matrix, $args{verbose} ) if ( $min_size >= 2 and $min_size != $max_size );
# epsilon is the stepsize and auction algorithm terminates with a feasible assignment if the problem data are integer and epsilon < 1/min(N,M).
# There is a trade-off between runtime and the chosen stepsize. Using the largest possible increment accelerates the algorithm.
$inicial_price = $args{inicial_price};
$price_object{$_} = $inicial_price for ( 0 .. $max_size - 1 );
#$price_object{$_} = sprintf( "%.0f", rand($max_matrix_value) ) for ( 0 .. $max_size - 1 ); # random values for initial prices
# inicial epsilon value
my $epsilon = int (($max_matrix_value/2) * exp ( -1 * $max_size/$min_size )); # exp (1) = e = 2.71828182845905
$epsilon = $args{inicial_stepsize} if ( defined $args{inicial_stepsize} );
$epsilon = 1/(1+$min_size) if ($epsilon < 1/$min_size);
# (inicial epsilon value) / factor ** k = 1/$min_size --> epsilon * $min_size = factor ** k
# log(epsilon * $min_size) = k * log(factor) --> log(factor) = (log(epsilon * $min_size)) / k --> factor = exp ( (log(epsilon * $min_size)) / k )
#$max_epsilon_scaling = 20;
#my $factor = exp ( (log( $epsilon * (1+$min_size) )) / ($max_epsilon_scaling-1) ); print "\n \$max_epsilon_scaling = $max_epsilon_scaling ; \$factor = $factor \n";
lib/Algorithm/Bertsekas.pm view on Meta::CPAN
$epsilon_scaling++;
$iter_count_local = 0;
%assignned_object = ();
%assignned_person = ();
%seen_person = ();
$seen_assignned_objects{$_} = 0 for ( 0 .. $max_size - 1 );
while ( (scalar keys %assignned_person) < $max_size ){ # while there is at least one element not assigned.
$iter_count_global++;
$iter_count_local++;
auctionRound( \@matrix, $epsilon, $args{verbose} );
if ( $args{verbose} >= 10 ){
for my $i ( -1 .. $#matrix ) {
if ($i >= 0){ printf $output " %2s [", $i; } else{ printf $output "object"; }
for my $j ( 0 .. $#{$matrix[$i]} ) {
if ($i >= 0){ printf $output (" %${matrix_spaces}.${decimals}f", $matrix[$i]->[$j]); } else{ printf $output (" %${matrix_spaces}.0f", $j); }
if ( defined $assignned_person{$i} and $j == $assignned_person{$i} ){ print $output "**"; } else{ print $output " "; }
}
if ($i >= 0){ print $output "]\n"; } else{ print $output "\n\n"; }
}
}
}
$epsilon = $epsilon / $factor ; # (1/2): smooth convergence
if ( not $feasible_assignment_condition and $epsilon < 1/$min_size ){
$epsilon = 1/(1+$min_size);
$feasible_assignment_condition = 1;
}
}
$epsilon = $factor * $epsilon; # correcting information for printing
my %seeN;
my %seeM;
foreach my $person ( sort { $a <=> $b } keys %assignned_person ){
my $object = $assignned_person{$person};
$matrix_index[$person] = $object;
#print " \$need_transpose = $need_transpose ; \$matrix_index[$person] = $object ; \$index_i = $person ; \$index_j = $object --> $index_correlation{$object} ;";
my $index_i = $need_transpose ? $index_correlation{$object} // $object : $person;
my $index_j = $need_transpose ? $person : $index_correlation{$object} // $object;
$output_index[$index_i] = $index_j;
$seeN{$index_i}++;
$seeM{$index_j}++;
#print " \$output_index[$index_i] = $index_j \n";
next unless ( defined $matrix_input[$index_i] and defined $matrix_input[$index_i]->[$index_j] );
$assignment_hash{ $index_i } = $index_j;
$optimal_benefit += $matrix_input[$index_i]->[$index_j];
}
for my $i ( 0 .. $original_max_size - 1 ) {
for my $j ( 0 .. $original_max_size - 1 ) {
next if ($seeN{$i} or $seeM{$j});
$output_index[$i] = $j;
$seeN{$i}++;
$seeM{$j}++;
last;
}}
if ( $args{verbose} >= 8 ){
printf $output "\n\$optimal_benefit = $optimal_benefit ; \$iter_count_global = $iter_count_global ; \$epsilon = %.4g ; \@output_index = (@output_index) \n", $epsilon;
}
print_screen_messages( \@matrix, \@matrix_index, \@matrix_input, \@output_index, $optimal_benefit, $args{verbose}, $epsilon ) ;
return ( $optimal_benefit, \%assignment_hash, \@output_index ) ;
}
sub transpose {
my $matrix_ref = shift;
my @transpose;
for my $i ( 0 .. $#{$matrix_ref} ) {
for my $j ( 0 .. $#{$matrix_ref->[$i]} ) {
$transpose[$j]->[$i] = $matrix_ref->[$i]->[$j];
}
}
return \@transpose;
}
sub delete_multiple_columns { # if the column elements do not change the final result
my ( $matrix_ref, $verbose ) = @_;
my %lower_values;
my %intersection_columns;
my $number_of_columns_deleted = 0;
for my $i ( 0 .. $#{$matrix_ref} ) {
for my $j ( 0 .. $#{$matrix_ref->[$i]} ) {
$lower_values{ $i }{ $matrix_ref->[$i]->[$j] }{ $j }++;
}
}
# consider N rows < M columns
# remove the matching columns whose elements are never among the N largest elements in each row
foreach my $index_i ( sort { $a <=> $b } keys %lower_values ){
my $num_higher_values = 0;
foreach my $matrix_value ( sort { $b <=> $a } keys %{$lower_values{$index_i}} ){
foreach my $index_j ( sort { $b <=> $a } keys %{$lower_values{$index_i}{$matrix_value}} ){
$intersection_columns{$index_j}++ if ( $num_higher_values++ >= $min_size );
$number_of_columns_deleted++ if ( defined $intersection_columns{$index_j} and $intersection_columns{$index_j} >= $min_size );
}
}
}
if ( $verbose >= 5 ){
print "\n";
for my $i ( 0 .. $#{$matrix_ref} ) {
print " [";
for my $j ( 0 .. $#{$matrix_ref->[$i]} ) {
printf (" %${matrix_spaces}.${decimals}f", $matrix_ref->[$i]->[$j] );
if ( defined $intersection_columns{$j} and $intersection_columns{$j} == $min_size ){ print "**"; } else{ print " "; }
}
lib/Algorithm/Bertsekas.pm view on Meta::CPAN
print " [";
for my $j ( 0 .. $#{$matrix_ref->[$i]} ) {
printf (" %${matrix_spaces}.${decimals}f ", $matrix_ref->[$i]->[$j] );
}
print "]\n";
}
print "\n";
}
$max_size = $max_size - $number_of_columns_deleted;
}
sub print_screen_messages {
my ( $matrix_ref, $matrix_index_ref, $matrix_input_ref, $output_index_ref, $optimal_benefit, $verbose, $epsilon ) = @_;
my @matrix = @$matrix_ref;
my @matrix_index = @$matrix_index_ref;
my @matrix_input = @$matrix_input_ref;
my @output_index = @$output_index_ref;
if ( $verbose >= 1 ){
print "\nObjective: ";
printf( $maximize_total_benefit ? "to Maximize the total benefit\n" : "to Minimize the total benefit\n" );
printf(" Number of left nodes: %u\n", $array1_size );
printf(" Number of right nodes: %u\n", $array2_size );
printf(" Number of edges: %u\n", $array1_size * $array2_size );
print "\nSolution:\n";
printf(" Optimal assignment: sum of values = %.${decimals}f \n", $optimal_benefit );
printf(" Feasible assignment condition: stepsize = %.4g < 1/$min_size = %.4g \n", $epsilon, 1/$min_size ) if ( $verbose >= 1 and $max_size >= 2 );
printf(" Number of iterations: %u \n", $iter_count_global ) if ( $verbose >= 1 );
print "\n row index = [";
for my $i ( 0 .. $#output_index ) {
printf("%${matrix_spaces}d ", $i);
}
print "]\n";
print " column index = [";
for my $index (@output_index) {
printf("%${matrix_spaces}d ", $index);
}
print "]\n";
print " matrix value = [";
for my $i ( 0 .. $#output_index ){
my $j = $output_index[$i];
last if not defined $j;
my $weight;
$weight = ( defined $matrix_input[$i] and defined $matrix_input[$i]->[$j] ) ? sprintf( "%${matrix_spaces}.${decimals}f ", $matrix_input[$i]->[$j] ) : ' ' x ($matrix_spaces+1) ;
print $weight;
}
print "]\n\n";
}
if ( $verbose >= 2 ){
my $index_length = length($original_max_size);
if ( $verbose >= 3 ){
printf " modified matrix %d x %d:\n", $#matrix + 1, $#{$matrix[0]} + 1;
for my $i ( 0 .. $#matrix ) {
print " [";
for my $j ( 0 .. $#{$matrix[$i]} ) {
printf(" %${matrix_spaces}.${decimals}f", $matrix[$i]->[$j] );
if ( $j == $matrix_index[$i] ){ print "**"; } else{ print " "; }
}
print "]\n";
}
print "\n";
}
print " original matrix $array1_size x $array2_size with solution:\n";
for my $i ( 0 .. $#matrix_input ) {
print " [";
for my $j ( 0 .. $#{$matrix_input[$i]} ) {
printf(" %${matrix_spaces}.${decimals}f", $matrix_input[$i]->[$j] );
if ( $j == $output_index[$i] ){ print "**"; } else{ print " "; }
}
print "]\n";
}
my %orderly_solution;
for my $i ( 0 .. $original_max_size - 1 ){
my $j = $output_index[$i];
my $weight = $max_matrix_value;
$weight = $matrix_input[$i]->[$j] if ( defined $matrix_input[$i] and defined $matrix_input[$i]->[$j] ); # condition for valid solution
$orderly_solution{ $weight } { $i } { 'index_array1' } = $i;
$orderly_solution{ $weight } { $i } { 'index_array2' } = $j;
}
print "\n Pairs (in ascending order of matrix values):\n";
my $sum_matrix_value = 0;
my $sum_spaces = 2 * $matrix_spaces - 1;
foreach my $matrix_value ( sort { $a <=> $b } keys %orderly_solution ){
foreach my $k ( sort { $a <=> $b } keys %{$orderly_solution{$matrix_value}} ){
my $index_array1 = $orderly_solution{ $matrix_value } { $k } { 'index_array1' };
my $index_array2 = $orderly_solution{ $matrix_value } { $k } { 'index_array2' };
$sum_matrix_value += $matrix_value if ( defined $matrix_input[$index_array1] and defined $matrix_input[$index_array1]->[$index_array2] );
my $weight = ( defined $matrix_input[$index_array1] and defined $matrix_input[$index_array1]->[$index_array2] ) ? sprintf( "%${matrix_spaces}.${decimals}f", $matrix_value ) : ' ' x $matrix_spaces ;
printf( " indexes ( %${index_length}d, %${index_length}d ), matrix value = $weight ; sum of values = %${sum_spaces}.${decimals}f \n", $index_array1, $index_array2, $sum_matrix_value );
}}
}
}
sub get_matrix_info {
my ( $matrix_ref, $verbose ) = @_;
my @matrix = @$matrix_ref;
my $min_matrix_value;
for my $i ( 0 .. $#matrix ) {
for my $j ( 0 .. $#{$matrix[$i]} ) {
my $char_number = length( $matrix[$i]->[$j] ); # count the number of characters
$matrix_spaces = $char_number if ( (not defined $matrix_spaces) || ($char_number > $matrix_spaces) );
$max_matrix_value = $matrix[$i]->[$j] if ( (not defined $max_matrix_value) || ($matrix[$i]->[$j] > $max_matrix_value) );
$min_matrix_value = $matrix[$i]->[$j] if ( (not defined $min_matrix_value) || ($matrix[$i]->[$j] < $min_matrix_value) );
}}
$decimals = length(($max_matrix_value =~ /[,.](\d+)/)[0]); # counting the number of digits after the decimal point
$decimals = 0 unless ( defined $decimals ); # for integers $decimals = 0
my $range = $max_matrix_value - $min_matrix_value; # $range >= 0
$range = 1 if ($range == 0);
if ( $verbose >= 4 ){
print "\n min_matrix_value = $min_matrix_value ; max_matrix_value = $max_matrix_value ; range = $range ; matrix_spaces = $matrix_spaces ; decimals = $decimals \n";
}
if ( $maximize_total_benefit ){
for my $i ( 0 .. $#matrix ) {
for my $j ( 0 .. $#{$matrix[$i]} ) {
$matrix[$i]->[$j] = $matrix[$i]->[$j] - $min_matrix_value ;
lib/Algorithm/Bertsekas.pm view on Meta::CPAN
printf " Auction Algorithm, (row, column) indexes --> \$i = %3d ; \$j = %3d ; \$value = $value ; \$sum = %8s \n", $i, $j, $sum;
}
### --- final --- ###
### --- simple and direct application --- ###
Example 1: Find the nearest neighbor, Minimize the total benefit.
my @array1 = ( 893, 401, 902, 576, 767, 917, 76, 464, 124, 207, 125, 530 );
my @array2 = ( 161, 559, 247, 478, 456 );
my @input_matrix;
for my $i ( 0 .. $#array1 ){
my @weight_function;
for my $j ( 0 .. $#array2 ){
my $weight = abs ($array1[$i] - $array2[$j]);
# $weight = ($array1[$i] - $array2[$j]) ** 2; # another option
push @weight_function, $weight;
}
push @input_matrix, \@weight_function;
}
161 559 247 478 456
893 [ 732 334 646 415 437 ]
401 [ 240 158 154 77 55 ]
902 [ 741 343 655 424 446 ]
576 [ 415 17 329 98 120 ]
767 [ 606 208 520 289 311 ]
917 [ 756 358 670 439 461 ]
76 [ 85 483 171 402 380 ]
464 [ 303 95 217 14 8 ]
124 [ 37 435 123 354 332 ]
207 [ 46 352 40 271 249 ]
125 [ 36 434 122 353 331 ]
530 [ 369 29 283 52 74 ]
my ( $optimal, $assignment_ref, $output_index_ref ) = auction( matrix_ref => \@input_matrix, maximize_total_benefit => 0, verbose => 5 );
Objective: to Minimize the total benefit
Number of left nodes: 12
Number of right nodes: 5
Number of edges: 60
Solution:
Optimal assignment: sum of values = 153
Feasible assignment condition: stepsize = 0.1667 < 1/5 = 0.2
Number of iterations: 50
row index = [ 0 1 2 3 4 5 6 7 8 9 10 11 ]
column index = [ 9 8 10 1 5 11 7 4 6 2 0 3 ]
matrix value = [ 17 8 40 36 52 ]
modified matrix 5 x 9:
[ 516 341 150 671 453 719 710 720** 387 ]
[ 598 739** 548 273 661 321 404 322 727 ]
[ 602 427 236 585 539 633 716** 634 473 ]
[ 679 658 467 354 742 402 485 403 704**]
[ 701 636 445 376 748** 424 507 425 682 ]
original matrix 12 x 5 with solution:
[ 732 334 646 415 437 ]
[ 240 158 154 77 55 ]
[ 741 343 655 424 446 ]
[ 415 17** 329 98 120 ]
[ 606 208 520 289 311 ]
[ 756 358 670 439 461 ]
[ 85 483 171 402 380 ]
[ 303 95 217 14 8**]
[ 37 435 123 354 332 ]
[ 46 352 40** 271 249 ]
[ 36** 434 122 353 331 ]
[ 369 29 283 52** 74 ]
Pairs (in ascending order of matrix values):
indexes ( 7, 4 ), matrix value = 8 ; sum of values = 8
indexes ( 3, 1 ), matrix value = 17 ; sum of values = 25
indexes ( 10, 0 ), matrix value = 36 ; sum of values = 61
indexes ( 9, 2 ), matrix value = 40 ; sum of values = 101
indexes ( 11, 3 ), matrix value = 52 ; sum of values = 153
indexes ( 0, 9 ), matrix value = ; sum of values = 153
indexes ( 1, 8 ), matrix value = ; sum of values = 153
indexes ( 2, 10 ), matrix value = ; sum of values = 153
indexes ( 4, 5 ), matrix value = ; sum of values = 153
indexes ( 5, 11 ), matrix value = ; sum of values = 153
indexes ( 6, 7 ), matrix value = ; sum of values = 153
indexes ( 8, 6 ), matrix value = ; sum of values = 153
Example 2: Maximize the total benefit.
my $N = 10;
my $M = 10;
my $r = 100;
my @input_matrix;
for my $i ( 0 .. $N - 1 ){
my @weight_function;
for my $j ( 0 .. $M - 1 ){
my $weight = sprintf( "%.0f", rand($r) );
push @weight_function, $weight;
}
push @input_matrix, \@weight_function;
}
Alternatively, we can define the matrix with its elements:
my @input_matrix = (
[ 84, 94, 75, 56, 66, 95, 39, 53, 73, 4 ],
[ 76, 71, 56, 49, 29, 1, 40, 40, 72, 72 ],
[ 85, 100, 71, 23, 47, 18, 82, 70, 30, 71 ],
[ 2, 95, 71, 89, 73, 73, 48, 52, 90, 51 ],
[ 65, 28, 77, 73, 24, 28, 75, 48, 8, 81 ],
[ 25, 27, 35, 89, 98, 10, 99, 3, 27, 4 ],
[ 58, 15, 99, 37, 92, 55, 52, 82, 73, 96 ],
[ 11, 75, 2, 1, 88, 43, 8, 28, 98, 20 ],
[ 52, 95, 10, 38, 41, 64, 20, 75, 1, 47 ],
[ 50, 80, 31, 90, 10, 83, 51, 55, 57, 40 ]
);
my ( $optimal, $assignment_ref, $output_index_ref ) = auction( matrix_ref => \@input_matrix, maximize_total_benefit => 1, verbose => 3 );
Objective: to Maximize the total benefit
Number of left nodes: 10
Number of right nodes: 10
Number of edges: 100
Solution:
Optimal assignment: sum of values = 893
Feasible assignment condition: stepsize = 0.09091 < 1/10 = 0.1
Number of iterations: 27
row index = [ 0 1 2 3 4 5 6 7 8 9 ]
column index = [ 5 0 1 8 9 6 2 4 7 3 ]
matrix value = [ 95 76 100 90 81 99 99 88 75 90 ]
original matrix 10 x 10 with solution:
[ 84 94 75 56 66 95** 39 53 73 4 ]
[ 76** 71 56 49 29 1 40 40 72 72 ]
[ 85 100** 71 23 47 18 82 70 30 71 ]
[ 2 95 71 89 73 73 48 52 90** 51 ]
[ 65 28 77 73 24 28 75 48 8 81**]
[ 25 27 35 89 98 10 99** 3 27 4 ]
[ 58 15 99** 37 92 55 52 82 73 96 ]
[ 11 75 2 1 88** 43 8 28 98 20 ]
[ 52 95 10 38 41 64 20 75** 1 47 ]
[ 50 80 31 90** 10 83 51 55 57 40 ]
Pairs (in ascending order of matrix values):
indexes ( 8, 7 ), matrix value = 75 ; sum of values = 75
indexes ( 1, 0 ), matrix value = 76 ; sum of values = 151
indexes ( 4, 9 ), matrix value = 81 ; sum of values = 232
indexes ( 7, 4 ), matrix value = 88 ; sum of values = 320
indexes ( 3, 8 ), matrix value = 90 ; sum of values = 410
indexes ( 9, 3 ), matrix value = 90 ; sum of values = 500
indexes ( 0, 5 ), matrix value = 95 ; sum of values = 595
indexes ( 5, 6 ), matrix value = 99 ; sum of values = 694
indexes ( 6, 2 ), matrix value = 99 ; sum of values = 793
indexes ( 2, 1 ), matrix value = 100 ; sum of values = 893
=head1 OPTIONS
matrix_ref => \@input_matrix, reference to array: matrix N x M.
maximize_total_benefit => 0, 0: minimize the total benefit ; 1: maximize the total benefit.
inicial_stepsize => 1, auction algorithm terminates with a feasible assignment if the problem data are integer and stepsize < 1/min(N,M).
inicial_price => 0,
verbose => 3, print messages on the screen, level of verbosity, 0: quiet; 1, 2, 3, 4, 5, 8, 9, 10: debug information.
=head1 EXPORT
"auction" function by default.
=head1 INPUT
The input matrix should be in a two dimensional array (array of array)
and the 'auction' subroutine expects a reference to this array.
=head1 OUTPUT
The $output_index_ref is the reference to the output_index array.
The $assignment_ref is the reference to the assignment hash.
The $optimal is the total benefit which can be a minimum or maximum value.
=head1 SEE ALSO
1. Network Optimization: Continuous and Discrete Models (1998).
Dimitri P. Bertsekas
http://web.mit.edu/dimitrib/www/netbook_Full_Book.pdf
2. Towards auction algorithms for large dense assignment problems (2008).
Libor Bus and Pavel Tvrdik
https://pdfs.semanticscholar.org/b759/b8fb205df73c810b483b5be2b1ded62309b4.pdf
3. https://github.com/EvanOman/AuctionAlgorithmCPP/blob/master/auction.cpp
This Perl algorithm started from this C++ implementation.
( run in 0.823 second using v1.01-cache-2.11-cpan-140bd7fdf52 )