dph-examples: Add handvectorised QuickHull.
authorGeorge Roldugin <groldugin@cse.unsw.edu.au>
Sun, 30 Sep 2012 08:32:41 +0000 (18:32 +1000)
committerGeorge Roldugin <groldugin@cse.unsw.edu.au>
Fri, 22 Feb 2013 05:12:13 +0000 (16:12 +1100)
dph-examples/dph-examples.cabal
dph-examples/examples/spectral/QuickHull/dph/Vectorised.copy.hcr
dph-examples/examples/spectral/QuickHull/handvec/Handvec.hs [new file with mode: 0644]
dph-examples/examples/spectral/QuickHull/handvec/HandvecWrp.hs [new file with mode: 0644]
dph-examples/examples/spectral/QuickHull/handvec/Main.hs [new file with mode: 0644]

index b05ea24..3fdac1e 100644 (file)
@@ -153,6 +153,16 @@ Executable dph-spectral-quickhull-vector
   ghc-options: -eventlog -rtsopts -threaded -fllvm -optlo-O3 -Odph -package dph-lifted-vseg -fcpr-off -fsimpl-tick-factor=1000
 
 
+Executable dph-spectral-quickhull-handvec
+  Build-depends:  base == 4.6.*, vector == 0.9.*, random == 1.0.*, old-time == 1.0.*, containers == 0.5.*, dph-base == 0.6.*, dph-prim-par == 0.6.*, dph-lifted-copy == 0.6.*, HUnit == 1.2.*
+  Main-is:        Main.hs
+  other-modules:  Handvec
+                  HandvecWrp
+                  Timing Points2D.Types SVG
+  hs-source-dirs: examples/spectral/QuickHull/handvec examples/spectral/QuickHull/lib lib
+  ghc-options:    -eventlog -rtsopts -threaded -fllvm -Odph -package dph-lifted-copy -fcpr-off -fsimpl-tick-factor=1000
+
+
 Executable dph-spectral-quicksort
   Build-depends: base == 4.6.*, vector == 0.9.*, random == 1.0.*, old-time == 1.1.*, containers == 0.5.*, dph-base == 0.8.*, dph-prim-par == 0.8.*, dph-lifted-vseg == 0.8.*, HUnit == 1.2.*, repa-flow == 3.2.*
   Main-is: Main.hs
index 7c8dbe2..8915267 100644 (file)
@@ -19,11 +19,11 @@ a1_r615 = PNested segd_r613 (a_r614 `cast` ...)
 lvl2_r616
 lvl2_r616 = $fDT(,) $fDTInt $fDTDouble
 
-lvl3_r617
-lvl3_r617 = sum_s $fNumInt $fEltInt
+sum_sIntInt
+sum_sIntInt = sum_s $fNumInt $fEltInt
 
-lvl4_r618
-lvl4_r618 = packByTag $fEltInt
+packByTagInt
+packByTagInt = packByTag $fEltInt
 
 packByTagDbl
 packByTagDbl = packByTag $fEltDouble
@@ -34,8 +34,8 @@ int1 = I# 1
 dbl0.0#
 dbl0.0# = D# 0.0
 
-a2_r61c
-a2_r61c =
+distAbsentErr
+distAbsentErr =
   absentError
     "ww_s5Pu{v} [lid] Dist{tc rE4}\n\
     \                   ((USegd{tc rLh}, Int{(w) tc 3J}), Int{(w) tc 3J})"#
@@ -154,11 +154,11 @@ $wlhsplit_r61h =
                 case elementsSel2_1 selGZ of _ { I# _ ->
                 case points_s_data `cast` ... of _ { P_2 xos yos ->
                 let {
-                  tagsGZ
-                  tagsGZ = tagsSel2 selGZ } in
+                  tagsZ
+                  tagsZ = tagsSel2 selGZ } in
                 (P_2
-                   ((packByTagDbl (xos `cast` ...) tagsGZ int1) `cast` ...)
-                   ((packByTagDbl (yos `cast` ...) tagsGZ int1) `cast` ...))
+                   ((packByTagDbl (xos `cast` ...) tagsZ int1) `cast` ...)
+                   ((packByTagDbl (yos `cast` ...) tagsZ int1) `cast` ...))
                 `cast` ...
                 }
                 }
@@ -183,240 +183,240 @@ $wlhsplit_r61h =
                 })
         of _ { UPSel2 usel2 upselrep2 ->
         case usel2
-        of ww2_X5PN { Usel22 usel2_tags usel2_idxs usel2_n0s usel2_n1s ->
+        of usel2' { Usel22 usel2_tags usel2_idxs usel2_n0s usel2_n1s ->
         let {
-          a5_X4Y6
-          a5_X4Y6 = UPSel2 ww2_X5PN upselrep2 } in
-        case elementsSel2_0 a5_X4Y6 of wild2_d3oN { I# n#_d3oP ->
+          upsel2_recr
+          upsel2_recr = UPSel2 usel2' upselrep2 } in -- upsel2 is recreated, check if it's optimised away
+        case elementsSel2_0 upsel2_recr of nGZ { I# n#nGZ ->
         case mkSegd
-               (replicate $fEltInt wild2_d3oN int1)
-               (enumFromStepLen int0 int1 wild2_d3oN)
-               wild2_d3oN
-        of _ { UPSegd ww7_s5Pr ww8_s5Ps ww9_s5Pt _ ->
+               (replicate $fEltInt nGZ int1)
+               (enumFromStepLen int0 int1 nGZ)
+               nGZ
+        of _ { UPSegd segdGZ_lens segdGZ_idxs segdGZ_nelts _ ->
         let {
-          a7_s56m
-          a7_s56m =
+          segdGZx1 -- one elt per segment in segd
+          segdGZx1 =
             mkSegd
-              (replicate $fEltInt wild2_d3oN int1)
-              (enumFromStepLen int0 int1 wild2_d3oN)
-              wild2_d3oN } in
+              (replicate $fEltInt nGZ int1)
+              (enumFromStepLen int0 int1 nGZ)
+              nGZ } in
         let {
-          a8_X57a
-          a8_X57a = UPSegd ww7_s5Pr ww8_s5Ps ww9_s5Pt a2_r61c } in
+          upsegdGZx1_recr -- recreated again
+          upsegdGZx1_recr = UPSegd segdGZ_lens segdGZ_idxs segdGZ_nelts distAbsentErr } in
         let {
-          segd1_a4Ml
-          segd1_a4Ml = plusSegd a8_X57a a7_s56m } in
-        case elementsSegd segd1_a4Ml of _ { I# k#nelts_segd1_a4Ml ->
+          segdGZx2
+          segdGZx2 = plusSegd upsegdGZx1_recr segdGZx1 } in
+        case elementsSegd segdGZx2 of _ { I# k#nelts_segdGZx2 ->
         let {
-          packed1_d3mB
-          packed1_d3mB =
+          packed_r
+          packed_r =
             let {
-              packed2_s3CI
-              packed2_s3CI =
-                case packed `cast` ... of _ { PNested segd2_a4Xt xs_a4Xu ->
+              packed2
+              packed2 =
+                case packed `cast` ... of _ { PNested packed_segd packed_s_data ->
                 let {
-                  tags_a4XC
-                  tags_a4XC = tagsSel2 a5_X4Y6 } in
+                  tagsZ
+                  tagsZ = tagsSel2 upsel2_recr } in
                 let {
-                  segd'_a4XG
-                  segd'_a4XG =
+                  segd'segdGZ'
+                  segd'segdGZ' =
                     lengthsToSegd
-                      (lvl4_r618 (lengthsSegd segd2_a4Xt) tags_a4XC int0) } in
+                      (packByTagInt (lengthsSegd packed_segd) tagsZ int0) } in
                 (PNested
-                   segd'_a4XG
-                   (case segd'_a4XG
+                   segd'segdGZ'
+                   (case segd'segdGZ'
                     of _ { UPSegd ipv_s5oB ipv1_s5oC ipv2_s5oD ipv3_s5oE ->
-                    case xs_a4Xu `cast` ... of _ { P_2 x1s y1s ->
+                    case packed_s_data `cast` ... of _ { P_2 x1s y1s ->
                     let {
-                      tags1_a4VN
-                      tags1_a4VN = replicate_s $fEltInt segd2_a4Xt tags_a4XC } in
+                      tagsZ_r
+                      tagsZ_r = replicate_s $fEltInt packed_segd tagsZ } in
                     (P_2
-                       ((packByTagDbl (x1s `cast` ...) tags1_a4VN int0) `cast` ...)
-                       ((packByTagDbl (y1s `cast` ...) tags1_a4VN int0) `cast` ...))
+                       ((packByTagDbl (x1s `cast` ...) tagsZ_r int0) `cast` ...)
+                       ((packByTagDbl (y1s `cast` ...) tagsZ_r int0) `cast` ...))
                     `cast` ...
                     }
                     }))
                 `cast` ...
                 } } in
-            case packed2_s3CI `cast` ... of _ { PNested segd2_a4Xt xs_a4Xu ->
+            case packed2 `cast` ... of _ { PNested packed_segd packed_s_data ->
             let {
-              ns_a58G
-              ns_a58G =
-                replicate_s $fEltInt segd1_a4Ml (lengthsSegd segd2_a4Xt) } in
+              doubleLens
+              doubleLens =
+                replicate_s $fEltInt segdGZx2 (lengthsSegd packed_segd) } in
             let {
-              xsegd'_a58K
-              xsegd'_a58K = lengthsToSegd ns_a58G } in
+              segd_r'
+              segd_r' = lengthsToSegd doubleLens } in
             (PNested
-               xsegd'_a58K
-               (case xsegd'_a58K
-                of segd3_a5nC { UPSegd ipv_s5oS ipv1_s5oT ipv2_s5oU ipv3_s5oV ->
-                case xs_a4Xu `cast` ... of _ { P_2 x1s y1s ->
+               segd_r'
+               (case segd_r'
+                of segd_r'' { UPSegd ipv_s5oS ipv1_s5oT ipv2_s5oU ipv3_s5oV ->
+                case packed_s_data `cast` ... of _ { P_2 x1s y1s ->
                 let {
-                  ixs_a54J
-                  ixs_a54J =
+                  ids_to_take
+                  ids_to_take =
                     enumFromStepLenEach
-                      (elementsSegd segd3_a5nC)
-                      (replicate_s $fEltInt segd1_a4Ml (indicesSegd segd2_a4Xt))
-                      (replicate $fEltInt (elementsSegd segd1_a4Ml) int1)
-                      ns_a58G } in
+                      (elementsSegd segd_r'')
+                      (replicate_s $fEltInt segdGZx2 (indicesSegd packed_segd))
+                      (replicate $fEltInt (elementsSegd segdGZx2) int1)
+                      doubleLens } in
                 (P_2
-                   ((bpermute $fEltDouble (x1s `cast` ...) ixs_a54J) `cast` ...)
-                   ((bpermute $fEltDouble (y1s `cast` ...) ixs_a54J) `cast` ...))
+                   ((bpermute $fEltDouble (x1s `cast` ...) ids_to_take) `cast` ...)
+                   ((bpermute $fEltDouble (y1s `cast` ...) ids_to_take) `cast` ...))
                 `cast` ...
                 }
                 }))
             `cast` ...
             } } in
         let {
-          packed2_a3iy
-          packed2_a3iy = packed1_d3mB } in
+          packed_r'
+          packed_r' = packed_r } in
         case ($wlhsplit_r61h
-                k#nelts_segd1_a4Ml
-                packed2_a3iy
+                k#nelts_segdGZx2
+                packed_r'
                 (let {
-                   pm_s3CH
-                   pm_s3CH =
-                     case points_s `cast` ... of _ { PNested segd2_a4O9 xs_a4Oa ->
-                     case xs_a4Oa `cast` ... of _ { P_2 x1s y1s ->
+                   pm_s
+                   pm_s =
+                     case points_s `cast` ... of _ { PNested segd_pts pts_s_dat ->
+                     case pts_s_dat `cast` ... of _ { P_2 x1s y1s ->
                      let {
-                       ixs_a54J
-                       ixs_a54J =
-                         case cross `cast` ... of _ { PNested segd3_a54a xs1_a54b ->
+                       ids_to_take
+                       ids_to_take =
+                         case cross `cast` ... of _ { PNested segd_crs crs_s_dat ->
                          zipWith
                            $fEltInt
                            $fEltInt
                            $fEltInt
                            $fNumInt_$c+
-                           (indicesSegd segd2_a4O9)
+                           (indicesSegd segd_pts)
                            (fsts
                               $fEltInt
                               $fEltDouble
                               (fold1_s
                                  lvl10_r61g
                                  ($vmaxIndexP_$smax')
-                                 segd3_a54a
+                                 segd_crs
                                  (zip
                                     $fEltInt
                                     $fEltDouble
-                                    (indices_s segd3_a54a)
-                                    (xs1_a54b `cast` ...))))
+                                    (indices_s segd_crs)
+                                    (crs_s_dat `cast` ...))))
                          } } in
                      let {
-                       a9_s5QL
-                       a9_s5QL = bpermute $fEltDouble (y1s `cast` ...) ixs_a54J } in
+                       maxys
+                       maxys = bpermute $fEltDouble (y1s `cast` ...) ids_to_take } in
                      let {
-                       a10_s5QH
-                       a10_s5QH = bpermute $fEltDouble (x1s `cast` ...) ixs_a54J } in
+                       maxxs
+                       maxxs = bpermute $fEltDouble (x1s `cast` ...) ids_to_take } in
                      let {
-                       a11_s5QN
-                       a11_s5QN = P_2 (a10_s5QH `cast` ...) (a9_s5QL `cast` ...) } in
+                       maxxys
+                       maxxys = P_2 (maxxs `cast` ...) (maxys `cast` ...) } in
                      let {
-                       pm1_s3CJ
-                       pm1_s3CJ = a11_s5QN `cast` ... } in
-                     case pm1_s3CJ `cast` ... of _ { P_2 as1_s5Ue bs1_s5Uf ->
+                       maxxys'
+                       maxxys' = maxxys `cast` ... } in
+                     case maxxys' `cast` ... of _ { P_2 maxxs maxys ->
                      let {
-                       tagsGZ
-                       tagsGZ = tagsSel2 a5_X4Y6 } in
+                       tagsZ
+                       tagsZ = tagsSel2 upsel2_recr } in
                      (P_2
-                        ((packByTagDbl (as1_s5Ue `cast` ...) tagsGZ int0) `cast` ...)
-                        ((packByTagDbl (bs1_s5Uf `cast` ...) tagsGZ int0) `cast` ...))
+                        ((packByTagDbl (maxxs `cast` ...) tagsZ int0) `cast` ...)
+                        ((packByTagDbl (maxys `cast` ...) tagsZ int0) `cast` ...))
                      `cast` ...
                      }
                      }
                      } } in
                  (P_2
                     (case line_starts `cast` ... of _ { P_2 x1s y1s ->
-                     case pm_s3CH `cast` ... of _ { P_2 as1_X6OI bs1_X6OL ->
+                     case pm_s `cast` ... of _ { P_2 pm_xs pm_ys ->
                      let {
-                       tagsGZ
-                       tagsGZ = tagsSel2 a5_X4Y6 } in
+                       tagsZ
+                       tagsZ = tagsSel2 upsel2_recr } in
                      (P_2
                         ((append_s
                             $fEltDouble
-                            segd1_a4Ml
-                            a8_X57a
-                            (packByTagDbl (x1s `cast` ...) tagsGZ int0)
-                            a7_s56m
-                            (as1_X6OI `cast` ...))
+                            segdGZx2
+                            upsegdGZx1_recr
+                            (packByTagDbl (x1s `cast` ...) tagsZ int0)
+                            segdGZx1
+                            (pm_xs `cast` ...))
                          `cast` ...)
                         ((append_s
                             $fEltDouble
-                            segd1_a4Ml
-                            a8_X57a
-                            (packByTagDbl (y1s `cast` ...) tagsGZ int0)
-                            a7_s56m
-                            (bs1_X6OL `cast` ...))
+                            segdGZx2
+                            upsegdGZx1_recr
+                            (packByTagDbl (y1s `cast` ...) tagsZ int0)
+                            segdGZx1
+                            (pm_ys `cast` ...))
                          `cast` ...))
                      `cast` ...
                      }
                      })
-                    (case pm_s3CH `cast` ... of _ { P_2 as_X5Rg bs_X5Ri ->
-                     case line_ends `cast` ... of _ { P_2 as1_a4TX bs1_a4TY ->
+                    (case pm_s `cast` ... of _ { P_2 pm_xs' pm_ys' ->
+                    case line_ends `cast` ... of _ { P_2 line_x2s line_y2s ->
                      let {
-                       tagsGZ
-                       tagsGZ = tagsSel2 a5_X4Y6 } in
+                       tagsZ
+                       tagsZ = tagsSel2 upsel2_recr } in
                      (P_2
                         ((append_s
                             $fEltDouble
-                            segd1_a4Ml
-                            a8_X57a
-                            (as_X5Rg `cast` ...)
-                            a7_s56m
-                            (packByTagDbl (as1_a4TX `cast` ...) tagsGZ int0))
+                            segdGZx2
+                            upsegdGZx1_recr
+                            (pm_xs' `cast` ...)
+                            segdGZx1
+                            (packByTagDbl (line_x2s `cast` ...) tagsZ int0))
                          `cast` ...)
                         ((append_s
                             $fEltDouble
-                            segd1_a4Ml
-                            a8_X57a
-                            (bs_X5Ri `cast` ...)
-                            a7_s56m
-                            (packByTagDbl (bs1_a4TY `cast` ...) tagsGZ int0))
+                            segdGZx2
+                            upsegdGZx1_recr
+                            (pm_ys' `cast` ...)
+                            segdGZx1
+                            (packByTagDbl (line_y2s `cast` ...) tagsZ int0))
                          `cast` ...))
                      `cast` ...
                      }
                      }))
                  `cast` ...))
              `cast` ...
-        of _ { PNested segd2_a4N4 xs_a4N5 ->
-        case elementsSel2_1 a5_X4Y6 of wild5_d3B8 { I# n#1_d3Ba ->
+        of _ { PNested points'_segd points'_data ->
+        case elementsSel2_1 upsel2_recr of nE { I# n#nE# ->
         let {
-          segd3_a5b7
-          segd3_a5b7 =
+          result_segd
+          result_segd =
             lengthsToSegd
               (combine2
                  $fEltInt
-                 (tagsSel2 a5_X4Y6)
-                 (repSel2 a5_X4Y6)
-                 (lvl3_r617 segd1_a4Ml (lengthsSegd segd2_a4N4))
-                 (replicate $fEltInt wild5_d3B8 int1)) } in
+                 (tagsSel2 upsel2_recr)
+                 (repSel2 upsel2_recr)
+                 (sum_sIntInt segdGZx2 (lengthsSegd points'_segd))
+                 (replicate $fEltInt nE int1)) } in
         (PNested
-           segd3_a5b7
-           (case segd3_a5b7
-            of segd4_a5nC { UPSegd ipv_s5ph ipv1_s5pi ipv2_s5pj ipv3_s5pk ->
-            case xs_a4N5 `cast` ... of _ { P_2 x1s y1s ->
-            case line_starts `cast` ... of _ { P_2 as1_X595 bs1_X597 ->
+           result_segd
+           (case result_segd
+            of result_segd' { UPSegd result_lens result_ids result_nelts ipv3_s5pk ->
+            case points'_data `cast` ... of _ { P_2 x1s y1s ->
+            case line_starts `cast` ... of _ { P_2 line_x1s line_x2s ->
             let {
-              tagsGZ
-              tagsGZ = tagsSel2 a5_X4Y6 } in
+              tagsZ
+              tagsZ = tagsSel2 upsel2_recr } in
             let {
-              sel_a5bV
-              sel_a5bV =
+              result_sel
+              result_sel =
                 tagsToSel2
-                  (replicate_s $fEltInt segd4_a5nC (tagsSel2 a5_X4Y6)) } in
+                  (replicate_s $fEltInt result_segd' (tagsSel2 upsel2_recr)) } in
             (P_2
                ((combine2
                    $fEltDouble
-                   (tagsSel2 sel_a5bV)
-                   (repSel2 sel_a5bV)
+                   (tagsSel2 result_sel)
+                   (repSel2 result_sel)
                    (x1s `cast` ...)
-                   (packByTagDbl (as1_X595 `cast` ...) tagsGZ int1))
+                   (packByTagDbl (line_x1s `cast` ...) tagsZ int1))
                 `cast` ...)
                ((combine2
                    $fEltDouble
-                   (tagsSel2 sel_a5bV)
-                   (repSel2 sel_a5bV)
+                   (tagsSel2 result_sel)
+                   (repSel2 result_sel)
                    (y1s `cast` ...)
-                   (packByTagDbl (bs1_X597 `cast` ...) tagsGZ int1))
+                   (packByTagDbl (line_x2s `cast` ...) tagsZ int1))
                 `cast` ...))
             `cast` ...
             }
@@ -636,7 +636,7 @@ $wlquickHull =
     let {
       a5_X5gr
       a5_X5gr = UPSel2 ww2_X5Qm ww1_s5PT } in
-    case elementsSel2_0 a5_X5gr of wild_d3oN { I# n#_d3oP ->
+    case elementsSel2_0 a5_X5gr of wild_d3oN { I# n#nGZ ->
     case mkSegd
            (replicate $fEltInt wild_d3oN int1)
            (enumFromStepLen int0 int1 wild_d3oN)
@@ -653,55 +653,55 @@ $wlquickHull =
       a8_X5iF
       a8_X5iF = UPSegd ww7_s5PW ww8_s5PX ww9_s5PY a3_r61x } in
     let {
-      segd1_a4Ml
-      segd1_a4Ml = plusSegd a8_X5iF a7_s5hP } in
-    case elementsSegd segd1_a4Ml of wild3_a2RP { I# k#nelts_segd1_a4Ml ->
+      segdGZx2
+      segdGZx2 = plusSegd a8_X5iF a7_s5hP } in
+    case elementsSegd segdGZx2 of wild3_a2RP { I# k#nelts_segdGZx2 ->
     case ($wlhsplit_r61h
-            k#nelts_segd1_a4Ml
-            (case w1_s5Q4 `cast` ... of _ { PNested segd2_a4Xt xs_a4Xu ->
+            k#nelts_segdGZx2
+            (case w1_s5Q4 `cast` ... of _ { PNested packed_segd packed_s_data ->
              let {
-               tags_a4XC
-               tags_a4XC = tagsSel2 a5_X5gr } in
+               tagsZ
+               tagsZ = tagsSel2 a5_X5gr } in
              let {
-               segd'_a4XG
-               segd'_a4XG =
+               segd'segdGZ'
+               segd'segdGZ' =
                  lengthsToSegd
-                   (lvl4_r618 (lengthsSegd segd2_a4Xt) tags_a4XC int0) } in
+                   (packByTagInt (lengthsSegd packed_segd) tagsZ int0) } in
              let {
-               ns_a58G
-               ns_a58G =
-                 replicate_s $fEltInt segd1_a4Ml (lengthsSegd segd'_a4XG) } in
+               doubleLens
+               doubleLens =
+                 replicate_s $fEltInt segdGZx2 (lengthsSegd segd'segdGZ') } in
              let {
-               xsegd'_a58K
-               xsegd'_a58K = lengthsToSegd ns_a58G } in
+               segd_r'
+               segd_r' = lengthsToSegd doubleLens } in
              (PNested
-                xsegd'_a58K
-                (case xsegd'_a58K
-                 of segd3_a5nC { UPSegd ipv_s5qk ipv1_s5ql ipv2_s5qm ipv3_s5qn ->
-                 case segd'_a4XG
+                segd_r'
+                (case segd_r'
+                 of segd_r'' { UPSegd ipv_s5qk ipv1_s5ql ipv2_s5qm ipv3_s5qn ->
+                 case segd'segdGZ'
                  of segd4_X5FV { UPSegd ipv4_s5q9 ipv5_s5qa ipv6_s5qb ipv7_s5qc ->
-                 case xs_a4Xu `cast` ... of _ { P_2 x1s y1s ->
+                 case packed_s_data `cast` ... of _ { P_2 x1s y1s ->
                  let {
-                   tags1_a4VN
-                   tags1_a4VN = replicate_s $fEltInt segd2_a4Xt tags_a4XC } in
+                   tagsZ_r
+                   tagsZ_r = replicate_s $fEltInt packed_segd tagsZ } in
                  let {
-                   ixs_a54J
-                   ixs_a54J =
+                   ids_to_take
+                   ids_to_take =
                      enumFromStepLenEach
-                       (elementsSegd segd3_a5nC)
-                       (replicate_s $fEltInt segd1_a4Ml (indicesSegd segd4_X5FV))
+                       (elementsSegd segd_r'')
+                       (replicate_s $fEltInt segdGZx2 (indicesSegd segd4_X5FV))
                        (replicate $fEltInt wild3_a2RP int1)
-                       ns_a58G } in
+                       doubleLens } in
                  (P_2
                     ((bpermute
                         $fEltDouble
-                        (packByTagDbl (x1s `cast` ...) tags1_a4VN int0)
-                        ixs_a54J)
+                        (packByTagDbl (x1s `cast` ...) tagsZ_r int0)
+                        ids_to_take)
                      `cast` ...)
                     ((bpermute
                         $fEltDouble
-                        (packByTagDbl (y1s `cast` ...) tags1_a4VN int0)
-                        ixs_a54J)
+                        (packByTagDbl (y1s `cast` ...) tagsZ_r int0)
+                        ids_to_take)
                      `cast` ...))
                  `cast` ...
                  }
@@ -723,37 +723,37 @@ $wlquickHull =
              let {
                minx_s3Co
                minx_s3Co =
-                 case w1_s5Q4 `cast` ... of _ { PNested segd2_a4O9 xs1_a4Oa ->
+                 case w1_s5Q4 `cast` ... of _ { PNested segd_pts xs1_a4Oa ->
                  case xs1_a4Oa `cast` ... of _ { P_2 x1s y1s ->
                  let {
-                   ixs_a54J
-                   ixs_a54J =
-                     case xs_s3Cr `cast` ... of _ { PNested segd3_a54a xs2_a54b ->
+                   ids_to_take
+                   ids_to_take =
+                     case xs_s3Cr `cast` ... of _ { PNested segd_crs xs2_a54b ->
                      zipWith
                        $fEltInt
                        $fEltInt
                        $fEltInt
                        $fNumInt_$c+
-                       (indicesSegd segd2_a4O9)
+                       (indicesSegd segd_pts)
                        (fsts
                           $fEltInt
                           $fEltDouble
                           (fold1_s
                              lvl24_r61w
                              ($vminIndexP_$smin')
-                             segd3_a54a
+                             segd_crs
                              (zip
                                 $fEltInt
                                 $fEltDouble
-                                (indices_s segd3_a54a)
+                                (indices_s segd_crs)
                                 (xs2_a54b `cast` ...))))
                      } } in
                  let {
                    a9_s5RY
-                   a9_s5RY = bpermute $fEltDouble (y1s `cast` ...) ixs_a54J } in
+                   a9_s5RY = bpermute $fEltDouble (y1s `cast` ...) ids_to_take } in
                  let {
                    a10_s5RU
-                   a10_s5RU = bpermute $fEltDouble (x1s `cast` ...) ixs_a54J } in
+                   a10_s5RU = bpermute $fEltDouble (x1s `cast` ...) ids_to_take } in
                  let {
                    a11_s5S0
                    a11_s5S0 = P_2 (a10_s5RU `cast` ...) (a9_s5RY `cast` ...) } in
@@ -762,11 +762,11 @@ $wlquickHull =
                    minx1_s3Cq = a11_s5S0 `cast` ... } in
                  case minx1_s3Cq `cast` ... of _ { P_2 as1_s5VM bs1_s5VN ->
                  let {
-                   tagsGZ
-                   tagsGZ = tagsSel2 a5_X5gr } in
+                   tagsZ
+                   tagsZ = tagsSel2 a5_X5gr } in
                  (P_2
-                    ((packByTagDbl (as1_s5VM `cast` ...) tagsGZ int0) `cast` ...)
-                    ((packByTagDbl (bs1_s5VN `cast` ...) tagsGZ int0) `cast` ...))
+                    ((packByTagDbl (as1_s5VM `cast` ...) tagsZ int0) `cast` ...)
+                    ((packByTagDbl (bs1_s5VN `cast` ...) tagsZ int0) `cast` ...))
                  `cast` ...
                  }
                  }
@@ -774,37 +774,37 @@ $wlquickHull =
              let {
                maxx_s3Cn
                maxx_s3Cn =
-                 case w1_s5Q4 `cast` ... of _ { PNested segd2_a4O9 xs1_a4Oa ->
+                 case w1_s5Q4 `cast` ... of _ { PNested segd_pts xs1_a4Oa ->
                  case xs1_a4Oa `cast` ... of _ { P_2 x1s y1s ->
                  let {
-                   ixs_a54J
-                   ixs_a54J =
-                     case xs_s3Cr `cast` ... of _ { PNested segd3_a54a xs2_a54b ->
+                   ids_to_take
+                   ids_to_take =
+                     case xs_s3Cr `cast` ... of _ { PNested segd_crs xs2_a54b ->
                      zipWith
                        $fEltInt
                        $fEltInt
                        $fEltInt
                        $fNumInt_$c+
-                       (indicesSegd segd2_a4O9)
+                       (indicesSegd segd_pts)
                        (fsts
                           $fEltInt
                           $fEltDouble
                           (fold1_s
                              lvl21_r61t
                              ($vmaxIndexP_$smax')
-                             segd3_a54a
+                             segd_crs
                              (zip
                                 $fEltInt
                                 $fEltDouble
-                                (indices_s segd3_a54a)
+                                (indices_s segd_crs)
                                 (xs2_a54b `cast` ...))))
                      } } in
                  let {
                    a9_s5RK
-                   a9_s5RK = bpermute $fEltDouble (y1s `cast` ...) ixs_a54J } in
+                   a9_s5RK = bpermute $fEltDouble (y1s `cast` ...) ids_to_take } in
                  let {
                    a10_s5RG
-                   a10_s5RG = bpermute $fEltDouble (x1s `cast` ...) ixs_a54J } in
+                   a10_s5RG = bpermute $fEltDouble (x1s `cast` ...) ids_to_take } in
                  let {
                    a11_s5RM
                    a11_s5RM = P_2 (a10_s5RG `cast` ...) (a9_s5RK `cast` ...) } in
@@ -813,11 +813,11 @@ $wlquickHull =
                    maxx1_s3Cp = a11_s5RM `cast` ... } in
                  case maxx1_s3Cp `cast` ... of _ { P_2 as1_s5Vw bs1_s5Vx ->
                  let {
-                   tagsGZ
-                   tagsGZ = tagsSel2 a5_X5gr } in
+                   tagsZ
+                   tagsZ = tagsSel2 a5_X5gr } in
                  (P_2
-                    ((packByTagDbl (as1_s5Vw `cast` ...) tagsGZ int0) `cast` ...)
-                    ((packByTagDbl (bs1_s5Vx `cast` ...) tagsGZ int0) `cast` ...))
+                    ((packByTagDbl (as1_s5Vw `cast` ...) tagsZ int0) `cast` ...)
+                    ((packByTagDbl (bs1_s5Vx `cast` ...) tagsZ int0) `cast` ...))
                  `cast` ...
                  }
                  }
@@ -828,7 +828,7 @@ $wlquickHull =
                  (P_2
                     ((append_s
                         $fEltDouble
-                        segd1_a4Ml
+                        segdGZx2
                         a8_X5iF
                         (as_X67i `cast` ...)
                         a7_s5hP
@@ -836,7 +836,7 @@ $wlquickHull =
                      `cast` ...)
                     ((append_s
                         $fEltDouble
-                        segd1_a4Ml
+                        segdGZx2
                         a8_X5iF
                         (bs_X67k `cast` ...)
                         a7_s5hP
@@ -850,7 +850,7 @@ $wlquickHull =
                  (P_2
                     ((append_s
                         $fEltDouble
-                        segd1_a4Ml
+                        segdGZx2
                         a8_X5iF
                         (as_X67i `cast` ...)
                         a7_s5hP
@@ -858,7 +858,7 @@ $wlquickHull =
                      `cast` ...)
                     ((append_s
                         $fEltDouble
-                        segd1_a4Ml
+                        segdGZx2
                         a8_X5iF
                         (bs_X67k `cast` ...)
                         a7_s5hP
@@ -869,56 +869,56 @@ $wlquickHull =
                  }))
              `cast` ...))
          `cast` ...
-    of _ { PNested segd2_a4N4 xs_a4N5 ->
+    of _ { PNested points'_segd points'_data ->
     case elementsSel2_1 a5_X5gr of _ { I# _ ->
     case w1_s5Q4 `cast` ... of _ { PNested segd3_a4Xt xs1_a4Xu ->
     let {
-      tags_a4XC
-      tags_a4XC = tagsSel2 a5_X5gr } in
+      tagsZ
+      tagsZ = tagsSel2 a5_X5gr } in
     let {
-      segd'_a4XG
-      segd'_a4XG =
+      segd'segdGZ'
+      segd'segdGZ' =
         lengthsToSegd
-          (lvl4_r618 (lengthsSegd segd3_a4Xt) tags_a4XC int1) } in
+          (packByTagInt (lengthsSegd segd3_a4Xt) tagsZ int1) } in
     let {
       segd4_a5b7
       segd4_a5b7 =
         lengthsToSegd
           (combine2
              $fEltInt
-             tags_a4XC
+             tagsZ
              (repSel2 a5_X5gr)
-             (lvl3_r617 segd1_a4Ml (lengthsSegd segd2_a4N4))
-             (lengthsSegd segd'_a4XG)) } in
+             (sum_sIntInt segdGZx2 (lengthsSegd points'_segd))
+             (lengthsSegd segd'segdGZ')) } in
     (PNested
        segd4_a5b7
        (case segd4_a5b7
         of segd5_a5nC { UPSegd ipv_s5qC ipv1_s5qD ipv2_s5qE ipv3_s5qF ->
-        case xs_a4N5 `cast` ... of _ { P_2 x1s y1s ->
-        case segd'_a4XG
+        case points'_data `cast` ... of _ { P_2 x1s y1s ->
+        case segd'segdGZ'
         of _ { UPSegd ipv4_s5qJ ipv5_s5qK ipv6_s5qL ipv7_s5qM ->
         case xs1_a4Xu `cast` ... of _ { P_2 as1_X5cv bs1_X5cx ->
         let {
-          tags1_a4VN
-          tags1_a4VN = replicate_s $fEltInt segd3_a4Xt tags_a4XC } in
+          tagsZ_r
+          tagsZ_r = replicate_s $fEltInt segd3_a4Xt tagsZ } in
         let {
-          sel_a5bV
-          sel_a5bV =
-            tagsToSel2 (replicate_s $fEltInt segd5_a5nC tags_a4XC) } in
+          result_sel
+          result_sel =
+            tagsToSel2 (replicate_s $fEltInt segd5_a5nC tagsZ) } in
         (P_2
            ((combine2
                $fEltDouble
-               (tagsSel2 sel_a5bV)
-               (repSel2 sel_a5bV)
+               (tagsSel2 result_sel)
+               (repSel2 result_sel)
                (x1s `cast` ...)
-               (packByTagDbl (as1_X5cv `cast` ...) tags1_a4VN int1))
+               (packByTagDbl (as1_X5cv `cast` ...) tagsZ_r int1))
             `cast` ...)
            ((combine2
                $fEltDouble
-               (tagsSel2 sel_a5bV)
-               (repSel2 sel_a5bV)
+               (tagsSel2 result_sel)
+               (repSel2 result_sel)
                (y1s `cast` ...)
-               (packByTagDbl (bs1_X5cx `cast` ...) tags1_a4VN int1))
+               (packByTagDbl (bs1_X5cx `cast` ...) tagsZ_r int1))
             `cast` ...))
         `cast` ...
         }
diff --git a/dph-examples/examples/spectral/QuickHull/handvec/Handvec.hs b/dph-examples/examples/spectral/QuickHull/handvec/Handvec.hs
new file mode 100644 (file)
index 0000000..6667e23
--- /dev/null
@@ -0,0 +1,245 @@
+{-# LANGUAGE MagicHash, BangPatterns, ParallelArrays, TypeOperators #-}
+
+module Handvec ( hsplit_v, Point, Line )
+where
+
+import Points2D.Types
+
+import Data.Array.Parallel as PA
+import Data.Array.Parallel.Base as B
+import Data.Array.Parallel.PArray ()
+import Data.Array.Parallel.PArray.PData
+import Data.Array.Parallel.PArray.PDataInstances
+import Data.Array.Parallel.PArray.ScalarInstances
+import Data.Array.Parallel.Lifted.Closure
+import Data.Array.Parallel.Prelude as P'
+import Data.Array.Parallel.Prelude.Double as D
+import Data.Array.Parallel.Unlifted as U
+import GHC.Exts
+
+import Debug.Trace as Debug
+
+import Prelude as P
+
+--type Tag = Bool -- from dph-prim-seq:D.A.P.Base
+
+-- not changed thanks to vectavoid (except tuples are dissolved)
+distance :: Double -> Double -- Point
+         -> Double -> Double -- Line start
+         -> Double -> Double -- Line end
+         -> Double          -- Distance
+distance xo yo x1 y1 x2 y2
+  = (x1 P.- xo) P.* (y2 P.- yo) P.- (y1 P.- yo) P.* (x2 P.- xo)
+{-# INLINE distance #-}
+
+hsplit_v :: PArray Point :-> Line :-> PArray Point
+hsplit_v = closure2 hsplit_s hsplit_l_wrp
+{-# NOINLINE hsplit_v #-}
+
+-- scalar version can be generalised to lifted hsplit:
+-- 1. construct a nested array out of $points$ array (single segment)
+-- 2. construct a singleton array containing $line$
+hsplit_s :: PArray Point -> Line -> PArray Point
+hsplit_s ps@(PArray n# pts) line
+  = let segd     = U.singletonSegd (I# n#)
+        pts_s    = PArray 1# (PNested segd pts)
+        lines    = PArray 1# (toPData line)
+        hull@(PArray _ (PNested segd' hull_pts))
+                 = hsplit_l_wrp pts_s lines
+        (I# n'#) = U.elementsSegd segd'
+    in Debug.trace "hsplit_s" $ PArray n'# hull_pts
+  where toPData ((x1,y1), (x2,y2)) = P_2 (P_2 (r x1) (r y1))
+                                         (P_2 (r x2) (r y2))
+        r x = PDouble (U.replicate 1 x)
+{-# NOINLINE hsplit_s #-}
+
+hsplit_l_wrp :: PArray (PArray Point) -> PArray Line -> PArray (PArray Point)
+hsplit_l_wrp (PArray n# pts) (PArray _ lns)
+  = let result@(PNested segd _) = hsplit_l_wrk n# pts lns
+        !(I# n'#) = U.elementsSegd segd
+    in  PArray n'# result
+{-# INLINE hsplit_l_wrp #-}
+
+hsplit_l_wrk :: Int# -> PData (PArray Point) -> PData Line -> PData (PArray Point)
+hsplit_l_wrk n# (PNested segd (P_2 (PDouble xos) _)) (P_2 (P_2 (PDouble x1s) _) _)
+  | Debug.trace ("hsplit_l_wrk: "++
+                 "nlines#: "     ++ show (I# n#) ++
+                 ", eltsSegd: "  ++ show (U.elementsSegd segd) ++
+                 ", lensSegd: "  ++ show (U.toList (U.lengthsSegd segd)) ++
+                 ", npts: "      ++ show (U.length xos) ++
+                 ", nlines arr: "++ show (U.length x1s)) False = undefined
+hsplit_l_wrk 0# _ _ = PNested U.emptySegd (P_2 (PDouble U.empty) (PDouble U.empty))
+hsplit_l_wrk nlines# points_s                                     -- nlines#: 6, npts >= 13
+         lines@(P_2 (P_2 (PDouble line_x1s) (PDouble line_y1s))
+                    (P_2 (PDouble line_x2s) (PDouble line_y2s)))
+  = let -- find distance-like measures between each point and its respective line
+        cross_s     = calc_cross_s  points_s lines
+        -- throw out all the points which are below the line (i.e. distance <= 0)
+        packed_s@(PNested segd (P_2 (PDouble xos) (PDouble yos)))
+                    = calc_packed_s points_s cross_s
+        -- for some of the lines there are no more points left (segment length 0)
+        -- the recursion for these lines has finished and we'll treat them separately
+        ( segd_lens, _, _) = getSegdParams segd             -- Segd ( [4, 5, 0, 2, 2,  0]
+                                                            --      , [0, 4, 9, 9, 11, 13]
+                                                            --      , 13 )
+        -- empty segments selector ((==0) -> 1; otherwise -> 0)
+        selE        = U.tagsToSel2                          -- Sel2 ( [0, 0, 1, 0, 0, 1]
+                    $ U.map B.fromBool                      --      , [0, 1, 0, 2, 3, 1]
+                    $ U.zipWith (P.==) segd_lens            --      , 4, 2 )
+                                    (U.replicate (I# nlines#) 0)
+        tagsE       = U.tagsSel2 selE                       -- [0, 0, 1, 0, 0, 1]
+        --- *** First take care of the lines for which the recursion hasn't finished *** ---
+        -- number of non-empty segments
+        !nNE@(I# nNE#) = U.elementsSel2_0 selE              -- 4
+        -- new segd with one element per non empty segments of original segd
+        segdNEx1    = U.mkSegd (U.replicate nNE 1)          -- Segd ( [1, 1, 1, 1]
+                               (U.enumFromStepLen 0 1 nNE)  --      , [0, 1, 3, 4]
+                               nNE                          --      , 4 )
+        segdNEx2    = segdNEx1 `U.plusSegd` segdNEx1        -- Segd ( [2, 2, 2, 2]
+                                                            --      , [0, 2, 4, 8]
+                                                            --      , 8 )
+        (lensx2,idsx2,neltsx2)
+                    = getSegdParams segdNEx2
+        -- recreate segdNE for some reason
+        segdNE      = U.lengthsToSegd                       -- Segd ( [4, 5, 2, 2]
+                    $ U.packByTag segd_lens tagsE 0         --      , [0, 4, 9, 11]
+                                                            --      , 13)
+        (lensNE,idsNE,_)
+                    = getSegdParams segdNE
+        tagsE_r     = U.replicate_s segd tagsE              -- [0, 0, 0, 0,     -- 4
+                                                            --  0, 0, 0, 0, 0,  -- 5
+                                                            --                  -- 0
+                                                            --  0, 0,           -- 2
+                                                            --  0, 0]           -- 2
+                                                            --                  -- 0
+        -- WATCH OUT: redundant packByTag. Basicaly remove all empty segments from segd
+        packedNE@(PNested _ (P_2 (PDouble xos') (PDouble yos')))
+                    = PNested segdNE (P_2 (PDouble $ U.packByTag xos tagsE_r 0)
+                                          (PDouble $ U.packByTag yos tagsE_r 0))
+        -- prepare segd to call hsplit anew
+        lens_r      = U.replicate_s segdNEx2 lensNE        -- Segd ( [4, 4, 5, 5, 2, 2, 2, 2]
+        segd_r      = U.lengthsToSegd lens_r               --      , [0, 4, 8,13,18,20,22,24]
+        (_,ixs_r,elts_r) = getSegdParams segd_r            --      ,  26
+        ids_to_take = U.enumFromStepLenEach elts_r                         -- 26
+                                            (U.replicate_s segdNEx2 idsNE) -- [0, 0, 4, 4, 9, 9,11,11]
+                                            (U.replicate   neltsx2  1)     -- [1, 1, 1, 1, 1, 1, 1, 1]
+                                            lens_r                         -- [4, 4, 5, 5, 2, 2, 2, 2]
+        -- ids_to_take: [0,1,2,3, 0,1,2,3, 4,5,6,7,8, 4,5,6,7,8, 9,10, 9,10, 11,12, 11,12]
+        packed_r@(PNested _ (P_2 (PDouble xos_r) (PDouble yos_r)))
+                    = PNested segd_r (P_2 (PDouble $ U.bpermute xos' ids_to_take)
+                                          (PDouble $ U.bpermute yos' ids_to_take))
+        pm_s@(P_2 (PDouble pm_xs) (PDouble pm_ys))
+                    = calc_pm_s points_s cross_s tagsE
+        line_1s'    = P_2 (PDouble $ U.append_s segdNEx2
+                                                segdNEx1
+                                                (U.packByTag line_x1s tagsE 0)
+                                                segdNEx1
+                                                pm_xs)
+                          (PDouble $ U.append_s segdNEx2
+                                                segdNEx1
+                                                (U.packByTag line_y1s tagsE 0)
+                                                segdNEx1
+                                                pm_ys)
+        line_2s'    = P_2 (PDouble $ U.append_s segdNEx2
+                                                segdNEx1
+                                                pm_xs
+                                                segdNEx1
+                                                (U.packByTag line_x2s tagsE 0))
+                          (PDouble $ U.append_s segdNEx2
+                                                segdNEx1
+                                                pm_ys
+                                                segdNEx1
+                                                (U.packByTag line_y2s tagsE 0))
+        lines'      = P_2 line_1s' line_2s'
+        -- finally make recursive call to compute convect hull from (when the recursion hasn't finished)
+        hullNE@(PNested hullNE_segd (P_2 (PDouble hullNE_xs) (PDouble hullNE_ys)))
+                    = hsplit_l_wrk nNE# packed_r lines'
+        (hullNE_lens, hullNE_ids, hullNE_nelts)
+                    = getSegdParams hullNE_segd
+        --- *** now deal with the empty segnments, i.e. lines which have no points above them *** ---
+        --- *** and will thus form part of the convex hull                                    *** ---
+        !nE@(I# nE#) = U.elementsSel2_1 selE                -- 2
+        result_segd = U.lengthsToSegd
+                    $ U.combine2 (U.tagsSel2 selE)
+                                 (U.repSel2  selE)
+                                 (U.sum_s segdNEx2 hullNE_lens) --
+                                 (U.replicate nE 1)        -- [1, 1]
+        result_sel  = U.tagsToSel2 (U.replicate_s result_segd tagsE)
+        result_tags = U.tagsSel2  result_sel
+        result_repsel = U.repSel2 result_sel
+        result_xs   = U.combine2 result_tags
+                                 result_repsel
+                                 hullNE_xs
+                                 (U.packByTag line_x1s tagsE 1)
+        result_ys   = U.combine2 result_tags
+                                 result_repsel
+                                 hullNE_ys
+                                 (U.packByTag line_y1s tagsE 1)
+    in  (PNested result_segd (P_2 (PDouble result_xs) (PDouble result_ys)))
+{-# INLINE hsplit_l_wrk #-}
+
+getSegdParams :: Segd -> (Array Int, Array Int, Int)
+getSegdParams segd = ( U.lengthsSegd  segd
+                     , U.indicesSegd segd
+                     , U.elementsSegd segd )
+{-# INLINE getSegdParams #-}
+
+calc_cross_s :: PData (PArray Point) -> PData Line -> PData (PArray Double)
+calc_cross_s (PNested segd (P_2 (PDouble xos) (PDouble yos)))
+             (P_2 (P_2 (PDouble x1s) (PDouble y1s))
+                  (P_2 (PDouble x2s) (PDouble y2s)))
+  = let distances = U.zipWith6 distance xos
+                                        yos
+                                        (U.replicate_s segd x1s)
+                                        (U.replicate_s segd y1s)
+                                        (U.replicate_s segd x2s)
+                                        (U.replicate_s segd y2s)
+    in PNested segd (PDouble distances)
+{-# INLINE calc_cross_s #-}
+
+calc_packed_s :: PData (PArray Point) -> PData (PArray Double) -> PData (PArray Point)
+calc_packed_s (PNested segd (P_2 (PDouble xos) (PDouble yos)))
+              (PNested _    (PDouble cross))
+  = let nelts  = U.elementsSegd segd
+        -- compute selector for positive elements ((>0) -> 1; otherwise -> 0)
+        selGZ  = U.tagsToSel2
+               $ U.map B.fromBool
+               $ U.zipWith (P.>) cross (U.replicate nelts 0.0)
+        -- compute segd for just the positive elements
+        segdGZ = U.lengthsToSegd (U.count_s segd (tagsSel2 selGZ) 1)
+        -- get the actual points corresponding to positive elements
+        tagsGZ = U.tagsSel2 selGZ
+        xosGZ  = U.packByTag xos tagsGZ 1
+        yosGZ  = U.packByTag yos tagsGZ 1
+    in  (PNested segdGZ (P_2 (PDouble xosGZ) (PDouble yosGZ)))
+{-# INLINE calc_packed_s #-}
+
+calc_pm_s :: PData (PArray Point) -> PData (PArray Double) -> Array Tag -> PData Point
+calc_pm_s (PNested segd_pts (P_2 (PDouble x1s) (PDouble y1s))) -- points_s
+          (PNested segd_crs (PDouble cross))                   -- cross_s
+          tags
+  = let idexed = U.zip (U.indices_s segd_crs) cross
+        max_ids = U.fsts (U.fold1_s max' segd_crs idexed) 
+        -- find indices to take from the points array
+        ids     = U.zipWith (P.+) (U.indicesSegd segd_pts)
+                                max_ids
+        max_xs  = U.bpermute x1s ids
+        max_ys  = U.bpermute y1s ids
+        -- however we are only interested in the ones which are from the non-empty segments
+        -- (TODO revise comment)
+    in  P_2 (PDouble $ U.packByTag max_xs tags 0)
+            (PDouble $ U.packByTag max_ys tags 0)
+{-# INLINE calc_pm_s #-}
+
+-- from dph-lifted-copy:D.A.P.Prelude.Int
+max' :: P.Ord b => (a, b) -> (a, b) -> (a, b)
+max' (i,x) (j,y) | x P.>= y    = (i,x)
+                 | P.otherwise = (j,y)
+{-# INLINE max' #-}
+
+min' :: P.Ord b => (a, b) -> (a, b) -> (a, b)
+min' (i,x) (j,y) | x P.<= y    = (i,x)
+                 | P.otherwise = (j,y)
+{-# INLINE min' #-}
+
+
diff --git a/dph-examples/examples/spectral/QuickHull/handvec/HandvecWrp.hs b/dph-examples/examples/spectral/QuickHull/handvec/HandvecWrp.hs
new file mode 100644 (file)
index 0000000..347554a
--- /dev/null
@@ -0,0 +1,37 @@
+{-# LANGUAGE MagicHash, BangPatterns, ParallelArrays, TypeOperators #-}
+{-# OPTIONS_GHC -fvectorise #-}
+
+module HandvecWrp ( quickhullPA )
+where
+
+import Points2D.Types
+
+import Data.Array.Parallel as PA
+import Data.Array.Parallel.Prelude.Double as D
+import qualified Data.Array.Parallel.Prelude.Int as I
+
+import Prelude as P
+
+import Handvec
+
+-- Manually override vectorisation of hsplit
+-- hsplit :: [:Point:] -> Line -> [:Point:]
+hsplit :: PArray Point -> Line -> PArray Point
+hsplit ps _ = ps -- should never be called
+{-# VECTORISE hsplit = hsplit_v #-}
+-- vectorisation NOT overridden without NOINLINE pragma:
+{-# NOINLINE hsplit #-}
+
+quickHull :: [:Point:] -> [:Point:]
+quickHull points
+  | lengthP points I.== 0 = points
+  | otherwise
+  = concatP [: fromPArrayP (hsplit (toPArrayP points) ends) | ends <- [: (minx, maxx), (maxx, minx) :] :]
+  where
+    xs   = [: x | (x, y) <- points :]
+    minx = points PA.!: minIndexP xs
+    maxx = points PA.!: maxIndexP xs
+
+quickhullPA :: PArray Point -> PArray Point
+quickhullPA ps = toPArrayP (quickHull (fromPArrayP ps))
+{-# NOINLINE quickhullPA #-}
diff --git a/dph-examples/examples/spectral/QuickHull/handvec/Main.hs b/dph-examples/examples/spectral/QuickHull/handvec/Main.hs
new file mode 100644 (file)
index 0000000..3ca732a
--- /dev/null
@@ -0,0 +1,59 @@
+
+import HandvecWrp
+import SVG
+import Timing
+import Points2D.Generate
+import Points2D.Types
+import System.Environment
+import Data.Array.Parallel.PArray      as P
+
+main :: IO ()
+main 
+ = do  args    <- getArgs
+       case args of
+         [pointCount]  
+           -> run (read pointCount) Nothing
+       
+         [pointCount, fileSVG]
+           -> run (read pointCount) (Just fileSVG)
+
+         _ -> do
+               putStr usage
+               return ()
+
+
+-- | Command line usage information.
+usage :: String
+usage  = unlines
+       [ "Usage: quickhull <points> [out.svg]" ]
+
+
+-- | Run the benchmark.
+run    :: Int                  -- ^ How many points to use.
+       -> Maybe String         -- ^ File name to dump an SVG of the output to.
+       -> IO ()
+       
+run pointCount mFileSVG
+ = do
+       vPoints <- pointsPArrayOfUArray
+               $ genPointsDisc pointCount (400, 400) 350 
+
+       -- Force points to create the input vector.
+       vPoints `seq` return ()
+
+       -- Compute the convex hull.
+       (vHull, tElapsed)
+               <- time 
+               $  let  vHull   = quickhullPA vPoints
+                  in   vHull `seq` return vHull
+                                       
+       -- Print how long it took.
+       putStr $ prettyTime tElapsed
+       
+       -- If we were asked for an SVG then write it out to file.
+       maybe   (return ())
+               (\fileSVG -> 
+                       writeFile fileSVG
+                        $ makeSVG      (roundPoints $ P.toList vPoints) 
+                                       (roundPoints $ P.toList vHull))
+               mFileSVG