Home | History | Annotate | Line # | Download | only in fpe
fpu_rem.c revision 1.1
      1  1.1  briggs /*	$NetBSD: fpu_rem.c,v 1.1 1995/11/03 04:47:17 briggs 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.1  briggs 
     34  1.1  briggs #include <sys/types.h>
     35  1.1  briggs #include <sys/signal.h>
     36  1.1  briggs #include <machine/frame.h>
     37  1.1  briggs 
     38  1.1  briggs #include "fpu_emulate.h"
     39  1.1  briggs 
     40  1.1  briggs /*
     41  1.1  briggs  *       ALGORITHM
     42  1.1  briggs  *
     43  1.1  briggs  *       Step 1.  Save and strip signs of X and Y: signX := sign(X),
     44  1.1  briggs  *                signY := sign(Y), X := *X*, Y := *Y*,
     45  1.1  briggs  *                signQ := signX EOR signY. Record whether MOD or REM
     46  1.1  briggs  *                is requested.
     47  1.1  briggs  *
     48  1.1  briggs  *       Step 2.  Set L := expo(X)-expo(Y), k := 0, Q := 0.
     49  1.1  briggs  *                If (L < 0) then
     50  1.1  briggs  *                   R := X, go to Step 4.
     51  1.1  briggs  *                else
     52  1.1  briggs  *                   R := 2^(-L)X, j := L.
     53  1.1  briggs  *                endif
     54  1.1  briggs  *
     55  1.1  briggs  *       Step 3.  Perform MOD(X,Y)
     56  1.1  briggs  *            3.1 If R = Y, go to Step 9.
     57  1.1  briggs  *            3.2 If R > Y, then { R := R - Y, Q := Q + 1}
     58  1.1  briggs  *            3.3 If j = 0, go to Step 4.
     59  1.1  briggs  *            3.4 k := k + 1, j := j - 1, Q := 2Q, R := 2R. Go to
     60  1.1  briggs  *                Step 3.1.
     61  1.1  briggs  *
     62  1.1  briggs  *       Step 4.  At this point, R = X - QY = MOD(X,Y). Set
     63  1.1  briggs  *                Last_Subtract := false (used in Step 7 below). If
     64  1.1  briggs  *                MOD is requested, go to Step 6.
     65  1.1  briggs  *
     66  1.1  briggs  *       Step 5.  R = MOD(X,Y), but REM(X,Y) is requested.
     67  1.1  briggs  *            5.1 If R < Y/2, then R = MOD(X,Y) = REM(X,Y). Go to
     68  1.1  briggs  *                Step 6.
     69  1.1  briggs  *            5.2 If R > Y/2, then { set Last_Subtract := true,
     70  1.1  briggs  *                Q := Q + 1, Y := signY*Y }. Go to Step 6.
     71  1.1  briggs  *            5.3 This is the tricky case of R = Y/2. If Q is odd,
     72  1.1  briggs  *                then { Q := Q + 1, signX := -signX }.
     73  1.1  briggs  *
     74  1.1  briggs  *       Step 6.  R := signX*R.
     75  1.1  briggs  *
     76  1.1  briggs  *       Step 7.  If Last_Subtract = true, R := R - Y.
     77  1.1  briggs  *
     78  1.1  briggs  *       Step 8.  Return signQ, last 7 bits of Q, and R as required.
     79  1.1  briggs  *
     80  1.1  briggs  *       Step 9.  At this point, R = 2^(-j)*X - Q Y = Y. Thus,
     81  1.1  briggs  *                X = 2^(j)*(Q+1)Y. set Q := 2^(j)*(Q+1),
     82  1.1  briggs  *                R := 0. Return signQ, last 7 bits of Q, and R.
     83  1.1  briggs  */
     84  1.1  briggs 
     85  1.1  briggs static struct fpn * __fpu_modrem __P((struct fpemu *fe, int modrem));
     86  1.1  briggs 
     87  1.1  briggs static struct fpn *
     88  1.1  briggs __fpu_modrem(fe, modrem)
     89  1.1  briggs      struct fpemu *fe;
     90  1.1  briggs      int modrem;
     91  1.1  briggs {
     92  1.1  briggs     static struct fpn X, Y;
     93  1.1  briggs     struct fpn *x, *y, *r;
     94  1.1  briggs     u_int signX, signY, signQ;
     95  1.1  briggs     int i, j, k, l, q;
     96  1.1  briggs     int Last_Subtract;
     97  1.1  briggs 
     98  1.1  briggs     CPYFPN(&X, &fe->fe_f1);
     99  1.1  briggs     CPYFPN(&Y, &fe->fe_f2);
    100  1.1  briggs     x = &X;
    101  1.1  briggs     y = &Y;
    102  1.1  briggs     r = &fe->fe_f2;
    103  1.1  briggs 
    104  1.1  briggs     /*
    105  1.1  briggs      * Step 1
    106  1.1  briggs      */
    107  1.1  briggs     signX = x->fp_sign;
    108  1.1  briggs     signY = y->fp_sign;
    109  1.1  briggs     signQ = (signX ^ signY);
    110  1.1  briggs     x->fp_sign = y->fp_sign = 0;
    111  1.1  briggs 
    112  1.1  briggs     /*
    113  1.1  briggs      * Step 2
    114  1.1  briggs      */
    115  1.1  briggs     l = x->fp_exp - y->fp_exp;
    116  1.1  briggs     k = 0;
    117  1.1  briggs     q = 0;
    118  1.1  briggs     if (l < 0) {
    119  1.1  briggs 	goto Step4;
    120  1.1  briggs     } else {
    121  1.1  briggs 	CPYFPN(r, x);
    122  1.1  briggs 	r->fp_exp -= l;
    123  1.1  briggs 	j = l;
    124  1.1  briggs 
    125  1.1  briggs 	/*
    126  1.1  briggs 	 * Step 3
    127  1.1  briggs 	 */
    128  1.1  briggs 	while (y->fp_exp != r->fp_exp || y->fp_mant[0] != r->fp_mant[0] ||
    129  1.1  briggs 	       y->fp_mant[1] != r->fp_mant[1] ||
    130  1.1  briggs 	       y->fp_mant[2] != r->fp_mant[2] ||
    131  1.1  briggs 	       y->fp_mant[3] != r->fp_mant[3]) {
    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.1  briggs 		y->fp_mant[2] < r->fp_mant[2] ||
    137  1.1  briggs 		y->fp_mant[3] < r->fp_mant[3]) {
    138  1.1  briggs 		CPYFPN(&fe->fe_f1, r);
    139  1.1  briggs 		CPYFPN(&fe->fe_f2, y);
    140  1.1  briggs 		fe->fe_f2.fp_sign = 1;
    141  1.1  briggs 		r = fpu_add(fe);
    142  1.1  briggs 		q++;
    143  1.1  briggs 	    }
    144  1.1  briggs 
    145  1.1  briggs 	    /* Step 3.3 */
    146  1.1  briggs 	    if (j == 0)
    147  1.1  briggs 		goto Step4;
    148  1.1  briggs 
    149  1.1  briggs 	    /* Step 3.4 */
    150  1.1  briggs 	    k++;
    151  1.1  briggs 	    j--;
    152  1.1  briggs 	    q += q;
    153  1.1  briggs 	    r->fp_exp++;
    154  1.1  briggs 	}
    155  1.1  briggs 	/* Step 9 */
    156  1.1  briggs 
    157  1.1  briggs     }
    158  1.1  briggs  Step4:
    159  1.1  briggs     Last_Subtract = 0;
    160  1.1  briggs     if (modrem == 0)
    161  1.1  briggs 	goto Step6;
    162  1.1  briggs 
    163  1.1  briggs     /*
    164  1.1  briggs      * Step 5
    165  1.1  briggs      */
    166  1.1  briggs     /* Step 5.1 */
    167  1.1  briggs     if (r->fp_exp + 1 < y->fp_exp ||
    168  1.1  briggs 	r->fp_exp + 1 == y->fp_exp &&
    169  1.1  briggs 	(r->fp_mant[0] < y->fp_mant[0] || r->fp_mant[1] < y->fp_mant[1] ||
    170  1.1  briggs 	 r->fp_mant[2] < y->fp_mant[3] || r->fp_mant[4] < y->fp_mant[4]))
    171  1.1  briggs 	/* if r < y/2 */
    172  1.1  briggs 	goto Step6;
    173  1.1  briggs     /* Step 5.2 */
    174  1.1  briggs     if (r->fp_exp + 1 != y->fp_exp ||
    175  1.1  briggs 	r->fp_mant[0] != y->fp_mant[0] || r->fp_mant[1] != y->fp_mant[1] ||
    176  1.1  briggs 	r->fp_mant[2] != y->fp_mant[2] || r->fp_mant[3] != y->fp_mant[3]) {
    177  1.1  briggs 	/* if (!(r < y/2) && !(r == y/2)) */
    178  1.1  briggs 	Last_Subtract = 1;
    179  1.1  briggs 	q++;
    180  1.1  briggs 	y->fp_sign = signY;
    181  1.1  briggs     } else {
    182  1.1  briggs 	/* Step 5.3 */
    183  1.1  briggs 	/* r == y/2 */
    184  1.1  briggs 	if (q % 2) {
    185  1.1  briggs 	    q++;
    186  1.1  briggs 	    signX = !signX;
    187  1.1  briggs 	}
    188  1.1  briggs     }
    189  1.1  briggs 
    190  1.1  briggs  Step6:
    191  1.1  briggs     r->fp_sign = signX;
    192  1.1  briggs 
    193  1.1  briggs     /*
    194  1.1  briggs      * Step 7
    195  1.1  briggs      */
    196  1.1  briggs     if (Last_Subtract) {
    197  1.1  briggs 	CPYFPN(&fe->fe_f1, r);
    198  1.1  briggs 	CPYFPN(&fe->fe_f2, y);
    199  1.1  briggs 	fe->fe_f2.fp_sign = !y->fp_sign;
    200  1.1  briggs 	r = fpu_add(fe);
    201  1.1  briggs     }
    202  1.1  briggs     /*
    203  1.1  briggs      * Step 8
    204  1.1  briggs      */
    205  1.1  briggs     q &= 0x7f;
    206  1.1  briggs     q |= (signQ << 7);
    207  1.1  briggs     fe->fe_fpframe->fpf_fpsr =
    208  1.1  briggs 	fe->fe_fpsr =
    209  1.1  briggs 	    (fe->fe_fpsr & ~FPSR_QTT) | (q << 16);
    210  1.1  briggs     return r;
    211  1.1  briggs 
    212  1.1  briggs  Step9:
    213  1.1  briggs     fe->fe_f1.fp_class = FPC_ZERO;
    214  1.1  briggs     q++;
    215  1.1  briggs     q &= 0x7f;
    216  1.1  briggs     q |= (signQ << 7);
    217  1.1  briggs     fe->fe_fpframe->fpf_fpsr =
    218  1.1  briggs 	fe->fe_fpsr =
    219  1.1  briggs 	    (fe->fe_fpsr & ~FPSR_QTT) | (q << 16);
    220  1.1  briggs     return &fe->fe_f1;
    221  1.1  briggs }
    222  1.1  briggs 
    223  1.1  briggs struct fpn *
    224  1.1  briggs fpu_rem(fe)
    225  1.1  briggs      struct fpemu *fe;
    226  1.1  briggs {
    227  1.1  briggs   return __fpu_modrem(fe, 1);
    228  1.1  briggs }
    229  1.1  briggs 
    230  1.1  briggs struct fpn *
    231  1.1  briggs fpu_mod(fe)
    232  1.1  briggs      struct fpemu *fe;
    233  1.1  briggs {
    234  1.1  briggs   return __fpu_modrem(fe, 0);
    235  1.1  briggs }
    236