xref: /openbmc/linux/lib/math/div64.c (revision 3dc167ba)
12c64e9cbSAndy Shevchenko // SPDX-License-Identifier: GPL-2.0
22c64e9cbSAndy Shevchenko /*
32c64e9cbSAndy Shevchenko  * Copyright (C) 2003 Bernardo Innocenti <bernie@develer.com>
42c64e9cbSAndy Shevchenko  *
52c64e9cbSAndy Shevchenko  * Based on former do_div() implementation from asm-parisc/div64.h:
62c64e9cbSAndy Shevchenko  *	Copyright (C) 1999 Hewlett-Packard Co
72c64e9cbSAndy Shevchenko  *	Copyright (C) 1999 David Mosberger-Tang <davidm@hpl.hp.com>
82c64e9cbSAndy Shevchenko  *
92c64e9cbSAndy Shevchenko  *
102c64e9cbSAndy Shevchenko  * Generic C version of 64bit/32bit division and modulo, with
112c64e9cbSAndy Shevchenko  * 64bit result and 32bit remainder.
122c64e9cbSAndy Shevchenko  *
132c64e9cbSAndy Shevchenko  * The fast case for (n>>32 == 0) is handled inline by do_div().
142c64e9cbSAndy Shevchenko  *
152c64e9cbSAndy Shevchenko  * Code generated for this function might be very inefficient
162c64e9cbSAndy Shevchenko  * for some CPUs. __div64_32() can be overridden by linking arch-specific
172c64e9cbSAndy Shevchenko  * assembly versions such as arch/ppc/lib/div64.S and arch/sh/lib/div64.S
182c64e9cbSAndy Shevchenko  * or by defining a preprocessor macro in arch/include/asm/div64.h.
192c64e9cbSAndy Shevchenko  */
202c64e9cbSAndy Shevchenko 
212c64e9cbSAndy Shevchenko #include <linux/export.h>
222c64e9cbSAndy Shevchenko #include <linux/kernel.h>
232c64e9cbSAndy Shevchenko #include <linux/math64.h>
242c64e9cbSAndy Shevchenko 
252c64e9cbSAndy Shevchenko /* Not needed on 64bit architectures */
262c64e9cbSAndy Shevchenko #if BITS_PER_LONG == 32
272c64e9cbSAndy Shevchenko 
282c64e9cbSAndy Shevchenko #ifndef __div64_32
292c64e9cbSAndy Shevchenko uint32_t __attribute__((weak)) __div64_32(uint64_t *n, uint32_t base)
302c64e9cbSAndy Shevchenko {
312c64e9cbSAndy Shevchenko 	uint64_t rem = *n;
322c64e9cbSAndy Shevchenko 	uint64_t b = base;
332c64e9cbSAndy Shevchenko 	uint64_t res, d = 1;
342c64e9cbSAndy Shevchenko 	uint32_t high = rem >> 32;
352c64e9cbSAndy Shevchenko 
362c64e9cbSAndy Shevchenko 	/* Reduce the thing a bit first */
372c64e9cbSAndy Shevchenko 	res = 0;
382c64e9cbSAndy Shevchenko 	if (high >= base) {
392c64e9cbSAndy Shevchenko 		high /= base;
402c64e9cbSAndy Shevchenko 		res = (uint64_t) high << 32;
412c64e9cbSAndy Shevchenko 		rem -= (uint64_t) (high*base) << 32;
422c64e9cbSAndy Shevchenko 	}
432c64e9cbSAndy Shevchenko 
442c64e9cbSAndy Shevchenko 	while ((int64_t)b > 0 && b < rem) {
452c64e9cbSAndy Shevchenko 		b = b+b;
462c64e9cbSAndy Shevchenko 		d = d+d;
472c64e9cbSAndy Shevchenko 	}
482c64e9cbSAndy Shevchenko 
492c64e9cbSAndy Shevchenko 	do {
502c64e9cbSAndy Shevchenko 		if (rem >= b) {
512c64e9cbSAndy Shevchenko 			rem -= b;
522c64e9cbSAndy Shevchenko 			res += d;
532c64e9cbSAndy Shevchenko 		}
542c64e9cbSAndy Shevchenko 		b >>= 1;
552c64e9cbSAndy Shevchenko 		d >>= 1;
562c64e9cbSAndy Shevchenko 	} while (d);
572c64e9cbSAndy Shevchenko 
582c64e9cbSAndy Shevchenko 	*n = res;
592c64e9cbSAndy Shevchenko 	return rem;
602c64e9cbSAndy Shevchenko }
612c64e9cbSAndy Shevchenko EXPORT_SYMBOL(__div64_32);
622c64e9cbSAndy Shevchenko #endif
632c64e9cbSAndy Shevchenko 
642c64e9cbSAndy Shevchenko /**
652c64e9cbSAndy Shevchenko  * div_s64_rem - signed 64bit divide with 64bit divisor and remainder
662c64e9cbSAndy Shevchenko  * @dividend:	64bit dividend
672c64e9cbSAndy Shevchenko  * @divisor:	64bit divisor
682c64e9cbSAndy Shevchenko  * @remainder:  64bit remainder
692c64e9cbSAndy Shevchenko  */
702c64e9cbSAndy Shevchenko #ifndef div_s64_rem
712c64e9cbSAndy Shevchenko s64 div_s64_rem(s64 dividend, s32 divisor, s32 *remainder)
722c64e9cbSAndy Shevchenko {
732c64e9cbSAndy Shevchenko 	u64 quotient;
742c64e9cbSAndy Shevchenko 
752c64e9cbSAndy Shevchenko 	if (dividend < 0) {
762c64e9cbSAndy Shevchenko 		quotient = div_u64_rem(-dividend, abs(divisor), (u32 *)remainder);
772c64e9cbSAndy Shevchenko 		*remainder = -*remainder;
782c64e9cbSAndy Shevchenko 		if (divisor > 0)
792c64e9cbSAndy Shevchenko 			quotient = -quotient;
802c64e9cbSAndy Shevchenko 	} else {
812c64e9cbSAndy Shevchenko 		quotient = div_u64_rem(dividend, abs(divisor), (u32 *)remainder);
822c64e9cbSAndy Shevchenko 		if (divisor < 0)
832c64e9cbSAndy Shevchenko 			quotient = -quotient;
842c64e9cbSAndy Shevchenko 	}
852c64e9cbSAndy Shevchenko 	return quotient;
862c64e9cbSAndy Shevchenko }
872c64e9cbSAndy Shevchenko EXPORT_SYMBOL(div_s64_rem);
882c64e9cbSAndy Shevchenko #endif
892c64e9cbSAndy Shevchenko 
902c64e9cbSAndy Shevchenko /**
912c64e9cbSAndy Shevchenko  * div64_u64_rem - unsigned 64bit divide with 64bit divisor and remainder
922c64e9cbSAndy Shevchenko  * @dividend:	64bit dividend
932c64e9cbSAndy Shevchenko  * @divisor:	64bit divisor
942c64e9cbSAndy Shevchenko  * @remainder:  64bit remainder
952c64e9cbSAndy Shevchenko  *
962c64e9cbSAndy Shevchenko  * This implementation is a comparable to algorithm used by div64_u64.
972c64e9cbSAndy Shevchenko  * But this operation, which includes math for calculating the remainder,
982c64e9cbSAndy Shevchenko  * is kept distinct to avoid slowing down the div64_u64 operation on 32bit
992c64e9cbSAndy Shevchenko  * systems.
1002c64e9cbSAndy Shevchenko  */
1012c64e9cbSAndy Shevchenko #ifndef div64_u64_rem
1022c64e9cbSAndy Shevchenko u64 div64_u64_rem(u64 dividend, u64 divisor, u64 *remainder)
1032c64e9cbSAndy Shevchenko {
1042c64e9cbSAndy Shevchenko 	u32 high = divisor >> 32;
1052c64e9cbSAndy Shevchenko 	u64 quot;
1062c64e9cbSAndy Shevchenko 
1072c64e9cbSAndy Shevchenko 	if (high == 0) {
1082c64e9cbSAndy Shevchenko 		u32 rem32;
1092c64e9cbSAndy Shevchenko 		quot = div_u64_rem(dividend, divisor, &rem32);
1102c64e9cbSAndy Shevchenko 		*remainder = rem32;
1112c64e9cbSAndy Shevchenko 	} else {
1122c64e9cbSAndy Shevchenko 		int n = fls(high);
1132c64e9cbSAndy Shevchenko 		quot = div_u64(dividend >> n, divisor >> n);
1142c64e9cbSAndy Shevchenko 
1152c64e9cbSAndy Shevchenko 		if (quot != 0)
1162c64e9cbSAndy Shevchenko 			quot--;
1172c64e9cbSAndy Shevchenko 
1182c64e9cbSAndy Shevchenko 		*remainder = dividend - quot * divisor;
1192c64e9cbSAndy Shevchenko 		if (*remainder >= divisor) {
1202c64e9cbSAndy Shevchenko 			quot++;
1212c64e9cbSAndy Shevchenko 			*remainder -= divisor;
1222c64e9cbSAndy Shevchenko 		}
1232c64e9cbSAndy Shevchenko 	}
1242c64e9cbSAndy Shevchenko 
1252c64e9cbSAndy Shevchenko 	return quot;
1262c64e9cbSAndy Shevchenko }
1272c64e9cbSAndy Shevchenko EXPORT_SYMBOL(div64_u64_rem);
1282c64e9cbSAndy Shevchenko #endif
1292c64e9cbSAndy Shevchenko 
1302c64e9cbSAndy Shevchenko /**
1312c64e9cbSAndy Shevchenko  * div64_u64 - unsigned 64bit divide with 64bit divisor
1322c64e9cbSAndy Shevchenko  * @dividend:	64bit dividend
1332c64e9cbSAndy Shevchenko  * @divisor:	64bit divisor
1342c64e9cbSAndy Shevchenko  *
1352c64e9cbSAndy Shevchenko  * This implementation is a modified version of the algorithm proposed
1362c64e9cbSAndy Shevchenko  * by the book 'Hacker's Delight'.  The original source and full proof
1372c64e9cbSAndy Shevchenko  * can be found here and is available for use without restriction.
1382c64e9cbSAndy Shevchenko  *
1392c64e9cbSAndy Shevchenko  * 'http://www.hackersdelight.org/hdcodetxt/divDouble.c.txt'
1402c64e9cbSAndy Shevchenko  */
1412c64e9cbSAndy Shevchenko #ifndef div64_u64
1422c64e9cbSAndy Shevchenko u64 div64_u64(u64 dividend, u64 divisor)
1432c64e9cbSAndy Shevchenko {
1442c64e9cbSAndy Shevchenko 	u32 high = divisor >> 32;
1452c64e9cbSAndy Shevchenko 	u64 quot;
1462c64e9cbSAndy Shevchenko 
1472c64e9cbSAndy Shevchenko 	if (high == 0) {
1482c64e9cbSAndy Shevchenko 		quot = div_u64(dividend, divisor);
1492c64e9cbSAndy Shevchenko 	} else {
1502c64e9cbSAndy Shevchenko 		int n = fls(high);
1512c64e9cbSAndy Shevchenko 		quot = div_u64(dividend >> n, divisor >> n);
1522c64e9cbSAndy Shevchenko 
1532c64e9cbSAndy Shevchenko 		if (quot != 0)
1542c64e9cbSAndy Shevchenko 			quot--;
1552c64e9cbSAndy Shevchenko 		if ((dividend - quot * divisor) >= divisor)
1562c64e9cbSAndy Shevchenko 			quot++;
1572c64e9cbSAndy Shevchenko 	}
1582c64e9cbSAndy Shevchenko 
1592c64e9cbSAndy Shevchenko 	return quot;
1602c64e9cbSAndy Shevchenko }
1612c64e9cbSAndy Shevchenko EXPORT_SYMBOL(div64_u64);
1622c64e9cbSAndy Shevchenko #endif
1632c64e9cbSAndy Shevchenko 
1642c64e9cbSAndy Shevchenko /**
1652c64e9cbSAndy Shevchenko  * div64_s64 - signed 64bit divide with 64bit divisor
1662c64e9cbSAndy Shevchenko  * @dividend:	64bit dividend
1672c64e9cbSAndy Shevchenko  * @divisor:	64bit divisor
1682c64e9cbSAndy Shevchenko  */
1692c64e9cbSAndy Shevchenko #ifndef div64_s64
1702c64e9cbSAndy Shevchenko s64 div64_s64(s64 dividend, s64 divisor)
1712c64e9cbSAndy Shevchenko {
1722c64e9cbSAndy Shevchenko 	s64 quot, t;
1732c64e9cbSAndy Shevchenko 
1742c64e9cbSAndy Shevchenko 	quot = div64_u64(abs(dividend), abs(divisor));
1752c64e9cbSAndy Shevchenko 	t = (dividend ^ divisor) >> 63;
1762c64e9cbSAndy Shevchenko 
1772c64e9cbSAndy Shevchenko 	return (quot ^ t) - t;
1782c64e9cbSAndy Shevchenko }
1792c64e9cbSAndy Shevchenko EXPORT_SYMBOL(div64_s64);
1802c64e9cbSAndy Shevchenko #endif
1812c64e9cbSAndy Shevchenko 
1822c64e9cbSAndy Shevchenko #endif /* BITS_PER_LONG == 32 */
1832c64e9cbSAndy Shevchenko 
1842c64e9cbSAndy Shevchenko /*
1852c64e9cbSAndy Shevchenko  * Iterative div/mod for use when dividend is not expected to be much
1862c64e9cbSAndy Shevchenko  * bigger than divisor.
1872c64e9cbSAndy Shevchenko  */
1882c64e9cbSAndy Shevchenko u32 iter_div_u64_rem(u64 dividend, u32 divisor, u64 *remainder)
1892c64e9cbSAndy Shevchenko {
1902c64e9cbSAndy Shevchenko 	return __iter_div_u64_rem(dividend, divisor, remainder);
1912c64e9cbSAndy Shevchenko }
1922c64e9cbSAndy Shevchenko EXPORT_SYMBOL(iter_div_u64_rem);
1933dc167baSOleg Nesterov 
1943dc167baSOleg Nesterov #ifndef mul_u64_u64_div_u64
1953dc167baSOleg Nesterov u64 mul_u64_u64_div_u64(u64 a, u64 b, u64 c)
1963dc167baSOleg Nesterov {
1973dc167baSOleg Nesterov 	u64 res = 0, div, rem;
1983dc167baSOleg Nesterov 	int shift;
1993dc167baSOleg Nesterov 
2003dc167baSOleg Nesterov 	/* can a * b overflow ? */
2013dc167baSOleg Nesterov 	if (ilog2(a) + ilog2(b) > 62) {
2023dc167baSOleg Nesterov 		/*
2033dc167baSOleg Nesterov 		 * (b * a) / c is equal to
2043dc167baSOleg Nesterov 		 *
2053dc167baSOleg Nesterov 		 *      (b / c) * a +
2063dc167baSOleg Nesterov 		 *      (b % c) * a / c
2073dc167baSOleg Nesterov 		 *
2083dc167baSOleg Nesterov 		 * if nothing overflows. Can the 1st multiplication
2093dc167baSOleg Nesterov 		 * overflow? Yes, but we do not care: this can only
2103dc167baSOleg Nesterov 		 * happen if the end result can't fit in u64 anyway.
2113dc167baSOleg Nesterov 		 *
2123dc167baSOleg Nesterov 		 * So the code below does
2133dc167baSOleg Nesterov 		 *
2143dc167baSOleg Nesterov 		 *      res = (b / c) * a;
2153dc167baSOleg Nesterov 		 *      b = b % c;
2163dc167baSOleg Nesterov 		 */
2173dc167baSOleg Nesterov 		div = div64_u64_rem(b, c, &rem);
2183dc167baSOleg Nesterov 		res = div * a;
2193dc167baSOleg Nesterov 		b = rem;
2203dc167baSOleg Nesterov 
2213dc167baSOleg Nesterov 		shift = ilog2(a) + ilog2(b) - 62;
2223dc167baSOleg Nesterov 		if (shift > 0) {
2233dc167baSOleg Nesterov 			/* drop precision */
2243dc167baSOleg Nesterov 			b >>= shift;
2253dc167baSOleg Nesterov 			c >>= shift;
2263dc167baSOleg Nesterov 			if (!c)
2273dc167baSOleg Nesterov 				return res;
2283dc167baSOleg Nesterov 		}
2293dc167baSOleg Nesterov 	}
2303dc167baSOleg Nesterov 
2313dc167baSOleg Nesterov 	return res + div64_u64(a * b, c);
2323dc167baSOleg Nesterov }
2333dc167baSOleg Nesterov #endif
234