Home | History | Annotate | Line # | Download | only in fpe
fpu_rem.c revision 1.6.78.1
      1  1.6.78.1    yamt /*	$NetBSD: fpu_rem.c,v 1.6.78.1 2009/05/04 08:11:25 yamt 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.6.78.1    yamt __KERNEL_RCSID(0, "$NetBSD: fpu_rem.c,v 1.6.78.1 2009/05/04 08:11:25 yamt 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.6.78.1    yamt static struct fpn * __fpu_modrem(struct fpemu *fe, int modrem);
     89       1.1  briggs 
     90       1.1  briggs static struct fpn *
     91  1.6.78.1    yamt __fpu_modrem(struct fpemu *fe, int modrem)
     92       1.1  briggs {
     93       1.1  briggs     static struct fpn X, Y;
     94       1.1  briggs     struct fpn *x, *y, *r;
     95       1.1  briggs     u_int signX, signY, signQ;
     96       1.2  briggs     int j, k, l, q;
     97       1.1  briggs     int Last_Subtract;
     98       1.1  briggs 
     99       1.1  briggs     CPYFPN(&X, &fe->fe_f1);
    100       1.1  briggs     CPYFPN(&Y, &fe->fe_f2);
    101       1.1  briggs     x = &X;
    102       1.1  briggs     y = &Y;
    103       1.1  briggs     r = &fe->fe_f2;
    104       1.1  briggs 
    105       1.1  briggs     /*
    106       1.1  briggs      * Step 1
    107       1.1  briggs      */
    108       1.1  briggs     signX = x->fp_sign;
    109       1.1  briggs     signY = y->fp_sign;
    110       1.1  briggs     signQ = (signX ^ signY);
    111       1.1  briggs     x->fp_sign = y->fp_sign = 0;
    112       1.1  briggs 
    113       1.1  briggs     /*
    114       1.1  briggs      * Step 2
    115       1.1  briggs      */
    116       1.1  briggs     l = x->fp_exp - y->fp_exp;
    117       1.1  briggs     k = 0;
    118       1.1  briggs     q = 0;
    119       1.3  briggs     if (l >= 0) {
    120       1.1  briggs 	CPYFPN(r, x);
    121       1.1  briggs 	r->fp_exp -= l;
    122       1.1  briggs 	j = l;
    123       1.1  briggs 
    124       1.1  briggs 	/*
    125       1.1  briggs 	 * Step 3
    126       1.1  briggs 	 */
    127       1.1  briggs 	while (y->fp_exp != r->fp_exp || y->fp_mant[0] != r->fp_mant[0] ||
    128       1.1  briggs 	       y->fp_mant[1] != r->fp_mant[1] ||
    129       1.4  briggs 	       y->fp_mant[2] != r->fp_mant[2]) {
    130       1.1  briggs 
    131       1.1  briggs 	    /* Step 3.2 */
    132       1.1  briggs 	    if (y->fp_exp < r->fp_exp || y->fp_mant[0] < r->fp_mant[0] ||
    133       1.1  briggs 		y->fp_mant[1] < r->fp_mant[1] ||
    134       1.4  briggs 		y->fp_mant[2] < r->fp_mant[2]) {
    135       1.1  briggs 		CPYFPN(&fe->fe_f1, r);
    136       1.1  briggs 		CPYFPN(&fe->fe_f2, y);
    137       1.1  briggs 		fe->fe_f2.fp_sign = 1;
    138       1.1  briggs 		r = fpu_add(fe);
    139       1.1  briggs 		q++;
    140       1.1  briggs 	    }
    141       1.1  briggs 
    142       1.1  briggs 	    /* Step 3.3 */
    143       1.1  briggs 	    if (j == 0)
    144       1.1  briggs 		goto Step4;
    145       1.1  briggs 
    146       1.1  briggs 	    /* Step 3.4 */
    147       1.1  briggs 	    k++;
    148       1.1  briggs 	    j--;
    149       1.1  briggs 	    q += q;
    150       1.1  briggs 	    r->fp_exp++;
    151       1.1  briggs 	}
    152       1.3  briggs 	/* Step 9 */
    153       1.3  briggs 	goto Step9;
    154       1.1  briggs     }
    155       1.1  briggs  Step4:
    156       1.1  briggs     Last_Subtract = 0;
    157       1.1  briggs     if (modrem == 0)
    158       1.1  briggs 	goto Step6;
    159       1.1  briggs 
    160       1.1  briggs     /*
    161       1.1  briggs      * Step 5
    162       1.1  briggs      */
    163       1.1  briggs     /* Step 5.1 */
    164       1.1  briggs     if (r->fp_exp + 1 < y->fp_exp ||
    165       1.2  briggs 	(r->fp_exp + 1 == y->fp_exp &&
    166       1.2  briggs 	 (r->fp_mant[0] < y->fp_mant[0] || r->fp_mant[1] < y->fp_mant[1] ||
    167       1.4  briggs 	  r->fp_mant[2] < y->fp_mant[2])))
    168       1.1  briggs 	/* if r < y/2 */
    169       1.1  briggs 	goto Step6;
    170       1.1  briggs     /* Step 5.2 */
    171       1.1  briggs     if (r->fp_exp + 1 != y->fp_exp ||
    172       1.1  briggs 	r->fp_mant[0] != y->fp_mant[0] || r->fp_mant[1] != y->fp_mant[1] ||
    173       1.4  briggs 	r->fp_mant[2] != y->fp_mant[2]) {
    174       1.1  briggs 	/* if (!(r < y/2) && !(r == y/2)) */
    175       1.1  briggs 	Last_Subtract = 1;
    176       1.1  briggs 	q++;
    177       1.1  briggs 	y->fp_sign = signY;
    178       1.1  briggs     } else {
    179       1.1  briggs 	/* Step 5.3 */
    180       1.1  briggs 	/* r == y/2 */
    181       1.1  briggs 	if (q % 2) {
    182       1.1  briggs 	    q++;
    183       1.1  briggs 	    signX = !signX;
    184       1.1  briggs 	}
    185       1.1  briggs     }
    186       1.1  briggs 
    187       1.1  briggs  Step6:
    188       1.1  briggs     r->fp_sign = signX;
    189       1.1  briggs 
    190       1.1  briggs     /*
    191       1.1  briggs      * Step 7
    192       1.1  briggs      */
    193       1.1  briggs     if (Last_Subtract) {
    194       1.1  briggs 	CPYFPN(&fe->fe_f1, r);
    195       1.1  briggs 	CPYFPN(&fe->fe_f2, y);
    196       1.1  briggs 	fe->fe_f2.fp_sign = !y->fp_sign;
    197       1.1  briggs 	r = fpu_add(fe);
    198       1.1  briggs     }
    199       1.1  briggs     /*
    200       1.1  briggs      * Step 8
    201       1.1  briggs      */
    202       1.1  briggs     q &= 0x7f;
    203       1.1  briggs     q |= (signQ << 7);
    204       1.1  briggs     fe->fe_fpframe->fpf_fpsr =
    205       1.1  briggs 	fe->fe_fpsr =
    206       1.1  briggs 	    (fe->fe_fpsr & ~FPSR_QTT) | (q << 16);
    207       1.1  briggs     return r;
    208       1.3  briggs 
    209       1.3  briggs  Step9:
    210       1.3  briggs     fe->fe_f1.fp_class = FPC_ZERO;
    211       1.3  briggs     q++;
    212       1.3  briggs     q &= 0x7f;
    213       1.3  briggs     q |= (signQ << 7);
    214       1.3  briggs     fe->fe_fpframe->fpf_fpsr =
    215       1.3  briggs 	fe->fe_fpsr =
    216       1.3  briggs 	    (fe->fe_fpsr & ~FPSR_QTT) | (q << 16);
    217       1.3  briggs     return &fe->fe_f1;
    218       1.1  briggs }
    219       1.1  briggs 
    220       1.1  briggs struct fpn *
    221  1.6.78.1    yamt fpu_rem(struct fpemu *fe)
    222       1.1  briggs {
    223       1.1  briggs   return __fpu_modrem(fe, 1);
    224       1.1  briggs }
    225       1.1  briggs 
    226       1.1  briggs struct fpn *
    227  1.6.78.1    yamt fpu_mod(struct fpemu *fe)
    228       1.1  briggs {
    229       1.1  briggs   return __fpu_modrem(fe, 0);
    230       1.1  briggs }
    231