REDUCE 3.6, 15-Jul-95, patched to 6 Mar 96 ...
if lisp !*rounded then rounded_was_on := t
else rounded_was_on := nil;
mat1 := mat((1,2,3,4,5),(2,3,4,5,6),(3,4,5,6,7),(4,5,6,7,8),(5,6,7,8,9));
[1 2 3 4 5]
[ ]
[2 3 4 5 6]
[ ]
mat1 := [3 4 5 6 7]
[ ]
[4 5 6 7 8]
[ ]
[5 6 7 8 9]
mat2 := mat((1,1,1,1),(2,2,2,2),(3,3,3,3),(4,4,4,4));
[1 1 1 1]
[ ]
[2 2 2 2]
mat2 := [ ]
[3 3 3 3]
[ ]
[4 4 4 4]
mat3 := mat((x),(x),(x),(x));
[x]
[ ]
[x]
mat3 := [ ]
[x]
[ ]
[x]
mat4 := mat((3,3),(4,4),(5,5),(6,6));
[3 3]
[ ]
[4 4]
mat4 := [ ]
[5 5]
[ ]
[6 6]
mat5 := mat((1,2,1,1),(1,2,3,1),(4,5,1,2),(3,4,5,6));
[1 2 1 1]
[ ]
[1 2 3 1]
mat5 := [ ]
[4 5 1 2]
[ ]
[3 4 5 6]
mat6 := mat((i+1,i+2,i+3),(4,5,2),(1,i,0));
[i + 1 i + 2 i + 3]
[ ]
mat6 := [ 4 5 2 ]
[ ]
[ 1 i 0 ]
mat7 := mat((1,1,0),(1,3,1),(0,1,1));
[1 1 0]
[ ]
mat7 := [1 3 1]
[ ]
[0 1 1]
mat8 := mat((1,3),(-4,3));
[1 3]
mat8 := [ ]
[-4 3]
mat9 := mat((1,2,3,4),(9,8,7,6));
[1 2 3 4]
mat9 := [ ]
[9 8 7 6]
poly := x^7+x^5+4*x^4+5*x^3+12;
7 5 4 3
poly := x + x + 4*x + 5*x + 12
poly1 := x^2+x*y^3+x*y*z^3+y*x+2+y*3;
2 3 3
poly1 := x + x*y + x*y*z + x*y + 3*y + 2
on errcont;
% Basis matrix manipulations.
add_columns(mat1,1,2,5*y);
[1 5*y + 2 3 4 5]
[ ]
[2 10*y + 3 4 5 6]
[ ]
[3 15*y + 4 5 6 7]
[ ]
[4 5*(4*y + 1) 6 7 8]
[ ]
[5 25*y + 6 7 8 9]
add_rows(mat1,1,2,x);
[ 1 2 3 4 5 ]
[ ]
[x + 2 2*x + 3 3*x + 4 4*x + 5 5*x + 6]
[ ]
[ 3 4 5 6 7 ]
[ ]
[ 4 5 6 7 8 ]
[ ]
[ 5 6 7 8 9 ]
add_to_columns(mat1,3,1000);
[1 2 1003 4 5]
[ ]
[2 3 1004 5 6]
[ ]
[3 4 1005 6 7]
[ ]
[4 5 1006 7 8]
[ ]
[5 6 1007 8 9]
add_to_columns(mat1,{1,2,3},y);
[y + 1 y + 2 y + 3 4 5]
[ ]
[y + 2 y + 3 y + 4 5 6]
[ ]
[y + 3 y + 4 y + 5 6 7]
[ ]
[y + 4 y + 5 y + 6 7 8]
[ ]
[y + 5 y + 6 y + 7 8 9]
add_to_rows(mat1,2,1000);
[ 1 2 3 4 5 ]
[ ]
[1002 1003 1004 1005 1006]
[ ]
[ 3 4 5 6 7 ]
[ ]
[ 4 5 6 7 8 ]
[ ]
[ 5 6 7 8 9 ]
add_to_rows(mat1,{1,2,3},x);
[x + 1 x + 2 x + 3 x + 4 x + 5]
[ ]
[x + 2 x + 3 x + 4 x + 5 x + 6]
[ ]
[x + 3 x + 4 x + 5 x + 6 x + 7]
[ ]
[ 4 5 6 7 8 ]
[ ]
[ 5 6 7 8 9 ]
augment_columns(mat1,2);
[2]
[ ]
[3]
[ ]
[4]
[ ]
[5]
[ ]
[6]
augment_columns(mat1,{1,2,5});
[1 2 5]
[ ]
[2 3 6]
[ ]
[3 4 7]
[ ]
[4 5 8]
[ ]
[5 6 9]
stack_rows(mat1,3);
[3 4 5 6 7]
stack_rows(mat1,{1,3,5});
[1 2 3 4 5]
[ ]
[3 4 5 6 7]
[ ]
[5 6 7 8 9]
char_poly(mat1,x);
3 2
x *(x - 25*x - 50)
column_dim(mat2);
4
row_dim(mat1);
5
copy_into(mat7,mat1,2,3);
[1 2 3 4 5]
[ ]
[2 3 1 1 0]
[ ]
[3 4 1 3 1]
[ ]
[4 5 0 1 1]
[ ]
[5 6 7 8 9]
copy_into(mat7,mat1,5,5);
***** Error in copy_into: the matrix
[1 1 0]
[ ]
[1 3 1]
[ ]
[0 1 1]
does not fit into
[1 2 3 4 5]
[ ]
[2 3 4 5 6]
[ ]
[3 4 5 6 7]
[ ]
[4 5 6 7 8]
[ ]
[5 6 7 8 9]
at position 5,5.
diagonal(3);
[3]
% diagonal can take both a list of arguments or just the arguments.
diagonal({mat2,mat6});
[1 1 1 1 0 0 0 ]
[ ]
[2 2 2 2 0 0 0 ]
[ ]
[3 3 3 3 0 0 0 ]
[ ]
[4 4 4 4 0 0 0 ]
[ ]
[0 0 0 0 i + 1 i + 2 i + 3]
[ ]
[0 0 0 0 4 5 2 ]
[ ]
[0 0 0 0 1 i 0 ]
diagonal(mat1,mat2,mat5);
[1 2 3 4 5 0 0 0 0 0 0 0 0]
[ ]
[2 3 4 5 6 0 0 0 0 0 0 0 0]
[ ]
[3 4 5 6 7 0 0 0 0 0 0 0 0]
[ ]
[4 5 6 7 8 0 0 0 0 0 0 0 0]
[ ]
[5 6 7 8 9 0 0 0 0 0 0 0 0]
[ ]
[0 0 0 0 0 1 1 1 1 0 0 0 0]
[ ]
[0 0 0 0 0 2 2 2 2 0 0 0 0]
[ ]
[0 0 0 0 0 3 3 3 3 0 0 0 0]
[ ]
[0 0 0 0 0 4 4 4 4 0 0 0 0]
[ ]
[0 0 0 0 0 0 0 0 0 1 2 1 1]
[ ]
[0 0 0 0 0 0 0 0 0 1 2 3 1]
[ ]
[0 0 0 0 0 0 0 0 0 4 5 1 2]
[ ]
[0 0 0 0 0 0 0 0 0 3 4 5 6]
extend(mat1,3,2,x);
[1 2 3 4 5 x x]
[ ]
[2 3 4 5 6 x x]
[ ]
[3 4 5 6 7 x x]
[ ]
[4 5 6 7 8 x x]
[ ]
[5 6 7 8 9 x x]
[ ]
[x x x x x x x]
[ ]
[x x x x x x x]
[ ]
[x x x x x x x]
find_companion(mat5,x);
2
x - 2*x - 2
get_columns(mat1,1);
{
[1]
[ ]
[2]
[ ]
[3]
[ ]
[4]
[ ]
[5]
}
get_columns(mat1,{1,2});
{
[1]
[ ]
[2]
[ ]
[3]
[ ]
[4]
[ ]
[5]
,
[2]
[ ]
[3]
[ ]
[4]
[ ]
[5]
[ ]
[6]
}
get_rows(mat1,3);
{
[3 4 5 6 7]
}
get_rows(mat1,{1,3});
{
[1 2 3 4 5]
,
[3 4 5 6 7]
}
hermitian_tp(mat6);
[ - i + 1 4 1 ]
[ ]
[ - i + 2 5 - i]
[ ]
[ - i + 3 2 0 ]
% matrix_augment and matrix_stack can take both a list of arguments
% or just the arguments.
matrix_augment({mat1,mat2});
***** Error in matrix_augment:
***** all input matrices must have the same row dimension.
matrix_augment(mat4,mat2,mat4);
[3 3 1 1 1 1 3 3]
[ ]
[4 4 2 2 2 2 4 4]
[ ]
[5 5 3 3 3 3 5 5]
[ ]
[6 6 4 4 4 4 6 6]
matrix_stack(mat1,mat2);
***** Error in matrix_stack:
***** all input matrices must have the same column dimension.
matrix_stack({mat6,mat((z,z,z)),mat7});
[i + 1 i + 2 i + 3]
[ ]
[ 4 5 2 ]
[ ]
[ 1 i 0 ]
[ ]
[ z z z ]
[ ]
[ 1 1 0 ]
[ ]
[ 1 3 1 ]
[ ]
[ 0 1 1 ]
minor(mat1,2,3);
[1 2 4 5]
[ ]
[3 4 6 7]
[ ]
[4 5 7 8]
[ ]
[5 6 8 9]
mult_columns(mat1,3,y);
[1 2 3*y 4 5]
[ ]
[2 3 4*y 5 6]
[ ]
[3 4 5*y 6 7]
[ ]
[4 5 6*y 7 8]
[ ]
[5 6 7*y 8 9]
mult_columns(mat1,{2,3,4},100);
[1 200 300 400 5]
[ ]
[2 300 400 500 6]
[ ]
[3 400 500 600 7]
[ ]
[4 500 600 700 8]
[ ]
[5 600 700 800 9]
mult_rows(mat1,2,x);
[ 1 2 3 4 5 ]
[ ]
[2*x 3*x 4*x 5*x 6*x]
[ ]
[ 3 4 5 6 7 ]
[ ]
[ 4 5 6 7 8 ]
[ ]
[ 5 6 7 8 9 ]
mult_rows(mat1,{1,3,5},10);
[10 20 30 40 50]
[ ]
[2 3 4 5 6 ]
[ ]
[30 40 50 60 70]
[ ]
[4 5 6 7 8 ]
[ ]
[50 60 70 80 90]
pivot(mat1,3,3);
[ - 4 - 2 2 4 ]
[------ ------ 0 --- --- ]
[ 5 5 5 5 ]
[ ]
[ - 2 - 1 1 2 ]
[------ ------ 0 --- --- ]
[ 5 5 5 5 ]
[ ]
[ 3 4 5 6 7 ]
[ ]
[ 2 1 - 1 - 2 ]
[ --- --- 0 ------ ------]
[ 5 5 5 5 ]
[ ]
[ 4 2 - 2 - 4 ]
[ --- --- 0 ------ ------]
[ 5 5 5 5 ]
rows_pivot(mat1,3,3,{1,5});
[ - 4 - 2 2 4 ]
[------ ------ 0 --- --- ]
[ 5 5 5 5 ]
[ ]
[ 2 3 4 5 6 ]
[ ]
[ 3 4 5 6 7 ]
[ ]
[ 4 5 6 7 8 ]
[ ]
[ 4 2 - 2 - 4 ]
[ --- --- 0 ------ ------]
[ 5 5 5 5 ]
remove_columns(mat1,3);
[1 2 4 5]
[ ]
[2 3 5 6]
[ ]
[3 4 6 7]
[ ]
[4 5 7 8]
[ ]
[5 6 8 9]
remove_columns(mat1,{2,3,4});
[1 5]
[ ]
[2 6]
[ ]
[3 7]
[ ]
[4 8]
[ ]
[5 9]
remove_rows(mat1,2);
[1 2 3 4 5]
[ ]
[3 4 5 6 7]
[ ]
[4 5 6 7 8]
[ ]
[5 6 7 8 9]
remove_rows(mat1,{1,3});
[2 3 4 5 6]
[ ]
[4 5 6 7 8]
[ ]
[5 6 7 8 9]
remove_rows(mat1,{1,2,3,4,5});
***** Warning in remove_rows:
all the rows have been removed. Returning nil.
swap_columns(mat1,2,4);
[1 4 3 2 5]
[ ]
[2 5 4 3 6]
[ ]
[3 6 5 4 7]
[ ]
[4 7 6 5 8]
[ ]
[5 8 7 6 9]
swap_rows(mat1,1,2);
[2 3 4 5 6]
[ ]
[1 2 3 4 5]
[ ]
[3 4 5 6 7]
[ ]
[4 5 6 7 8]
[ ]
[5 6 7 8 9]
swap_entries(mat1,{1,1},{5,5});
[9 2 3 4 5]
[ ]
[2 3 4 5 6]
[ ]
[3 4 5 6 7]
[ ]
[4 5 6 7 8]
[ ]
[5 6 7 8 1]
% Constructors - functions that create matrices.
band_matrix(x,5);
[x 0 0 0 0]
[ ]
[0 x 0 0 0]
[ ]
[0 0 x 0 0]
[ ]
[0 0 0 x 0]
[ ]
[0 0 0 0 x]
band_matrix({x,y,z},6);
[y z 0 0 0 0]
[ ]
[x y z 0 0 0]
[ ]
[0 x y z 0 0]
[ ]
[0 0 x y z 0]
[ ]
[0 0 0 x y z]
[ ]
[0 0 0 0 x y]
block_matrix(1,2,{mat1,mat2});
***** Error in block_matrix: row dimensions of
***** matrices into block_matrix are not compatible
block_matrix(2,3,{mat2,mat3,mat2,mat3,mat2,mat2});
[1 1 1 1 x 1 1 1 1]
[ ]
[2 2 2 2 x 2 2 2 2]
[ ]
[3 3 3 3 x 3 3 3 3]
[ ]
[4 4 4 4 x 4 4 4 4]
[ ]
[x 1 1 1 1 1 1 1 1]
[ ]
[x 2 2 2 2 2 2 2 2]
[ ]
[x 3 3 3 3 3 3 3 3]
[ ]
[x 4 4 4 4 4 4 4 4]
char_matrix(mat1,x);
[x - 1 -2 -3 -4 -5 ]
[ ]
[ -2 x - 3 -4 -5 -6 ]
[ ]
[ -3 -4 x - 5 -6 -7 ]
[ ]
[ -4 -5 -6 x - 7 -8 ]
[ ]
[ -5 -6 -7 -8 x - 9]
cfmat := coeff_matrix({x+y+4*z=10,y+x-z=20,x+y+4});
cfmat := {
[4 1 1]
[ ]
[-1 1 1]
[ ]
[0 1 1]
,
[z]
[ ]
[y]
[ ]
[x]
,
[10]
[ ]
[20]
[ ]
[-4]
}
first cfmat * second cfmat;
[x + y + 4*z]
[ ]
[ x + y - z ]
[ ]
[ x + y ]
third cfmat;
[10]
[ ]
[20]
[ ]
[-4]
companion(poly,x);
[0 0 0 0 0 0 -12]
[ ]
[1 0 0 0 0 0 0 ]
[ ]
[0 1 0 0 0 0 0 ]
[ ]
[0 0 1 0 0 0 -5 ]
[ ]
[0 0 0 1 0 0 -4 ]
[ ]
[0 0 0 0 1 0 -1 ]
[ ]
[0 0 0 0 0 1 0 ]
hessian(poly1,{w,x,y,z});
[0 0 0 0 ]
[ ]
[ 2 3 2 ]
[0 2 3*y + z + 1 3*y*z ]
[ ]
[ 2 3 2 ]
[0 3*y + z + 1 6*x*y 3*x*z ]
[ ]
[ 2 2 ]
[0 3*y*z 3*x*z 6*x*y*z]
hilbert(4,1);
[ 1 1 1 ]
[ 1 --- --- ---]
[ 2 3 4 ]
[ ]
[ 1 1 1 1 ]
[--- --- --- ---]
[ 2 3 4 5 ]
[ ]
[ 1 1 1 1 ]
[--- --- --- ---]
[ 3 4 5 6 ]
[ ]
[ 1 1 1 1 ]
[--- --- --- ---]
[ 4 5 6 7 ]
hilbert(3,y+x);
[ - 1 - 1 - 1 ]
[----------- ----------- -----------]
[ x + y - 2 x + y - 3 x + y - 4 ]
[ ]
[ - 1 - 1 - 1 ]
[----------- ----------- -----------]
[ x + y - 3 x + y - 4 x + y - 5 ]
[ ]
[ - 1 - 1 - 1 ]
[----------- ----------- -----------]
[ x + y - 4 x + y - 5 x + y - 6 ]
jacobian({x^4,x*y^2,x*y*z^3},{w,x,y,z});
[ 3 ]
[0 4*x 0 0 ]
[ ]
[ 2 ]
[0 y 2*x*y 0 ]
[ ]
[ 3 3 2]
[0 y*z x*z 3*x*y*z ]
jordan_block(x,5);
[x 1 0 0 0]
[ ]
[0 x 1 0 0]
[ ]
[0 0 x 1 0]
[ ]
[0 0 0 x 1]
[ ]
[0 0 0 0 x]
make_identity(11);
[1 0 0 0 0 0 0 0 0 0 0]
[ ]
[0 1 0 0 0 0 0 0 0 0 0]
[ ]
[0 0 1 0 0 0 0 0 0 0 0]
[ ]
[0 0 0 1 0 0 0 0 0 0 0]
[ ]
[0 0 0 0 1 0 0 0 0 0 0]
[ ]
[0 0 0 0 0 1 0 0 0 0 0]
[ ]
[0 0 0 0 0 0 1 0 0 0 0]
[ ]
[0 0 0 0 0 0 0 1 0 0 0]
[ ]
[0 0 0 0 0 0 0 0 1 0 0]
[ ]
[0 0 0 0 0 0 0 0 0 1 0]
[ ]
[0 0 0 0 0 0 0 0 0 0 1]
on rounded;
% makes things a bit easier to read.
random_matrix(3,3,100);
[ - 28.8708957912 - 83.4775707772 - 22.5350939803]
[ ]
[ - 83.4775707772 35.4030294062 5.76473742166 ]
[ ]
[ - 22.5350939803 5.76473742166 74.9835325987 ]
on not_negative;
random_matrix(3,3,100);
[52.1649704637 63.7238994877 0.886186275566]
[ ]
[63.7238994877 18.6365765396 77.194379299 ]
[ ]
[0.886186275566 77.194379299 4.09332974714 ]
on only_integer;
random_matrix(3,3,100);
[49 42 67]
[ ]
[42 69 84]
[ ]
[67 84 31]
on symmetric;
random_matrix(3,3,100);
[73 63 85]
[ ]
[63 5 81]
[ ]
[85 81 35]
off symmetric;
on upper_matrix;
random_matrix(3,3,100);
[4 70 28]
[ ]
[0 7 58]
[ ]
[0 0 19]
off upper_matrix;
on lower_matrix;
random_matrix(3,3,100);
[7 0 0 ]
[ ]
[46 49 0 ]
[ ]
[42 70 65]
off lower_matrix;
on imaginary;
off not_negative;
random_matrix(3,3,100);
[ 38*i - 14 - 72*i + 19 51*i - 30 ]
[ ]
[ - 99*i + 72 94*i - 59 3*i - 46 ]
[ ]
[ 47*i - 54 - 9*i - 73*i - 28]
off rounded;
% toeplitz and vandermonde can take both a list of arguments or just
% the arguments.
toeplitz({1,2,3,4,5});
[1 2 3 4 5]
[ ]
[2 1 2 3 4]
[ ]
[3 2 1 2 3]
[ ]
[4 3 2 1 2]
[ ]
[5 4 3 2 1]
toeplitz(x,y,z);
[x y z]
[ ]
[y x y]
[ ]
[z y x]
vandermonde({1,2,3,4,5});
[1 1 1 1 1 ]
[ ]
[1 2 4 8 16 ]
[ ]
[1 3 9 27 81 ]
[ ]
[1 4 16 64 256]
[ ]
[1 5 25 125 625]
vandermonde(x,y,z);
[ 2]
[1 x x ]
[ ]
[ 2]
[1 y y ]
[ ]
[ 2]
[1 z z ]
% kronecker_product
a1 := mat((1,2),(3,4),(5,6));
[1 2]
[ ]
a1 := [3 4]
[ ]
[5 6]
a2 := mat((1,x,1),(2,2,2),(3,3,3));
[1 x 1]
[ ]
a2 := [2 2 2]
[ ]
[3 3 3]
kronecker_product(a1,a2);
[1 x 1 2 2*x 2 ]
[ ]
[2 2 2 4 4 4 ]
[ ]
[3 3 3 6 6 6 ]
[ ]
[3 3*x 3 4 4*x 4 ]
[ ]
[6 6 6 8 8 8 ]
[ ]
[9 9 9 12 12 12]
[ ]
[5 5*x 5 6 6*x 6 ]
[ ]
[10 10 10 12 12 12]
[ ]
[15 15 15 18 18 18]
clear a1,a2;
% High level algorithms.
on rounded;
% makes output easier to read.
ch := cholesky(mat7);
ch := {
[1 0 0 ]
[ ]
[1 1.41421356237 0 ]
[ ]
[0 0.707106781187 0.707106781187]
,
[1 1 0 ]
[ ]
[0 1.41421356237 0.707106781187]
[ ]
[0 0 0.707106781187]
}
tp first ch - second ch;
[0 0 0]
[ ]
[0 0 0]
[ ]
[0 0 0]
tmp := first ch * second ch;
[1 1 0]
[ ]
tmp := [1 3.0 1]
[ ]
[0 1 1]
tmp - mat7;
[0 0 0]
[ ]
[0 0 0]
[ ]
[0 0 0]
off rounded;
gram_schmidt({1,0,0},{1,1,0},{1,1,1});
{{1,0,0},{0,1,0},{0,0,1}}
gram_schmidt({1,2},{3,4});
1 2 2*sqrt(5) - sqrt(5)
{{---------,---------},{-----------,------------}}
sqrt(5) sqrt(5) 5 5
on rounded;
% again, makes large quotients a bit more readable.
% The algorithm used for lu_decom sometimes swaps the rows of the input
% matrix so that (given matrix A, lu_decom(A) = {L,U,vec}), we find L*U
% does not equal A but a row equivalent of it. The call convert(A,vec)
% will return this row equivalent (ie: L*U = convert(A,vec)).
lu := lu_decom(mat5);
lu := {
[4 0 0 0 ]
[ ]
[1 0.75 0 0 ]
[ ]
[1 0.75 2.0 0 ]
[ ]
[3 0.25 4.0 4.33333333333]
,
[1 1.25 0.25 0.5 ]
[ ]
[0 1 1 0.666666666667]
[ ]
[0 0 1 0 ]
[ ]
[0 0 0 1 ]
,
[3,3,3,4]}
mat5;
[1 2 1 1]
[ ]
[1 2 3 1]
[ ]
[4 5 1 2]
[ ]
[3 4 5 6]
tmp := first lu * second lu;
[4 5.0 1 2.0]
[ ]
[1 2.0 1 1 ]
tmp := [ ]
[1 2.0 3.0 1 ]
[ ]
[3 4.0 5.0 6.0]
tmp1 := convert(mat5,third lu);
[4 5 1 2]
[ ]
[1 2 1 1]
tmp1 := [ ]
[1 2 3 1]
[ ]
[3 4 5 6]
tmp - tmp1;
[0 0 0 0]
[ ]
[0 0 0 0]
[ ]
[0 0 0 0]
[ ]
[0 0 0 0]
% and the complex case...
lu1 := lu_decom(mat6);
lu1 := {
[ 1 0 0 ]
[ ]
[ 4 - 4*i + 5 0 ]
[ ]
[i + 1 3 0.414634146341*i + 2.26829268293]
,
[1 i 0 ]
[ ]
[0 1 0.19512195122*i + 0.243902439024]
[ ]
[0 0 1 ]
,
[3,2,3]}
mat6;
[i + 1 i + 2 i + 3]
[ ]
[ 4 5 2 ]
[ ]
[ 1 i 0 ]
tmp := first lu1 * second lu1;
[ 1 i 0 ]
[ ]
tmp := [ 4 5 2.0 ]
[ ]
[i + 1 i + 2 i + 3.0]
tmp1 := convert(mat6,third lu1);
[ 1 i 0 ]
[ ]
tmp1 := [ 4 5 2 ]
[ ]
[i + 1 i + 2 i + 3]
tmp - tmp1;
[0 0 0]
[ ]
[0 0 0]
[ ]
[0 0 0]
mat9inv := pseudo_inverse(mat9);
[ - 0.199999999996 0.100000000013 ]
[ ]
[ - 0.0499999999988 0.0500000000037 ]
mat9inv := [ ]
[ 0.0999999999982 - 5.57816640101e-12]
[ ]
[ 0.249999999995 - 0.0500000000148 ]
mat9 * mat9inv;
[ 0.999999999982 - 0.0000000000557817125824]
[ ]
[5.54201129432e-12 1.00000000002 ]
simplex(min,2*x1+14*x2+36*x3,{-2*x1+x2+4*x3>=5,-x1-2*x2-3*x3<=2});
{45.0,{x1=0,x2=0,x3=1.25}}
simplex(max,10000 x1 + 1000 x2 + 100 x3 + 10 x4 + x5,{ x1 <= 1, 20 x1 +
x2 <= 100, 200 x1 + 20 x2 + x3 <= 10000, 2000 x1 + 200 x2 + 20 x3 + x4
<= 1000000, 20000 x1 + 2000 x2 + 200 x3 + 20 x4 + x5 <= 100000000});
{1.0e+8,{x1=0,x2=0,x3=0,x4=0,x5=100000000.0}}
simplex(max, 5 x1 + 4 x2 + 3 x3,
{ 2 x1 + 3 x2 + x3 <= 5,
4 x1 + x2 + 2 x3 <= 11,
3 x1 + 4 x2 + 2 x3 <= 8 });
{13.0,{x1=2.0,x2=0,x3=1.0}}
simplex(min,3 x1 + 5 x2,{ x1 + 2 x2 >= 2, 22 x1 + x2 >= 3});
{5.04651162791,{x1=0.093023255813953,x2=0.95348837209302}}
simplex(max,10x+5y+5.5z,{5x+3z<=200,0.2x+0.1y+0.5z<=12,0.1x+0.2y+0.3z<=9,
30x+10y+50z<=1500});
{525.0,{x=40.0,y=25.0,z=0}}
%example of extra variables (>=0) being added.
simplex(min,x-y,{x>=-3});
*** Warning: variable y not defined in input. Has been defined as >=0.
***** Error in simplex: The problem is unbounded.
% unfeasible as simplex algorithm implies all x>=0.
simplex(min,x,{x<=-100});
***** Error in simplex: Problem has no feasible solution.
% three error examples.
simplex(maxx,x,{x>=5});
***** Error in simplex(first argument): must be either max or min.
simplex(max,x,x>=5);
***** Error in simplex(third argument}: must be a list.
simplex(max,x,{x<=y});
***** Error in simplex(third argument):
***** right hand side of each inequality must be a number
simplex(max, 346 X11 + 346 X12 + 248 X21 + 248 X22 + 399 X31 + 399 X32 +
200 Y11 + 200 Y12 + 75 Y21 + 75 Y22 + 2.35 Z1 + 3.5 Z2,
{
4 X11 + 4 X12 + 2 X21 + 2 X22 + X31 + X32 + 250 Y11 + 250 Y12 + 125 Y21 +
125 Y22 <= 25000,
X11 + X12 + X21 + X22 + X31 + X32 + 2 Y11 + 2 Y12 + Y21 + Y22 <= 300,
20 X11 + 15 X12 + 30 Y11 + 20 Y21 + Z1 <= 1500,
40 X12 + 35 X22 + 50 X32 + 15 Y12 + 10 Y22 + Z2 = 5000,
X31 = 0,
Y11 + Y12 <= 50,
Y21 + Y22 <= 100
});
{99250.0,
{y21=0,
y22=0,
x31=0,
x11=75.0,
z1=0,
x21=225.0,
z2=5000.0,
x32=0,
x22=0,
x12=0,
y12=0,
y11=0}}
% from Marc van Dongen. Finding the first feasible solution for the
% solution of systems of linear diophantine inequalities.
simplex(max,0,{
3*X259+4*X261+3*X262+2*X263+X269+2*X270+3*X271+4*X272+5*X273+X229=2,
7*X259+11*X261+8*X262+5*X263+3*X269+6*X270+9*X271+12*X272+15*X273+X229=4,
2*X259+5*X261+4*X262+3*X263+3*X268+4*X269+5*X270+6*X271+7*X272+8*X273=1,
X262+2*X263+5*X268+4*X269+3*X270+2*X271+X272+2*X229=1,
X259+X262+2*X263+4*X268+3*X269+2*X270+X271-X273+3*X229=2,
X259+2*X261+2*X262+2*X263+3*X268+3*X269+3*X270+3*X271+3*X272+3*X273+X229=1,
X259+X261+X262+X263+X268+X269+X270+X271+X272+X273+X229=1});
{0,
{x229=0.5,
x259=0.5,
x261=0,
x262=0,
x263=0,
x268=0,
x269=0,
x270=0,
x271=0,
x272=0,
x273=0}}
svd_ans := svd(mat8);
svd_ans := {
[ 0.289784137735 0.957092029805]
[ ]
[ - 0.957092029805 0.289784137735]
,
[5.1491628629 0 ]
[ ]
[ 0 2.9130948854]
,
[ - 0.687215403194 0.726453707825 ]
[ ]
[ - 0.726453707825 - 0.687215403194]
}
tmp := tp first svd_ans * second svd_ans * third svd_ans;
[ 0.99999998509 2.9999999859 ]
tmp := [ ]
[ - 4.00000004924 2.99999995342]
tmp - mat8;
[ - 0.0000000149096008872 - 0.0000000141042817425]
[ ]
[ - 0.0000000492430629606 - 0.0000000465832750152]
mat9inv := pseudo_inverse(mat9);
[ - 0.199999999996 0.100000000013 ]
[ ]
[ - 0.0499999999988 0.0500000000037 ]
mat9inv := [ ]
[ 0.0999999999982 - 5.57816640101e-12]
[ ]
[ 0.249999999995 - 0.0500000000148 ]
mat9 * mat9inv;
[ 0.999999999982 - 0.0000000000557817125824]
[ ]
[5.54201129432e-12 1.00000000002 ]
% Predicates.
matrixp(mat1);
t
matrixp(poly);
squarep(mat2);
t
squarep(mat3);
symmetricp(mat1);
t
symmetricp(mat3);
if not rounded_was_on then off rounded;
END;
(TIME: linalg 3359 3469)