Given Y is nxk A is mxn B is mxk To solve AY=B for Y From http://math.stackexchange.com/questions/22616/solving-non-square-matrix-equations We only need to gaussian eliminate the augmented matrix (A|B) But I am at lost in finding a good way to solve XA=B where X is mxn A is nxk B is mxk How to solve this in a systematic way like solving AY=B??? Background of this question: THIS IS NOT A HOMEWORK QUESTION It concerns about the proof of the variation of parameters method http://tutorial.math.lamar.edu/Classes/DE/VariationofParameters.aspx When I try to do the proof myself, I always end up having second derivative terms of A" and B" i.e. Start with yp=Ay1+By2 (where A,B,y1,y2,p,q,r are functions of x and y1,y2 solves the homgeonous equation y"+py'+qy=r) yp'=Ay1'+A'y1+By2'+B'y2 yp"=Ay1"+2A'y1'+A"y1+By2"+2B'y2'+B"y2 I am also confused about the assumption A'y1+B'y2=0 http://www.physicsforums.com/showthread.php?t=646671 and this, being the most relevant of the google results, does not adequately address it So if I don't put in a single assumption (other than y1 and y2 will solve the corresponding homogenous Diff eqt hence theri corredongin term vanishes), this is what I end up View attachment 6563 Which is basically a matrix equation of the form XA=B Since these are not square matrices, I cannot invert them in the usual sense (by using determinants) I then attempt to solve this by trying to find a matrix X mutlipleid both sides on the left hand of the attached equation hoping to reduce the 2x3 matrix into I3 However I don't know how to do this systematically P.S. As a test of the idea, I tried to apply it to dot products of vectors in R2, the result (by actually checking the 4 systems of equations that results) is that they are not invertable even with this idea. I then tried to check the 3x3 case but solving 9 equations in 6 unknowns is too tedious and I am not sure whether I can interpret the gaussian elimination correctly...)
A nice example: \( \begin{pmatrix} ? & ? \\ ? & ? \\ ? & ? \end{pmatrix} \begin{pmatrix} -180400 & 596397 & 73596 \\ -151575 & 74496 & -506972 \end{pmatrix} = \begin{pmatrix} 274325 & 372909 & 1594512 \\ 360800 & -1192794 & -147192 \\ 36275 & 2013108 & 2829244 \end{pmatrix} \) Because \(A = \begin{pmatrix} -180400 & 596397 & 73596 \\ -151575 & 74496 & -506972 \end{pmatrix} = \begin{pmatrix} 7 & 24 \\ -24 & 7 \end{pmatrix} \begin{pmatrix} 19 & 0 & 0 \\ 0 & 23 & 0 \end{pmatrix} \begin{pmatrix} 200 & 201 & 1068 \\ -375 & 1032 & -124 \\ -1020 & -340 & 255 \end{pmatrix}\) it follows that A has a pseudo-inverse given by: \(\tilde{A} = \begin{pmatrix} 200 & 201 & 1068 \\ -375 & 1032 & -124 \\ -1020 & -340 & 255 \end{pmatrix}^{-1} \begin{pmatrix} \frac{1}{19} & 0 \\ 0 & \frac{1}{23} \\ 0 & 0 \end{pmatrix} \begin{pmatrix} 7 & 24 \\ -24 & 7 \end{pmatrix}^{-1} = \left( \frac{1}{1221025} \begin{pmatrix} 200 & 201 & 1068 \\ -375 & 1032 & -124 \\ -1020 & -340 & 255 \end{pmatrix}^{T} \right) \begin{pmatrix} \frac{1}{19} & 0 \\ 0 & \frac{1}{23} \\ 0 & 0 \end{pmatrix} \left( \frac{1}{625} \begin{pmatrix} 7 & 24 \\ -24 & 7 \end{pmatrix}^{T} \right) \\ \quad \quad \quad \quad = \left( \frac{1}{1221025} \begin{pmatrix} 200 & -375 & -1020 \\ 201 & 1032 & -340 \\ 1068 & -124 & 255 \end{pmatrix} \right) \begin{pmatrix} \frac{1}{19} & 0 \\ 0 & \frac{1}{23} \\ 0 & 0 \end{pmatrix} \left( \frac{1}{625} \begin{pmatrix} 7 & -24 \\ 24 & 7 \end{pmatrix} \right) = \frac{1}{333492453125} \begin{pmatrix} -138800, & -160275 \\ 502953 & 26304 \\ 115404 & -606028 \end{pmatrix}\) Thus \(A \tilde{A} = \frac{1}{333492453125} \begin{pmatrix} -180400 & 596397 & 73596 \\ -151575 & 74496 & -506972 \end{pmatrix} . \begin{pmatrix}-138800 & -160275 \\ 502953 & 26304 \\ 115404 & -606028 \end{pmatrix} = \begin{pmatrix} 1 & 0 \\ 0 & 1\end{pmatrix}\) for the win. But \(\tilde{A} A = \frac{1}{333492453125} \begin{pmatrix} -138800 & -160275 \\ 502953 & 26304 \\ 115404 & -606028 \end{pmatrix} . \begin{pmatrix}-180400 & 596397 & 73596 \\ -151575 & 74496 & -506972 \end{pmatrix} = \frac{1}{169} \begin{pmatrix}25 & -48 & 36 \\ -48 & 153 & 12 \\ 36 & 12 & 160\end{pmatrix} \) The latter is disappointing but what are you going to do. You can't recover a rank 3 matrix from a series of multiplications which include rank 2 matrices. In a certain sense, it is the best we can do. But in \(X A = B\) we multiply both sides on the right by \(\tilde{A}\) and get: \( X = X A \tilde{A} = B \tilde{A}\) and for this choice of \(A\) this is exact. \( X = \frac{1}{333492453125} \begin{pmatrix} 274325 & 372909 & 1594512 \\ 360800 & -1192794 & -147192 \\ 36275 & 2013108 & 2829244 \end{pmatrix} \begin{pmatrix} -138800 & -160275 \\ 502953 & 26304 \\ 115404 & -606028 \end{pmatrix}\) which I promise is nice. The route to getting a pseudo-inverse (or at least a numeric approximation of it) can be gotten from the singular value decomposition of a square or rectangular matrix which decomposes into a triple product analogous to my triple product which is equal to A above.
The above may or may not have and exact solution but the Penrose pseudoinverse always provides the minimum norm solution. In other words, \(||BA^{pseudo}A-B||=min\)
Bingo. I will point out that when doing this for (if memory serves) rank deficient problems you will still not be able to get a good answer. For example, it is often used to do linear best fit where the noise is presumed to be Gaussian. It is not unheard of to find that running the same best fit procedure on subsets of the data leads to wildly different answers, all of which seem to do the job quite well. A wall of avoiding this fluctuation in solutions is to regularise the problem by adding something like \(\lambda \mathbb{I}_{n}\) to the matrix before inverting, where \(\lambda\) is a small positive number. It basically just nudges the matrix to be inverted away from the problematic configuration. If you make \(\lambda\) huge then you deform the answer too much, but for small \(\lambda\) you can get an answer extremely close to what you expected (ie suppose you made the data yourself so you know the answer, as a means of debugging your algorithm) but without having problems of instability. Once took me a day of algorithm deconstruction and "print(EVERYTHING)" commands to realise that was the problem in one of my algorithms....
You do not need to resort to this trick if you use Finsim Math. Finsim Math takes care of the instability automatically, without any external intervention. Yes, I can believe that. There is an excellent book, by Bingulac, tha shows all these algorithms in detail. He has his own algorithm in computing the pseudo-inverse, it is one of the best in terms of speed and stability.
