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