xref: /illumos-gate/usr/src/lib/libm/i386/src/exp10.S (revision 55fea89dcaa64928bed4327112404dcb3e07b79f)
1/*
2 * CDDL HEADER START
3 *
4 * The contents of this file are subject to the terms of the
5 * Common Development and Distribution License (the "License").
6 * You may not use this file except in compliance with the License.
7 *
8 * You can obtain a copy of the license at usr/src/OPENSOLARIS.LICENSE
9 * or http://www.opensolaris.org/os/licensing.
10 * See the License for the specific language governing permissions
11 * and limitations under the License.
12 *
13 * When distributing Covered Code, include this CDDL HEADER in each
14 * file and include the License file at usr/src/OPENSOLARIS.LICENSE.
15 * If applicable, add the following below this CDDL HEADER, with the
16 * fields enclosed by brackets "[]" replaced with your own identifying
17 * information: Portions Copyright [yyyy] [name of copyright owner]
18 *
19 * CDDL HEADER END
20 */
21/*
22 * Copyright 2011 Nexenta Systems, Inc.  All rights reserved.
23 */
24/*
25 * Copyright 2005 Sun Microsystems, Inc.  All rights reserved.
26 * Use is subject to license terms.
27 */
28
29        .file "exp10.s"
30
31#include "libm.h"
32
33	ENTRY(exp10)
34	movl	8(%esp),%ecx		/ ecx <-- hi_32(x)
35	andl	$0x7fffffff,%ecx	/ ecx <-- hi_32(|x|)
36	cmpl	$0x3fd34413,%ecx	/ Is |x| < log10(2)?
37	jb	.shortcut		/ If so, take a shortcut.
38	je	.check_tail		/ maybe |x| only slightly < log10(2)
39	cmpl	$0x7ff00000,%ecx	/ hi_32(|x|) >= hi_32(INF)?
40	jae	.not_finite		/ if so, x is not finite
41.finite_non_special:			/ Here, log10(2) < |x| < INF
42	fldl	4(%esp)			/ push x (=arg)
43
44	subl	$8,%esp			/ save RP and set round-to-64-bits
45	fstcw	(%esp)
46	movw	(%esp),%ax
47	movw	%ax,4(%esp)
48	orw	$0x0300,%ax
49	movw	%ax,(%esp)
50	fldcw	(%esp)
51
52	fldl2t				/ push log2(10)  }NOT for xtndd_dbl
53	fmulp	%st,%st(1)		/ z = x*log2(10) }NOT for xtndd_dbl
54	fld	%st(0)			/ duplicate stack top
55	frndint				/ [z],z
56	fucom				/ z integral?
57	fstsw  %ax
58	sahf
59	je      .z_integral		/ branch if z integral
60	fxch				/ z, [z]
61	fsub	%st(1),%st		/ z-[z], [z]
62	f2xm1				/ 2**(z-[z])-1, [z]
63	fld1				/ 1,2**(z-[z])-1, [z]
64	faddp	%st,%st(1)		/ 2**(z-[z]), [z]
65	fscale				/ 2**z = 10**(arg), [z]
66	fstp	%st(1)
67
68	fstcw	(%esp)			/ restore old RP
69	movw	(%esp),%dx
70	andw	$0xfcff,%dx
71	movw	4(%esp),%cx
72	andw	$0x0300,%cx
73	orw	%dx,%cx
74	movw	%cx,(%esp)
75	fldcw	(%esp)
76	add	$8,%esp
77
78	ret
79
80.z_integral:				/ here, z is integral
81	fstp	%st(0)			/ ,z
82	fld1				/ 1 = 2**0, z
83	fscale				/ 2**(0 + z) = 2**z = 10**(arg), z
84	fstp	%st(1)			/ 10**(arg)
85
86	fstcw	(%esp)			/ restore old RP
87	movw	(%esp),%dx
88	andw	$0xfcff,%dx
89	movw	4(%esp),%cx
90	andw	$0x0300,%cx
91	orw	%dx,%cx
92	movw	%cx,(%esp)
93	fldcw	(%esp)
94	add	$8,%esp
95
96	ret
97
98.check_tail:
99	movl	4(%esp),%edx		/ edx <-- lo_32(x)
100	cmpl	$0x509f79fe,%edx	/ Is |x| slightly > log10(2)?
101	ja	.finite_non_special	/ branch if |x| slightly > log10(2)
102.shortcut:
103	/ Here, |x| < log10(2), so |z| = |x*log2(10)| < 1
104	/ whence z is in f2xm1's domain.
105	fldl	4(%esp)			/ push x (=arg)
106	fldl2t				/ push log2(10)  }NOT for xtndd_dbl
107	fmulp	%st,%st(1)		/ z = x*log2(10) }NOT for xtndd_dbl
108	f2xm1				/ 2**z - 1
109	fld1				/ 1,2**z - 1
110	faddp	%st,%st(1)		/ 2**z = 10**x
111	ret
112
113.not_finite:
114	cmpl	$0x7ff00000,%ecx	/ hi_32(|x|) > hi_32(INF)?
115	ja	.NaN_or_pinf		/ if so, x is NaN
116	movl	4(%esp),%edx		/ edx <-- lo_32(x)
117	cmpl	$0,%edx			/ lo_32(x) = 0?
118	jne	.NaN_or_pinf		/ if not, x is NaN
119	movl	8(%esp),%eax		/ eax <-- hi_32(x)
120	andl	$0x80000000,%eax	/ here, x is infinite, but +/-?
121	jz	.NaN_or_pinf		/ branch if x = +INF
122	fldz				/ Here, x = -inf, so return 0
123	ret
124
125.NaN_or_pinf:
126	/ Here, x = NaN or +inf, so load x and return immediately.
127	fldl	4(%esp)
128	fwait
129	ret
130	.align	4
131	SET_SIZE(exp10)
132