m68k/math-emu/fp_arith.c

31 struct fp_ext *fp_fabs(struct fp_ext *dest, struct fp_ext *src)
35 	fp_monadic_check(dest, src);
37 	dest->sign = 0;
39 	return dest;
42 struct fp_ext *fp_fneg(struct fp_ext *dest, struct fp_ext *src)
46 	fp_monadic_check(dest, src);
48 	dest->sign = !dest->sign;
50 	return dest;
58 struct fp_ext *fp_fadd(struct fp_ext *dest, struct fp_ext *src)
64 	fp_dyadic_check(dest, src);
66 	if (IS_INF(dest)) {
68 		if (IS_INF(src) && (src->sign != dest->sign))
69 			fp_set_nan(dest);
70 		return dest;
73 		fp_copy_ext(dest, src);
74 		return dest;
77 	if (IS_ZERO(dest)) {
79 			if (src->sign != dest->sign) {
81 					dest->sign = 1;
83 					dest->sign = 0;
86 			fp_copy_ext(dest, src);
87 		return dest;
90 	dest->lowmant = src->lowmant = 0;
92 	if ((diff = dest->exp - src->exp) > 0)
95 		fp_denormalize(dest, diff);
97 	if (dest->sign == src->sign) {
98 		if (fp_addmant(dest, src))
99 			if (!fp_addcarry(dest))
100 				return dest;
102 		if (dest->mant.m64 < src->mant.m64) {
103 			fp_submant(dest, src, dest);
104 			dest->sign = !dest->sign;
106 			fp_submant(dest, dest, src);
109 	return dest;
117 struct fp_ext *fp_fsub(struct fp_ext *dest, struct fp_ext *src)
122 	return fp_fadd(dest, src);
126 struct fp_ext *fp_fcmp(struct fp_ext *dest, struct fp_ext *src)
130 	FPDATA->temp[1] = *dest;
135 struct fp_ext *fp_ftst(struct fp_ext *dest, struct fp_ext *src)
139 	(void)dest;
144 struct fp_ext *fp_fmul(struct fp_ext *dest, struct fp_ext *src)
151 	fp_dyadic_check(dest, src);
154 	dest->sign = src->sign ^ dest->sign;
157 	if (IS_INF(dest)) {
159 			fp_set_nan(dest);
160 		return dest;
163 		if (IS_ZERO(dest))
164 			fp_set_nan(dest);
166 			fp_copy_ext(dest, src);
167 		return dest;
173 	if (IS_ZERO(dest) || IS_ZERO(src)) {
174 		dest->exp = 0;
175 		dest->mant.m64 = 0;
176 		dest->lowmant = 0;
178 		return dest;
181 	exp = dest->exp + src->exp - 0x3ffe;
186 	if ((long)dest->mant.m32[0] >= 0)
187 		exp -= fp_overnormalize(dest);
192 	fp_multiplymant(&temp, dest, src);
198 		fp_putmant128(dest, &temp, 1);
200 		fp_putmant128(dest, &temp, 0);
203 		fp_set_ovrflw(dest);
204 		return dest;
206 	dest->exp = exp;
209 		fp_denormalize(dest, -exp);
212 	return dest;
219    of src / dest, the result is actually dest / src. */
221 struct fp_ext *fp_fdiv(struct fp_ext *dest, struct fp_ext *src)
228 	fp_dyadic_check(dest, src);
231 	dest->sign = src->sign ^ dest->sign;
234 	if (IS_INF(dest)) {
237 			fp_set_nan(dest);
239 		return dest;
243 		dest->exp = 0;
244 		dest->mant.m64 = 0;
245 		dest->lowmant = 0;
247 		return dest;
251 	if (IS_ZERO(dest)) {
254 			fp_set_nan(dest);
256 		return dest;
261 		dest->exp = 0x7fff;
262 		dest->mant.m64 = 0;
264 		return dest;
267 	exp = dest->exp - src->exp + 0x3fff;
272 	if ((long)dest->mant.m32[0] >= 0)
273 		exp -= fp_overnormalize(dest);
278 	fp_dividemant(&temp, dest, src);
284 		fp_putmant128(dest, &temp, 32);
286 		fp_putmant128(dest, &temp, 31);
289 		fp_set_ovrflw(dest);
290 		return dest;
292 	dest->exp = exp;
295 		fp_denormalize(dest, -exp);
298 	return dest;
301 struct fp_ext *fp_fsglmul(struct fp_ext *dest, struct fp_ext *src)
307 	fp_dyadic_check(dest, src);
310 	dest->sign = src->sign ^ dest->sign;
313 	if (IS_INF(dest)) {
315 			fp_set_nan(dest);
316 		return dest;
319 		if (IS_ZERO(dest))
320 			fp_set_nan(dest);
322 			fp_copy_ext(dest, src);
323 		return dest;
329 	if (IS_ZERO(dest) || IS_ZERO(src)) {
330 		dest->exp = 0;
331 		dest->mant.m64 = 0;
332 		dest->lowmant = 0;
334 		return dest;
337 	exp = dest->exp + src->exp - 0x3ffe;
340 	fp_mul64(dest->mant.m32[0], dest->mant.m32[1],
341 		 dest->mant.m32[0] & 0xffffff00,
345 		fp_set_ovrflw(dest);
346 		return dest;
348 	dest->exp = exp;
351 		fp_denormalize(dest, -exp);
354 	return dest;
357 struct fp_ext *fp_fsgldiv(struct fp_ext *dest, struct fp_ext *src)
364 	fp_dyadic_check(dest, src);
367 	dest->sign = src->sign ^ dest->sign;
370 	if (IS_INF(dest)) {
373 			fp_set_nan(dest);
375 		return dest;
379 		dest->exp = 0;
380 		dest->mant.m64 = 0;
381 		dest->lowmant = 0;
383 		return dest;
387 	if (IS_ZERO(dest)) {
390 			fp_set_nan(dest);
392 		return dest;
397 		dest->exp = 0x7fff;
398 		dest->mant.m64 = 0;
400 		return dest;
403 	exp = dest->exp - src->exp + 0x3fff;
405 	dest->mant.m32[0] &= 0xffffff00;
409 	if (dest->mant.m32[0] >= src->mant.m32[0]) {
410 		fp_sub64(dest->mant, src->mant);
411 		fp_div64(quot, rem, dest->mant.m32[0], 0, src->mant.m32[0]);
412 		dest->mant.m32[0] = 0x80000000 | (quot >> 1);
413 		dest->mant.m32[1] = (quot & 1) | rem;	/* only for rounding */
415 		fp_div64(quot, rem, dest->mant.m32[0], 0, src->mant.m32[0]);
416 		dest->mant.m32[0] = quot;
417 		dest->mant.m32[1] = rem;		/* only for rounding */
422 		fp_set_ovrflw(dest);
423 		return dest;
425 	dest->exp = exp;
428 		fp_denormalize(dest, -exp);
431 	return dest;
440 static void fp_roundint(struct fp_ext *dest, int mode)
445 	if (!fp_normalize_ext(dest))
449 	if (IS_INF(dest) || IS_ZERO(dest))
453 	oldmant = dest->mant;
454 	switch (dest->exp) {
456 		dest->mant.m64 = 0;
459 		dest->mant.m32[0] &= 0xffffffffU << (0x401e - dest->exp);
460 		dest->mant.m32[1] = 0;
461 		if (oldmant.m64 == dest->mant.m64)
465 		dest->mant.m32[1] &= 0xffffffffU << (0x403e - dest->exp);
466 		if (oldmant.m32[1] == dest->mant.m32[1])
488 		switch (dest->exp) {
500 			mask = 1 << (0x401d - dest->exp);
505 			if (!(oldmant.m32[0] << (dest->exp - 0x3ffd)) &&
518 			mask = 1 << (0x403d - dest->exp);
523 			if (!(oldmant.m32[1] << (dest->exp - 0x401d)))
533 		if (dest->sign ^ (mode - FPCR_ROUND_RM))
538 	switch (dest->exp) {
540 		dest->exp = 0x3fff;
541 		dest->mant.m64 = 1ULL << 63;
544 		mask = 1 << (0x401e - dest->exp);
545 		if (dest->mant.m32[0] += mask)
547 		dest->mant.m32[0] = 0x80000000;
548 		dest->exp++;
551 		mask = 1 << (0x403e - dest->exp);
552 		if (dest->mant.m32[1] += mask)
554 		if (dest->mant.m32[0] += 1)
556 		dest->mant.m32[0] = 0x80000000;
557                 dest->exp++;
566 static struct fp_ext *modrem_kernel(struct fp_ext *dest, struct fp_ext *src,
571 	fp_dyadic_check(dest, src);
574 	if (IS_INF(dest) || IS_ZERO(src)) {
575 		fp_set_nan(dest);
576 		return dest;
578 	if (IS_ZERO(dest) || IS_INF(src))
579 		return dest;
582 	fp_copy_ext(&tmp, dest);
586 	fp_fsub(dest, &tmp);
589 	fp_set_quotient((dest->mant.m64 & 0x7f) | (dest->sign << 7));
590 	return dest;
598    fmod(src,dest) = (dest - (src * floor(dest / src))) */
600 struct fp_ext *fp_fmod(struct fp_ext *dest, struct fp_ext *src)
603 	return modrem_kernel(dest, src, FPCR_ROUND_RZ);
608    frem(src,dest) = (dest - (src * round(dest / src)))
611 struct fp_ext *fp_frem(struct fp_ext *dest, struct fp_ext *src)
614 	return modrem_kernel(dest, src, FPCR_ROUND_RN);
617 struct fp_ext *fp_fint(struct fp_ext *dest, struct fp_ext *src)
621 	fp_copy_ext(dest, src);
623 	fp_roundint(dest, FPDATA->rnd);
625 	return dest;
628 struct fp_ext *fp_fintrz(struct fp_ext *dest, struct fp_ext *src)
632 	fp_copy_ext(dest, src);
634 	fp_roundint(dest, FPCR_ROUND_RZ);
636 	return dest;
639 struct fp_ext *fp_fscale(struct fp_ext *dest, struct fp_ext *src)
645 	fp_dyadic_check(dest, src);
649 		fp_set_nan(dest);
650 		return dest;
652 	if (IS_INF(dest))
653 		return dest;
656 	if (IS_ZERO(src) || IS_ZERO(dest))
657 		return dest;
661 		fp_set_ovrflw(dest);
662 		return dest;
672 	scale += dest->exp;
675 		fp_set_ovrflw(dest);
678 		fp_denormalize(dest, -scale);
680 		dest->exp = scale;
682 	return dest;