1// polynomial for approximating erfc(x)*exp(x*x)
2//
3// Copyright (c) 2022-2023, Arm Limited.
4// SPDX-License-Identifier: MIT OR Apache-2.0 WITH LLVM-exception
5
6deg = 12; // poly degree
7
8itv = parse(__argv[0]);
9
10bounds = [|3.725290298461914e-9,
11           0.18920711500272103,
12           0.41421356237309515,
13           0.681792830507429,
14           1,
15           1.378414230005442,
16           1.8284271247461903,
17           2.363585661014858,
18           3,
19           3.756828460010884,
20           4.656854249492381,
21           5.727171322029716,
22           7,
23           8.513656920021768,
24           10.313708498984761,
25           12.454342644059432,
26           15,
27           18.027313840043536,
28           21.627416997969522,
29           25.908685288118864,
30           31|];
31
32a = bounds[itv];
33b = bounds[itv + 1];
34
35f = proc(y) {
36  t = y + a;
37  return erfc(t) * exp(t*t);
38};
39
40poly = fpminimax(f(x), deg, [|double ...|], [0;b-a]);
41
42display = hexadecimal;
43print("rel error:", accurateinfnorm(1-poly(x)/f(x), [a;b], 30));
44print("in [",a,b,"]");
45print("coeffs:");
46for i from 0 to deg do coeff(poly, i);
47