Home | History | Annotate | Download | only in common
      1 /*
      2  * CDDL HEADER START
      3  *
      4  * The contents of this file are subject to the terms of the
      5  * Common Development and Distribution License, Version 1.0 only
      6  * (the "License").  You may not use this file except in compliance
      7  * with the License.
      8  *
      9  * You can obtain a copy of the license at usr/src/OPENSOLARIS.LICENSE
     10  * or http://www.opensolaris.org/os/licensing.
     11  * See the License for the specific language governing permissions
     12  * and limitations under the License.
     13  *
     14  * When distributing Covered Code, include this CDDL HEADER in each
     15  * file and include the License file at usr/src/OPENSOLARIS.LICENSE.
     16  * If applicable, add the following below this CDDL HEADER, with the
     17  * fields enclosed by brackets "[]" replaced with your own identifying
     18  * information: Portions Copyright [yyyy] [name of copyright owner]
     19  *
     20  * CDDL HEADER END
     21  */
     22 #pragma ident	"%Z%%M%	%I%	%E% SMI"
     23 
     24 /*
     25  * Copyright (c) 1988 by Sun Microsystems, Inc.
     26  */
     27 
     28 #include "_Qquad.h"
     29 #include "_Qglobals.h"
     30 
     31 
     32 void
     33 _fp_mul(px, py, pz)
     34 	unpacked       *px, *py, *pz;
     35 
     36 {
     37 	unpacked       *pt;
     38 	unsigned       acc[4];		/* Product accumulator. */
     39 	unsigned       i,j,y,*x,s,r,c;
     40 
     41 	if ((int) px->fpclass < (int) py->fpclass) {
     42 		pt = px;
     43 		px = py;
     44 		py = pt;
     45 	}
     46 	/* Now class(x) >= class(y).  */
     47 
     48 	*pz = *px;
     49 	pz->sign = px->sign ^ py->sign;
     50 
     51 	switch (px->fpclass) {
     52 	case fp_quiet:
     53 	case fp_signaling:
     54 	case fp_zero:
     55 		return;
     56 	case fp_infinity:
     57 		if (py->fpclass == fp_zero) {
     58 			fpu_error_nan(pz);
     59 			pz->fpclass = fp_quiet;
     60 		}
     61 		return;
     62 	case fp_normal:
     63 		if (py->fpclass == fp_zero) {
     64 			pz->fpclass = fp_zero;
     65 			return;
     66 		}
     67 	}
     68 
     69 	/* Now x and y are both normal or subnormal. */
     70 
     71 	x = px->significand;	/* save typing */
     72 
     73 	s=r=acc[0]=acc[1]=acc[2]=acc[3]=0;	/* intialize acc to zero */
     74 
     75 	y = py->significand[3];		/* py->significand[3] * x */
     76 	if(y!=0) {
     77 	    j=1;
     78 	    do {
     79 		s |= r;		/* shift acc right one bit */
     80 		r  = acc[3]&1;
     81 		acc[3] = ((acc[2]&1)<<31)|(acc[3]>>1);
     82 		acc[2] = ((acc[1]&1)<<31)|(acc[2]>>1);
     83 		acc[1] = ((acc[0]&1)<<31)|(acc[1]>>1);
     84 		acc[0] = (acc[0]>>1);
     85 		if(j&y) {		/* bit i of y != 0, add x to acc */
     86 			c = 0;
     87 			c = fpu_add3wc(&acc[3],acc[3],x[3],c);
     88 			c = fpu_add3wc(&acc[2],acc[2],x[2],c);
     89 			c = fpu_add3wc(&acc[1],acc[1],x[1],c);
     90 			c = fpu_add3wc(&acc[0],acc[0],x[0],c);
     91 		}
     92 		j += j;
     93 	    } while (j!=0);
     94 	}
     95 
     96 	y = py->significand[2];		/* py->significand[2] * x */
     97 	if(y!=0) {
     98 	    j=1;
     99 	    do {
    100 		s |= r;		/* shift acc right one bit */
    101 		r  = acc[3]&1;
    102 		acc[3] = ((acc[2]&1)<<31)|(acc[3]>>1);
    103 		acc[2] = ((acc[1]&1)<<31)|(acc[2]>>1);
    104 		acc[1] = ((acc[0]&1)<<31)|(acc[1]>>1);
    105 		acc[0] = (acc[0]>>1);
    106 		if(j&y) {		/* bit i of y != 0, add x to acc */
    107 			c = 0;
    108 			c = fpu_add3wc(&acc[3],acc[3],x[3],c);
    109 			c = fpu_add3wc(&acc[2],acc[2],x[2],c);
    110 			c = fpu_add3wc(&acc[1],acc[1],x[1],c);
    111 			c = fpu_add3wc(&acc[0],acc[0],x[0],c);
    112 		}
    113 		j += j;
    114 	    } while (j!=0);
    115 	} else {
    116 		s |= r|(acc[3]&0x7fffffff);
    117 		r  = (acc[3]&0x80000000)>>31;
    118 		acc[3]=acc[2];acc[2]=acc[1];acc[1]=acc[0];acc[0]=0;
    119 	}
    120 
    121 	y = py->significand[1];		/* py->significand[1] * x */
    122 	if(y!=0) {
    123 	    j=1;
    124 	    do {
    125 		s |= r;		/* shift acc right one bit */
    126 		r  = acc[3]&1;
    127 		acc[3] = ((acc[2]&1)<<31)|(acc[3]>>1);
    128 		acc[2] = ((acc[1]&1)<<31)|(acc[2]>>1);
    129 		acc[1] = ((acc[0]&1)<<31)|(acc[1]>>1);
    130 		acc[0] = (acc[0]>>1);
    131 		if(j&y) {		/* bit i of y != 0, add x to acc */
    132 			c = 0;
    133 			c = fpu_add3wc(&acc[3],acc[3],x[3],c);
    134 			c = fpu_add3wc(&acc[2],acc[2],x[2],c);
    135 			c = fpu_add3wc(&acc[1],acc[1],x[1],c);
    136 			c = fpu_add3wc(&acc[0],acc[0],x[0],c);
    137 		}
    138 		j += j;
    139 	    } while (j!=0);
    140 	} else {
    141 		s |= r|(acc[3]&0x7fffffff);
    142 		r  = (acc[3]&0x80000000)>>31;
    143 		acc[3]=acc[2];acc[2]=acc[1];acc[1]=acc[0];acc[0]=0;
    144 	}
    145 
    146     					/* py->significand[0] * x */
    147 	y = py->significand[0];		/* y is of form 0x0001???? */
    148 	j=1;
    149 	do {
    150 		s |= r;		/* shift acc right one bit */
    151 		r  = acc[3]&1;
    152 		acc[3] = ((acc[2]&1)<<31)|(acc[3]>>1);
    153 		acc[2] = ((acc[1]&1)<<31)|(acc[2]>>1);
    154 		acc[1] = ((acc[0]&1)<<31)|(acc[1]>>1);
    155 		acc[0] = (acc[0]>>1);
    156 		if(j&y) {		/* bit i of y != 0, add x to acc */
    157 			c = 0;
    158 			c = fpu_add3wc(&acc[3],acc[3],x[3],c);
    159 			c = fpu_add3wc(&acc[2],acc[2],x[2],c);
    160 			c = fpu_add3wc(&acc[1],acc[1],x[1],c);
    161 			c = fpu_add3wc(&acc[0],acc[0],x[0],c);
    162 		}
    163 		j += j;
    164 	} while (j<=y);
    165 
    166 	if(acc[0]>=0x20000) {	/* right shift one bit to normalize */
    167 		pz->exponent = px->exponent + py->exponent + 1;
    168 		pz->sticky = s|r;
    169 		pz->rounded = acc[3]&1;
    170 		pz->significand[3]=((acc[2]&1)<<31)|(acc[3]>>1);
    171 		pz->significand[2]=((acc[1]&1)<<31)|(acc[2]>>1);
    172 		pz->significand[1]=((acc[0]&1)<<31)|(acc[1]>>1);
    173 		pz->significand[0]=(acc[0]>>1);
    174 	} else {
    175 		pz->exponent = px->exponent + py->exponent;
    176 		pz->sticky = s;
    177 		pz->rounded = r;
    178 		pz->significand[3]=acc[3];
    179 		pz->significand[2]=acc[2];
    180 		pz->significand[1]=acc[1];
    181 		pz->significand[0]=acc[0];
    182 	}
    183 }
    184