Take the 2-minute tour ×
Mathematica Stack Exchange is a question and answer site for users of Mathematica. It's 100% free, no registration required.

Is there a simpler way of performing Gaussian Elimination other than using RowReduce? Such as a single built in function?

Edit:
Look at the example from our simulation class. Not too difficult, but using this method to solve problems of the sort is new to most of us. We are solving for P# of course.

Example

Also, to those asking... While I can see what is going on below, I really don't understand what it all means. Asking why not RowReduce? I guess not everyone is at that level of use yet, and I don't like just cutting and pasting internet code without understanding it. I simply wondered if there was a function that would do what the code did, but be built in.

  GaussianElimination[m_?MatrixQ, v_?VectorQ] :=
    Last /@ RowReduce[Flatten /@ Transpose[{m, v}]]
share|improve this question
5  
What's wrong with RowReduce? How is it not satisfactory and what do you mean by "simpler" solution? –  rm -rf Feb 20 '12 at 7:36
1  
Do you mean LinearSolve? –  user21 Feb 20 '12 at 8:31
2  
If you just need to solve an equation, use LinearSolve, no need for implementing a specific method. If you need the specific method, then what R.M. said. –  Szabolcs Feb 20 '12 at 8:31
    
Simpler, like a function versus having to device my own. I will look at LinearSolve, but might just use RowReduce. –  FossilizedCarlos Feb 20 '12 at 8:42
1  
Could you fill in some of the mathematical details in you question, then we could better direct you to the function you should use? –  rcollyer Feb 21 '12 at 2:26
show 2 more comments

1 Answer

up vote 14 down vote accepted

Based upon your update, you are trying to solve the system

$$\mathbf{A}\vec{x} = \vec{b}$$

for $\vec{x}$, so LinearSolve is exactly what you want. Also, it has the exact form

LinearSolve[A, b]

that you're asking for. Internally it uses a form of Gaussian elimination to solve such systems; this is most likely a variant of LU decomposition, but other methods are available. If you have more than one $\vec{b}$, you can use the form

solv = LinearSolve[A]

which returns a LinearSolveFunction which you can apply to each $\vec{x}$ in turn via

solv[b]

Edit: In the case of your example, RowReduce will return the identity matrix as your matrix is invertible (non-singular), so it would not be immediately useful. You could make it "useful" and create an augmented matrix, via

 augA = ArrayFlatten[{{#, IdentityMatrix[Length@#]}}]& @ A

which creates $$\left(\mathbf{A}\, |\, \mathbf{I} \right).$$ Then,

 redAugA = RowReduce[augA]

gives a matrix of the form $$\left(\mathbf{I}\, |\, \mathbf{A}^{-1} \right),$$ and the inverse is extractable via

redAugA[[All, Length@A + 1 ;; ]]

which uses the shorthand form of Part and Span to extract only the columns you want. But, if your going to go to the trouble of getting the inverse, you might as well use Inverse[A] directly.

However, if your matrix is singular, i.e. MatrixRank[A] < Length[A], then you need to use LeastSquares which returns the vector, $\vec{x}$, that minimizes $\lVert\mathbf{A}\vec{x} - \vec{b}\rVert_2$ where $\lVert\cdot\rVert_2$ refers to the standard Euclidean norm. Which has the same calling convention

LeastSquares[A, b]

but it lacks the pre-calculation capabilities of LinearSolve. If you need those, then you would first decompose the matrix using QRDecomposition and then LinearSolve is used, as follows

{q,r} = QRDecomposition[A];
LinearSolve[r, q.b]

Or, if you want a single function that operates like the second form of LinearSolve but with the least squares minimization,

savedLeastSquares[m_?MatrixQ]:= 
 Module[{q,r}, 
  {q,r} = QRDecomposition[m];
  LinearSolve[r, q.#]&
 ]
share|improve this answer
    
One could also use SingularValueDecomposition[] instead of QRDecomposition[]... –  J. M. May 4 '12 at 5:17
    
Your same idea of augmentation can be used to solve for A.x=b: instead of augmenting by an identity matrix, put b in a new (rightmost) column. Example: In[50]:= A = RandomInteger[{-10, 10}, {4, 4}]; b = RandomInteger[{-10, 10}, 4]; LinearSolve[A, b] - RowReduce[ArrayFlatten[{{#, Transpose[{b}]}}] &@A][[All, -1]] Out[52]= {0, 0, 0, 0} –  Daniel Lichtblau May 31 '13 at 16:59
    
@DanielLichtblau absolutely, where do you think I got the idea? –  rcollyer May 31 '13 at 17:07
add comment

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

Not the answer you're looking for? Browse other questions tagged or ask your own question.