From f9cb6cb7305b6e62ad0d2e3f21815fd2b8676bff Mon Sep 17 00:00:00 2001 From: Matt Pagan Date: Mon, 31 Mar 2014 16:03:31 -0500 Subject: [PATCH 1/9] Added an 'elementary' function. --- numpy/lib/twodim_base.py | 116 +++++++++++++++++++++++++++++++++++++++ 1 file changed, 116 insertions(+) diff --git a/numpy/lib/twodim_base.py b/numpy/lib/twodim_base.py index 20b5cdd67a7b..a79b2af711bf 100644 --- a/numpy/lib/twodim_base.py +++ b/numpy/lib/twodim_base.py @@ -998,3 +998,119 @@ def triu_indices_from(arr, k=0): if arr.ndim != 2: raise ValueError("input array must be 2-d") return triu_indices(arr.shape[-2], k=k, m=arr.shape[-1]) + + +def elementary(N, i, j=None, t=None, dtype=float): + """ + Return an array arranged like an elementary matrix. + + Parameters + ---------- + N : int + The size of the NxN array to be returned. Elementary matrices + are be square. + i : int + The index of the first row on which operations are to be + performed. + j : int + If set, the index of the second row on which operations are to + be performed. + t : scalar + If set, the factor by which a given row will be multiplied. + + Returns + ------- + m : ndarray of shape (NxN) + An identity array that has had a single row operation performed + on it. + + See also + -------- + eye, identity + + Examples + ------- + To swap the the first and third rows of a 4x4 identity matirx: + + >>> L = elementary(4, 0, 2) + >>> L + array([[ 0., 0., 1., 0.], + [ 0., 1., 0., 0.], + [ 1., 0., 0., 0.], + [ 0., 0., 0., 1.]]) + + When the elementary array is multiplied by a given matrix, the + result is the matrix with its first and third rows swapped. + + >>> H = np.matrix([[ 2, 3, 5, 7], + [11, 13, 17, 19], + [23, 29, 31, 37], + [41, 43, 47, 53]]) + >>> L*H + matrix([[ 23., 29., 31., 37.], + [ 11., 13., 17., 19.], + [ 2., 3., 5., 7.], + [ 41., 43., 47., 53.]]) + + If the matrix is multiplied by the elementary array (i.e., the + multiplication takes place in reverse order), the result is the + given matrix with its first and third columns swapped. + + >>> H*L + matrix([[ 5., 3., 2., 7.], + [ 17., 13., 11., 19.], + [ 31., 29., 23., 37.], + [ 47., 43., 41., 53.]]) + + To add a multiple of a matrix's second row to the same matrix's + fourth row: + + >>> H + matrix([[ 2, 3, 5, 7], + [11, 13, 17, 19], + [23, 29, 31, 37], + [41, 43, 47, 53]]) + >> M=elementary(4, 1, 3, 10000, int) + >> M + array([[ 1, 0, 0, 0], + [ 0, 1, 0, 0], + [ 0, 0, 1, 0], + [ 0, 10000, 0, 1]]) + >> M*H + matrix([[ 2, 3, 5, 7], + [ 11, 13, 17, 19], + [ 23, 29, 31, 37], + [110041, 130043, 170047, 190053]]) + + To multiply only the last column of a matrix by a scalar: + + >>> H + matrix([[ 2, 3, 5, 7], + [11, 13, 17, 19], + [23, 29, 31, 37], + [41, 43, 47, 53]]) + >>> H*elementary(len(H), len(H)-1, t=2) + matrix([[ 2., 3., 5., 14.], + [ 11., 13., 17., 38.], + [ 23., 29., 31., 74.], + [ 41., 43., 47., 106.]]) + + Any nonsingular matrix can be formed by taking the product of an array + of elementary matrices. + + """ + m=eye(N, dtype=dtype) + if j==None and t==None: + raise ValueError("One or more of {0} and {1} must be set.".format( \ + 'j', 't')) + elif t==None: + swap = array(m[i]) + m[i] = m[j] + m[j] = swap + return m + elif j==None: + m[i] *= t + return m + else: + m[j] += (t * m[i]) + return m From 159b9ae42ddd1087d2965e4c10efbe983b25111e Mon Sep 17 00:00:00 2001 From: Matt Pagan Date: Tue, 1 Apr 2014 01:21:17 +0000 Subject: [PATCH 2/9] Update twodim_base.py --- numpy/lib/twodim_base.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/numpy/lib/twodim_base.py b/numpy/lib/twodim_base.py index a79b2af711bf..323134589ce1 100644 --- a/numpy/lib/twodim_base.py +++ b/numpy/lib/twodim_base.py @@ -1002,7 +1002,7 @@ def triu_indices_from(arr, k=0): def elementary(N, i, j=None, t=None, dtype=float): """ - Return an array arranged like an elementary matrix. + Return a 2-D array arranged as an elementary matrix. Parameters ---------- From 0e10a449ee2c483594d4d1303f0313eba91b9ae5 Mon Sep 17 00:00:00 2001 From: Matt Pagan Date: Wed, 9 Apr 2014 16:12:00 -0500 Subject: [PATCH 3/9] Tightened up code per github suggestions. --- numpy/lib/twodim_base.py | 41 ++++++++++++++++++++-------------------- 1 file changed, 21 insertions(+), 20 deletions(-) diff --git a/numpy/lib/twodim_base.py b/numpy/lib/twodim_base.py index a79b2af711bf..cf3b6683ce85 100644 --- a/numpy/lib/twodim_base.py +++ b/numpy/lib/twodim_base.py @@ -1000,7 +1000,7 @@ def triu_indices_from(arr, k=0): return triu_indices(arr.shape[-2], k=k, m=arr.shape[-1]) -def elementary(N, i, j=None, t=None, dtype=float): +def elementary(N, i, j=None, multiplier=None, dtype=float): """ Return an array arranged like an elementary matrix. @@ -1008,14 +1008,14 @@ def elementary(N, i, j=None, t=None, dtype=float): ---------- N : int The size of the NxN array to be returned. Elementary matrices - are be square. + are square by definition. i : int The index of the first row on which operations are to be performed. j : int If set, the index of the second row on which operations are to be performed. - t : scalar + multiplier : scalar If set, the factor by which a given row will be multiplied. Returns @@ -1028,8 +1028,12 @@ def elementary(N, i, j=None, t=None, dtype=float): -------- eye, identity + Notes + ----- + This function's primary utility is in matrix multiplication. + Examples - ------- + -------- To swap the the first and third rows of a 4x4 identity matirx: >>> L = elementary(4, 0, 2) @@ -1046,7 +1050,7 @@ def elementary(N, i, j=None, t=None, dtype=float): [11, 13, 17, 19], [23, 29, 31, 37], [41, 43, 47, 53]]) - >>> L*H + >>> L * H matrix([[ 23., 29., 31., 37.], [ 11., 13., 17., 19.], [ 2., 3., 5., 7.], @@ -1056,7 +1060,7 @@ def elementary(N, i, j=None, t=None, dtype=float): multiplication takes place in reverse order), the result is the given matrix with its first and third columns swapped. - >>> H*L + >>> H * L matrix([[ 5., 3., 2., 7.], [ 17., 13., 11., 19.], [ 31., 29., 23., 37.], @@ -1070,13 +1074,13 @@ def elementary(N, i, j=None, t=None, dtype=float): [11, 13, 17, 19], [23, 29, 31, 37], [41, 43, 47, 53]]) - >> M=elementary(4, 1, 3, 10000, int) + >> M = elementary(4, 1, 3, 10000, int) >> M array([[ 1, 0, 0, 0], [ 0, 1, 0, 0], [ 0, 0, 1, 0], [ 0, 10000, 0, 1]]) - >> M*H + >> M * H matrix([[ 2, 3, 5, 7], [ 11, 13, 17, 19], [ 23, 29, 31, 37], @@ -1089,7 +1093,7 @@ def elementary(N, i, j=None, t=None, dtype=float): [11, 13, 17, 19], [23, 29, 31, 37], [41, 43, 47, 53]]) - >>> H*elementary(len(H), len(H)-1, t=2) + >>> H * elementary(len(H), len(H)-1, multiplier=2) matrix([[ 2., 3., 5., 14.], [ 11., 13., 17., 38.], [ 23., 29., 31., 74.], @@ -1099,18 +1103,15 @@ def elementary(N, i, j=None, t=None, dtype=float): of elementary matrices. """ - m=eye(N, dtype=dtype) - if j==None and t==None: - raise ValueError("One or more of {0} and {1} must be set.".format( \ - 'j', 't')) - elif t==None: - swap = array(m[i]) - m[i] = m[j] - m[j] = swap + m = eye(N, dtype=dtype) + if j is None and multiplier is None: + raise ValueError("One or more of 'j' and 'multiplier' must be set.") + elif multiplier is None: + m[i], m[j] = m[j], m[i] return m - elif j==None: - m[i] *= t + elif j is None: + m[i] *= multiplier return m else: - m[j] += (t * m[i]) + m[j] += (multiplier * m[i]) return m From c3f45d369adc5374ca1656e419bb6da9d1ac3d39 Mon Sep 17 00:00:00 2001 From: Matt Pagan Date: Thu, 10 Apr 2014 01:41:56 -0500 Subject: [PATCH 4/9] Stopped using np.matrix, brought back the swap variable. --- numpy/lib/twodim_base.py | 113 +++++++++++++++++++-------------------- 1 file changed, 55 insertions(+), 58 deletions(-) diff --git a/numpy/lib/twodim_base.py b/numpy/lib/twodim_base.py index 48f28508208c..9cc7063952e3 100644 --- a/numpy/lib/twodim_base.py +++ b/numpy/lib/twodim_base.py @@ -1007,34 +1007,28 @@ def elementary(N, i, j=None, multiplier=None, dtype=float): Parameters ---------- N : int - The size of the NxN array to be returned. Elementary matrices - are square by definition. + The size of the NxN array to be returned. Elementary matrices are + square by definition. i : int - The index of the first row on which operations are to be - performed. + The index of the first row on which operations are to be performed. j : int - If set, the index of the second row on which operations are to - be performed. + If set, the index of the second row on which operations are to be + performed. multiplier : scalar If set, the factor by which a given row will be multiplied. Returns ------- m : ndarray of shape (NxN) - An identity array that has had a single row operation performed - on it. + An identity array that has had a single row operation performed on it. See also -------- eye, identity - Notes - ----- - This function's primary utility is in matrix multiplication. - Examples -------- - To swap the the first and third rows of a 4x4 identity matirx: + To swap the the first and third rows of a 4x4 identity array: >>> L = elementary(4, 0, 2) >>> L @@ -1043,71 +1037,73 @@ def elementary(N, i, j=None, multiplier=None, dtype=float): [ 1., 0., 0., 0.], [ 0., 0., 0., 1.]]) - When the elementary array is multiplied by a given matrix, the - result is the matrix with its first and third rows swapped. + When the elementary array is matrix-multiplied by a given array, the + result is the given array with its first and third rows swapped. - >>> H = np.matrix([[ 2, 3, 5, 7], - [11, 13, 17, 19], - [23, 29, 31, 37], - [41, 43, 47, 53]]) - >>> L * H - matrix([[ 23., 29., 31., 37.], - [ 11., 13., 17., 19.], - [ 2., 3., 5., 7.], - [ 41., 43., 47., 53.]]) + >>> H = np.array([[ 2, 3, 5, 7], + [11, 13, 17, 19], + [23, 29, 31, 37], + [41, 43, 47, 53]]) + >>> L.dot(H) + array([[ 23., 29., 31., 37.], + [ 11., 13., 17., 19.], + [ 2., 3., 5., 7.], + [ 41., 43., 47., 53.]]) - If the matrix is multiplied by the elementary array (i.e., the - multiplication takes place in reverse order), the result is the - given matrix with its first and third columns swapped. + If the array is matrix-multiplied by the elementary array (i.e., the + multiplication takes place in reverse order), the result is the given + array with its first and third columns swapped. - >>> H * L - matrix([[ 5., 3., 2., 7.], - [ 17., 13., 11., 19.], - [ 31., 29., 23., 37.], - [ 47., 43., 41., 53.]]) + >>> H.dot(L) + array([[ 5., 3., 2., 7.], + [ 17., 13., 11., 19.], + [ 31., 29., 23., 37.], + [ 47., 43., 41., 53.]]) - To add a multiple of a matrix's second row to the same matrix's - fourth row: + To add a multiple of an array's second row to the same array's fourth row: >>> H - matrix([[ 2, 3, 5, 7], - [11, 13, 17, 19], - [23, 29, 31, 37], - [41, 43, 47, 53]]) + array([[ 2, 3, 5, 7], + [11, 13, 17, 19], + [23, 29, 31, 37], + [41, 43, 47, 53]]) >> M = elementary(4, 1, 3, 10000, int) >> M array([[ 1, 0, 0, 0], [ 0, 1, 0, 0], [ 0, 0, 1, 0], [ 0, 10000, 0, 1]]) - >> M * H - matrix([[ 2, 3, 5, 7], - [ 11, 13, 17, 19], - [ 23, 29, 31, 37], - [110041, 130043, 170047, 190053]]) + >> M.dot(H) + array([[ 2, 3, 5, 7], + [ 11, 13, 17, 19], + [ 23, 29, 31, 37], + [110041, 130043, 170047, 190053]]) - To multiply only the last column of a matrix by a scalar: + To multiply only the last column of an array by a scalar: >>> H - matrix([[ 2, 3, 5, 7], - [11, 13, 17, 19], - [23, 29, 31, 37], - [41, 43, 47, 53]]) - >>> H * elementary(len(H), len(H)-1, multiplier=2) - matrix([[ 2., 3., 5., 14.], - [ 11., 13., 17., 38.], - [ 23., 29., 31., 74.], - [ 41., 43., 47., 106.]]) - - Any nonsingular matrix can be formed by taking the product of an array - of elementary matrices. + array([[ 2, 3, 5, 7], + [11, 13, 17, 19], + [23, 29, 31, 37], + [41, 43, 47, 53]]) + >>> G = elementary(len(H), len(H)-1, multiplier=2) + >>> H.dot(G) + array([[ 2., 3., 5., 14.], + [ 11., 13., 17., 38.], + [ 23., 29., 31., 74.], + [ 41., 43., 47., 106.]]) + + Any nonsingular array can be formed by taking the product of an array of + elementary arrays. """ m = eye(N, dtype=dtype) if j is None and multiplier is None: - raise ValueError("One or more of 'j' and 'multiplier' must be set.") + raise ValueError("'j' and 'multiplier' are both set to None.") elif multiplier is None: - m[i], m[j] = m[j], m[i] + swap = array(m[i]) + m[i] = m[j] + m[j] = swap return m elif j is None: m[i] *= multiplier @@ -1115,3 +1111,4 @@ def elementary(N, i, j=None, multiplier=None, dtype=float): else: m[j] += (multiplier * m[i]) return m + From 9dba2618220d00aa22da9322f202c7e15e8c6bae Mon Sep 17 00:00:00 2001 From: Matt Pagan Date: Sun, 13 Apr 2014 20:19:35 -0500 Subject: [PATCH 5/9] Added tests for the elementary function. --- numpy/lib/tests/test_twodim_base.py | 24 ++++++++++++++++++++++++ 1 file changed, 24 insertions(+) diff --git a/numpy/lib/tests/test_twodim_base.py b/numpy/lib/tests/test_twodim_base.py index 9e81cfe4b26a..5371c6f400c4 100644 --- a/numpy/lib/tests/test_twodim_base.py +++ b/numpy/lib/tests/test_twodim_base.py @@ -465,6 +465,30 @@ def test_dtypes(self): # assert_array_almost_equal). yield (assert_array_equal, v, expected) +class TestElementary(object): + def test_switch_rows(self): + G = elementary(4, 0, 3) + assert_array_equal(G, np.array([[0., 0., 0., 1.], + [0., 1., 0., 0.], + [0., 0., 1., 0.], + [1., 0., 0., 0.]])) + + + def test_add_row_multiple(self): + H = elementary(5, 1, 3, 7, 'i') + assert_array_equal(H, np.array([[1, 0, 0, 0, 0], + [0, 1, 0, 0, 0], + [0, 0, 1, 0, 0], + [0, 7, 0, 1, 0], + [0, 0, 0, 0, 1]])) + + + def test_multiply_row(self): + I = elementary(3, 0, multiplier=6, dtype=int) + assert_array_equal(I, np.array([[6, 0, 0], + [0, 1, 0], + [0, 0, 1]])) + if __name__ == "__main__": run_module_suite() From 1da28fd6c9a77803089df2a6a0bdb005c20cfb0e Mon Sep 17 00:00:00 2001 From: Matt Pagan Date: Sun, 13 Apr 2014 20:52:06 -0500 Subject: [PATCH 6/9] Added 'elementary' to __all__ --- numpy/lib/twodim_base.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/numpy/lib/twodim_base.py b/numpy/lib/twodim_base.py index 9cc7063952e3..ce8c7b514c6f 100644 --- a/numpy/lib/twodim_base.py +++ b/numpy/lib/twodim_base.py @@ -6,7 +6,7 @@ __all__ = ['diag', 'diagflat', 'eye', 'fliplr', 'flipud', 'rot90', 'tri', 'triu', 'tril', 'vander', 'histogram2d', 'mask_indices', 'tril_indices', 'tril_indices_from', 'triu_indices', - 'triu_indices_from', + 'triu_indices_from', 'elementary', ] from numpy.core.numeric import ( From 4de601b664c828be59d8f41c8d772fe9f2c6a9c7 Mon Sep 17 00:00:00 2001 From: Matt Pagan Date: Sun, 13 Apr 2014 21:09:56 -0500 Subject: [PATCH 7/9] Tests for 'elementary' should inherit TestCase. --- numpy/lib/tests/test_twodim_base.py | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/numpy/lib/tests/test_twodim_base.py b/numpy/lib/tests/test_twodim_base.py index 5371c6f400c4..0919eb621f98 100644 --- a/numpy/lib/tests/test_twodim_base.py +++ b/numpy/lib/tests/test_twodim_base.py @@ -465,10 +465,10 @@ def test_dtypes(self): # assert_array_almost_equal). yield (assert_array_equal, v, expected) -class TestElementary(object): +class TestElementary(TestCase): def test_switch_rows(self): G = elementary(4, 0, 3) - assert_array_equal(G, np.array([[0., 0., 0., 1.], + assert_equal(G, np.array([[0., 0., 0., 1.], [0., 1., 0., 0.], [0., 0., 1., 0.], [1., 0., 0., 0.]])) @@ -476,7 +476,7 @@ def test_switch_rows(self): def test_add_row_multiple(self): H = elementary(5, 1, 3, 7, 'i') - assert_array_equal(H, np.array([[1, 0, 0, 0, 0], + assert_equal(H, np.array([[1, 0, 0, 0, 0], [0, 1, 0, 0, 0], [0, 0, 1, 0, 0], [0, 7, 0, 1, 0], @@ -485,7 +485,7 @@ def test_add_row_multiple(self): def test_multiply_row(self): I = elementary(3, 0, multiplier=6, dtype=int) - assert_array_equal(I, np.array([[6, 0, 0], + assert_equal(I, np.array([[6, 0, 0], [0, 1, 0], [0, 0, 1]])) From 08965f56aab992654ed7664df70a1f8e034673d2 Mon Sep 17 00:00:00 2001 From: Matt Pagan Date: Sun, 13 Apr 2014 21:23:38 -0500 Subject: [PATCH 8/9] Eliminated the swap variable correctly. --- numpy/lib/twodim_base.py | 4 +--- 1 file changed, 1 insertion(+), 3 deletions(-) diff --git a/numpy/lib/twodim_base.py b/numpy/lib/twodim_base.py index ce8c7b514c6f..3639c5f0caef 100644 --- a/numpy/lib/twodim_base.py +++ b/numpy/lib/twodim_base.py @@ -1101,9 +1101,7 @@ def elementary(N, i, j=None, multiplier=None, dtype=float): if j is None and multiplier is None: raise ValueError("'j' and 'multiplier' are both set to None.") elif multiplier is None: - swap = array(m[i]) - m[i] = m[j] - m[j] = swap + m[[i, j]] = m[[j, i]] return m elif j is None: m[i] *= multiplier From 480c51b17c608a2cbcca65b8a633dbe1b1654070 Mon Sep 17 00:00:00 2001 From: Matt Pagan Date: Mon, 14 Apr 2014 14:25:45 -0500 Subject: [PATCH 9/9] Tests need to call np.elementary rather than elementary. --- numpy/lib/tests/test_twodim_base.py | 7 ++++--- 1 file changed, 4 insertions(+), 3 deletions(-) diff --git a/numpy/lib/tests/test_twodim_base.py b/numpy/lib/tests/test_twodim_base.py index 0919eb621f98..9548cad30069 100644 --- a/numpy/lib/tests/test_twodim_base.py +++ b/numpy/lib/tests/test_twodim_base.py @@ -465,9 +465,10 @@ def test_dtypes(self): # assert_array_almost_equal). yield (assert_array_equal, v, expected) + class TestElementary(TestCase): def test_switch_rows(self): - G = elementary(4, 0, 3) + G = np.elementary(4, 0, 3) assert_equal(G, np.array([[0., 0., 0., 1.], [0., 1., 0., 0.], [0., 0., 1., 0.], @@ -475,7 +476,7 @@ def test_switch_rows(self): def test_add_row_multiple(self): - H = elementary(5, 1, 3, 7, 'i') + H = np.elementary(5, 1, 3, 7, 'i') assert_equal(H, np.array([[1, 0, 0, 0, 0], [0, 1, 0, 0, 0], [0, 0, 1, 0, 0], @@ -484,7 +485,7 @@ def test_add_row_multiple(self): def test_multiply_row(self): - I = elementary(3, 0, multiplier=6, dtype=int) + I = np.elementary(3, 0, multiplier=6, dtype=int) assert_equal(I, np.array([[6, 0, 0], [0, 1, 0], [0, 0, 1]]))