|
| 1 | + |
| 2 | +using LoopVectorization, Test |
| 3 | + |
| 4 | +function reference_mul4!(target_arr, src, range_a, range_b, padded_axis_a, padded_axis_b) |
| 5 | + @inbounds @fastmath for a1i ∈ eachindex(range_a), |
| 6 | + a2i ∈ eachindex(range_a), |
| 7 | + b1i ∈ eachindex(range_b), |
| 8 | + b2i ∈ eachindex(range_b) |
| 9 | + a1 = range_a[a1i]; |
| 10 | + a2 = range_a[a2i] |
| 11 | + b1 = range_b[b1i]; b2 = range_b[b2i] |
| 12 | + contribution = zero(eltype(target_arr)) |
| 13 | + for i_a ∈ padded_axis_a, i_b ∈ padded_axis_b |
| 14 | + contribution += src[i_a, i_b] * src[i_a+a1,i_b+b1] * src[i_a+a2,i_b+b2] |
| 15 | + end |
| 16 | + target_arr[b1i, b2i] += contribution |
| 17 | + end |
| 18 | +end |
| 19 | +function mul4_turbo_v1!(target_arr, src, range_a, range_b, padded_axis_a, padded_axis_b) |
| 20 | + @turbo for a1i ∈ eachindex(range_a), |
| 21 | + a2i ∈ eachindex(range_a), |
| 22 | + b1i ∈ eachindex(range_b), |
| 23 | + b2i ∈ eachindex(range_b) |
| 24 | + a1 = range_a[a1i]; |
| 25 | + a2 = range_a[a2i] |
| 26 | + b1 = range_b[b1i]; b2 = range_b[b2i] |
| 27 | + contribution = zero(eltype(target_arr)) |
| 28 | + for i_a ∈ padded_axis_a, i_b ∈ padded_axis_b |
| 29 | + contribution += src[i_a, i_b] * src[i_a+a1,i_b+b1] * src[i_a+a2,i_b+b2] |
| 30 | + end |
| 31 | + target_arr[b1i, b2i] += contribution |
| 32 | + end |
| 33 | + target_arr |
| 34 | +end |
| 35 | +function mul4_turbo_v2!(target_arr, src, range_a, range_b, padded_axis_a, padded_axis_b) |
| 36 | + @turbo for a1i ∈ eachindex(range_a), |
| 37 | + a2i ∈ eachindex(range_a), |
| 38 | + b1i ∈ eachindex(range_b), |
| 39 | + b2i ∈ eachindex(range_b) |
| 40 | + a1 = range_a[a1i]; |
| 41 | + a2 = range_a[a2i] |
| 42 | + b1 = range_b[b1i]; b2 = range_b[b2i] |
| 43 | + contribution = target_arr[b1i, b2i] |
| 44 | + for i_a ∈ padded_axis_a, i_b ∈ padded_axis_b |
| 45 | + contribution += src[i_a, i_b] * src[i_a+a1,i_b+b1] * src[i_a+a2,i_b+b2] |
| 46 | + end |
| 47 | + target_arr[b1i, b2i] = contribution |
| 48 | + end |
| 49 | +end |
| 50 | +function mul4_turbo_v3!(target_arr, src, range_a, range_b, padded_axis_a, padded_axis_b) |
| 51 | + @turbo for a1i ∈ eachindex(range_a), |
| 52 | + a2i ∈ eachindex(range_a), |
| 53 | + b1i ∈ eachindex(range_b), |
| 54 | + b2i ∈ eachindex(range_b), |
| 55 | + i_a ∈ padded_axis_a, |
| 56 | + i_b ∈ padded_axis_b |
| 57 | + a1 = range_a[a1i]; |
| 58 | + a2 = range_a[a2i] |
| 59 | + b1 = range_b[b1i]; b2 = range_b[b2i] |
| 60 | + target_arr[b1i, b2i] += src[i_a, i_b] * src[i_a+a1,i_b+b1] * src[i_a+a2,i_b+b2] |
| 61 | + end |
| 62 | +end |
| 63 | +function mul4_turbo_v4!(target_arr, src, range_a, range_b, padded_axis_a, padded_axis_b) |
| 64 | + @turbo for b1i ∈ eachindex(range_b), b2i ∈ eachindex(range_b) |
| 65 | + b1 = range_b[b1i]; b2 = range_b[b2i] |
| 66 | + contribution = zero(eltype(target_arr)) |
| 67 | + for i_a ∈ padded_axis_a, i_b ∈ padded_axis_b, a1i ∈ eachindex(range_a), a2i ∈ eachindex(range_a) |
| 68 | + a1 = range_a[a1i]; a2 = range_a[a2i] |
| 69 | + contribution += src[i_a, i_b] * src[i_a+a1,i_b+b1] * src[i_a+a2,i_b+b2] |
| 70 | + end |
| 71 | + target_arr[b1i, b2i] += contribution |
| 72 | + end |
| 73 | +end |
| 74 | + |
| 75 | +@testset "Inner reductions" begin |
| 76 | + src = ones(19, 101) |
| 77 | + max_a = 7; max_b = 9 |
| 78 | + range_a = -max_a:max_a |
| 79 | + range_b = -max_b:max_b |
| 80 | + |
| 81 | + target_arr = zeros(Float64, length(range_b), length(range_b)) |
| 82 | + |
| 83 | + padded_axis_b = (first(axes(src,2)) .+ max_b):(last(axes(src,2)) - max_b) |
| 84 | + padded_axis_a = (first(axes(src,1)) .+ max_a):(last(axes(src,1)) - max_a) |
| 85 | + |
| 86 | + target_ref = zero(target_arr); |
| 87 | + reference_mul4!(target_ref, src, range_a, range_b, padded_axis_a, padded_axis_b) |
| 88 | + |
| 89 | + mul4_turbo_v1!(target_arr, src, range_a, range_b, padded_axis_a, padded_axis_b) |
| 90 | + @test target_arr ≈ target_ref |
| 91 | + target_arr .= 0; |
| 92 | + mul4_turbo_v2!(target_arr, src, range_a, range_b, padded_axis_a, padded_axis_b) |
| 93 | + @test target_arr ≈ target_ref |
| 94 | + target_arr .= 0; |
| 95 | + mul4_turbo_v3!(target_arr, src, range_a, range_b, padded_axis_a, padded_axis_b) |
| 96 | + @test target_arr ≈ target_ref |
| 97 | + target_arr .= 0; |
| 98 | + mul4_turbo_v4!(target_arr, src, range_a, range_b, padded_axis_a, padded_axis_b) |
| 99 | + @test target_arr ≈ target_ref |
| 100 | + |
| 101 | +end |
| 102 | + |
0 commit comments