Update test_BLAS.jl

JuliaSparse · May 14, 2024 · 80ec792 · 80ec792
1 parent 16f294a
commit 80ec792
Showing 1 changed file with 84 additions and 55 deletions.
diff --git a/test/test_BLAS.jl b/test/test_BLAS.jl
@@ -25,74 +25,103 @@ macro test_blas(ex)
 end
 
 @testset "matrix-vector multiplication (non-square)" begin
-    for i = 1:5
-        a = sprand(10, 5, 0.5)
-        b = rand(5)
-        @test_blas maximum(abs.(a*b - Array(a)*b)) < 100*eps()
-        b = rand(5, 5)
-        @test_blas maximum(abs.(a*b - Array(a)*b)) < 100*eps()
-        b = rand(10)
-        @test_blas maximum(abs.(a'*b - Array(a)'*b)) < 100*eps()
-        @test_blas maximum(abs.(transpose(a)*b - Array(a)'*b)) < 100*eps()
-        b = rand(10,10)
-        @test_blas maximum(abs.(a'*b - Array(a)'*b)) < 100*eps()
-        @test_blas maximum(abs.(transpose(a)*b - Array(a)'*b)) < 100*eps()
+    T = Float64
+    @testset "interface = $interface" for interface in ("LP64", "ILP64")
+        S = interface == "LP64" ? Int32 : Int64
+        for i = 1:5
+            a = sprand(T, 10, 5, 0.5)
+            a = SparseMatrixCSC{T, S}(a)
+            b = rand(T, 5)
+            @test_blas maximum(abs.(a*b - Array(a)*b)) < 100*eps(T)
+            b = rand(T, 5, 5)
+            @test_blas maximum(abs.(a*b - Array(a)*b)) < 100*eps(T)
+            b = rand(T, 10)
+            @test_blas maximum(abs.(a'*b - Array(a)'*b)) < 100*eps(T)
+            @test_blas maximum(abs.(transpose(a)*b - Array(a)'*b)) < 100*eps(T)
+            b = rand(T, 10, 10)
+            @test_blas maximum(abs.(a'*b - Array(a)'*b)) < 100*eps()
+            @test_blas maximum(abs.(transpose(a)*b - Array(a)'*b)) < 100*eps(T)
+        end
     end
 end
 
 #?
 @testset "complex matrix-vector multiplication" begin
-    for i = 1:5
-        a = I + im * 0.1*sprandn(5, 5, 0.2)
-        b = randn(5,3) + im*randn(5,3)
-        c = randn(5) + im*randn(5)
-        d = randn(5) + im*randn(5)
-        α = rand(ComplexF64)
-        β = rand(ComplexF64)
-        @test_blas (maximum(abs.(a*b - Array(a)*b)) < 100*eps())
-        @test_blas (maximum(abs.(a'*b - Array(a)'*b)) < 100*eps())
-        @test_blas (maximum(abs.(transpose(a)*b - transpose(Array(a))*b)) < 100*eps())
-        @test_blas (maximum(abs.(mul!(similar(b), a, b) - Array(a)*b)) < 100*eps())
-        @test_blas (maximum(abs.(mul!(similar(c), a, c) - Array(a)*c)) < 100*eps())
-        @test_blas (maximum(abs.(mul!(similar(b), transpose(a), b) - transpose(Array(a))*b)) < 100*eps())
-        @test_blas (maximum(abs.(mul!(similar(c), transpose(a), c) - transpose(Array(a))*c)) < 100*eps())
-        @test_blas (maximum(abs.(mul!(copy(b), a, b, α, β) - (α*(Array(a)*b) + β*b))) < 100*eps())
-        @test_blas (maximum(abs.(mul!(copy(b), transpose(a), b, α, β) - (α*(transpose(Array(a))*b) + β*b))) < 100*eps())
-        @test_blas (maximum(abs.(mul!(copy(c), transpose(a), c, α, β) - (α*(transpose(Array(a))*c) + β*c))) < 100*eps())
-        α = β = 1 # test conversion to float
-        @test_blas (maximum(abs.(mul!(copy(b), a, b, α, β) - (α*(Array(a)*b) + β*b))) < 100*eps())
-        @test_blas (maximum(abs.(mul!(copy(b), transpose(a), b, α, β) - (α*(transpose(Array(a))*b) + β*b))) < 100*eps())
-        @test_blas (maximum(abs.(mul!(copy(c), transpose(a), c, α, β) - (α*(transpose(Array(a))*c) + β*c))) < 100*eps())
+    T = ComplexF64
+    R = Float64
+    @testset "interface = $interface" for interface in ("LP64", "ILP64")
+        S = interface == "LP64" ? Int32 : Int64
+        for i = 1:5
+            a = I + im * 0.1 * sprandn(T, 5, 5, 0.2)
+            a = SparseMatrixCSC{T,S}(a)
+            b = randn(T, 5, 3) + im * randn(T, 5, 3)
+            c = randn(T, 5) + im * randn(T, 5)
+            d = randn(T, 5) + im * randn(T, 5)
+            α = rand(T)
+            β = rand(T)
+            @test_blas (maximum(abs.(a*b - Array(a)*b)) < 100*eps(R))
+            @test_blas (maximum(abs.(a'*b - Array(a)'*b)) < 100*eps(R))
+            @test_blas (maximum(abs.(transpose(a)*b - transpose(Array(a))*b)) < 100*eps(R))
+            @test_blas (maximum(abs.(mul!(similar(b), a, b) - Array(a)*b)) < 100*eps(R))
+            @test_blas (maximum(abs.(mul!(similar(c), a, c) - Array(a)*c)) < 100*eps(R))
+            @test_blas (maximum(abs.(mul!(similar(b), transpose(a), b) - transpose(Array(a))*b)) < 100*eps(R))
+            @test_blas (maximum(abs.(mul!(similar(c), transpose(a), c) - transpose(Array(a))*c)) < 100*eps(R))
+            @test_blas (maximum(abs.(mul!(copy(b), a, b, α, β) - (α*(Array(a)*b) + β*b))) < 100*eps(R))
+            @test_blas (maximum(abs.(mul!(copy(b), transpose(a), b, α, β) - (α*(transpose(Array(a))*b) + β*b))) < 100*eps(R))
+            @test_blas (maximum(abs.(mul!(copy(c), transpose(a), c, α, β) - (α*(transpose(Array(a))*c) + β*c))) < 100*eps(R))
+            α = β = 1 # test conversion to float
+            @test_blas (maximum(abs.(mul!(copy(b), a, b, α, β) - (α*(Array(a)*b) + β*b))) < 100*eps(R))
+            @test_blas (maximum(abs.(mul!(copy(b), transpose(a), b, α, β) - (α*(transpose(Array(a))*b) + β*b))) < 100*eps(R))
+            @test_blas (maximum(abs.(mul!(copy(c), transpose(a), c, α, β) - (α*(transpose(Array(a))*c) + β*c))) < 100*eps(R))
 
-        c = randn(6) + im*randn(6)
-        @test_throws DimensionMismatch transpose(a)*c
-        @test_throws DimensionMismatch a.*c
-        @test_throws DimensionMismatch a.*c
+            c = randn(T, 6) + im * randn(T, 6)
+            @test_throws DimensionMismatch transpose(a)*c
+            @test_throws DimensionMismatch a.*c
+            @test_throws DimensionMismatch a.*c
+        end
     end
 end
 
 @testset "triangular" begin
-    n = 100
-    A = sprandn(n, n, 0.5) + sqrt(n)*I
-    b = rand(n)
-    symA = A + transpose(A)
-    trilA = tril(A)
-    triuA = triu(A)
-    trilUA = tril(A, -1) + I
-    triuUA = triu(A, 1)  + I
+    @testset "T = $T" for T in (Float32, Float64, ComplexF32, ComplexF64)
+        @testset "interface = $interface" for interface in ("LP64", "ILP64")
+            S = interface == "LP64" ? Int32 : Int64
+            n = 100
+            A = sprandn(T, n, n, 0.5) + sqrt(n)*I
+            A = SparseMatrixCSC{T,S}(A)
+            b = rand(T, n)
+            trilA = tril(A)
+            triuA = triu(A)
+            trilUA = tril(A, -1) + I
+            triuUA = triu(A, 1)  + I
 
-    @test_blas LowerTriangular(trilA) \ b ≈ Array(LowerTriangular(trilA)) \ b
-    @test_blas LowerTriangular(trilA) * b ≈ Array(LowerTriangular(trilA)) * b
+            @test_blas LowerTriangular(trilA) \ b ≈ Array(LowerTriangular(trilA)) \ b
+            @test_blas LowerTriangular(trilA) * b ≈ Array(LowerTriangular(trilA)) * b
 
-    @test_blas UpperTriangular(triuA) \ b ≈ Array(UpperTriangular(triuA)) \ b
-    @test_blas UpperTriangular(triuA) * b ≈ Array(UpperTriangular(triuA)) * b
+            @test_blas UpperTriangular(triuA) \ b ≈ Array(UpperTriangular(triuA)) \ b
+            @test_blas UpperTriangular(triuA) * b ≈ Array(UpperTriangular(triuA)) * b
 
-    @test_blas UnitLowerTriangular(trilUA) \ b ≈ Array(UnitLowerTriangular(trilUA)) \ b
-    @test_blas UnitLowerTriangular(trilUA) * b ≈ Array(UnitLowerTriangular(trilUA)) * b
+            @test_blas UnitLowerTriangular(trilUA) \ b ≈ Array(UnitLowerTriangular(trilUA)) \ b
+            @test_blas UnitLowerTriangular(trilUA) * b ≈ Array(UnitLowerTriangular(trilUA)) * b
 
-    @test_blas UnitUpperTriangular(triuUA) \ b ≈ Array(UnitUpperTriangular(triuUA)) \ b
-    @test_blas UnitUpperTriangular(triuUA) * b ≈ Array(UnitUpperTriangular(triuUA)) * b
+            @test_blas UnitUpperTriangular(triuUA) \ b ≈ Array(UnitUpperTriangular(triuUA)) \ b
+            @test_blas UnitUpperTriangular(triuUA) * b ≈ Array(UnitUpperTriangular(triuUA)) * b
+        end
+    end
+end
 
-    @test_blas Symmetric(symA) * b ≈ Array(Symmetric(symA)) * b
-    @test_blas Hermitian(symA) * b ≈ Array(Hermitian(symA)) * b
+@testset "Symmetric -- Hermitian" begin
+    @testset "T = $T" for T in (Float32, Float64, ComplexF32, ComplexF64)
+        @testset "interface = $interface" for interface in ("LP64", "ILP64")
+            S = interface == "LP64" ? Int32 : Int64
+            n = 100
+            A = sprandn(T, n, n, 0.5) + sqrt(n)*I
+            b = rand(T, n)
+            symA = A + transpose(A)
+            hermA = A + adjoint(A)
+
+            @test_blas Symmetric(symA) * b ≈ Array(Symmetric(symA)) * b
+            @test_blas Hermitian(hermA) * b ≈ Array(Hermitian(hermA)) * b
+        end
+    end
 end