1 /*-
2  * SPDX-License-Identifier: BSD-2-Clause
4  * Copyright (c) 2012 Stephen Montgomery-Smith <stephen@FreeBSD.ORG>
32  * Fairgrieve, and Ping Tak Peter Tang, published in ACM Transactions on
33  * Mathematical Software, Volume 23 Issue 3, 1997, Pages 299-335,
38  * XXX comments were removed automatically, and even short ones on the right
60 FOUR_SQRT_MIN =		0x1p-61,
62 m_e =			2.7182818285e0,		/*  0xadf854.0p-22 */
63 m_ln2 =			6.9314718056e-1,	/*  0xb17218.0p-24 */
64 pio2_hi =		1.5707962513e0,		/*  0xc90fda.0p-23 */
66 SQRT_3_EPSILON =	5.9801995673e-4,	/*  0x9cc471.0p-34 */
67 SQRT_6_EPSILON =	8.4572793338e-4,	/*  0xddb3d7.0p-34 */
68 SQRT_MIN =		0x1p-63;
71 pio2_lo =		7.5497899549e-8,	/*  0xa22169.0p-47 */
72 tiny =			0x1p-100;
80 		return ((hypot_a_b - b) / 2);  in f()
94 	S = hypotf(x, y - 1);  in do_hard_work()
103 		} else if (x >= FLT_EPSILON * fabsf(y - 1)) {  in do_hard_work()
104 			Am1 = f(x, 1 + y, R) + f(x, 1 - y, S);  in do_hard_work()
107 			*rx = x / sqrtf((1 - y) * (1 + y));  in do_hard_work()
109 			*rx = log1pf((y - 1) + sqrtf((y - 1) * (y + 1)));  in do_hard_work()
112 		*rx = logf(A + sqrtf(A * A - 1));  in do_hard_work()
131 		} else if (x >= FLT_EPSILON * fabsf(y - 1)) {  in do_hard_work()
132 			Amy = f(x, y + 1, R) + f(x, y - 1, S);  in do_hard_work()
136 			    sqrtf((y + 1) * (y - 1));  in do_hard_work()
139 			*sqrt_A2my2 = sqrtf((1 - y) * (1 + y));  in do_hard_work()
170 			w = clog_for_large_values(-z) + m_ln2;  in casinhf()
216 			return (CMPLXF(y + y, -INFINITY));  in cacosf()
218 			return (CMPLXF(x + x, -y));  in cacosf()
229 			ry = -ry;  in cacosf()
234 		return (CMPLXF(0, -y));  in cacosf()
239 		return (CMPLXF(pio2_hi - (x - pio2_lo), -y));  in cacosf()
246 			rx = acosf(-B);  in cacosf()
251 			rx = atan2f(sqrt_A2mx2, -new_x);  in cacosf()
254 		ry = -ry;  in cacosf()
323 #define	BIAS	(FLT_MAX_EXP - 1)  in real_part_reciprocal()
325 	if (ix - iy >= CUTOFF << 23 || isinf(x))  in real_part_reciprocal()
327 	if (iy - ix >= CUTOFF << 23)  in real_part_reciprocal()
329 	if (ix <= (BIAS + FLT_MAX_EXP / 2 - CUTOFF) << 23)  in real_part_reciprocal()
331 	SET_FLOAT_WORD(scale, 0x7f800000 - ix);  in real_part_reciprocal()
372 		rx = (m_ln2 - logf(ay)) / 2;  in catanhf()
374 		rx = log1pf(4 * ax / sum_squares(ax - 1, ay)) / 4;  in catanhf()
377 		ry = atan2f(2, -ay) / 2;  in catanhf()
379 		ry = atan2f(2 * ay, (1 - ax) * (1 + ax)) / 2;  in catanhf()
381 		ry = atan2f(2 * ay, (1 - ax) * (1 + ax) - ay * ay) / 2;  in catanhf()