1/* ix87 specific implementation of arcsinh.
2   Copyright (C) 1996-2022 Free Software Foundation, Inc.
3   This file is part of the GNU C Library.
4
5   The GNU C Library is free software; you can redistribute it and/or
6   modify it under the terms of the GNU Lesser General Public
7   License as published by the Free Software Foundation; either
8   version 2.1 of the License, or (at your option) any later version.
9
10   The GNU C Library is distributed in the hope that it will be useful,
11   but WITHOUT ANY WARRANTY; without even the implied warranty of
12   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
13   Lesser General Public License for more details.
14
15   You should have received a copy of the GNU Lesser General Public
16   License along with the GNU C Library; if not, see
17   <https://www.gnu.org/licenses/>.  */
18
19#include <libm-alias-ldouble.h>
20#include <machine/asm.h>
21
22	.section .rodata
23
24	.align ALIGNARG(4)
25	.type huge,@object
26huge:	.tfloat 1e+4930
27	ASM_SIZE_DIRECTIVE(huge)
28	.align ALIGNARG(4)
29	/* Please note that we use double value for 1.0.  This number
30	   has an exact representation and so we don't get accuracy
31	   problems.  The advantage is that the code is simpler.  */
32	.type one,@object
33one:	.double 1.0
34	ASM_SIZE_DIRECTIVE(one)
35	/* It is not important that this constant is precise.  It is only
36	   a value which is known to be on the safe side for using the
37	   fyl2xp1 instruction.  */
38	.type limit,@object
39limit:	.double 0.29
40	ASM_SIZE_DIRECTIVE(limit)
41
42#ifdef PIC
43#define MO(op) op##@GOTOFF(%edx)
44#else
45#define MO(op) op
46#endif
47
48	.text
49ENTRY(__asinhl)
50	movl	12(%esp), %ecx
51	movl	$0x7fff, %eax
52	andl	%ecx, %eax
53	andl	$0x8000, %ecx
54	movl	%eax, %edx
55	orl	$0xffff8000, %edx
56	incl	%edx
57	jz	7f			// x in �Inf or NaN
58	xorl	%ecx, 12(%esp)
59	fldt	4(%esp)			// |x|
60	cmpl	$0x3fde, %eax
61	jb	2f			// |x| < 2^-34
62	fldln2				// log(2) : |x|
63	cmpl	$0x4020, %eax
64	fxch				// |x| : log(2)
65	ja	3f			// |x| > 2^34
66#ifdef	PIC
67	LOAD_PIC_REG (dx)
68#endif
69	cmpl	$0x4000, %eax
70	ja	5f			// |x| > 2
71
72	// 2^-34 <= |x| <= 2 => y = sign(x)*log1p(|x|+|x|^2/(1+sqrt(1+|x|^2)))
73	fld	%st			// |x| : |x| : log(2)
74	fmul	%st(1)			// |x|^2 : |x| : log(2)
75	fld	%st			// |x|^2 : |x|^2 : |x| : log(2)
76	faddl	MO(one)			// 1+|x|^2 : |x|^2 : |x| : log(2)
77	fsqrt				// sqrt(1+|x|^2) : |x|^2 : |x| : log(2)
78	faddl	MO(one)			// 1+sqrt(1+|x|^2) : |x|^2 : |x| : log(2)
79	fdivrp				// |x|^2/(1+sqrt(1+|x|^2)) : |x| : log(2)
80	faddp				// |x|+|x|^2/(1+sqrt(1+|x|^2)) : log(2)
81	fcoml	MO(limit)
82	fnstsw
83	sahf
84	ja	6f
85	fyl2xp1
86	jecxz	4f
87	fchs
884:	ret
89
907:	fldt	4(%esp)
91	fadd	%st
92	ret
93
946:	faddl	MO(one)
95	fyl2x
96	jecxz	4f
97	fchs
984:	ret
99
100	// |x| < 2^-34 => y = x (inexact iff |x| != 0.0)
101	.align ALIGNARG(4)
1022:
103#ifdef	PIC
104	LOAD_PIC_REG (dx)
105#endif
106	jecxz	4f
107	fchs				// x
1084:	fld	%st			// x : x
109	fldt	MO(huge)		// huge : x : x
110	faddp				// huge+x : x
111	fstp	%st(0)			// x
112	cmpl	$0x0001, %eax
113	jae	8f
114	fld	%st(0)
115	fmul	%st(0)
116	fstp	%st(0)
1178:	ret
118
119	// |x| > 2^34 => y = sign(x) * (log(|x|) + log(2))
120	.align ALIGNARG(4)
1213:	fyl2x				// log(|x|)
122	fldln2				// log(2) : log(|x|)
123	faddp				// log(|x|)+log(2)
124	jecxz	4f
125	fchs
1264:	ret
127
128	// |x| > 2 => y = sign(x) * log(2*|x| + 1/(|x|+sqrt(x*x+1)))
129	.align ALIGNARG(4)
1305:	fld	%st			// |x| : |x| : log(2)
131	fadd	%st, %st(1)		// |x| : 2*|x| : log(2)
132	fld	%st			// |x| : |x| : 2*|x| : log(2)
133	fmul	%st(1)			// |x|^2 : |x| : 2*|x| : log(2)
134	faddl	MO(one)			// 1+|x|^2 : |x| : 2*|x| : log(2)
135	fsqrt				// sqrt(1+|x|^2) : |x| : 2*|x| : log(2)
136	faddp				// |x|+sqrt(1+|x|^2) : 2*|x| : log(2)
137	fdivrl	MO(one)			// 1/(|x|+sqrt(1+|x|^2)) : 2*|x| : log(2)
138	faddp				// 2*|x|+1/(|x|+sqrt(1+|x|^2)) : log(2)
139	fyl2x				// log(2*|x|+1/(|x|+sqrt(1+|x|^2)))
140	jecxz	4f
141	fchs
1424:	ret
143END(__asinhl)
144libm_alias_ldouble (__asinh, asinh)
145