Fix Tullio issue #131.

Author: chriselrod

Project.toml

@@ -1,7 +1,7 @@
 name = "LoopVectorization"
 uuid = "bdcacae8-1622-11e9-2a5c-532679323890"
 authors = ["Chris Elrod <[email protected]>"]
-version = "0.12.97"
+version = "0.12.98"
 
 [deps]
 ArrayInterface = "4fba245c-0d91-5ea0-9b3e-6abc04ee57a9"

src/codegen/lower_load.jl

@@ -456,8 +456,8 @@ function lower_load_collection!(
   offset_dummy_loop = Loop(first(opindices), MaybeKnown(1), MaybeKnown(1024), MaybeKnown(1), Symbol(""), Symbol(""))
   unrollcurl₂ = unrolled_curly(op, nouter, offset_dummy_loop, vloop, mask, 1) # interleave always 1 here
   inds = mem_offset_u(op, ua, inds_calc_by_ptr_offset, false, 0, ls, false)
+  # @show op suffix, inds
   falseexpr = Expr(:call, lv(:False)); rs = staticexpr(reg_size(ls));
-
   opu₁, opu₂ = isunrolled_sym(op, u₁loopsym, u₂loopsym, vloopsym, ls)
   manualunrollu₁ = if opu₁ && u₁ > 1 # both unrolled
     if isknown(step(u₁loop)) && sum(Base.Fix2(===,u₁loopsym), getindicesonly(op)) == 1

src/codegen/lower_memory_common.jl

@@ -307,12 +307,14 @@ function mem_offset_u(
       ind_by_offset = inds_calc_by_ptr_offset[n] | (ind === CONSTANTZEROINDEX)
       offset = convert(Int, offsets[n])
       stride = convert(Int, strides[n])
+      # @show stride, ind, u₂loopsym
       indvectorized = ind === vloopsym
       indvectorizedmm = _mm & indvectorized
       if ind === u₁loopsym
+        # TODO: should it be stride*incr₁ ???
         addvectoroffset!(ret, indvectorizedmm, incr₁, u₁step, vstep, stride, ind, offset, ind_by_offset, indvectorized) # 9 arg
       elseif ind === u₂loopsym
-        addvectoroffset!(ret, indvectorizedmm, incr₂, u₂step, vstep, stride, ind, offset, ind_by_offset, indvectorized) # 9 arg
+        addvectoroffset!(ret, indvectorizedmm, stride * incr₂, u₂step, vstep, stride, ind, offset, ind_by_offset, indvectorized) # 9 arg
       elseif loopedindex[n]
         addoffset!(ret, indvectorizedmm, vstep, stride, ind, offset, ind_by_offset) # 7 arg
       else

src/modeling/determinestrategy.jl

@@ -771,40 +771,45 @@ end
 #     mno, id
 # end
 function maxnegativeoffset(ls::LoopSet, op::Operation, u::Symbol)
-    mno::Int = typemin(Int)
-    id = 0
-    isknown(step(getloop(ls, u))) || return mno, id
-    omop = offsetloadcollection(ls)
-    collectionid, opind = omop.opidcollectionmap[identifier(op)]
-    collectionid == 0 && return mno, id
-    @unpack opids = omop
+  mno::Int = typemin(Int)
+  id = 0
+  isknown(step(getloop(ls, u))) || return mno, id
+  omop = offsetloadcollection(ls)
+  collectionid, opind = omop.opidcollectionmap[identifier(op)]
+  collectionid == 0 && return mno, id
+  @unpack opids = omop
 
-    # offsetcol = offsets[collectionid]
-    opidcol = opids[collectionid]
-    opid = identifier(op)
-    # opoffs = offsetcol[opind]
-    opoffs = getoffsets(op)
-    ops = operations(ls)
-    opindices = getindicesonly(op)
-    for (i,oppid) ∈ enumerate(opidcol)
-        opid == oppid && continue
-        opp = ops[oppid]
-        oppoffs = getoffsets(opp)
-        mnonew::Int = typemin(Int)
-        for i ∈ eachindex(opindices)
-            if opindices[i] === u
-                mnonew = ((opoffs[i] % Int) - (oppoffs[i] % Int))
-            elseif opoffs[i] != oppoffs[i]
-                mnonew = 1
-                break
-            end
-        end
-        if mno < mnonew < 0
-            mno = mnonew
-            id = identifier(opp)
-        end
+  opidcol = opids[collectionid]
+  opid = identifier(op)
+  opoffs = getoffsets(op)
+  opstrd = getstrides(op)
+  ops = operations(ls)
+  opindices = getindicesonly(op)
+  for oppid ∈ opidcol
+    opid == oppid && continue
+    opp = ops[oppid]
+    oppoffs = getoffsets(opp)
+    oppstrd = getstrides(opp)
+    mnonew::Int = typemin(Int)
+    for i ∈ eachindex(opindices)
+      strd = opstrd[i]
+      strd == oppstrd[i] == 1 || continue
+      if opindices[i] === u
+        mnonew = (opoffs[i] % Int) - (oppoffs[i] % Int)
+        # mnonew_t, mnonew_rem = divrem((opoffs[i] % Int) - (oppoffs[i] % Int), strd % Int)
+        # mnonew_rem == 0 || continue
+        # mnonew = mnonew_t
+      elseif opoffs[i] != oppoffs[i]
+        mnonew = 1
+        break
+      end
     end
-    mno, id
+    if mno < mnonew < 0
+      mno = mnonew
+      id = identifier(opp)
+    end
+  end
+  mno, id
 end
 function maxnegativeoffset(ls::LoopSet, op::Operation, unrollsyms::UnrollSymbols)
     @unpack u₁loopsym, u₂loopsym, vloopsym = unrollsyms

test/shuffleloadstores.jl

@@ -401,6 +401,31 @@ function reverse_part_ref(n1,n2)
   return A
 end
 
+
+function tullio_issue_131_ref(arr)
+  M, N = size(arr)
+  out = zeros(M >>> 1, N >>> 1)
+  @inbounds @fastmath for j in axes(out,2)
+    for i in axes(out,1)
+      out[i, j] = arr[2i, 2j] +  arr[2i - 1, 2j] + arr[2i - 1, 2j - 1] + arr[2i, 2j - 1]
+    end
+  end
+  out
+end
+
+
+function tullio_issue_131(arr)
+  M, N = size(arr)
+  out = zeros(M >>> 1, N >>> 1)
+  @turbo for j in axes(out,2)
+    for i in axes(out,1)
+      out[i, j] = arr[2i, 2j] +  arr[2i - 1, 2j] + arr[2i - 1, 2j - 1] + arr[2i, 2j - 1]
+    end
+  end
+  out
+end
+
+
 @testset "shuffles load/stores" begin
   @show @__LINE__
   for i ∈ 1:128
@@ -425,8 +450,10 @@ end
     end
     @test qsimd ≈ Base.vect(qdot_affine(xqv, yqv)...) ≈ Base.vect(qdot_stride(xqv, yqv)...)
 
-    if VERSION ≥ v"1.6.0-rc1"
-      for j ∈ max(1,i-5):i+5, k ∈ max(1,i-5,i+5)
+    for j ∈ max(1,i-5):i+5, k ∈ max(1,i-5,i+5)
+      A = rand(j+1, k);
+      @test tullio_issue_131(A) ≈ tullio_issue_131_ref(A)
+      if VERSION ≥ v"1.6.0-rc1"
         Ac = rand(Complex{Float64}, j, i);
         Bc = rand(Complex{Float64}, i, k);
         Cc1 = Ac*Bc;