Skip to content

Commit 39cb467

Browse files
committed
share internal_remquo for difference address space implementations
1 parent 30373c2 commit 39cb467

File tree

2 files changed

+173
-173
lines changed

2 files changed

+173
-173
lines changed

libclc/clc/lib/generic/math/clc_remquo.cl

Lines changed: 173 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -18,6 +18,179 @@
1818
#include <clc/math/math.h>
1919
#include <clc/shared/clc_max.h>
2020

21+
#define _sHighMask 0xfffff000u
22+
#define _iMaxQExp 0xbu
23+
// To prevent YLow to be denormal it should be checked
24+
// that Exp(Y) <= -127+23 (worst case when only last bit is non zero)
25+
// Exp(Y) < -103 -> Y < 0x0C000000
26+
// That value is used to construct _iYSub by setting up first bit to 1.
27+
// _iYCmp is get from max acceptable value 0x797fffff:
28+
// 0x797fffff - 0x8C000000 = 0x(1)ED7FFFFF
29+
#define _iYSub 0x8C000000u
30+
#define _iYCmp 0xED7FFFFFu
31+
#define _iOne 0x00000001u
32+
33+
static _CLC_INLINE int internal_remquo(float x, float y, private float *r,
34+
private uint *q) {
35+
uint signif_x, signif_y, rem_bit, quo_bit, tmp_x, tmp_y;
36+
int exp_x, exp_y, i, j;
37+
uint expabs_diff, special_op = 0, signbit_x, signbit_y, sign = 1;
38+
float result, abs_x, abs_y;
39+
float zero = 0.0f;
40+
int nRet = 0;
41+
// Remove sign bits
42+
tmp_x = ((*(int *)&x)) & EXSIGNBIT_SP32;
43+
tmp_y = ((*(int *)&y)) & EXSIGNBIT_SP32;
44+
signbit_x = (uint)((*(int *)&x) >> 31);
45+
signbit_y = (uint)((*(int *)&y) >> 31);
46+
if (signbit_x ^ signbit_y)
47+
sign = (-sign);
48+
// Get float absolute values
49+
abs_x = *(float *)&tmp_x;
50+
abs_y = *(float *)&tmp_y;
51+
// Remove exponent bias
52+
exp_x = (int)((tmp_x & (0x7F800000)) >> 23) - 127;
53+
exp_y = (int)((tmp_y & (0x7F800000)) >> 23) - 127;
54+
// Test for NaNs, Infs, and Zeros
55+
if ((exp_x == 0x00000080) || (exp_y == 0x00000080) || (tmp_x == 0) ||
56+
(tmp_y == 0))
57+
special_op++;
58+
// Get significands
59+
signif_x = (tmp_x & MANTBITS_SP32);
60+
signif_y = (tmp_y & MANTBITS_SP32);
61+
// Process NaNs, Infs, and Zeros
62+
if (special_op) {
63+
(*q) = 0;
64+
// x is NaN
65+
if ((signif_x != 0) && (exp_x == 0x00000080))
66+
result = x * 1.7f;
67+
// y is NaN
68+
else if ((signif_y != 0) && (exp_y == 0x00000080))
69+
result = y * 1.7f;
70+
// y is zero
71+
else if (abs_y == zero) {
72+
result = zero / zero;
73+
nRet = 1;
74+
}
75+
// x is zero
76+
else if (abs_x == zero)
77+
result = x;
78+
// x is Inf
79+
else if ((signif_x == 0) && (exp_x == 0x00000080))
80+
result = zero / zero;
81+
// y is Inf
82+
else
83+
result = x;
84+
(*r) = (result);
85+
return nRet;
86+
}
87+
// If x < y, fast return
88+
if (abs_x <= abs_y) {
89+
(*q) = 1 * sign;
90+
if (abs_x == abs_y) {
91+
(*r) = (zero * x);
92+
return nRet;
93+
}
94+
// Is x too big to scale up by 2.0f?
95+
if (exp_x != 127) {
96+
if ((2.0f * abs_x) <= abs_y) {
97+
(*q) = 0;
98+
(*r) = x;
99+
return nRet;
100+
}
101+
}
102+
result = abs_x - abs_y;
103+
if (signbit_x) {
104+
result = -result;
105+
}
106+
(*r) = (result);
107+
return nRet;
108+
}
109+
// Check for denormal x and y, adjust and normalize
110+
if ((exp_x == -127) && (signif_x != 0)) {
111+
exp_x = -126;
112+
while (signif_x <= MANTBITS_SP32) {
113+
exp_x--;
114+
signif_x <<= 1;
115+
};
116+
} else
117+
signif_x = (signif_x | (0x00800000L));
118+
if ((exp_y == -127) && (signif_y != 0)) {
119+
exp_y = -126;
120+
while (signif_y <= MANTBITS_SP32) {
121+
exp_y--;
122+
signif_y <<= 1;
123+
};
124+
} else
125+
signif_y = (signif_y | (0x00800000L));
126+
//
127+
// Main computational path
128+
//
129+
// Calculate exponent difference
130+
expabs_diff = (exp_x - exp_y) + 1;
131+
rem_bit = signif_x;
132+
quo_bit = 0;
133+
for (i = 0; i < expabs_diff; i++) {
134+
quo_bit = quo_bit << 1;
135+
if (rem_bit >= signif_y) {
136+
rem_bit -= signif_y;
137+
quo_bit++;
138+
}
139+
rem_bit <<= 1;
140+
}
141+
// Zero remquo ... return immediately with sign of x
142+
if (rem_bit == 0) {
143+
(*q) = ((uint)(0x7FFFFFFFL & quo_bit)) * sign;
144+
(*r) = (zero * x);
145+
return nRet;
146+
}
147+
// Adjust remquo
148+
rem_bit >>= 1;
149+
// Set exponent base, unbiased
150+
j = exp_y;
151+
// Calculate normalization shift
152+
while (rem_bit <= MANTBITS_SP32) {
153+
j--;
154+
rem_bit <<= 1;
155+
};
156+
// Prepare normal results
157+
if (j >= -126) {
158+
// Remove explicit 1
159+
rem_bit &= MANTBITS_SP32;
160+
// Set final exponent ... add exponent bias
161+
j = j + 127;
162+
}
163+
// Prepare denormal results
164+
else {
165+
// Determine denormalization shift count
166+
j = -j - 126;
167+
// Denormalization
168+
rem_bit >>= j;
169+
// Set final exponent ... denorms are 0
170+
j = 0;
171+
}
172+
rem_bit = (((uint)(j)) << 23) | rem_bit;
173+
// Create float result and adjust if >= .5 * divisor
174+
result = *(float *)&rem_bit;
175+
if ((2.0f * result) >= abs_y) {
176+
if ((2.0f * result) == abs_y) {
177+
if (quo_bit & 0x01) {
178+
result = -result;
179+
quo_bit++;
180+
}
181+
} else {
182+
result = result - abs_y;
183+
quo_bit++;
184+
}
185+
}
186+
// Final adjust for sign of input
187+
(*q) = ((uint)(0x7FFFFFFF & (quo_bit))) * sign;
188+
if (signbit_x)
189+
result = -result;
190+
(*r) = (result);
191+
return nRet;
192+
}
193+
21194
#define __CLC_ADDRESS_SPACE private
22195
#include <clc_remquo.inc>
23196
#undef __CLC_ADDRESS_SPACE

libclc/clc/lib/generic/math/clc_remquo.inc

Lines changed: 0 additions & 173 deletions
Original file line numberDiff line numberDiff line change
@@ -6,179 +6,6 @@
66
//
77
//===----------------------------------------------------------------------===//
88

9-
#define _sHighMask 0xfffff000u
10-
#define _iMaxQExp 0xbu
11-
// To prevent YLow to be denormal it should be checked
12-
// that Exp(Y) <= -127+23 (worst case when only last bit is non zero)
13-
// Exp(Y) < -103 -> Y < 0x0C000000
14-
// That value is used to construct _iYSub by setting up first bit to 1.
15-
// _iYCmp is get from max acceptable value 0x797fffff:
16-
// 0x797fffff - 0x8C000000 = 0x(1)ED7FFFFF
17-
#define _iYSub 0x8C000000u
18-
#define _iYCmp 0xED7FFFFFu
19-
#define _iOne 0x00000001u
20-
21-
static _CLC_INLINE _CLC_OVERLOAD int
22-
internal_remquo(float x, float y, float *r, __CLC_ADDRESS_SPACE uint *q) {
23-
uint signif_x, signif_y, rem_bit, quo_bit, tmp_x, tmp_y;
24-
int exp_x, exp_y, i, j;
25-
uint expabs_diff, special_op = 0, signbit_x, signbit_y, sign = 1;
26-
float result, abs_x, abs_y;
27-
float zero = 0.0f;
28-
int nRet = 0;
29-
// Remove sign bits
30-
tmp_x = ((*(int *)&x)) & EXSIGNBIT_SP32;
31-
tmp_y = ((*(int *)&y)) & EXSIGNBIT_SP32;
32-
signbit_x = (uint)((*(int *)&x) >> 31);
33-
signbit_y = (uint)((*(int *)&y) >> 31);
34-
if (signbit_x ^ signbit_y)
35-
sign = (-sign);
36-
// Get float absolute values
37-
abs_x = *(float *)&tmp_x;
38-
abs_y = *(float *)&tmp_y;
39-
// Remove exponent bias
40-
exp_x = (int)((tmp_x & (0x7F800000)) >> 23) - 127;
41-
exp_y = (int)((tmp_y & (0x7F800000)) >> 23) - 127;
42-
// Test for NaNs, Infs, and Zeros
43-
if ((exp_x == 0x00000080) || (exp_y == 0x00000080) || (tmp_x == 0) ||
44-
(tmp_y == 0))
45-
special_op++;
46-
// Get significands
47-
signif_x = (tmp_x & MANTBITS_SP32);
48-
signif_y = (tmp_y & MANTBITS_SP32);
49-
// Process NaNs, Infs, and Zeros
50-
if (special_op) {
51-
(*q) = 0;
52-
// x is NaN
53-
if ((signif_x != 0) && (exp_x == 0x00000080))
54-
result = x * 1.7f;
55-
// y is NaN
56-
else if ((signif_y != 0) && (exp_y == 0x00000080))
57-
result = y * 1.7f;
58-
// y is zero
59-
else if (abs_y == zero) {
60-
result = zero / zero;
61-
nRet = 1;
62-
}
63-
// x is zero
64-
else if (abs_x == zero)
65-
result = x;
66-
// x is Inf
67-
else if ((signif_x == 0) && (exp_x == 0x00000080))
68-
result = zero / zero;
69-
// y is Inf
70-
else
71-
result = x;
72-
(*r) = (result);
73-
return nRet;
74-
}
75-
// If x < y, fast return
76-
if (abs_x <= abs_y) {
77-
(*q) = 1 * sign;
78-
if (abs_x == abs_y) {
79-
(*r) = (zero * x);
80-
return nRet;
81-
}
82-
// Is x too big to scale up by 2.0f?
83-
if (exp_x != 127) {
84-
if ((2.0f * abs_x) <= abs_y) {
85-
(*q) = 0;
86-
(*r) = x;
87-
return nRet;
88-
}
89-
}
90-
result = abs_x - abs_y;
91-
if (signbit_x) {
92-
result = -result;
93-
}
94-
(*r) = (result);
95-
return nRet;
96-
}
97-
// Check for denormal x and y, adjust and normalize
98-
if ((exp_x == -127) && (signif_x != 0)) {
99-
exp_x = -126;
100-
while (signif_x <= MANTBITS_SP32) {
101-
exp_x--;
102-
signif_x <<= 1;
103-
};
104-
} else
105-
signif_x = (signif_x | (0x00800000L));
106-
if ((exp_y == -127) && (signif_y != 0)) {
107-
exp_y = -126;
108-
while (signif_y <= MANTBITS_SP32) {
109-
exp_y--;
110-
signif_y <<= 1;
111-
};
112-
} else
113-
signif_y = (signif_y | (0x00800000L));
114-
//
115-
// Main computational path
116-
//
117-
// Calculate exponent difference
118-
expabs_diff = (exp_x - exp_y) + 1;
119-
rem_bit = signif_x;
120-
quo_bit = 0;
121-
for (i = 0; i < expabs_diff; i++) {
122-
quo_bit = quo_bit << 1;
123-
if (rem_bit >= signif_y) {
124-
rem_bit -= signif_y;
125-
quo_bit++;
126-
}
127-
rem_bit <<= 1;
128-
}
129-
// Zero remquo ... return immediately with sign of x
130-
if (rem_bit == 0) {
131-
(*q) = ((uint)(0x7FFFFFFFL & quo_bit)) * sign;
132-
(*r) = (zero * x);
133-
return nRet;
134-
}
135-
// Adjust remquo
136-
rem_bit >>= 1;
137-
// Set exponent base, unbiased
138-
j = exp_y;
139-
// Calculate normalization shift
140-
while (rem_bit <= MANTBITS_SP32) {
141-
j--;
142-
rem_bit <<= 1;
143-
};
144-
// Prepare normal results
145-
if (j >= -126) {
146-
// Remove explicit 1
147-
rem_bit &= MANTBITS_SP32;
148-
// Set final exponent ... add exponent bias
149-
j = j + 127;
150-
}
151-
// Prepare denormal results
152-
else {
153-
// Determine denormalization shift count
154-
j = -j - 126;
155-
// Denormalization
156-
rem_bit >>= j;
157-
// Set final exponent ... denorms are 0
158-
j = 0;
159-
}
160-
rem_bit = (((uint)(j)) << 23) | rem_bit;
161-
// Create float result and adjust if >= .5 * divisor
162-
result = *(float *)&rem_bit;
163-
if ((2.0f * result) >= abs_y) {
164-
if ((2.0f * result) == abs_y) {
165-
if (quo_bit & 0x01) {
166-
result = -result;
167-
quo_bit++;
168-
}
169-
} else {
170-
result = result - abs_y;
171-
quo_bit++;
172-
}
173-
}
174-
// Final adjust for sign of input
175-
(*q) = ((uint)(0x7FFFFFFF & (quo_bit))) * sign;
176-
if (signbit_x)
177-
result = -result;
178-
(*r) = (result);
179-
return nRet;
180-
}
181-
1829
_CLC_DEF _CLC_OVERLOAD float __clc_remquo(float x, float y,
18310
__CLC_ADDRESS_SPACE int *quo) {
18411
float vr1;

0 commit comments

Comments
 (0)