Home | History | Annotate | Line # | Download | only in fpe
fpu_rem.c revision 1.8
      1  1.8     dsl /*	$NetBSD: fpu_rem.c,v 1.8 2009/03/14 15:36:09 dsl Exp $	*/
      2  1.1  briggs 
      3  1.1  briggs /*
      4  1.1  briggs  * Copyright (c) 1995  Ken Nakata
      5  1.1  briggs  *	All rights reserved.
      6  1.1  briggs  *
      7  1.1  briggs  * Redistribution and use in source and binary forms, with or without
      8  1.1  briggs  * modification, are permitted provided that the following conditions
      9  1.1  briggs  * are met:
     10  1.1  briggs  * 1. Redistributions of source code must retain the above copyright
     11  1.1  briggs  *    notice, this list of conditions and the following disclaimer.
     12  1.1  briggs  * 2. Redistributions in binary form must reproduce the above copyright
     13  1.1  briggs  *    notice, this list of conditions and the following disclaimer in the
     14  1.1  briggs  *    documentation and/or other materials provided with the distribution.
     15  1.1  briggs  * 3. Neither the name of the author nor the names of its contributors
     16  1.1  briggs  *    may be used to endorse or promote products derived from this software
     17  1.1  briggs  *    without specific prior written permission.
     18  1.1  briggs  *
     19  1.1  briggs  * THIS SOFTWARE IS PROVIDED BY THE AUTHOR AND CONTRIBUTORS ``AS IS'' AND
     20  1.1  briggs  * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
     21  1.1  briggs  * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
     22  1.1  briggs  * ARE DISCLAIMED.  IN NO EVENT SHALL THE AUTHOR OR CONTRIBUTORS BE LIABLE
     23  1.1  briggs  * FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
     24  1.1  briggs  * DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS
     25  1.1  briggs  * OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
     26  1.1  briggs  * HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
     27  1.1  briggs  * LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY
     28  1.1  briggs  * OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
     29  1.1  briggs  * SUCH DAMAGE.
     30  1.1  briggs  *
     31  1.1  briggs  *	@(#)fpu_rem.c	10/24/95
     32  1.1  briggs  */
     33  1.5   lukem 
     34  1.5   lukem #include <sys/cdefs.h>
     35  1.8     dsl __KERNEL_RCSID(0, "$NetBSD: fpu_rem.c,v 1.8 2009/03/14 15:36:09 dsl Exp $");
     36  1.1  briggs 
     37  1.1  briggs #include <sys/types.h>
     38  1.1  briggs #include <sys/signal.h>
     39  1.1  briggs #include <machine/frame.h>
     40  1.1  briggs 
     41  1.1  briggs #include "fpu_emulate.h"
     42  1.1  briggs 
     43  1.1  briggs /*
     44  1.1  briggs  *       ALGORITHM
     45  1.1  briggs  *
     46  1.1  briggs  *       Step 1.  Save and strip signs of X and Y: signX := sign(X),
     47  1.1  briggs  *                signY := sign(Y), X := *X*, Y := *Y*,
     48  1.1  briggs  *                signQ := signX EOR signY. Record whether MOD or REM
     49  1.1  briggs  *                is requested.
     50  1.1  briggs  *
     51  1.1  briggs  *       Step 2.  Set L := expo(X)-expo(Y), k := 0, Q := 0.
     52  1.1  briggs  *                If (L < 0) then
     53  1.1  briggs  *                   R := X, go to Step 4.
     54  1.1  briggs  *                else
     55  1.1  briggs  *                   R := 2^(-L)X, j := L.
     56  1.1  briggs  *                endif
     57  1.1  briggs  *
     58  1.1  briggs  *       Step 3.  Perform MOD(X,Y)
     59  1.1  briggs  *            3.1 If R = Y, go to Step 9.
     60  1.1  briggs  *            3.2 If R > Y, then { R := R - Y, Q := Q + 1}
     61  1.1  briggs  *            3.3 If j = 0, go to Step 4.
     62  1.1  briggs  *            3.4 k := k + 1, j := j - 1, Q := 2Q, R := 2R. Go to
     63  1.1  briggs  *                Step 3.1.
     64  1.1  briggs  *
     65  1.1  briggs  *       Step 4.  At this point, R = X - QY = MOD(X,Y). Set
     66  1.1  briggs  *                Last_Subtract := false (used in Step 7 below). If
     67  1.1  briggs  *                MOD is requested, go to Step 6.
     68  1.1  briggs  *
     69  1.1  briggs  *       Step 5.  R = MOD(X,Y), but REM(X,Y) is requested.
     70  1.1  briggs  *            5.1 If R < Y/2, then R = MOD(X,Y) = REM(X,Y). Go to
     71  1.1  briggs  *                Step 6.
     72  1.1  briggs  *            5.2 If R > Y/2, then { set Last_Subtract := true,
     73  1.1  briggs  *                Q := Q + 1, Y := signY*Y }. Go to Step 6.
     74  1.1  briggs  *            5.3 This is the tricky case of R = Y/2. If Q is odd,
     75  1.1  briggs  *                then { Q := Q + 1, signX := -signX }.
     76  1.1  briggs  *
     77  1.1  briggs  *       Step 6.  R := signX*R.
     78  1.1  briggs  *
     79  1.1  briggs  *       Step 7.  If Last_Subtract = true, R := R - Y.
     80  1.1  briggs  *
     81  1.1  briggs  *       Step 8.  Return signQ, last 7 bits of Q, and R as required.
     82  1.1  briggs  *
     83  1.1  briggs  *       Step 9.  At this point, R = 2^(-j)*X - Q Y = Y. Thus,
     84  1.1  briggs  *                X = 2^(j)*(Q+1)Y. set Q := 2^(j)*(Q+1),
     85  1.1  briggs  *                R := 0. Return signQ, last 7 bits of Q, and R.
     86  1.1  briggs  */
     87  1.1  briggs 
     88  1.7     dsl static struct fpn * __fpu_modrem(struct fpemu *fe, int modrem);
     89  1.1  briggs 
     90  1.1  briggs static struct fpn *
     91  1.1  briggs __fpu_modrem(fe, modrem)
     92  1.1  briggs      struct fpemu *fe;
     93  1.1  briggs      int modrem;
     94  1.1  briggs {
     95  1.1  briggs     static struct fpn X, Y;
     96  1.1  briggs     struct fpn *x, *y, *r;
     97  1.1  briggs     u_int signX, signY, signQ;
     98  1.2  briggs     int j, k, l, q;
     99  1.1  briggs     int Last_Subtract;
    100  1.1  briggs 
    101  1.1  briggs     CPYFPN(&X, &fe->fe_f1);
    102  1.1  briggs     CPYFPN(&Y, &fe->fe_f2);
    103  1.1  briggs     x = &X;
    104  1.1  briggs     y = &Y;
    105  1.1  briggs     r = &fe->fe_f2;
    106  1.1  briggs 
    107  1.1  briggs     /*
    108  1.1  briggs      * Step 1
    109  1.1  briggs      */
    110  1.1  briggs     signX = x->fp_sign;
    111  1.1  briggs     signY = y->fp_sign;
    112  1.1  briggs     signQ = (signX ^ signY);
    113  1.1  briggs     x->fp_sign = y->fp_sign = 0;
    114  1.1  briggs 
    115  1.1  briggs     /*
    116  1.1  briggs      * Step 2
    117  1.1  briggs      */
    118  1.1  briggs     l = x->fp_exp - y->fp_exp;
    119  1.1  briggs     k = 0;
    120  1.1  briggs     q = 0;
    121  1.3  briggs     if (l >= 0) {
    122  1.1  briggs 	CPYFPN(r, x);
    123  1.1  briggs 	r->fp_exp -= l;
    124  1.1  briggs 	j = l;
    125  1.1  briggs 
    126  1.1  briggs 	/*
    127  1.1  briggs 	 * Step 3
    128  1.1  briggs 	 */
    129  1.1  briggs 	while (y->fp_exp != r->fp_exp || y->fp_mant[0] != r->fp_mant[0] ||
    130  1.1  briggs 	       y->fp_mant[1] != r->fp_mant[1] ||
    131  1.4  briggs 	       y->fp_mant[2] != r->fp_mant[2]) {
    132  1.1  briggs 
    133  1.1  briggs 	    /* Step 3.2 */
    134  1.1  briggs 	    if (y->fp_exp < r->fp_exp || y->fp_mant[0] < r->fp_mant[0] ||
    135  1.1  briggs 		y->fp_mant[1] < r->fp_mant[1] ||
    136  1.4  briggs 		y->fp_mant[2] < r->fp_mant[2]) {
    137  1.1  briggs 		CPYFPN(&fe->fe_f1, r);
    138  1.1  briggs 		CPYFPN(&fe->fe_f2, y);
    139  1.1  briggs 		fe->fe_f2.fp_sign = 1;
    140  1.1  briggs 		r = fpu_add(fe);
    141  1.1  briggs 		q++;
    142  1.1  briggs 	    }
    143  1.1  briggs 
    144  1.1  briggs 	    /* Step 3.3 */
    145  1.1  briggs 	    if (j == 0)
    146  1.1  briggs 		goto Step4;
    147  1.1  briggs 
    148  1.1  briggs 	    /* Step 3.4 */
    149  1.1  briggs 	    k++;
    150  1.1  briggs 	    j--;
    151  1.1  briggs 	    q += q;
    152  1.1  briggs 	    r->fp_exp++;
    153  1.1  briggs 	}
    154  1.3  briggs 	/* Step 9 */
    155  1.3  briggs 	goto Step9;
    156  1.1  briggs     }
    157  1.1  briggs  Step4:
    158  1.1  briggs     Last_Subtract = 0;
    159  1.1  briggs     if (modrem == 0)
    160  1.1  briggs 	goto Step6;
    161  1.1  briggs 
    162  1.1  briggs     /*
    163  1.1  briggs      * Step 5
    164  1.1  briggs      */
    165  1.1  briggs     /* Step 5.1 */
    166  1.1  briggs     if (r->fp_exp + 1 < y->fp_exp ||
    167  1.2  briggs 	(r->fp_exp + 1 == y->fp_exp &&
    168  1.2  briggs 	 (r->fp_mant[0] < y->fp_mant[0] || r->fp_mant[1] < y->fp_mant[1] ||
    169  1.4  briggs 	  r->fp_mant[2] < y->fp_mant[2])))
    170  1.1  briggs 	/* if r < y/2 */
    171  1.1  briggs 	goto Step6;
    172  1.1  briggs     /* Step 5.2 */
    173  1.1  briggs     if (r->fp_exp + 1 != y->fp_exp ||
    174  1.1  briggs 	r->fp_mant[0] != y->fp_mant[0] || r->fp_mant[1] != y->fp_mant[1] ||
    175  1.4  briggs 	r->fp_mant[2] != y->fp_mant[2]) {
    176  1.1  briggs 	/* if (!(r < y/2) && !(r == y/2)) */
    177  1.1  briggs 	Last_Subtract = 1;
    178  1.1  briggs 	q++;
    179  1.1  briggs 	y->fp_sign = signY;
    180  1.1  briggs     } else {
    181  1.1  briggs 	/* Step 5.3 */
    182  1.1  briggs 	/* r == y/2 */
    183  1.1  briggs 	if (q % 2) {
    184  1.1  briggs 	    q++;
    185  1.1  briggs 	    signX = !signX;
    186  1.1  briggs 	}
    187  1.1  briggs     }
    188  1.1  briggs 
    189  1.1  briggs  Step6:
    190  1.1  briggs     r->fp_sign = signX;
    191  1.1  briggs 
    192  1.1  briggs     /*
    193  1.1  briggs      * Step 7
    194  1.1  briggs      */
    195  1.1  briggs     if (Last_Subtract) {
    196  1.1  briggs 	CPYFPN(&fe->fe_f1, r);
    197  1.1  briggs 	CPYFPN(&fe->fe_f2, y);
    198  1.1  briggs 	fe->fe_f2.fp_sign = !y->fp_sign;
    199  1.1  briggs 	r = fpu_add(fe);
    200  1.1  briggs     }
    201  1.1  briggs     /*
    202  1.1  briggs      * Step 8
    203  1.1  briggs      */
    204  1.1  briggs     q &= 0x7f;
    205  1.1  briggs     q |= (signQ << 7);
    206  1.1  briggs     fe->fe_fpframe->fpf_fpsr =
    207  1.1  briggs 	fe->fe_fpsr =
    208  1.1  briggs 	    (fe->fe_fpsr & ~FPSR_QTT) | (q << 16);
    209  1.1  briggs     return r;
    210  1.3  briggs 
    211  1.3  briggs  Step9:
    212  1.3  briggs     fe->fe_f1.fp_class = FPC_ZERO;
    213  1.3  briggs     q++;
    214  1.3  briggs     q &= 0x7f;
    215  1.3  briggs     q |= (signQ << 7);
    216  1.3  briggs     fe->fe_fpframe->fpf_fpsr =
    217  1.3  briggs 	fe->fe_fpsr =
    218  1.3  briggs 	    (fe->fe_fpsr & ~FPSR_QTT) | (q << 16);
    219  1.3  briggs     return &fe->fe_f1;
    220  1.1  briggs }
    221  1.1  briggs 
    222  1.1  briggs struct fpn *
    223  1.8     dsl fpu_rem(struct fpemu *fe)
    224  1.1  briggs {
    225  1.1  briggs   return __fpu_modrem(fe, 1);
    226  1.1  briggs }
    227  1.1  briggs 
    228  1.1  briggs struct fpn *
    229  1.8     dsl fpu_mod(struct fpemu *fe)
    230  1.1  briggs {
    231  1.1  briggs   return __fpu_modrem(fe, 0);
    232  1.1  briggs }
    233