@@ -12,8 +12,8 @@ struct HelmholtzHodge{T}
12
12
end
13
13
14
14
function HelmholtzHodge (:: Type{T} , N:: Int ) where T
15
- Q = Vector {BandedQ{T}} (N)
16
- R = Vector {BandedMatrix{T, Matrix{T}}} (N)
15
+ Q = Vector {BandedQ{T}} (undef, N)
16
+ R = Vector {BandedMatrix{T, Matrix{T}}} (undef, N)
17
17
for m = 1 : N
18
18
Q[m], R[m] = qr (helmholtzhodgeconversion (T, N, m))
19
19
end
@@ -74,7 +74,7 @@ function solvehelmholtzhodge!(HH::HelmholtzHodge{T}, m::Int) where T
74
74
75
75
# Step 2: backsolve with (square) R
76
76
77
- BandedMatrices . tbsv ! (' U' , ' N' , ' N' , size (R. data, 2 ), R. u, pointer (R. data), size (R. data, 1 ), pointer (X), 1 )
77
+ hhdtbsv ! (' U' , ' N' , ' N' , size (R. data, 2 ), R. u, pointer (R. data), size (R. data, 1 ), pointer (X), 1 )
78
78
79
79
X
80
80
end
@@ -170,3 +170,29 @@ function writeout2!(X, U1, U2, N, m)
170
170
end
171
171
U2[N- m, 2 m] = - X[2 N- 2 m]
172
172
end
173
+
174
+ using LinearAlgebra
175
+ import LinearAlgebra: BlasInt
176
+ import LinearAlgebra. BLAS: libblas, @blasfunc
177
+
178
+ for (fname, elty) in ((:dtbsv_ ,:Float64 ),
179
+ (:stbsv_ ,:Float32 ),
180
+ (:ztbsv_ ,:ComplexF64 ),
181
+ (:ctbsv_ ,:ComplexF32 ))
182
+ @eval begin
183
+ function hhdtbsv! (uplo:: Char , trans:: Char , diag:: Char ,
184
+ n:: Int , k:: Int , A:: Ptr{$elty} , lda:: Int ,
185
+ x:: Ptr{$elty} , incx:: Int )
186
+ ccall ((@blasfunc ($ fname), libblas), Nothing,
187
+ (Ref{UInt8}, Ref{UInt8}, Ref{UInt8},
188
+ Ref{BlasInt}, Ref{BlasInt},
189
+ Ptr{$ elty}, Ref{BlasInt},
190
+ Ptr{$ elty}, Ref{BlasInt}),
191
+ uplo, trans, diag,
192
+ n, k,
193
+ A, lda,
194
+ x, incx)
195
+ x
196
+ end
197
+ end
198
+ end
0 commit comments