Mathematica 9 is now available
Services & Resources / Wolfram Forums
-----
 /
MathGroup Archive
2000
*January
*February
*March
*April
*May
*June
*July
*August
*September
*October
*November
*December
*Archive Index
*Ask about this page
*Print this page
*Give us feedback
*Sign up for the Wolfram Insider

MathGroup Archive 2000

[Date Index] [Thread Index] [Author Index]

Search the Archive

Re: A few questions about RowReduce - bis

  • To: mathgroup at smc.vnet.net
  • Subject: [mg24103] Re: [mg24074] A few questions about RowReduce - bis
  • From: "Preston Nichols" <pnichols at wittenberg.edu>
  • Date: Tue, 27 Jun 2000 00:51:57 -0400 (EDT)
  • Sender: owner-wri-mathgroup at wolfram.com

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
>



  • Prev by Date: Re: "misbehaving" Union function
  • Next by Date: TableForm
  • Previous by thread: A few questions about RowReduce - bis
  • Next by thread: Looking for Multidimensional Scaling (MDS) Notebooks...