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