macros/tableau.pl
DONE: change find_next_basis_from_pivot to next_basis_from_pivot
DONE: add phase2 to some match phase1 for some of the pivots -- added as main:: subroutine
DONE: add a generic _solve that solves tableau once it has a filled basis -- needs to be tested
add something that will fill the basis.
regularize the use of m and n -- should they be the number of
constraints and the number of decision(problem) variables or should
they be the size of the complete tableau.
we probably need both -- need names for everything
current tableau returns the complete tableau minus the objective function
row.
(do we need a second objective function row? -- think we can skip this for now.)
# We're going to have several types
# MathObject Matrices Value::Matrix
# tableaus form John Jones macros
# MathObject tableaus
# Containing an matrix $A coefficients for constraint
# A vertical vector $b for constants for constraints
# A horizontal vector $c for coefficients for objective function
# A vertical vector corresponding to the value $z of the objective function
# dimensions $n problem vectors, $m constraints = $m slack variables
# A basis Value::Set -- positions for columns which are independent and
# whose associated variables can be determined
# uniquely from the parameter variables.
# The non-basis (parameter) variables are set to zero.
#
# state variables (assuming parameter variables are zero or when given parameter variables)
# create the methods for updating the various containers
# create the method for printing the tableau with all its decorations
# possibly with switches to turn the decorations on and off.
The structure of the tableau is:
-----------------------------------------------------------
| | | | |
| A | S | 0 | b |
| | | | |
----------------------------------------------
| -c | 0 | 1 | 0 |
----------------------------------------------
Matrix A, the initial constraint matrix is n=num_problem_vars by m=num_slack_vars
Matrix S, the initial slack variables is m by m
Matrix b, the initial constraint constants is n by 1 (Matrix or ColumnVector)
Matrix c, the objective function coefficients matrix is 1 by n
Matrix which changes with state:
Matrix upper_tableau: m by n+2 (A|S|0|b)
Matrix tableau m+1 by n+2 (A/(-c)) | S/0 | 0/1 |b/z
Matrix current_obj_coeff = 1 by n+m matrix (negative of current coeff for obj_function)
Matrix current_last_row = 1 by n+m+2 matrix tableau is upper_tableau over current_last_row
The next to the last column holds z or objective value
z(...x^i...) = c_i* x^i (Einstein summation convention)
FIXME: ?? allow c to be a 2 by n matrix so that you can do phase1 calculations easily
ANS( $tableau->cmp(checker=>tableauEquivalence()) );
# Note: it is important to include the () at the end of tableauEquivalence
# tableauEquivalence compares two matrices up to
# reshuffling the rows and multiplying each row by a constant.
# It is equivalent up to multiplying on the left by a permuation matrix
# or a (non-uniformly constant) diagonal matrix.
# It is appropriate for comparing augmented matrices representing a system of equations
# since the order of the equations is unimportant. This applies to tableaus for
# Linear Optimization Problems being solved using the simplex method.
Returns the solution variables to the tableau assuming that the parameter (non-basis) variables have been set to zero. It returns a list in array context and a reference to an array in scalar context.
Parameters: ()
Return:
Useage:
ANS($constraints->cmp()->withPostFilter(
linebreak_at_commas()
));
Replaces commas with line breaks in the latex presentations of the answer checker. Used most often when $constraints is a LinearInequality math object.
MathObject Matrix methods: http://webwork.maa.org/wiki/Matrix_(MathObject_Class) MathObject Contexts: http://webwork.maa.org/wiki/Common_Contexts CPAN RealMatrix docs: http://search.cpan.org/~leto/Math-MatrixReal-2.09/lib/Math/MatrixReal.pm
More references: "Matrix.pm" in lib
Tableau->new(A=>Matrix, b=>Vector or Matrix, c=>Vector or Matrix)
A => undef, # original constraint matrix MathObjectMatrix
b => undef, # constraint constants ColumnVector or MathObjectMatrix n by 1
c => undef, # coefficients for objective function Vector or MathObjectMatrix 1 by n
obj_row => undef, # contains the negative of the coefficients of the objective function.
z => undef, # value for objective function Real
n => undef, # dimension of problem variables (columns in A)
m => undef, # dimension of slack variables (rows in A)
S => undef, # square m by m MathObjectMatrix for slack variables. default is the identity
M => undef, # matrix (m by m+n+1+1) consisting of all original columns and all
rows except for the objective function row. The m+n+1 column and
is the objective_value column. It is all zero with a pivot in the objective row.
The current version of this accessed by Tableau->upper_tableau (A | S |0 | b)
#FIXME
obj_col_num => undef, # have obj_col on the left or on the right? FIXME? obj_column_position
# perhaps not store this column at all and add it when items are printed?
basis => List | Set, # unordered list describing the current basis columns corresponding
to determined variables. With a basis argument this sets a new state defined by that basis.
current_constraint_matrix=>(m by n matrix), # the current version of [A | S]
current_b=> (1 by m matrix or Column vector) # the current version of the constraint constants b
current_basis_matrix => (m by m invertible matrix) a square invertible matrix
# corresponding to the current basis columns
# flag indicating the column (1 or n+m+1) for the objective value
constraint_labels => undef, (not sure if this remains relevant after pivots)
problem_var_labels => undef,
slack_var_labels => undef,
Notes: (1 by m MathObjectMatrix) <= Value::Matrix->new($self->b->transpose->value )
Useage:
lop_display($tableau, align=>'cccc|cc|c|c', toplevel=>[qw(x1,x2,x3,x4,s1,s2,P,b)])
Pretty prints the output of a matrix as a LOP with separating labels and variable labels.
[complementary_basis_set] = $self->primal_basis_to_dual(primal_basis_set)
[complementary_basis_set] = $self->dual_basis_to_primal(dual_basis_set)
ARRAY reference = $self->basis_columns()
[3,4] = $self->basis_columns([3,4])
Sets or returns the basis_columns as an ARRAY reference
$self->objective_row
Parameters: ()
Returns:
$self->current_tableau
Parameters: () or (list)
Returns: A MathObjectMatrix
Useage:
$MathObjectmatrix = $self->current_tableau
$MathObjectmatrix = $self->current_tableau(3,4) #updates basis to (3,4)
Returns the current constraint matrix as a MathObjectMatrix, including the constraint constants, problem variable coefficients, slack variable coefficients AND the row containing the objective function coefficients. ------------- |A | S |0| b| ------------- | -c |z|z*| -------------
If a list of basis columns is passed as an argument then $self->basis() is called to switch the tableau to the new basis before returning the tableau.
[x1,.....xn,]= $self->statevars()
MathObjectList = $self->basis
MathObjectList = $self->basis(3,4)
MathObjectList = $self->basis([3,4])
MathObjectList = $self->basis(Set(3,4))
MathObjectList = $self->basis(List(3,4))
to obtain ARRAY reference use
[3,4]== $self->basis(Set3,4)->value
Returns a MathObjectList containing the current basis columns. If basis columns are provided as arguments it updates all elements of the tableau to present the view corresponding to the new choice of basis columns.
($col1,$col2,..,$flag) = $self->find_next_basis (max/min, obj_row_number)
In phase 2 of the simplex method calculates the next basis. $optimum or $unbounded is set if the process has found on the optimum value, or the column $col gives a certificate of unboundedness.
$flag can be either 'optimum' or 'unbounded' in which case the basis returned is the current basis. is a list of column numbers.
FIXME Should we change this so that ($basis,$flag) is returned instead? $basis and $flag are very different things. $basis could be a set or list type but in that case it can't have undef as an entry. It probably can have '' (an empty string)
($row, $col,$optimum,$unbounded) = $self->find_next_pivot (max/minm obj_row_number)
This is used in phase2 so the possible outcomes are only $optimum and $unbounded. $infeasible is not possible. Use the lowest index strategy to find the next pivot point. This calls find_pivot_row and find_pivot_column. $row and $col are undefined if either $optimum or $unbounded is set.
List(basis) = $self->find_next_basis_from_pivot (pivot_row, pivot_column)
Calculate the next basis from the current basis given the pivot position.
($index, $value, $optimum) = $self->find_pivot_column (max/min, obj_row_number)
This finds the left most obj function coefficient that is negative (for maximizing) or positive (for minimizing) and returns the value and the index. Only the index is really needed for this method. The row number is included because there might be more than one objective function in the table (for example when using the Auxiliary method in phase1 of the simplex method.) If there is no coefficient of the appropriate sign then the $optimum flag is set and $index and $value are undefined.
($index, $value, $unbounded) = $self->find_pivot_row(col_number)
Compares the ratio $b[$j]/a[$j, $col_number] and chooses the smallest non-negative entry. It assumes that we are in phase2 of simplex methods so that $b[j]>0; If all entries are negative (or infinity) then the $unbounded flag is set and returned and the $index and $value quantities are undefined.
($index, $value) = $self->find_leaving_column(obj_row_number)
Finds the non-basis column with a non-zero entry in the given row. When called with the pivot row number this index gives the column which will be removed from the basis while the pivot col number gives the basis column which will become a parameter column.
($row, $col, $feasible, $infeasible) = $self->next_short_cut_pivot
Following the short-cut algorithm this chooses the next pivot by choosing the row with the most negative constraint constant entry (top most first in case of tie) and then the left most negative entry in that row.
The process stops with either $feasible=1 (state variables give a feasible point for the constraints) or $infeasible=1 (a row in the tableau shows that the LOP has empty domain.)
($basis->value, $flag) = $self->next_short_cut_basis()
In phase 1 of the simplex method calculates the next basis for the short cut method. $flag is set to 'feasible_point' if the basis and its corresponding tableau is associated with a basic feasible point (a point on a corner of the domain of the LOP). The tableau is ready for phase 2 processing. $flag is set to 'infeasible_lop' which means that the tableau has a row which demonstrates that the LOP constraints are inconsistent and the domain is empty. In these cases the basis returned is the current basis of the tableau object.
Otherwise the $basis->value returned is the next basis that should be used in the short_cut method and $flag contains undef.
($index, $value, $feasible)=$self->find_short_cut_row
Find the most negative entry in the constraint column vector $b. If all entries are positive then the tableau represents a feasible point, $feasible is set to 1 and $index and $value are undefined.
($index, $value, $infeasible) = $self->find_short_cut_column(row_index)
Find the left most negative entry in the specified row. If all coefficients are positive then the tableau represents an infeasible LOP, the $infeasible flag is set, and the $index and $value are undefined.
(or tableau pivot???)
Tableau = $self->row_reduce(3,4)
MathObjectMatrix = $self->row_reduce(3,4)->current_tableau
Row reduce matrix so that column 4 is a basis column. Used in pivoting for simplex method. Returns tableau object.
TableauObject = $self->dual_lop
Creates the tableau of the LOP which is dual to the linear optimization problem represented by the current tableau.
These are specialized routines used in the simplex method
@array = $self->primal2dual(2,3,4)
Maps LOP column indices to dual LOP indicies (basis of complementary slack property)
@array = $self->dual2primal(2,3,4)
Maps dual LOP column indices to primal LOP indicies (basis of complementary slack property). Inverse of primal2dual method.
$self->isOptimal('min'| 'max')
Returns 1 or 0
This checks to see if the state is a local minimum or maximum for the objective function -- it does not check whether the stateis feasible.
Checks to see if the current state is feasible or whether it requires further phase 1 processing.
These are generic matrix routines. Perhaps some or all of these should be added to the file Value::Matrix?
$self->row_slice
Parameter: @slice or \@slice
Return: MathObject matrix
MathObjectMatrix = $self->row_slice(3,4)
MathObjectMatrix = $self->row_slice([3,4])
Similar to $self->extract_rows (or $self->rows) but returns a MathObjectmatrix
$self->extract_rows
Parameter: @slice or \@slice
Return: two dimensional array ref
ARRAY reference = $self->extract_rows(@slice)
ARRAY reference = $self->extract_rows([@slice])
$self->column_slice()
Parameter: @slice or \@slice
Return: two dimensional array ref
ARRAY reference = $self->extract_rows(@slice)
ARRAY reference = $self->extract_rows([@slice])
$self->extract_columns
Parameter: @slice or \@slice
Return: two dimensional array ref
ARRAY reference = $self->extract_columns(@slice)
ARRAY reference = $self->extract_columns([@slice])
Parameter: @slice or \@slice
Return: MathObject List of row references
MathObjectList = $self->extract_rows_to_list(@slice)
MathObjectList = $self->extract_rows_to_list([@slice])
$self->extract_columns_to_list
Parameter: @slice or \@slice
Return: MathObject List of Matrix references ?
ARRAY reference = $self->extract_columns_to_list(@slice)
ARRAY reference = $self->extract_columns_to_list([@slice])
$self->submatrix
Parameter:(rows=>\@row_slice,columns=>\@column_slice)
Return: MathObject matrix
MathObjectMatrix = $self->submatrix([[1,2,3],[2,4,5]])
Extracts a submatrix from a Matrix and returns it as MathObjectMatrix.
Indices for MathObjectMatrices start at 1.
$Matrix->change_matrix_entry([i,j,k],$value)
Taken from MatrixReduce.pl. Written by Davide Cervone.
perhaps "assign" would be a better name for this?