src/graph.bend
bend-visgraph@1.0.0.0 · 517 lines · 0 definitions
1# Natural and horizontal visibility graphs, as ts2vg 1.2.4 builds them, in
2# F32. Each def below that has a ts2vg twin does the same F32 operations in
3# the same order as the twin does in float64.
4import Base
5
6type Graph is Data:
7 Natural{}
8 Horizontal{}
9
10type Direction is Data:
11 Undirected{}
12 LeftToRight{}
13 TopToBottom{}
14
15type Weight is Data:
16 Unweighted{}
17 Distance{}
18 SqDistance{}
19 VDistance{}
20 AbsVDistance{}
21 HDistance{}
22 AbsHDistance{}
23 Slope{}
24 AbsSlope{}
25 Angle{}
26 AbsAngle{}
27 NumPenetrations{}
28
29# The x coordinate of each value: its index, or the value at the same
30# index of an array.
31type Xs is Type:
32 Index{}
33 At{xs: Array<F32>}
34
35# The edges, as a rope: ECon adds an edge and EJoin joins two ropes, so the
36# two halves of a split join in one step.
37type Edges is Data:
38 ENil{}
39 ECon{src: U32, dst: U32, w: F32, next: Edges}
40 EJoin{a: Edges, b: Edges}
41
42# F32 constants, from divisions by zero: the JS lane builds an F32 from its
43# bits with two typed arrays, too slow for a sweep that passes NaN at every
44# point it visits.
45
46def nan() -> F32:
47 F32.div(0.0, 0.0)
48
49def inf() -> F32:
50 F32.div(1.0, 0.0)
51
52def neg_inf() -> F32:
53 F32.div(F32.neg(1.0), 0.0)
54
55def is_nan(+a: F32) -> Bool:
56 F32.is_ne(a, a)
57
58def pick(c: Bool, a: F32, b: F32) -> F32:
59 match c:
60 case True{}:
61 a
62 case False{}:
63 b
64
65def pick_u32(c: Bool, a: U32, b: U32) -> U32:
66 match c:
67 case True{}:
68 a
69 case False{}:
70 b
71
72# Weights: ts2vg's _weight_* functions. `s` is the slope the sweep measured,
73# or NaN when it measured none.
74
75def sq_distance(+x1: F32, +x2: F32, +y1: F32, +y2: F32) -> F32:
76 F32.add(F32.mul(F32.sub(x2, x1), F32.sub(x2, x1)), F32.mul(F32.sub(y2, y1), F32.sub(y2, y1)))
77
78def known_slope(+s: F32, +x1: F32, +x2: F32, +y1: F32, +y2: F32) -> F32:
79 pick(is_nan(s), F32.div(F32.sub(y2, y1), F32.sub(x2, x1)), s)
80
81def weight(+k: Weight, +x1: F32, +x2: F32, +y1: F32, +y2: F32, +s: F32) -> F32:
82 match k:
83 case Unweighted{}:
84 0.0
85 case Distance{}:
86 F32.sqrt(sq_distance(x1, x2, y1, y2))
87 case SqDistance{}:
88 sq_distance(x1, x2, y1, y2)
89 case VDistance{}:
90 F32.sub(y2, y1)
91 case AbsVDistance{}:
92 F32.abs(F32.sub(y2, y1))
93 case HDistance{}:
94 F32.sub(x2, x1)
95 case AbsHDistance{}:
96 F32.abs(F32.sub(x2, x1))
97 case Slope{}:
98 known_slope(s, x1, x2, y1, y2)
99 case AbsSlope{}:
100 F32.abs(known_slope(s, x1, x2, y1, y2))
101 case Angle{}:
102 F32.atan(known_slope(s, x1, x2, y1, y2))
103 case AbsAngle{}:
104 F32.atan(F32.abs(known_slope(s, x1, x2, y1, y2)))
105 case NumPenetrations{}:
106 0.0
107
108# Edges: ts2vg's add_edge, which drops an edge whose weight is at most lo or
109# at least hi.
110
111def put(drop: Bool, i1: U32, i2: U32, w: F32, es: Edges) -> Edges:
112 match drop:
113 case True{}:
114 es
115 case False{}:
116 ECon{i1, i2, w, es}
117
118def add_edge(+k: Weight, +lo: F32, +hi: F32, i1: U32, i2: U32, +x1: F32, +x2: F32, +y1: F32, +y2: F32, +s: F32, es: Edges) -> Edges:
119 +w = weight(k, x1, x2, y1, y2, s)
120 put(F32.is_le(w, lo) || F32.is_ge(w, hi), i1, i2, w, es)
121
122# A visible point right of the pivot a: the edge (a, b).
123def add_right(seen: Bool, +k: Weight, +lo: F32, +hi: F32, ia: U32, ib: U32, xa: F32, xb: F32, ya: F32, yb: F32, s: F32, es: Edges) -> Edges:
124 match seen:
125 case True{}:
126 add_edge(k, lo, hi, ia, ib, xa, xb, ya, yb, s, es)
127 case False{}:
128 es
129
130# A visible point left of the pivot a: the edge (a, b) for top_to_bottom,
131# else (b, a).
132def add_left_dir(dir: Direction, +k: Weight, +lo: F32, +hi: F32, ia: U32, ib: U32, xa: F32, xb: F32, ya: F32, yb: F32, s: F32, es: Edges) -> Edges:
133 match dir:
134 case TopToBottom{}:
135 add_edge(k, lo, hi, ia, ib, xa, xb, ya, yb, s, es)
136 case _:
137 add_edge(k, lo, hi, ib, ia, xb, xa, yb, ya, s, es)
138
139def add_left(seen: Bool, +dir: Direction, +k: Weight, +lo: F32, +hi: F32, ia: U32, ib: U32, xa: F32, xb: F32, ya: F32, yb: F32, s: F32, es: Edges) -> Edges:
140 match seen:
141 case True{}:
142 add_left_dir(dir, k, lo, hi, ia, ib, xa, xb, ya, yb, s, es)
143 case False{}:
144 es
145
146# Sweeps from the pivot a over the points of one side, nearest first.
147
148def abs_tol() -> F32:
149 {0.00000000000001 : F32}
150
151def rel_tol() -> F32:
152 {0.00000000000001 : F32}
153
154# max(|xa|, |ya|), which a sweep computes once.
155def scale(+xa: F32, +ya: F32) -> F32:
156 F32.max(F32.abs(xa), F32.abs(ya))
157
158# ts2vg's tolerance for the points a and b.
159# F32.max only picks one of its words, so grouping the four maxima this way
160# gives ts2vg's word.
161def tol(+scale_a: F32, +xb: F32, +yb: F32) -> F32:
162 F32.max(abs_tol(), F32.mul(rel_tol(), F32.max(scale_a, F32.max(F32.abs(xb), F32.abs(yb)))))
163
164def greater(+a: F32, +b: F32, +t: F32) -> Bool:
165 F32.is_gt(F32.sub(a, b), t)
166
167# greater(s, best, tol(scale_a, xb, yb)), with the tolerance computed only
168# when s is above best: the tolerance is positive, so otherwise it fails
169# anyway.
170def clears(up: Bool, +s: F32, +best: F32, +scale_a: F32, +xb: F32, +yb: F32) -> Bool:
171 match up:
172 case True{}:
173 greater(s, best, tol(scale_a, xb, yb))
174 case False{}:
175 False{}
176
177# Divide and conquer: ts2vg's _compute_graph takes an interval's first
178# highest point as its pivot, sweeps the interval's points left and right of
179# it, and recurses on both sides. Every point is the pivot of one interval,
180# so once each point knows how many points its interval holds on each side,
181# the sweeps no longer depend on each other.
182
183# The stack holds the points no later point has beaten yet. A point pops the
184# entries it beats; the entry left on top bounds its interval on the left,
185# and the point that pops an entry bounds that entry's interval on the right.
186type Stack is Data:
187 KNil{}
188 KCon{i: U32, key: F32, next: Stack}
189
190# What a later point must beat to pop y. An equal point beats nothing, so an
191# interval's first highest point is its pivot. A NaN beats nothing either,
192# and ts2vg's argmax takes it only as an interval's first point, so a NaN at
193# i is the pivot of the right side of i - 1: it takes the key of the entry
194# below it and leaves the stack with that entry.
195def key(+y: F32, +below: F32) -> F32:
196 pick(is_nan(y), below, y)
197
198# The work only balances the split, so F32 rounding does no harm.
199type Popped is Type:
200 Popped{first: U32, below: F32, st: Stack, rights: Array<U32>, visits: F32}
201
202def pops(rest: Stack, beat: Bool, +y: F32, +i: U32, +j: U32, +kj: F32, rights: Array<U32>, +visits: F32) -> Popped:
203 match rest beat:
204 case KCon{+j2, +k2, rest2} True{}:
205 +r = (i - 1 - j : U32)
206 pops(rest2, F32.is_gt(y, k2), y, i, j2, k2, Array.set(U32, rights, j, r), F32.add(visits, U32.to_f32(r)))
207 case KNil{} True{}:
208 +r = (i - 1 - j : U32)
209 Popped{0, inf(), KNil{}, Array.set(U32, rights, j, r), F32.add(visits, U32.to_f32(r))}
210 case _ False{}:
211 Popped{(j + 1 : U32), kj, KCon{j, kj, rest}, rights, visits}
212
213def flush(st: Stack, +n: U32, rights: Array<U32>, +visits: F32) -> Array<U32> & F32:
214 match st:
215 case KCon{+j, k, rest}:
216 +r = (n - 1 - j : U32)
217 flush(rest, n, Array.set(U32, rights, j, r), F32.add(visits, U32.to_f32(r)))
218 case KNil{}:
219 (rights, visits)
220
221# lefts and rights count the points each point's left and right sweeps visit.
222type Cols is Type:
223 Cols{ys: Array<F32>, lefts: Array<U32>, rights: Array<U32>}
224
225type Arrived is Type:
226 OpensWithNaN{rights: Array<U32>, visits: F32}
227 PopsAll{rights: Array<U32>, visits: F32}
228 Bounded{first: U32, below: F32, st: Stack, rights: Array<U32>, visits: F32}
229
230def lpops(rest: Stack, beat: Bool, +y: F32, +i: U32, +j: U32, +kj: F32, rights: Array<U32>, +visits: F32) -> Arrived:
231 match rest beat:
232 case KCon{+j2, +k2, rest2} True{}:
233 +r = (i - 1 - j : U32)
234 lpops(rest2, F32.is_gt(y, k2), y, i, j2, k2, Array.set(U32, rights, j, r), F32.add(visits, U32.to_f32(r)))
235 case KNil{} True{}:
236 +r = (i - 1 - j : U32)
237 PopsAll{Array.set(U32, rights, j, r), F32.add(visits, U32.to_f32(r))}
238 case _ False{}:
239 Bounded{(j + 1 : U32), kj, KCon{j, kj, rest}, rights, visits}
240
241def arrive(st: Stack, nan: Bool, +y: F32, +i: U32, rights: Array<U32>, +visits: F32) -> Arrived:
242 match st nan:
243 case KNil{} True{}:
244 OpensWithNaN{rights, visits}
245 case KNil{} False{}:
246 PopsAll{rights, visits}
247 case KCon{+j, +kj, rest} True{}:
248 Bounded{(j + 1 : U32), kj, KCon{j, kj, rest}, rights, visits}
249 case KCon{+j, +kj, rest} False{}:
250 lpops(rest, F32.is_gt(y, kj), y, i, j, kj, rights, visits)
251
252# Latest first.
253type Unbounded is Data:
254 UNil{}
255 UCon{i: U32, y: F32, next: Unbounded}
256
257# A block of points runs the stack on its own points. Its points that pop
258# its whole stack, and the NaNs that open it, need the points before it.
259type Block is Data:
260 Block{start: U32, opening_nans: U32, unbounded: Unbounded, st: Stack, visits: F32}
261
262type Local is Type:
263 Local{b: Block, c: Cols}
264
265def local(count: Nat, +i: U32, +y: F32, a: Arrived, yr: Array<F32> & F32, lefts: Array<U32>, +start: U32, +opening_nans: U32, unbounded: Unbounded) -> Local:
266 match count a:
267 case 1n+c OpensWithNaN{rights, +visits}:
268 (ys, +y2) = yr
269 +i2 = (i + 1 : U32)
270 local(c, i2, y2, arrive(KNil{}, is_nan(y2), y2, i2, rights, F32.add(visits, 1.0)), Array.get(F32, ys, (i2 + 1 : U32)), Array.set(U32, lefts, i, 0), start, (opening_nans + 1 : U32), unbounded)
271 case 1n+c PopsAll{rights, +visits}:
272 (ys, +y2) = yr
273 +i2 = (i + 1 : U32)
274 local(c, i2, y2, arrive(KCon{i, y, KNil{}}, is_nan(y2), y2, i2, rights, visits), Array.get(F32, ys, (i2 + 1 : U32)), lefts, start, opening_nans, UCon{i, y, unbounded})
275 case 1n+c Bounded{+first, +below, st, rights, +visits}:
276 (ys, +y2) = yr
277 +i2 = (i + 1 : U32)
278 +l = (i - first : U32)
279 local(c, i2, y2, arrive(KCon{i, key(y, below), st}, is_nan(y2), y2, i2, rights, F32.add(visits, U32.to_f32((l + 1 : U32)))), Array.get(F32, ys, (i2 + 1 : U32)), Array.set(U32, lefts, i, l), start, opening_nans, unbounded)
280 case 0n OpensWithNaN{rights, +visits}:
281 (ys, y2) = yr
282 Local{Block{start, (opening_nans + 1 : U32), unbounded, KNil{}, F32.add(visits, 1.0)}, Cols{ys, Array.set(U32, lefts, i, 0), rights}}
283 case 0n PopsAll{rights, +visits}:
284 (ys, y2) = yr
285 Local{Block{start, opening_nans, UCon{i, y, unbounded}, KCon{i, y, KNil{}}, visits}, Cols{ys, lefts, rights}}
286 case 0n Bounded{+first, +below, st, rights, +visits}:
287 (ys, y2) = yr
288 +l = (i - first : U32)
289 Local{Block{start, opening_nans, unbounded, KCon{i, key(y, below), st}, F32.add(visits, U32.to_f32((l + 1 : U32)))}, Cols{ys, Array.set(U32, lefts, i, l), rights}}
290
291def block_pass(count: Nat, +s: U32, yr: Array<F32> & F32, lefts: Array<U32>, rights: Array<U32>) -> Local:
292 match count:
293 case 0n:
294 (ys, y) = yr
295 Local{Block{s, 0, UNil{}, KNil{}, 0.0}, Cols{ys, lefts, rights}}
296 case 1n+c:
297 (ys, +y) = yr
298 local(c, s, y, arrive(KNil{}, is_nan(y), y, s, rights, 0.0), Array.get(F32, ys, (s + 1 : U32)), lefts, s, 0, UNil{})
299
300type Merged is Type:
301 Merged{st: Stack, lefts: Array<U32>, rights: Array<U32>, visits: F32}
302
303def push_with_key(count: Nat, +i: U32, +k: F32, st: Stack) -> Stack:
304 match count:
305 case 0n:
306 st
307 case 1n+c:
308 push_with_key(c, (i + 1 : U32), k, KCon{i, k, st})
309
310def push_opening_nans(count: Nat, +i: U32, st: Stack) -> Stack:
311 match st:
312 case KCon{+j, +kj, rest}:
313 push_with_key(count, i, kj, KCon{j, kj, rest})
314 case KNil{}:
315 push_with_key(count, i, inf(), KNil{})
316
317def pops_from(st: Stack, +y: F32, +i: U32, rights: Array<U32>, +visits: F32) -> Popped:
318 match st:
319 case KCon{+j, +kj, rest}:
320 pops(rest, F32.is_gt(y, kj), y, i, j, kj, rights, visits)
321 case KNil{}:
322 Popped{0, inf(), KNil{}, rights, visits}
323
324def set_left(+i: U32, lefts: Array<U32>, p: Popped) -> Merged:
325 Popped{+first, below, st, rights, +visits} = p
326 +l = (i - first : U32)
327 Merged{st, Array.set(U32, lefts, i, l), rights, F32.add(visits, U32.to_f32((l + 1 : U32)))}
328
329def in_order(u: Unbounded, acc: Unbounded) -> Unbounded:
330 match u:
331 case UNil{}:
332 acc
333 case UCon{i, y, next}:
334 in_order(next, UCon{i, y, acc})
335
336def bound_unbounded(u: Unbounded, m: Merged) -> Merged:
337 match u:
338 case UNil{}:
339 m
340 case UCon{+i, +y, next}:
341 Merged{st, lefts, rights, +visits} = m
342 bound_unbounded(next, set_left(i, lefts, pops_from(st, y, i, rights, visits)))
343
344def onto(st: Stack, acc: Stack) -> Stack:
345 match st:
346 case KNil{}:
347 acc
348 case KCon{j, k, rest}:
349 onto(rest, KCon{j, k, acc})
350
351def push_block_stack(st: Stack, m: Merged) -> Merged:
352 Merged{g, lefts, rights, visits} = m
353 Merged{onto(onto(st, KNil{}), g), lefts, rights, visits}
354
355# The merge replays the stack over every point before the block. The points
356# the block's own stack could not bound pop it in order; the others sit on
357# the block's stack, which goes on top.
358def merge_block(b: Block, m: Merged) -> Merged:
359 Block{+start, +opening_nans, unbounded, st, +visits} = b
360 Merged{g, lefts, rights, +v} = m
361 push_block_stack(st, bound_unbounded(in_order(unbounded, UNil{}), Merged{push_opening_nans(U32.to_nat(opening_nans), start, g), lefts, rights, F32.add(v, visits)}))
362
363# Bend matches only parameters, so a read's result must arrive as one: a
364# sweep reads y and x at `ahead`, where the next sweep starts, and hands them
365# over. The sweeps read x through `x_at`, so x as the index costs no read.
366type Walked is Type:
367 Walked{yr: Array<F32> & F32, xr: Array<F32> & F32, es: Edges}
368
369def at_index(xs: Array<F32>, +j: U32) -> Array<F32> & F32:
370 (xs, U32.to_f32(j))
371
372def natural_right(~x_at: Array<F32> -> U32 -> Array<F32> & F32, count: Nat, +j: U32, yr: Array<F32> & F32, xr: Array<F32> & F32, +k: Weight, +lo: F32, +hi: F32, +ia: U32, +xa: F32, +ya: F32, +scale_a: F32, +best: F32, es: Edges, +ahead: U32) -> Walked:
373 match count:
374 case 0n:
375 (ys, yb) = yr
376 (xs, xb) = xr
377 Walked{Array.get(F32, ys, ahead), Array.get(F32, xs, ahead), es}
378 case 1n+c:
379 (ys, +yb) = yr
380 (xs, +xb) = xr
381 +s = F32.div(F32.sub(yb, ya), F32.sub(xb, xa))
382 +seen = clears(F32.is_gt(s, best), s, best, scale_a, xb, yb)
383 +j2 = (j + 1 : U32)
384 natural_right(~x_at, c, j2, Array.get(F32, ys, j2), x_at(xs, j2), k, lo, hi, ia, xa, ya, scale_a, pick(seen, s, best), add_right(seen, k, lo, hi, ia, j, xa, xb, ya, yb, s, es), ahead)
385
386def natural_left(~x_at: Array<F32> -> U32 -> Array<F32> & F32, count: Nat, +j: U32, yr: Array<F32> & F32, xr: Array<F32> & F32, +dir: Direction, +k: Weight, +lo: F32, +hi: F32, +ia: U32, +xa: F32, +ya: F32, +scale_a: F32, +best: F32, es: Edges, +ahead: U32) -> Walked:
387 match count:
388 case 0n:
389 (ys, yb) = yr
390 (xs, xb) = xr
391 Walked{Array.get(F32, ys, ahead), Array.get(F32, xs, ahead), es}
392 case 1n+c:
393 (ys, +yb) = yr
394 (xs, +xb) = xr
395 +s = F32.div(F32.sub(yb, ya), F32.neg(F32.sub(xb, xa)))
396 +seen = clears(F32.is_gt(s, best), s, best, scale_a, xb, yb)
397 +j2 = (j - 1 : U32)
398 natural_left(~x_at, c, j2, Array.get(F32, ys, j2), x_at(xs, j2), dir, k, lo, hi, ia, xa, ya, scale_a, pick(seen, s, best), add_left(seen, dir, k, lo, hi, ia, j, xa, xb, ya, yb, F32.neg(s), es), ahead)
399
400def horizontal_right(~x_at: Array<F32> -> U32 -> Array<F32> & F32, count: Nat, +j: U32, yr: Array<F32> & F32, xr: Array<F32> & F32, +k: Weight, +lo: F32, +hi: F32, +ia: U32, +xa: F32, +ya: F32, +best: F32, es: Edges, +ahead: U32) -> Walked:
401 match count:
402 case 0n:
403 (ys, yb) = yr
404 (xs, xb) = xr
405 Walked{Array.get(F32, ys, ahead), Array.get(F32, xs, ahead), es}
406 case 1n+c:
407 (ys, +yb) = yr
408 (xs, xb) = xr
409 +seen = F32.is_gt(yb, best)
410 +j2 = (j + 1 : U32)
411 horizontal_right(~x_at, c, j2, Array.get(F32, ys, j2), x_at(xs, j2), k, lo, hi, ia, xa, ya, pick(seen, yb, best), add_right(seen, k, lo, hi, ia, j, xa, xb, ya, yb, nan(), es), ahead)
412
413def horizontal_left(~x_at: Array<F32> -> U32 -> Array<F32> & F32, count: Nat, +j: U32, yr: Array<F32> & F32, xr: Array<F32> & F32, +dir: Direction, +k: Weight, +lo: F32, +hi: F32, +ia: U32, +xa: F32, +ya: F32, +best: F32, es: Edges, +ahead: U32) -> Walked:
414 match count:
415 case 0n:
416 (ys, yb) = yr
417 (xs, xb) = xr
418 Walked{Array.get(F32, ys, ahead), Array.get(F32, xs, ahead), es}
419 case 1n+c:
420 (ys, +yb) = yr
421 (xs, xb) = xr
422 +seen = F32.is_gt(yb, best)
423 +j2 = (j - 1 : U32)
424 horizontal_left(~x_at, c, j2, Array.get(F32, ys, j2), x_at(xs, j2), dir, k, lo, hi, ia, xa, ya, pick(seen, yb, best), add_left(seen, dir, k, lo, hi, ia, j, xa, xb, ya, yb, nan(), es), ahead)
425
426# Each sweep, with x as the index (`index`) or read from `xs`.
427def right_sweep(g: Graph, index: Bool, count: Nat, yr: Array<F32> & F32, xs: Array<F32>, +k: Weight, +lo: F32, +hi: F32, +ia: U32, +xa: F32, +ya: F32, es: Edges) -> Walked:
428 match g index:
429 case Natural{} True{}:
430 natural_right(~(a => +j => at_index(a, j)), count, (ia + 1 : U32), yr, at_index(xs, (ia + 1 : U32)), k, lo, hi, ia, xa, ya, scale(xa, ya), neg_inf(), es, (ia - 1 : U32))
431 case Natural{} False{}:
432 natural_right(~(a => +j => Array.get(F32, a, j)), count, (ia + 1 : U32), yr, Array.get(F32, xs, (ia + 1 : U32)), k, lo, hi, ia, xa, ya, scale(xa, ya), neg_inf(), es, (ia - 1 : U32))
433 case Horizontal{} True{}:
434 horizontal_right(~(a => +j => at_index(a, j)), count, (ia + 1 : U32), yr, at_index(xs, (ia + 1 : U32)), k, lo, hi, ia, xa, ya, neg_inf(), es, (ia - 1 : U32))
435 case Horizontal{} False{}:
436 horizontal_right(~(a => +j => Array.get(F32, a, j)), count, (ia + 1 : U32), yr, Array.get(F32, xs, (ia + 1 : U32)), k, lo, hi, ia, xa, ya, neg_inf(), es, (ia - 1 : U32))
437
438def left_sweep(g: Graph, index: Bool, count: Nat, +dir: Direction, +k: Weight, +lo: F32, +hi: F32, +ia: U32, +xa: F32, +ya: F32, w: Walked) -> Walked:
439 match g index:
440 case Natural{} True{}:
441 Walked{yr, xr, es} = w
442 (xs, xv) = xr
443 natural_left(~(a => +j => at_index(a, j)), count, (ia - 1 : U32), yr, at_index(xs, (ia - 1 : U32)), dir, k, lo, hi, ia, xa, ya, scale(xa, ya), neg_inf(), es, (ia + 1 : U32))
444 case Natural{} False{}:
445 Walked{yr, xr, es} = w
446 (xs, xv) = xr
447 natural_left(~(a => +j => Array.get(F32, a, j)), count, (ia - 1 : U32), yr, Array.get(F32, xs, (ia - 1 : U32)), dir, k, lo, hi, ia, xa, ya, scale(xa, ya), neg_inf(), es, (ia + 1 : U32))
448 case Horizontal{} True{}:
449 Walked{yr, xr, es} = w
450 (xs, xv) = xr
451 horizontal_left(~(a => +j => at_index(a, j)), count, (ia - 1 : U32), yr, at_index(xs, (ia - 1 : U32)), dir, k, lo, hi, ia, xa, ya, neg_inf(), es, (ia + 1 : U32))
452 case Horizontal{} False{}:
453 Walked{yr, xr, es} = w
454 (xs, xv) = xr
455 horizontal_left(~(a => +j => Array.get(F32, a, j)), count, (ia - 1 : U32), yr, Array.get(F32, xs, (ia - 1 : U32)), dir, k, lo, hi, ia, xa, ya, neg_inf(), es, (ia + 1 : U32))
456
457type Out is Type:
458 Out{cols: Cols, xs: Array<F32>, es: Edges}
459
460def points(count: Nat, +i: U32, lr: Array<U32> & U32, rr: Array<U32> & U32, w: Walked, +index: Bool, +g: Graph, +dir: Direction, +k: Weight, +lo: F32, +hi: F32) -> Out:
461 match count:
462 case 0n:
463 (lefts, l) = lr
464 (rights, r) = rr
465 Walked{yr, xr, es} = w
466 (ys, y) = yr
467 (xs, x) = xr
468 Out{Cols{ys, lefts, rights}, xs, es}
469 case 1n+c:
470 (lefts, +l) = lr
471 (rights, +r) = rr
472 Walked{yr, xr, es} = w
473 (ys, +y) = yr
474 (xs, +xv) = xr
475 +x = pick(index, U32.to_f32(i), xv)
476 +i2 = (i + 1 : U32)
477 points(c, i2, Array.get(U32, lefts, i2), Array.get(U32, rights, i2), left_sweep(g, index, U32.to_nat(l), dir, k, lo, hi, i, x, y, right_sweep(g, index, U32.to_nat(r), Array.get(F32, ys, i2), xs, k, lo, hi, i, x, y, es)), index, g, dir, k, lo, hi)
478
479def edges_of(o: Out) -> Edges:
480 Out{cols, xs, es} = o
481 es
482
483def swept(+index: Bool, +g: Graph, +dir: Direction, +k: Weight, +lo: F32, +hi: F32, +n: U32, ys: Array<F32>, xs: Array<F32>, lefts: Array<U32>, r: Array<U32> & F32) -> Edges:
484 (rights, visits) = r
485 edges_of(points(U32.to_nat(n), 0, Array.get(U32, lefts, 0), Array.get(U32, rights, 0), Walked{Array.get(F32, ys, 0), Array.get(F32, xs, 0), ENil{}}, index, g, dir, k, lo, hi))
486
487def bounded(+index: Bool, +g: Graph, +dir: Direction, +k: Weight, +lo: F32, +hi: F32, +n: U32, ys: Array<F32>, xs: Array<F32>, m: Merged) -> Edges:
488 Merged{st, lefts, rights, +visits} = m
489 swept(index, g, dir, k, lo, hi, n, ys, xs, lefts, flush(st, n, rights, visits))
490
491def stacked(+index: Bool, +g: Graph, +dir: Direction, +k: Weight, +lo: F32, +hi: F32, +n: U32, xs: Array<F32>, l: Local) -> Edges:
492 Local{b, c} = l
493 Cols{ys, lefts, rights} = c
494 bounded(index, g, dir, k, lo, hi, n, ys, xs, merge_block(b, Merged{KNil{}, lefts, rights, 0.0}))
495
496def counts(+n: U32) -> Array<U32>:
497 Array.new(U32, Nat.add(1n, U32.log2(n)), 0)
498
499# A stand-in for `xs` when x is the index, which no sweep reads.
500def no_xs() -> Array<F32>:
501 Array.new(F32, 0n, 0.0)
502
503def built(+index: Bool, +g: Graph, +dir: Direction, +k: Weight, +lo: F32, +hi: F32, ys: Array<F32>, xs: Array<F32>, +n: U32) -> Edges:
504 stacked(index, g, dir, k, lo, hi, n, xs, block_pass(U32.to_nat(n), 0, Array.get(F32, ys, 0), counts(n), counts(n)))
505
506# The graph of the `n` values `ys` at `xs`, on one thread, as one block.
507# `lo` and `hi` bound the weights, as ts2vg's min_weight and max_weight;
508# pass neg_inf() and inf() for no bound. The x of `Index{}` is the index as
509# an F32, exact for up to 2^24 values; the x of `At{xs}` must increase.
510# src/parallel.bend builds the same graph on every thread, and relies on
511# Base's @unsafe Array.fork.
512def build(+g: Graph, +dir: Direction, +k: Weight, +lo: F32, +hi: F32, ys: Array<F32>, xs: Xs, +n: U32) -> Edges:
513 match xs:
514 case Index{}:
515 built(True{}, g, dir, k, lo, hi, ys, no_xs(), n)
516 case At{a}:
517 built(False{}, g, dir, k, lo, hi, ys, a, n)