Fixed mistake in prefix sum

Author: artemgl

src/App/Program.fs

@@ -10,56 +10,69 @@ open GraphBLAS.FSharp.Predefined
 
 [<EntryPoint>]
 let main argv =
-    let fstMatrix = COOMatrix(100, 100, [|0;1;1;2;5|], [|4;1;50;2;5|], [|1.0;2.0;76.0;3.0;6.0|])
-    let sndMatrix = COOMatrix(100, 100, [|0;1;1;2;5|], [|4;2;50;3;5|], [|-0.8;2.0;-76.0;3.0;6.0|])
-    let workflow =
-        opencl {
-            let! result = fstMatrix.EWiseAdd sndMatrix None FloatSemiring.addMult
-            let! _ = result.ToHost ()
-            return result
-        }
-    let res: COOMatrix<float> = downcast oclContext.RunSync workflow
+    // // let fstMatrix = COOMatrix(1, 2, [|0|], [|1|], [|true|])
+    // // let sndMatrix = COOMatrix(1, 2, [|0|], [|0|], [|true|])
+    // let fstMatrix = COOMatrix(100, 100, [|0;1;1;2;5|], [|4;1;50;2;5|], [|1.0;2.0;76.0;3.0;6.0|])
+    // let sndMatrix = COOMatrix(100, 100, [|0;1;1;2;5|], [|4;2;50;3;5|], [|-0.8;2.0;-76.0;3.0;6.0|])
+    // let workflow =
+    //     opencl {
+    //         let! result = fstMatrix.EWiseAdd sndMatrix None FloatSemiring.addMult
+    //         let! _ = result.ToHost ()
+    //         return result
+    //     }
+    // let res: COOMatrix<float> = downcast oclContext.RunSync workflow
 
-    let indices = res.Indices
-    // let rows = res.Rows
-    // let columns = res.Columns
-    let values = res.Values
+    // let indices = res.Indices
+    // // let rows = res.Rows
+    // // let columns = res.Columns
+    // let values = res.Values
 
-    // for i in 0 .. rows.Length - 1 do
-    //     printfn "(%i, %i, %A)" rows.[i] columns.[i] values.[i]
+    // // for i in 0 .. rows.Length - 1 do
+    // //     printfn "(%i, %i, %A)" rows.[i] columns.[i] values.[i]
 
-    for i in 0 .. indices.Length - 1 do
-        let index = indices.[i]
-        let a = int (index >>> 32)
-        let b = int index
-        printfn "(%i, %i, %A)" a b values.[i]
+    // for i in 0 .. indices.Length - 1 do
+    //     let index = indices.[i]
+    //     let a = int (index >>> 32)
+    //     let b = int index
+    //     printfn "(%i, %i, %A)" a b values.[i]
 
     // let allIndicesLength = 42
-    // let workGroupSize = 256
-    // let longSide = 72
-    // let shortSide = 27
+    let workGroupSize = 256
+    let longSide = 72
+    let shortSide = 27
+    let defaultIdx = Unchecked.defaultof<int>
 
-    // let command =
-    //     <@
-    //         fun (ndRange: _1D)
-    //             (arr: int[]) ->
+    // let arrLength = 10
 
-    //             arr.[ndRange.GlobalID0] <-
-    //                 if true then 41 else 42
-    //     @>
+    let command =
+        <@
+            fun (ndRange: _1D)
+                (buf: int[]) ->
 
-    // let translate2opencl (provider: ComputeProvider) (command: Quotations.Expr<(_1D -> int[] -> unit)>) : string =
-    //     let options = ComputeProvider.DefaultOptions_p
-    //     let tOptions = []
-    //     provider.SetCompileOptions options
+                let i = ndRange.GlobalID0
+                let f a b =
+                    if a > b
+                    then
+                        let x = a * a
+                        if x > 42
+                        then x + b
+                        else a + x
+                    else a - b
+                buf.[i] <- f 2 3
+        @>
 
-    //     let kernel = System.Activator.CreateInstance<Kernel<_1D>>()
-    //     CLCodeGenerator.GenerateKernel(command, provider, kernel, tOptions) |> ignore
-    //     let code = (kernel :> ICLKernel).Source.ToString()
-    //     code
+    let translate2opencl (provider: ComputeProvider) (command: Quotations.Expr<(_1D -> int[] -> unit)>) : string =
+        let options = ComputeProvider.DefaultOptions_p
+        let tOptions = []
+        provider.SetCompileOptions options
 
-    // let code = translate2opencl oclContext.Provider command
+        let kernel = System.Activator.CreateInstance<Kernel<_1D>>()
+        CLCodeGenerator.GenerateKernel(command, provider, kernel, tOptions) |> ignore
+        let code = (kernel :> ICLKernel).Source.ToString()
+        code
 
-    // printf "%s" code
+    let code = translate2opencl oclContext.Provider command
+
+    printf "%s" code
 
     0

src/GraphBLAS-sharp/Implementations.fs

@@ -143,6 +143,14 @@ and COOMatrix<'a when 'a : struct and 'a : equality>(rowCount: int, columnCount:
                 yield (uint64 rows.[i]) <<< 32 ||| (uint64 columns.[i]) |]
         COOMatrix(rowCount, columnCount, indices, values)
 
+    override this.ToString() =
+        [
+            sprintf "COO Matrix %ix%i \n" rowCount columnCount
+            sprintf "Indices: %A \n" indices
+            sprintf "Values: %A \n" values
+        ]
+        |> String.concat ""
+
     override this.Clear () = failwith "Not Implemented"
     override this.Copy () = failwith "Not Implemented"
     override this.Resize a b = failwith "Not Implemented"
@@ -246,7 +254,7 @@ and COOMatrix<'a when 'a : struct and 'a : equality>(rowCount: int, columnCount:
 
         let allIndicesLength = allIndices.Length
 
-        let createSortedConcatenation =
+        let merge =
             <@
                 fun (ndRange: _1D)
                     (firstIndicesBuffer: uint64[])
@@ -259,7 +267,8 @@ and COOMatrix<'a when 'a : struct and 'a : equality>(rowCount: int, columnCount:
                     let i = ndRange.GlobalID0
 
                     if i < allIndicesLength then
-                        let knots = localArray<int> 2
+                        let mutable beginIdxLocal = local ()
+                        let mutable endIdxLocal = local ()
                         let localID = ndRange.LocalID0
                         if localID < 2 then
                             let mutable x = localID * (workGroupSize - 1) + i - 1
@@ -278,11 +287,12 @@ and COOMatrix<'a when 'a : struct and 'a : equality>(rowCount: int, columnCount:
                                 let secondIndex = secondIndicesBuffer.[diagonalNumber - middleIdx]
                                 if firstIndex < secondIndex then leftEdge <- middleIdx + 1 else rightEdge <- middleIdx - 1
 
-                            knots.[localID] <- leftEdge
+                            // Here localID equals either 0 or 1
+                            if localID = 0 then beginIdxLocal <- leftEdge else endIdxLocal <- leftEdge
                         barrier ()
 
-                        let beginIdx = knots.[0] // BANK CONFLICTS?
-                        let endIdx = knots.[1]
+                        let beginIdx = beginIdxLocal
+                        let endIdx = endIdxLocal
                         let firstLocalLength = endIdx - beginIdx
                         let mutable x = workGroupSize - firstLocalLength
                         if endIdx = longSide then x <- shortSide - i + localID + beginIdx
@@ -331,7 +341,7 @@ and COOMatrix<'a when 'a : struct and 'a : equality>(rowCount: int, columnCount:
                             allValuesBuffer.[i] <- firstValuesBuffer.[beginIdx + boundaryX]
             @>
 
-        let createSortedConcatenation =
+        let merge =
             opencl {
                 let binder kernelP =
                     let ndRange = _1D(workSize allIndices.Length, workGroupSize)
@@ -343,7 +353,7 @@ and COOMatrix<'a when 'a : struct and 'a : equality>(rowCount: int, columnCount:
                         secondValues
                         allIndices
                         allValues
-                do! RunCommand createSortedConcatenation binder
+                do! RunCommand merge binder
             }
 
         let auxiliaryArray = Array.create allIndices.Length 1
@@ -400,29 +410,25 @@ and COOMatrix<'a when 'a : struct and 'a : equality>(rowCount: int, columnCount:
                         resultValuesBuffer.[index] <- allValuesBuffer.[i]
             @>
 
-        let resultIndices = Array.zeroCreate allIndices.Length
-        let resultValues = Array.create allValues.Length zero
-
-        let createUnion =
-            opencl {
-                let! prefixSumArray = prefixSum3 auxiliaryArray
-                let binder kernelP =
-                    let ndRange = _1D(workSize auxiliaryArray.Length, workGroupSize)
-                    kernelP
-                        ndRange
-                        allIndices
-                        allValues
-                        prefixSumArray
-                        resultIndices
-                        resultValues
-                do! RunCommand createUnion binder
-            }
-
         opencl {
-            do! createSortedConcatenation
+            do! merge
             do! filterThroughMask
             do! fillAuxiliaryArray
-            do! createUnion
+            let! prefixSumArray, resultLength = prefixSum auxiliaryArray
+            let! _ = ToHost resultLength
+            let resultLength = resultLength.[0]
+            let resultIndices = Array.zeroCreate resultLength
+            let resultValues = Array.create resultLength zero
+            let binder kernelP =
+                let ndRange = _1D(workSize prefixSumArray.Length, workGroupSize)
+                kernelP
+                    ndRange
+                    allIndices
+                    allValues
+                    prefixSumArray
+                    resultIndices
+                    resultValues
+            do! RunCommand createUnion binder
 
             return upcast COOMatrix<'a>(this.RowCount, this.ColumnCount, resultIndices, resultValues)
         }
@@ -638,7 +644,7 @@ and SparseVector<'a when 'a : struct and 'a : equality>(size: int, indices: int[
 
         let createUnion =
             opencl {
-                let! prefixSumArray = Toolbox.prefixSum2 auxiliaryArray
+                let! prefixSumArray, _ = prefixSum auxiliaryArray
                 let binder kernelP =
                     let ndRange = _1D(workSize auxiliaryArray.Length, workGroupSize)
                     kernelP