Re: A few questions about RowReduce - bis

*To*: mathgroup at smc.vnet.net*Subject*: [mg24256] Re: [mg24074] A few questions about RowReduce - bis*From*: Murray Eisenberg <murray at math.umass.edu>*Date*: Wed, 5 Jul 2000 23:10:37 -0400 (EDT)*Organization*: Mathematics & Statistics, Univ. of Mass./Amherst*References*: <8j9dsq$51v@smc.vnet.net>*Sender*: owner-wri-mathgroup at wolfram.com

A reason to find a row-echelon form that is not the canonical "reduced" row-echelon form, is that you can answer questions about the nature of solutions of linear systems, or linear independence of vectors, etc., without having to do -- or asking Mathematica to do -- all the extra work to get the reduced row-echelon form. Preston Nichols wrote: > > Jack-Michel, > > Here are some ideas which I hope will be helpful, though they do not > address your questions about a "non-reduced" or "partial" row-echelon form. > (To be honest, I can't think of why that specific form would be > essential.) My discussion focuses on equations for the image of a matrix. > > Here's your matrix: > > In[1]:= n = 5; > m = Table[(i + j - 1)^2, {i, n}, {j, n}]; > > A vector of variables (note I changed yi from you message to y[i]): > > In[2]:= v = Table[y[i], {i, n}]; > > The same structure can be built this way: > > In[3]:= v = Array[y,n]; > > The vector v appended to m as the last column: > > In[4]:= << "LinearAlgebra`MatrixManipulation`"; > mv = AppendRows[m, Transpose[{v}]] > > The same operation, without requiring the MatrixManipulation package: > > In[5]:= mv = Transpose[Append[Transpose[m], v]] > > RowReduce[mv] permits division by polynomials in the y-variables, even > though such expression may have the value zero. Row-reduction of mv, > treating the y-variables as generic, leads to (from your own Mathematica > function): > > {{1, 4, 9, 16, 25, y[1]}, > {0, -7, -20, -39, -64, -4y[1] + y[2]}, > {0, 0, 8/7, 24/7, 48/7, (17*y[1])/7 - (20*y[2])/7 + y[3]}, > {0, 0, 0, 0, 0, -y[1] + 3*y[2] - 3*y[3] + y[4]}, > {0, 0, 0, 0, 0, -3*y[1] + 8*y[2] - 6*y[3] + y[5]}} > > from which can be extracted equations for the image of the matrix: > > -y[1] + 3 y[2] - 3 y[3] + y[4]==0, > -3 y[1] + 8 y[2] - 6 y[3] + y[5]==0 > > One of your questions is how to get Mathematica to perform such a > calculation. I found three ways, each having its own interest. > > ==================================================== > > Mathematica solution 1: The "Left Nullspace" > > In[6]:= leftnullspace = NullSpace[Transpose[m]] > > Out[6]= {{-3, 8, -6, 0, 1}, {-1, 3, -3, 1, 0}} > > These vectors span the orthogonal complement of the image of m, so by > setting dot products with v to zero, you can recover the desired equations. > > In[7]:= Thread[leftnullspace.v == 0] > > Out[7]= {-3 y[1] + 8 y[2] - 6 y[3] + y[5] == 0, > -y[1] + 3 y[2] - 3 y[3] + y[4] == 0} > > (Thread changes an "equation of lists" into a "list of equations".) > > This strikes me as "the right way", because (1) it is brief, using > high-level built-in Mathematica functions, and (2) it relies on a > fundamental theorem of finite-dimensional linear algebra. An aphorism: > Sometimes Mathematica does mathematics for you; sometimes you must do > mathematics for Mathematica. > > ==================================================== > > Mathematica solution 2: Use Eliminate > > In[8]:= ImageEquations[mat_?MatrixQ, y_Symbol] := Module[{m, n, x}, > {m, n} = Dimensions[mat]; > Eliminate[Thread[mat.Array[x, n] == Array[y, m]], Array[x, n]] > ] > > In[9]:= ImageEquations[m, y] > > Out[9]= y[1] == 3 y[2] - 3 y[3] + y[4] && y[2] == 3 y[3] - 3 y[4] + y[5] > > Note that this method returns equations, rather than vectors (of > coefficients for the equations) as my other solutions do. > > To understand how this function works, examine the step-by-step output from > input lines such as > > u = Array[x, n] > v = Array[y, n] > > m.u == v > > Thread[%] (* "equation of lists" to "list of equations" again *) > > Eliminate[%, u] > > This may not always give exactly the same equations as the "left nullspace" > method, but the difference is due only to a change of basis for the > orthogonal space to the image of m. Mathematically, the answers are > equivalent. > > ==================================================== > > Mathematica solution 3: Use RowReduce after all! > > Here's a function which appends an IdentityMatrix of appropriate size to > the right of the given matrix, and then applies RowReduce to the whole > thing. (This is what one does when using row operations to find the > inverse of a square matrix, but it can be done, in principle, on any matrix.) > > Here's a version that makes use of LinearAlgebra`MatrixManipulation`: > > In[10]:= doubleRowReduce[mat_?MatrixQ] := > RowReduce[ > AppendRows[mat, IdentityMatrix[First[Dimensions[mat]]]] > ] > > Here's a version of the same function that doesn't require the package: > > In[11]:= doubleRowReduce[mat_?MatrixQ] := > RowReduce[ > MapThread[Join, {mat, IdentityMatrix[First[Dimensions[mat]]]}] > ] > > The coefficients for the equations of the image of the matrix m appear at > the right in the bottom rows of doubleRowReduce[m]. The rows to use are > those which reduce to zeroes in m. (I won't try to reproduce the output > here; try it and see!) You can think of the columns of the appended > identity matrix as holding coefficients of the y-variables. The original > problem required a cummulative "record" of the arithmetical operations > performed on each of the y[i]'s, so giving each y[i] a column of its own > makes sense. > > Here's a function to extract the desired coefficients from the results of > doubleRowReduce, using TakeMatrix from LinearAlgebra`MatrixManipulation. I > suspect this could be done in a better way, but it works: > > In[12]:= imagecoefficients[mat_?MatrixQ] := > With[{m = Length[mat], n = Length[First[mat]]}, > TakeMatrix[mat, > {-Count[TakeColumns[mat, n - m], _?(Union[#] == {0} &)], > n - m + 1}, > {-1, -1} > ] > ] > > In[13]:= imagecoefficients[doubleRowReduce[m]] > > Out[13]= {{1, 0, -6, 8, -3}, {0, 1, -3, 3, -1}} > > Again, this may not always give exactly the same equations as the other > methods, but the difference is due only to a change of basis in the > y-variables, leaving the image and its complementary space invariant. > > ==================================================== > > I hope this is of some use for you. > > Preston Nichols > Mathematics and Computer Science > Wittenberg University > > At 02:26 AM 06/23/2000 -400, jmcornil wrote: > >First question : >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> > > > >The function RowReduce perform a Gauss-Jordan reduction putting > >1 on the "diagonal" and 0 elsewhere so long it is possible. But is there > >an other function (or an option for RowReduce) performing a simple > >Gauss reduction i.e. putting 0 ONLY UNDER the first diagonal ? > > > >By the way I saw there is an option "Method" for the function. > > > >In[22]:= Options[RowReduce] > >Out[22]= {Method -> Automatic, Modulus -> 0, ZeroTest -> Automatic} > > > >and I did not found any help on it. How to use it ? > > > >Second question : >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> > > > >I did not foud any function to perform a "limited" Gauss reduction, i.e. > >a function performing > >a reduction up to a column number . Such an option woud be useful in a > >case as follows > >to get equations of the image of a matrice. > > > >In the following exampe, if you let RowReduce working as it will, you > >also reduce the last column > >(since formally the polynoms in yi are not ZERO) and you cannot get the > >equations of image of the matrice. > > > >In[23]:=n = 5; > > m = Table[(i + j - 1)^2, {i, n}, {j, n}]; > > v = Table[yi, {i, n}]; > > << "LinearAlgebra`MatrixManipulation`"; > > mv = AppendRows[m, Transpose[{v}]] > >Out[27]=... > > > >In[30]:=RowReduce[mv] > > > >Out[31]=... > > > >Third and (last) question : >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> > > > >So I qickly wrote some lines to get my aim. > >But I am just beginning to program with Mathematica > >and it is written as with Pascal and other classical languages. > >Is ther a mean to write such a function using more pattern-matching ? > > > >gaussElim[m_, jmax_] := > > Module[{}, > > mm = m; > > nl = Dimensions[mm][[1]]; > > nc = Dimensions[mm][[2]]; > > For[i0 = 1; j0 = 1, > > i0 <= nl; j0 <= Min[nc, jmax], > > i0++ ; j0++ , > > ip = i0; > > While[ip < nl && mm[[ip, j0]] == 0, ip++]; > > If[mm[[ip, j0]] == 0, i0--; Continue[]]; > > > > {mm[[i0]], mm[[ip]]} = {mm[[ip]], mm[[i0]]}; > > > > For[i = i0 + 1, > > i <= nl, > > i++, > > mm[[i]] = Simplify[mm[[i]] - mm[[i, j0]]/mm[[i0, > >j0]]*mm[[i0]]]]; > > ]; > > mm] > > > >In[29]:= gaussElim[mv, 5] > >Out[29]= ... > > > > > >Thanks very much for answer, > > > >Jack-Michel CORNIL > > > >VERSAILLES-FRANCE > > -- Murray Eisenberg murray at math.umass.edu Mathematics & Statistics Dept. phone 413 549-1020 (H) Univ. of Massachusetts 413 545-2859 (W) Amherst, MA 01003-4515