util.bend

unified-nsga-iii@0.2.2.0 · 344 lines · 48 definitions

1# Shared helpers. Parallelize independent list splits; keep these total.
2# Population maps use mid-split + `a b = f(lo) f(hi)` (bend guide Parallelism).
3# Fuel is |xs| so the Nat shrinks; halves stay observationally ordered.
4
5import Base
6
7def list_len(-A: Data, xs: List<&2, A>) -> Nat:
8  List.length(&2, A, xs)
9
10# |take(xs ++ replicate(n, fill), n)| = n. Identity when |xs| = n
11# (take never reads the pad). Used to lock Run / select to the law size.
12def list_take_pad(-A: Data, +fill: A, xs: List<&2, A>, +n: Nat) -> List<&2, A>:
13  List.take(&2, A, List.append(&2, A, xs, List.replicate(A, n, fill)), n)
14
15# IEEE +0.0 is the all-false Word(32). Definitional on the 0.0 literal,
16# so p=0 operator laws skip F32.is_gt / F32.is_eq (opaque Base laws).
17def word_all_false(n: Nat, w: Word(n)) -> Bool:
18  match n:
19    case 0n:
20      True{}
21    case 1n+p:
22      match w:
23        case WCon{False{}, t}:
24          word_all_false(p, t)
25        case WCon{True{}, t}:
26          False{}
27
28def f32_is_zero(x: F32) -> Bool:
29  match x:
30    case F32{w}:
31      word_all_false(32n, w)
32
33def nat_min.fin(+a: Nat, +b: Nat, c: Bool) -> Nat:
34  match c:
35    case True{}:
36      a
37    case False{}:
38      b
39
40def nat_min(+a: Nat, +b: Nat) -> Nat:
41  nat_min.fin(a, b, Nat.is_le(a, b))
42
43def nat_max.fin(+a: Nat, +b: Nat, c: Bool) -> Nat:
44  match c:
45    case True{}:
46      a
47    case False{}:
48      b
49
50def nat_max(+a: Nat, +b: Nat) -> Nat:
51  nat_max.fin(a, b, Nat.is_ge(a, b))
52
53def nat_pow(+b: Nat, e: Nat) -> Nat:
54  match e:
55    case 0n:
56      1n
57    case 1n+p:
58      Nat.mul(b, nat_pow(b, p))
59
60# C(n, k). Das–Dennis count = C(p + M - 1, M - 1).
61def binomial(+n: Nat, +k: Nat) -> Nat:
62  match n k:
63    case x 0n:
64      1n
65    case 0n 1n+pk:
66      0n
67    case 1n+pn 1n+pk:
68      Nat.add(binomial(pn, pk), binomial(pn, 1n+pk))
69
70def in_nats(+x: Nat, xs: List<&2, Nat>) -> Bool:
71  List.contains(~Nat, ~Nat.is_eq, xs, x)
72
73def keep_if(-A: Data, +h: A, rest: List<&2, A>, keep: Bool) -> List<&2, A>:
74  match keep:
75    case True{}:
76      h <> rest
77    case False{}:
78      rest
79
80def minus.put(h: Nat, rest: List<&2, Nat>, drop: Bool) -> List<&2, Nat>:
81  match drop:
82    case True{}:
83      rest
84    case False{}:
85      h <> rest
86
87def minus(+xs: List<&2, Nat>, +ys: List<&2, Nat>) -> List<&2, Nat>:
88  match xs:
89    case Nil{}:
90      Nil{}
91    case h <> t:
92      minus.put(h, minus(t, ys), in_nats(h, ys))
93
94def get_nat(xs: List<&2, Nat>, i: Nat) -> Nat:
95  Maybe.default(&2, Nat, List.get(&2, Nat, xs, i), 0n)
96
97def inc_at(+xs: List<&2, Nat>, +i: Nat) -> List<&2, Nat>:
98  List.set(&2, Nat, xs, i, 1n+get_nat(xs, i))
99
100def zeros_nat(n: Nat) -> List<&2, Nat>:
101  List.replicate(Nat, n, 0n)
102
103# Sum of a Nat histogram. Used to lock mid-split count conservation.
104def nats_sum(xs: List<&2, Nat>) -> Nat:
105  match xs:
106    case Nil{}:
107      0n
108    case +x <> rest:
109      Nat.add(x, nats_sum(rest))
110
111# 1 when i is a valid bin of a length-n histogram, else 0.
112# List.set (and therefore inc_at) is a no-op at i >= n.
113def bit(b: Bool) -> Nat:
114  match b:
115    case True{}:
116      1n
117    case False{}:
118      0n
119
120def bin_hit(+i: Nat, +n: Nat) -> Nat:
121  bit(Nat.is_lt(i, n))
122
123def nats_add(a: List<&2, Nat>, b: List<&2, Nat>) -> List<&2, Nat>:
124  match a b:
125    case Nil{} Nil{}:
126      Nil{}
127    case x <> xs y <> ys:
128      Nat.add(x, y) <> nats_add(xs, ys)
129    case Nil{} ys:
130      ys
131    case xs Nil{}:
132      xs
133
134def vec_min(a: List<&2, F32>, b: List<&2, F32>) -> List<&2, F32>:
135  match a b:
136    case Nil{} Nil{}:
137      Nil{}
138    case x <> xs y <> ys:
139      F32.min(x, y) <> vec_min(xs, ys)
140    case Nil{} y <> ys:
141      y <> ys
142    case x <> xs Nil{}:
143      x <> xs
144
145def vec_max(a: List<&2, F32>, b: List<&2, F32>) -> List<&2, F32>:
146  match a b:
147    case Nil{} Nil{}:
148      Nil{}
149    case x <> xs y <> ys:
150      F32.max(x, y) <> vec_max(xs, ys)
151    case Nil{} y <> ys:
152      y <> ys
153    case x <> xs Nil{}:
154      x <> xs
155
156def vec_sub(a: List<&2, F32>, b: List<&2, F32>) -> List<&2, F32>:
157  match a b:
158    case Nil{} Nil{}:
159      Nil{}
160    case x <> xs y <> ys:
161      F32.sub(x, y) <> vec_sub(xs, ys)
162    case Nil{} y <> ys:
163      Nil{}
164    case x <> xs Nil{}:
165      Nil{}
166
167def vec_scale(+s: F32, v: List<&2, F32>) -> List<&2, F32>:
168  match v:
169    case Nil{}:
170      Nil{}
171    case x <> xs:
172      F32.mul(s, x) <> vec_scale(s, xs)
173
174def vec_dot(a: List<&2, F32>, b: List<&2, F32>) -> F32:
175  match a b:
176    case Nil{} Nil{}:
177      0.0
178    case x <> xs y <> ys:
179      F32.add(F32.mul(x, y), vec_dot(xs, ys))
180    case Nil{} y <> ys:
181      0.0
182    case x <> xs Nil{}:
183      0.0
184
185def vec_sqnorm(+v: List<&2, F32>) -> F32:
186  vec_dot(v, v)
187
188def ones_f32(n: Nat) -> List<&2, F32>:
189  List.replicate(F32, n, 1.0)
190
191def zeros_f32(n: Nat) -> List<&2, F32>:
192  List.replicate(F32, n, 0.0)
193
194def get_f32(xs: List<&2, F32>, i: Nat) -> F32:
195  Maybe.default(&2, F32, List.get(&2, F32, xs, i), 0.0)
196
197def set_f32(xs: List<&2, F32>, i: Nat, x: F32) -> List<&2, F32>:
198  List.set(&2, F32, xs, i, x)
199
200def f32_or_one(+span: F32, tiny: Bool) -> F32:
201  match tiny:
202    case True{}:
203      1.0
204    case False{}:
205      span
206
207# C# Normalize: denom < 1e-12 → 1.0 (intercept span uses the same floor).
208def denom(+lo: F32, +hi: F32) -> F32:
209  +span = F32.sub(hi, lo)
210  f32_or_one(span, F32.is_lt(span, 1.0e-12))
211
212def vec_norm_div(f: List<&2, F32>, ideal: List<&2, F32>, nadir: List<&2, F32>) -> List<&2, F32>:
213  match f ideal nadir:
214    case Nil{} los his:
215      Nil{}
216    case x <> xs +lo <> los +hi <> his:
217      F32.div(F32.sub(x, lo), denom(lo, hi)) <> vec_norm_div(xs, los, his)
218    case x <> xs Nil{} his:
219      Nil{}
220    case x <> xs lo <> los Nil{}:
221      Nil{}
222
223def csv_step(x: F32, acc: String) -> String:
224  match acc:
225    case SNil{}:
226      F32.show(x)
227    case SCon{h, t}:
228      F32.show(x) ++ "," ++ SCon{h, t}
229
230def show_f32s(xs: List<&2, F32>) -> String:
231  List.foldr(~&2, ~F32, ~String, ~csv_step, xs, "")
232
233def f32s_eq(a: List<&2, F32>, b: List<&2, F32>) -> Bool:
234  match a b:
235    case Nil{} Nil{}:
236      True{}
237    case x <> xs y <> ys:
238      Bool.and(F32.is_eq(x, y), f32s_eq(xs, ys))
239    case Nil{} y <> ys:
240      False{}
241    case x <> xs Nil{}:
242      False{}
243
244def has_f32s(+needle: List<&2, F32>, xs: List<&2, List<&2, F32>>) -> Bool:
245  match xs:
246    case Nil{}:
247      False{}
248    case h <> t:
249      Bool.or(f32s_eq(needle, h), has_f32s(needle, t))
250
251# Decision-variable duplicate key: 12 decimal places, round(x*1e12)/1e12.
252# C# DecisionKey is double.ToString("G12") (12 significant digits). Those
253# are not the same format. Near magnitude 1, F32 ULP is coarser than
254# 1e-12, so 12-dp rounding often does nothing. Below about 1e-5 a 12-dp
255# bin is wider than the F32 ULP and the round changes the value
256# (1.234567e-8 → 1.2346e-8). This tree keeps 12-dp keys. It does not
257# match C# G12. See src/g12_key.bend.
258def g12_scale() -> F32:
259  1000000000000.0
260
261def g12_round_abs(+x: F32) -> F32:
262  +s = g12_scale()
263  F32.div(F32.round(F32.mul(x, s)), s)
264
265def g12_one.signed(+x: F32, neg: Bool) -> F32:
266  match neg:
267    case True{}:
268      F32.neg(g12_round_abs(F32.neg(x)))
269    case False{}:
270      g12_round_abs(x)
271
272def g12_one(+x: F32) -> F32:
273  g12_one.signed(x, F32.is_lt(x, 0.0))
274
275# Balanced fork-join map. Singleton / empty stay sequential. Fuel 0 is a
276# total fallback (public callers pass |xs|). `~f` is a closed top-level.
277# Mid take/drop is inlined: a computed pair cannot be matched.
278def map_halves.go(
279  ~A: Data,
280  ~B: Data,
281  ~C: Data,
282  ~f: @+c:C -> @+x:A -> B,
283  +ctx: C,
284  +fuel: Nat,
285  xs: List<&2, A>
286) -> List<&2, B>:
287  match fuel xs:
288    case 0n Nil{}:
289      Nil{}
290    case 0n h <> t:
291      f(ctx, h) <> map_halves.go(~A, ~B, ~C, ~f, ctx, 0n, t)
292    case 1n+p Nil{}:
293      Nil{}
294    case 1n+p h <> Nil{}:
295      f(ctx, h) <> Nil{}
296    case 1n+p h <> t:
297      +all = {h <> t : List<&2, A>}
298      +mid = Nat.div(list_len(A, all), 2n)
299      a b = map_halves.go(~A, ~B, ~C, ~f, ctx, p, List.take(&2, A, all, mid)) map_halves.go(~A, ~B, ~C, ~f, ctx, p, List.drop(&2, A, all, mid))
300      List.append(&2, B, a, b)
301
302def map_halves(~A: Data, ~B: Data, ~C: Data, ~f: @+c:C -> @+x:A -> B, +ctx: C, +xs: List<&2, A>) -> List<&2, B>:
303  map_halves.go(~A, ~B, ~C, ~f, ctx, list_len(A, xs), xs)
304
305# Independent keep/drop. Mid-split; `take ++ drop` matches the old
306# tail-first cons_if walk (original relative order).
307def filter_halves.go(
308  ~A: Data,
309  ~C: Data,
310  ~p: @+c:C -> @+x:A -> Bool,
311  +ctx: C,
312  +fuel: Nat,
313  xs: List<&2, A>
314) -> List<&2, A>:
315  match fuel xs:
316    case 0n Nil{}:
317      Nil{}
318    case 0n +h <> t:
319      keep_if(A, h, filter_halves.go(~A, ~C, ~p, ctx, 0n, t), p(ctx, h))
320    case 1n+p0 Nil{}:
321      Nil{}
322    case 1n+p0 +h <> Nil{}:
323      keep_if(A, h, Nil{}, p(ctx, h))
324    case 1n+p0 h <> t:
325      +all = {h <> t : List<&2, A>}
326      +mid = Nat.div(list_len(A, all), 2n)
327      a b = filter_halves.go(~A, ~C, ~p, ctx, p0, List.take(&2, A, all, mid)) filter_halves.go(~A, ~C, ~p, ctx, p0, List.drop(&2, A, all, mid))
328      List.append(&2, A, a, b)
329
330def filter_halves(~A: Data, ~C: Data, ~p: @+c:C -> @+x:A -> Bool, +ctx: C, +xs: List<&2, A>) -> List<&2, A>:
331  filter_halves.go(~A, ~C, ~p, ctx, list_len(A, xs), xs)
332
333# Sequential cons walk (guide: the list is the shape a flat loop walks).
334# n_var is 12–30 on the oracles; mid-split take/drop of that spine is
335# more expensive than a tail walk. Same rounding as the old map_halves.
336def g12_vec.go(xs: List<&2, F32>, acc: List<&2, F32>) -> List<&2, F32>:
337  match xs:
338    case Nil{}:
339      List.reverse(&2, F32, acc)
340    case +x <> rest:
341      g12_vec.go(rest, g12_one(x) <> acc)
342
343def g12_vec(xs: List<&2, F32>) -> List<&2, F32>:
344  g12_vec.go(xs, Nil{})