s [ Pobierz całość w formacie PDF ]

kind of situation that might result:
îø ùø
1 0 0 a b c
ïø úø
0 1 0 d e f
ïø úø
ïø úø
ïø 0 0 1 g h k úø
. (3.6.1)
ïø úø
. . . .
ïø úø
. . . .
ðø ûø
. . · · · · · · . .
3 5 2 1 4 "
The matrix above represents schematically the reduced row echelon form in a problem
where there are five unknowns (n = 5), the pseudorank r = 3, just one right-hand side
vector is given (p = 1), and the permuations that were carried out on the columns are
recorded in the array Ä : [3, 5, 2, 1, 4] shown in the last row of the matrix as it would be
storded in a computation.
The question now is, how do we express the general solution of the given set of equations?
To find the answer, let s go back to the set of equations that (3.6.1) stands for. The first of
these is
x3 = c - ax1 - bx4 (3.6.2)
106 Numerical linear algebra
because the numbering of the unknowns is as shown in the Ä array. The next two equations
are
x5 = f - dx1 - ex4
(3.6.3)
x2 = k - gx1 - hx4.
If we add the two trivial equations x1 = x1 and x4 = x4, then we get the whole solution
vector which, after re-ordering the equations, can be written as
îø ùø îø ùø îø ùø îø ùø
x1 0 -1 0
ïø úø ïø úø ïø úø ïø úø
x2 k g h
ïø úø ïø úø ïø úø ïø úø
ïø úø ïø úø ïø úø ïø úø
ïø x3 úø = ïø c úø +(-x1) " ïø a úø +(-x4) " ïø b úø . (3.6.4)
ïø úø ïø úø ïø úø ïø úø
ðø ûø ðø ûø ðø ûø ðø -1
ûø
x4 0 0
x5 f d e
Now we are looking at a display of the output as we would like our subroutine to give
it. The three vectors on the right side of (3.6.4) are, respectively, a particular solution of
the given system of equations, and the two vectors of a basis for the kernel of the coefficient
matrix.
The question can now be rephrased: exactly what operations must be done to the matrix
shown in (3.6.1) that represents the situation at the end of the back solution, in order to
obtain the three vectors in (3.6.4)?
The first things to do are, as we have previously noted, to append the negative of a 2× 2
identity matrix to the bottom of the fourth and fifth columns of (3.6.1), and to lengthen
the last column on the right by appending two more zeros. That brings us to the matrix
îø ùø
1 0 0 a b c
ïø úø
0 1 0 d e f
ïø úø
ïø úø
0 0 1 g h k
ïø úø
ïø úø . (3.6.5)
ïø -1 0 0 úø
ïø úø
ðø ûø
0 -1 0
[3 5 2 1 4] "
The first two of the three long columns above will be the basis for the kernel, and the last
column above will be the particular solution, but only after we do the right rearrangement.
Now here is the punch line: the right rearrangement to do is to permute the rows of
those three long columns as described by the permuation Ä.
That means that the first row becomes the third, the second row becomes the fifth, the
third row becomes the second, the fourth row is the new first, and the old fifth row is the
new fourth. The reader is invited to carry out on the rows the interchanges just described,
and to compare the result with what we want, namely with (3.6.4). It will be seen that we
have gotten the desired result.
The point that is just a little surprising is that to undo the column interchanges that
are recorded by Ä, we do row interchanges. Just roughly, the reason for this is that we
begin by wanting to solve Ax = b, and instead we end up solving (AE)y = b, where E is
3.6 To unscramble the eggs 107
a matrix obtained from the identity by elementary column operations. Evidently, x = Ey,
which means that we must perform row operations on y to recover the answers in the right
order.
Now we can leave the example above, and state the rule in general. We are given p
systems of m simultaneous equations each, all having a commom m × n coefficient matrix
A, in n unknowns. At the end of the back solution we will have before us a matrix of the
form
îø ùø
ïø úø
I(r, r) B(r, n - r) P (r, p) (3.6.6)
ðø ûø
where I(r, r) is the r×r identity matrix, r is the pseudorank of A, and B and P are matrices
of the sizes shown.
We adjoin under B the negative of the (n - r) × (n - r) identity matrix, and under P
we adjoin an (n - r) × p block of zeros. Next, we forget the identity matrix on the left,
and we consider the entire remaining n × (n - r + p) matrix as a whole, call it T , say. Now
we exchange the rows of T according to the permuation array Ä. Preciesly, row 1 of T will
be row Ä1 of the new T , row 2 will be row Ä2, . . . . Conceptually, we should regard the
old T and the new T as occupying different areas of storage, so that the new T is just a
rearrangement of the rows of the old.
Now the first n - r columns of the newT are a basis for the kernel of A, and should be
output as such, and the jth one of the last p columns of the new T is a particular solution
of the jth one of the input systems of equations, and should be output as such.
Although conceptually we should think of the old T and the new T as occupying distinct
arrays in memory, in fact it is perfectly possible to carry out the whole row interchange
procedure described above in just one array, the one that holds T , without ever  stepping
on our own toes, so let s consider that problem.
Suppose a linear array a = [a1, . . . , an] is given, along with a permutation array Ä =
[Ä1, . . . , Än]. We want to rearrange the entries of the array a according to the permuation Ä
without using any additional array storage. Thus the present array a1 will end up as the
output aÄ1, the initial a2 will end up as aÄ2, etc.
To do this with no extra array storage, let s first pick up the element a1 and move it to
aÄ1, being careful to store the original aÄ1 in a temporary location t so it won t be destroyed.
Next we move the contents of t to its destination, and so forth. After a certain number of
steps (maybe only 1), we will be back to a1.
Her s an example to help clarify the situation. Suppose the arrays a and Ä at input time
were:
a =[5, 7, 13, 9, 2, 8]
(3.6.7)
Ä =[3, 4, 5, 2, 1, 6].
So we move the 5 in position a1 to position a3 (after putting the 13 into a safe place), and [ Pobierz całość w formacie PDF ]