|
|
1.1 root 1: /* $NetBSD: fpu_hyperb.c,v 1.17 2016/12/05 15:31:01 isaki Exp $ */
2:
3: /*
4: * Copyright (c) 1995 Ken Nakata
5: * All rights reserved.
6: *
7: * Redistribution and use in source and binary forms, with or without
8: * modification, are permitted provided that the following conditions
9: * are met:
10: * 1. Redistributions of source code must retain the above copyright
11: * notice, this list of conditions and the following disclaimer.
12: * 2. Redistributions in binary form must reproduce the above copyright
13: * notice, this list of conditions and the following disclaimer in the
14: * documentation and/or other materials provided with the distribution.
15: * 3. Neither the name of the author nor the names of its contributors
16: * may be used to endorse or promote products derived from this software
17: * without specific prior written permission.
18: *
19: * THIS SOFTWARE IS PROVIDED BY THE AUTHOR AND CONTRIBUTORS ``AS IS'' AND
20: * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
21: * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
22: * ARE DISCLAIMED. IN NO EVENT SHALL THE AUTHOR OR CONTRIBUTORS BE LIABLE
23: * FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
24: * DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS
25: * OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
26: * HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
27: * LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY
28: * OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
29: * SUCH DAMAGE.
30: *
31: * @(#)fpu_hyperb.c 10/24/95
32: */
33:
34: /*
35: * Copyright (c) 2011 Tetsuya Isaki. All rights reserved.
36: *
37: * Redistribution and use in source and binary forms, with or without
38: * modification, are permitted provided that the following conditions
39: * are met:
40: * 1. Redistributions of source code must retain the above copyright
41: * notice, this list of conditions and the following disclaimer.
42: * 2. Redistributions in binary form must reproduce the above copyright
43: * notice, this list of conditions and the following disclaimer in the
44: * documentation and/or other materials provided with the distribution.
45: *
46: * THIS SOFTWARE IS PROVIDED BY THE AUTHOR ``AS IS'' AND ANY EXPRESS OR
47: * IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES
48: * OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED.
49: * IN NO EVENT SHALL THE AUTHOR BE LIABLE FOR ANY DIRECT, INDIRECT,
50: * INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING,
51: * BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
52: * LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED
53: * AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY,
54: * OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY
55: * OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
56: * SUCH DAMAGE.
57: */
58:
59: #include "fpu_emulate.h"
60:
61: /*
62: * fpu_hyperb.c: defines the following functions
63: *
64: * fpu_atanh(), fpu_cosh(), fpu_sinh(), and fpu_tanh()
65: */
66:
67: /*
68: * 1 1 + x
69: * atanh(x) = ---*log(-------)
70: * 2 1 - x
71: */
72: struct fpn *
73: fpu_atanh(struct fpemu *fe)
74: {
75: struct fpn x;
76: struct fpn t;
77: struct fpn *r;
78:
79: if (ISNAN(&fe->fe_f2))
80: return &fe->fe_f2;
81: if (ISINF(&fe->fe_f2))
82: return fpu_newnan(fe);
83:
84: /* if x is +0/-0, return +0/-0 */
85: if (ISZERO(&fe->fe_f2))
86: return &fe->fe_f2;
87:
88: /*
89: * -INF if x == -1
90: * +INF if x == 1
91: */
92: r = &fe->fe_f2;
93: if (r->fp_exp == 0 && r->fp_mant[0] == FP_1 &&
94: r->fp_mant[1] == 0 && r->fp_mant[2] == 0) {
95: r->fp_class = FPC_INF;
96: return r;
97: }
98:
99: /* NAN if |x| > 1 */
100: if (fe->fe_f2.fp_exp >= 0)
101: return fpu_newnan(fe);
102:
103: CPYFPN(&x, &fe->fe_f2);
104:
105: /* t := 1 - x */
106: fpu_const(&fe->fe_f1, FPU_CONST_1);
107: fe->fe_f2.fp_sign = !fe->fe_f2.fp_sign;
108: r = fpu_add(fe);
109: CPYFPN(&t, r);
110:
111: /* r := 1 + x */
112: fpu_const(&fe->fe_f1, FPU_CONST_1);
113: CPYFPN(&fe->fe_f2, &x);
114: r = fpu_add(fe);
115:
116: /* (1-x)/(1+x) */
117: CPYFPN(&fe->fe_f1, r);
118: CPYFPN(&fe->fe_f2, &t);
119: r = fpu_div(fe);
120:
121: /* log((1-x)/(1+x)) */
122: CPYFPN(&fe->fe_f2, r);
123: r = fpu_logn(fe);
124:
125: /* r /= 2 */
126: r->fp_exp--;
127:
128: return r;
129: }
130:
131: /*
132: * exp(x) + exp(-x)
133: * cosh(x) = ------------------
134: * 2
135: */
136: struct fpn *
137: fpu_cosh(struct fpemu *fe)
138: {
139: struct fpn x, *fp;
140:
141: if (ISNAN(&fe->fe_f2))
142: return &fe->fe_f2;
143:
144: if (ISINF(&fe->fe_f2)) {
145: fe->fe_f2.fp_sign = 0;
146: return &fe->fe_f2;
147: }
148:
149: fp = fpu_etox(fe);
150: CPYFPN(&x, fp);
151:
152: fpu_const(&fe->fe_f1, FPU_CONST_1);
153: CPYFPN(&fe->fe_f2, fp);
154: fp = fpu_div(fe);
155:
156: CPYFPN(&fe->fe_f1, fp);
157: CPYFPN(&fe->fe_f2, &x);
158: fp = fpu_add(fe);
159:
160: fp->fp_exp--;
161:
162: return fp;
163: }
164:
165: /*
166: * exp(x) - exp(-x)
167: * sinh(x) = ------------------
168: * 2
169: */
170: struct fpn *
171: fpu_sinh(struct fpemu *fe)
172: {
173: struct fpn x, *fp;
174:
175: if (ISNAN(&fe->fe_f2))
176: return &fe->fe_f2;
177: if (ISINF(&fe->fe_f2))
178: return &fe->fe_f2;
179:
180: /* if x is +0/-0, return +0/-0 */
181: if (ISZERO(&fe->fe_f2))
182: return &fe->fe_f2;
183:
184: fp = fpu_etox(fe);
185: CPYFPN(&x, fp);
186:
187: fpu_const(&fe->fe_f1, FPU_CONST_1);
188: CPYFPN(&fe->fe_f2, fp);
189: fp = fpu_div(fe);
190:
191: fp->fp_sign = 1;
192: CPYFPN(&fe->fe_f1, fp);
193: CPYFPN(&fe->fe_f2, &x);
194: fp = fpu_add(fe);
195:
196: fp->fp_exp--;
197:
198: return fp;
199: }
200:
201: /*
202: * sinh(x)
203: * tanh(x) = ---------
204: * cosh(x)
205: */
206: struct fpn *
207: fpu_tanh(struct fpemu *fe)
208: {
209: struct fpn x;
210: struct fpn s;
211: struct fpn *r;
212: struct fpn epx;
213: struct fpn emx;
214:
215: if (ISNAN(&fe->fe_f2))
216: return &fe->fe_f2;
217:
218: /* if x is +0/-0, return +0/-0 */
219: if (ISZERO(&fe->fe_f2))
220: return &fe->fe_f2;
221:
222: CPYFPN(&x, &fe->fe_f2);
223:
224: if (ISINF(&fe->fe_f2)) {
225: fpu_const(&fe->fe_f2, FPU_CONST_1);
226: fe->fe_f2.fp_sign = x.fp_sign;
227: return &fe->fe_f2;
228: }
229: if (fe->fe_f2.fp_exp >= 8) {
230: fpu_const(&fe->fe_f2, FPU_CONST_1);
231: fe->fe_f2.fp_sign = x.fp_sign;
232: fe->fe_fpsr |= FPSR_INEX2; /* inexact */
233: return &fe->fe_f2;
234: }
235:
236: /* exp(x) */
237: CPYFPN(&fe->fe_f2, &x);
238: r = fpu_etox(fe);
239: CPYFPN(&epx, r);
240:
241: /* exp(-x) */
242: CPYFPN(&fe->fe_f2, &epx);
243: fpu_const(&fe->fe_f1, FPU_CONST_1);
244: r = fpu_div(fe);
245: CPYFPN(&emx, r);
246:
247: /* exp(x) - exp(-x) */
248: CPYFPN(&fe->fe_f1, &epx);
249: CPYFPN(&fe->fe_f2, &emx);
250: fe->fe_f2.fp_sign = !fe->fe_f2.fp_sign;
251: r = fpu_add(fe);
252: CPYFPN(&s, r);
253:
254: /* exp(x) + exp(-x) */
255: CPYFPN(&fe->fe_f1, &epx);
256: CPYFPN(&fe->fe_f2, &emx);
257: r = fpu_add(fe);
258:
259: CPYFPN(&fe->fe_f2, r);
260: CPYFPN(&fe->fe_f1, &s);
261: r = fpu_div(fe);
262:
263: return r;
264: }
This archive runs on limited infrastructure. Preserving old code on modern bandwidth. Automated agents are requested to crawl responsibly.