xref: /freebsd-src/contrib/arm-optimized-routines/math/tools/asin.sollya (revision f3087bef11543b42e0d69b708f367097a4118d24)
1*f3087befSAndrew Turner// polynomial for approximating asin(x)
2*f3087befSAndrew Turner//
3*f3087befSAndrew Turner// Copyright (c) 2023-2024, Arm Limited.
4*f3087befSAndrew Turner// SPDX-License-Identifier: MIT OR Apache-2.0 WITH LLVM-exception
5*f3087befSAndrew Turner
6*f3087befSAndrew Turnerf = asin(x);
7*f3087befSAndrew Turnerdtype = double;
8*f3087befSAndrew Turner
9*f3087befSAndrew Turnerprec=256;
10*f3087befSAndrew Turner
11*f3087befSAndrew Turnera = 0x1p-106;
12*f3087befSAndrew Turnerb = 0.25;
13*f3087befSAndrew Turner
14*f3087befSAndrew Turnerdeg = 11;
15*f3087befSAndrew Turner
16*f3087befSAndrew Turnerbackward = proc(poly, d) {
17*f3087befSAndrew Turner  return d + d ^ 3 * poly(d * d);
18*f3087befSAndrew Turner};
19*f3087befSAndrew Turner
20*f3087befSAndrew Turnerforward = proc(f, d) {
21*f3087befSAndrew Turner  return (f(sqrt(d))-sqrt(d))/(d*sqrt(d));
22*f3087befSAndrew Turner};
23*f3087befSAndrew Turner
24*f3087befSAndrew Turnerpoly = fpminimax(forward(f, x), [|0,...,deg|], [|dtype ...|], [a;b], relative, floating);
25*f3087befSAndrew Turner
26*f3087befSAndrew Turnerdisplay = hexadecimal!;
27*f3087befSAndrew Turnerprint("rel error:", dirtyinfnorm(1-backward(poly, x)/f(x), [a;b]));
28*f3087befSAndrew Turnerprint("in [", a, b, "]");
29*f3087befSAndrew Turnerfor i from 0 to deg do print(coeff(poly, i));
30