11using Test
22using BlockSparseArrays
3- using BlockSparseArrays:
4- BlockSparseArray, svd, notrunc, truncbelow, truncdim, BlockDiagonal
3+ using BlockSparseArrays: BlockSparseArray, svd, notrunc, truncbelow, truncdim, BlockDiagonal
54using BlockArrays
65using LinearAlgebra: LinearAlgebra, Diagonal, svdvals
6+ using Random
77
88function test_svd (a, usv)
9- U, S, V = usv
9+ U, S, V = usv
1010
11- @test U * Diagonal (S) * V' ≈ a
12- @test U' * U ≈ LinearAlgebra. I
13- @test V' * V ≈ LinearAlgebra. I
11+ @test U * Diagonal (S) * V' ≈ a
12+ @test U' * U ≈ LinearAlgebra. I
13+ @test V' * V ≈ LinearAlgebra. I
1414end
1515
1616# regular matrix
1717# --------------
1818sizes = ((3 , 3 ), (4 , 3 ), (3 , 4 ))
1919eltypes = (Float32, Float64, ComplexF64)
2020@testset " ($m , $n ) Matrix{$T }" for ((m, n), T) in Iterators. product (sizes, eltypes)
21- a = rand (m, n)
22- usv = @inferred svd (a)
23- test_svd (a, usv)
21+ a = rand (m, n)
22+ usv = @inferred svd (a)
23+ test_svd (a, usv)
2424end
2525
2626# block matrix
2727# ------------
2828blockszs = (([2 , 2 ], [2 , 2 ]), ([2 , 2 ], [2 , 3 ]), ([2 , 2 , 1 ], [2 , 3 ]), ([2 , 3 ], [2 ]))
2929@testset " ($m , $n ) BlockMatrix{$T }" for ((m, n), T) in Iterators. product (blockszs, eltypes)
30- a = mortar ([rand (T, i, j) for i in m, j in n])
31- usv = svd (a)
32- test_svd (a, usv)
33- @test usv. U isa BlockedMatrix
34- @test usv. Vt isa BlockedMatrix
35- @test usv. S isa BlockedVector
30+ a = mortar ([rand (T, i, j) for i in m, j in n])
31+ usv = svd (a)
32+ test_svd (a, usv)
33+ @test usv. U isa BlockedMatrix
34+ @test usv. Vt isa BlockedMatrix
35+ @test usv. S isa BlockedVector
3636end
3737
3838# Block-Diagonal matrices
3939# -----------------------
4040@testset " ($m , $n ) BlockDiagonal{$T }" for ((m, n), T) in
4141 Iterators. product (blockszs, eltypes)
42- a = BlockDiagonal ([rand (T, i, j) for (i, j) in zip (m, n)])
43- usv = svd (a)
44- # TODO : `BlockDiagonal * Adjoint` errors
45- test_svd (a, usv)
46- @test usv. U isa BlockDiagonal
47- @test usv. Vt isa BlockDiagonal
48- @test usv. S isa BlockVector
42+ a = BlockDiagonal ([rand (T, i, j) for (i, j) in zip (m, n)])
43+ usv = svd (a)
44+ # TODO : `BlockDiagonal * Adjoint` errors
45+ test_svd (a, usv)
46+ @test usv. U isa BlockDiagonal
47+ @test usv. Vt isa BlockDiagonal
48+ @test usv. S isa BlockVector
4949end
5050
5151a = mortar ([rand (2 , 2 ) for i in 1 : 2 , j in 1 : 3 ])
@@ -60,24 +60,24 @@ test_svd(a, usv)
6060# -----------
6161@testset " ($m , $n ) BlockDiagonal{$T }" for ((m, n), T) in
6262 Iterators. product (blockszs, eltypes)
63- a = BlockSparseArray {T} (m, n)
64- for i in LinearAlgebra. diagind (blocks (a))
65- I = CartesianIndices (blocks (a))[i]
66- a[Block (I. I... )] = rand (T, size (blocks (a)[i]))
67- end
68- perm = Random. randperm (length (m))
69- a = a[Block .(perm), Block .(1 : length (n))]
63+ a = BlockSparseArray {T} (m, n)
64+ for i in LinearAlgebra. diagind (blocks (a))
65+ I = CartesianIndices (blocks (a))[i]
66+ a[Block (I. I... )] = rand (T, size (blocks (a)[i]))
67+ end
68+ perm = Random. randperm (length (m))
69+ a = a[Block .(perm), Block .(1 : length (n))]
7070
71- # errors because `blocks(a)[CartesianIndex.(...)]` is not implemented
72- usv = svd (a)
73- # TODO : `BlockDiagonal * Adjoint` errors
74- test_svd (a, usv)
75- @test usv. U isa BlockDiagonal
76- @test usv. Vt isa BlockDiagonal
77- @test usv. S isa BlockVector
71+ # errors because `blocks(a)[CartesianIndex.(...)]` is not implemented
72+ usv = svd (a)
73+ # TODO : `BlockDiagonal * Adjoint` errors
74+ test_svd (a, usv)
75+ @test usv. U isa BlockDiagonal
76+ @test usv. Vt isa BlockDiagonal
77+ @test usv. S isa BlockVector
7878
79- test_svd (a, usv2)
80- @test usv. U isa BlockDiagonal
81- @test usv. Vt isa BlockDiagonal
82- @test usv. S isa BlockVector
79+ test_svd (a, usv2)
80+ @test usv. U isa BlockDiagonal
81+ @test usv. Vt isa BlockDiagonal
82+ @test usv. S isa BlockVector
8383end
0 commit comments