Move KomaCore to use KernelAbstractions.jl (to support CUDA, AMD, Metal, etc)

[cncastillo]

I started doing the first changes to support KernelAbstractions.jl.

At the same time, I also changed how the GPU drivers are loaded using Package extensions, to increase the loading speed.

These changes need CUDA 4.4. I included it in the compat, but while testing I saw that KomaMRI loads CUDA 3.0 sometimes. Maybe KomaMRI also needs an extension for the compat.

[codecov]

Codecov Report

Merging #179 (b13b9c7) into master (c7c3f24) will decrease coverage by 0.80%.
Report is 8 commits behind head on master.
The diff coverage is 35.08%.

Additional details and impacted files

@@            Coverage Diff             @@
##           master     #179      +/-   ##
==========================================
- Coverage   92.45%   91.66%   -0.80%     
==========================================
  Files          39       40       +1     
  Lines        2135     2160      +25     
==========================================
+ Hits         1974     1980       +6     
- Misses        161      180      +19     

Flag	Coverage Δ
core	90.37% <35.08%> (-1.17%)	⬇️
komamri	93.84% <ø> (ø)
plots	95.00% <ø> (ø)

Flags with carried forward coverage won't be shown. Click here to find out more.

Files Changed	Coverage Δ
KomaMRICore/src/simulation/BlochKA/KAtest.jl	0.00% <0.00%> (ø)
KomaMRICore/src/simulation/GPUFunctions.jl	38.63% <41.02%> (-15.65%)	⬇️
KomaMRICore/src/KomaMRICore.jl	100.00% <100.00%> (ø)
KomaMRICore/src/simulation/SimulatorCore.jl	87.09% <100.00%> (+6.04%)	⬆️

... and 1 file with indirect coverage changes

[cncastillo]

I think the way that get_flip_angles calculates the flip-angle is inaccurate. This comes from the fact that the RFs are interpreted as samples with zero-order hold. Thus the integral is calculated as $\alpha \propto \sum_i A_i \Delta_i$ instead of $\alpha \propto \sum_i \frac{A_i + A_{i+1}}{2} \Delta_i$ (trapezoidal rule). This generates errors in the PulseDesigner. RF_sinc test as $\alpha \neq 90$.

We could take this refactoring as an opportunity to remove the difference between RFs and GRs. This implies various things:

• Increased accuracy in RF blocks, needing fewer samples (equivalently, we could increase the default Δt_rf).
• The removal of the difference between $RF(A_N, T_N)$ and and $Grad(A_N, T_{N-1})$
• Removal of the stupid function ⏢ to do zero-order hold in the RFs. That does not even work for non-uniform RFs (src/datatypes/Sequence.jl).
• Plotting fewer points for RFs, affecting get_theo_A(r::RF) and get_theo_A(r::RF) (src/simulation/KeyValuesCalculation.jl).

Maybe it would be better to do:

perform_step!(p::Phantom{T}, seq::DiscreteSequence{T}, sig::AbstractArray{Complex{T}}, M::Mag{T}, sim_method::Bloch) where {T<:Real}
        #Motion
        xt = p.x .+ p.ux(p.x, p.y, p.z, s.t)
        yt = p.y .+ p.uy(p.x, p.y, p.z, s.t)
        zt = p.z .+ p.uz(p.x, p.y, p.z, s.t)
        #Effective field
        Bxy = s.B1  # <----  TODO: $B{xy} = (B_{xy}(t_i) + B_{xy}(t_{i+1}) )/2$
        ΔBz = p.Δw ./ T(2π * γ) .- s.Δf ./ T(γ) # Off-resonance ΔB_0 = (B_0 - ω_rf/γ)
        Bz = (s.Gx .* xt .+ s.Gy .* yt .+ s.Gz .* zt) .+ ΔBz # <---- TODO: $B{z} = (B_{z}(t_i) + B_{z}(t_{i+1}) )/2$
        #Spinor Rotation
        B = sqrt.(abs.(Bxy) .^ 2 .+ abs.(Bz) .^ 2)
        φ = T(-2π * γ) * (B .* s.Δt) 
        mul!( Q(φ, s.B1 ./ B, Bz ./ B), M )
        #Relaxation
        M.xy .= M.xy .* exp.(-s.Δt ./ p.T2)
        M.z  .= M.z  .* exp.(-s.Δt ./ p.T1) .+ p.ρ .* (1 .- exp.(-s.Δt ./ p.T1))
        #Saving results
        if s.ADC || sim_method.save_everystep
                signal .= sum(M.xy)
        end
        return nothing
end

[cncastillo]

This was done in #351.