Compatibility with CUDA.jl about quadgk.jl HOT 2 OPEN

It seems that this issue didn't attention from the devs. @weiyuanlou were you able to make this work?

from quadgk.jl.

Hi,
I don't think f_int.(Zcu), which is solving multiple integrals in parallel, can be currently solved with QuadGK.jl. The adaptation logic in the code will not port well to the GPU.

You could attempt using quadgk! to evaluate an inplace integrand that evaluates all of your different integrands on the GPU in an embarrasingly parallel fashion.

using CUDA
using QuadGK

Z = Array([z for z in 1:1:4]);
Zcu = CuArray([z for z in 1:1:4]);

function f_int(z_ob) 

    iii(z) =  z_ob - z  
    
    return quadgk(z -> iii(z), 0, 1, rtol=1e-4)[1]
    
end

f_int.(Z)  # This works 

quadgk!((y, z) -> y .= Zcu .- z, similar(Zcu), 0, 1, rtol=1e-4) # try this

I haven't tested this (I don't have a gpu), and I think it should parallelize well for large Zcu.

Now, with version 2.9 you will be able to use a BatchIntegrand to parallelize the evaluation of a single integral. That pr (#80) should allow, for example BatchIntegrand((y, x) -> y .= f.(x), CuArray{Float64}(undef, 1000), CuVector{Float64}(undef, 1000); max_batch=1000). This is not embarrasingly parallel since it only parallelizes within each adaptive step

QuadGK.jl is missing support for batched+inplace integrands, which would allow for more parallelism, but it would be worth trying to see how much mileage you can get out of either of the above approaches first.

from quadgk.jl.