Skip to content

Commit

Permalink
added __setitem__ to hmat
Browse files Browse the repository at this point in the history
  • Loading branch information
maekke97 committed May 4, 2017
1 parent 9d569f5 commit 38a562c
Show file tree
Hide file tree
Showing 4 changed files with 124 additions and 15 deletions.
14 changes: 9 additions & 5 deletions HierMat/hmat.py
Original file line number Diff line number Diff line change
Expand Up @@ -64,6 +64,14 @@ def __getitem__(self, item):
structured_blocks = {block.parent_index: block for block in self.blocks}
return structured_blocks[item]

def __setitem__(self, key, value):
"""Set block at index key to value"""
if isinstance(key, int):
self.blocks[key] = value
else:
index = self.blocks.index(self[key])
self.blocks[index] = value

def __repr__(self):
return '<HMat with {content}>'.format(content=self.blocks if self.blocks else self.content)

Expand Down Expand Up @@ -608,11 +616,7 @@ def transpose(self):
return HMat(content=out_content, blocks=out_blocks, shape=out_shape, parent_index=self.parent_index)

def inv(self):
"""Invert the matrix according to chapter 3.7 in :cite:`hackbusch2015hierarchical`
.. caution::
Inversion will change the block structure. The result will always have 4 blocks
"""Invert the matrix according to chapter 7.5 in :cite:`hackbusch2015hierarchical`
:raises: :class:`numpy.LinAlgError` if matrix is singular
"""
Expand Down
2 changes: 1 addition & 1 deletion setup.py
Original file line number Diff line number Diff line change
Expand Up @@ -35,7 +35,7 @@

setup(
name='HierMat',
version='0.6.5',
version='0.6.6',
packages=['HierMat'],
url='http://hierarchicalmatrices.readthedocs.io/en/latest/index.html',
download_url='https://github.com/maekke97/HierarchicalMatrices',
Expand Down
11 changes: 11 additions & 0 deletions tests/test_hmat.py
Original file line number Diff line number Diff line change
Expand Up @@ -72,6 +72,17 @@ def test_get_item(self):
self.assertEqual(self.consistent2[4], self.cmat5T)
self.assertEqual(self.consistent2[2, 3], self.cmat5T)

def test_set_item(self):
self.consistent1[0, 0] = self.cmat1T
self.assertEqual(self.consistent1[0, 0], self.cmat1T)
self.consistent1[0, 2] = self.cmat2T
self.consistent1[0, 3] = self.cmat3T
self.consistent1[3, 0] = self.cmat4T
self.consistent1[3, 2] = self.cmat5T
self.consistent1[3, 3] = self.cmat6T
self.consistent1.shape = (6, 5)
self.assertEqual(self.consistent1, self.consistent2)

def test_determine_block_structure(self):
check = {(0, 0): (3, 4), (0, 4): (3, 2), (3, 0): (4, 2), (3, 2): (4, 4)}
self.assertEqual(check, self.hmat.block_structure())
Expand Down
112 changes: 103 additions & 9 deletions tests/test_math.py
Original file line number Diff line number Diff line change
Expand Up @@ -4,21 +4,23 @@


class TestMath(TestCase):
precision = 14

def test_rmat_vec_mul(self):
left = numpy.matrix(numpy.random.rand(3, 1))
right = numpy.matrix(numpy.random.rand(3, 1))
rmat = RMat(left, right, 1)
x = numpy.random.rand(3, 1)
check1 = rmat.to_matrix() * x
check2 = rmat * x
self.assertAlmostEqual(numpy.linalg.norm(check1 - check2), 0, places=14)
self.assertAlmostEqual(numpy.linalg.norm(check1 - check2), 0, places=self.precision)
left = numpy.matrix(numpy.random.rand(5, 2))
right = numpy.matrix(numpy.random.rand(4, 2))
rmat = RMat(left, right, 2)
x = numpy.random.rand(4, 1)
check1 = rmat.to_matrix() * x
check2 = rmat * x
self.assertAlmostEqual(numpy.linalg.norm(check1 - check2), 0, places=14)
self.assertAlmostEqual(numpy.linalg.norm(check1 - check2), 0, places=self.precision)

def test_rmat_add_exact(self):
left1 = numpy.matrix(numpy.random.rand(3, 1))
Expand All @@ -29,7 +31,7 @@ def test_rmat_add_exact(self):
rmat2 = RMat(left2, right2)
check1 = rmat1.to_matrix() + rmat2.to_matrix()
check2 = (rmat1 + rmat2).to_matrix()
self.assertAlmostEqual(numpy.linalg.norm(check1 - check2), 0, places=14)
self.assertAlmostEqual(numpy.linalg.norm(check1 - check2), 0, places=self.precision)
left1 = numpy.matrix(numpy.random.rand(5, 1))
right1 = numpy.matrix(numpy.random.rand(4, 1))
rmat1 = RMat(left1, right1)
Expand All @@ -38,7 +40,7 @@ def test_rmat_add_exact(self):
rmat2 = RMat(left2, right2)
check1 = rmat1.to_matrix() + rmat2.to_matrix()
check2 = (rmat1 + rmat2).to_matrix()
self.assertAlmostEqual(numpy.linalg.norm(check1 - check2), 0, places=14)
self.assertAlmostEqual(numpy.linalg.norm(check1 - check2), 0, places=self.precision)

def test_rmat_add_subtract_exact(self):
left1 = numpy.matrix(numpy.random.rand(3, 1))
Expand All @@ -50,8 +52,8 @@ def test_rmat_add_subtract_exact(self):
step1 = rmat1 + rmat2
check1 = step1 - rmat2
check2 = step1 - rmat1
self.assertAlmostEqual(check1.norm(), rmat1.norm(), places=14)
self.assertAlmostEqual(check2.norm(), rmat2.norm(), places=14)
self.assertAlmostEqual(check1.norm(), rmat1.norm(), places=self.precision)
self.assertAlmostEqual(check2.norm(), rmat2.norm(), places=self.precision)
left1 = numpy.matrix(numpy.random.rand(5, 1))
right1 = numpy.matrix(numpy.random.rand(4, 1))
rmat1 = RMat(left1, right1)
Expand All @@ -61,8 +63,10 @@ def test_rmat_add_subtract_exact(self):
step1 = rmat1 + rmat2
check1 = step1 - rmat2
check2 = step1 - rmat1
self.assertAlmostEqual(check1.norm(), rmat1.norm(), places=14)
self.assertAlmostEqual(check2.norm(), rmat2.norm(), places=14)
self.assertAlmostEqual(check1.norm(), rmat1.norm(), places=self.precision)
self.assertAlmostEqual(check2.norm(), rmat2.norm(), places=self.precision)
res = rmat1 + rmat1 - 2 * rmat1
self.assertAlmostEqual(res.norm(), 0, places=self.precision)

def test_hmat_vec_mul(self):
blocks = []
Expand All @@ -78,4 +82,94 @@ def test_hmat_vec_mul(self):
check1 = hmat * x
hmat_full = hmat.to_matrix()
check2 = hmat_full * x
self.assertAlmostEqual(numpy.linalg.norm(check1 - check2), 0, places=14)
self.assertAlmostEqual(numpy.linalg.norm(check1 - check2), 0, places=self.precision)

def test_hmat_add_subtract(self):
blocks = []
for i in xrange(3):
for j in xrange(3):
blocks.append(HMat(content=numpy.matrix(numpy.random.rand(2, 2)),
shape=(2, 2),
parent_index=(2 * i, 2 * j)
)
)
hmat = HMat(blocks=blocks, shape=(6, 6), parent_index=(0, 0))
res = hmat - hmat
self.assertAlmostEqual(res.norm(), 0, places=self.precision)
res = hmat + hmat - 2 * hmat
self.assertAlmostEqual(res.norm(), 0, places=self.precision)
x = numpy.random.rand(6, 1)
check1 = (hmat + hmat) * x
check2 = hmat * x + hmat * x
self.assertAlmostEqual(numpy.linalg.norm(check1 - check2), 0, places=self.precision)

def test_hmat_multiplication(self):
blocks1 = []
for i in xrange(2):
for j in xrange(2):
blocks1.append(HMat(content=numpy.matrix(numpy.random.rand(2, 2)),
shape=(2, 2),
parent_index=(2 * i, 2 * j)
)
)
hmat1 = HMat(blocks=blocks1, shape=(4, 4), parent_index=(0, 0))
x = numpy.random.rand(4, 1)
res1 = (hmat1 * hmat1) * x
res2 = hmat1 * (hmat1 * x)
self.assertAlmostEqual(numpy.linalg.norm(res1 - res2), 0, places=self.precision)

def test_hmat_inversion(self):
blocks1 = []
for i in xrange(2):
for j in xrange(2):
blocks1.append(HMat(content=numpy.matrix(numpy.random.rand(2, 2)),
shape=(2, 2),
parent_index=(2 * i, 2 * j)
)
)
hmat1 = HMat(blocks=blocks1, shape=(4, 4), parent_index=(0, 0))
res = hmat1 * hmat1.inv()
check = numpy.matrix(numpy.eye(4))
self.assertAlmostEqual(numpy.linalg.norm(res.to_matrix() - check), 0, places=12)
blocks2 = []
for i in xrange(2):
for j in xrange(2):
blocks2.append(HMat(content=numpy.matrix(numpy.random.rand(3, 3)),
shape=(3, 3),
parent_index=(3 * i, 3 * j)
)
)
hmat2 = HMat(blocks=blocks2, shape=(6, 6), parent_index=(0, 0))
res = hmat2 * hmat2.inv()
check = numpy.matrix(numpy.eye(6))
self.assertAlmostEqual(numpy.linalg.norm(res.to_matrix() - check), 0, places=12)
blocks3 = []
for i in xrange(3):
for j in xrange(3):
blocks3.append(HMat(content=numpy.matrix(numpy.random.rand(3, 3)),
shape=(3, 3),
parent_index=(3 * i, 3 * j)
)
)
hmat3 = HMat(blocks=blocks3, shape=(9, 9), parent_index=(0, 0))
x = numpy.random.rand(9, 1)
y = hmat3 * x
z = hmat3.inv() * y
self.assertAlmostEqual(numpy.linalg.norm(x - z), 0, places=11)
# blocks4 = []
# for outer_i in xrange(2):
# for outer_j in xrange(2):
# inner_blocks = []
# for i in xrange(3):
# for j in xrange(3):
# inner_blocks.append(HMat(content=numpy.matrix(numpy.random.rand(3, 3)),
# shape=(3, 3),
# parent_index=(3 * i, 3 * j)
# )
# )
# blocks4.append(HMat(blocks=inner_blocks, shape=(9, 9), parent_index=(9 * outer_i, 9 * outer_j)))
# hmat4 = HMat(blocks=blocks4, shape=(18, 18), parent_index=(0, 0))
# x = numpy.random.rand(18, 1)
# y = hmat4 * x
# z = hmat4.inv() * y
# self.assertAlmostEqual(numpy.linalg.norm(x - z), 0, places=12)

0 comments on commit 38a562c

Please sign in to comment.