Book a Demo!
CoCalc Logo Icon
StoreFeaturesDocsShareSupportNewsAboutPoliciesSign UpSign In
torvalds
GitHub Repository: torvalds/linux
 Path: blob/master/arch/mips/math-emu/sp_mul.c
26442 views
1
// SPDX-License-Identifier: GPL-2.0-only

2
/* IEEE754 floating point arithmetic

3
 * single precision

4
 */

5
/*

6
 * MIPS floating point support

7
 * Copyright (C) 1994-2000 Algorithmics Ltd.

8
 */

9


10
#include "ieee754sp.h"

11


12
union ieee754sp ieee754sp_mul(union ieee754sp x, union ieee754sp y)

13
{

14
	int re;

15
	int rs;

16
	unsigned int rm;

17
	unsigned short lxm;

18
	unsigned short hxm;

19
	unsigned short lym;

20
	unsigned short hym;

21
	unsigned int lrm;

22
	unsigned int hrm;

23
	unsigned int t;

24
	unsigned int at;

25


26
	COMPXSP;

27
	COMPYSP;

28


29
	EXPLODEXSP;

30
	EXPLODEYSP;

31


32
	ieee754_clearcx();

33


34
	FLUSHXSP;

35
	FLUSHYSP;

36


37
	switch (CLPAIR(xc, yc)) {

38
	case CLPAIR(IEEE754_CLASS_QNAN, IEEE754_CLASS_SNAN):

39
	case CLPAIR(IEEE754_CLASS_ZERO, IEEE754_CLASS_SNAN):

40
	case CLPAIR(IEEE754_CLASS_NORM, IEEE754_CLASS_SNAN):

41
	case CLPAIR(IEEE754_CLASS_DNORM, IEEE754_CLASS_SNAN):

42
	case CLPAIR(IEEE754_CLASS_INF, IEEE754_CLASS_SNAN):

43
		return ieee754sp_nanxcpt(y);

44


45
	case CLPAIR(IEEE754_CLASS_SNAN, IEEE754_CLASS_SNAN):

46
	case CLPAIR(IEEE754_CLASS_SNAN, IEEE754_CLASS_QNAN):

47
	case CLPAIR(IEEE754_CLASS_SNAN, IEEE754_CLASS_ZERO):

48
	case CLPAIR(IEEE754_CLASS_SNAN, IEEE754_CLASS_NORM):

49
	case CLPAIR(IEEE754_CLASS_SNAN, IEEE754_CLASS_DNORM):

50
	case CLPAIR(IEEE754_CLASS_SNAN, IEEE754_CLASS_INF):

51
		return ieee754sp_nanxcpt(x);

52


53
	case CLPAIR(IEEE754_CLASS_ZERO, IEEE754_CLASS_QNAN):

54
	case CLPAIR(IEEE754_CLASS_NORM, IEEE754_CLASS_QNAN):

55
	case CLPAIR(IEEE754_CLASS_DNORM, IEEE754_CLASS_QNAN):

56
	case CLPAIR(IEEE754_CLASS_INF, IEEE754_CLASS_QNAN):

57
		return y;

58


59
	case CLPAIR(IEEE754_CLASS_QNAN, IEEE754_CLASS_QNAN):

60
	case CLPAIR(IEEE754_CLASS_QNAN, IEEE754_CLASS_ZERO):

61
	case CLPAIR(IEEE754_CLASS_QNAN, IEEE754_CLASS_NORM):

62
	case CLPAIR(IEEE754_CLASS_QNAN, IEEE754_CLASS_DNORM):

63
	case CLPAIR(IEEE754_CLASS_QNAN, IEEE754_CLASS_INF):

64
		return x;

65


66


67
	/*

68
	 * Infinity handling

69
	 */

70
	case CLPAIR(IEEE754_CLASS_INF, IEEE754_CLASS_ZERO):

71
	case CLPAIR(IEEE754_CLASS_ZERO, IEEE754_CLASS_INF):

72
		ieee754_setcx(IEEE754_INVALID_OPERATION);

73
		return ieee754sp_indef();

74


75
	case CLPAIR(IEEE754_CLASS_NORM, IEEE754_CLASS_INF):

76
	case CLPAIR(IEEE754_CLASS_DNORM, IEEE754_CLASS_INF):

77
	case CLPAIR(IEEE754_CLASS_INF, IEEE754_CLASS_NORM):

78
	case CLPAIR(IEEE754_CLASS_INF, IEEE754_CLASS_DNORM):

79
	case CLPAIR(IEEE754_CLASS_INF, IEEE754_CLASS_INF):

80
		return ieee754sp_inf(xs ^ ys);

81


82
	case CLPAIR(IEEE754_CLASS_ZERO, IEEE754_CLASS_ZERO):

83
	case CLPAIR(IEEE754_CLASS_ZERO, IEEE754_CLASS_NORM):

84
	case CLPAIR(IEEE754_CLASS_ZERO, IEEE754_CLASS_DNORM):

85
	case CLPAIR(IEEE754_CLASS_NORM, IEEE754_CLASS_ZERO):

86
	case CLPAIR(IEEE754_CLASS_DNORM, IEEE754_CLASS_ZERO):

87
		return ieee754sp_zero(xs ^ ys);

88


89


90
	case CLPAIR(IEEE754_CLASS_DNORM, IEEE754_CLASS_DNORM):

91
		SPDNORMX;

92
		fallthrough;

93
	case CLPAIR(IEEE754_CLASS_NORM, IEEE754_CLASS_DNORM):

94
		SPDNORMY;

95
		break;

96


97
	case CLPAIR(IEEE754_CLASS_DNORM, IEEE754_CLASS_NORM):

98
		SPDNORMX;

99
		break;

100


101
	case CLPAIR(IEEE754_CLASS_NORM, IEEE754_CLASS_NORM):

102
		break;

103
	}

104
	/* rm = xm * ym, re = xe+ye basically */

105
	assert(xm & SP_HIDDEN_BIT);

106
	assert(ym & SP_HIDDEN_BIT);

107


108
	re = xe + ye;

109
	rs = xs ^ ys;

110


111
	/* shunt to top of word */

112
	xm <<= 32 - (SP_FBITS + 1);

113
	ym <<= 32 - (SP_FBITS + 1);

114


115
	/*

116
	 * Multiply 32 bits xm, ym to give high 32 bits rm with stickness.

117
	 */

118
	lxm = xm & 0xffff;

119
	hxm = xm >> 16;

120
	lym = ym & 0xffff;

121
	hym = ym >> 16;

122


123
	lrm = lxm * lym;	/* 16 * 16 => 32 */

124
	hrm = hxm * hym;	/* 16 * 16 => 32 */

125


126
	t = lxm * hym; /* 16 * 16 => 32 */

127
	at = lrm + (t << 16);

128
	hrm += at < lrm;

129
	lrm = at;

130
	hrm = hrm + (t >> 16);

131


132
	t = hxm * lym; /* 16 * 16 => 32 */

133
	at = lrm + (t << 16);

134
	hrm += at < lrm;

135
	lrm = at;

136
	hrm = hrm + (t >> 16);

137


138
	rm = hrm | (lrm != 0);

139


140
	/*

141
	 * Sticky shift down to normal rounding precision.

142
	 */

143
	if ((int) rm < 0) {

144
		rm = (rm >> (32 - (SP_FBITS + 1 + 3))) |

145
		    ((rm << (SP_FBITS + 1 + 3)) != 0);

146
		re++;

147
	} else {

148
		rm = (rm >> (32 - (SP_FBITS + 1 + 3 + 1))) |

149
		     ((rm << (SP_FBITS + 1 + 3 + 1)) != 0);

150
	}

151
	assert(rm & (SP_HIDDEN_BIT << 3));

152


153
	return ieee754sp_format(rs, re, rm);

154
}

155


156