I have found a code that pivotize a square matrix for LU decomposition, but I can't understand some of them.
def pivotize(m):
"""Creates the pivoting matrix for m."""
n = len(m)
ID = [[float(i == j) for i in xrange(n)] for j in xrange(n)]
for j in xrange(n):
row = max(xrange(j, n), key=lambda i: abs(m[i][j]))
if j != row:
ID[j], ID[row] = ID[row], ID[j]
return ID
First, isn't the line for ID simply the identity matrix? any advantage doing this?
Second, I can't really understand the line for row. I know lambda is used to define a function in the text, and it simply returns the value of M_ij once the value of i is supplied (and the value of j depends on the for loop), but what is i?
And isn't xrange similar to range? But what does it return here?
And when combined with the function max, what happens? I just don't know what are the things inside the max function that are being compared.
Sorry if this question sounds stupid. I'm pretty new to programming
First, isn't the line for ID simply the identity matrix?
Yes.
Second, I can't really understand the line for row....
See this for a discussion about the max/key/lambda interaction. To answer "what is i
?", its the argument to the lambda function, i
could equivalently be x
for foo
. (For clarity, yielding abs(m[x][j])
and abs(m[foo][j])
, respectively).
And isn't xrange similar to range?
Yes. In Python 2 xrange
returns a sequence object that lazily computes the next value, only when needed. See this for more info. When looped through in entirety, range
and xrange
will produce the same results, with different implementations.
But what does it return here?
Line 5's xrange(n)
will return the integer values from 0 to (n-1) while Line 6's xrange(j, n)
will return the integer values from j to (n-1).
EDIT
More on lambda:
Consider how you might take a given sequence of numbers, and double each of them. First define a function that doubles a number x
and returns that value. Then map that function to each element of a sequence.
# Traditional
def double(x): return x*2
print map(double, [1,2,3,4]) # [2, 4, 6, 8]
You could have alternatively used an anonymous (lambda) function to do the same thing:
# Lambda
print map(lambda i: i*2, [1,2,3,4]) # [2, 4, 6, 8]
Note that everything is the same, except the definition of the double
function is gone and the reference to the function in the call to map
is replaced by the "inline" lambda function.
EDIT 2
And when combined with the function max, what happens?
This line row = max(xrange(j, n), key=lambda i: abs(m[i][j]))
can be broken down as follows:
xrange(j,n)
generates a sequence of integers from j
(inclusive) to n
(exclusive).i
in the lambda function. The lambda function "returns" the absolute value of the i
th row and j
th column.[1]max
function then finds the greatest value of these "lambda outputs" and sets row
equal to it.This could have alternatively been written as a max of a list comprehension:
row = max( [abs(m[i][j]) for i in xrange(j,n)] )
Or as Dan D. points out in his comment, written as a generator expression (without creating an intermediary list) as simply:
row = max( abs(m[i][j]) for i in xrange(j,n) )
Notes:
[1] Some assumptions here, but row-column is a standard way of expressing matrices.