From 2269b4f0eb675209d680027227cfb555929b83da Mon Sep 17 00:00:00 2001 From: Alexis Montoison Date: Sun, 7 Sep 2025 12:32:51 -0500 Subject: [PATCH 01/18] Add an extension for oneAPI.jl --- lib/MadNLPGPU/Project.toml | 5 +- .../MadNLPGPUOneAPIExt/MadNLPGPUOneAPIExt.jl | 25 + .../ext/MadNLPGPUOneAPIExt/oneapi.jl | 47 ++ .../ext/MadNLPGPUOneAPIExt/oneapi_dense.jl | 155 ++++++ .../ext/MadNLPGPUOneAPIExt/oneapi_sparse.jl | 495 ++++++++++++++++++ .../ext/MadNLPGPUOneAPIExt/onemkl.jl | 321 ++++++++++++ 6 files changed, 1047 insertions(+), 1 deletion(-) create mode 100644 lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/MadNLPGPUOneAPIExt.jl create mode 100644 lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/oneapi.jl create mode 100644 lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/oneapi_dense.jl create mode 100644 lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/oneapi_sparse.jl create mode 100644 lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/onemkl.jl diff --git a/lib/MadNLPGPU/Project.toml b/lib/MadNLPGPU/Project.toml index c5c1ddabb..7c937e5c2 100644 --- a/lib/MadNLPGPU/Project.toml +++ b/lib/MadNLPGPU/Project.toml @@ -14,9 +14,11 @@ SparseArrays = "2f01184e-e22b-5df5-ae63-d93ebab69eaf" [weakdeps] AMDGPU = "21141c5a-9bdb-4563-92ae-f87d6854732e" +oneAPI = "8f75cd03-7ff8-4ecb-9b8f-daf728133b1b" [extensions] MadNLPGPUAMDGPUExt = "AMDGPU" +MadNLPGPUOneAPIExt = "oneAPI" [compat] AMD = "0.5" @@ -27,6 +29,7 @@ KernelAbstractions = "0.9" MadNLP = "0.8.12" MadNLPTests = "0.5.3" Metis = "1" +oneAPI = "2.2.0" julia = "1.10" [extras] @@ -34,4 +37,4 @@ MadNLPTests = "b52a2a03-04ab-4a5f-9698-6a2deff93217" Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40" [targets] -test = ["Test", "MadNLPTests", "AMDGPU"] +test = ["Test", "MadNLPTests", "AMDGPU", "oneAPI"] diff --git a/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/MadNLPGPUOneAPIExt.jl b/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/MadNLPGPUOneAPIExt.jl new file mode 100644 index 000000000..d0f9d5bc9 --- /dev/null +++ b/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/MadNLPGPUOneAPIExt.jl @@ -0,0 +1,25 @@ +module MadNLPGPUAMDGPUExt + +import LinearAlgebra +import SparseArrays: SparseMatrixCSC, nonzeros, nnz +import LinearAlgebra: Symmetric + +import MadNLP +import MadNLPGPU + +import KernelAbstractions: synchronize + +using oneAPI +using oneAPI.oneMKL, onAPI.Support + +function __init__() + setglobal!(MadNLPGPU, :LapackOneMKLSolver, LapackOneMKLSolver) + return +end + +include("oneapi_dense.jl") +include("oneapi_sparse.jl") +include("onemkl.jl") +include("oneapi.jl") + +end diff --git a/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/oneapi.jl b/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/oneapi.jl new file mode 100644 index 000000000..cfeaf2d91 --- /dev/null +++ b/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/oneapi.jl @@ -0,0 +1,47 @@ +#= + MadNLP.MadNLPOptions +=# + +function MadNLP.MadNLPOptions{T}( + nlp::MadNLP.AbstractNLPModel{T,VT}; + dense_callback = MadNLP.is_dense_callback(nlp), + callback = dense_callback ? MadNLP.DenseCallback : MadNLP.SparseCallback, + kkt_system = dense_callback ? MadNLP.DenseCondensedKKTSystem : MadNLP.SparseCondensedKKTSystem, + linear_solver = MadNLPGPU.LapackOneMKLSolver, + tol = MadNLP.get_tolerance(T,kkt_system), + bound_relax_factor = tol, +) where {T, VT <: oneVector{T}} + return MadNLP.MadNLPOptions{T}( + tol = tol, + callback = callback, + kkt_system = kkt_system, + linear_solver = linear_solver, + bound_relax_factor = bound_relax_factor, + ) +end + +#= + SparseMatrixCSC to oneSparseMatrixCSC +=# + +function oneMKL.oneSparseMatrixCSC{Tv,Ti}(A::SparseMatrixCSC{Tv,Ti}) where {Tv,Ti} + return oneMKL.oneSparseMatrixCSC{Tv,Ti}( + oneVector(A.colptr), + oneVector(A.rowval), + oneVector(A.nzval), + size(A), + ) +end + +#= + oneSparseMatrixCSC to oneMatrix +=# + +function MadNLPGPU.gpu_transfer!(y::oneMatrix{T}, x::oneMKL.oneSparseMatrixCSC{T}) where {T} + n = size(y, 2) + fill!(y, zero(T)) + backend = OneAPIBackend() + MadNLPGPU._csc_to_dense_kernel!(backend)(y, x.colPtr, x.rowVal, x.nzVal, ndrange = n) + synchronize(backend) + return +end diff --git a/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/oneapi_dense.jl b/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/oneapi_dense.jl new file mode 100644 index 000000000..61d908e36 --- /dev/null +++ b/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/oneapi_dense.jl @@ -0,0 +1,155 @@ +######################################################################## +##### oneAPI wrappers for DenseKKTSystem / DenseCondensedKKTSystem ##### +######################################################################## + +#= + MadNLP.symul! +=# + +MadNLP.symul!(y, A, x::oneVector{T}, α = one(T), β = zero(T)) where T = oneMKL.symv!('L', T(α), A, x, T(β), y) + +#= + MadNLP._ger! +=# + +MadNLP._ger!(alpha::Number, x::oneVector{T}, y::oneVector{T}, A::oneMatrix{T}) where T = oneMKL.ger!(alpha, x, y, A) + +#= + MadNLP._madnlp_unsafe_wrap +=# + +function MadNLP._madnlp_unsafe_wrap(vec::VT, n, shift=1) where {T, VT <: oneVector{T}} + return view(vec,shift:shift+n-1) +end + +#= + MadNLP.diag! +=# + +function MadNLP.diag!(dest::oneVector{T}, src::oneMatrix{T}) where {T} + @assert length(dest) == size(src, 1) + backend = OneAPIBackend() + MadNLPGPU._copy_diag_kernel!(backend)(dest, src, ndrange = length(dest)) + synchronize(backend) + return +end + +#= + MadNLP.diag_add! +=# + +function MadNLP.diag_add!(dest::oneMatrix, src1::oneVector, src2::oneVector) + backend = OneAPIBackend() + MadNLPGPU._add_diagonal_kernel!(backend)(dest, src1, src2, ndrange = size(dest, 1)) + synchronize(backend) + return +end + +#= + MadNLP._set_diag! +=# + +function MadNLP._set_diag!(A::oneMatrix, inds, a) + if !isempty(inds) + backend = OneAPIBackend() + MadNLPGPU._set_diag_kernel!(backend)(A, inds, a; ndrange = length(inds)) + synchronize(backend) + end + return +end + +#= + MadNLP._build_dense_kkt_system! +=# + +function MadNLP._build_dense_kkt_system!( + dest::oneMatrix, + hess::oneMatrix, + jac::oneMatrix, + pr_diag::oneVector, + du_diag::oneVector, + diag_hess::oneVector, + ind_ineq::AbstractVector, + n, + m, + ns, +) + ind_ineq_gpu = oneVector(ind_ineq) + ndrange = (n + m + ns, n) + backend = OneAPIBackend() + MadNLPGPU._build_dense_kkt_system_kernel!(backend)( + dest, + hess, + jac, + pr_diag, + du_diag, + diag_hess, + ind_ineq_gpu, + n, + m, + ns, + ndrange = ndrange, + ) + synchronize(backend) + return +end + +#= + MadNLP._build_ineq_jac! +=# + +function MadNLP._build_ineq_jac!( + dest::oneMatrix, + jac::oneMatrix, + diag_buffer::oneVector, + ind_ineq::AbstractVector, + n, + m_ineq, +) + (m_ineq == 0) && return # nothing to do if no ineq. constraints + ind_ineq_gpu = oneVector(ind_ineq) + ndrange = (m_ineq, n) + backend = OneAPIBackend() + MadNLPGPU._build_jacobian_condensed_kernel!(backend)( + dest, + jac, + diag_buffer, + ind_ineq_gpu, + m_ineq, + ndrange = ndrange, + ) + synchronize(backend) + return +end + +#= + MadNLP._build_condensed_kkt_system! +=# + +function MadNLP._build_condensed_kkt_system!( + dest::oneMatrix, + hess::oneMatrix, + jac::oneMatrix, + pr_diag::oneVector, + du_diag::oneVector, + ind_eq::AbstractVector, + n, + m_eq, +) + ind_eq_gpu = oneVector(ind_eq) + ndrange = (n + m_eq, n) + backend = OneAPIBackend() + MadNLPGPU._build_condensed_kkt_system_kernel!(backend)( + dest, + hess, + jac, + pr_diag, + du_diag, + ind_eq_gpu, + n, + m_eq, + ndrange = ndrange, + ) + synchronize(backend) + return +end diff --git a/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/oneapi_sparse.jl b/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/oneapi_sparse.jl new file mode 100644 index 000000000..6fb494ee9 --- /dev/null +++ b/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/oneapi_sparse.jl @@ -0,0 +1,495 @@ +######################################################## +##### oneAPI wrappers for SparseCondensedKKTSystem ##### +######################################################## + +#= + SparseMatrixCOO to oneMKLMatrixCSC +=# + +function MadNLP.transfer!( + dest::oneMKL.oneMKLMatrixCSC, + src::MadNLP.SparseMatrixCOO, + map, +) + return copyto!(view(dest.nzVal, map), src.V) +end + +#= + MadNLP.mul! with SparseCondensedKKTSystem + oneAPI +=# + +function MadNLP.mul!( + w::MadNLP.AbstractKKTVector{T,VT}, + kkt::MadNLP.SparseCondensedKKTSystem, + x::MadNLP.AbstractKKTVector, + alpha = one(T), + beta = zero(T), +) where {T,VT<:oneVector{T}} + n = size(kkt.hess_com, 1) + m = size(kkt.jt_csc, 2) + + # Decompose results + xx = view(MadNLP.full(x), 1:n) + xs = view(MadNLP.full(x), n+1:n+m) + xz = view(MadNLP.full(x), n+m+1:n+2*m) + + # Decompose buffers + wx = view(MadNLP.full(w), 1:n) + ws = view(MadNLP.full(w), n+1:n+m) + wz = view(MadNLP.full(w), n+m+1:n+2*m) + + MadNLP.mul!(wx, kkt.hess_com, xx, alpha, beta) + MadNLP.mul!(wx, kkt.hess_com', xx, alpha, one(T)) + MadNLP.mul!(wx, kkt.jt_csc, xz, alpha, beta) + if !isempty(kkt.ext.diag_map_to) + backend = OneAPIBackend() + MadNLPGPU._diag_operation_kernel!(backend)( + wx, + kkt.hess_com.nzVal, + xx, + alpha, + kkt.ext.diag_map_to, + kkt.ext.diag_map_fr; + ndrange = length(kkt.ext.diag_map_to), + ) + synchronize(backend) + end + + MadNLP.mul!(wz, kkt.jt_csc', xx, alpha, one(T)) + MadNLP.axpy!(-alpha, xz, ws) + MadNLP.axpy!(-alpha, xs, wz) + return MadNLP._kktmul!( + w, + x, + kkt.reg, + kkt.du_diag, + kkt.l_lower, + kkt.u_lower, + kkt.l_diag, + kkt.u_diag, + alpha, + beta, + ) +end + +function MadNLP.mul_hess_blk!( + wx::VT, + kkt::Union{MadNLP.SparseKKTSystem,MadNLP.SparseCondensedKKTSystem}, + t, +) where {T,VT<:oneVector{T}} + n = size(kkt.hess_com, 1) + wxx = @view(wx[1:n]) + tx = @view(t[1:n]) + + MadNLP.mul!(wxx, kkt.hess_com, tx, one(T), zero(T)) + MadNLP.mul!(wxx, kkt.hess_com', tx, one(T), one(T)) + if !isempty(kkt.ext.diag_map_to) + backend = OneAPIBackend() + MadNLPGPU._diag_operation_kernel!(backend)( + wxx, + kkt.hess_com.nzVal, + tx, + one(T), + kkt.ext.diag_map_to, + kkt.ext.diag_map_fr; + ndrange = length(kkt.ext.diag_map_to), + ) + synchronize(backend) + end + + fill!(@view(wx[n+1:end]), 0) + wx .+= t .* kkt.pr_diag + return +end + +function MadNLP.get_tril_to_full(csc::oneMKL.oneMKLMatrixCSC{Tv,Ti}) where {Tv,Ti} + cscind = MadNLP.SparseMatrixCSC{Int,Ti}( + Symmetric( + MadNLP.SparseMatrixCSC{Int,Ti}( + size(csc)..., + Array(csc.colPtr), + Array(csc.rowVal), + collect(1:MadNLP.nnz(csc)), + ), + :L, + ), + ) + return oneMKL.oneMKLMatrixCSC{Tv,Ti}( + oneArray(cscind.colptr), + oneArray(cscind.rowval), + oneVector{Tv}(undef, MadNLP.nnz(cscind)), + size(csc), + ), + view(csc.nzVal, oneArray(cscind.nzval)) +end + +function MadNLP.get_sparse_condensed_ext( + ::Type{VT}, + hess_com, + jptr, + jt_map, + hess_map, +) where {T,VT<:oneVector{T}} + zvals = oneVector{Int}(1:length(hess_map)) + hess_com_ptr = map((i, j) -> (i, j), hess_map, zvals) + if length(hess_com_ptr) > 0 # otherwise error is thrown + sort!(hess_com_ptr) + end + + jvals = oneVector{Int}(1:length(jt_map)) + jt_csc_ptr = map((i, j) -> (i, j), jt_map, jvals) + if length(jt_csc_ptr) > 0 # otherwise error is thrown + sort!(jt_csc_ptr) + end + + by = (i, j) -> i[1] != j[1] + jptrptr = MadNLP.getptr(jptr, by = by) + hess_com_ptrptr = MadNLP.getptr(hess_com_ptr, by = by) + jt_csc_ptrptr = MadNLP.getptr(jt_csc_ptr, by = by) + + diag_map_to, diag_map_fr = get_diagonal_mapping(hess_com.colPtr, hess_com.rowVal) + + return ( + jptrptr = jptrptr, + hess_com_ptr = hess_com_ptr, + hess_com_ptrptr = hess_com_ptrptr, + jt_csc_ptr = jt_csc_ptr, + jt_csc_ptrptr = jt_csc_ptrptr, + diag_map_to = diag_map_to, + diag_map_fr = diag_map_fr, + ) +end + +function get_diagonal_mapping(colptr, rowval) + nnz = length(rowval) + if nnz == 0 + return similar(colptr, 0), similar(colptr, 0) + end + inds1 = findall( + map( + (x, y) -> ((x <= nnz) && (x != y)), + @view(colptr[1:end-1]), + @view(colptr[2:end]) + ), + ) + if length(inds1) == 0 + return similar(rows, 0), similar(ptrs, 0) + end + ptrs = colptr[inds1] + rows = rowval[ptrs] + inds2 = findall(inds1 .== rows) + if length(inds2) == 0 + return similar(rows, 0), similar(ptrs, 0) + end + + return rows[inds2], ptrs[inds2] +end + +function MadNLP._sym_length(Jt::oneMKL.oneMKLMatrixCSC) + return mapreduce( + (x, y) -> begin + z = x - y + div(z^2 + z, 2) + end, + +, + @view(Jt.colPtr[2:end]), + @view(Jt.colPtr[1:end-1]) + ) +end + +function MadNLP._first_and_last_col(sym2::oneVector, ptr2) + AMDGPU.@allowscalar begin + first = sym2[1][2] + last = sym2[ptr2[end]][2] + end + return (first, last) +end + +MadNLP.nzval(H::oneMKL.oneMKLMatrixCSC) = H.nzVal + +function MadNLP._get_sparse_csc(dims, colptr::oneVector, rowval, nzval) + return oneMKL.oneMKLMatrixCSC(colptr, rowval, nzval, dims) +end + +function getij(idx, n) + j = ceil(Int, ((2n + 1) - sqrt((2n + 1)^2 - 8 * idx)) / 2) + i = idx - div((j - 1) * (2n - j), 2) + return (i, j) +end + +#= + MadNLP._set_colptr! +=# + +function MadNLP._set_colptr!(colptr::oneVector, ptr2, sym2, guide) + if length(ptr2) > 1 # otherwise error is thrown + backend = OneAPIBackend() + MadNLPGPU._set_colptr_kernel!(backend)( + colptr, + sym2, + ptr2, + guide; + ndrange = length(ptr2) - 1, + ) + synchronize(backend) + end + return +end + + +#= + MadNLP.tril_to_full! +=# + +function MadNLP.tril_to_full!(dense::oneMatrix{T}) where {T} + n = size(dense, 1) + backend = OneAPIBackend() + MadNLPGPU._tril_to_full_kernel!(backend)(dense; ndrange = div(n^2 + n, 2)) + synchronize(backend) + return +end + +#= + MadNLP.force_lower_triangular! +=# + +function MadNLP.force_lower_triangular!(I::oneVector{T}, J) where {T} + if !isempty(I) + backend = OneAPIBackend() + MadNLPGPU._force_lower_triangular_kernel!(backend)(I, J; ndrange = length(I)) + synchronize(backend) + end + return +end + +#= + MadNLP.coo_to_csc +=# + +function MadNLP.coo_to_csc( + coo::MadNLP.SparseMatrixCOO{T,I,VT,VI}, +) where {T,I,VT<:oneArray,VI<:oneArray} + zvals = oneVector{Int}(1:length(coo.I)) + coord = map((i, j, k) -> ((i, j), k), coo.I, coo.J, zvals) + if length(coord) > 0 + sort!(coord, lt = (((i, j), k), ((n, m), l)) -> (j, i) < (m, n)) + end + + mapptr = MadNLP.getptr(coord; by = ((x1, x2), (y1, y2)) -> x1 != y1) + + colptr = similar(coo.I, size(coo, 2) + 1) + + coord_csc = coord[@view(mapptr[1:end-1])] + + backend = OneAPIBackend() + if length(coord_csc) > 0 + MadNLPGPU._set_coo_to_colptr_kernel!(backend)( + colptr, + coord_csc, + ndrange = length(coord_csc), + ) + synchronize(backend) + else + fill!(colptr, one(Int)) + end + + rowval = map(x -> x[1][1], coord_csc) + nzval = similar(rowval, T) + + csc = oneMKL.oneMKLMatrixCSC(colptr, rowval, nzval, size(coo)) + + cscmap = similar(coo.I, Int) + if length(mapptr) > 1 + MadNLPGPU._set_coo_to_csc_map_kernel!(backend)( + cscmap, + mapptr, + coord, + ndrange = length(mapptr) - 1, + ) + synchronize(backend) + end + + return csc, cscmap +end + +#= + MadNLP.build_condensed_aug_coord! +=# + +function MadNLP.build_condensed_aug_coord!( + kkt::MadNLP.AbstractCondensedKKTSystem{T,VT,MT}, +) where {T,VT,MT<:oneMKL.oneMKLMatrixCSC{T}} + fill!(kkt.aug_com.nzVal, zero(T)) + backend = OneAPIBackend() + if length(kkt.hptr) > 0 + MadNLPGPU._transfer_hessian_kernel!(backend)( + kkt.aug_com.nzVal, + kkt.hptr, + kkt.hess_com.nzVal; + ndrange = length(kkt.hptr), + ) + synchronize(backend) + end + if length(kkt.dptr) > 0 + MadNLPGPU._transfer_hessian_kernel!(backend)( + kkt.aug_com.nzVal, + kkt.dptr, + kkt.pr_diag; + ndrange = length(kkt.dptr), + ) + synchronize(backend) + end + if length(kkt.ext.jptrptr) > 1 # otherwise error is thrown + MadNLPGPU._transfer_jtsj_kernel!(backend)( + kkt.aug_com.nzVal, + kkt.jptr, + kkt.ext.jptrptr, + kkt.jt_csc.nzVal, + kkt.diag_buffer; + ndrange = length(kkt.ext.jptrptr) - 1, + ) + synchronize(backend) + end + return +end + +#= + MadNLP.compress_hessian! / MadNLP.compress_jacobian! +=# + +function MadNLP.compress_hessian!( + kkt::MadNLP.AbstractSparseKKTSystem{T,VT,MT}, +) where {T,VT,MT<:oneMKL.oneMKLMatrixCSC{T,Int32}} + fill!(kkt.hess_com.nzVal, zero(T)) + backend = OneAPIBackend() + if length(kkt.ext.hess_com_ptrptr) > 1 + MadNLPGPU._transfer_to_csc_kernel!(backend)( + kkt.hess_com.nzVal, + kkt.ext.hess_com_ptr, + kkt.ext.hess_com_ptrptr, + kkt.hess_raw.V; + ndrange = length(kkt.ext.hess_com_ptrptr) - 1, + ) + synchronize(backend) + end + return +end + +function MadNLP.compress_jacobian!( + kkt::MadNLP.SparseCondensedKKTSystem{T,VT,MT}, +) where {T,VT,MT<:oneMKL.oneMKLMatrixCSC{T,Int32}} + fill!(kkt.jt_csc.nzVal, zero(T)) + backend = OneAPIBackend() + if length(kkt.ext.jt_csc_ptrptr) > 1 # otherwise error is thrown + MadNLPGPU._transfer_to_csc_kernel!(backend)( + kkt.jt_csc.nzVal, + kkt.ext.jt_csc_ptr, + kkt.ext.jt_csc_ptrptr, + kkt.jt_coo.V; + ndrange = length(kkt.ext.jt_csc_ptrptr) - 1, + ) + synchronize(backend) + end + return +end + +#= + MadNLP._set_con_scale_sparse! +=# + +function MadNLP._set_con_scale_sparse!( + con_scale::VT, + jac_I, + jac_buffer, +) where {T,VT<:oneVector{T}} + ind_jac = oneVector{Int}(1:length(jac_I)) + inds = map((i, j) -> (i, j), jac_I, ind_jac) + !isempty(inds) && sort!(inds) + ptr = MadNLP.getptr(inds; by = ((x1, x2), (y1, y2)) -> x1 != y1) + if length(ptr) > 1 + backend = OneAPIBackend() + MadNLPGPU._set_con_scale_sparse_kernel!(backend)( + con_scale, + ptr, + inds, + jac_I, + jac_buffer; + ndrange = length(ptr) - 1, + ) + synchronize(backend) + end + return +end + +#= + MadNLP._build_condensed_aug_symbolic_hess +=# + +function MadNLP._build_condensed_aug_symbolic_hess( + H::oneMKL.oneMKLMatrixCSC{Tv,Ti}, + sym, + sym2, +) where {Tv,Ti} + if size(H, 2) > 0 + backend = OneAPIBackend() + MadNLPGPU._build_condensed_aug_symbolic_hess_kernel!(backend)( + sym, + sym2, + H.colPtr, + H.rowVal; + ndrange = size(H, 2), + ) + synchronize(backend) + end + return +end + +#= + MadNLP._build_condensed_aug_symbolic_jt +=# + +function MadNLP._build_condensed_aug_symbolic_jt( + Jt::oneMKL.oneMKLMatrixCSC{Tv,Ti}, + sym, + sym2, +) where {Tv,Ti} + if size(Jt, 2) > 0 + _offsets = map( + (i, j) -> div((j - i)^2 + (j - i), 2), + @view(Jt.colPtr[1:end-1]), + @view(Jt.colPtr[2:end]) + ) + offsets = cumsum(_offsets) + backend = OneAPIBackend() + MadNLPGPU._build_condensed_aug_symbolic_jt_kernel!(backend)( + sym, + sym2, + Jt.colPtr, + Jt.rowVal, + offsets; + ndrange = size(Jt, 2), + ) + synchronize(backend) + end + return +end + +#= + MadNLP._build_scale_augmented_system_coo! +=# + +function MadNLP._build_scale_augmented_system_coo!(dest, src, scaling::oneArray, n, m) + backend = OneAPIBackend() + MadNLPGPU._scale_augmented_system_coo_kernel!(backend)( + dest.V, + src.I, + src.J, + src.V, + scaling, + n, + m; + ndrange = nnz(src), + ) + synchronize(backend) + return +end diff --git a/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/onemkl.jl b/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/onemkl.jl new file mode 100644 index 000000000..d3555ac21 --- /dev/null +++ b/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/onemkl.jl @@ -0,0 +1,321 @@ +mutable struct LapackOneMKLSolver{T,MT} <: MadNLP.AbstractLinearSolver{T} + A::MT + fact::oneMatrix{T} + n::Int64 + sol::oneVector{T} + tau::oneVector{T} + Λ::oneVector{T} + info::oneVector{Cint} + ipiv::oneVector{Int64} + scratchpad::oneVector{T} + scratchpad_size::Int64 + opt::MadNLP.LapackOptions + logger::MadNLP.MadNLPLogger + + function LapackOneMKLSolver( + A::MT; + option_dict::Dict{Symbol,Any} = Dict{Symbol,Any}(), + opt = MadNLP.LapackOptions(), + logger = MadNLP.MadNLPLogger(), + kwargs..., + ) where {MT<:AbstractMatrix} + MadNLP.set_options!(opt, option_dict, kwargs...) + T = eltype(A) + m,n = size(A) + @assert m == n + fact = oneMatrix{T}(undef, m, n) + sol = oneVector{T}(undef, 0) + tau = oneVector{T}(undef, 0) + Λ = oneVector{T}(undef, 0) + info = oneVector{Cint}(undef, 1) + ipiv = oneVector{Int64}(undef, 0) + solver = new{T,MT}(A, fact, n, sol, tau, Λ, info, ipiv, opt, logger) + setup!(solver) + return solver + end +end + +MadNLP.improve!(M::LapackOneMKLSolver) = false +MadNLP.is_inertia(M::LapackOneMKLSolver) = (M.opt.lapack_algorithm == MadNLP.CHOLESKY) || (M.opt.lapack_algorithm == MadNLP.EVD) +function MadNLP.inertia(M::LapackOneMKLSolver) + if M.opt.lapack_algorithm == MadNLP.CHOLESKY + sum(M.info) == 0 ? (M.n, 0, 0) : (0, M.n, 0) + elseif M.opt.lapack_algorithm == MadNLP.EVD + numpos = count(λ -> λ > 0, M.Λ) + numneg = count(λ -> λ < 0, M.Λ) + numzero = M.n - numpos - numneg + (numpos, numzero, numneg) + else + error(M.logger, "Invalid lapack_algorithm") + end +end + +MadNLP.input_type(::Type{LapackOneMKLSolver}) = :dense +MadNLP.default_options(::Type{LapackOneMKLSolver}) = MadNLP.LapackOptions(MadNLP.EVD) +MadNLP.introduce(M::LapackOneMKLSolver) = "OneAPI -- ($(M.opt.lapack_algorithm))" +# MadNLP.introduce(M::LapackOneMKLSolver) = "OneAPI v$(oneAPI.version()) -- ($(M.opt.lapack_algorithm))" + +function setup!(M::LapackOneMKLSolver) + if M.opt.lapack_algorithm == MadNLP.LU + setup_lu!(M) + elseif M.opt.lapack_algorithm == MadNLP.QR + setup_qr!(M) + elseif M.opt.lapack_algorithm == MadNLP.CHOLESKY + setup_cholesky!(M) + elseif M.opt.lapack_algorithm == MadNLP.EVD + setup_evd!(M) + else + error(M.logger, "Invalid lapack_algorithm") + end +end + +function MadNLP.factorize!(M::LapackOneMKLSolver) + MadNLPGPU.gpu_transfer!(M.fact, M.A) + if M.opt.lapack_algorithm == MadNLP.LU + MadNLP.tril_to_full!(M.fact) + factorize_lu!(M) + elseif M.opt.lapack_algorithm == MadNLP.QR + MadNLP.tril_to_full!(M.fact) + factorize_qr!(M) + elseif M.opt.lapack_algorithm == MadNLP.CHOLESKY + factorize_cholesky!(M) + elseif M.opt.lapack_algorithm == MadNLP.EVD + factorize_evd!(M) + else + error(M.logger, "Invalid lapack_algorithm") + end +end + +for T in (:Float32, :Float64) + @eval begin + function MadNLP.solve!(M::LapackOneMKLSolver{$T}, x::oneVector{$T}) + if M.opt.lapack_algorithm == MadNLP.LU + solve_lu!(M, x) + elseif M.opt.lapack_algorithm == MadNLP.QR + solve_qr!(M, x) + elseif M.opt.lapack_algorithm == MadNLP.CHOLESKY + solve_cholesky!(M, x) + elseif M.opt.lapack_algorithm == MadNLP.EVD + solve_evd!(M, x) + else + error(M.logger, "Invalid lapack_algorithm") + end + end + + MadNLP.is_supported(::Type{LapackOneMKLSolver}, ::Type{$T}) = true + end +end + +function MadNLP.solve!(M::LapackOneMKLSolver, x::AbstractVector) + isempty(M.sol) && resize!(M.sol, M.n) + copyto!(M.sol, x) + MadNLP.solve!(M, M.sol) + copyto!(x, M.sol) + return x +end + +for (potrf, potrf_buffer, potrs, potrs_buffer, T) in + ((:onemklDpotrf, :onemklDpotrf_scratchpad_size, :onemklDpotrs, :onemklDpotrs_scratchpad_size, :Float64), + (:onemklSpotrf, :onemklSpotrf_scratchpad_size, :onemklSpotrs, :onemklSpotrs_scratchpad_size, :Float32)) + @eval begin + function setup_cholesky!(M::LapackOneMKLSolver{$T}) + Support.$potrf_buffer(M.device_queue, 'L', M.n, M.n) + Support.$potrs_buffer(M.device_queue, 'L', M.n, one(Int64), M.n, M.n) + return M + end + + function factorize_cholesky!(M::LapackOneMKLSolver{$T}) + Support.$potrf( + M.device_queue, + 'L', + M.n, + M.fact, + M.n, + M.scratchpad, + M.scratchpad_size, + ) + return M + end + + function solve_cholesky!(M::LapackOneMKLSolver{$T}, x::oneVector{$T}) + Support.$potrs( + M.device_queue, + 'L', + n, + one(Int64), + M.fact, + M.n, + x, + M.n, + M.scratchpad, + M.scratchpad_size, + ) + return x + end + end +end + +for (getrf, getrf_buffer, getrs, getrs_buffer, T) in + ((:onemklDgetrf, :onemklDgetrf_scratchpad_size, :onemklDgetrs, :onemklDgetrs_scratchpad_size, :Float64), + (:onemklSgetrf, :onemklSgetrf_scratchpad_size, :onemklSgetrs, :onemklSgetrs_scratchpad_size, :Float32)) + @eval begin + function setup_lu!(M::LapackOneMKLSolver{$T}) + resize!(M.ipiv, M.n) + Support.$getrf_buffer(M.device_queue, M.n, M.n, M.n) + Support.$getrs_buffer(M.device_queue, 'N', M.n, one(Int64), M.n, M.n) + return M + end + + function factorize_lu!(M::LapackOneMKLSolver{$T}) + Support.$getrf( + M.device_queue, + M.n, + M.n, + M.fact, + M.n, + M.ipiv, + M.scratchpad, + M.scratchpad_size, + ) + return M + end + + function solve_lu!(M::LapackOneMKLSolver{$T}, x::oneVector{$T}) + Support.$getrs( + M.device_queue, + 'N', + M.n, + one(Int64), + M.fact, + M.n, + M.ipiv, + x, + M.n, + M.scratchpad, + M.scratchpad_size, + ) + return x + end + end +end + +for (geqrf, geqrf_buffer, ormqr, ormqr_buffer, trsv, T) in + ((:onemklDgeqrf, :onemklDgeqrf_scratchpad_size, :onemklDormqr, :onemklDormqr_scratchpad_size, :onemklDtrsv, :Float64), + (:onemklSgeqrf, :onemklSgeqrf_scratchpad_size, :onemklSormqr, :onemklSormqr_scratchpad_size, :onemklStrsv, :Float32)) + @eval begin + function setup_qr!(M::LapackOneMKLSolver{$T}) + resize!(M.tau, M.n) + Support.$geqrf_buffer(M.device_queue, M.n, M.n, M.n) + Support.$ormqr_buffer(M.device_queue, side, trans, m, n, k, lda, ldc) + return M + end + + function factorize_qr!(M::LapackOneMKLSolver{$T}) + Support.$geqrf( + M.device_queue, + M.n, + M.n, + M.fact, + M.n, + M.tau, + M.scratchpad, + M.scratchpad_size, + ) + return M + end + + function solve_qr!(M::LapackOneMKLSolver{$T}, x::oneVector{$T}) + Support.$ormqr( + device_queue, + side, + trans, + M.n, + M.n, + M.n, + M.fact, + M.n, + M.tau, + c, + ldc, + M.scratchpad, + M.scratchpad_size, + ) + oneBLAS.$trsv( + oneBLAS.handle(), + oneBLAS.oneblas_side_left, + oneBLAS.oneblas_fill_upper, + oneBLAS.oneblas_operation_none, + oneBLAS.oneblas_diagonal_non_unit, + M.n, + one(Int64), + M.alpha, + M.fact, + M.n, + x, + M.n, + ) + return x + end + end +end + +for (syevd, syevd_buffer, gemv, T) in + ((:onemklDsyevd, :onemklDsyevd_scratchpad_size, :onemklDgemv, :Float64), + (:onemklSsyevd, :onemklSsyevd_scratchpad_size, :onemklSgemv, :Float32)) + @eval begin + function setup_evd!(M::LapackOneMKLSolver{$T}) + resize!(M.tau, M.n) + resize!(M.Λ, M.n) + Support.$syevd_buffer(M.device_queue, 'V', 'L', M.n, M.n) + return M + end + + function factorize_evd!(M::LapackOneMKLSolver{$T}) + Support.$syevd( + M.device_queue, + 'V', + 'L', + M.n, + M.fact, + M.n, + M.Λ, + M.scratchpad, + M.scratchpad_size, + ) + return M + end + + function solve_evd!(M::LapackOneMKLSolver{$T}, x::oneVector{$T}) + Support.$gemv( + M.device_queue, + 'T', + M.n, + M.n, + M.alpha, + M.fact, + M.n, + x, + one(Int64), + M.beta, + M.tau, + one(Int64), + ) + M.tau ./= M.Λ + Support.$gemv( + M.device_queue, + 'N', + M.n, + M.n, + M.alpha, + M.fact, + M.n, + M.tau, + one(Int64), + M.beta, + x, + one(Int64), + ) + return x + end + end +end From 4893842b2fbfe2f6f482a6130cfb4f2a3689aedd Mon Sep 17 00:00:00 2001 From: Michel Schanen Date: Mon, 8 Sep 2025 08:46:59 -0500 Subject: [PATCH 02/18] Fix module name --- lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/MadNLPGPUOneAPIExt.jl | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/MadNLPGPUOneAPIExt.jl b/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/MadNLPGPUOneAPIExt.jl index d0f9d5bc9..2060a25b9 100644 --- a/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/MadNLPGPUOneAPIExt.jl +++ b/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/MadNLPGPUOneAPIExt.jl @@ -1,4 +1,4 @@ -module MadNLPGPUAMDGPUExt +module MadNLPGPUOneAPIExt import LinearAlgebra import SparseArrays: SparseMatrixCSC, nonzeros, nnz @@ -10,7 +10,7 @@ import MadNLPGPU import KernelAbstractions: synchronize using oneAPI -using oneAPI.oneMKL, onAPI.Support +using oneAPI.oneMKL, oneAPI.Support function __init__() setglobal!(MadNLPGPU, :LapackOneMKLSolver, LapackOneMKLSolver) From d3fe153aaf6a29f0edc5925de196d31f28422dd6 Mon Sep 17 00:00:00 2001 From: Michel Schanen Date: Mon, 8 Sep 2025 09:15:45 -0500 Subject: [PATCH 03/18] Fix oneSparseMatrix --- .../ext/MadNLPGPUOneAPIExt/oneapi_sparse.jl | 26 +++++++++---------- 1 file changed, 13 insertions(+), 13 deletions(-) diff --git a/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/oneapi_sparse.jl b/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/oneapi_sparse.jl index 6fb494ee9..fa5075fd9 100644 --- a/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/oneapi_sparse.jl +++ b/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/oneapi_sparse.jl @@ -3,11 +3,11 @@ ######################################################## #= - SparseMatrixCOO to oneMKLMatrixCSC + SparseMatrixCOO to oneSparseMatrixCSC =# function MadNLP.transfer!( - dest::oneMKL.oneMKLMatrixCSC, + dest::oneMKL.oneSparseMatrixCSC, src::MadNLP.SparseMatrixCOO, map, ) @@ -102,7 +102,7 @@ function MadNLP.mul_hess_blk!( return end -function MadNLP.get_tril_to_full(csc::oneMKL.oneMKLMatrixCSC{Tv,Ti}) where {Tv,Ti} +function MadNLP.get_tril_to_full(csc::oneMKL.oneSparseMatrixCSC{Tv,Ti}) where {Tv,Ti} cscind = MadNLP.SparseMatrixCSC{Int,Ti}( Symmetric( MadNLP.SparseMatrixCSC{Int,Ti}( @@ -114,7 +114,7 @@ function MadNLP.get_tril_to_full(csc::oneMKL.oneMKLMatrixCSC{Tv,Ti}) where {Tv,T :L, ), ) - return oneMKL.oneMKLMatrixCSC{Tv,Ti}( + return oneMKL.oneSparseMatrixCSC{Tv,Ti}( oneArray(cscind.colptr), oneArray(cscind.rowval), oneVector{Tv}(undef, MadNLP.nnz(cscind)), @@ -185,7 +185,7 @@ function get_diagonal_mapping(colptr, rowval) return rows[inds2], ptrs[inds2] end -function MadNLP._sym_length(Jt::oneMKL.oneMKLMatrixCSC) +function MadNLP._sym_length(Jt::oneMKL.oneSparseMatrixCSC) return mapreduce( (x, y) -> begin z = x - y @@ -205,10 +205,10 @@ function MadNLP._first_and_last_col(sym2::oneVector, ptr2) return (first, last) end -MadNLP.nzval(H::oneMKL.oneMKLMatrixCSC) = H.nzVal +MadNLP.nzval(H::oneMKL.oneSparseMatrixCSC) = H.nzVal function MadNLP._get_sparse_csc(dims, colptr::oneVector, rowval, nzval) - return oneMKL.oneMKLMatrixCSC(colptr, rowval, nzval, dims) + return oneMKL.oneSparseMatrixCSC(colptr, rowval, nzval, dims) end function getij(idx, n) @@ -296,7 +296,7 @@ function MadNLP.coo_to_csc( rowval = map(x -> x[1][1], coord_csc) nzval = similar(rowval, T) - csc = oneMKL.oneMKLMatrixCSC(colptr, rowval, nzval, size(coo)) + csc = oneMKL.oneSparseMatrixCSC(colptr, rowval, nzval, size(coo)) cscmap = similar(coo.I, Int) if length(mapptr) > 1 @@ -318,7 +318,7 @@ end function MadNLP.build_condensed_aug_coord!( kkt::MadNLP.AbstractCondensedKKTSystem{T,VT,MT}, -) where {T,VT,MT<:oneMKL.oneMKLMatrixCSC{T}} +) where {T,VT,MT<:oneMKL.oneSparseMatrixCSC{T}} fill!(kkt.aug_com.nzVal, zero(T)) backend = OneAPIBackend() if length(kkt.hptr) > 0 @@ -359,7 +359,7 @@ end function MadNLP.compress_hessian!( kkt::MadNLP.AbstractSparseKKTSystem{T,VT,MT}, -) where {T,VT,MT<:oneMKL.oneMKLMatrixCSC{T,Int32}} +) where {T,VT,MT<:oneMKL.oneSparseMatrixCSC{T,Int32}} fill!(kkt.hess_com.nzVal, zero(T)) backend = OneAPIBackend() if length(kkt.ext.hess_com_ptrptr) > 1 @@ -377,7 +377,7 @@ end function MadNLP.compress_jacobian!( kkt::MadNLP.SparseCondensedKKTSystem{T,VT,MT}, -) where {T,VT,MT<:oneMKL.oneMKLMatrixCSC{T,Int32}} +) where {T,VT,MT<:oneMKL.oneSparseMatrixCSC{T,Int32}} fill!(kkt.jt_csc.nzVal, zero(T)) backend = OneAPIBackend() if length(kkt.ext.jt_csc_ptrptr) > 1 # otherwise error is thrown @@ -426,7 +426,7 @@ end =# function MadNLP._build_condensed_aug_symbolic_hess( - H::oneMKL.oneMKLMatrixCSC{Tv,Ti}, + H::oneMKL.oneSparseMatrixCSC{Tv,Ti}, sym, sym2, ) where {Tv,Ti} @@ -449,7 +449,7 @@ end =# function MadNLP._build_condensed_aug_symbolic_jt( - Jt::oneMKL.oneMKLMatrixCSC{Tv,Ti}, + Jt::oneMKL.oneSparseMatrixCSC{Tv,Ti}, sym, sym2, ) where {Tv,Ti} From ef72bf8271b7049a78d20c4c40e5754ad76dafd6 Mon Sep 17 00:00:00 2001 From: Michel Schanen Date: Mon, 15 Sep 2025 09:59:42 -0500 Subject: [PATCH 04/18] Add oneapi runner --- .github/workflows/test.yml | 2 +- lib/MadNLPGPU/test/runtests.jl | 2 +- 2 files changed, 2 insertions(+), 2 deletions(-) diff --git a/.github/workflows/test.yml b/.github/workflows/test.yml index a8ca3333e..703fc457d 100644 --- a/.github/workflows/test.yml +++ b/.github/workflows/test.yml @@ -108,7 +108,7 @@ jobs: fail-fast: false matrix: julia-version: ['1'] - gpu: [cuda, amdgpu] + gpu: [cuda, amdgpu, oneapi] runs-on: - self-hosted - ${{ matrix.gpu }} diff --git a/lib/MadNLPGPU/test/runtests.jl b/lib/MadNLPGPU/test/runtests.jl index 8f754f280..773d511b2 100644 --- a/lib/MadNLPGPU/test/runtests.jl +++ b/lib/MadNLPGPU/test/runtests.jl @@ -1,4 +1,4 @@ -using Test, CUDA, AMDGPU, MadNLP, MadNLPGPU, MadNLPTests +using Test, CUDA, AMDGPU, oneAPI, MadNLP, MadNLPGPU, MadNLPTests @testset "MadNLPGPU test" begin include("madnlpgpu_test.jl") From 8787f02e071c8215ff770d19b396f75437c9f309 Mon Sep 17 00:00:00 2001 From: Michel Schanen Date: Thu, 25 Sep 2025 11:18:16 -0500 Subject: [PATCH 05/18] Fix --- Project.toml | 1 + lib/MadNLPGPU/Project.toml | 2 ++ lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/MadNLPGPUOneAPIExt.jl | 1 + lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/oneapi_sparse.jl | 2 +- lib/MadNLPGPU/src/MadNLPGPU.jl | 3 ++- 5 files changed, 7 insertions(+), 2 deletions(-) diff --git a/Project.toml b/Project.toml index ff0ff68a8..1c0d9e3cd 100644 --- a/Project.toml +++ b/Project.toml @@ -23,6 +23,7 @@ MathOptInterface = "b8f27783-ece8-5eb3-8dc8-9495eed66fee" MadNLPMOI = "MathOptInterface" [compat] +GPUArraysCore = "0.2.0" LDLFactorizations = "0.10" MINLPTests = "~0.5" MadNLPTests = "0.5" diff --git a/lib/MadNLPGPU/Project.toml b/lib/MadNLPGPU/Project.toml index 7c937e5c2..7960dfb7c 100644 --- a/lib/MadNLPGPU/Project.toml +++ b/lib/MadNLPGPU/Project.toml @@ -6,6 +6,7 @@ version = "0.7.18" AMD = "14f7f29c-3bd6-536c-9a0b-7339e30b5a3e" CUDA = "052768ef-5323-5732-b1bb-66c8b64840ba" CUDSS = "45b445bb-4962-46a0-9369-b4df9d0f772e" +GPUArraysCore = "46192b85-c4d5-4398-a991-12ede77f4527" KernelAbstractions = "63c18a36-062a-441e-b654-da1e3ab1ce7c" LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e" MadNLP = "2621e9c9-9eb4-46b1-8089-e8c72242dfb6" @@ -25,6 +26,7 @@ AMD = "0.5" AMDGPU = "2" CUDA = "5.4.0" CUDSS = "0.6.4" +GPUArraysCore = "0.2" KernelAbstractions = "0.9" MadNLP = "0.8.12" MadNLPTests = "0.5.3" diff --git a/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/MadNLPGPUOneAPIExt.jl b/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/MadNLPGPUOneAPIExt.jl index 2060a25b9..e40cb9b60 100644 --- a/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/MadNLPGPUOneAPIExt.jl +++ b/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/MadNLPGPUOneAPIExt.jl @@ -8,6 +8,7 @@ import MadNLP import MadNLPGPU import KernelAbstractions: synchronize +import GPUArraysCore: @allowscalar using oneAPI using oneAPI.oneMKL, oneAPI.Support diff --git a/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/oneapi_sparse.jl b/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/oneapi_sparse.jl index fa5075fd9..a914075d3 100644 --- a/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/oneapi_sparse.jl +++ b/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/oneapi_sparse.jl @@ -198,7 +198,7 @@ function MadNLP._sym_length(Jt::oneMKL.oneSparseMatrixCSC) end function MadNLP._first_and_last_col(sym2::oneVector, ptr2) - AMDGPU.@allowscalar begin + @allowscalar begin first = sym2[1][2] last = sym2[ptr2[end]][2] end diff --git a/lib/MadNLPGPU/src/MadNLPGPU.jl b/lib/MadNLPGPU/src/MadNLPGPU.jl index 60667f0ba..6fbcf2540 100644 --- a/lib/MadNLPGPU/src/MadNLPGPU.jl +++ b/lib/MadNLPGPU/src/MadNLPGPU.jl @@ -39,7 +39,8 @@ include("LinearSolvers/cudss.jl") include("cuda.jl") global LapackROCmSolver -export LapackCUDASolver, CUDSSSolver, LapackROCmSolver +global LapackOneMKLSolver +export LapackCUDASolver, CUDSSSolver, LapackROCmSolver, LapackOneMKLSolver # re-export MadNLP, including deprecated names for name in names(MadNLP, all=true) From bbcd9b6a4d4eff6b66cd1971cb3a3360b21cbf44 Mon Sep 17 00:00:00 2001 From: Michel Schanen Date: Thu, 25 Sep 2025 11:21:56 -0500 Subject: [PATCH 06/18] GPUArraysCore fix --- Project.toml | 1 - 1 file changed, 1 deletion(-) diff --git a/Project.toml b/Project.toml index 1c0d9e3cd..ff0ff68a8 100644 --- a/Project.toml +++ b/Project.toml @@ -23,7 +23,6 @@ MathOptInterface = "b8f27783-ece8-5eb3-8dc8-9495eed66fee" MadNLPMOI = "MathOptInterface" [compat] -GPUArraysCore = "0.2.0" LDLFactorizations = "0.10" MINLPTests = "~0.5" MadNLPTests = "0.5" From 7f9b5af5a3e9226445406003bbd91afd1732392b Mon Sep 17 00:00:00 2001 From: Michel Schanen Date: Thu, 25 Sep 2025 13:57:55 -0500 Subject: [PATCH 07/18] WIP --- .../ext/MadNLPGPUOneAPIExt/onemkl.jl | 34 +++++++++++++----- lib/MadNLPGPU/test/madnlpgpu_test.jl | 36 +++++++++++++++++++ lib/MadNLPGPU/test/runtests.jl | 3 ++ 3 files changed, 65 insertions(+), 8 deletions(-) diff --git a/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/onemkl.jl b/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/onemkl.jl index d3555ac21..47e07f19d 100644 --- a/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/onemkl.jl +++ b/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/onemkl.jl @@ -9,6 +9,9 @@ mutable struct LapackOneMKLSolver{T,MT} <: MadNLP.AbstractLinearSolver{T} ipiv::oneVector{Int64} scratchpad::oneVector{T} scratchpad_size::Int64 + device_queue::SYCL.syclQueue_t + alpha::Base.RefValue{T} + beta::Base.RefValue{T} opt::MadNLP.LapackOptions logger::MadNLP.MadNLPLogger @@ -29,7 +32,14 @@ mutable struct LapackOneMKLSolver{T,MT} <: MadNLP.AbstractLinearSolver{T} Λ = oneVector{T}(undef, 0) info = oneVector{Cint}(undef, 1) ipiv = oneVector{Int64}(undef, 0) - solver = new{T,MT}(A, fact, n, sol, tau, Λ, info, ipiv, opt, logger) + scratchpad = oneVector{T}(undef, 0) + scratchpad_size = 0 + # Get the device queue from the oneAPI context + queue = oneAPI.global_queue(oneAPI.context(fact), oneAPI.device(fact)) + device_queue = oneAPI.sycl_queue(queue) + alpha = Ref{T}(1) + beta = Ref{T}(0) + solver = new{T,MT}(A, fact, n, sol, tau, Λ, info, ipiv, scratchpad, scratchpad_size, device_queue, alpha, beta, opt, logger) setup!(solver) return solver end @@ -119,8 +129,10 @@ for (potrf, potrf_buffer, potrs, potrs_buffer, T) in (:onemklSpotrf, :onemklSpotrf_scratchpad_size, :onemklSpotrs, :onemklSpotrs_scratchpad_size, :Float32)) @eval begin function setup_cholesky!(M::LapackOneMKLSolver{$T}) - Support.$potrf_buffer(M.device_queue, 'L', M.n, M.n) - Support.$potrs_buffer(M.device_queue, 'L', M.n, one(Int64), M.n, M.n) + potrf_scratchpad_size = Support.$potrf_buffer(M.device_queue, 'L', M.n, M.n) + potrs_scratchpad_size = Support.$potrs_buffer(M.device_queue, 'L', M.n, one(Int64), M.n, M.n) + M.scratchpad_size = max(potrf_scratchpad_size, potrs_scratchpad_size) + resize!(M.scratchpad, M.scratchpad_size) return M end @@ -161,8 +173,10 @@ for (getrf, getrf_buffer, getrs, getrs_buffer, T) in @eval begin function setup_lu!(M::LapackOneMKLSolver{$T}) resize!(M.ipiv, M.n) - Support.$getrf_buffer(M.device_queue, M.n, M.n, M.n) - Support.$getrs_buffer(M.device_queue, 'N', M.n, one(Int64), M.n, M.n) + getrf_scratchpad_size = Support.$getrf_buffer(M.device_queue, M.n, M.n, M.n) + getrs_scratchpad_size = Support.$getrs_buffer(M.device_queue, 'N', M.n, one(Int64), M.n, M.n) + M.scratchpad_size = max(getrf_scratchpad_size, getrs_scratchpad_size) + resize!(M.scratchpad, M.scratchpad_size) return M end @@ -205,8 +219,11 @@ for (geqrf, geqrf_buffer, ormqr, ormqr_buffer, trsv, T) in @eval begin function setup_qr!(M::LapackOneMKLSolver{$T}) resize!(M.tau, M.n) - Support.$geqrf_buffer(M.device_queue, M.n, M.n, M.n) - Support.$ormqr_buffer(M.device_queue, side, trans, m, n, k, lda, ldc) + geqrf_scratchpad_size = Support.$geqrf_buffer(M.device_queue, M.n, M.n, M.n) + # TODO: Fix ormqr buffer size calculation - undefined variables + # ormqr_scratchpad_size = Support.$ormqr_buffer(M.device_queue, side, trans, m, n, k, lda, ldc) + M.scratchpad_size = geqrf_scratchpad_size + resize!(M.scratchpad, M.scratchpad_size) return M end @@ -266,7 +283,8 @@ for (syevd, syevd_buffer, gemv, T) in function setup_evd!(M::LapackOneMKLSolver{$T}) resize!(M.tau, M.n) resize!(M.Λ, M.n) - Support.$syevd_buffer(M.device_queue, 'V', 'L', M.n, M.n) + M.scratchpad_size = Support.$syevd_buffer(M.device_queue, 'V', 'L', M.n, M.n) + resize!(M.scratchpad, M.scratchpad_size) return M end diff --git a/lib/MadNLPGPU/test/madnlpgpu_test.jl b/lib/MadNLPGPU/test/madnlpgpu_test.jl index b9573a83d..dd0d85925 100644 --- a/lib/MadNLPGPU/test/madnlpgpu_test.jl +++ b/lib/MadNLPGPU/test/madnlpgpu_test.jl @@ -175,6 +175,35 @@ rocm_testset = [ ], ] +oneapi_testset = [ + [ + "LapackOneMKLSolver-LU", + ()->MadNLP.Optimizer( + linear_solver=LapackOneMKLSolver, + lapack_algorithm=MadNLP.LU, + print_level=MadNLP.ERROR, + ), + [], + ], + [ + "LapackOneMKLSolver-QR", + ()->MadNLP.Optimizer( + linear_solver=LapackOneMKLSolver, + lapack_algorithm=MadNLP.QR, + print_level=MadNLP.ERROR, + ), + [], + ], + [ + "LapackOneMKLSolver-CHOLESKY", + ()->MadNLP.Optimizer( + linear_solver=LapackOneMKLSolver, + lapack_algorithm=MadNLP.CHOLESKY, + print_level=MadNLP.ERROR, + ), + ["infeasible", "lootsma", "eigmina", "lp_examodels_issue75"], # KKT system not PD + ], +] @testset "MadNLPGPU test" begin if CUDA.functional() MadNLPTests.test_linear_solver(LapackCUDASolver,Float32) @@ -191,4 +220,11 @@ rocm_testset = [ test_madnlp(name,optimizer_constructor,exclude; Arr=ROCArray) end end + if oneAPI.functional() + MadNLPTests.test_linear_solver(LapackOneMKLSolver,Float32) + MadNLPTests.test_linear_solver(LapackOneMKLSolver,Float64) + for (name,optimizer_constructor,exclude) in oneapi_testset + test_madnlp(name,optimizer_constructor,exclude; Arr=oneArray) + end + end end diff --git a/lib/MadNLPGPU/test/runtests.jl b/lib/MadNLPGPU/test/runtests.jl index 773d511b2..03e72528c 100644 --- a/lib/MadNLPGPU/test/runtests.jl +++ b/lib/MadNLPGPU/test/runtests.jl @@ -12,4 +12,7 @@ using Test, CUDA, AMDGPU, oneAPI, MadNLP, MadNLPGPU, MadNLPTests # Need to add support for CompactLBFGS in SparseCondensedKKTSystem (Issue #563) # include("sparsekkt_rocm.jl") end + # if oneAPI.functional() + # include("densekkt_oneapi.jl") + # end end From 4c53da9b3c1c54a3cbf2599afd0064340326e837 Mon Sep 17 00:00:00 2001 From: Michel Schanen Date: Tue, 30 Sep 2025 12:57:40 -0500 Subject: [PATCH 08/18] WIP --- lib/MadNLPGPU/Project.toml | 6 +- .../ext/MadNLPGPUOneAPIExt/oneapi.jl | 2 +- .../ext/MadNLPGPUOneAPIExt/oneapi_dense.jl | 12 +-- .../ext/MadNLPGPUOneAPIExt/oneapi_sparse.jl | 28 +++--- .../ext/MadNLPGPUOneAPIExt/onemkl.jl | 6 +- lib/MadNLPGPU/test/densekkt_oneapi.jl | 88 +++++++++++++++++++ lib/MadNLPGPU/test/madnlpgpu_test.jl | 40 ++++----- 7 files changed, 135 insertions(+), 47 deletions(-) create mode 100644 lib/MadNLPGPU/test/densekkt_oneapi.jl diff --git a/lib/MadNLPGPU/Project.toml b/lib/MadNLPGPU/Project.toml index 7960dfb7c..883110baf 100644 --- a/lib/MadNLPGPU/Project.toml +++ b/lib/MadNLPGPU/Project.toml @@ -12,10 +12,10 @@ LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e" MadNLP = "2621e9c9-9eb4-46b1-8089-e8c72242dfb6" Metis = "2679e427-3c69-5b7f-982b-ece356f1e94b" SparseArrays = "2f01184e-e22b-5df5-ae63-d93ebab69eaf" +oneAPI = "8f75cd03-7ff8-4ecb-9b8f-daf728133b1b" [weakdeps] AMDGPU = "21141c5a-9bdb-4563-92ae-f87d6854732e" -oneAPI = "8f75cd03-7ff8-4ecb-9b8f-daf728133b1b" [extensions] MadNLPGPUAMDGPUExt = "AMDGPU" @@ -31,12 +31,12 @@ KernelAbstractions = "0.9" MadNLP = "0.8.12" MadNLPTests = "0.5.3" Metis = "1" -oneAPI = "2.2.0" julia = "1.10" +oneAPI = "2.2.0" [extras] MadNLPTests = "b52a2a03-04ab-4a5f-9698-6a2deff93217" Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40" [targets] -test = ["Test", "MadNLPTests", "AMDGPU", "oneAPI"] +test = ["Test", "MadNLPTests", "AMDGPU"] diff --git a/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/oneapi.jl b/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/oneapi.jl index cfeaf2d91..6ad3066e0 100644 --- a/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/oneapi.jl +++ b/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/oneapi.jl @@ -40,7 +40,7 @@ end function MadNLPGPU.gpu_transfer!(y::oneMatrix{T}, x::oneMKL.oneSparseMatrixCSC{T}) where {T} n = size(y, 2) fill!(y, zero(T)) - backend = OneAPIBackend() + backend = oneAPIBackend() MadNLPGPU._csc_to_dense_kernel!(backend)(y, x.colPtr, x.rowVal, x.nzVal, ndrange = n) synchronize(backend) return diff --git a/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/oneapi_dense.jl b/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/oneapi_dense.jl index 61d908e36..6f3bf96bc 100644 --- a/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/oneapi_dense.jl +++ b/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/oneapi_dense.jl @@ -28,7 +28,7 @@ end function MadNLP.diag!(dest::oneVector{T}, src::oneMatrix{T}) where {T} @assert length(dest) == size(src, 1) - backend = OneAPIBackend() + backend = oneAPIBackend() MadNLPGPU._copy_diag_kernel!(backend)(dest, src, ndrange = length(dest)) synchronize(backend) return @@ -39,7 +39,7 @@ end =# function MadNLP.diag_add!(dest::oneMatrix, src1::oneVector, src2::oneVector) - backend = OneAPIBackend() + backend = oneAPIBackend() MadNLPGPU._add_diagonal_kernel!(backend)(dest, src1, src2, ndrange = size(dest, 1)) synchronize(backend) return @@ -51,7 +51,7 @@ end function MadNLP._set_diag!(A::oneMatrix, inds, a) if !isempty(inds) - backend = OneAPIBackend() + backend = oneAPIBackend() MadNLPGPU._set_diag_kernel!(backend)(A, inds, a; ndrange = length(inds)) synchronize(backend) end @@ -76,7 +76,7 @@ function MadNLP._build_dense_kkt_system!( ) ind_ineq_gpu = oneVector(ind_ineq) ndrange = (n + m + ns, n) - backend = OneAPIBackend() + backend = oneAPIBackend() MadNLPGPU._build_dense_kkt_system_kernel!(backend)( dest, hess, @@ -109,7 +109,7 @@ function MadNLP._build_ineq_jac!( (m_ineq == 0) && return # nothing to do if no ineq. constraints ind_ineq_gpu = oneVector(ind_ineq) ndrange = (m_ineq, n) - backend = OneAPIBackend() + backend = oneAPIBackend() MadNLPGPU._build_jacobian_condensed_kernel!(backend)( dest, jac, @@ -138,7 +138,7 @@ function MadNLP._build_condensed_kkt_system!( ) ind_eq_gpu = oneVector(ind_eq) ndrange = (n + m_eq, n) - backend = OneAPIBackend() + backend = oneAPIBackend() MadNLPGPU._build_condensed_kkt_system_kernel!(backend)( dest, hess, diff --git a/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/oneapi_sparse.jl b/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/oneapi_sparse.jl index a914075d3..dcc04db25 100644 --- a/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/oneapi_sparse.jl +++ b/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/oneapi_sparse.jl @@ -42,7 +42,7 @@ function MadNLP.mul!( MadNLP.mul!(wx, kkt.hess_com', xx, alpha, one(T)) MadNLP.mul!(wx, kkt.jt_csc, xz, alpha, beta) if !isempty(kkt.ext.diag_map_to) - backend = OneAPIBackend() + backend = oneAPIBackend() MadNLPGPU._diag_operation_kernel!(backend)( wx, kkt.hess_com.nzVal, @@ -84,7 +84,7 @@ function MadNLP.mul_hess_blk!( MadNLP.mul!(wxx, kkt.hess_com, tx, one(T), zero(T)) MadNLP.mul!(wxx, kkt.hess_com', tx, one(T), one(T)) if !isempty(kkt.ext.diag_map_to) - backend = OneAPIBackend() + backend = oneAPIBackend() MadNLPGPU._diag_operation_kernel!(backend)( wxx, kkt.hess_com.nzVal, @@ -223,7 +223,7 @@ end function MadNLP._set_colptr!(colptr::oneVector, ptr2, sym2, guide) if length(ptr2) > 1 # otherwise error is thrown - backend = OneAPIBackend() + backend = oneAPIBackend() MadNLPGPU._set_colptr_kernel!(backend)( colptr, sym2, @@ -243,7 +243,7 @@ end function MadNLP.tril_to_full!(dense::oneMatrix{T}) where {T} n = size(dense, 1) - backend = OneAPIBackend() + backend = oneAPIBackend() MadNLPGPU._tril_to_full_kernel!(backend)(dense; ndrange = div(n^2 + n, 2)) synchronize(backend) return @@ -255,7 +255,7 @@ end function MadNLP.force_lower_triangular!(I::oneVector{T}, J) where {T} if !isempty(I) - backend = OneAPIBackend() + backend = oneAPIBackend() MadNLPGPU._force_lower_triangular_kernel!(backend)(I, J; ndrange = length(I)) synchronize(backend) end @@ -281,7 +281,7 @@ function MadNLP.coo_to_csc( coord_csc = coord[@view(mapptr[1:end-1])] - backend = OneAPIBackend() + backend = oneAPIBackend() if length(coord_csc) > 0 MadNLPGPU._set_coo_to_colptr_kernel!(backend)( colptr, @@ -295,7 +295,7 @@ function MadNLP.coo_to_csc( rowval = map(x -> x[1][1], coord_csc) nzval = similar(rowval, T) - + @show size(coo) csc = oneMKL.oneSparseMatrixCSC(colptr, rowval, nzval, size(coo)) cscmap = similar(coo.I, Int) @@ -320,7 +320,7 @@ function MadNLP.build_condensed_aug_coord!( kkt::MadNLP.AbstractCondensedKKTSystem{T,VT,MT}, ) where {T,VT,MT<:oneMKL.oneSparseMatrixCSC{T}} fill!(kkt.aug_com.nzVal, zero(T)) - backend = OneAPIBackend() + backend = oneAPIBackend() if length(kkt.hptr) > 0 MadNLPGPU._transfer_hessian_kernel!(backend)( kkt.aug_com.nzVal, @@ -361,7 +361,7 @@ function MadNLP.compress_hessian!( kkt::MadNLP.AbstractSparseKKTSystem{T,VT,MT}, ) where {T,VT,MT<:oneMKL.oneSparseMatrixCSC{T,Int32}} fill!(kkt.hess_com.nzVal, zero(T)) - backend = OneAPIBackend() + backend = oneAPIBackend() if length(kkt.ext.hess_com_ptrptr) > 1 MadNLPGPU._transfer_to_csc_kernel!(backend)( kkt.hess_com.nzVal, @@ -379,7 +379,7 @@ function MadNLP.compress_jacobian!( kkt::MadNLP.SparseCondensedKKTSystem{T,VT,MT}, ) where {T,VT,MT<:oneMKL.oneSparseMatrixCSC{T,Int32}} fill!(kkt.jt_csc.nzVal, zero(T)) - backend = OneAPIBackend() + backend = oneAPIBackend() if length(kkt.ext.jt_csc_ptrptr) > 1 # otherwise error is thrown MadNLPGPU._transfer_to_csc_kernel!(backend)( kkt.jt_csc.nzVal, @@ -407,7 +407,7 @@ function MadNLP._set_con_scale_sparse!( !isempty(inds) && sort!(inds) ptr = MadNLP.getptr(inds; by = ((x1, x2), (y1, y2)) -> x1 != y1) if length(ptr) > 1 - backend = OneAPIBackend() + backend = oneAPIBackend() MadNLPGPU._set_con_scale_sparse_kernel!(backend)( con_scale, ptr, @@ -431,7 +431,7 @@ function MadNLP._build_condensed_aug_symbolic_hess( sym2, ) where {Tv,Ti} if size(H, 2) > 0 - backend = OneAPIBackend() + backend = oneAPIBackend() MadNLPGPU._build_condensed_aug_symbolic_hess_kernel!(backend)( sym, sym2, @@ -460,7 +460,7 @@ function MadNLP._build_condensed_aug_symbolic_jt( @view(Jt.colPtr[2:end]) ) offsets = cumsum(_offsets) - backend = OneAPIBackend() + backend = oneAPIBackend() MadNLPGPU._build_condensed_aug_symbolic_jt_kernel!(backend)( sym, sym2, @@ -479,7 +479,7 @@ end =# function MadNLP._build_scale_augmented_system_coo!(dest, src, scaling::oneArray, n, m) - backend = OneAPIBackend() + backend = oneAPIBackend() MadNLPGPU._scale_augmented_system_coo_kernel!(backend)( dest.V, src.I, diff --git a/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/onemkl.jl b/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/onemkl.jl index 47e07f19d..a8f4fc514 100644 --- a/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/onemkl.jl +++ b/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/onemkl.jl @@ -243,9 +243,9 @@ for (geqrf, geqrf_buffer, ormqr, ormqr_buffer, trsv, T) in function solve_qr!(M::LapackOneMKLSolver{$T}, x::oneVector{$T}) Support.$ormqr( - device_queue, - side, - trans, + M.device_queue, + M.side, + M.trans, M.n, M.n, M.n, diff --git a/lib/MadNLPGPU/test/densekkt_oneapi.jl b/lib/MadNLPGPU/test/densekkt_oneapi.jl new file mode 100644 index 000000000..896611201 --- /dev/null +++ b/lib/MadNLPGPU/test/densekkt_oneapi.jl @@ -0,0 +1,88 @@ +using oneAPI +using MadNLPTests + +function _compare_oneapi_with_cpu(KKTSystem, n, m, ind_fixed) + for (T,tol,atol) in [ + (Float32,1e-4,1e1), + (Float64,1e-8,1e-6) + ] + madnlp_options = Dict{Symbol, Any}( + :callback=>MadNLP.DenseCallback, + :kkt_system=>KKTSystem, + :linear_solver=>LapackOneMKLSolver, + :lapack_algorithm=>MadNLP.QR, + :print_level=>MadNLP.ERROR, + :tol=>tol + ) + + # Host evaluator + nlph = MadNLPTests.DenseDummyQP(zeros(T,n); m=m, fixed_variables=ind_fixed) + # Device evaluator + nlpd = MadNLPTests.DenseDummyQP(oneAPI.zeros(T,n); m=m, fixed_variables=oneArray(ind_fixed)) + + # Solve on CPU + h_solver = MadNLPSolver(nlph; madnlp_options...) + results_cpu = MadNLP.solve!(h_solver) + + # Solve on GPU + d_solver = MadNLPSolver(nlpd; madnlp_options...) + results_gpu = MadNLP.solve!(d_solver) + + @test isa(d_solver.kkt, KKTSystem{T}) + # Check that both results match exactly + if T == Float64 + @test h_solver.cnt.k == d_solver.cnt.k + @test results_cpu.objective ≈ results_gpu.objective + @test results_cpu.solution ≈ Array(results_gpu.solution) atol=atol + @test results_cpu.multipliers ≈ Array(results_gpu.multipliers) atol=atol + end + end +end + +@testset "MadNLPGPU -- LapackOneMKLSolver -- ($(kkt_system))" for kkt_system in [ + MadNLP.DenseKKTSystem, + MadNLP.DenseCondensedKKTSystem, + ] + @testset "Size: ($n, $m)" for (n, m) in [(10, 0), (10, 5), (50, 10)] + _compare_oneapi_with_cpu(kkt_system, n, m, Int[]) + end + @testset "Fixed variables" for (n,m) in [(10, 0), (10, 5), (50, 10)] + _compare_oneapi_with_cpu(kkt_system, n, m, Int[1, 2]) + end +end + +@testset "MadNLP -- LapackOneMKLSolver: $QN + $KKT" for QN in [ + MadNLP.BFGS, + MadNLP.DampedBFGS, +], KKT in [ + MadNLP.DenseKKTSystem, + MadNLP.DenseCondensedKKTSystem, +] + @testset "Size: ($n, $m)" for (n, m) in [(10, 0), (10, 5), (50, 10)] + nlp = MadNLPTests.DenseDummyQP(zeros(Float64, n); m=m) + solver_exact = MadNLPSolver( + nlp; + callback=MadNLP.DenseCallback, + print_level=MadNLP.ERROR, + kkt_system=KKT, + linear_solver=LapackOneMKLSolver, + ) + results_ref = MadNLP.solve!(solver_exact) + + nlp = MadNLPTests.DenseDummyQP(oneAPI.zeros(Float64, n); m=m) + solver_qn = MadNLPSolver( + nlp; + callback=MadNLP.DenseCallback, + print_level=MadNLP.ERROR, + kkt_system=KKT, + hessian_approximation=QN, + linear_solver=LapackOneMKLSolver, + ) + results_qn = MadNLP.solve!(solver_qn) + + @test results_qn.status == MadNLP.SOLVE_SUCCEEDED + @test results_qn.objective ≈ results_ref.objective atol=1e-6 + @test Array(results_qn.solution) ≈ Array(results_ref.solution) atol=1e-6 + @test solver_qn.cnt.lag_hess_cnt == 0 + end +end diff --git a/lib/MadNLPGPU/test/madnlpgpu_test.jl b/lib/MadNLPGPU/test/madnlpgpu_test.jl index dd0d85925..d0772e2ef 100644 --- a/lib/MadNLPGPU/test/madnlpgpu_test.jl +++ b/lib/MadNLPGPU/test/madnlpgpu_test.jl @@ -185,24 +185,24 @@ oneapi_testset = [ ), [], ], - [ - "LapackOneMKLSolver-QR", - ()->MadNLP.Optimizer( - linear_solver=LapackOneMKLSolver, - lapack_algorithm=MadNLP.QR, - print_level=MadNLP.ERROR, - ), - [], - ], - [ - "LapackOneMKLSolver-CHOLESKY", - ()->MadNLP.Optimizer( - linear_solver=LapackOneMKLSolver, - lapack_algorithm=MadNLP.CHOLESKY, - print_level=MadNLP.ERROR, - ), - ["infeasible", "lootsma", "eigmina", "lp_examodels_issue75"], # KKT system not PD - ], + # [ + # "LapackOneMKLSolver-QR", + # ()->MadNLP.Optimizer( + # linear_solver=LapackOneMKLSolver, + # lapack_algorithm=MadNLP.QR, + # print_level=MadNLP.ERROR, + # ), + # [], + # ], + # [ + # "LapackOneMKLSolver-CHOLESKY", + # ()->MadNLP.Optimizer( + # linear_solver=LapackOneMKLSolver, + # lapack_algorithm=MadNLP.CHOLESKY, + # print_level=MadNLP.ERROR, + # ), + # ["infeasible", "lootsma", "eigmina", "lp_examodels_issue75"], # KKT system not PD + # ], ] @testset "MadNLPGPU test" begin if CUDA.functional() @@ -214,8 +214,8 @@ oneapi_testset = [ end end if AMDGPU.functional() - MadNLPTests.test_linear_solver(LapackROCmSolver,Float32) - MadNLPTests.test_linear_solver(LapackROCmSolver,Float64) + MadNLPTests.test_linear_solver(LapackROCSolver,Float32) + MadNLPTests.test_linear_solver(LapackROCSolver,Float64) for (name,optimizer_constructor,exclude) in rocm_testset test_madnlp(name,optimizer_constructor,exclude; Arr=ROCArray) end From abfbc928c810d43bcf22ccf8619687e9d3756d53 Mon Sep 17 00:00:00 2001 From: Github action runner Date: Wed, 4 Feb 2026 10:20:45 -0800 Subject: [PATCH 09/18] Rebase fixes --- .../ext/MadNLPGPUOneAPIExt/oneapi.jl | 30 +++++++++++++++++++ .../ext/MadNLPGPUOneAPIExt/oneapi_dense.jl | 8 +---- lib/MadNLPGPU/test/runtests.jl | 6 ++-- 3 files changed, 34 insertions(+), 10 deletions(-) diff --git a/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/oneapi.jl b/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/oneapi.jl index 6ad3066e0..5669c2867 100644 --- a/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/oneapi.jl +++ b/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/oneapi.jl @@ -45,3 +45,33 @@ function MadNLPGPU.gpu_transfer!(y::oneMatrix{T}, x::oneMKL.oneSparseMatrixCSC{T synchronize(backend) return end + +#= + MadNLP._syr! +=# + +MadNLP._syr!(uplo::Char, alpha::T, x::oneVector{T}, A::oneMatrix{T}) where T = oneMKL.syr!(uplo, alpha, x, A) + +#= + MadNLP._symv! +=# + +MadNLP._symv!(uplo::Char, alpha::T, A::oneMatrix{T}, x::oneVector{T}, beta::T, y::oneVector{T}) where T = oneMKL.symv!(uplo, alpha, A, x, beta, y) + +#= + MadNLP._syrk! +=# + +MadNLP._syrk!(uplo::Char, trans::Char, alpha::T, A::oneMatrix{T}, beta::T, C::oneMatrix{T}) where T = oneMKL.syrk!(uplo, trans, alpha, A, beta, C) + +#= + MadNLP._trsm! +=# + +MadNLP._trsm!(side::Char, uplo::Char, transa::Char, diag::Char, alpha::T, A::oneMatrix{T}, B::oneMatrix{T}) where T = oneMKL.trsm!(side, uplo, transa, diag, alpha, A, B) + +#= + MadNLP._dgmm! +=# + +MadNLP._dgmm!(side::Char, A::oneMatrix{T}, x::oneVector{T}, B::oneMatrix{T}) where T = oneMKL.dgmm!(side, A, x, B) diff --git a/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/oneapi_dense.jl b/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/oneapi_dense.jl index 6f3bf96bc..2f2578e3a 100644 --- a/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/oneapi_dense.jl +++ b/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/oneapi_dense.jl @@ -2,17 +2,11 @@ ##### oneAPI wrappers for DenseKKTSystem / DenseCondensedKKTSystem ##### ######################################################################## -#= - MadNLP.symul! -=# - -MadNLP.symul!(y, A, x::oneVector{T}, α = one(T), β = zero(T)) where T = oneMKL.symv!('L', T(α), A, x, T(β), y) - #= MadNLP._ger! =# -MadNLP._ger!(alpha::Number, x::oneVector{T}, y::oneVector{T}, A::oneMatrix{T}) where T = oneMKL.ger!(alpha, x, y, A) +MadNLP._ger!(alpha::T, x::oneVector{T}, y::oneVector{T}, A::oneMatrix{T}) where T = oneMKL.ger!(alpha, x, y, A) #= MadNLP._madnlp_unsafe_wrap diff --git a/lib/MadNLPGPU/test/runtests.jl b/lib/MadNLPGPU/test/runtests.jl index 03e72528c..0a6cbe337 100644 --- a/lib/MadNLPGPU/test/runtests.jl +++ b/lib/MadNLPGPU/test/runtests.jl @@ -12,7 +12,7 @@ using Test, CUDA, AMDGPU, oneAPI, MadNLP, MadNLPGPU, MadNLPTests # Need to add support for CompactLBFGS in SparseCondensedKKTSystem (Issue #563) # include("sparsekkt_rocm.jl") end - # if oneAPI.functional() - # include("densekkt_oneapi.jl") - # end + if oneAPI.functional() + include("densekkt_oneapi.jl") + end end From 190f57692f5458672418d98875a1b6343312c1d8 Mon Sep 17 00:00:00 2001 From: Github action runner Date: Wed, 4 Feb 2026 11:23:50 -0800 Subject: [PATCH 10/18] oneAPI fixes --- .../ext/MadNLPGPUOneAPIExt/oneapi.jl | 22 ++++++++++ .../ext/MadNLPGPUOneAPIExt/oneapi_sparse.jl | 1 - .../ext/MadNLPGPUOneAPIExt/onemkl.jl | 42 +++++++------------ 3 files changed, 37 insertions(+), 28 deletions(-) diff --git a/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/oneapi.jl b/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/oneapi.jl index 5669c2867..a1012b054 100644 --- a/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/oneapi.jl +++ b/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/oneapi.jl @@ -75,3 +75,25 @@ MadNLP._trsm!(side::Char, uplo::Char, transa::Char, diag::Char, alpha::T, A::one =# MadNLP._dgmm!(side::Char, A::oneMatrix{T}, x::oneVector{T}, B::oneMatrix{T}) where T = oneMKL.dgmm!(side, A, x, B) + +#= + MadNLP.get_sd / MadNLP.get_sc + Workaround for norm() on views of oneArray triggering scalar indexing in Julia 1.12+ + See https://github.com/JuliaGPU/CUDA.jl/issues/2811 for similar issue in CUDA.jl +=# + +if VERSION > v"1.11" + function MadNLP.get_sd(l::oneVector{T}, zl_r, zu_r, s_max) where T + return max( + s_max, + (my1norm(l)+my1norm(zl_r)+my1norm(zu_r)) / max(1, (length(l)+length(zl_r)+length(zu_r))), + ) / s_max + end + function MadNLP.get_sc(zl_r::SubArray{T,1,VT}, zu_r, s_max) where {T, VT <: oneVector{T}} + return max( + s_max, + (my1norm(zl_r)+my1norm(zu_r)) / max(1,length(zl_r)+length(zu_r)), + ) / s_max + end + my1norm(x) = mapreduce(abs, +, x) +end diff --git a/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/oneapi_sparse.jl b/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/oneapi_sparse.jl index dcc04db25..1069ad8f2 100644 --- a/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/oneapi_sparse.jl +++ b/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/oneapi_sparse.jl @@ -295,7 +295,6 @@ function MadNLP.coo_to_csc( rowval = map(x -> x[1][1], coord_csc) nzval = similar(rowval, T) - @show size(coo) csc = oneMKL.oneSparseMatrixCSC(colptr, rowval, nzval, size(coo)) cscmap = similar(coo.I, Int) diff --git a/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/onemkl.jl b/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/onemkl.jl index a8f4fc514..ec1bb29a7 100644 --- a/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/onemkl.jl +++ b/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/onemkl.jl @@ -220,9 +220,8 @@ for (geqrf, geqrf_buffer, ormqr, ormqr_buffer, trsv, T) in function setup_qr!(M::LapackOneMKLSolver{$T}) resize!(M.tau, M.n) geqrf_scratchpad_size = Support.$geqrf_buffer(M.device_queue, M.n, M.n, M.n) - # TODO: Fix ormqr buffer size calculation - undefined variables - # ormqr_scratchpad_size = Support.$ormqr_buffer(M.device_queue, side, trans, m, n, k, lda, ldc) - M.scratchpad_size = geqrf_scratchpad_size + ormqr_scratchpad_size = Support.$ormqr_buffer(M.device_queue, 'L', 'T', M.n, one(Int64), M.n, M.n, M.n) + M.scratchpad_size = max(geqrf_scratchpad_size, ormqr_scratchpad_size) resize!(M.scratchpad, M.scratchpad_size) return M end @@ -242,35 +241,24 @@ for (geqrf, geqrf_buffer, ormqr, ormqr_buffer, trsv, T) in end function solve_qr!(M::LapackOneMKLSolver{$T}, x::oneVector{$T}) + # Apply Q^T to x: x = Q^T * x Support.$ormqr( M.device_queue, - M.side, - M.trans, - M.n, - M.n, - M.n, - M.fact, - M.n, - M.tau, - c, - ldc, + 'L', # side (left multiplication) + 'T', # trans (transpose) + M.n, # m + one(Int64), # n (single RHS) + M.n, # k + M.fact, # A + M.n, # lda + M.tau, # tau + x, # c (the RHS vector) + M.n, # ldc M.scratchpad, M.scratchpad_size, ) - oneBLAS.$trsv( - oneBLAS.handle(), - oneBLAS.oneblas_side_left, - oneBLAS.oneblas_fill_upper, - oneBLAS.oneblas_operation_none, - oneBLAS.oneblas_diagonal_non_unit, - M.n, - one(Int64), - M.alpha, - M.fact, - M.n, - x, - M.n, - ) + # Solve R*x = Q^T*b using triangular solve + oneMKL.trsv!('U', 'N', 'N', M.fact, x) # upper, no-trans, non-unit diagonal return x end end From 8dc137470bcf8c4adac8b4430bead6be8432ea5f Mon Sep 17 00:00:00 2001 From: Github action runner Date: Wed, 4 Feb 2026 11:32:50 -0800 Subject: [PATCH 11/18] Format --- .../ext/MadNLPGPUOneAPIExt/oneapi.jl | 38 +++--- .../ext/MadNLPGPUOneAPIExt/oneapi_dense.jl | 60 ++++----- .../ext/MadNLPGPUOneAPIExt/oneapi_sparse.jl | 118 +++++++++--------- .../ext/MadNLPGPUOneAPIExt/onemkl.jl | 50 ++++---- lib/MadNLPGPU/test/densekkt_oneapi.jl | 69 +++++----- lib/MadNLPGPU/test/madnlpgpu_test.jl | 20 +-- 6 files changed, 182 insertions(+), 173 deletions(-) diff --git a/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/oneapi.jl b/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/oneapi.jl index a1012b054..b09a84a47 100644 --- a/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/oneapi.jl +++ b/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/oneapi.jl @@ -3,14 +3,14 @@ =# function MadNLP.MadNLPOptions{T}( - nlp::MadNLP.AbstractNLPModel{T,VT}; - dense_callback = MadNLP.is_dense_callback(nlp), - callback = dense_callback ? MadNLP.DenseCallback : MadNLP.SparseCallback, - kkt_system = dense_callback ? MadNLP.DenseCondensedKKTSystem : MadNLP.SparseCondensedKKTSystem, - linear_solver = MadNLPGPU.LapackOneMKLSolver, - tol = MadNLP.get_tolerance(T,kkt_system), - bound_relax_factor = tol, -) where {T, VT <: oneVector{T}} + nlp::MadNLP.AbstractNLPModel{T, VT}; + dense_callback = MadNLP.is_dense_callback(nlp), + callback = dense_callback ? MadNLP.DenseCallback : MadNLP.SparseCallback, + kkt_system = dense_callback ? MadNLP.DenseCondensedKKTSystem : MadNLP.SparseCondensedKKTSystem, + linear_solver = MadNLPGPU.LapackOneMKLSolver, + tol = MadNLP.get_tolerance(T, kkt_system), + bound_relax_factor = tol, + ) where {T, VT <: oneVector{T}} return MadNLP.MadNLPOptions{T}( tol = tol, callback = callback, @@ -24,8 +24,8 @@ end SparseMatrixCSC to oneSparseMatrixCSC =# -function oneMKL.oneSparseMatrixCSC{Tv,Ti}(A::SparseMatrixCSC{Tv,Ti}) where {Tv,Ti} - return oneMKL.oneSparseMatrixCSC{Tv,Ti}( +function oneMKL.oneSparseMatrixCSC{Tv, Ti}(A::SparseMatrixCSC{Tv, Ti}) where {Tv, Ti} + return oneMKL.oneSparseMatrixCSC{Tv, Ti}( oneVector(A.colptr), oneVector(A.rowval), oneVector(A.nzval), @@ -50,31 +50,31 @@ end MadNLP._syr! =# -MadNLP._syr!(uplo::Char, alpha::T, x::oneVector{T}, A::oneMatrix{T}) where T = oneMKL.syr!(uplo, alpha, x, A) +MadNLP._syr!(uplo::Char, alpha::T, x::oneVector{T}, A::oneMatrix{T}) where {T} = oneMKL.syr!(uplo, alpha, x, A) #= MadNLP._symv! =# -MadNLP._symv!(uplo::Char, alpha::T, A::oneMatrix{T}, x::oneVector{T}, beta::T, y::oneVector{T}) where T = oneMKL.symv!(uplo, alpha, A, x, beta, y) +MadNLP._symv!(uplo::Char, alpha::T, A::oneMatrix{T}, x::oneVector{T}, beta::T, y::oneVector{T}) where {T} = oneMKL.symv!(uplo, alpha, A, x, beta, y) #= MadNLP._syrk! =# -MadNLP._syrk!(uplo::Char, trans::Char, alpha::T, A::oneMatrix{T}, beta::T, C::oneMatrix{T}) where T = oneMKL.syrk!(uplo, trans, alpha, A, beta, C) +MadNLP._syrk!(uplo::Char, trans::Char, alpha::T, A::oneMatrix{T}, beta::T, C::oneMatrix{T}) where {T} = oneMKL.syrk!(uplo, trans, alpha, A, beta, C) #= MadNLP._trsm! =# -MadNLP._trsm!(side::Char, uplo::Char, transa::Char, diag::Char, alpha::T, A::oneMatrix{T}, B::oneMatrix{T}) where T = oneMKL.trsm!(side, uplo, transa, diag, alpha, A, B) +MadNLP._trsm!(side::Char, uplo::Char, transa::Char, diag::Char, alpha::T, A::oneMatrix{T}, B::oneMatrix{T}) where {T} = oneMKL.trsm!(side, uplo, transa, diag, alpha, A, B) #= MadNLP._dgmm! =# -MadNLP._dgmm!(side::Char, A::oneMatrix{T}, x::oneVector{T}, B::oneMatrix{T}) where T = oneMKL.dgmm!(side, A, x, B) +MadNLP._dgmm!(side::Char, A::oneMatrix{T}, x::oneVector{T}, B::oneMatrix{T}) where {T} = oneMKL.dgmm!(side, A, x, B) #= MadNLP.get_sd / MadNLP.get_sc @@ -83,16 +83,16 @@ MadNLP._dgmm!(side::Char, A::oneMatrix{T}, x::oneVector{T}, B::oneMatrix{T}) whe =# if VERSION > v"1.11" - function MadNLP.get_sd(l::oneVector{T}, zl_r, zu_r, s_max) where T + function MadNLP.get_sd(l::oneVector{T}, zl_r, zu_r, s_max) where {T} return max( s_max, - (my1norm(l)+my1norm(zl_r)+my1norm(zu_r)) / max(1, (length(l)+length(zl_r)+length(zu_r))), + (my1norm(l) + my1norm(zl_r) + my1norm(zu_r)) / max(1, (length(l) + length(zl_r) + length(zu_r))), ) / s_max end - function MadNLP.get_sc(zl_r::SubArray{T,1,VT}, zu_r, s_max) where {T, VT <: oneVector{T}} + function MadNLP.get_sc(zl_r::SubArray{T, 1, VT}, zu_r, s_max) where {T, VT <: oneVector{T}} return max( s_max, - (my1norm(zl_r)+my1norm(zu_r)) / max(1,length(zl_r)+length(zu_r)), + (my1norm(zl_r) + my1norm(zu_r)) / max(1, length(zl_r) + length(zu_r)), ) / s_max end my1norm(x) = mapreduce(abs, +, x) diff --git a/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/oneapi_dense.jl b/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/oneapi_dense.jl index 2f2578e3a..555e730ce 100644 --- a/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/oneapi_dense.jl +++ b/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/oneapi_dense.jl @@ -6,14 +6,14 @@ MadNLP._ger! =# -MadNLP._ger!(alpha::T, x::oneVector{T}, y::oneVector{T}, A::oneMatrix{T}) where T = oneMKL.ger!(alpha, x, y, A) +MadNLP._ger!(alpha::T, x::oneVector{T}, y::oneVector{T}, A::oneMatrix{T}) where {T} = oneMKL.ger!(alpha, x, y, A) #= MadNLP._madnlp_unsafe_wrap =# -function MadNLP._madnlp_unsafe_wrap(vec::VT, n, shift=1) where {T, VT <: oneVector{T}} - return view(vec,shift:shift+n-1) +function MadNLP._madnlp_unsafe_wrap(vec::VT, n, shift = 1) where {T, VT <: oneVector{T}} + return view(vec, shift:(shift + n - 1)) end #= @@ -57,17 +57,17 @@ end =# function MadNLP._build_dense_kkt_system!( - dest::oneMatrix, - hess::oneMatrix, - jac::oneMatrix, - pr_diag::oneVector, - du_diag::oneVector, - diag_hess::oneVector, - ind_ineq::AbstractVector, - n, - m, - ns, -) + dest::oneMatrix, + hess::oneMatrix, + jac::oneMatrix, + pr_diag::oneVector, + du_diag::oneVector, + diag_hess::oneVector, + ind_ineq::AbstractVector, + n, + m, + ns, + ) ind_ineq_gpu = oneVector(ind_ineq) ndrange = (n + m + ns, n) backend = oneAPIBackend() @@ -93,13 +93,13 @@ end =# function MadNLP._build_ineq_jac!( - dest::oneMatrix, - jac::oneMatrix, - diag_buffer::oneVector, - ind_ineq::AbstractVector, - n, - m_ineq, -) + dest::oneMatrix, + jac::oneMatrix, + diag_buffer::oneVector, + ind_ineq::AbstractVector, + n, + m_ineq, + ) (m_ineq == 0) && return # nothing to do if no ineq. constraints ind_ineq_gpu = oneVector(ind_ineq) ndrange = (m_ineq, n) @@ -121,15 +121,15 @@ end =# function MadNLP._build_condensed_kkt_system!( - dest::oneMatrix, - hess::oneMatrix, - jac::oneMatrix, - pr_diag::oneVector, - du_diag::oneVector, - ind_eq::AbstractVector, - n, - m_eq, -) + dest::oneMatrix, + hess::oneMatrix, + jac::oneMatrix, + pr_diag::oneVector, + du_diag::oneVector, + ind_eq::AbstractVector, + n, + m_eq, + ) ind_eq_gpu = oneVector(ind_eq) ndrange = (n + m_eq, n) backend = oneAPIBackend() diff --git a/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/oneapi_sparse.jl b/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/oneapi_sparse.jl index 1069ad8f2..aaa1dffea 100644 --- a/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/oneapi_sparse.jl +++ b/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/oneapi_sparse.jl @@ -7,10 +7,10 @@ =# function MadNLP.transfer!( - dest::oneMKL.oneSparseMatrixCSC, - src::MadNLP.SparseMatrixCOO, - map, -) + dest::oneMKL.oneSparseMatrixCSC, + src::MadNLP.SparseMatrixCOO, + map, + ) return copyto!(view(dest.nzVal, map), src.V) end @@ -19,24 +19,24 @@ end =# function MadNLP.mul!( - w::MadNLP.AbstractKKTVector{T,VT}, - kkt::MadNLP.SparseCondensedKKTSystem, - x::MadNLP.AbstractKKTVector, - alpha = one(T), - beta = zero(T), -) where {T,VT<:oneVector{T}} + w::MadNLP.AbstractKKTVector{T, VT}, + kkt::MadNLP.SparseCondensedKKTSystem, + x::MadNLP.AbstractKKTVector, + alpha = one(T), + beta = zero(T), + ) where {T, VT <: oneVector{T}} n = size(kkt.hess_com, 1) m = size(kkt.jt_csc, 2) # Decompose results xx = view(MadNLP.full(x), 1:n) - xs = view(MadNLP.full(x), n+1:n+m) - xz = view(MadNLP.full(x), n+m+1:n+2*m) + xs = view(MadNLP.full(x), (n + 1):(n + m)) + xz = view(MadNLP.full(x), (n + m + 1):(n + 2 * m)) # Decompose buffers wx = view(MadNLP.full(w), 1:n) - ws = view(MadNLP.full(w), n+1:n+m) - wz = view(MadNLP.full(w), n+m+1:n+2*m) + ws = view(MadNLP.full(w), (n + 1):(n + m)) + wz = view(MadNLP.full(w), (n + m + 1):(n + 2 * m)) MadNLP.mul!(wx, kkt.hess_com, xx, alpha, beta) MadNLP.mul!(wx, kkt.hess_com', xx, alpha, one(T)) @@ -73,10 +73,10 @@ function MadNLP.mul!( end function MadNLP.mul_hess_blk!( - wx::VT, - kkt::Union{MadNLP.SparseKKTSystem,MadNLP.SparseCondensedKKTSystem}, - t, -) where {T,VT<:oneVector{T}} + wx::VT, + kkt::Union{MadNLP.SparseKKTSystem, MadNLP.SparseCondensedKKTSystem}, + t, + ) where {T, VT <: oneVector{T}} n = size(kkt.hess_com, 1) wxx = @view(wx[1:n]) tx = @view(t[1:n]) @@ -97,15 +97,15 @@ function MadNLP.mul_hess_blk!( synchronize(backend) end - fill!(@view(wx[n+1:end]), 0) + fill!(@view(wx[(n + 1):end]), 0) wx .+= t .* kkt.pr_diag return end -function MadNLP.get_tril_to_full(csc::oneMKL.oneSparseMatrixCSC{Tv,Ti}) where {Tv,Ti} - cscind = MadNLP.SparseMatrixCSC{Int,Ti}( +function MadNLP.get_tril_to_full(csc::oneMKL.oneSparseMatrixCSC{Tv, Ti}) where {Tv, Ti} + cscind = MadNLP.SparseMatrixCSC{Int, Ti}( Symmetric( - MadNLP.SparseMatrixCSC{Int,Ti}( + MadNLP.SparseMatrixCSC{Int, Ti}( size(csc)..., Array(csc.colPtr), Array(csc.rowVal), @@ -114,22 +114,22 @@ function MadNLP.get_tril_to_full(csc::oneMKL.oneSparseMatrixCSC{Tv,Ti}) where {T :L, ), ) - return oneMKL.oneSparseMatrixCSC{Tv,Ti}( - oneArray(cscind.colptr), - oneArray(cscind.rowval), - oneVector{Tv}(undef, MadNLP.nnz(cscind)), - size(csc), - ), - view(csc.nzVal, oneArray(cscind.nzval)) + return oneMKL.oneSparseMatrixCSC{Tv, Ti}( + oneArray(cscind.colptr), + oneArray(cscind.rowval), + oneVector{Tv}(undef, MadNLP.nnz(cscind)), + size(csc), + ), + view(csc.nzVal, oneArray(cscind.nzval)) end function MadNLP.get_sparse_condensed_ext( - ::Type{VT}, - hess_com, - jptr, - jt_map, - hess_map, -) where {T,VT<:oneVector{T}} + ::Type{VT}, + hess_com, + jptr, + jt_map, + hess_map, + ) where {T, VT <: oneVector{T}} zvals = oneVector{Int}(1:length(hess_map)) hess_com_ptr = map((i, j) -> (i, j), hess_map, zvals) if length(hess_com_ptr) > 0 # otherwise error is thrown @@ -168,7 +168,7 @@ function get_diagonal_mapping(colptr, rowval) inds1 = findall( map( (x, y) -> ((x <= nnz) && (x != y)), - @view(colptr[1:end-1]), + @view(colptr[1:(end - 1)]), @view(colptr[2:end]) ), ) @@ -193,7 +193,7 @@ function MadNLP._sym_length(Jt::oneMKL.oneSparseMatrixCSC) end, +, @view(Jt.colPtr[2:end]), - @view(Jt.colPtr[1:end-1]) + @view(Jt.colPtr[1:(end - 1)]) ) end @@ -267,8 +267,8 @@ end =# function MadNLP.coo_to_csc( - coo::MadNLP.SparseMatrixCOO{T,I,VT,VI}, -) where {T,I,VT<:oneArray,VI<:oneArray} + coo::MadNLP.SparseMatrixCOO{T, I, VT, VI}, + ) where {T, I, VT <: oneArray, VI <: oneArray} zvals = oneVector{Int}(1:length(coo.I)) coord = map((i, j, k) -> ((i, j), k), coo.I, coo.J, zvals) if length(coord) > 0 @@ -279,7 +279,7 @@ function MadNLP.coo_to_csc( colptr = similar(coo.I, size(coo, 2) + 1) - coord_csc = coord[@view(mapptr[1:end-1])] + coord_csc = coord[@view(mapptr[1:(end - 1)])] backend = oneAPIBackend() if length(coord_csc) > 0 @@ -316,8 +316,8 @@ end =# function MadNLP.build_condensed_aug_coord!( - kkt::MadNLP.AbstractCondensedKKTSystem{T,VT,MT}, -) where {T,VT,MT<:oneMKL.oneSparseMatrixCSC{T}} + kkt::MadNLP.AbstractCondensedKKTSystem{T, VT, MT}, + ) where {T, VT, MT <: oneMKL.oneSparseMatrixCSC{T}} fill!(kkt.aug_com.nzVal, zero(T)) backend = oneAPIBackend() if length(kkt.hptr) > 0 @@ -357,8 +357,8 @@ end =# function MadNLP.compress_hessian!( - kkt::MadNLP.AbstractSparseKKTSystem{T,VT,MT}, -) where {T,VT,MT<:oneMKL.oneSparseMatrixCSC{T,Int32}} + kkt::MadNLP.AbstractSparseKKTSystem{T, VT, MT}, + ) where {T, VT, MT <: oneMKL.oneSparseMatrixCSC{T, Int32}} fill!(kkt.hess_com.nzVal, zero(T)) backend = oneAPIBackend() if length(kkt.ext.hess_com_ptrptr) > 1 @@ -375,8 +375,8 @@ function MadNLP.compress_hessian!( end function MadNLP.compress_jacobian!( - kkt::MadNLP.SparseCondensedKKTSystem{T,VT,MT}, -) where {T,VT,MT<:oneMKL.oneSparseMatrixCSC{T,Int32}} + kkt::MadNLP.SparseCondensedKKTSystem{T, VT, MT}, + ) where {T, VT, MT <: oneMKL.oneSparseMatrixCSC{T, Int32}} fill!(kkt.jt_csc.nzVal, zero(T)) backend = oneAPIBackend() if length(kkt.ext.jt_csc_ptrptr) > 1 # otherwise error is thrown @@ -397,10 +397,10 @@ end =# function MadNLP._set_con_scale_sparse!( - con_scale::VT, - jac_I, - jac_buffer, -) where {T,VT<:oneVector{T}} + con_scale::VT, + jac_I, + jac_buffer, + ) where {T, VT <: oneVector{T}} ind_jac = oneVector{Int}(1:length(jac_I)) inds = map((i, j) -> (i, j), jac_I, ind_jac) !isempty(inds) && sort!(inds) @@ -425,10 +425,10 @@ end =# function MadNLP._build_condensed_aug_symbolic_hess( - H::oneMKL.oneSparseMatrixCSC{Tv,Ti}, - sym, - sym2, -) where {Tv,Ti} + H::oneMKL.oneSparseMatrixCSC{Tv, Ti}, + sym, + sym2, + ) where {Tv, Ti} if size(H, 2) > 0 backend = oneAPIBackend() MadNLPGPU._build_condensed_aug_symbolic_hess_kernel!(backend)( @@ -448,14 +448,14 @@ end =# function MadNLP._build_condensed_aug_symbolic_jt( - Jt::oneMKL.oneSparseMatrixCSC{Tv,Ti}, - sym, - sym2, -) where {Tv,Ti} + Jt::oneMKL.oneSparseMatrixCSC{Tv, Ti}, + sym, + sym2, + ) where {Tv, Ti} if size(Jt, 2) > 0 _offsets = map( (i, j) -> div((j - i)^2 + (j - i), 2), - @view(Jt.colPtr[1:end-1]), + @view(Jt.colPtr[1:(end - 1)]), @view(Jt.colPtr[2:end]) ) offsets = cumsum(_offsets) diff --git a/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/onemkl.jl b/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/onemkl.jl index ec1bb29a7..74676c64f 100644 --- a/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/onemkl.jl +++ b/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/onemkl.jl @@ -1,4 +1,4 @@ -mutable struct LapackOneMKLSolver{T,MT} <: MadNLP.AbstractLinearSolver{T} +mutable struct LapackOneMKLSolver{T, MT} <: MadNLP.AbstractLinearSolver{T} A::MT fact::oneMatrix{T} n::Int64 @@ -16,15 +16,15 @@ mutable struct LapackOneMKLSolver{T,MT} <: MadNLP.AbstractLinearSolver{T} logger::MadNLP.MadNLPLogger function LapackOneMKLSolver( - A::MT; - option_dict::Dict{Symbol,Any} = Dict{Symbol,Any}(), - opt = MadNLP.LapackOptions(), - logger = MadNLP.MadNLPLogger(), - kwargs..., - ) where {MT<:AbstractMatrix} + A::MT; + option_dict::Dict{Symbol, Any} = Dict{Symbol, Any}(), + opt = MadNLP.LapackOptions(), + logger = MadNLP.MadNLPLogger(), + kwargs..., + ) where {MT <: AbstractMatrix} MadNLP.set_options!(opt, option_dict, kwargs...) T = eltype(A) - m,n = size(A) + m, n = size(A) @assert m == n fact = oneMatrix{T}(undef, m, n) sol = oneVector{T}(undef, 0) @@ -39,7 +39,7 @@ mutable struct LapackOneMKLSolver{T,MT} <: MadNLP.AbstractLinearSolver{T} device_queue = oneAPI.sycl_queue(queue) alpha = Ref{T}(1) beta = Ref{T}(0) - solver = new{T,MT}(A, fact, n, sol, tau, Λ, info, ipiv, scratchpad, scratchpad_size, device_queue, alpha, beta, opt, logger) + solver = new{T, MT}(A, fact, n, sol, tau, Λ, info, ipiv, scratchpad, scratchpad_size, device_queue, alpha, beta, opt, logger) setup!(solver) return solver end @@ -48,7 +48,7 @@ end MadNLP.improve!(M::LapackOneMKLSolver) = false MadNLP.is_inertia(M::LapackOneMKLSolver) = (M.opt.lapack_algorithm == MadNLP.CHOLESKY) || (M.opt.lapack_algorithm == MadNLP.EVD) function MadNLP.inertia(M::LapackOneMKLSolver) - if M.opt.lapack_algorithm == MadNLP.CHOLESKY + return if M.opt.lapack_algorithm == MadNLP.CHOLESKY sum(M.info) == 0 ? (M.n, 0, 0) : (0, M.n, 0) elseif M.opt.lapack_algorithm == MadNLP.EVD numpos = count(λ -> λ > 0, M.Λ) @@ -66,7 +66,7 @@ MadNLP.introduce(M::LapackOneMKLSolver) = "OneAPI -- ($(M.opt.lapack_algorithm)) # MadNLP.introduce(M::LapackOneMKLSolver) = "OneAPI v$(oneAPI.version()) -- ($(M.opt.lapack_algorithm))" function setup!(M::LapackOneMKLSolver) - if M.opt.lapack_algorithm == MadNLP.LU + return if M.opt.lapack_algorithm == MadNLP.LU setup_lu!(M) elseif M.opt.lapack_algorithm == MadNLP.QR setup_qr!(M) @@ -81,7 +81,7 @@ end function MadNLP.factorize!(M::LapackOneMKLSolver) MadNLPGPU.gpu_transfer!(M.fact, M.A) - if M.opt.lapack_algorithm == MadNLP.LU + return if M.opt.lapack_algorithm == MadNLP.LU MadNLP.tril_to_full!(M.fact) factorize_lu!(M) elseif M.opt.lapack_algorithm == MadNLP.QR @@ -99,7 +99,7 @@ end for T in (:Float32, :Float64) @eval begin function MadNLP.solve!(M::LapackOneMKLSolver{$T}, x::oneVector{$T}) - if M.opt.lapack_algorithm == MadNLP.LU + return if M.opt.lapack_algorithm == MadNLP.LU solve_lu!(M, x) elseif M.opt.lapack_algorithm == MadNLP.QR solve_qr!(M, x) @@ -125,8 +125,10 @@ function MadNLP.solve!(M::LapackOneMKLSolver, x::AbstractVector) end for (potrf, potrf_buffer, potrs, potrs_buffer, T) in - ((:onemklDpotrf, :onemklDpotrf_scratchpad_size, :onemklDpotrs, :onemklDpotrs_scratchpad_size, :Float64), - (:onemklSpotrf, :onemklSpotrf_scratchpad_size, :onemklSpotrs, :onemklSpotrs_scratchpad_size, :Float32)) + ( + (:onemklDpotrf, :onemklDpotrf_scratchpad_size, :onemklDpotrs, :onemklDpotrs_scratchpad_size, :Float64), + (:onemklSpotrf, :onemklSpotrf_scratchpad_size, :onemklSpotrs, :onemklSpotrs_scratchpad_size, :Float32), + ) @eval begin function setup_cholesky!(M::LapackOneMKLSolver{$T}) potrf_scratchpad_size = Support.$potrf_buffer(M.device_queue, 'L', M.n, M.n) @@ -168,8 +170,10 @@ for (potrf, potrf_buffer, potrs, potrs_buffer, T) in end for (getrf, getrf_buffer, getrs, getrs_buffer, T) in - ((:onemklDgetrf, :onemklDgetrf_scratchpad_size, :onemklDgetrs, :onemklDgetrs_scratchpad_size, :Float64), - (:onemklSgetrf, :onemklSgetrf_scratchpad_size, :onemklSgetrs, :onemklSgetrs_scratchpad_size, :Float32)) + ( + (:onemklDgetrf, :onemklDgetrf_scratchpad_size, :onemklDgetrs, :onemklDgetrs_scratchpad_size, :Float64), + (:onemklSgetrf, :onemklSgetrf_scratchpad_size, :onemklSgetrs, :onemklSgetrs_scratchpad_size, :Float32), + ) @eval begin function setup_lu!(M::LapackOneMKLSolver{$T}) resize!(M.ipiv, M.n) @@ -214,8 +218,10 @@ for (getrf, getrf_buffer, getrs, getrs_buffer, T) in end for (geqrf, geqrf_buffer, ormqr, ormqr_buffer, trsv, T) in - ((:onemklDgeqrf, :onemklDgeqrf_scratchpad_size, :onemklDormqr, :onemklDormqr_scratchpad_size, :onemklDtrsv, :Float64), - (:onemklSgeqrf, :onemklSgeqrf_scratchpad_size, :onemklSormqr, :onemklSormqr_scratchpad_size, :onemklStrsv, :Float32)) + ( + (:onemklDgeqrf, :onemklDgeqrf_scratchpad_size, :onemklDormqr, :onemklDormqr_scratchpad_size, :onemklDtrsv, :Float64), + (:onemklSgeqrf, :onemklSgeqrf_scratchpad_size, :onemklSormqr, :onemklSormqr_scratchpad_size, :onemklStrsv, :Float32), + ) @eval begin function setup_qr!(M::LapackOneMKLSolver{$T}) resize!(M.tau, M.n) @@ -265,8 +271,10 @@ for (geqrf, geqrf_buffer, ormqr, ormqr_buffer, trsv, T) in end for (syevd, syevd_buffer, gemv, T) in - ((:onemklDsyevd, :onemklDsyevd_scratchpad_size, :onemklDgemv, :Float64), - (:onemklSsyevd, :onemklSsyevd_scratchpad_size, :onemklSgemv, :Float32)) + ( + (:onemklDsyevd, :onemklDsyevd_scratchpad_size, :onemklDgemv, :Float64), + (:onemklSsyevd, :onemklSsyevd_scratchpad_size, :onemklSgemv, :Float32), + ) @eval begin function setup_evd!(M::LapackOneMKLSolver{$T}) resize!(M.tau, M.n) diff --git a/lib/MadNLPGPU/test/densekkt_oneapi.jl b/lib/MadNLPGPU/test/densekkt_oneapi.jl index 896611201..1c137dc5c 100644 --- a/lib/MadNLPGPU/test/densekkt_oneapi.jl +++ b/lib/MadNLPGPU/test/densekkt_oneapi.jl @@ -2,23 +2,23 @@ using oneAPI using MadNLPTests function _compare_oneapi_with_cpu(KKTSystem, n, m, ind_fixed) - for (T,tol,atol) in [ - (Float32,1e-4,1e1), - (Float64,1e-8,1e-6) + for (T, tol, atol) in [ + (Float32, 1.0e-4, 1.0e1), + (Float64, 1.0e-8, 1.0e-6), ] madnlp_options = Dict{Symbol, Any}( - :callback=>MadNLP.DenseCallback, - :kkt_system=>KKTSystem, - :linear_solver=>LapackOneMKLSolver, - :lapack_algorithm=>MadNLP.QR, - :print_level=>MadNLP.ERROR, - :tol=>tol + :callback => MadNLP.DenseCallback, + :kkt_system => KKTSystem, + :linear_solver => LapackOneMKLSolver, + :lapack_algorithm => MadNLP.QR, + :print_level => MadNLP.ERROR, + :tol => tol ) # Host evaluator - nlph = MadNLPTests.DenseDummyQP(zeros(T,n); m=m, fixed_variables=ind_fixed) + nlph = MadNLPTests.DenseDummyQP(zeros(T, n); m = m, fixed_variables = ind_fixed) # Device evaluator - nlpd = MadNLPTests.DenseDummyQP(oneAPI.zeros(T,n); m=m, fixed_variables=oneArray(ind_fixed)) + nlpd = MadNLPTests.DenseDummyQP(oneAPI.zeros(T, n); m = m, fixed_variables = oneArray(ind_fixed)) # Solve on CPU h_solver = MadNLPSolver(nlph; madnlp_options...) @@ -33,10 +33,11 @@ function _compare_oneapi_with_cpu(KKTSystem, n, m, ind_fixed) if T == Float64 @test h_solver.cnt.k == d_solver.cnt.k @test results_cpu.objective ≈ results_gpu.objective - @test results_cpu.solution ≈ Array(results_gpu.solution) atol=atol - @test results_cpu.multipliers ≈ Array(results_gpu.multipliers) atol=atol - end + @test results_cpu.solution ≈ Array(results_gpu.solution) atol = atol + @test results_cpu.multipliers ≈ Array(results_gpu.multipliers) atol = atol + end end + return end @testset "MadNLPGPU -- LapackOneMKLSolver -- ($(kkt_system))" for kkt_system in [ @@ -46,43 +47,43 @@ end @testset "Size: ($n, $m)" for (n, m) in [(10, 0), (10, 5), (50, 10)] _compare_oneapi_with_cpu(kkt_system, n, m, Int[]) end - @testset "Fixed variables" for (n,m) in [(10, 0), (10, 5), (50, 10)] + @testset "Fixed variables" for (n, m) in [(10, 0), (10, 5), (50, 10)] _compare_oneapi_with_cpu(kkt_system, n, m, Int[1, 2]) end end @testset "MadNLP -- LapackOneMKLSolver: $QN + $KKT" for QN in [ - MadNLP.BFGS, - MadNLP.DampedBFGS, -], KKT in [ - MadNLP.DenseKKTSystem, - MadNLP.DenseCondensedKKTSystem, -] + MadNLP.BFGS, + MadNLP.DampedBFGS, + ], KKT in [ + MadNLP.DenseKKTSystem, + MadNLP.DenseCondensedKKTSystem, + ] @testset "Size: ($n, $m)" for (n, m) in [(10, 0), (10, 5), (50, 10)] - nlp = MadNLPTests.DenseDummyQP(zeros(Float64, n); m=m) + nlp = MadNLPTests.DenseDummyQP(zeros(Float64, n); m = m) solver_exact = MadNLPSolver( nlp; - callback=MadNLP.DenseCallback, - print_level=MadNLP.ERROR, - kkt_system=KKT, - linear_solver=LapackOneMKLSolver, + callback = MadNLP.DenseCallback, + print_level = MadNLP.ERROR, + kkt_system = KKT, + linear_solver = LapackOneMKLSolver, ) results_ref = MadNLP.solve!(solver_exact) - nlp = MadNLPTests.DenseDummyQP(oneAPI.zeros(Float64, n); m=m) + nlp = MadNLPTests.DenseDummyQP(oneAPI.zeros(Float64, n); m = m) solver_qn = MadNLPSolver( nlp; - callback=MadNLP.DenseCallback, - print_level=MadNLP.ERROR, - kkt_system=KKT, - hessian_approximation=QN, - linear_solver=LapackOneMKLSolver, + callback = MadNLP.DenseCallback, + print_level = MadNLP.ERROR, + kkt_system = KKT, + hessian_approximation = QN, + linear_solver = LapackOneMKLSolver, ) results_qn = MadNLP.solve!(solver_qn) @test results_qn.status == MadNLP.SOLVE_SUCCEEDED - @test results_qn.objective ≈ results_ref.objective atol=1e-6 - @test Array(results_qn.solution) ≈ Array(results_ref.solution) atol=1e-6 + @test results_qn.objective ≈ results_ref.objective atol = 1.0e-6 + @test Array(results_qn.solution) ≈ Array(results_ref.solution) atol = 1.0e-6 @test solver_qn.cnt.lag_hess_cnt == 0 end end diff --git a/lib/MadNLPGPU/test/madnlpgpu_test.jl b/lib/MadNLPGPU/test/madnlpgpu_test.jl index d0772e2ef..32ece9d22 100644 --- a/lib/MadNLPGPU/test/madnlpgpu_test.jl +++ b/lib/MadNLPGPU/test/madnlpgpu_test.jl @@ -178,10 +178,10 @@ rocm_testset = [ oneapi_testset = [ [ "LapackOneMKLSolver-LU", - ()->MadNLP.Optimizer( - linear_solver=LapackOneMKLSolver, - lapack_algorithm=MadNLP.LU, - print_level=MadNLP.ERROR, + () -> MadNLP.Optimizer( + linear_solver = LapackOneMKLSolver, + lapack_algorithm = MadNLP.LU, + print_level = MadNLP.ERROR, ), [], ], @@ -214,17 +214,17 @@ oneapi_testset = [ end end if AMDGPU.functional() - MadNLPTests.test_linear_solver(LapackROCSolver,Float32) - MadNLPTests.test_linear_solver(LapackROCSolver,Float64) + MadNLPTests.test_linear_solver(LapackROCSolver, Float32) + MadNLPTests.test_linear_solver(LapackROCSolver, Float64) for (name,optimizer_constructor,exclude) in rocm_testset test_madnlp(name,optimizer_constructor,exclude; Arr=ROCArray) end end if oneAPI.functional() - MadNLPTests.test_linear_solver(LapackOneMKLSolver,Float32) - MadNLPTests.test_linear_solver(LapackOneMKLSolver,Float64) - for (name,optimizer_constructor,exclude) in oneapi_testset - test_madnlp(name,optimizer_constructor,exclude; Arr=oneArray) + MadNLPTests.test_linear_solver(LapackOneMKLSolver, Float32) + MadNLPTests.test_linear_solver(LapackOneMKLSolver, Float64) + for (name, optimizer_constructor, exclude) in oneapi_testset + test_madnlp(name, optimizer_constructor, exclude; Arr = oneArray) end end end From 4727a17c26d09534bec2b2e859bc4d407c4e8bfc Mon Sep 17 00:00:00 2001 From: Github action runner Date: Wed, 4 Feb 2026 11:36:02 -0800 Subject: [PATCH 12/18] Fix Project.toml --- lib/MadNLPGPU/Project.toml | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/lib/MadNLPGPU/Project.toml b/lib/MadNLPGPU/Project.toml index 883110baf..90eb5e0a0 100644 --- a/lib/MadNLPGPU/Project.toml +++ b/lib/MadNLPGPU/Project.toml @@ -12,10 +12,10 @@ LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e" MadNLP = "2621e9c9-9eb4-46b1-8089-e8c72242dfb6" Metis = "2679e427-3c69-5b7f-982b-ece356f1e94b" SparseArrays = "2f01184e-e22b-5df5-ae63-d93ebab69eaf" -oneAPI = "8f75cd03-7ff8-4ecb-9b8f-daf728133b1b" [weakdeps] AMDGPU = "21141c5a-9bdb-4563-92ae-f87d6854732e" +oneAPI = "8f75cd03-7ff8-4ecb-9b8f-daf728133b1b" [extensions] MadNLPGPUAMDGPUExt = "AMDGPU" @@ -39,4 +39,4 @@ MadNLPTests = "b52a2a03-04ab-4a5f-9698-6a2deff93217" Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40" [targets] -test = ["Test", "MadNLPTests", "AMDGPU"] +test = ["Test", "MadNLPTests", "AMDGPU", "oneAPI"] From 050a2ae1f0591d59fb6255268c2da1e4492a2e9c Mon Sep 17 00:00:00 2001 From: Github action runner Date: Wed, 4 Feb 2026 12:22:04 -0800 Subject: [PATCH 13/18] Split oneAPI and AMDGPU package loading in tests --- .github/workflows/test.yml | 14 ++++++++++++++ lib/MadNLPGPU/Project.toml | 9 +-------- lib/MadNLPGPU/test/Project.toml | 6 ++++++ lib/MadNLPGPU/test/madnlpgpu_test.jl | 26 +++++++++++++------------- lib/MadNLPGPU/test/runtests.jl | 18 ++++++++++++++---- 5 files changed, 48 insertions(+), 25 deletions(-) create mode 100644 lib/MadNLPGPU/test/Project.toml diff --git a/.github/workflows/test.yml b/.github/workflows/test.yml index 703fc457d..45648e0d9 100644 --- a/.github/workflows/test.yml +++ b/.github/workflows/test.yml @@ -118,8 +118,22 @@ jobs: with: channel: ${{ matrix.julia-version }} - uses: julia-actions/cache@v1 + - name: Add AMDGPU + if: matrix.gpu == 'amdgpu' + shell: julia --project=./lib/MadNLPGPU --color=yes {0} + run: | + using Pkg + Pkg.add("AMDGPU") + - name: Add oneAPI + if: matrix.gpu == 'oneapi' + shell: julia --project=./lib/MadNLPGPU --color=yes {0} + run: | + using Pkg + Pkg.add("oneAPI") - name: test MadNLPGPU shell: julia --project=./lib/MadNLPGPU --color=yes {0} + env: + MADNLP_GPU_BACKEND: ${{ matrix.gpu }} run: | using Pkg Pkg.Registry.update() diff --git a/lib/MadNLPGPU/Project.toml b/lib/MadNLPGPU/Project.toml index 90eb5e0a0..39e207c18 100644 --- a/lib/MadNLPGPU/Project.toml +++ b/lib/MadNLPGPU/Project.toml @@ -32,11 +32,4 @@ MadNLP = "0.8.12" MadNLPTests = "0.5.3" Metis = "1" julia = "1.10" -oneAPI = "2.2.0" - -[extras] -MadNLPTests = "b52a2a03-04ab-4a5f-9698-6a2deff93217" -Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40" - -[targets] -test = ["Test", "MadNLPTests", "AMDGPU", "oneAPI"] +oneAPI = "2.6.0" diff --git a/lib/MadNLPGPU/test/Project.toml b/lib/MadNLPGPU/test/Project.toml new file mode 100644 index 000000000..24bfd4e1e --- /dev/null +++ b/lib/MadNLPGPU/test/Project.toml @@ -0,0 +1,6 @@ +[deps] +CUDA = "052768ef-5323-5732-b1bb-66c8b64840ba" +MadNLP = "2621e9c9-9eb4-46b1-8089-e8c72242dfb6" +MadNLPGPU = "d72a61cc-809d-412f-99be-fd81f4b8a598" +MadNLPTests = "b52a2a03-04ab-4a5f-9698-6a2deff93217" +Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40" diff --git a/lib/MadNLPGPU/test/madnlpgpu_test.jl b/lib/MadNLPGPU/test/madnlpgpu_test.jl index 32ece9d22..b5fff98e6 100644 --- a/lib/MadNLPGPU/test/madnlpgpu_test.jl +++ b/lib/MadNLPGPU/test/madnlpgpu_test.jl @@ -178,10 +178,10 @@ rocm_testset = [ oneapi_testset = [ [ "LapackOneMKLSolver-LU", - () -> MadNLP.Optimizer( - linear_solver = LapackOneMKLSolver, - lapack_algorithm = MadNLP.LU, - print_level = MadNLP.ERROR, + ()->MadNLP.Optimizer( + linear_solver=LapackOneMKLSolver, + lapack_algorithm=MadNLP.LU, + print_level=MadNLP.ERROR, ), [], ], @@ -205,7 +205,7 @@ oneapi_testset = [ # ], ] @testset "MadNLPGPU test" begin - if CUDA.functional() + if CUDA.functional() && GPU_BACKEND == "cuda" MadNLPTests.test_linear_solver(LapackCUDASolver,Float32) MadNLPTests.test_linear_solver(LapackCUDASolver,Float64) # Test LapackGPU wrapper @@ -213,18 +213,18 @@ oneapi_testset = [ test_madnlp(name,optimizer_constructor,exclude; Arr=CuArray) end end - if AMDGPU.functional() - MadNLPTests.test_linear_solver(LapackROCSolver, Float32) - MadNLPTests.test_linear_solver(LapackROCSolver, Float64) + if GPU_BACKEND == "amdgpu" && AMDGPU.functional() + MadNLPTests.test_linear_solver(LapackROCmSolver,Float32) + MadNLPTests.test_linear_solver(LapackROCmSolver,Float64) for (name,optimizer_constructor,exclude) in rocm_testset test_madnlp(name,optimizer_constructor,exclude; Arr=ROCArray) end end - if oneAPI.functional() - MadNLPTests.test_linear_solver(LapackOneMKLSolver, Float32) - MadNLPTests.test_linear_solver(LapackOneMKLSolver, Float64) - for (name, optimizer_constructor, exclude) in oneapi_testset - test_madnlp(name, optimizer_constructor, exclude; Arr = oneArray) + if GPU_BACKEND == "oneapi" && oneAPI.functional() + MadNLPTests.test_linear_solver(LapackOneMKLSolver,Float32) + MadNLPTests.test_linear_solver(LapackOneMKLSolver,Float64) + for (name,optimizer_constructor,exclude) in oneapi_testset + test_madnlp(name,optimizer_constructor,exclude; Arr=oneArray) end end end diff --git a/lib/MadNLPGPU/test/runtests.jl b/lib/MadNLPGPU/test/runtests.jl index 0a6cbe337..683321385 100644 --- a/lib/MadNLPGPU/test/runtests.jl +++ b/lib/MadNLPGPU/test/runtests.jl @@ -1,18 +1,28 @@ -using Test, CUDA, AMDGPU, oneAPI, MadNLP, MadNLPGPU, MadNLPTests +using Test, CUDA, MadNLP, MadNLPGPU, MadNLPTests + +# Get backend from environment +const GPU_BACKEND = get(ENV, "MADNLP_GPU_BACKEND", "cuda") + +# Conditionally load GPU backends +if GPU_BACKEND == "amdgpu" + using AMDGPU +elseif GPU_BACKEND == "oneapi" + using oneAPI +end @testset "MadNLPGPU test" begin include("madnlpgpu_test.jl") - if CUDA.functional() + if GPU_BACKEND == "cuda" && CUDA.functional() include("densekkt_cuda.jl") # Need to add support for CompactLBFGS in SparseCondensedKKTSystem (Issue #563) # include("sparsekkt_cuda.jl") end - if AMDGPU.functional() + if GPU_BACKEND == "amdgpu" && AMDGPU.functional() include("densekkt_rocm.jl") # Need to add support for CompactLBFGS in SparseCondensedKKTSystem (Issue #563) # include("sparsekkt_rocm.jl") end - if oneAPI.functional() + if GPU_BACKEND == "oneapi" && oneAPI.functional() include("densekkt_oneapi.jl") end end From 064b2b2ffacab68fc295be1a0f46ea1ee002c0e0 Mon Sep 17 00:00:00 2001 From: Github action runner Date: Wed, 4 Feb 2026 12:26:05 -0800 Subject: [PATCH 14/18] Fix --- lib/MadNLPGPU/Project.toml | 1 - lib/MadNLPGPU/test/Project.toml | 3 +++ 2 files changed, 3 insertions(+), 1 deletion(-) diff --git a/lib/MadNLPGPU/Project.toml b/lib/MadNLPGPU/Project.toml index 39e207c18..10a52c7c8 100644 --- a/lib/MadNLPGPU/Project.toml +++ b/lib/MadNLPGPU/Project.toml @@ -29,7 +29,6 @@ CUDSS = "0.6.4" GPUArraysCore = "0.2" KernelAbstractions = "0.9" MadNLP = "0.8.12" -MadNLPTests = "0.5.3" Metis = "1" julia = "1.10" oneAPI = "2.6.0" diff --git a/lib/MadNLPGPU/test/Project.toml b/lib/MadNLPGPU/test/Project.toml index 24bfd4e1e..c1c2a45b0 100644 --- a/lib/MadNLPGPU/test/Project.toml +++ b/lib/MadNLPGPU/test/Project.toml @@ -4,3 +4,6 @@ MadNLP = "2621e9c9-9eb4-46b1-8089-e8c72242dfb6" MadNLPGPU = "d72a61cc-809d-412f-99be-fd81f4b8a598" MadNLPTests = "b52a2a03-04ab-4a5f-9698-6a2deff93217" Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40" + +[compat] +MadNLPTests = "0.5.3" From 815bf063e34ed2b0af10e550768c23165bde4ad7 Mon Sep 17 00:00:00 2001 From: Github action runner Date: Wed, 4 Feb 2026 12:26:21 -0800 Subject: [PATCH 15/18] Format --- lib/MadNLPGPU/test/madnlpgpu_test.jl | 16 ++++++++-------- 1 file changed, 8 insertions(+), 8 deletions(-) diff --git a/lib/MadNLPGPU/test/madnlpgpu_test.jl b/lib/MadNLPGPU/test/madnlpgpu_test.jl index b5fff98e6..21e701507 100644 --- a/lib/MadNLPGPU/test/madnlpgpu_test.jl +++ b/lib/MadNLPGPU/test/madnlpgpu_test.jl @@ -178,10 +178,10 @@ rocm_testset = [ oneapi_testset = [ [ "LapackOneMKLSolver-LU", - ()->MadNLP.Optimizer( - linear_solver=LapackOneMKLSolver, - lapack_algorithm=MadNLP.LU, - print_level=MadNLP.ERROR, + () -> MadNLP.Optimizer( + linear_solver = LapackOneMKLSolver, + lapack_algorithm = MadNLP.LU, + print_level = MadNLP.ERROR, ), [], ], @@ -221,10 +221,10 @@ oneapi_testset = [ end end if GPU_BACKEND == "oneapi" && oneAPI.functional() - MadNLPTests.test_linear_solver(LapackOneMKLSolver,Float32) - MadNLPTests.test_linear_solver(LapackOneMKLSolver,Float64) - for (name,optimizer_constructor,exclude) in oneapi_testset - test_madnlp(name,optimizer_constructor,exclude; Arr=oneArray) + MadNLPTests.test_linear_solver(LapackOneMKLSolver, Float32) + MadNLPTests.test_linear_solver(LapackOneMKLSolver, Float64) + for (name, optimizer_constructor, exclude) in oneapi_testset + test_madnlp(name, optimizer_constructor, exclude; Arr = oneArray) end end end From f867c6856b1a308d05a9c62c0cdd2b3c3cc9f841 Mon Sep 17 00:00:00 2001 From: Github action runner Date: Wed, 4 Feb 2026 12:33:43 -0800 Subject: [PATCH 16/18] Fix tests --- .github/workflows/test.yml | 19 ++++++++++++------- 1 file changed, 12 insertions(+), 7 deletions(-) diff --git a/.github/workflows/test.yml b/.github/workflows/test.yml index 45648e0d9..375d6a569 100644 --- a/.github/workflows/test.yml +++ b/.github/workflows/test.yml @@ -118,27 +118,32 @@ jobs: with: channel: ${{ matrix.julia-version }} - uses: julia-actions/cache@v1 + - name: Setup test environment + shell: julia --project=./lib/MadNLPGPU/test --color=yes {0} + run: | + using Pkg + Pkg.Registry.update() + Pkg.develop(path=".") + Pkg.develop(path="./lib/MadNLPGPU") + Pkg.develop(path="./lib/MadNLPTests") - name: Add AMDGPU if: matrix.gpu == 'amdgpu' - shell: julia --project=./lib/MadNLPGPU --color=yes {0} + shell: julia --project=./lib/MadNLPGPU/test --color=yes {0} run: | using Pkg Pkg.add("AMDGPU") - name: Add oneAPI if: matrix.gpu == 'oneapi' - shell: julia --project=./lib/MadNLPGPU --color=yes {0} + shell: julia --project=./lib/MadNLPGPU/test --color=yes {0} run: | using Pkg Pkg.add("oneAPI") - name: test MadNLPGPU - shell: julia --project=./lib/MadNLPGPU --color=yes {0} + shell: julia --project=./lib/MadNLPGPU/test --color=yes {0} env: MADNLP_GPU_BACKEND: ${{ matrix.gpu }} run: | - using Pkg - Pkg.Registry.update() - Pkg.develop(path=".") - Pkg.test("MadNLPGPU", coverage=true) + include("lib/MadNLPGPU/test/runtests.jl") - uses: julia-actions/julia-processcoverage@v1 with: directories: lib/MadNLPGPU/src From 910553620002ac98331a813d4ffb33659f190da8 Mon Sep 17 00:00:00 2001 From: Github action runner Date: Wed, 4 Feb 2026 12:36:24 -0800 Subject: [PATCH 17/18] Fix tests --- .github/workflows/test.yml | 4 +--- 1 file changed, 1 insertion(+), 3 deletions(-) diff --git a/.github/workflows/test.yml b/.github/workflows/test.yml index 375d6a569..7fc474e04 100644 --- a/.github/workflows/test.yml +++ b/.github/workflows/test.yml @@ -139,11 +139,9 @@ jobs: using Pkg Pkg.add("oneAPI") - name: test MadNLPGPU - shell: julia --project=./lib/MadNLPGPU/test --color=yes {0} env: MADNLP_GPU_BACKEND: ${{ matrix.gpu }} - run: | - include("lib/MadNLPGPU/test/runtests.jl") + run: julia --project=./lib/MadNLPGPU/test --color=yes -e 'include("lib/MadNLPGPU/test/runtests.jl")' - uses: julia-actions/julia-processcoverage@v1 with: directories: lib/MadNLPGPU/src From 96689aaba6bc5d88c98711d6b567eb1ddc127c9c Mon Sep 17 00:00:00 2001 From: Michel Schanen Date: Wed, 11 Feb 2026 14:16:46 -0600 Subject: [PATCH 18/18] potrf fix, memory reclaim fix --- .../ext/MadNLPGPUOneAPIExt/oneapi.jl | 39 +++++++++++-------- .../ext/MadNLPGPUOneAPIExt/onemkl.jl | 12 +++++- 2 files changed, 34 insertions(+), 17 deletions(-) diff --git a/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/oneapi.jl b/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/oneapi.jl index b09a84a47..5b4fc7e28 100644 --- a/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/oneapi.jl +++ b/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/oneapi.jl @@ -77,23 +77,30 @@ MadNLP._trsm!(side::Char, uplo::Char, transa::Char, diag::Char, alpha::T, A::one MadNLP._dgmm!(side::Char, A::oneMatrix{T}, x::oneVector{T}, B::oneMatrix{T}) where {T} = oneMKL.dgmm!(side, A, x, B) #= - MadNLP.get_sd / MadNLP.get_sc - Workaround for norm() on views of oneArray triggering scalar indexing in Julia 1.12+ - See https://github.com/JuliaGPU/CUDA.jl/issues/2811 for similar issue in CUDA.jl + LinearAlgebra.norm for oneAPI arrays. + GPUArrays' norm uses LinearAlgebra.norm as a map function in mapreduce, + which fails on oneAPI because norm on scalars generates jl_f_throw_methoderror + calls that can't be compiled to SPIRV. Replace with abs-based implementation. + The SubArray method also avoids scalar indexing on views. =# -if VERSION > v"1.11" - function MadNLP.get_sd(l::oneVector{T}, zl_r, zu_r, s_max) where {T} - return max( - s_max, - (my1norm(l) + my1norm(zl_r) + my1norm(zu_r)) / max(1, (length(l) + length(zl_r) + length(zu_r))), - ) / s_max +function _onenorm(v, p::Real) + isempty(v) && return float(zero(eltype(v))) + if p == Inf + return maximum(abs, v) + elseif p == -Inf + return minimum(abs, v) + elseif p == 1 + return mapreduce(abs, +, v) + elseif p == 2 + return sqrt(mapreduce(abs2, +, v)) + elseif p == 0 + return float(count(!iszero, v)) + else + spp = float(p) + return mapreduce(x -> abs(x)^spp, +, v)^inv(spp) end - function MadNLP.get_sc(zl_r::SubArray{T, 1, VT}, zu_r, s_max) where {T, VT <: oneVector{T}} - return max( - s_max, - (my1norm(zl_r) + my1norm(zu_r)) / max(1, length(zl_r) + length(zu_r)), - ) / s_max - end - my1norm(x) = mapreduce(abs, +, x) end + +LinearAlgebra.norm(v::oneArray{<:Number}, p::Real = 2) = _onenorm(v, p) +LinearAlgebra.norm(v::SubArray{<:Number, <:Any, <:oneArray}, p::Real = 2) = _onenorm(v, p) diff --git a/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/onemkl.jl b/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/onemkl.jl index 74676c64f..fd76e665f 100644 --- a/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/onemkl.jl +++ b/lib/MadNLPGPU/ext/MadNLPGPUOneAPIExt/onemkl.jl @@ -22,6 +22,16 @@ mutable struct LapackOneMKLSolver{T, MT} <: MadNLP.AbstractLinearSolver{T} logger = MadNLP.MadNLPLogger(), kwargs..., ) where {MT <: AbstractMatrix} + # Reclaim GPU resources from previous solvers. Julia's GC doesn't see GPU + # memory pressure, so old oneArray/MKL handles accumulate. + # Synchronize FIRST to ensure all pending GPU operations complete, making + # it safe for GC finalizers to free GPU buffers. Then GC to collect stale + # objects, flush deferred MKL sparse handles, and GC again. + oneAPI.synchronize() + GC.gc(true) + oneAPI.oneL0._run_reclaim_callbacks() + GC.gc(true) + oneAPI.synchronize() MadNLP.set_options!(opt, option_dict, kwargs...) T = eltype(A) m, n = size(A) @@ -155,7 +165,7 @@ for (potrf, potrf_buffer, potrs, potrs_buffer, T) in Support.$potrs( M.device_queue, 'L', - n, + M.n, one(Int64), M.fact, M.n,