1 /*
2  * Copyright 2015 Advanced Micro Devices, Inc.
3  *
4  * Permission is hereby granted, free of charge, to any person obtaining a
5  * copy of this software and associated documentation files (the "Software"),
6  * to deal in the Software without restriction, including without limitation
7  * the rights to use, copy, modify, merge, publish, distribute, sublicense,
8  * and/or sell copies of the Software, and to permit persons to whom the
9  * Software is furnished to do so, subject to the following conditions:
10  *
11  * The above copyright notice and this permission notice shall be included in
12  * all copies or substantial portions of the Software.
13  *
14  * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
15  * IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
16  * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT.  IN NO EVENT SHALL
17  * THE COPYRIGHT HOLDER(S) OR AUTHOR(S) BE LIABLE FOR ANY CLAIM, DAMAGES OR
18  * OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE,
19  * ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR
20  * OTHER DEALINGS IN THE SOFTWARE.
21  *
22  */
23 #include <asm/div64.h>
24 
25 enum ppevvmath_constants {
26 	/* We multiply all original integers with 2^SHIFT_AMOUNT to get the fInt representation */
27 	SHIFT_AMOUNT	= 16,
28 
29 	/* Change this value to change the number of decimal places in the final output - 5 is a good default */
30 	PRECISION	=  5,
31 
32 	SHIFTED_2	= (2 << SHIFT_AMOUNT),
33 
34 	/* 32767 - Might change in the future */
35 	MAX		= (1 << (SHIFT_AMOUNT - 1)) - 1,
36 };
37 
38 /* -------------------------------------------------------------------------------
39  * NEW TYPE - fINT
40  * -------------------------------------------------------------------------------
41  * A variable of type fInt can be accessed in 3 ways using the dot (.) operator
42  * fInt A;
43  * A.full => The full number as it is. Generally not easy to read
44  * A.partial.real => Only the integer portion
45  * A.partial.decimal => Only the fractional portion
46  */
47 typedef union _fInt {
48     int full;
49     struct _partial {
50         unsigned int decimal: SHIFT_AMOUNT; /*Needs to always be unsigned*/
51         int real: 32 - SHIFT_AMOUNT;
52     } partial;
53 } fInt;
54 
55 /* -------------------------------------------------------------------------------
56  * Function Declarations
57  *  -------------------------------------------------------------------------------
58  */
59 static fInt ConvertToFraction(int);                       /* Use this to convert an INT to a FINT */
60 static fInt Convert_ULONG_ToFraction(uint32_t);           /* Use this to convert an uint32_t to a FINT */
61 static fInt GetScaledFraction(int, int);                  /* Use this to convert an INT to a FINT after scaling it by a factor */
62 static int ConvertBackToInteger(fInt);                    /* Convert a FINT back to an INT that is scaled by 1000 (i.e. last 3 digits are the decimal digits) */
63 
64 static fInt fNegate(fInt);                                /* Returns -1 * input fInt value */
65 static fInt fAdd (fInt, fInt);                            /* Returns the sum of two fInt numbers */
66 static fInt fSubtract (fInt A, fInt B);                   /* Returns A-B - Sometimes easier than Adding negative numbers */
67 static fInt fMultiply (fInt, fInt);                       /* Returns the product of two fInt numbers */
68 static fInt fDivide (fInt A, fInt B);                     /* Returns A/B */
69 static fInt fGetSquare(fInt);                             /* Returns the square of a fInt number */
70 static fInt fSqrt(fInt);                                  /* Returns the Square Root of a fInt number */
71 
72 static int uAbs(int);                                     /* Returns the Absolute value of the Int */
73 static int uPow(int base, int exponent);                  /* Returns base^exponent an INT */
74 
75 static void SolveQuadracticEqn(fInt, fInt, fInt, fInt[]); /* Returns the 2 roots via the array */
76 static bool Equal(fInt, fInt);                            /* Returns true if two fInts are equal to each other */
77 static bool GreaterThan(fInt A, fInt B);                  /* Returns true if A > B */
78 
79 static fInt fExponential(fInt exponent);                  /* Can be used to calculate e^exponent */
80 static fInt fNaturalLog(fInt value);                      /* Can be used to calculate ln(value) */
81 
82 /* Fuse decoding functions
83  * -------------------------------------------------------------------------------------
84  */
85 static fInt fDecodeLinearFuse(uint32_t fuse_value, fInt f_min, fInt f_range, uint32_t bitlength);
86 static fInt fDecodeLogisticFuse(uint32_t fuse_value, fInt f_average, fInt f_range, uint32_t bitlength);
87 static fInt fDecodeLeakageID (uint32_t leakageID_fuse, fInt ln_max_div_min, fInt f_min, uint32_t bitlength);
88 
89 /* Internal Support Functions - Use these ONLY for testing or adding to internal functions
90  * -------------------------------------------------------------------------------------
91  * Some of the following functions take two INTs as their input - This is unsafe for a variety of reasons.
92  */
93 static fInt Divide (int, int);                            /* Divide two INTs and return result as FINT */
94 static fInt fNegate(fInt);
95 
96 static int uGetScaledDecimal (fInt);                      /* Internal function */
97 static int GetReal (fInt A);                              /* Internal function */
98 
99 /* -------------------------------------------------------------------------------------
100  * TROUBLESHOOTING INFORMATION
101  * -------------------------------------------------------------------------------------
102  * 1) ConvertToFraction - InputOutOfRangeException: Only accepts numbers smaller than MAX (default: 32767)
103  * 2) fAdd - OutputOutOfRangeException: Output bigger than MAX (default: 32767)
104  * 3) fMultiply - OutputOutOfRangeException:
105  * 4) fGetSquare - OutputOutOfRangeException:
106  * 5) fDivide - DivideByZeroException
107  * 6) fSqrt - NegativeSquareRootException: Input cannot be a negative number
108  */
109 
110 /* -------------------------------------------------------------------------------------
111  * START OF CODE
112  * -------------------------------------------------------------------------------------
113  */
fExponential(fInt exponent)114 static fInt fExponential(fInt exponent)        /*Can be used to calculate e^exponent*/
115 {
116 	uint32_t i;
117 	bool bNegated = false;
118 
119 	fInt fPositiveOne = ConvertToFraction(1);
120 	fInt fZERO = ConvertToFraction(0);
121 
122 	fInt lower_bound = Divide(78, 10000);
123 	fInt solution = fPositiveOne; /*Starting off with baseline of 1 */
124 	fInt error_term;
125 
126 	static const uint32_t k_array[11] = {55452, 27726, 13863, 6931, 4055, 2231, 1178, 606, 308, 155, 78};
127 	static const uint32_t expk_array[11] = {2560000, 160000, 40000, 20000, 15000, 12500, 11250, 10625, 10313, 10156, 10078};
128 
129 	if (GreaterThan(fZERO, exponent)) {
130 		exponent = fNegate(exponent);
131 		bNegated = true;
132 	}
133 
134 	while (GreaterThan(exponent, lower_bound)) {
135 		for (i = 0; i < 11; i++) {
136 			if (GreaterThan(exponent, GetScaledFraction(k_array[i], 10000))) {
137 				exponent = fSubtract(exponent, GetScaledFraction(k_array[i], 10000));
138 				solution = fMultiply(solution, GetScaledFraction(expk_array[i], 10000));
139 			}
140 		}
141 	}
142 
143 	error_term = fAdd(fPositiveOne, exponent);
144 
145 	solution = fMultiply(solution, error_term);
146 
147 	if (bNegated)
148 		solution = fDivide(fPositiveOne, solution);
149 
150 	return solution;
151 }
152 
fNaturalLog(fInt value)153 static fInt fNaturalLog(fInt value)
154 {
155 	uint32_t i;
156 	fInt upper_bound = Divide(8, 1000);
157 	fInt fNegativeOne = ConvertToFraction(-1);
158 	fInt solution = ConvertToFraction(0); /*Starting off with baseline of 0 */
159 	fInt error_term;
160 
161 	static const uint32_t k_array[10] = {160000, 40000, 20000, 15000, 12500, 11250, 10625, 10313, 10156, 10078};
162 	static const uint32_t logk_array[10] = {27726, 13863, 6931, 4055, 2231, 1178, 606, 308, 155, 78};
163 
164 	while (GreaterThan(fAdd(value, fNegativeOne), upper_bound)) {
165 		for (i = 0; i < 10; i++) {
166 			if (GreaterThan(value, GetScaledFraction(k_array[i], 10000))) {
167 				value = fDivide(value, GetScaledFraction(k_array[i], 10000));
168 				solution = fAdd(solution, GetScaledFraction(logk_array[i], 10000));
169 			}
170 		}
171 	}
172 
173 	error_term = fAdd(fNegativeOne, value);
174 
175 	return fAdd(solution, error_term);
176 }
177 
fDecodeLinearFuse(uint32_t fuse_value,fInt f_min,fInt f_range,uint32_t bitlength)178 static fInt fDecodeLinearFuse(uint32_t fuse_value, fInt f_min, fInt f_range, uint32_t bitlength)
179 {
180 	fInt f_fuse_value = Convert_ULONG_ToFraction(fuse_value);
181 	fInt f_bit_max_value = Convert_ULONG_ToFraction((uPow(2, bitlength)) - 1);
182 
183 	fInt f_decoded_value;
184 
185 	f_decoded_value = fDivide(f_fuse_value, f_bit_max_value);
186 	f_decoded_value = fMultiply(f_decoded_value, f_range);
187 	f_decoded_value = fAdd(f_decoded_value, f_min);
188 
189 	return f_decoded_value;
190 }
191 
192 
fDecodeLogisticFuse(uint32_t fuse_value,fInt f_average,fInt f_range,uint32_t bitlength)193 static fInt fDecodeLogisticFuse(uint32_t fuse_value, fInt f_average, fInt f_range, uint32_t bitlength)
194 {
195 	fInt f_fuse_value = Convert_ULONG_ToFraction(fuse_value);
196 	fInt f_bit_max_value = Convert_ULONG_ToFraction((uPow(2, bitlength)) - 1);
197 
198 	fInt f_CONSTANT_NEG13 = ConvertToFraction(-13);
199 	fInt f_CONSTANT1 = ConvertToFraction(1);
200 
201 	fInt f_decoded_value;
202 
203 	f_decoded_value = fSubtract(fDivide(f_bit_max_value, f_fuse_value), f_CONSTANT1);
204 	f_decoded_value = fNaturalLog(f_decoded_value);
205 	f_decoded_value = fMultiply(f_decoded_value, fDivide(f_range, f_CONSTANT_NEG13));
206 	f_decoded_value = fAdd(f_decoded_value, f_average);
207 
208 	return f_decoded_value;
209 }
210 
fDecodeLeakageID(uint32_t leakageID_fuse,fInt ln_max_div_min,fInt f_min,uint32_t bitlength)211 static fInt fDecodeLeakageID (uint32_t leakageID_fuse, fInt ln_max_div_min, fInt f_min, uint32_t bitlength)
212 {
213 	fInt fLeakage;
214 	fInt f_bit_max_value = Convert_ULONG_ToFraction((uPow(2, bitlength)) - 1);
215 
216 	fLeakage = fMultiply(ln_max_div_min, Convert_ULONG_ToFraction(leakageID_fuse));
217 	fLeakage = fDivide(fLeakage, f_bit_max_value);
218 	fLeakage = fExponential(fLeakage);
219 	fLeakage = fMultiply(fLeakage, f_min);
220 
221 	return fLeakage;
222 }
223 
ConvertToFraction(int X)224 static fInt ConvertToFraction(int X) /*Add all range checking here. Is it possible to make fInt a private declaration? */
225 {
226 	fInt temp;
227 
228 	if (X <= MAX)
229 		temp.full = (X << SHIFT_AMOUNT);
230 	else
231 		temp.full = 0;
232 
233 	return temp;
234 }
235 
fNegate(fInt X)236 static fInt fNegate(fInt X)
237 {
238 	fInt CONSTANT_NEGONE = ConvertToFraction(-1);
239 	return fMultiply(X, CONSTANT_NEGONE);
240 }
241 
Convert_ULONG_ToFraction(uint32_t X)242 static fInt Convert_ULONG_ToFraction(uint32_t X)
243 {
244 	fInt temp;
245 
246 	if (X <= MAX)
247 		temp.full = (X << SHIFT_AMOUNT);
248 	else
249 		temp.full = 0;
250 
251 	return temp;
252 }
253 
GetScaledFraction(int X,int factor)254 static fInt GetScaledFraction(int X, int factor)
255 {
256 	int times_shifted, factor_shifted;
257 	bool bNEGATED;
258 	fInt fValue;
259 
260 	times_shifted = 0;
261 	factor_shifted = 0;
262 	bNEGATED = false;
263 
264 	if (X < 0) {
265 		X = -1*X;
266 		bNEGATED = true;
267 	}
268 
269 	if (factor < 0) {
270 		factor = -1*factor;
271 		bNEGATED = !bNEGATED; /*If bNEGATED = true due to X < 0, this will cover the case of negative cancelling negative */
272 	}
273 
274 	if ((X > MAX) || factor > MAX) {
275 		if ((X/factor) <= MAX) {
276 			while (X > MAX) {
277 				X = X >> 1;
278 				times_shifted++;
279 			}
280 
281 			while (factor > MAX) {
282 				factor = factor >> 1;
283 				factor_shifted++;
284 			}
285 		} else {
286 			fValue.full = 0;
287 			return fValue;
288 		}
289 	}
290 
291 	if (factor == 1)
292 		return ConvertToFraction(X);
293 
294 	fValue = fDivide(ConvertToFraction(X * uPow(-1, bNEGATED)), ConvertToFraction(factor));
295 
296 	fValue.full = fValue.full << times_shifted;
297 	fValue.full = fValue.full >> factor_shifted;
298 
299 	return fValue;
300 }
301 
302 /* Addition using two fInts */
fAdd(fInt X,fInt Y)303 static fInt fAdd (fInt X, fInt Y)
304 {
305 	fInt Sum;
306 
307 	Sum.full = X.full + Y.full;
308 
309 	return Sum;
310 }
311 
312 /* Addition using two fInts */
fSubtract(fInt X,fInt Y)313 static fInt fSubtract (fInt X, fInt Y)
314 {
315 	fInt Difference;
316 
317 	Difference.full = X.full - Y.full;
318 
319 	return Difference;
320 }
321 
Equal(fInt A,fInt B)322 static bool Equal(fInt A, fInt B)
323 {
324 	if (A.full == B.full)
325 		return true;
326 	else
327 		return false;
328 }
329 
GreaterThan(fInt A,fInt B)330 static bool GreaterThan(fInt A, fInt B)
331 {
332 	if (A.full > B.full)
333 		return true;
334 	else
335 		return false;
336 }
337 
fMultiply(fInt X,fInt Y)338 static fInt fMultiply (fInt X, fInt Y) /* Uses 64-bit integers (int64_t) */
339 {
340 	fInt Product;
341 	int64_t tempProduct;
342 
343 	/*The following is for a very specific common case: Non-zero number with ONLY fractional portion*/
344 	/* TEMPORARILY DISABLED - CAN BE USED TO IMPROVE PRECISION
345 	bool X_LessThanOne, Y_LessThanOne;
346 
347 	X_LessThanOne = (X.partial.real == 0 && X.partial.decimal != 0 && X.full >= 0);
348 	Y_LessThanOne = (Y.partial.real == 0 && Y.partial.decimal != 0 && Y.full >= 0);
349 
350 	if (X_LessThanOne && Y_LessThanOne) {
351 		Product.full = X.full * Y.full;
352 		return Product
353 	}*/
354 
355 	tempProduct = ((int64_t)X.full) * ((int64_t)Y.full); /*Q(16,16)*Q(16,16) = Q(32, 32) - Might become a negative number! */
356 	tempProduct = tempProduct >> 16; /*Remove lagging 16 bits - Will lose some precision from decimal; */
357 	Product.full = (int)tempProduct; /*The int64_t will lose the leading 16 bits that were part of the integer portion */
358 
359 	return Product;
360 }
361 
fDivide(fInt X,fInt Y)362 static fInt fDivide (fInt X, fInt Y)
363 {
364 	fInt fZERO, fQuotient;
365 	int64_t longlongX, longlongY;
366 
367 	fZERO = ConvertToFraction(0);
368 
369 	if (Equal(Y, fZERO))
370 		return fZERO;
371 
372 	longlongX = (int64_t)X.full;
373 	longlongY = (int64_t)Y.full;
374 
375 	longlongX = longlongX << 16; /*Q(16,16) -> Q(32,32) */
376 
377 	div64_s64(longlongX, longlongY); /*Q(32,32) divided by Q(16,16) = Q(16,16) Back to original format */
378 
379 	fQuotient.full = (int)longlongX;
380 	return fQuotient;
381 }
382 
ConvertBackToInteger(fInt A)383 static int ConvertBackToInteger (fInt A) /*THIS is the function that will be used to check with the Golden settings table*/
384 {
385 	fInt fullNumber, scaledDecimal, scaledReal;
386 
387 	scaledReal.full = GetReal(A) * uPow(10, PRECISION-1); /* DOUBLE CHECK THISSSS!!! */
388 
389 	scaledDecimal.full = uGetScaledDecimal(A);
390 
391 	fullNumber = fAdd(scaledDecimal, scaledReal);
392 
393 	return fullNumber.full;
394 }
395 
fGetSquare(fInt A)396 static fInt fGetSquare(fInt A)
397 {
398 	return fMultiply(A, A);
399 }
400 
401 /* x_new = x_old - (x_old^2 - C) / (2 * x_old) */
fSqrt(fInt num)402 static fInt fSqrt(fInt num)
403 {
404 	fInt F_divide_Fprime, Fprime;
405 	fInt test;
406 	fInt twoShifted;
407 	int seed, counter, error;
408 	fInt x_new, x_old, C, y;
409 
410 	fInt fZERO = ConvertToFraction(0);
411 
412 	/* (0 > num) is the same as (num < 0), i.e., num is negative */
413 
414 	if (GreaterThan(fZERO, num) || Equal(fZERO, num))
415 		return fZERO;
416 
417 	C = num;
418 
419 	if (num.partial.real > 3000)
420 		seed = 60;
421 	else if (num.partial.real > 1000)
422 		seed = 30;
423 	else if (num.partial.real > 100)
424 		seed = 10;
425 	else
426 		seed = 2;
427 
428 	counter = 0;
429 
430 	if (Equal(num, fZERO)) /*Square Root of Zero is zero */
431 		return fZERO;
432 
433 	twoShifted = ConvertToFraction(2);
434 	x_new = ConvertToFraction(seed);
435 
436 	do {
437 		counter++;
438 
439 		x_old.full = x_new.full;
440 
441 		test = fGetSquare(x_old); /*1.75*1.75 is reverting back to 1 when shifted down */
442 		y = fSubtract(test, C); /*y = f(x) = x^2 - C; */
443 
444 		Fprime = fMultiply(twoShifted, x_old);
445 		F_divide_Fprime = fDivide(y, Fprime);
446 
447 		x_new = fSubtract(x_old, F_divide_Fprime);
448 
449 		error = ConvertBackToInteger(x_new) - ConvertBackToInteger(x_old);
450 
451 		if (counter > 20) /*20 is already way too many iterations. If we dont have an answer by then, we never will*/
452 			return x_new;
453 
454 	} while (uAbs(error) > 0);
455 
456 	return x_new;
457 }
458 
SolveQuadracticEqn(fInt A,fInt B,fInt C,fInt Roots[])459 static void SolveQuadracticEqn(fInt A, fInt B, fInt C, fInt Roots[])
460 {
461 	fInt *pRoots = &Roots[0];
462 	fInt temp, root_first, root_second;
463 	fInt f_CONSTANT10, f_CONSTANT100;
464 
465 	f_CONSTANT100 = ConvertToFraction(100);
466 	f_CONSTANT10 = ConvertToFraction(10);
467 
468 	while (GreaterThan(A, f_CONSTANT100) || GreaterThan(B, f_CONSTANT100) || GreaterThan(C, f_CONSTANT100)) {
469 		A = fDivide(A, f_CONSTANT10);
470 		B = fDivide(B, f_CONSTANT10);
471 		C = fDivide(C, f_CONSTANT10);
472 	}
473 
474 	temp = fMultiply(ConvertToFraction(4), A); /* root = 4*A */
475 	temp = fMultiply(temp, C); /* root = 4*A*C */
476 	temp = fSubtract(fGetSquare(B), temp); /* root = b^2 - 4AC */
477 	temp = fSqrt(temp); /*root = Sqrt (b^2 - 4AC); */
478 
479 	root_first = fSubtract(fNegate(B), temp); /* b - Sqrt(b^2 - 4AC) */
480 	root_second = fAdd(fNegate(B), temp); /* b + Sqrt(b^2 - 4AC) */
481 
482 	root_first = fDivide(root_first, ConvertToFraction(2)); /* [b +- Sqrt(b^2 - 4AC)]/[2] */
483 	root_first = fDivide(root_first, A); /*[b +- Sqrt(b^2 - 4AC)]/[2*A] */
484 
485 	root_second = fDivide(root_second, ConvertToFraction(2)); /* [b +- Sqrt(b^2 - 4AC)]/[2] */
486 	root_second = fDivide(root_second, A); /*[b +- Sqrt(b^2 - 4AC)]/[2*A] */
487 
488 	*(pRoots + 0) = root_first;
489 	*(pRoots + 1) = root_second;
490 }
491 
492 /* -----------------------------------------------------------------------------
493  * SUPPORT FUNCTIONS
494  * -----------------------------------------------------------------------------
495  */
496 
497 /* Conversion Functions */
GetReal(fInt A)498 static int GetReal (fInt A)
499 {
500 	return (A.full >> SHIFT_AMOUNT);
501 }
502 
Divide(int X,int Y)503 static fInt Divide (int X, int Y)
504 {
505 	fInt A, B, Quotient;
506 
507 	A.full = X << SHIFT_AMOUNT;
508 	B.full = Y << SHIFT_AMOUNT;
509 
510 	Quotient = fDivide(A, B);
511 
512 	return Quotient;
513 }
514 
uGetScaledDecimal(fInt A)515 static int uGetScaledDecimal (fInt A) /*Converts the fractional portion to whole integers - Costly function */
516 {
517 	int dec[PRECISION];
518 	int i, scaledDecimal = 0, tmp = A.partial.decimal;
519 
520 	for (i = 0; i < PRECISION; i++) {
521 		dec[i] = tmp / (1 << SHIFT_AMOUNT);
522 		tmp = tmp - ((1 << SHIFT_AMOUNT)*dec[i]);
523 		tmp *= 10;
524 		scaledDecimal = scaledDecimal + dec[i]*uPow(10, PRECISION - 1 - i);
525 	}
526 
527 	return scaledDecimal;
528 }
529 
uPow(int base,int power)530 static int uPow(int base, int power)
531 {
532 	if (power == 0)
533 		return 1;
534 	else
535 		return (base)*uPow(base, power - 1);
536 }
537 
uAbs(int X)538 static int uAbs(int X)
539 {
540 	if (X < 0)
541 		return (X * -1);
542 	else
543 		return X;
544 }
545 
fRoundUpByStepSize(fInt A,fInt fStepSize,bool error_term)546 static fInt fRoundUpByStepSize(fInt A, fInt fStepSize, bool error_term)
547 {
548 	fInt solution;
549 
550 	solution = fDivide(A, fStepSize);
551 	solution.partial.decimal = 0; /*All fractional digits changes to 0 */
552 
553 	if (error_term)
554 		solution.partial.real += 1; /*Error term of 1 added */
555 
556 	solution = fMultiply(solution, fStepSize);
557 	solution = fAdd(solution, fStepSize);
558 
559 	return solution;
560 }
561 
562