Linear Elimination

I have an undertermined system 
	A.x=b
The deal is to emulate the output from linear solve using basic linear algebra.

Example

We are going to convert a constrained minimization problem (with linear constraints) into an unconstrained problem over a lower dimensional space and see how well things work.

Here is our little example from class.  It looks like A.x=b where

A={
{1.0,   1.0, 1.0},
{1.0,-3.0,  1.0}
};
b={1.0,4.0};

This is two equations (the rows are linearly independent) in three variables (there are three columns).  Two of the columns are Linearly independent and following instructions in the book we have permuted the variables to get them to the front of the matrix.

The front 2×2 of the matrix is invertible because of our assumptions and the book wants to call it B and the rest of the matrix N.  I can not use N in Mathematica so I am going to write BMat and NMat

BMat=A[[1;;2,1;;2]];
NMat=A[[1;;2,3;;3]];
Inverse[BMat].A
Inverse[BMat].b

And we can assemble this using the instructions in the text to give an expression for the solution space for the constraints.   Another way is to use the symbolics in Mathematica

{x1,x2,x3}/.Solve[A.{x1,x2,x3}==b,{x1,x2,x3}]

Here is another way which is actually better behaved for large problems.  In some sense it trys to pick the best values and scalings

NS=NullSpace[A];
xP=LinearSolve[A,b];
xP+Array[α,Length[NS]].NS

On our little problem all of them will work equally well.

Using the Example

Quad[{x1_,x2_,x3_}]= (1/2)(((x1+x2-x3))^2)+(1/2)(((x3-2 x2+x1))^2)+(1/2)(((x1+3 x2-3x3))^2) ;
FindMinimum[{
Quad[{x1,x2,x3}],
{
x1+x2+x3==1.0,
x1+-3x2+x3==4.0
}
},
{x1,x2,x3}
]

Version 1

Elim1[x3_]={x1,x2,x3}/.Solve[{
x1+x2+x3==1.0,
x1+-3x2+x3==4.0
}][[1]]

FindMinimum[
Quad[Elim1[x3]],
{x3}
]

Version 2

A={
{1.0,  1.0,1.0},
{1.0,-3.0, 1.0}
};
b={1.0,4.0};
xp=LinearSolve[A,b];
ns=NullSpace[A];
Elim2[α_]=xp+α ns[[1]]

FindMinimum[
Quad[Elim2[α]],
{α}
]

We are going to convert a constrained minimization problem (with linear constraints) into an unconstrained problem over a lower dimensional space and see how well things work.

Here is our little example from class.  It looks like A.x=b where

A={
{1.0,   1.0, 1.0},
{1.0,-3.0,  1.0}
};
b={1.0,4.0};

This is two equations (the rows are linearly independent) in three variables (there are three columns).  Two of the columns are Linearly independent and following instructions in the book we have permuted the variables to get them to the front of the matrix.

The front 2×2 of the matrix is invertible because of our assumptions and the book wants to call it B and the rest of the matrix N.  I can not use N in Mathematica so I am going to write BMat and NMat

BMat=A[[1;;2,1;;2]];
NMat=A[[1;;2,3;;3]];
Inverse[BMat].A
Inverse[BMat].b

And we can assemble this using the instructions in the text to give an expression for the solution space for the constraints.   Another way is to use the symbolics in Mathematica

{x1,x2,x3}/.Solve[A.{x1,x2,x3}==b,{x1,x2,x3}]

Here is another way which is actually better behaved for large problems.  In some sense it trys to pick the best values and scalings

NS=NullSpace[A];
xP=LinearSolve[A,b];
xP+Array[α,Length[NS]].NS

On our little problem all of them will work equally well.