| 123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471472473474475476477478479480481482483484485486487488489490491492493494495496497498499500501502503504505506507508509510511512513514515516517518519520521522523524525526527528529530531532533534535536537538539540541542543544545546547548549550551552553554555556557558559560 | /* *  linux/arch/arm/vfp/vfpdouble.c * * This code is derived in part from John R. Housers softfloat library, which * carries the following notice: * * =========================================================================== * This C source file is part of the SoftFloat IEC/IEEE Floating-point * Arithmetic Package, Release 2. * * Written by John R. Hauser.  This work was made possible in part by the * International Computer Science Institute, located at Suite 600, 1947 Center * Street, Berkeley, California 94704.  Funding was partially provided by the * National Science Foundation under grant MIP-9311980.  The original version * of this code was written as part of a project to build a fixed-point vector * processor in collaboration with the University of California at Berkeley, * overseen by Profs. Nelson Morgan and John Wawrzynek.  More information * is available through the web page `http://HTTP.CS.Berkeley.EDU/~jhauser/ * arithmetic/softfloat.html'. * * THIS SOFTWARE IS DISTRIBUTED AS IS, FOR FREE.  Although reasonable effort * has been made to avoid it, THIS SOFTWARE MAY CONTAIN FAULTS THAT WILL AT * TIMES RESULT IN INCORRECT BEHAVIOR.  USE OF THIS SOFTWARE IS RESTRICTED TO * PERSONS AND ORGANIZATIONS WHO CAN AND WILL TAKE FULL RESPONSIBILITY FOR ANY * AND ALL LOSSES, COSTS, OR OTHER PROBLEMS ARISING FROM ITS USE. * * Derivative works are acceptable, even for commercial purposes, so long as * (1) they include prominent notice that the work is derivative, and (2) they * include prominent notice akin to these three paragraphs for those parts of * this code that are retained. * =========================================================================== */#include <linux/kernel.h>#include <linux/bitops.h>#include <asm/div64.h>#include <asm/vfp.h>#include "vfpinstr.h"#include "vfp.h"static struct vfp_double vfp_double_default_qnan = {	.exponent	= 2047,	.sign		= 0,	.significand	= VFP_DOUBLE_SIGNIFICAND_QNAN,};static void vfp_double_dump(const char *str, struct vfp_double *d){	pr_debug("VFP: %s: sign=%d exponent=%d significand=%016llx\n",		 str, d->sign != 0, d->exponent, d->significand);}static void vfp_double_normalise_denormal(struct vfp_double *vd){	int bits = 31 - fls(vd->significand >> 32);	if (bits == 31)		bits = 63 - fls(vd->significand);	vfp_double_dump("normalise_denormal: in", vd);	if (bits) {		vd->exponent -= bits - 1;		vd->significand <<= bits;	}	vfp_double_dump("normalise_denormal: out", vd);}u32 vfp_double_normaliseround(int dd, struct vfp_double *vd, u32 fpscr, u32 exceptions, const char *func){	u64 significand, incr;	int exponent, shift, underflow;	u32 rmode;	vfp_double_dump("pack: in", vd);	/*	 * Infinities and NaNs are a special case.	 */	if (vd->exponent == 2047 && (vd->significand == 0 || exceptions))		goto pack;	/*	 * Special-case zero.	 */	if (vd->significand == 0) {		vd->exponent = 0;		goto pack;	}	exponent = vd->exponent;	significand = vd->significand;	shift = 32 - fls(significand >> 32);	if (shift == 32)		shift = 64 - fls(significand);	if (shift) {		exponent -= shift;		significand <<= shift;	}#ifdef DEBUG	vd->exponent = exponent;	vd->significand = significand;	vfp_double_dump("pack: normalised", vd);#endif	/*	 * Tiny number?	 */	underflow = exponent < 0;	if (underflow) {		significand = vfp_shiftright64jamming(significand, -exponent);		exponent = 0;#ifdef DEBUG		vd->exponent = exponent;		vd->significand = significand;		vfp_double_dump("pack: tiny number", vd);#endif		if (!(significand & ((1ULL << (VFP_DOUBLE_LOW_BITS + 1)) - 1)))			underflow = 0;	}	/*	 * Select rounding increment.	 */	incr = 0;	rmode = fpscr & FPSCR_RMODE_MASK;	if (rmode == FPSCR_ROUND_NEAREST) {		incr = 1ULL << VFP_DOUBLE_LOW_BITS;		if ((significand & (1ULL << (VFP_DOUBLE_LOW_BITS + 1))) == 0)			incr -= 1;	} else if (rmode == FPSCR_ROUND_TOZERO) {		incr = 0;	} else if ((rmode == FPSCR_ROUND_PLUSINF) ^ (vd->sign != 0))		incr = (1ULL << (VFP_DOUBLE_LOW_BITS + 1)) - 1;	pr_debug("VFP: rounding increment = 0x%08llx\n", incr);	/*	 * Is our rounding going to overflow?	 */	if ((significand + incr) < significand) {		exponent += 1;		significand = (significand >> 1) | (significand & 1);		incr >>= 1;#ifdef DEBUG		vd->exponent = exponent;		vd->significand = significand;		vfp_double_dump("pack: overflow", vd);#endif	}	/*	 * If any of the low bits (which will be shifted out of the	 * number) are non-zero, the result is inexact.	 */	if (significand & ((1 << (VFP_DOUBLE_LOW_BITS + 1)) - 1))		exceptions |= FPSCR_IXC;	/*	 * Do our rounding.	 */	significand += incr;	/*	 * Infinity?	 */	if (exponent >= 2046) {		exceptions |= FPSCR_OFC | FPSCR_IXC;		if (incr == 0) {			vd->exponent = 2045;			vd->significand = 0x7fffffffffffffffULL;		} else {			vd->exponent = 2047;		/* infinity */			vd->significand = 0;		}	} else {		if (significand >> (VFP_DOUBLE_LOW_BITS + 1) == 0)			exponent = 0;		if (exponent || significand > 0x8000000000000000ULL)			underflow = 0;		if (underflow)			exceptions |= FPSCR_UFC;		vd->exponent = exponent;		vd->significand = significand >> 1;	} pack:	vfp_double_dump("pack: final", vd);	{		s64 d = vfp_double_pack(vd);		pr_debug("VFP: %s: d(d%d)=%016llx exceptions=%08x\n", func,			 dd, d, exceptions);		vfp_put_double(d, dd);	}	return exceptions;}/* * Propagate the NaN, setting exceptions if it is signalling. * 'n' is always a NaN.  'm' may be a number, NaN or infinity. */static u32vfp_propagate_nan(struct vfp_double *vdd, struct vfp_double *vdn,		  struct vfp_double *vdm, u32 fpscr){	struct vfp_double *nan;	int tn, tm = 0;	tn = vfp_double_type(vdn);	if (vdm)		tm = vfp_double_type(vdm);	if (fpscr & FPSCR_DEFAULT_NAN)		/*		 * Default NaN mode - always returns a quiet NaN		 */		nan = &vfp_double_default_qnan;	else {		/*		 * Contemporary mode - select the first signalling		 * NAN, or if neither are signalling, the first		 * quiet NAN.		 */		if (tn == VFP_SNAN || (tm != VFP_SNAN && tn == VFP_QNAN))			nan = vdn;		else			nan = vdm;		/*		 * Make the NaN quiet.		 */		nan->significand |= VFP_DOUBLE_SIGNIFICAND_QNAN;	}	*vdd = *nan;	/*	 * If one was a signalling NAN, raise invalid operation.	 */	return tn == VFP_SNAN || tm == VFP_SNAN ? FPSCR_IOC : VFP_NAN_FLAG;}/* * Extended operations */static u32 vfp_double_fabs(int dd, int unused, int dm, u32 fpscr){	vfp_put_double(vfp_double_packed_abs(vfp_get_double(dm)), dd);	return 0;}static u32 vfp_double_fcpy(int dd, int unused, int dm, u32 fpscr){	vfp_put_double(vfp_get_double(dm), dd);	return 0;}static u32 vfp_double_fneg(int dd, int unused, int dm, u32 fpscr){	vfp_put_double(vfp_double_packed_negate(vfp_get_double(dm)), dd);	return 0;}static u32 vfp_double_fsqrt(int dd, int unused, int dm, u32 fpscr){	struct vfp_double vdm, vdd;	int ret, tm;	vfp_double_unpack(&vdm, vfp_get_double(dm));	tm = vfp_double_type(&vdm);	if (tm & (VFP_NAN|VFP_INFINITY)) {		struct vfp_double *vdp = &vdd;		if (tm & VFP_NAN)			ret = vfp_propagate_nan(vdp, &vdm, NULL, fpscr);		else if (vdm.sign == 0) { sqrt_copy:			vdp = &vdm;			ret = 0;		} else { sqrt_invalid:			vdp = &vfp_double_default_qnan;			ret = FPSCR_IOC;		}		vfp_put_double(vfp_double_pack(vdp), dd);		return ret;	}	/*	 * sqrt(+/- 0) == +/- 0	 */	if (tm & VFP_ZERO)		goto sqrt_copy;	/*	 * Normalise a denormalised number	 */	if (tm & VFP_DENORMAL)		vfp_double_normalise_denormal(&vdm);	/*	 * sqrt(<0) = invalid	 */	if (vdm.sign)		goto sqrt_invalid;	vfp_double_dump("sqrt", &vdm);	/*	 * Estimate the square root.	 */	vdd.sign = 0;	vdd.exponent = ((vdm.exponent - 1023) >> 1) + 1023;	vdd.significand = (u64)vfp_estimate_sqrt_significand(vdm.exponent, vdm.significand >> 32) << 31;	vfp_double_dump("sqrt estimate1", &vdd);	vdm.significand >>= 1 + (vdm.exponent & 1);	vdd.significand += 2 + vfp_estimate_div128to64(vdm.significand, 0, vdd.significand);	vfp_double_dump("sqrt estimate2", &vdd);	/*	 * And now adjust.	 */	if ((vdd.significand & VFP_DOUBLE_LOW_BITS_MASK) <= 5) {		if (vdd.significand < 2) {			vdd.significand = ~0ULL;		} else {			u64 termh, terml, remh, reml;			vdm.significand <<= 2;			mul64to128(&termh, &terml, vdd.significand, vdd.significand);			sub128(&remh, &reml, vdm.significand, 0, termh, terml);			while ((s64)remh < 0) {				vdd.significand -= 1;				shift64left(&termh, &terml, vdd.significand);				terml |= 1;				add128(&remh, &reml, remh, reml, termh, terml);			}			vdd.significand |= (remh | reml) != 0;		}	}	vdd.significand = vfp_shiftright64jamming(vdd.significand, 1);	return vfp_double_normaliseround(dd, &vdd, fpscr, 0, "fsqrt");}/* * Equal	:= ZC * Less than	:= N * Greater than	:= C * Unordered	:= CV */static u32 vfp_compare(int dd, int signal_on_qnan, int dm, u32 fpscr){	s64 d, m;	u32 ret = 0;	m = vfp_get_double(dm);	if (vfp_double_packed_exponent(m) == 2047 && vfp_double_packed_mantissa(m)) {		ret |= FPSCR_C | FPSCR_V;		if (signal_on_qnan || !(vfp_double_packed_mantissa(m) & (1ULL << (VFP_DOUBLE_MANTISSA_BITS - 1))))			/*			 * Signalling NaN, or signalling on quiet NaN			 */			ret |= FPSCR_IOC;	}	d = vfp_get_double(dd);	if (vfp_double_packed_exponent(d) == 2047 && vfp_double_packed_mantissa(d)) {		ret |= FPSCR_C | FPSCR_V;		if (signal_on_qnan || !(vfp_double_packed_mantissa(d) & (1ULL << (VFP_DOUBLE_MANTISSA_BITS - 1))))			/*			 * Signalling NaN, or signalling on quiet NaN			 */			ret |= FPSCR_IOC;	}	if (ret == 0) {		if (d == m || vfp_double_packed_abs(d | m) == 0) {			/*			 * equal			 */			ret |= FPSCR_Z | FPSCR_C;		} else if (vfp_double_packed_sign(d ^ m)) {			/*			 * different signs			 */			if (vfp_double_packed_sign(d))				/*				 * d is negative, so d < m				 */				ret |= FPSCR_N;			else				/*				 * d is positive, so d > m				 */				ret |= FPSCR_C;		} else if ((vfp_double_packed_sign(d) != 0) ^ (d < m)) {			/*			 * d < m			 */			ret |= FPSCR_N;		} else if ((vfp_double_packed_sign(d) != 0) ^ (d > m)) {			/*			 * d > m			 */			ret |= FPSCR_C;		}	}	return ret;}static u32 vfp_double_fcmp(int dd, int unused, int dm, u32 fpscr){	return vfp_compare(dd, 0, dm, fpscr);}static u32 vfp_double_fcmpe(int dd, int unused, int dm, u32 fpscr){	return vfp_compare(dd, 1, dm, fpscr);}static u32 vfp_double_fcmpz(int dd, int unused, int dm, u32 fpscr){	return vfp_compare(dd, 0, VFP_REG_ZERO, fpscr);}static u32 vfp_double_fcmpez(int dd, int unused, int dm, u32 fpscr){	return vfp_compare(dd, 1, VFP_REG_ZERO, fpscr);}static u32 vfp_double_fcvts(int sd, int unused, int dm, u32 fpscr){	struct vfp_double vdm;	struct vfp_single vsd;	int tm;	u32 exceptions = 0;	vfp_double_unpack(&vdm, vfp_get_double(dm));	tm = vfp_double_type(&vdm);	/*	 * If we have a signalling NaN, signal invalid operation.	 */	if (tm == VFP_SNAN)		exceptions = FPSCR_IOC;	if (tm & VFP_DENORMAL)		vfp_double_normalise_denormal(&vdm);	vsd.sign = vdm.sign;	vsd.significand = vfp_hi64to32jamming(vdm.significand);	/*	 * If we have an infinity or a NaN, the exponent must be 255	 */	if (tm & (VFP_INFINITY|VFP_NAN)) {		vsd.exponent = 255;		if (tm == VFP_QNAN)			vsd.significand |= VFP_SINGLE_SIGNIFICAND_QNAN;		goto pack_nan;	} else if (tm & VFP_ZERO)		vsd.exponent = 0;	else		vsd.exponent = vdm.exponent - (1023 - 127);	return vfp_single_normaliseround(sd, &vsd, fpscr, exceptions, "fcvts"); pack_nan:	vfp_put_float(vfp_single_pack(&vsd), sd);	return exceptions;}static u32 vfp_double_fuito(int dd, int unused, int dm, u32 fpscr){	struct vfp_double vdm;	u32 m = vfp_get_float(dm);	vdm.sign = 0;	vdm.exponent = 1023 + 63 - 1;	vdm.significand = (u64)m;	return vfp_double_normaliseround(dd, &vdm, fpscr, 0, "fuito");}static u32 vfp_double_fsito(int dd, int unused, int dm, u32 fpscr){	struct vfp_double vdm;	u32 m = vfp_get_float(dm);	vdm.sign = (m & 0x80000000) >> 16;	vdm.exponent = 1023 + 63 - 1;	vdm.significand = vdm.sign ? -m : m;	return vfp_double_normaliseround(dd, &vdm, fpscr, 0, "fsito");}static u32 vfp_double_ftoui(int sd, int unused, int dm, u32 fpscr){	struct vfp_double vdm;	u32 d, exceptions = 0;	int rmode = fpscr & FPSCR_RMODE_MASK;	int tm;	vfp_double_unpack(&vdm, vfp_get_double(dm));	/*	 * Do we have a denormalised number?	 */	tm = vfp_double_type(&vdm);	if (tm & VFP_DENORMAL)		exceptions |= FPSCR_IDC;	if (tm & VFP_NAN)		vdm.sign = 0;	if (vdm.exponent >= 1023 + 32) {		d = vdm.sign ? 0 : 0xffffffff;		exceptions = FPSCR_IOC;	} else if (vdm.exponent >= 1023 - 1) {		int shift = 1023 + 63 - vdm.exponent;		u64 rem, incr = 0;		/*		 * 2^0 <= m < 2^32-2^8		 */		d = (vdm.significand << 1) >> shift;		rem = vdm.significand << (65 - shift);		if (rmode == FPSCR_ROUND_NEAREST) {			incr = 0x8000000000000000ULL;			if ((d & 1) == 0)				incr -= 1;		} else if (rmode == FPSCR_ROUND_TOZERO) {			incr = 0;		} else if ((rmode == FPSCR_ROUND_PLUSINF) ^ (vdm.sign != 0)) {			incr = ~0ULL;		}		if ((rem + incr) < rem) {			if (d < 0xffffffff)				d += 1;			else				exceptions |= FPSCR_IOC;		}		if (d && vdm.sign) {			d = 0;			exceptions |= FPSCR_IOC;		} else if (rem)			exceptions |= FPSCR_IXC;	} else {
 |