LinearSolver.mesa, a definition of a Cedar module for solving linear equations and inequations by gaussian elimination and the simplex method. Part of the Juno 2 graphics programming system.
Coded by Greg Nelson, November 1983.
Last Edited by: Gnelson, December 11, 1983 10:45 pm
Last Edited by: Beach, December 11, 1984 7:55:06 pm PST
LinearSolver: CEDAR DEFINITIONS = BEGIN
RowMax: INT = 100;
Row: TYPE = [0 .. RowMax);
Col:
TYPE = [0 .. ColMax);
Tableau: TYPE = REF TableauRec;
TableauRec:
TYPE =
RECORD [
n: [0 .. RowMax],
m: [0 .. ColMax],
maxN: [0 .. RowMax],
maxM: [0 .. ColMax],
a:
REF
ARRAY Row
OF
ARRAY Col
OF
REAL,
a[i][j] is only relevant if i is in [0, n) and j is in [0, m). Thus we declare a matrix of n rows by m columns, which can grow as large as RowMax by ColMax.
x: REF ARRAY Col OF Unknown,
y: REF ARRAY Row OF Unknown,
unity: Unknown,
satisfiable: BOOL,
nameCount: NAT
];
Unknown: TYPE = REF UnknownRec;
UnknownRec:
TYPE =
RECORD [r, name:
REF, i, j:
INT, restricted:
BOOL];
r is for the client: LinearSolver does not use it.
name is for the client: LinearSolver does not use it.
i is the index of the row owned by this variable, or -1 if it doesn't own a row
j is the index of the column owned by this variable, or -1 if it doesn't own a column
NewUnknown: PROCEDURE [t: Tableau, name: REF] RETURNS [Unknown];
The client of LinearSolver creates (with NewUnknown) objects of type LinearSolver.Unknown, which are references to records. A LinearSolver.Unknown represents an unknown to be solved for. Clients should ignore the components of
Push, Pop: PROCEDURE;
C: PROCEDURE [x: LinearMonomial] RETURNS [r: REF];
AssertZero:
PROCEDURE [t: Tableau, c: Constraint, propEQ: PropEQProc ←
NIL];
Asserts that the linear monomial c = 0
AssertZero2:
PROCEDURE [t: Tableau, c1:
REAL, v1: Unknown, c2:
REAL, v2: Unknown, propEQ: PropEQProc ←
NIL];
Asserts that c1*v1 + c2*v2 = 0
AssertZero3:
PROCEDURE [t: Tableau, c1:
REAL, v1: Unknown, c2:
REAL, v2: Unknown, c3:
REAL, v3: Unknown, propEQ: PropEQProc ←
NIL];
Asserts that c1*v1 + c2*v2 + c3*v3 = 0
AssertGEZero:
PROCEDURE [t: Tableau, c: Constraint, propEQ: PropEQProc ←
NIL];
Asserts that the linear monomial c >= 0
AssertGEZero2:
PROCEDURE [t: Tableau, c1:
REAL, v1: Unknown, c2:
REAL, v2: Unknown, propEQ: PropEQProc ←
NIL];
Asserts that c1*v1 + c2*v2 >= 0
AssertGEZero3:
PROCEDURE [t: Tableau, c1:
REAL, v1: Unknown, c2:
REAL, v2: Unknown, c3:
REAL, v3: Unknown, propEQ: PropEQProc ←
NIL];
Asserts that c1*v1 + c2*v2 + c3*v3 >= 0
Restrict:
PROCEDURE [t: Tableau, v: Unknown, propEQ: PropEQProc ←
NIL];
Asserts that the variable v is >= 0
PropEQProc: TYPE = PROCEDURE [t: Tableau, u, v: Unknown];
Maximize:
PROCEDURE [t: Tableau, v: Unknown]
RETURNS [unbounded:
BOOLEAN];
Arranges the tableau so that the solution maximizes variable v
Satisfiable: PROCEDURE [t: Tableau] RETURNS [BOOL];
Solution: PROCEDURE [t: Tableau, v: Unknown] -- to be called only if Satisfiable[] -- RETURNS [REAL];
The module represents a stack of constraints, initially empty. The client pushes a constraint onto the stack with Assert. The effect of Push is to mark the current depth of the stack; the effect of Pop is to remove constraints from the stack until the last marked depth is restored. Satisfiable[] returns true if the current constraints can be satisfied. In this case, a solution to the constraints can be obtained by calling Solution[v] for each Unknown v. Note that this solution may not be unique. If you are interested in getting some particular solution instead of any solution, then you have to read about pivoting, below.
Constraint: TYPE = LIST OF REF -- to LinearMonomial --;
LinearMonomial: TYPE = RECORD [coefficient: REAL, unknown: Unknown];
"unity" is an Unknown representing one unit. It exists before any Unknowns are created. Solution[unity] = 1 invariantly. The constraint represented by a list of linear monomials is that the some of them is zero. Thus the constraint x+x=4 would be asserted Assert[[[2, x], [-4, unity]]], which asserts that 2*x+(-4)*unity = 0, or that x+x=4.
Example. A madmen obsessed by a fear that 2 + 2 = 4 will fail to be true in the near future decides to write a program to keep checking and alert him as soon as it fails. Since efficiency is not crucial, he chooses not to compute 2 + 2 and compare it with 4, but to solve x + x = 4 and compare the solution value of x with 2. He could do so with the following code.
-- Madman's code.
-- BEGIN OPEN LinearSolver;
-- x: Unknown;
-- WHILE TRUE DO
-- Push[];
-- x := NewUnknown[];
-- Assert[LIST[[2, x], [-4, unity]]];
-- IF NOT Satisfiable[] OR Solution[x] # 2 THEN SIGNAL Madman;
-- Pop[]
-- ENDLOOP
-- END.
Init: PROCEDURE RETURNS [t: Tableau]; -- Client must call it once before using the solver.
END. -- of LinearSolver.mesa