Home | History | Annotate | Line # | Download | only in fpe
fpu_exp.c revision 1.7 
      1 /*	$NetBSD: fpu_exp.c,v 1.7 2013/04/20 04:38:51 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_exp.c	10/24/95
     32  */
     33 
     34 #include <sys/cdefs.h>
     35 __KERNEL_RCSID(0, "$NetBSD: fpu_exp.c,v 1.7 2013/04/20 04:38:51 isaki Exp $");
     36 
     37 #include "fpu_emulate.h"
     38 
     39 /* The number of items to terminate the Taylor expansion */
     40 #define MAX_ITEMS	(2000)
     41 
     42 /*
     43  * fpu_exp.c: defines fpu_etox(), fpu_etoxm1(), fpu_tentox(), and fpu_twotox();
     44  */
     45 
     46 /*
     47  *                  x^2   x^3   x^4
     48  * exp(x) = 1 + x + --- + --- + --- + ...
     49  *                   2!    3!    4!
     50  */
     51 static struct fpn *
     52 fpu_etox_taylor(struct fpemu *fe)
     53 {
     54 	struct fpn res;
     55 	struct fpn x;
     56 	struct fpn s0;
     57 	struct fpn *s1;
     58 	struct fpn *r;
     59 	uint32_t k;
     60 
     61 	CPYFPN(&x, &fe->fe_f2);
     62 	CPYFPN(&s0, &fe->fe_f2);
     63 
     64 	/* res := 1 + x */
     65 	fpu_const(&fe->fe_f1, FPU_CONST_1);
     66 	r = fpu_add(fe);
     67 	CPYFPN(&res, r);
     68 
     69 	k = 2;
     70 	for (; k < MAX_ITEMS; k++) {
     71 		/* s1 = s0 * x / k */
     72 		CPYFPN(&fe->fe_f1, &s0);
     73 		CPYFPN(&fe->fe_f2, &x);
     74 		r = fpu_mul(fe);
     75 
     76 		CPYFPN(&fe->fe_f1, r);
     77 		fpu_explode(fe, &fe->fe_f2, FTYPE_LNG, &k);
     78 		s1 = fpu_div(fe);
     79 
     80 		/* break if s1 is enough small */
     81 		if (ISZERO(s1))
     82 			break;
     83 		if (res.fp_exp - s1->fp_exp >= FP_NMANT)
     84 			break;
     85 
     86 		/* s0 := s1 for next loop */
     87 		CPYFPN(&s0, s1);
     88 
     89 		/* res += s1 */
     90 		CPYFPN(&fe->fe_f2, s1);
     91 		CPYFPN(&fe->fe_f1, &res);
     92 		r = fpu_add(fe);
     93 		CPYFPN(&res, r);
     94 	}
     95 
     96 	CPYFPN(&fe->fe_f2, &res);
     97 	return &fe->fe_f2;
     98 }
     99 
    100 /*
    101  * exp(x)
    102  */
    103 struct fpn *
    104 fpu_etox(struct fpemu *fe)
    105 {
    106 	struct fpn *fp;
    107 
    108 	if (ISNAN(&fe->fe_f2))
    109 		return &fe->fe_f2;
    110 	if (ISINF(&fe->fe_f2)) {
    111 		if (fe->fe_f2.fp_sign)
    112 			fpu_const(&fe->fe_f2, FPU_CONST_0);
    113 		return &fe->fe_f2;
    114 	}
    115 
    116 	if (fe->fe_f2.fp_sign == 0) {
    117 		/* exp(x) */
    118 		fp = fpu_etox_taylor(fe);
    119 	} else {
    120 		/* 1/exp(-x) */
    121 		fe->fe_f2.fp_sign = 0;
    122 		fp = fpu_etox_taylor(fe);
    123 
    124 		CPYFPN(&fe->fe_f2, fp);
    125 		fpu_const(&fe->fe_f1, FPU_CONST_1);
    126 		fp = fpu_div(fe);
    127 	}
    128 
    129 	return fp;
    130 }
    131 
    132 /*
    133  * exp(x) - 1
    134  */
    135 struct fpn *
    136 fpu_etoxm1(struct fpemu *fe)
    137 {
    138 	struct fpn *fp;
    139 
    140 	fp = fpu_etox(fe);
    141 
    142 	CPYFPN(&fe->fe_f1, fp);
    143 	/* build a 1.0 */
    144 	fp = fpu_const(&fe->fe_f2, FPU_CONST_1);
    145 	fe->fe_f2.fp_sign = !fe->fe_f2.fp_sign;
    146 	/* fp = f2 - 1.0 */
    147 	fp = fpu_add(fe);
    148 
    149 	return fp;
    150 }
    151 
    152 /*
    153  * 10^x = exp(x * ln10)
    154  */
    155 struct fpn *
    156 fpu_tentox(struct fpemu *fe)
    157 {
    158 	struct fpn *fp;
    159 
    160 	/* build a ln10 */
    161 	fp = fpu_const(&fe->fe_f1, FPU_CONST_LN_10);
    162 	/* fp = ln10 * f2 */
    163 	fp = fpu_mul(fe);
    164 
    165 	/* copy the result to the src opr */
    166 	CPYFPN(&fe->fe_f2, fp);
    167 
    168 	return fpu_etox(fe);
    169 }
    170 
    171 /*
    172  * 2^x = exp(x * ln2)
    173  */
    174 struct fpn *
    175 fpu_twotox(struct fpemu *fe)
    176 {
    177 	struct fpn *fp;
    178 
    179 	/* build a ln2 */
    180 	fp = fpu_const(&fe->fe_f1, FPU_CONST_LN_2);
    181 	/* fp = ln2 * f2 */
    182 	fp = fpu_mul(fe);
    183 
    184 	/* copy the result to the src opr */
    185 	CPYFPN(&fe->fe_f2, fp);
    186 
    187 	return fpu_etox(fe);
    188 }
    189