From cb822c23569d9a243ea0ce17c4e772992f7f9ab7 Mon Sep 17 00:00:00 2001 From: "Viral B. Shah" Date: Fri, 16 Oct 2026 03:52:58 +0000 Subject: [PATCH 1/4] Precompile `show`, `\`, products, broadcast and concatenation of sparse Float64 arrays Co-Authored-By: Claude Opus 5.5 --- src/SparseArrays.jl | 2 ++ src/precompile.jl | 24 ++++++++++++++++++++++++ test/solvers/threads.jl | 18 ++++++++++++++++++ 3 files changed, 44 insertions(+) create mode 100644 src/precompile.jl diff --git a/src/SparseArrays.jl b/src/SparseArrays.jl index 22733cd9..77e456df 100644 --- a/src/SparseArrays.jl +++ b/src/SparseArrays.jl @@ -486,4 +486,6 @@ zero(a::AbstractSparseArray{Tv,Ti}) where {Tv,Ti} = spzeros(Tv, Ti, size(a)...) LinearAlgebra.diagzero(D::Diagonal{<:AbstractSparseMatrix{T}},i,j) where {T} = spzeros(T, size(D.diag[i], 1), size(D.diag[j], 2)) +include("precompile.jl") + end diff --git a/src/precompile.jl b/src/precompile.jl new file mode 100644 index 00000000..4e18e66e --- /dev/null +++ b/src/precompile.jl @@ -0,0 +1,24 @@ +# This file is a part of Julia. License is MIT: https://julialang.org/license + +# Runs the common operations on Float64 arrays while the package image is built, so that their +# native code is cached; `show` and the `\` polyalgorithm otherwise take seconds on first call. +if Base.generating_output() + let + A = sparse([1, 2, 3, 1, 3], [1, 2, 3, 3, 1], [4.0, 5.0, 6.0, 1.0, 1.0]) + x = [1.0, 2.0, 3.0] + v = sparsevec([1, 3], [1.0, 2.0], 3) + A[1, 2]; A[1:2, 2:3]; A[:, 2]; A[2, 2] = 5.0 + A * x; A * A; A' * x; A * v; A + A; A - A; 2 * A + A .* 2; A .+ A; A .* A + copy(A'); sum(A); findnz(A); dropzeros(A); diag(A); kron(A, A) + hcat(A, A); vcat(A, A); [A A; A A] + v + v; dot(v, v) + sprint(show, MIME"text/plain"(), A); sprint(show, A) + sprint(show, MIME"text/plain"(), v); sprint(show, v) + @static if Base.USE_GPL_LIBS + A \ x; lu(A) \ x; qr(A) \ x + H = Hermitian(A' * A + 10 * LinearAlgebra.I) + cholesky(H) \ x; ldlt(H) \ x + end + end +end diff --git a/test/solvers/threads.jl b/test/solvers/threads.jl index 0adc73f3..897ac07c 100644 --- a/test/solvers/threads.jl +++ b/test/solvers/threads.jl @@ -41,6 +41,24 @@ end end end +@static if COMPREHENSIVE +@testset "the package image holds native code for `show` and `\\`" begin + script = """ + using LinearAlgebra + function cached(f, types) + mi = Base.method_instance(f, types) + return mi !== nothing && isdefined(mi, :cache) && mi.cache.invoke != C_NULL + end + S = SparseMatrixCSC{Float64,Int} + if Base.JLOptions().use_pkgimages != 0 + @test cached(show, (typeof(IOBuffer()), MIME"text/plain", S)) + @test cached(\\, (S, Vector{Float64})) + end + """ + @test success(pipeline(testprocess(script); stdout, stderr)) +end +end + @testset "SuiteSparse library directory override (#250)" begin L = SparseArrays.LibSuiteSparse original_dir = L.libdir() From b3941bf041e65d455cd6ce92cb487c5a4cd83074 Mon Sep 17 00:00:00 2001 From: "Viral B. Shah" Date: Tue, 27 Oct 2026 00:42:26 -0400 Subject: [PATCH 2/4] Precompile `sprand` of a Float64 matrix Co-Authored-By: Claude Opus 5.5 --- src/precompile.jl | 1 + test/solvers/threads.jl | 1 + 2 files changed, 2 insertions(+) diff --git a/src/precompile.jl b/src/precompile.jl index 4e18e66e..4ea9f3de 100644 --- a/src/precompile.jl +++ b/src/precompile.jl @@ -13,6 +13,7 @@ if Base.generating_output() copy(A'); sum(A); findnz(A); dropzeros(A); diag(A); kron(A, A) hcat(A, A); vcat(A, A); [A A; A A] v + v; dot(v, v) + sprand(3, 3, 0.5) sprint(show, MIME"text/plain"(), A); sprint(show, A) sprint(show, MIME"text/plain"(), v); sprint(show, v) @static if Base.USE_GPL_LIBS diff --git a/test/solvers/threads.jl b/test/solvers/threads.jl index 897ac07c..4052add6 100644 --- a/test/solvers/threads.jl +++ b/test/solvers/threads.jl @@ -53,6 +53,7 @@ end if Base.JLOptions().use_pkgimages != 0 @test cached(show, (typeof(IOBuffer()), MIME"text/plain", S)) @test cached(\\, (S, Vector{Float64})) + @test cached(sprand, (Int, Int, Float64)) end """ @test success(pipeline(testprocess(script); stdout, stderr)) From 799f843e46d54f031403c11746d3b4fd688c5d7c Mon Sep 17 00:00:00 2001 From: "Viral B. Shah" Date: Tue, 3 Nov 2026 15:55:14 +0000 Subject: [PATCH 3/4] Broaden the precompile workload to construction, slicing, reductions and factorizations Covers `spdiagm`, `spzeros`, `sparse(I, m, n)`, `mul!`, row and logical indexing, reductions with `dims`, `norm`, `findall`/`findmax`, `triu`/`tril`, `blockdiag`, `permute`, mixed sparse-dense `hcat`, products with `Symmetric` and `Diagonal`, triangular solves, and reuse of `lu`/`cholesky` factorizations, still on Float64 with `Int` indices only. Co-Authored-By: Claude Opus 5.5 --- src/precompile.jl | 28 ++++++++++++++++++++-------- 1 file changed, 20 insertions(+), 8 deletions(-) diff --git a/src/precompile.jl b/src/precompile.jl index 4ea9f3de..1cf82c41 100644 --- a/src/precompile.jl +++ b/src/precompile.jl @@ -6,20 +6,32 @@ if Base.generating_output() let A = sparse([1, 2, 3, 1, 3], [1, 2, 3, 3, 1], [4.0, 5.0, 6.0, 1.0, 1.0]) x = [1.0, 2.0, 3.0] + y = similar(x) + X = [1.0 2.0; 3.0 4.0; 5.0 6.0] v = sparsevec([1, 3], [1.0, 2.0], 3) - A[1, 2]; A[1:2, 2:3]; A[:, 2]; A[2, 2] = 5.0 - A * x; A * A; A' * x; A * v; A + A; A - A; 2 * A - A .* 2; A .+ A; A .* A - copy(A'); sum(A); findnz(A); dropzeros(A); diag(A); kron(A, A) - hcat(A, A); vcat(A, A); [A A; A A] - v + v; dot(v, v) + spzeros(3, 3); sparse(1.0I, 3, 3); sparse(X); sparse(x) + spdiagm(-1 => [1.0, 1.0], 0 => x, 1 => [1.0, 1.0]) + A[1, 2]; A[1:2, 2:3]; A[:, 2]; A[2, :]; A[[3, 1], [1, 3]]; A[A .> 2.5]; A[2, 2] = 5.0 + view(A, :, 2:3) * [1.0, 1.0] + A * x; A * X; A * A; A' * x; transpose(A) * x; A * v; A + A; A - A; 2 * A + mul!(y, A, x); mul!(y, A, x, 2.0, 1.0) + A .* 2; A .+ A; A .* A; A .+ x; abs.(A) + copy(A'); sparse(A'); sum(A); sum(A; dims = 1); sum(A; dims = 2); maximum(A) + norm(A); opnorm(A, 1); tr(A); findnz(A); findall(!iszero, A); findmax(A) + dropzeros(A); droptol!(copy(A), 1.5); diag(A); triu(A); tril(A); issymmetric(A) + kron(A, A); blockdiag(A, A); permute(A, [3, 2, 1], [1, 2, 3]) + hcat(A, A); vcat(A, A); [A A; A A]; hcat(A, X) + Symmetric(A) * x; Diagonal(x) * A; A * Diagonal(x) + LowerTriangular(A) \ x; UpperTriangular(A) \ x + v + v; v .* v; dot(v, v); norm(v); v[2:3]; Vector(v) sprand(3, 3, 0.5) sprint(show, MIME"text/plain"(), A); sprint(show, A) sprint(show, MIME"text/plain"(), v); sprint(show, v) @static if Base.USE_GPL_LIBS - A \ x; lu(A) \ x; qr(A) \ x + A \ x; qr(A) \ x + F = lu(A); F \ x; lu!(F, A); det(F) H = Hermitian(A' * A + 10 * LinearAlgebra.I) - cholesky(H) \ x; ldlt(H) \ x + C = cholesky(H); C \ x; logdet(C); ldlt(H) \ x end end end From 281766173890daad1d006027d8b839e4fa0c24eb Mon Sep 17 00:00:00 2001 From: "Viral B. Shah" Date: Tue, 3 Nov 2026 22:26:41 +0000 Subject: [PATCH 4/4] Precompile the same workload for ComplexF64 Each eltype's block is `@eval`ed rather than run in a loop: a top-level loop is compiled, which inlines small methods such as `sprand` and leaves no native code for them in the package image. Co-Authored-By: Claude Opus 5.5 --- src/precompile.jl | 69 ++++++++++++++++++++++++++--------------------- 1 file changed, 38 insertions(+), 31 deletions(-) diff --git a/src/precompile.jl b/src/precompile.jl index 1cf82c41..9d1e22f1 100644 --- a/src/precompile.jl +++ b/src/precompile.jl @@ -1,37 +1,44 @@ # This file is a part of Julia. License is MIT: https://julialang.org/license -# Runs the common operations on Float64 arrays while the package image is built, so that their -# native code is cached; `show` and the `\` polyalgorithm otherwise take seconds on first call. +# Runs the common operations on Float64 and ComplexF64 arrays while the package image is built, +# so that their native code is cached; `show` and the `\` polyalgorithm otherwise take seconds +# on first call. Each eltype's block runs in the interpreter, without a loop, so that it +# calls every method dynamically as the REPL does: compiled code would inline small methods +# such as `sprand`, which then get no native code of their own. if Base.generating_output() - let - A = sparse([1, 2, 3, 1, 3], [1, 2, 3, 3, 1], [4.0, 5.0, 6.0, 1.0, 1.0]) - x = [1.0, 2.0, 3.0] - y = similar(x) - X = [1.0 2.0; 3.0 4.0; 5.0 6.0] - v = sparsevec([1, 3], [1.0, 2.0], 3) - spzeros(3, 3); sparse(1.0I, 3, 3); sparse(X); sparse(x) - spdiagm(-1 => [1.0, 1.0], 0 => x, 1 => [1.0, 1.0]) - A[1, 2]; A[1:2, 2:3]; A[:, 2]; A[2, :]; A[[3, 1], [1, 3]]; A[A .> 2.5]; A[2, 2] = 5.0 - view(A, :, 2:3) * [1.0, 1.0] - A * x; A * X; A * A; A' * x; transpose(A) * x; A * v; A + A; A - A; 2 * A - mul!(y, A, x); mul!(y, A, x, 2.0, 1.0) - A .* 2; A .+ A; A .* A; A .+ x; abs.(A) - copy(A'); sparse(A'); sum(A); sum(A; dims = 1); sum(A; dims = 2); maximum(A) - norm(A); opnorm(A, 1); tr(A); findnz(A); findall(!iszero, A); findmax(A) - dropzeros(A); droptol!(copy(A), 1.5); diag(A); triu(A); tril(A); issymmetric(A) - kron(A, A); blockdiag(A, A); permute(A, [3, 2, 1], [1, 2, 3]) - hcat(A, A); vcat(A, A); [A A; A A]; hcat(A, X) - Symmetric(A) * x; Diagonal(x) * A; A * Diagonal(x) - LowerTriangular(A) \ x; UpperTriangular(A) \ x - v + v; v .* v; dot(v, v); norm(v); v[2:3]; Vector(v) - sprand(3, 3, 0.5) - sprint(show, MIME"text/plain"(), A); sprint(show, A) - sprint(show, MIME"text/plain"(), v); sprint(show, v) - @static if Base.USE_GPL_LIBS - A \ x; qr(A) \ x - F = lu(A); F \ x; lu!(F, A); det(F) - H = Hermitian(A' * A + 10 * LinearAlgebra.I) - C = cholesky(H); C \ x; logdet(C); ldlt(H) \ x + for T in (Float64, ComplexF64) + @eval let T = $T + A = sparse([1, 2, 3, 1, 3], [1, 2, 3, 3, 1], T[4, 5, 6, 1, 1]) + x = T[1, 2, 3] + y = similar(x) + X = T[1 2; 3 4; 5 6] + v = sparsevec([1, 3], T[1, 2], 3) + spzeros(T, 3, 3); sparse(one(T) * I, 3, 3); sparse(X); sparse(x) + spdiagm(-1 => T[1, 1], 0 => x, 1 => T[1, 1]) + A[1, 2]; A[1:2, 2:3]; A[:, 2]; A[2, :]; A[[3, 1], [1, 3]]; A[2, 2] = 5 + view(A, :, 2:3) * T[1, 1] + A * x; A * X; A * A; A' * x; transpose(A) * x; A * v; A + A; A - A; 2 * A + mul!(y, A, x); mul!(y, A, x, T(2), T(1)) + A .* 2; A .+ A; A .* A; A .+ x; abs.(A) + copy(A'); sparse(A'); sum(A); sum(A; dims = 1); sum(A; dims = 2) + norm(A); opnorm(A, 1); tr(A); findnz(A); findall(!iszero, A) + dropzeros(A); droptol!(copy(A), 1.5); diag(A); triu(A); tril(A); issymmetric(A) + kron(A, A); blockdiag(A, A); permute(A, [3, 2, 1], [1, 2, 3]) + hcat(A, A); vcat(A, A); [A A; A A]; hcat(A, X) + Symmetric(A) * x; Diagonal(x) * A; A * Diagonal(x) + LowerTriangular(A) \ x; UpperTriangular(A) \ x + v + v; v .* v; dot(v, v); norm(v); v[2:3]; Vector(v) + sprint(show, MIME"text/plain"(), A); sprint(show, A) + sprint(show, MIME"text/plain"(), v); sprint(show, v) + if T <: Real + A[A .> 2.5]; maximum(A); findmax(A); sprand(3, 3, 0.5) + end + @static if Base.USE_GPL_LIBS + A \ x; qr(A) \ x + F = lu(A); F \ x; lu!(F, A); det(F) + H = Hermitian(A' * A + 10 * LinearAlgebra.I) + C = cholesky(H); C \ x; logdet(C); ldlt(H) \ x + end end end end