1 /*
2  * Copyright (C) 2015 The Android Open Source Project
3  *
4  * Licensed under the Apache License, Version 2.0 (the "License");
5  * you may not use this file except in compliance with the License.
6  * You may obtain a copy of the License at
7  *
8  *      http://www.apache.org/licenses/LICENSE-2.0
9  *
10  * Unless required by applicable law or agreed to in writing, software
11  * distributed under the License is distributed on an "AS IS" BASIS,
12  * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
13  * See the License for the specific language governing permissions and
14  * limitations under the License.
15  */
16 
17 #include <stdlib.h>
18 #include <stdio.h>
19 #include <math.h>
20 
21 // This utility generates a Java file containing data used in the CTS test for
22 // FP16 arithmetic.  The Java file containes a class with the following fields:
23 //     * A 1D array of length 'n' containing various constants used in the test
24 //     * Four n x n x 2 arrays, each containing the reference output for
25 //     pair-wise addition, subtraction, multiplication and division.  The
26 //     reference output is a range accounting for tolerable error.  The
27 //     acceptable error is 3 x ULP for division and 1 x ULP for other
28 //     operations.
29 
30 typedef __fp16 half;
31 
32 // Macros for names of the package, class and fields in the generated java file
33 #define PACKAGE_NAME "android.renderscript.cts"
34 #define CLASS_NAME "Float16TestData"
35 #define INPUT_ARRAY "input"
36 #define OUTPUT_ARRAY_ADD "ReferenceOutputForAdd"
37 #define OUTPUT_ARRAY_SUB "ReferenceOutputForSub"
38 #define OUTPUT_ARRAY_MUL "ReferenceOutputForMul"
39 #define OUTPUT_ARRAY_DIV "ReferenceOutputForDiv"
40 
41 // Structure to hold an FP16 constant and its human-readable description, to be
42 // added as a comment, in the generated Java file
43 typedef struct {
44   unsigned short value;
45   const char* description;
46 } FP16Constant;
47 
48 FP16Constant input[] = {
49     { 0b0011110000000000, "one" },
50     { 0b0100000000000000, "two" },
51     { 0b0000000000000001, "smallest subnormal" },
52     { 0b0000001111111111, "largest subnormal" },
53     { 0b0000010000000000, "smallest normal" },
54     { 0b0111101111111111, "largest normal" },
55     { 0x3880, "0.562500" },
56     { 0x3e80, "1.625000" },
57     { 0x5140, "42.000000" },
58     { 0x5ac0, "216.000000" },
59     { 0x6c75, "4564.000000" },
60     { 0x7b53, "60000.000000" },
61     { 0b1011110000000000, "negative one" },
62     { 0b1100000000000000, "negative two" },
63     { 0b1000000000000001, "negative (smallest subnormal)" },
64     { 0b1000001111111111, "negative (largest subnormal)" },
65     { 0b1000010000000000, "negative (smallest normal)" },
66     { 0b1111101111111111, "negative (largest normal)" },
67     { 0xb880, "-0.562500" },
68     { 0xbe80, "-1.625000" },
69     { 0xd140, "-42.000000" },
70     { 0xdac0, "-216.000000" },
71     { 0xec75, "-4564.000000" },
72     { 0xfb53, "-60000.000000" },
73     { 0b0000000000000000, "zero" },
74     { 0b0111110000000000, "infinity" },
75     { 0b1000000000000000, "negative zero" },
76     { 0b1111110000000000, "negative infinity" },
77     { 0b0111110000000001, "nan" },
78 };
79 
80 const int numInputs = sizeof(input) / sizeof(FP16Constant);
81 
82 // 16-bit masks for extracting sign, exponent and mantissa bits
83 static unsigned short SIGN_MASK     = 0x8000;
84 static unsigned short EXPONENT_MASK = 0x7C00;
85 static unsigned short MANTISSA_MASK = 0x03FF;
86 
87 // NaN has all exponent bits set to 1 and a non-zero mantissa
isFloat16NaN(unsigned short val)88 int isFloat16NaN(unsigned short val) {
89   return (val & EXPONENT_MASK) == EXPONENT_MASK &&
90          (val & MANTISSA_MASK) != 0;
91 }
92 
93 // Infinity has all exponent bits set to 1 and zeroes in mantissa
isFloat16Infinite(unsigned short val)94 int isFloat16Infinite(unsigned short val) {
95   return (val & EXPONENT_MASK) == EXPONENT_MASK &&
96          (val & MANTISSA_MASK) == 0;
97 }
98 
99 // Subnormal numbers have exponent bits set to 0 and a non-zero mantissa
isFloat16SubNormal(unsigned short val)100 int isFloat16SubNormal(unsigned short val) {
101     return (val & EXPONENT_MASK) == 0 && (val & MANTISSA_MASK) != 0;
102 }
103 
104 // Negativity test checks the sign bit
isFloat16Negative(unsigned short val)105 int isFloat16Negative(unsigned short val) {
106     return (val & SIGN_MASK) != 0;
107 }
108 
109 // Interpret a short as a FP16 value and convert to float
half2float(unsigned short s)110 float half2float(unsigned short s) {
111   half h = *(half *) &s;
112   return (float) h;
113 }
114 
115 // Return the short value representing a float value in FP16
float2half(float f)116 unsigned short float2half(float f) {
117   half h = (half) f;
118   return *(unsigned short *) &h;
119 }
120 
121 // Compute ULP for 'value' and store value +/- tolerance * ULP in bounds sarray
getErrorBar(unsigned short value,int tolerance,unsigned short bounds[2])122 void getErrorBar(unsigned short value, int tolerance, unsigned short bounds[2]) {
123   // Validate 'tolerance' parameter
124   if (tolerance != 1 && tolerance != 3) {
125     fprintf(stderr, "Allowed ULP error should either be 1 or 3, and not %d\n",
126             tolerance);
127     exit(0);
128   }
129 
130   half hValue = *(half *) &value;
131   half ulp;
132 
133   // For Infinity and NaN, bounds are equal to 'value'
134   if (isFloat16Infinite(value) || isFloat16NaN(value)) {
135     bounds[0] = value;
136     bounds[1] = value;
137     return;
138   }
139 
140   // Compute ULP
141   if (isFloat16SubNormal(value)) {
142     // 1 ulp for a subnormal number is the smallest possible subnormal
143     unsigned short ulpInShort = 0b0000000000000001;
144     ulp = *(half *) &ulpInShort;
145   }
146   else {
147     // 1 ulp for a non-subnormal number is (b - a) where
148     //   - a has same exponent as 'value', zeroes for sign and mantissa
149     //   - b has same exponent and sign as 'a', and has '1' in the mantissa
150     // (b - a) gives the ULP by getting rid of the implied '1' at the front of
151     // the mantissa
152     unsigned short a = (value & EXPONENT_MASK);
153     unsigned short b = (a | 1);
154     half hA = *(half *) &a;
155     half hB = *(half *) &b;
156     ulp = hB - hA;
157   }
158 
159   // Compute error bar based on error tolerance
160   half lb = hValue - tolerance * ulp;
161   half ub = hValue + tolerance * ulp;
162   if (lb > ub) {
163     fprintf(stderr, "Warning! inconsistency in bounds\n");
164     fprintf(stderr, "Value: %f, ulp: %f\n", (float) hValue, (float) ulp);
165     fprintf(stderr, "lb: %f ub: %f\n", (float) lb, (float) ub);
166     fprintf(stderr, "lb: %x ub: %x\n", *(unsigned short *) &lb, *(unsigned short *) &ub);
167   }
168 
169   // Set the bounds
170   bounds[0] = *(unsigned short *) &lb;
171   bounds[1] = *(unsigned short *) &ub;
172 
173   // RS allows flush-to-zero for sub-normal results in relaxed precision.
174   // Flush lower bound of a positive sub-normal result to zero.
175   if (!isFloat16Negative(bounds[0]) && isFloat16SubNormal(bounds[0]))
176     bounds[0] = 0x0;
177   // Flush upper bound of a negative sub-normal result to negative zero.
178   if (isFloat16Negative(bounds[1]) && isFloat16SubNormal(bounds[1]))
179     bounds[1] = 0x0 | SIGN_MASK;
180 
181 }
182 
183 // Utilities that take 'unsigned short' representations of two fp16 values and
184 // return the result of an arithmetic operation as an 'unsigned short'.
185 typedef unsigned short operation_t(unsigned short, unsigned short);
186 
add(unsigned short a,unsigned short b)187 unsigned short add(unsigned short a, unsigned short b) {
188   float op1 = half2float(a);
189   float op2 = half2float(b);
190   return float2half(op1 + op2);
191 }
192 
subtract(unsigned short a,unsigned short b)193 unsigned short subtract(unsigned short a, unsigned short b) {
194   float op1 = half2float(a);
195   float op2 = half2float(b);
196   return float2half(op1 - op2);
197 }
198 
multiply(unsigned short a,unsigned short b)199 unsigned short multiply(unsigned short a, unsigned short b) {
200   float op1 = half2float(a);
201   float op2 = half2float(b);
202   return float2half(op1 * op2);
203 }
204 
divide(unsigned short a,unsigned short b)205 unsigned short divide(unsigned short a, unsigned short b) {
206   float op1 = half2float(a);
207   float op2 = half2float(b);
208   return float2half(op1 / op2);
209 }
210 
211 // Print Java code that initializes the input array (along with the description
212 // of the constant as a comment)
printInput()213 void printInput() {
214   printf("static short[] %s = {\n", INPUT_ARRAY);
215 
216   for (int x = 0; x < numInputs; x ++)
217     printf("(short) 0x%04x, // %s\n", input[x].value, input[x].description);
218 
219   printf("};\n\n");
220 }
221 
222 // Print Java code that initializes the output array with the acceptable bounds
223 // on the output.  For each pair of inputs, bounds are calculated on the result
224 // from applying 'operation' on the pair.
printReferenceOutput(const char * fieldName,operation_t operation,int tolerance)225 void printReferenceOutput(const char *fieldName, operation_t operation,
226                           int tolerance) {
227   unsigned short result;
228   unsigned short resultBounds[2];
229 
230   printf("static short[][][] %s = {\n", fieldName);
231 
232   for (int x = 0; x < numInputs; x ++) {
233     printf("{");
234     for (int y = 0; y < numInputs; y ++) {
235       // Apply 'operation' and compute error bounds for the result.
236       result = operation(input[x].value, input[y].value);
237       getErrorBar(result, tolerance, resultBounds);
238 
239       printf("{ (short) 0x%04x, (short) 0x%04x},", resultBounds[0],
240                                                    resultBounds[1]);
241     }
242     printf("},\n");
243   }
244 
245   printf("};\n\n");
246 }
247 
248 const char *preamble = "/*\n"
249 " * Copyright (C) 2015 The Android Open Source Project\n"
250 " *\n"
251 " * Licensed under the Apache License, Version 2.0 (the \"License\");\n"
252 " * you may not use this file except in compliance with the License.\n"
253 " * You may obtain a copy of the License at\n"
254 " *\n"
255 " *      http://www.apache.org/licenses/LICENSE-2.0\n"
256 " *\n"
257 " * Unless required by applicable law or agreed to in writing, software\n"
258 " * distributed under the License is distributed on an \"AS IS\" BASIS,\n"
259 " * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.\n"
260 " * See the License for the specific language governing permissions and\n"
261 " * limitations under the License.\n"
262 " */\n"
263 "\n"
264 "/* Don't edit this file!  It is auto-generated by float16_gen.sh */\n\n"
265 "package "PACKAGE_NAME";\n\n"
266 "public class "CLASS_NAME" {\n";
267 
main()268 int main() {
269   // Print a preamble with copyright and class declaration, followed by the
270   // input FP16 array, and reference outputs for pair-wise arithmetic
271   // operations.
272   printf("%s", preamble);
273   printInput();
274 
275   printReferenceOutput(OUTPUT_ARRAY_ADD, add, 1);
276   printReferenceOutput(OUTPUT_ARRAY_SUB, subtract, 1);
277   printReferenceOutput(OUTPUT_ARRAY_MUL, multiply, 1);
278   printReferenceOutput(OUTPUT_ARRAY_DIV, divide, 3);
279 
280   printf("}");
281 }
282