1
1
using BlockArrays: blocksizes
2
- using LinearAlgebra: LinearAlgebra
2
+ using DiagonalArrays: diagonal
3
+ using LinearAlgebra: LinearAlgebra, Diagonal
3
4
using MatrixAlgebraKit:
4
5
MatrixAlgebraKit,
6
+ TruncationStrategy,
5
7
default_eig_algorithm,
6
8
default_eigh_algorithm,
9
+ diagview,
7
10
eig_full!,
11
+ eig_trunc!,
8
12
eig_vals!,
9
13
eigh_full!,
10
- eigh_vals!
14
+ eigh_trunc!,
15
+ eigh_vals!,
16
+ findtruncated
11
17
12
18
for f in [:default_eig_algorithm , :default_eigh_algorithm ]
13
19
@eval begin
@@ -37,7 +43,7 @@ for f in [:eig_full!, :eigh_full!]
37
43
for I in eachunstoredblockdiagindex (A)
38
44
# TODO : Support setting `LinearAlgebra.I` directly, and/or
39
45
# using `FillArrays.Eye`.
40
- V[I] = LinearAlgebra. I (first ( blocksizes (A)[ Int .( Tuple (I)) ... ] ))
46
+ V[I] = LinearAlgebra. I (size ( @view (V[I]), 1 ))
41
47
end
42
48
return (D, V)
43
49
end
@@ -61,3 +67,24 @@ for f in [:eig_vals!, :eigh_vals!]
61
67
end
62
68
end
63
69
end
70
+
71
+ const TBlockDV = Tuple{AbstractBlockSparseMatrix,AbstractBlockSparseMatrix}
72
+
73
+ for f in [:eig_trunc! , :eigh_trunc! ]
74
+ @eval begin
75
+ function MatrixAlgebraKit. truncate! (
76
+ :: typeof ($ f), (D, V):: TBlockDV , strategy:: TruncationStrategy
77
+ )
78
+ return MatrixAlgebraKit. truncate! (
79
+ $ f, (D, V), BlockPermutedDiagonalTruncationStrategy (strategy)
80
+ )
81
+ end
82
+ function MatrixAlgebraKit. truncate! (
83
+ :: typeof ($ f), (D, V):: TBlockDV , strategy:: BlockPermutedDiagonalTruncationStrategy
84
+ )
85
+ d = diagview (D)
86
+ ind = findtruncated (d, strategy)
87
+ return diagonal (d[ind]), V[:, ind]
88
+ end
89
+ end
90
+ end
0 commit comments