diff --git a/functional_algorithms/algorithms.py b/functional_algorithms/algorithms.py index 109e321..d85a81c 100644 --- a/functional_algorithms/algorithms.py +++ b/functional_algorithms/algorithms.py @@ -338,12 +338,6 @@ def real_asinh(ctx, x: float): This algorithm is based on the StableHLO v1.1.4 function CHLO_AsinhOp. - To avoid overflow in x * x, we use - - asinh(x) = log(2) + log(x) - - when abs(x) > sqrt(max), - To avoid underflow in 1 + x * x, we'll define z = hypot(1, x) and write @@ -360,6 +354,12 @@ def real_asinh(ctx, x: float): It turns out, this is accurate for all abs(x) < sqrt(max). + To avoid overflow in x ** 2, we'll use + + asinh(x) = log(2) + log(x) + + when abs(x) > sqrt(max), + For x < 0, we'll use asinh(x) = -asinh(-x) @@ -512,3 +512,58 @@ def acos(ctx, z: complex | float): if z.is_complex: return complex_acos(ctx, z) return real_acos(ctx, z) + + +def complex_acosh(ctx, z: complex): + """Inverse hyperbolic cosine on complex input: + + acosh(z) = sqrt(z - 1) / sqrt(1 - z) * acos(z) + = I * acos(z) # when z.imag >= 0 + = -I * acos(z) # otherwise + """ + w = complex_acos(ctx, z) + r = ctx.complex(-w.imag, w.real) + return ctx.select(z.imag < 0, -r, r) + + +def real_acosh(ctx, x: float): + """Inverse hyperbolic cosine on real input: + + acosh(x) = log(x + sqrt(x * x - 1)) + = log(x + sqrt(x+1)*sqrt(x-1))) + = log(1 + x-1 + sqrt(x+1)*sqrt(x-1))) + = log1p(sqrt(x-1) * (sqrt(x+1) + sqrt(x-1))) + + The last expression avoids errors from cancellations when x is + close to one. This also ensures the nan result when x < 1 because + sqrt(x') returns nan when x' < 0. + + To avoid overflow in multiplication for large x (x > max / 2), + we'll use + + acosh(x) = log(2) + log(x) + + """ + one = ctx.constant(1, x) + two = ctx.constant(2, x) + sqxm1 = ctx.sqrt(x - one) + sqxp1 = ctx.sqrt(x + one) + a0 = ctx.log(two) + ctx.log(x) + a1 = ctx.log1p(sqxm1 * (sqxp1 + sqxm1)) + + safe_max_limit_coefficient = ctx.parameters.get("safe_max_limit_coefficient", None) + if safe_max_limit_coefficient is None: + safe_max_limit = ctx.constant("largest", x) / 2 + else: + safe_max_limit = ctx.constant("largest", x) * safe_max_limit_coefficient + return ctx.select(x >= safe_max_limit, a0, a1) + + +def acosh(ctx, z: complex | float): + """Inverse hyperbolic cosine on complex and real inputs. + + See complex_acosh and real_acosh for more information. + """ + if z.is_complex: + return complex_acosh(ctx, z) + return real_acosh(ctx, z) diff --git a/functional_algorithms/context.py b/functional_algorithms/context.py index 2805290..3fa78c6 100644 --- a/functional_algorithms/context.py +++ b/functional_algorithms/context.py @@ -388,8 +388,8 @@ def ceil(self, x): def floor(self, x): return Expr(self, "floor", (x,)) - def copysign(self, x): - return Expr(self, "copysign", (x,)) + def copysign(self, x, y): + return Expr(self, "copysign", (x, y)) def round(self, x): return Expr(self, "round", (x,)) diff --git a/functional_algorithms/expr.py b/functional_algorithms/expr.py index dda35c1..b242a71 100644 --- a/functional_algorithms/expr.py +++ b/functional_algorithms/expr.py @@ -526,6 +526,9 @@ def real(self): def imag(self): return self.context.imag(self) + def conj(self): + return self.context.conj(self) + def __lt__(self, other): return self.context.lt(self, other) @@ -653,6 +656,8 @@ def get_type(self): "floor", "logical_not", "sign", + "copysign", + "conj", }: return self.operands[0].get_type() elif self.kind in { diff --git a/functional_algorithms/rewrite.py b/functional_algorithms/rewrite.py index 07b5d96..fd26eb1 100644 --- a/functional_algorithms/rewrite.py +++ b/functional_algorithms/rewrite.py @@ -225,6 +225,24 @@ def negative(self, expr): if x.kind == "negative": return x.operands[0] + def conj(self, expr): + + (x,) = expr.operands + + if x.kind == "constant": + value, like = x.operands + if isinstance(value, (int, float)): + return x + if isinstance(value, complex): + return x.context.constant(value.conjugate(), like) + + if x.kind == "complex": + real, imag = x.operands + return x.context.complex(real, -imag) + + if x.kind == "conj": + return x + def real(self, expr): (x,) = expr.operands diff --git a/functional_algorithms/targets/cpp.py b/functional_algorithms/targets/cpp.py index e557ea8..c6b51b1 100644 --- a/functional_algorithms/targets/cpp.py +++ b/functional_algorithms/targets/cpp.py @@ -19,6 +19,7 @@ absolute=[(":float32",), (":float64",), (":complex64",), (":complex128",)], hypot=[(":float32", ":float32"), (":float64", ":float64")], acos=[(":float32",), (":float64",), (":complex64",), (":complex128",)], + acosh=[(":float32",), (":float64",), (":complex64",), (":complex128",)], asin=[(":float32",), (":float64",), (":complex64",), (":complex128",)], asinh=[(":float32",), (":float64",), (":complex64",), (":complex128",)], # complex_asin=[(":complex", ":complex")], diff --git a/functional_algorithms/targets/numpy.py b/functional_algorithms/targets/numpy.py index 19faa86..aea07d7 100644 --- a/functional_algorithms/targets/numpy.py +++ b/functional_algorithms/targets/numpy.py @@ -27,6 +27,7 @@ def make_complex(r, i): trace_arguments = dict( absolute=[(":complex128",), (":complex64",)], acos=[(":complex128",), (":complex64",), (":float64",), (":float32",)], + acosh=[(":complex128",), (":complex64",), (":float64",), (":float32",)], asin=[(":complex128",), (":complex64",), (":float64",), (":float32",)], asinh=[(":complex128",), (":complex64",), (":float64",), (":float32",)], hypot=[(":float32", ":float32"), (":float64", ":float64")], @@ -86,7 +87,7 @@ def make_comment(message): log10="numpy.log10({0})", ceil="numpy.ceil({0})", floor="numpy.floor({0})", - copysign="numpy.copysign({0})", + copysign="numpy.copysign({0}, {1})", round=NotImplemented, sign="numpy.sign({0})", trunc="numpy.trunc({0})", @@ -111,7 +112,8 @@ def make_comment(message): largest="numpy.finfo({type}).max", posinf="{type}(numpy.inf)", neginf="-{type}(numpy.inf)", - pi="numpy.pi", + pi="{type}(numpy.pi)", + nan="{type}(numpy.nan)", ) type_to_target = dict( diff --git a/functional_algorithms/targets/python.py b/functional_algorithms/targets/python.py index 5ad92e8..284d514 100644 --- a/functional_algorithms/targets/python.py +++ b/functional_algorithms/targets/python.py @@ -18,6 +18,7 @@ trace_arguments = dict( absolute=[(":complex",)], acos=[(":complex",), (":float",)], + acosh=[(":complex",), (":float",)], asin=[(":complex",), (":float",)], asinh=[(":complex",), (":float",)], hypot=[(":float", ":float")], diff --git a/functional_algorithms/targets/stablehlo.py b/functional_algorithms/targets/stablehlo.py index 54d9cbc..a4cbe44 100644 --- a/functional_algorithms/targets/stablehlo.py +++ b/functional_algorithms/targets/stablehlo.py @@ -11,6 +11,7 @@ trace_arguments = dict( absolute=[(":complex",)], acos=[(":float",), (":complex",)], + acosh=[(":float",), (":complex",)], asin=[(":float",), (":complex",)], asinh=[(":float",), (":complex",)], hypot=[(":float", ":float")], @@ -55,10 +56,10 @@ def make_comment(message): bitwise_right_shift=NotImplemented, maximum="StableHLO_MaxOp", minimum="StableHLO_MinOp", - acos=NotImplemented, - acosh=NotImplemented, + acos="CHLO_AcosOp", + acosh="CHLO_AcoshOp", asin="CHLO_AsinOp", - asinh=NotImplemented, + asinh="CHLO_AsinhOp", atan=NotImplemented, atanh=NotImplemented, atan2="StableHLO_Atan2Op", diff --git a/functional_algorithms/targets/xla_client.py b/functional_algorithms/targets/xla_client.py index 1da83cd..44bc9aa 100644 --- a/functional_algorithms/targets/xla_client.py +++ b/functional_algorithms/targets/xla_client.py @@ -16,12 +16,15 @@ hypot=[(":float", ":float")], complex_acos=[(":complex", ":complex")], real_acos=[(":float", ":float")], + complex_acosh=[(":complex", ":complex")], + real_acosh=[(":float", ":float")], complex_asin=[(":complex", ":complex")], real_asin=[(":float", ":float")], complex_asinh=[(":complex", ":complex")], real_asinh=[(":float", ":float")], square=[(":float", ":float"), (":complex", ":complex")], acos=[(":float", ":float"), (":complex", ":complex")], + acosh=[(":float", ":float"), (":complex", ":complex")], asin=[(":float", ":float"), (":complex", ":complex")], asinh=[(":float", ":float"), (":complex", ":complex")], ) diff --git a/functional_algorithms/tests/test_algorithms.py b/functional_algorithms/tests/test_algorithms.py index d0b3941..52876ce 100644 --- a/functional_algorithms/tests/test_algorithms.py +++ b/functional_algorithms/tests/test_algorithms.py @@ -16,12 +16,12 @@ def dtype_name(request): return request.param -@pytest.fixture(scope="function", params=["absolute", "acos", "asin", "asinh", "hypot", "square"]) +@pytest.fixture(scope="function", params=["absolute", "acos", "acosh", "asin", "asinh", "hypot", "square"]) def func_name(request): return request.param -@pytest.fixture(scope="function", params=["absolute", "acos", "asin", "asinh", "square"]) +@pytest.fixture(scope="function", params=["absolute", "acos", "acosh", "asin", "asinh", "square"]) def unary_func_name(request): return request.param @@ -45,7 +45,7 @@ def test_unary(dtype_name, unary_func_name, flush_subnormals): func = targets.numpy.as_function(graph2, debug=0) - if unary_func_name in {"acos", "asin", "asinh"}: + if unary_func_name in {"acos", "asin", "asinh", "acosh"}: extra_prec_multiplier = 20 else: extra_prec_multiplier = 1 @@ -54,10 +54,13 @@ def test_unary(dtype_name, unary_func_name, flush_subnormals): unary_func_name, ) - size = 31 + size = 51 if dtype in {numpy.complex64, numpy.complex128}: samples = utils.complex_samples( - (size, size), dtype=dtype, include_huge=True, include_subnormal=not flush_subnormals + (size, size), + dtype=dtype, + include_huge=True, + include_subnormal=not flush_subnormals, ).flatten() else: samples = utils.real_samples(size * size, dtype=dtype, include_subnormal=not flush_subnormals).flatten() diff --git a/functional_algorithms/utils.py b/functional_algorithms/utils.py index 716c410..7054cb2 100644 --- a/functional_algorithms/utils.py +++ b/functional_algorithms/utils.py @@ -470,6 +470,33 @@ def arccos(self, x): return ctx.nan return ctx.acos(x) + def arccosh(self, x): + ctx = x.context + + if isinstance(x, ctx.mpc): + # Workaround mpmath 1.3 bug in acosh(+-inf+-infj) evaluation + # (see mpmath/mpmath#749). + pi = ctx.pi + inf = ctx.inf + zero = ctx.zero + if ctx.isinf(x.real): + sign_imag = -1 if x.imag < 0 else 1 + imag = ( + (3 if x.real < 0 else 1) * sign_imag * pi / 4 + if ctx.isinf(x.imag) + else (sign_imag * pi if x.real < 0 else zero) + ) + return ctx.make_mpc((inf._mpf_, imag._mpf_)) + elif ctx.isinf(x.imag): + sign_imag = -1 if x.imag < 0 else 1 + imag = sign_imag * pi / 2 + return ctx.make_mpc((inf._mpf_, imag._mpf_)) + else: + if x < 1: + # otherwise, mpmath.acosh would return complex value + return ctx.nan + return ctx.acosh(x) + class numpy_with_mpmath: """Namespace of universal functions on numpy arrays that use mpmath @@ -482,7 +509,7 @@ def __init__(self, **params): self.params = params def __getattr__(self, name): - name = dict(asinh="arcsinh", acos="arccos", asin="arcsin").get(name, name) + name = dict(asinh="arcsinh", acos="arccos", asin="arcsin", acosh="arccosh").get(name, name) if name in self._vfunc_cache: return self._vfunc_cache[name] if hasattr(mpmath_array_api, name): @@ -847,7 +874,7 @@ def validate_function( else: v1 = func(sample) v2 = reference_results[index][()] - assert v1.dtype == v2.dtype, (v1, v2) + assert v1.dtype == v2.dtype, (sample, v1, v2) ulp = diff_ulp(v1, v2, flush_subnormals=flush_subnormals) ulp_stats[ulp] += 1 if ulp > 2 and verbose: diff --git a/results/README.md b/results/README.md index 840934f..07e8e1a 100644 --- a/results/README.md +++ b/results/README.md @@ -16,6 +16,10 @@ MPMath functions using multi-precision arithmetic. | acos | float64 | 986273 | 12357 | 1366 | 5 | - | - | | acos | complex64 | 804838 | 196477 | 680 | 6 | - | - | | acos | complex128 | 701506 | 300255 | 238 | 2 | - | - | +| acosh | float32 | 989143 | 10829 | 29 | - | - | - | +| acosh | float64 | 947718 | 52275 | 8 | - | - | - | +| acosh | complex64 | 804838 | 196477 | 680 | 6 | - | - | +| acosh | complex128 | 701506 | 300255 | 238 | 2 | - | - | | asin | float32 | 937353 | 61812 | 810 | 26 | - | - | | asin | float64 | 983317 | 16662 | 22 | - | - | - | | asin | complex64 | 808487 | 191878 | 1592 | 44 | - | - | diff --git a/results/cpp/acosh.cpp b/results/cpp/acosh.cpp new file mode 100644 index 0000000..b7366fa --- /dev/null +++ b/results/cpp/acosh.cpp @@ -0,0 +1,210 @@ +// This file is generated using functional_algorithms tool (0.4.0), see +// https://github.com/pearu/functional_algorithms +// for more information. + + + +#include +#include +#include +#include +#include + + +float acosh_0(float z) { + float constant_2 = 2; + float constant_1 = 1; + float sqrt_subtract_z_constant_1 = std::sqrt((z) - (constant_1)); + return (((z) >= ((std::numeric_limits::max()) / (constant_2))) + ? ((std::log(constant_2)) + (std::log(z))) + : (std::log1p((sqrt_subtract_z_constant_1) * + ((std::sqrt((z) + (constant_1))) + + (sqrt_subtract_z_constant_1))))); +} + +double acosh_1(double z) { + double constant_2 = 2; + double constant_1 = 1; + double sqrt_subtract_z_constant_1 = std::sqrt((z) - (constant_1)); + return (((z) >= ((std::numeric_limits::max()) / (constant_2))) + ? ((std::log(constant_2)) + (std::log(z))) + : (std::log1p((sqrt_subtract_z_constant_1) * + ((std::sqrt((z) + (constant_1))) + + (sqrt_subtract_z_constant_1))))); +} + +std::complex acosh_2(std::complex z) { + float signed_y = (z).imag(); + float zero = 0; + float y = std::abs(signed_y); + float signed_x = (z).real(); + float x = std::abs(signed_x); + float safe_max = (std::sqrt(std::numeric_limits::max())) / (8); + float safe_max_opt = + (((x) < ((safe_max) * (1000000000000.0))) ? ((safe_max) * (1e-06)) + : ((safe_max) * (100.0))); + bool y_gt_safe_max_opt = (y) >= (safe_max_opt); + float mx = ((y_gt_safe_max_opt) ? (y) : (x)); + float two = 2; + float half = 0.5; + float xoy = (((y_gt_safe_max_opt) && + (!((y) == (std::numeric_limits::infinity())))) + ? ((x) / (y)) + : (zero)); + float one = 1; + bool logical_and_lt_y_safe_min_lt_x_one = + ((y) < ((std::sqrt(std::numeric_limits::min())) * (4))) && + ((x) < (one)); + float xp1 = (x) + (one); + float xm1 = (x) - (one); + float abs_xp1 = std::abs(xp1); + float _hypot_1_mx = std::max(abs_xp1, y); + float mn = std::min(abs_xp1, y); + float sqrt_two = std::sqrt(two); + float mn_over_mx = (mn) / (_hypot_1_mx); + float _hypot_1_r = (mn_over_mx) * (mn_over_mx); + float sqa = std::sqrt((one) + (_hypot_1_r)); + float r = + (((_hypot_1_mx) == (mn)) + ? ((sqrt_two) * (_hypot_1_mx)) + : (((((sqa) == (one)) && ((_hypot_1_r) > (zero))) + ? ((_hypot_1_mx) + (((_hypot_1_mx) * (_hypot_1_r)) / (two))) + : ((_hypot_1_mx) * (sqa))))); + float abs_xm1 = std::abs(xm1); + float _hypot_2_mx = std::max(abs_xm1, y); + float _hypot_2_mn = std::min(abs_xm1, y); + float _hypot_2_mn_over_mx = (_hypot_2_mn) / (_hypot_2_mx); + float _hypot_2_r = (_hypot_2_mn_over_mx) * (_hypot_2_mn_over_mx); + float _hypot_2_sqa = std::sqrt((one) + (_hypot_2_r)); + float s = + (((_hypot_2_mx) == (_hypot_2_mn)) + ? ((sqrt_two) * (_hypot_2_mx)) + : (((((_hypot_2_sqa) == (one)) && ((_hypot_2_r) > (zero))) + ? ((_hypot_2_mx) + (((_hypot_2_mx) * (_hypot_2_r)) / (two))) + : ((_hypot_2_mx) * (_hypot_2_sqa))))); + float a = (half) * ((r) + (s)); + float ap1 = (a) + (one); + float yy = (y) * (y); + float half_yy = (half) * (yy); + float rpxp1 = (r) + (xp1); + float divide_half_yy_rpxp1 = (half_yy) / (rpxp1); + float spxm1 = (s) + (xm1); + float smxm1 = (s) - (xm1); + float x_ge_1_or_not = + (((x) >= (one)) + ? ((divide_half_yy_rpxp1) + ((half) * (spxm1))) + : ((((a) <= (1.5)) ? ((divide_half_yy_rpxp1) + ((half_yy) / (smxm1))) + : ((a) - (one))))); + float am1 = + ((logical_and_lt_y_safe_min_lt_x_one) ? (-(((xp1) * (xm1)) / (ap1))) + : (x_ge_1_or_not)); + float sq = std::sqrt((am1) * (ap1)); + float imag = (((mx) >= (((y_gt_safe_max_opt) ? (safe_max_opt) : (safe_max)))) + ? (((std::log(two)) + (std::log(mx))) + + ((half) * (std::log1p((xoy) * (xoy))))) + : (((logical_and_lt_y_safe_min_lt_x_one) + ? ((y) / (sq)) + : (std::log1p((am1) + (sq)))))); + float half_apx = (half) * ((a) + (x)); + float acos_real = std::atan2( + (((std::max(x, y)) >= (safe_max)) + ? (y) + : ((((x) <= (one)) + ? (std::sqrt((half_apx) * (((yy) / (rpxp1)) + (smxm1)))) + : ((y) * (std::sqrt(((half_apx) / (rpxp1)) + + ((half_apx) / (spxm1)))))))), + signed_x); + std::complex complex_negative_acos_signed_imag_acos_real = + std::complex(-((((signed_y) < (zero)) ? (imag) : (-(imag)))), + acos_real); + return (((signed_y) < (0)) ? (-(complex_negative_acos_signed_imag_acos_real)) + : (complex_negative_acos_signed_imag_acos_real)); +} + +std::complex acosh_3(std::complex z) { + double signed_y = (z).imag(); + double zero = 0; + double y = std::abs(signed_y); + double signed_x = (z).real(); + double x = std::abs(signed_x); + double safe_max = (std::sqrt(std::numeric_limits::max())) / (8); + double safe_max_opt = + (((x) < ((safe_max) * (1000000000000.0))) ? ((safe_max) * (1e-06)) + : ((safe_max) * (100.0))); + bool y_gt_safe_max_opt = (y) >= (safe_max_opt); + double mx = ((y_gt_safe_max_opt) ? (y) : (x)); + double two = 2; + double half = 0.5; + double xoy = (((y_gt_safe_max_opt) && + (!((y) == (std::numeric_limits::infinity())))) + ? ((x) / (y)) + : (zero)); + double one = 1; + bool logical_and_lt_y_safe_min_lt_x_one = + ((y) < ((std::sqrt(std::numeric_limits::min())) * (4))) && + ((x) < (one)); + double xp1 = (x) + (one); + double xm1 = (x) - (one); + double abs_xp1 = std::abs(xp1); + double _hypot_1_mx = std::max(abs_xp1, y); + double mn = std::min(abs_xp1, y); + double sqrt_two = std::sqrt(two); + double mn_over_mx = (mn) / (_hypot_1_mx); + double _hypot_1_r = (mn_over_mx) * (mn_over_mx); + double sqa = std::sqrt((one) + (_hypot_1_r)); + double r = + (((_hypot_1_mx) == (mn)) + ? ((sqrt_two) * (_hypot_1_mx)) + : (((((sqa) == (one)) && ((_hypot_1_r) > (zero))) + ? ((_hypot_1_mx) + (((_hypot_1_mx) * (_hypot_1_r)) / (two))) + : ((_hypot_1_mx) * (sqa))))); + double abs_xm1 = std::abs(xm1); + double _hypot_2_mx = std::max(abs_xm1, y); + double _hypot_2_mn = std::min(abs_xm1, y); + double _hypot_2_mn_over_mx = (_hypot_2_mn) / (_hypot_2_mx); + double _hypot_2_r = (_hypot_2_mn_over_mx) * (_hypot_2_mn_over_mx); + double _hypot_2_sqa = std::sqrt((one) + (_hypot_2_r)); + double s = + (((_hypot_2_mx) == (_hypot_2_mn)) + ? ((sqrt_two) * (_hypot_2_mx)) + : (((((_hypot_2_sqa) == (one)) && ((_hypot_2_r) > (zero))) + ? ((_hypot_2_mx) + (((_hypot_2_mx) * (_hypot_2_r)) / (two))) + : ((_hypot_2_mx) * (_hypot_2_sqa))))); + double a = (half) * ((r) + (s)); + double ap1 = (a) + (one); + double yy = (y) * (y); + double half_yy = (half) * (yy); + double rpxp1 = (r) + (xp1); + double divide_half_yy_rpxp1 = (half_yy) / (rpxp1); + double spxm1 = (s) + (xm1); + double smxm1 = (s) - (xm1); + double x_ge_1_or_not = + (((x) >= (one)) + ? ((divide_half_yy_rpxp1) + ((half) * (spxm1))) + : ((((a) <= (1.5)) ? ((divide_half_yy_rpxp1) + ((half_yy) / (smxm1))) + : ((a) - (one))))); + double am1 = + ((logical_and_lt_y_safe_min_lt_x_one) ? (-(((xp1) * (xm1)) / (ap1))) + : (x_ge_1_or_not)); + double sq = std::sqrt((am1) * (ap1)); + double imag = (((mx) >= (((y_gt_safe_max_opt) ? (safe_max_opt) : (safe_max)))) + ? (((std::log(two)) + (std::log(mx))) + + ((half) * (std::log1p((xoy) * (xoy))))) + : (((logical_and_lt_y_safe_min_lt_x_one) + ? ((y) / (sq)) + : (std::log1p((am1) + (sq)))))); + double half_apx = (half) * ((a) + (x)); + double acos_real = std::atan2( + (((std::max(x, y)) >= (safe_max)) + ? (y) + : ((((x) <= (one)) + ? (std::sqrt((half_apx) * (((yy) / (rpxp1)) + (smxm1)))) + : ((y) * (std::sqrt(((half_apx) / (rpxp1)) + + ((half_apx) / (spxm1)))))))), + signed_x); + std::complex complex_negative_acos_signed_imag_acos_real = + std::complex(-((((signed_y) < (zero)) ? (imag) : (-(imag)))), + acos_real); + return (((signed_y) < (0)) ? (-(complex_negative_acos_signed_imag_acos_real)) + : (complex_negative_acos_signed_imag_acos_real)); +} \ No newline at end of file diff --git a/results/estimate_accuracy.py b/results/estimate_accuracy.py index 37accfa..9e28c2e 100644 --- a/results/estimate_accuracy.py +++ b/results/estimate_accuracy.py @@ -13,6 +13,12 @@ def get_inputs(): ("acos", np.float64, {}), ("acos", np.complex64, {}), ("acos", np.complex128, {}), + ("acosh", np.float32, {}), + # *(("acosh", np.float32, dict(safe_max_limit_coefficient=v)) for v in [0.51, 0.5]), + # *(("acosh", np.float64, dict(safe_max_limit_coefficient=v)) for v in [0.51, 0.5]), + ("acosh", np.float64, {}), + ("acosh", np.complex64, {}), + ("acosh", np.complex128, {}), ("asin", np.float32, {}), ("asin", np.float64, {}), ("asin", np.complex64, {}), @@ -79,7 +85,7 @@ def main(): graph2 = graph.implement_missing(fa.targets.numpy).simplify() func = fa.targets.numpy.as_function(graph2, debug=0) - if func_name in {"asin", "asinh", "acos"}: + if func_name in {"asin", "asinh", "acos", "acosh"}: extra_prec_multiplier = 20 else: extra_prec_multiplier = 1 @@ -126,6 +132,7 @@ def main(): else: cols.append(f"-") print(row_prefix + " | ".join(cols) + " |", file=f) + print(row_prefix + " | ".join(cols) + " |") f.flush() print(f"Created {target_file}") diff --git a/results/numpy/acosh.py b/results/numpy/acosh.py new file mode 100644 index 0000000..fc2a004 --- /dev/null +++ b/results/numpy/acosh.py @@ -0,0 +1,255 @@ +# This file is generated using functional_algorithms tool (0.4.0), see +# https://github.com/pearu/functional_algorithms +# for more information. + + +import numpy +import warnings + + +def make_complex(r, i): + if r.dtype == numpy.float32 and i.dtype == numpy.float32: + return numpy.array([r, i]).view(numpy.complex64)[0] + elif i.dtype == numpy.float64 and i.dtype == numpy.float64: + return numpy.array([r, i]).view(numpy.complex128)[0] + raise NotImplementedError((r.dtype, i.dtype)) + + +def acosh_0(z: numpy.complex128) -> numpy.complex128: + with warnings.catch_warnings(action="ignore"): + z = numpy.complex128(z) + signed_y: numpy.float64 = (z).imag + zero: numpy.float64 = numpy.float64(0) + y: numpy.float64 = numpy.abs(signed_y) + signed_x: numpy.float64 = (z).real + x: numpy.float64 = numpy.abs(signed_x) + safe_max: numpy.float64 = (numpy.sqrt(numpy.float64(numpy.finfo(numpy.float64).max))) / (numpy.float64(8)) + safe_max_opt: numpy.float64 = ( + ((safe_max) * (numpy.float64(1e-06))) + if ((x) < ((safe_max) * (numpy.float64(1000000000000.0)))) + else ((safe_max) * (numpy.float64(100.0))) + ) + y_gt_safe_max_opt: numpy.bool_ = (y) >= (safe_max_opt) + mx: numpy.float64 = (y) if (y_gt_safe_max_opt) else (x) + two: numpy.float64 = numpy.float64(2) + half: numpy.float64 = numpy.float64(0.5) + xoy: numpy.float64 = ( + ((x) / (y)) + if ((y_gt_safe_max_opt) and (not (numpy.equal(y, numpy.float64(numpy.float64(numpy.inf)), dtype=numpy.bool_)))) + else (zero) + ) + one: numpy.float64 = numpy.float64(1) + logical_and_lt_y_safe_min_lt_x_one: numpy.bool_ = ( + (y) < ((numpy.sqrt(numpy.float64(numpy.finfo(numpy.float64).tiny))) * (numpy.float64(4))) + ) and ((x) < (one)) + xp1: numpy.float64 = (x) + (one) + xm1: numpy.float64 = (x) - (one) + abs_xp1: numpy.float64 = numpy.abs(xp1) + _hypot_1_mx: numpy.float64 = max(abs_xp1, y) + mn: numpy.float64 = min(abs_xp1, y) + sqrt_two: numpy.float64 = numpy.sqrt(two) + mn_over_mx: numpy.float64 = (mn) / (_hypot_1_mx) + _hypot_1_r: numpy.float64 = (mn_over_mx) * (mn_over_mx) + sqa: numpy.float64 = numpy.sqrt((one) + (_hypot_1_r)) + r: numpy.float64 = ( + ((sqrt_two) * (_hypot_1_mx)) + if (numpy.equal(_hypot_1_mx, mn, dtype=numpy.bool_)) + else ( + ((_hypot_1_mx) + (((_hypot_1_mx) * (_hypot_1_r)) / (two))) + if ((numpy.equal(sqa, one, dtype=numpy.bool_)) and ((_hypot_1_r) > (zero))) + else ((_hypot_1_mx) * (sqa)) + ) + ) + abs_xm1: numpy.float64 = numpy.abs(xm1) + _hypot_2_mx: numpy.float64 = max(abs_xm1, y) + _hypot_2_mn: numpy.float64 = min(abs_xm1, y) + _hypot_2_mn_over_mx: numpy.float64 = (_hypot_2_mn) / (_hypot_2_mx) + _hypot_2_r: numpy.float64 = (_hypot_2_mn_over_mx) * (_hypot_2_mn_over_mx) + _hypot_2_sqa: numpy.float64 = numpy.sqrt((one) + (_hypot_2_r)) + s: numpy.float64 = ( + ((sqrt_two) * (_hypot_2_mx)) + if (numpy.equal(_hypot_2_mx, _hypot_2_mn, dtype=numpy.bool_)) + else ( + ((_hypot_2_mx) + (((_hypot_2_mx) * (_hypot_2_r)) / (two))) + if ((numpy.equal(_hypot_2_sqa, one, dtype=numpy.bool_)) and ((_hypot_2_r) > (zero))) + else ((_hypot_2_mx) * (_hypot_2_sqa)) + ) + ) + a: numpy.float64 = (half) * ((r) + (s)) + ap1: numpy.float64 = (a) + (one) + yy: numpy.float64 = (y) * (y) + half_yy: numpy.float64 = (half) * (yy) + rpxp1: numpy.float64 = (r) + (xp1) + divide_half_yy_rpxp1: numpy.float64 = (half_yy) / (rpxp1) + spxm1: numpy.float64 = (s) + (xm1) + smxm1: numpy.float64 = (s) - (xm1) + x_ge_1_or_not: numpy.float64 = ( + ((divide_half_yy_rpxp1) + ((half) * (spxm1))) + if ((x) >= (one)) + else (((divide_half_yy_rpxp1) + ((half_yy) / (smxm1))) if ((a) <= (numpy.float64(1.5))) else ((a) - (one))) + ) + am1: numpy.float64 = (-(((xp1) * (xm1)) / (ap1))) if (logical_and_lt_y_safe_min_lt_x_one) else (x_ge_1_or_not) + sq: numpy.float64 = numpy.sqrt((am1) * (ap1)) + imag: numpy.float64 = ( + (((numpy.log(two)) + (numpy.log(mx))) + ((half) * (numpy.log1p((xoy) * (xoy))))) + if ((mx) >= ((safe_max_opt) if (y_gt_safe_max_opt) else (safe_max))) + else (((y) / (sq)) if (logical_and_lt_y_safe_min_lt_x_one) else (numpy.log1p((am1) + (sq)))) + ) + half_apx: numpy.float64 = (half) * ((a) + (x)) + acos_real: numpy.float64 = numpy.arctan2( + ( + (y) + if ((max(x, y)) >= (safe_max)) + else ( + (numpy.sqrt((half_apx) * (((yy) / (rpxp1)) + (smxm1)))) + if ((x) <= (one)) + else ((y) * (numpy.sqrt(((half_apx) / (rpxp1)) + ((half_apx) / (spxm1))))) + ) + ), + signed_x, + ) + complex_negative_acos_signed_imag_acos_real: numpy.complex128 = make_complex( + -((imag) if ((signed_y) < (zero)) else (-(imag))), acos_real + ) + result = ( + (-(complex_negative_acos_signed_imag_acos_real)) + if ((signed_y) < (numpy.float64(0))) + else (complex_negative_acos_signed_imag_acos_real) + ) + return result + + +def acosh_1(z: numpy.complex64) -> numpy.complex64: + with warnings.catch_warnings(action="ignore"): + z = numpy.complex64(z) + signed_y: numpy.float32 = (z).imag + zero: numpy.float32 = numpy.float32(0) + y: numpy.float32 = numpy.abs(signed_y) + signed_x: numpy.float32 = (z).real + x: numpy.float32 = numpy.abs(signed_x) + safe_max: numpy.float32 = (numpy.sqrt(numpy.float32(numpy.finfo(numpy.float32).max))) / (numpy.float32(8)) + safe_max_opt: numpy.float32 = ( + ((safe_max) * (numpy.float32(1e-06))) + if ((x) < ((safe_max) * (numpy.float32(1000000000000.0)))) + else ((safe_max) * (numpy.float32(100.0))) + ) + y_gt_safe_max_opt: numpy.bool_ = (y) >= (safe_max_opt) + mx: numpy.float32 = (y) if (y_gt_safe_max_opt) else (x) + two: numpy.float32 = numpy.float32(2) + half: numpy.float32 = numpy.float32(0.5) + xoy: numpy.float32 = ( + ((x) / (y)) + if ((y_gt_safe_max_opt) and (not (numpy.equal(y, numpy.float32(numpy.float32(numpy.inf)), dtype=numpy.bool_)))) + else (zero) + ) + one: numpy.float32 = numpy.float32(1) + logical_and_lt_y_safe_min_lt_x_one: numpy.bool_ = ( + (y) < ((numpy.sqrt(numpy.float32(numpy.finfo(numpy.float32).tiny))) * (numpy.float32(4))) + ) and ((x) < (one)) + xp1: numpy.float32 = (x) + (one) + xm1: numpy.float32 = (x) - (one) + abs_xp1: numpy.float32 = numpy.abs(xp1) + _hypot_1_mx: numpy.float32 = max(abs_xp1, y) + mn: numpy.float32 = min(abs_xp1, y) + sqrt_two: numpy.float32 = numpy.sqrt(two) + mn_over_mx: numpy.float32 = (mn) / (_hypot_1_mx) + _hypot_1_r: numpy.float32 = (mn_over_mx) * (mn_over_mx) + sqa: numpy.float32 = numpy.sqrt((one) + (_hypot_1_r)) + r: numpy.float32 = ( + ((sqrt_two) * (_hypot_1_mx)) + if (numpy.equal(_hypot_1_mx, mn, dtype=numpy.bool_)) + else ( + ((_hypot_1_mx) + (((_hypot_1_mx) * (_hypot_1_r)) / (two))) + if ((numpy.equal(sqa, one, dtype=numpy.bool_)) and ((_hypot_1_r) > (zero))) + else ((_hypot_1_mx) * (sqa)) + ) + ) + abs_xm1: numpy.float32 = numpy.abs(xm1) + _hypot_2_mx: numpy.float32 = max(abs_xm1, y) + _hypot_2_mn: numpy.float32 = min(abs_xm1, y) + _hypot_2_mn_over_mx: numpy.float32 = (_hypot_2_mn) / (_hypot_2_mx) + _hypot_2_r: numpy.float32 = (_hypot_2_mn_over_mx) * (_hypot_2_mn_over_mx) + _hypot_2_sqa: numpy.float32 = numpy.sqrt((one) + (_hypot_2_r)) + s: numpy.float32 = ( + ((sqrt_two) * (_hypot_2_mx)) + if (numpy.equal(_hypot_2_mx, _hypot_2_mn, dtype=numpy.bool_)) + else ( + ((_hypot_2_mx) + (((_hypot_2_mx) * (_hypot_2_r)) / (two))) + if ((numpy.equal(_hypot_2_sqa, one, dtype=numpy.bool_)) and ((_hypot_2_r) > (zero))) + else ((_hypot_2_mx) * (_hypot_2_sqa)) + ) + ) + a: numpy.float32 = (half) * ((r) + (s)) + ap1: numpy.float32 = (a) + (one) + yy: numpy.float32 = (y) * (y) + half_yy: numpy.float32 = (half) * (yy) + rpxp1: numpy.float32 = (r) + (xp1) + divide_half_yy_rpxp1: numpy.float32 = (half_yy) / (rpxp1) + spxm1: numpy.float32 = (s) + (xm1) + smxm1: numpy.float32 = (s) - (xm1) + x_ge_1_or_not: numpy.float32 = ( + ((divide_half_yy_rpxp1) + ((half) * (spxm1))) + if ((x) >= (one)) + else (((divide_half_yy_rpxp1) + ((half_yy) / (smxm1))) if ((a) <= (numpy.float32(1.5))) else ((a) - (one))) + ) + am1: numpy.float32 = (-(((xp1) * (xm1)) / (ap1))) if (logical_and_lt_y_safe_min_lt_x_one) else (x_ge_1_or_not) + sq: numpy.float32 = numpy.sqrt((am1) * (ap1)) + imag: numpy.float32 = ( + (((numpy.log(two)) + (numpy.log(mx))) + ((half) * (numpy.log1p((xoy) * (xoy))))) + if ((mx) >= ((safe_max_opt) if (y_gt_safe_max_opt) else (safe_max))) + else (((y) / (sq)) if (logical_and_lt_y_safe_min_lt_x_one) else (numpy.log1p((am1) + (sq)))) + ) + half_apx: numpy.float32 = (half) * ((a) + (x)) + acos_real: numpy.float32 = numpy.arctan2( + ( + (y) + if ((max(x, y)) >= (safe_max)) + else ( + (numpy.sqrt((half_apx) * (((yy) / (rpxp1)) + (smxm1)))) + if ((x) <= (one)) + else ((y) * (numpy.sqrt(((half_apx) / (rpxp1)) + ((half_apx) / (spxm1))))) + ) + ), + signed_x, + ) + complex_negative_acos_signed_imag_acos_real: numpy.complex64 = make_complex( + -((imag) if ((signed_y) < (zero)) else (-(imag))), acos_real + ) + result = ( + (-(complex_negative_acos_signed_imag_acos_real)) + if ((signed_y) < (numpy.float32(0))) + else (complex_negative_acos_signed_imag_acos_real) + ) + return result + + +def acosh_2(z: numpy.float64) -> numpy.float64: + with warnings.catch_warnings(action="ignore"): + z = numpy.float64(z) + constant_2: numpy.float64 = numpy.float64(2) + constant_1: numpy.float64 = numpy.float64(1) + sqrt_subtract_z_constant_1: numpy.float64 = numpy.sqrt((z) - (constant_1)) + result = ( + ((numpy.log(constant_2)) + (numpy.log(z))) + if ((z) >= ((numpy.float64(numpy.finfo(numpy.float64).max)) / (constant_2))) + else ( + numpy.log1p((sqrt_subtract_z_constant_1) * ((numpy.sqrt((z) + (constant_1))) + (sqrt_subtract_z_constant_1))) + ) + ) + return result + + +def acosh_3(z: numpy.float32) -> numpy.float32: + with warnings.catch_warnings(action="ignore"): + z = numpy.float32(z) + constant_2: numpy.float32 = numpy.float32(2) + constant_1: numpy.float32 = numpy.float32(1) + sqrt_subtract_z_constant_1: numpy.float32 = numpy.sqrt((z) - (constant_1)) + result = ( + ((numpy.log(constant_2)) + (numpy.log(z))) + if ((z) >= ((numpy.float32(numpy.finfo(numpy.float32).max)) / (constant_2))) + else ( + numpy.log1p((sqrt_subtract_z_constant_1) * ((numpy.sqrt((z) + (constant_1))) + (sqrt_subtract_z_constant_1))) + ) + ) + return result diff --git a/results/python/acosh.py b/results/python/acosh.py new file mode 100644 index 0000000..a93186b --- /dev/null +++ b/results/python/acosh.py @@ -0,0 +1,109 @@ +# This file is generated using functional_algorithms tool (0.4.0), see +# https://github.com/pearu/functional_algorithms +# for more information. + + +import math +import sys + + +def acosh_0(z: complex) -> complex: + signed_y: float = (z).imag + zero: float = 0 + y: float = abs(signed_y) + signed_x: float = (z).real + x: float = abs(signed_x) + safe_max: float = (math.sqrt(sys.float_info.max)) / (8) + safe_max_opt: float = ((safe_max) * (1e-06)) if ((x) < ((safe_max) * (1000000000000.0))) else ((safe_max) * (100.0)) + y_gt_safe_max_opt: bool = (y) >= (safe_max_opt) + mx: float = (y) if (y_gt_safe_max_opt) else (x) + two: float = 2 + half: float = 0.5 + xoy: float = ((x) / (y)) if ((y_gt_safe_max_opt) and (not ((y) == (math.inf)))) else (zero) + one: float = 1 + logical_and_lt_y_safe_min_lt_x_one: bool = ((y) < ((math.sqrt(sys.float_info.min)) * (4))) and ((x) < (one)) + xp1: float = (x) + (one) + xm1: float = (x) - (one) + abs_xp1: float = abs(xp1) + _hypot_1_mx: float = max(abs_xp1, y) + mn: float = min(abs_xp1, y) + sqrt_two: float = math.sqrt(two) + mn_over_mx: float = (mn) / (_hypot_1_mx) + _hypot_1_r: float = (mn_over_mx) * (mn_over_mx) + sqa: float = math.sqrt((one) + (_hypot_1_r)) + r: float = ( + ((sqrt_two) * (_hypot_1_mx)) + if ((_hypot_1_mx) == (mn)) + else ( + ((_hypot_1_mx) + (((_hypot_1_mx) * (_hypot_1_r)) / (two))) + if (((sqa) == (one)) and ((_hypot_1_r) > (zero))) + else ((_hypot_1_mx) * (sqa)) + ) + ) + abs_xm1: float = abs(xm1) + _hypot_2_mx: float = max(abs_xm1, y) + _hypot_2_mn: float = min(abs_xm1, y) + _hypot_2_mn_over_mx: float = (_hypot_2_mn) / (_hypot_2_mx) + _hypot_2_r: float = (_hypot_2_mn_over_mx) * (_hypot_2_mn_over_mx) + _hypot_2_sqa: float = math.sqrt((one) + (_hypot_2_r)) + s: float = ( + ((sqrt_two) * (_hypot_2_mx)) + if ((_hypot_2_mx) == (_hypot_2_mn)) + else ( + ((_hypot_2_mx) + (((_hypot_2_mx) * (_hypot_2_r)) / (two))) + if (((_hypot_2_sqa) == (one)) and ((_hypot_2_r) > (zero))) + else ((_hypot_2_mx) * (_hypot_2_sqa)) + ) + ) + a: float = (half) * ((r) + (s)) + ap1: float = (a) + (one) + yy: float = (y) * (y) + half_yy: float = (half) * (yy) + rpxp1: float = (r) + (xp1) + divide_half_yy_rpxp1: float = (half_yy) / (rpxp1) + spxm1: float = (s) + (xm1) + smxm1: float = (s) - (xm1) + x_ge_1_or_not: float = ( + ((divide_half_yy_rpxp1) + ((half) * (spxm1))) + if ((x) >= (one)) + else (((divide_half_yy_rpxp1) + ((half_yy) / (smxm1))) if ((a) <= (1.5)) else ((a) - (one))) + ) + am1: float = (-(((xp1) * (xm1)) / (ap1))) if (logical_and_lt_y_safe_min_lt_x_one) else (x_ge_1_or_not) + sq: float = math.sqrt((am1) * (ap1)) + imag: float = ( + (((math.log(two)) + (math.log(mx))) + ((half) * (math.log1p((xoy) * (xoy))))) + if ((mx) >= ((safe_max_opt) if (y_gt_safe_max_opt) else (safe_max))) + else (((y) / (sq)) if (logical_and_lt_y_safe_min_lt_x_one) else (math.log1p((am1) + (sq)))) + ) + half_apx: float = (half) * ((a) + (x)) + acos_real: float = math.atan2( + ( + (y) + if ((max(x, y)) >= (safe_max)) + else ( + (math.sqrt((half_apx) * (((yy) / (rpxp1)) + (smxm1)))) + if ((x) <= (one)) + else ((y) * (math.sqrt(((half_apx) / (rpxp1)) + ((half_apx) / (spxm1))))) + ) + ), + signed_x, + ) + complex_negative_acos_signed_imag_acos_real: complex = complex( + -((imag) if ((signed_y) < (zero)) else (-(imag))), acos_real + ) + return ( + (-(complex_negative_acos_signed_imag_acos_real)) + if ((signed_y) < (0)) + else (complex_negative_acos_signed_imag_acos_real) + ) + + +def acosh_1(z: float) -> float: + constant_2: float = 2 + constant_1: float = 1 + sqrt_subtract_z_constant_1: float = math.sqrt((z) - (constant_1)) + return ( + ((math.log(constant_2)) + (math.log(z))) + if ((z) >= ((sys.float_info.max) / (constant_2))) + else (math.log1p((sqrt_subtract_z_constant_1) * ((math.sqrt((z) + (constant_1))) + (sqrt_subtract_z_constant_1)))) + ) diff --git a/results/stablehlo/acosh.td b/results/stablehlo/acosh.td new file mode 100644 index 0000000..91c543c --- /dev/null +++ b/results/stablehlo/acosh.td @@ -0,0 +1,254 @@ +// This file is generated using functional_algorithms tool (0.4.0), see +// https://github.com/pearu/functional_algorithms +// for more information. + + + + +def : Pat<(acosh_0 NonComplexElementType:$z), + (StableHLO_SelectOp + (StableHLO_CompareOp + $z, + (StableHLO_DivOp + (StableHLO_ConstantLikeMaxFiniteValue $z), + (StableHLO_ConstantLike<"2">:$constant_2 $z)), + StableHLO_ComparisonDirectionValue<"GE">, + (STABLEHLO_DEFAULT_COMPARISON_TYPE)), + (StableHLO_AddOp + (StableHLO_LogOp $constant_2), + (StableHLO_LogOp $z)), + (StableHLO_Log1pOp + (StableHLO_MulOp + (StableHLO_SqrtOp:$sqrt_subtract_z_constant_1 + (StableHLO_SubtractOp + $z, + (StableHLO_ConstantLike<"1">:$constant_1 $z))), + (StableHLO_AddOp + (StableHLO_SqrtOp + (StableHLO_AddOp $z, $constant_1)), + $sqrt_subtract_z_constant_1))))>; + +def : Pat<(acosh_1 ComplexElementType:$z), + (StableHLO_SelectOp + (StableHLO_CompareOp + (StableHLO_ImagOp:$signed_y $z), + (StableHLO_ConstantLike<"0"> $signed_y), + StableHLO_ComparisonDirectionValue<"LT">, + (STABLEHLO_DEFAULT_COMPARISON_TYPE)), + (StableHLO_NegOp + (StableHLO_ComplexOp:$complex_negative_acos_signed_imag_acos_real + (StableHLO_NegOp + (StableHLO_SelectOp + (StableHLO_CompareOp + $signed_y, + (StableHLO_ConstantLike<"0">:$zero (StableHLO_RealOp:$signed_x $z)), + StableHLO_ComparisonDirectionValue<"LT">, + (STABLEHLO_DEFAULT_COMPARISON_TYPE)), + (StableHLO_SelectOp:$imag + (StableHLO_CompareOp + (StableHLO_SelectOp:$mx + (StableHLO_CompareOp:$y_gt_safe_max_opt + (StableHLO_AbsOp:$y $signed_y), + (StableHLO_SelectOp:$safe_max_opt + (StableHLO_CompareOp + (StableHLO_AbsOp:$x $signed_x), + (StableHLO_MulOp + (StableHLO_DivOp:$safe_max + (StableHLO_SqrtOp + (StableHLO_ConstantLikeMaxFiniteValue $signed_x)), + (StableHLO_ConstantLike<"8"> $signed_x)), + (StableHLO_ConstantLike<"1000000000000.0"> $signed_x)), + StableHLO_ComparisonDirectionValue<"LT">, + (STABLEHLO_DEFAULT_COMPARISON_TYPE)), + (StableHLO_MulOp + $safe_max, + (StableHLO_ConstantLike<"1e-06"> $signed_x)), + (StableHLO_MulOp + $safe_max, + (StableHLO_ConstantLike<"100.0"> $signed_x))), + StableHLO_ComparisonDirectionValue<"GE">, + (STABLEHLO_DEFAULT_COMPARISON_TYPE)), + $y, + $x), + (StableHLO_SelectOp $y_gt_safe_max_opt, $safe_max_opt, $safe_max), + StableHLO_ComparisonDirectionValue<"GE">, + (STABLEHLO_DEFAULT_COMPARISON_TYPE)), + (StableHLO_AddOp + (StableHLO_AddOp + (StableHLO_LogOp + (StableHLO_ConstantLike<"2">:$two $signed_x)), + (StableHLO_LogOp $mx)), + (StableHLO_MulOp + (StableHLO_ConstantLike<"0.5">:$half $signed_x), + (StableHLO_Log1pOp + (StableHLO_MulOp + (StableHLO_SelectOp:$xoy + (StableHLO_AndOp + $y_gt_safe_max_opt, + (StableHLO_NotOp + (StableHLO_CompareOp + $y, + (StableHLO_ConstantLikePosInfValue $signed_y), + StableHLO_ComparisonDirectionValue<"EQ">, + (STABLEHLO_DEFAULT_COMPARISON_TYPE)))), + (StableHLO_DivOp $x, $y), + $zero), + $xoy)))), + (StableHLO_SelectOp + (StableHLO_AndOp:$logical_and_lt_y_safe_min_lt_x_one + (StableHLO_CompareOp + $y, + (StableHLO_MulOp + (StableHLO_SqrtOp + (StableHLO_ConstantLikeSmallestNormalizedValue $signed_x)), + (StableHLO_ConstantLike<"4"> $signed_x)), + StableHLO_ComparisonDirectionValue<"LT">, + (STABLEHLO_DEFAULT_COMPARISON_TYPE)), + (StableHLO_CompareOp + $x, + (StableHLO_ConstantLike<"1">:$one $signed_x), + StableHLO_ComparisonDirectionValue<"LT">, + (STABLEHLO_DEFAULT_COMPARISON_TYPE))), + (StableHLO_DivOp + $y, + (StableHLO_SqrtOp:$sq + (StableHLO_MulOp + (StableHLO_SelectOp:$am1 + $logical_and_lt_y_safe_min_lt_x_one, + (StableHLO_NegOp + (StableHLO_DivOp + (StableHLO_MulOp + (StableHLO_AddOp:$xp1 $x, $one), + (StableHLO_SubtractOp:$xm1 $x, $one)), + (StableHLO_AddOp:$ap1 + (StableHLO_MulOp:$a + $half, + (StableHLO_AddOp + (StableHLO_SelectOp:$r + (StableHLO_CompareOp + (StableHLO_MaxOp:$_hypot_1_mx + (StableHLO_AbsOp:$abs_xp1 $xp1), + $y), + (StableHLO_MinOp:$mn $abs_xp1, $y), + StableHLO_ComparisonDirectionValue<"EQ">, + (STABLEHLO_DEFAULT_COMPARISON_TYPE)), + (StableHLO_MulOp + (StableHLO_SqrtOp:$sqrt_two $two), + $_hypot_1_mx), + (StableHLO_SelectOp + (StableHLO_AndOp + (StableHLO_CompareOp + (StableHLO_SqrtOp:$sqa + (StableHLO_AddOp + $one, + (StableHLO_MulOp:$_hypot_1_r + (StableHLO_DivOp:$mn_over_mx $mn, $_hypot_1_mx), + $mn_over_mx))), + $one, + StableHLO_ComparisonDirectionValue<"EQ">, + (STABLEHLO_DEFAULT_COMPARISON_TYPE)), + (StableHLO_CompareOp + $_hypot_1_r, + $zero, + StableHLO_ComparisonDirectionValue<"GT">, + (STABLEHLO_DEFAULT_COMPARISON_TYPE))), + (StableHLO_AddOp + $_hypot_1_mx, + (StableHLO_DivOp + (StableHLO_MulOp $_hypot_1_mx, $_hypot_1_r), + $two)), + (StableHLO_MulOp $_hypot_1_mx, $sqa))), + (StableHLO_SelectOp:$s + (StableHLO_CompareOp + (StableHLO_MaxOp:$_hypot_2_mx + (StableHLO_AbsOp:$abs_xm1 $xm1), + $y), + (StableHLO_MinOp:$_hypot_2_mn $abs_xm1, $y), + StableHLO_ComparisonDirectionValue<"EQ">, + (STABLEHLO_DEFAULT_COMPARISON_TYPE)), + (StableHLO_MulOp $sqrt_two, $_hypot_2_mx), + (StableHLO_SelectOp + (StableHLO_AndOp + (StableHLO_CompareOp + (StableHLO_SqrtOp:$_hypot_2_sqa + (StableHLO_AddOp + $one, + (StableHLO_MulOp:$_hypot_2_r + (StableHLO_DivOp:$_hypot_2_mn_over_mx $_hypot_2_mn, $_hypot_2_mx), + $_hypot_2_mn_over_mx))), + $one, + StableHLO_ComparisonDirectionValue<"EQ">, + (STABLEHLO_DEFAULT_COMPARISON_TYPE)), + (StableHLO_CompareOp + $_hypot_2_r, + $zero, + StableHLO_ComparisonDirectionValue<"GT">, + (STABLEHLO_DEFAULT_COMPARISON_TYPE))), + (StableHLO_AddOp + $_hypot_2_mx, + (StableHLO_DivOp + (StableHLO_MulOp $_hypot_2_mx, $_hypot_2_r), + $two)), + (StableHLO_MulOp $_hypot_2_mx, $_hypot_2_sqa))))), + $one))), + (StableHLO_SelectOp:$x_ge_1_or_not + (StableHLO_CompareOp + $x, + $one, + StableHLO_ComparisonDirectionValue<"GE">, + (STABLEHLO_DEFAULT_COMPARISON_TYPE)), + (StableHLO_AddOp + (StableHLO_DivOp:$divide_half_yy_rpxp1 + (StableHLO_MulOp:$half_yy + $half, + (StableHLO_MulOp:$yy $y, $y)), + (StableHLO_AddOp:$rpxp1 $r, $xp1)), + (StableHLO_MulOp + $half, + (StableHLO_AddOp:$spxm1 $s, $xm1))), + (StableHLO_SelectOp + (StableHLO_CompareOp + $a, + (StableHLO_ConstantLike<"1.5"> $signed_x), + StableHLO_ComparisonDirectionValue<"LE">, + (STABLEHLO_DEFAULT_COMPARISON_TYPE)), + (StableHLO_AddOp + $divide_half_yy_rpxp1, + (StableHLO_DivOp + $half_yy, + (StableHLO_SubtractOp:$smxm1 $s, $xm1))), + (StableHLO_SubtractOp $a, $one)))), + $ap1))), + (StableHLO_Log1pOp + (StableHLO_AddOp $am1, $sq)))), + (StableHLO_NegOp $imag))), + (StableHLO_Atan2Op:$acos_real + (StableHLO_SelectOp + (StableHLO_CompareOp + (StableHLO_MaxOp $x, $y), + $safe_max, + StableHLO_ComparisonDirectionValue<"GE">, + (STABLEHLO_DEFAULT_COMPARISON_TYPE)), + $y, + (StableHLO_SelectOp + (StableHLO_CompareOp + $x, + $one, + StableHLO_ComparisonDirectionValue<"LE">, + (STABLEHLO_DEFAULT_COMPARISON_TYPE)), + (StableHLO_SqrtOp + (StableHLO_MulOp + (StableHLO_MulOp:$half_apx + $half, + (StableHLO_AddOp $a, $x)), + (StableHLO_AddOp + (StableHLO_DivOp $yy, $rpxp1), + $smxm1))), + (StableHLO_MulOp + $y, + (StableHLO_SqrtOp + (StableHLO_AddOp + (StableHLO_DivOp $half_apx, $rpxp1), + (StableHLO_DivOp $half_apx, $spxm1)))))), + $signed_x))), + $complex_negative_acos_signed_imag_acos_real)>; \ No newline at end of file diff --git a/results/xla_client/acosh.cc b/results/xla_client/acosh.cc new file mode 100644 index 0000000..5d55690 --- /dev/null +++ b/results/xla_client/acosh.cc @@ -0,0 +1,111 @@ +// This file is generated using functional_algorithms tool (0.4.0), see +// https://github.com/pearu/functional_algorithms +// for more information. + + + +#include "xla/client/lib/constants.h" +#include "xla/client/xla_builder.h" +#include + + +template +XlaOp acosh_0(XlaOp z) { + FloatType constant_2 = 2; + XlaOp constant_constant_1 = ScalarLike(z, 1); + XlaOp sqrt_subtract_z_constant_constant_1 = Sqrt(Sub(z, constant_constant_1)); + return Select(Ge(z, ScalarLike(z, (std::numeric_limits::max()) / + (constant_2))), + Add(ScalarLike(z, std::log(constant_2)), Log(z)), + Log1p(Mul(sqrt_subtract_z_constant_constant_1, + Add(Sqrt(Add(z, constant_constant_1)), + sqrt_subtract_z_constant_constant_1)))); +} + +template +XlaOp acosh_1(XlaOp z) { + XlaOp signed_y = Imag(z); + FloatType zero_ = 0; + XlaOp zero = ScalarLike(signed_x, zero_); + XlaOp y = Abs(signed_y); + XlaOp signed_x = Real(z); + XlaOp x = Abs(signed_x); + FloatType safe_max_ = + (std::sqrt(std::numeric_limits::max())) / (8); + XlaOp safe_max_opt = + Select(Lt(x, ScalarLike(signed_x, (safe_max_) * (1000000000000.0))), + ScalarLike(signed_x, (safe_max_) * (1e-06)), + ScalarLike(signed_x, (safe_max_) * (100.0))); + XlaOp y_gt_safe_max_opt = Ge(y, safe_max_opt); + XlaOp mx = Select(y_gt_safe_max_opt, y, x); + XlaOp safe_max = ScalarLike(signed_x, safe_max_); + FloatType two_ = 2; + XlaOp half = ScalarLike(signed_x, 0.5); + XlaOp xoy = Select( + And(y_gt_safe_max_opt, + Not(Eq(y, ScalarLike(signed_y, + std::numeric_limits::infinity())))), + Div(x, y), zero); + XlaOp one = ScalarLike(signed_x, 1); + XlaOp logical_and_lt_y_safe_min_lt_x_one = And( + Lt(y, + ScalarLike(signed_x, + (std::sqrt(std::numeric_limits::min())) * (4))), + Lt(x, one)); + XlaOp xp1 = Add(x, one); + XlaOp xm1 = Sub(x, one); + XlaOp abs_xp1 = Abs(xp1); + XlaOp _hypot_1_mx = Max(abs_xp1, y); + XlaOp mn = Min(abs_xp1, y); + XlaOp sqrt_two = ScalarLike(signed_x, std::sqrt(two_)); + XlaOp _hypot_1_r = Square(Div(mn, _hypot_1_mx)); + XlaOp sqa = Sqrt(Add(one, _hypot_1_r)); + XlaOp two = ScalarLike(signed_x, two_); + XlaOp r = + Select(Eq(_hypot_1_mx, mn), Mul(sqrt_two, _hypot_1_mx), + Select(And(Eq(sqa, one), Gt(_hypot_1_r, zero)), + Add(_hypot_1_mx, Div(Mul(_hypot_1_mx, _hypot_1_r), two)), + Mul(_hypot_1_mx, sqa))); + XlaOp abs_xm1 = Abs(xm1); + XlaOp _hypot_2_mx = Max(abs_xm1, y); + XlaOp _hypot_2_mn = Min(abs_xm1, y); + XlaOp _hypot_2_r = Square(Div(_hypot_2_mn, _hypot_2_mx)); + XlaOp _hypot_2_sqa = Sqrt(Add(one, _hypot_2_r)); + XlaOp s = + Select(Eq(_hypot_2_mx, _hypot_2_mn), Mul(sqrt_two, _hypot_2_mx), + Select(And(Eq(_hypot_2_sqa, one), Gt(_hypot_2_r, zero)), + Add(_hypot_2_mx, Div(Mul(_hypot_2_mx, _hypot_2_r), two)), + Mul(_hypot_2_mx, _hypot_2_sqa))); + XlaOp a = Mul(half, Add(r, s)); + XlaOp ap1 = Add(a, one); + XlaOp yy = Mul(y, y); + XlaOp half_yy = Mul(half, yy); + XlaOp rpxp1 = Add(r, xp1); + XlaOp divide_half_yy_rpxp1 = Div(half_yy, rpxp1); + XlaOp spxm1 = Add(s, xm1); + XlaOp smxm1 = Sub(s, xm1); + XlaOp x_ge_1_or_not = Select( + Ge(x, one), Add(divide_half_yy_rpxp1, Mul(half, spxm1)), + Select(Le(a, ScalarLike(signed_x, 1.5)), + Add(divide_half_yy_rpxp1, Div(half_yy, smxm1)), Sub(a, one))); + XlaOp am1 = Select(logical_and_lt_y_safe_min_lt_x_one, + Neg(Div(Mul(xp1, xm1), ap1)), x_ge_1_or_not); + XlaOp sq = Sqrt(Mul(am1, ap1)); + XlaOp imag = Select(Ge(mx, Select(y_gt_safe_max_opt, safe_max_opt, safe_max)), + Add(Add(ScalarLike(signed_x, std::log(two_)), Log(mx)), + Mul(half, Log1p(Mul(xoy, xoy)))), + Select(logical_and_lt_y_safe_min_lt_x_one, Div(y, sq), + Log1p(Add(am1, sq)))); + XlaOp half_apx = Mul(half, Add(a, x)); + XlaOp acos_real = Atan2( + Select(Ge(Max(x, y), safe_max), y, + Select(Le(x, one), Sqrt(Mul(half_apx, Add(Div(yy, rpxp1), smxm1))), + Mul(y, Sqrt(Add(Div(half_apx, rpxp1), + Div(half_apx, spxm1)))))), + signed_x); + XlaOp complex_negative_acos_signed_imag_acos_real = + Complex(Neg(Select(Lt(signed_y, zero), imag, Neg(imag))), acos_real); + return Select(Lt(signed_y, ScalarLike(signed_y, zero_)), + Neg(complex_negative_acos_signed_imag_acos_real), + complex_negative_acos_signed_imag_acos_real); +} \ No newline at end of file diff --git a/results/xla_client/complex_acosh.cc b/results/xla_client/complex_acosh.cc new file mode 100644 index 0000000..e099950 --- /dev/null +++ b/results/xla_client/complex_acosh.cc @@ -0,0 +1,98 @@ +// This file is generated using functional_algorithms tool (0.4.0), see +// https://github.com/pearu/functional_algorithms +// for more information. + + + +#include "xla/client/lib/constants.h" +#include "xla/client/xla_builder.h" +#include + + +template +XlaOp complex_acosh_0(XlaOp z) { + XlaOp signed_y = Imag(z); + FloatType zero_ = 0; + XlaOp zero = ScalarLike(signed_x, zero_); + XlaOp y = Abs(signed_y); + XlaOp signed_x = Real(z); + XlaOp x = Abs(signed_x); + FloatType safe_max_ = + (std::sqrt(std::numeric_limits::max())) / (8); + XlaOp safe_max_opt = + Select(Lt(x, ScalarLike(signed_x, (safe_max_) * (1000000000000.0))), + ScalarLike(signed_x, (safe_max_) * (1e-06)), + ScalarLike(signed_x, (safe_max_) * (100.0))); + XlaOp y_gt_safe_max_opt = Ge(y, safe_max_opt); + XlaOp mx = Select(y_gt_safe_max_opt, y, x); + XlaOp safe_max = ScalarLike(signed_x, safe_max_); + FloatType two_ = 2; + XlaOp half = ScalarLike(signed_x, 0.5); + XlaOp xoy = Select( + And(y_gt_safe_max_opt, + Not(Eq(y, ScalarLike(signed_y, + std::numeric_limits::infinity())))), + Div(x, y), zero); + XlaOp one = ScalarLike(signed_x, 1); + XlaOp logical_and_lt_y_safe_min_lt_x_one = And( + Lt(y, + ScalarLike(signed_x, + (std::sqrt(std::numeric_limits::min())) * (4))), + Lt(x, one)); + XlaOp xp1 = Add(x, one); + XlaOp xm1 = Sub(x, one); + XlaOp abs_xp1 = Abs(xp1); + XlaOp _hypot_1_mx = Max(abs_xp1, y); + XlaOp mn = Min(abs_xp1, y); + XlaOp sqrt_two = ScalarLike(signed_x, std::sqrt(two_)); + XlaOp _hypot_1_r = Square(Div(mn, _hypot_1_mx)); + XlaOp sqa = Sqrt(Add(one, _hypot_1_r)); + XlaOp two = ScalarLike(signed_x, two_); + XlaOp r = + Select(Eq(_hypot_1_mx, mn), Mul(sqrt_two, _hypot_1_mx), + Select(And(Eq(sqa, one), Gt(_hypot_1_r, zero)), + Add(_hypot_1_mx, Div(Mul(_hypot_1_mx, _hypot_1_r), two)), + Mul(_hypot_1_mx, sqa))); + XlaOp abs_xm1 = Abs(xm1); + XlaOp _hypot_2_mx = Max(abs_xm1, y); + XlaOp _hypot_2_mn = Min(abs_xm1, y); + XlaOp _hypot_2_r = Square(Div(_hypot_2_mn, _hypot_2_mx)); + XlaOp _hypot_2_sqa = Sqrt(Add(one, _hypot_2_r)); + XlaOp s = + Select(Eq(_hypot_2_mx, _hypot_2_mn), Mul(sqrt_two, _hypot_2_mx), + Select(And(Eq(_hypot_2_sqa, one), Gt(_hypot_2_r, zero)), + Add(_hypot_2_mx, Div(Mul(_hypot_2_mx, _hypot_2_r), two)), + Mul(_hypot_2_mx, _hypot_2_sqa))); + XlaOp a = Mul(half, Add(r, s)); + XlaOp ap1 = Add(a, one); + XlaOp yy = Mul(y, y); + XlaOp half_yy = Mul(half, yy); + XlaOp rpxp1 = Add(r, xp1); + XlaOp divide_half_yy_rpxp1 = Div(half_yy, rpxp1); + XlaOp spxm1 = Add(s, xm1); + XlaOp smxm1 = Sub(s, xm1); + XlaOp x_ge_1_or_not = Select( + Ge(x, one), Add(divide_half_yy_rpxp1, Mul(half, spxm1)), + Select(Le(a, ScalarLike(signed_x, 1.5)), + Add(divide_half_yy_rpxp1, Div(half_yy, smxm1)), Sub(a, one))); + XlaOp am1 = Select(logical_and_lt_y_safe_min_lt_x_one, + Neg(Div(Mul(xp1, xm1), ap1)), x_ge_1_or_not); + XlaOp sq = Sqrt(Mul(am1, ap1)); + XlaOp imag = Select(Ge(mx, Select(y_gt_safe_max_opt, safe_max_opt, safe_max)), + Add(Add(ScalarLike(signed_x, std::log(two_)), Log(mx)), + Mul(half, Log1p(Mul(xoy, xoy)))), + Select(logical_and_lt_y_safe_min_lt_x_one, Div(y, sq), + Log1p(Add(am1, sq)))); + XlaOp half_apx = Mul(half, Add(a, x)); + XlaOp acos_real = Atan2( + Select(Ge(Max(x, y), safe_max), y, + Select(Le(x, one), Sqrt(Mul(half_apx, Add(Div(yy, rpxp1), smxm1))), + Mul(y, Sqrt(Add(Div(half_apx, rpxp1), + Div(half_apx, spxm1)))))), + signed_x); + XlaOp complex_negative_acos_signed_imag_acos_real = + Complex(Neg(Select(Lt(signed_y, zero), imag, Neg(imag))), acos_real); + return Select(Lt(signed_y, ScalarLike(signed_y, zero_)), + Neg(complex_negative_acos_signed_imag_acos_real), + complex_negative_acos_signed_imag_acos_real); +} \ No newline at end of file diff --git a/results/xla_client/real_acosh.cc b/results/xla_client/real_acosh.cc new file mode 100644 index 0000000..27fa6c6 --- /dev/null +++ b/results/xla_client/real_acosh.cc @@ -0,0 +1,23 @@ +// This file is generated using functional_algorithms tool (0.4.0), see +// https://github.com/pearu/functional_algorithms +// for more information. + + + +#include "xla/client/lib/constants.h" +#include "xla/client/xla_builder.h" +#include + + +template +XlaOp real_acosh_0(XlaOp x) { + FloatType constant_2 = 2; + XlaOp constant_constant_1 = ScalarLike(x, 1); + XlaOp sqrt_subtract_x_constant_constant_1 = Sqrt(Sub(x, constant_constant_1)); + return Select(Ge(x, ScalarLike(x, (std::numeric_limits::max()) / + (constant_2))), + Add(ScalarLike(x, std::log(constant_2)), Log(x)), + Log1p(Mul(sqrt_subtract_x_constant_constant_1, + Add(Sqrt(Add(x, constant_constant_1)), + sqrt_subtract_x_constant_constant_1)))); +} \ No newline at end of file