2024-03-03 09:25:50

by Uwe Kleine-König

[permalink] [raw]
Subject: [PATCH v2] mul_u64_u64_div_u64: Increase precision by conditionally swapping a and b

As indicated in the added comment, the algorithm works better if b is
big. As multiplication is commutative, a and b can be swapped. Do this
if a is bigger than b.

Signed-off-by: Uwe Kleine-König <[email protected]>
---
Changes since v1:
- Make use of swap() (Thanks Marc)
- Fix a typo in a code comment (Thanks Randy)
- Fix a typo in the commit log (s/If/if/; noticed myself)

v1 got a Tested-by from Biju; I didn't add it here as the patch changed.
I'm optimistic that this v2 would pass his tests, too, but I don't wanna
assume stuff when adding tags.

Best regards
Uwe

lib/math/div64.c | 15 +++++++++++++++
1 file changed, 15 insertions(+)

diff --git a/lib/math/div64.c b/lib/math/div64.c
index 55a81782e271..191761b1b623 100644
--- a/lib/math/div64.c
+++ b/lib/math/div64.c
@@ -22,6 +22,7 @@
#include <linux/export.h>
#include <linux/math.h>
#include <linux/math64.h>
+#include <linux/minmax.h>
#include <linux/log2.h>

/* Not needed on 64bit architectures */
@@ -190,6 +191,20 @@ u64 mul_u64_u64_div_u64(u64 a, u64 b, u64 c)

/* can a * b overflow ? */
if (ilog2(a) + ilog2(b) > 62) {
+ /*
+ * Note that the algorithm after the if block below might lose
+ * some precision and the result is more exact for b > a. So
+ * exchange a and b if a is bigger than b.
+ *
+ * For example with a = 43980465100800, b = 100000000, c = 1000000000
+ * the below calculation doesn't modify b at all because div == 0
+ * and then shift becomes 45 + 26 - 62 = 9 and so the result
+ * becomes 4398035251080. However with a and b swapped the exact
+ * result is calculated (i.e. 4398046510080).
+ */
+ if (a > b)
+ swap(a, b);
+
/*
* (b * a) / c is equal to
*

base-commit: 1870cdc0e8dee32e3c221704a2977898ba4c10e8
--
2.43.0



2024-03-03 09:55:39

by Biju Das

[permalink] [raw]
Subject: RE: [PATCH v2] mul_u64_u64_div_u64: Increase precision by conditionally swapping a and b

Hi Uwe Kleine-König,

> -----Original Message-----
> From: Uwe Kleine-König <[email protected]>
> Sent: Sunday, March 3, 2024 9:24 AM
> Subject: [PATCH v2] mul_u64_u64_div_u64: Increase precision by conditionally swapping a and b
>
> As indicated in the added comment, the algorithm works better if b is big. As multiplication is
> commutative, a and b can be swapped. Do this if a is bigger than b.
>
> Signed-off-by: Uwe Kleine-König <[email protected]>
> ---
> Changes since v1:
> - Make use of swap() (Thanks Marc)
> - Fix a typo in a code comment (Thanks Randy)
> - Fix a typo in the commit log (s/If/if/; noticed myself)
>
> v1 got a Tested-by from Biju; I didn't add it here as the patch changed.
> I'm optimistic that this v2 would pass his tests, too, but I don't wanna assume stuff when adding tags.

Tested-by: Biju Das <[email protected]>

Done the testing by dropping rzg2l_gpt_mul_u64_u64_div_u64() in patch[1]

[1]
https://patchwork.kernel.org/project/linux-renesas-soc/patch/[email protected]/

Cheers,
Biju
>
> lib/math/div64.c | 15 +++++++++++++++
> 1 file changed, 15 insertions(+)
>
> diff --git a/lib/math/div64.c b/lib/math/div64.c index 55a81782e271..191761b1b623 100644
> --- a/lib/math/div64.c
> +++ b/lib/math/div64.c
> @@ -22,6 +22,7 @@
> #include <linux/export.h>
> #include <linux/math.h>
> #include <linux/math64.h>
> +#include <linux/minmax.h>
> #include <linux/log2.h>
>
> /* Not needed on 64bit architectures */ @@ -190,6 +191,20 @@ u64 mul_u64_u64_div_u64(u64 a, u64 b, u64
> c)
>
> /* can a * b overflow ? */
> if (ilog2(a) + ilog2(b) > 62) {
> + /*
> + * Note that the algorithm after the if block below might lose
> + * some precision and the result is more exact for b > a. So
> + * exchange a and b if a is bigger than b.
> + *
> + * For example with a = 43980465100800, b = 100000000, c = 1000000000
> + * the below calculation doesn't modify b at all because div == 0
> + * and then shift becomes 45 + 26 - 62 = 9 and so the result
> + * becomes 4398035251080. However with a and b swapped the exact
> + * result is calculated (i.e. 4398046510080).
> + */
> + if (a > b)
> + swap(a, b);
> +
> /*
> * (b * a) / c is equal to
> *
>
> base-commit: 1870cdc0e8dee32e3c221704a2977898ba4c10e8
> --
> 2.43.0