From 5c09099d73ae24392ee06627f562080364085ed7 Mon Sep 17 00:00:00 2001 From: Kent Overstreet Date: Tue, 11 Jul 2023 20:35:06 -0400 Subject: mean and variance: Promote to lib/math Small statistics library, for taking in a series of value and computing mean, weighted mean, standard deviation and weighted deviation. The main use case is for statistics on latency measurements. Signed-off-by: Kent Overstreet Cc: Daniel Hill Cc: Darrick J. Wong --- MAINTAINERS | 9 ++ fs/bcachefs/Kconfig | 10 +- fs/bcachefs/Makefile | 3 - fs/bcachefs/mean_and_variance.c | 165 ------------------------ fs/bcachefs/mean_and_variance.h | 201 ----------------------------- fs/bcachefs/mean_and_variance_test.c | 240 ----------------------------------- fs/bcachefs/util.c | 2 +- fs/bcachefs/util.h | 3 +- include/linux/mean_and_variance.h | 201 +++++++++++++++++++++++++++++ lib/Kconfig.debug | 9 ++ lib/math/Kconfig | 3 + lib/math/Makefile | 2 + lib/math/mean_and_variance.c | 164 ++++++++++++++++++++++++ lib/math/mean_and_variance_test.c | 239 ++++++++++++++++++++++++++++++++++ 14 files changed, 630 insertions(+), 621 deletions(-) delete mode 100644 fs/bcachefs/mean_and_variance.c delete mode 100644 fs/bcachefs/mean_and_variance.h delete mode 100644 fs/bcachefs/mean_and_variance_test.c create mode 100644 include/linux/mean_and_variance.h create mode 100644 lib/math/mean_and_variance.c create mode 100644 lib/math/mean_and_variance_test.c diff --git a/MAINTAINERS b/MAINTAINERS index 8d1052fa6a69..de635cfd354d 100644 --- a/MAINTAINERS +++ b/MAINTAINERS @@ -13379,6 +13379,15 @@ S: Maintained F: drivers/net/mdio/mdio-regmap.c F: include/linux/mdio/mdio-regmap.h +MEAN AND VARIANCE LIBRARY +M: Daniel B. Hill +M: Kent Overstreet +S: Maintained +T: git https://github.com/YellowOnion/linux/ +F: include/linux/mean_and_variance.h +F: lib/math/mean_and_variance.c +F: lib/math/mean_and_variance_test.c + MEASUREMENT COMPUTING CIO-DAC IIO DRIVER M: William Breathitt Gray L: linux-iio@vger.kernel.org diff --git a/fs/bcachefs/Kconfig b/fs/bcachefs/Kconfig index 5cdfef3b551a..72d1179262b3 100644 --- a/fs/bcachefs/Kconfig +++ b/fs/bcachefs/Kconfig @@ -24,6 +24,7 @@ config BCACHEFS_FS select XXHASH select SRCU select SYMBOLIC_ERRNAME + select MEAN_AND_VARIANCE help The bcachefs filesystem - a modern, copy on write filesystem, with support for multiple devices, compression, checksumming, etc. @@ -86,12 +87,3 @@ config BCACHEFS_SIX_OPTIMISTIC_SPIN Instead of immediately sleeping when attempting to take a six lock that is held by another thread, spin for a short while, as long as the thread owning the lock is running. - -config MEAN_AND_VARIANCE_UNIT_TEST - tristate "mean_and_variance unit tests" if !KUNIT_ALL_TESTS - depends on KUNIT - depends on BCACHEFS_FS - default KUNIT_ALL_TESTS - help - This option enables the kunit tests for mean_and_variance module. - If unsure, say N. diff --git a/fs/bcachefs/Makefile b/fs/bcachefs/Makefile index 1a05cecda7cc..b11ba74b8ad4 100644 --- a/fs/bcachefs/Makefile +++ b/fs/bcachefs/Makefile @@ -57,7 +57,6 @@ bcachefs-y := \ keylist.o \ logged_ops.o \ lru.o \ - mean_and_variance.o \ migrate.o \ move.o \ movinggc.o \ @@ -88,5 +87,3 @@ bcachefs-y := \ util.o \ varint.o \ xattr.o - -obj-$(CONFIG_MEAN_AND_VARIANCE_UNIT_TEST) += mean_and_variance_test.o diff --git a/fs/bcachefs/mean_and_variance.c b/fs/bcachefs/mean_and_variance.c deleted file mode 100644 index bf0ef668fd38..000000000000 --- a/fs/bcachefs/mean_and_variance.c +++ /dev/null @@ -1,165 +0,0 @@ -// SPDX-License-Identifier: GPL-2.0 -/* - * Functions for incremental mean and variance. - * - * This program is free software; you can redistribute it and/or modify it - * under the terms of the GNU General Public License version 2 as published by - * the Free Software Foundation. - * - * This program is distributed in the hope that it will be useful, but WITHOUT - * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or - * FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for - * more details. - * - * Copyright © 2022 Daniel B. Hill - * - * Author: Daniel B. Hill - * - * Description: - * - * This is includes some incremental algorithms for mean and variance calculation - * - * Derived from the paper: https://fanf2.user.srcf.net/hermes/doc/antiforgery/stats.pdf - * - * Create a struct and if it's the weighted variant set the w field (weight = 2^k). - * - * Use mean_and_variance[_weighted]_update() on the struct to update it's state. - * - * Use the mean_and_variance[_weighted]_get_* functions to calculate the mean and variance, some computation - * is deferred to these functions for performance reasons. - * - * see lib/math/mean_and_variance_test.c for examples of usage. - * - * DO NOT access the mean and variance fields of the weighted variants directly. - * DO NOT change the weight after calling update. - */ - -#include -#include -#include -#include -#include -#include -#include - -#include "mean_and_variance.h" - -u128_u u128_div(u128_u n, u64 d) -{ - u128_u r; - u64 rem; - u64 hi = u128_hi(n); - u64 lo = u128_lo(n); - u64 h = hi & ((u64) U32_MAX << 32); - u64 l = (hi & (u64) U32_MAX) << 32; - - r = u128_shl(u64_to_u128(div64_u64_rem(h, d, &rem)), 64); - r = u128_add(r, u128_shl(u64_to_u128(div64_u64_rem(l + (rem << 32), d, &rem)), 32)); - r = u128_add(r, u64_to_u128(div64_u64_rem(lo + (rem << 32), d, &rem))); - return r; -} -EXPORT_SYMBOL_GPL(u128_div); - -/** - * mean_and_variance_get_mean() - get mean from @s - * @s: mean and variance number of samples and their sums - */ -s64 mean_and_variance_get_mean(struct mean_and_variance s) -{ - return s.n ? div64_u64(s.sum, s.n) : 0; -} -EXPORT_SYMBOL_GPL(mean_and_variance_get_mean); - -/** - * mean_and_variance_get_variance() - get variance from @s1 - * @s1: mean and variance number of samples and sums - * - * see linked pdf equation 12. - */ -u64 mean_and_variance_get_variance(struct mean_and_variance s1) -{ - if (s1.n) { - u128_u s2 = u128_div(s1.sum_squares, s1.n); - u64 s3 = abs(mean_and_variance_get_mean(s1)); - - return u128_lo(u128_sub(s2, u128_square(s3))); - } else { - return 0; - } -} -EXPORT_SYMBOL_GPL(mean_and_variance_get_variance); - -/** - * mean_and_variance_get_stddev() - get standard deviation from @s - * @s: mean and variance number of samples and their sums - */ -u32 mean_and_variance_get_stddev(struct mean_and_variance s) -{ - return int_sqrt64(mean_and_variance_get_variance(s)); -} -EXPORT_SYMBOL_GPL(mean_and_variance_get_stddev); - -/** - * mean_and_variance_weighted_update() - exponentially weighted variant of mean_and_variance_update() - * @s: mean and variance number of samples and their sums - * @x: new value to include in the &mean_and_variance_weighted - * - * see linked pdf: function derived from equations 140-143 where alpha = 2^w. - * values are stored bitshifted for performance and added precision. - */ -void mean_and_variance_weighted_update(struct mean_and_variance_weighted *s, s64 x) -{ - // previous weighted variance. - u8 w = s->weight; - u64 var_w0 = s->variance; - // new value weighted. - s64 x_w = x << w; - s64 diff_w = x_w - s->mean; - s64 diff = fast_divpow2(diff_w, w); - // new mean weighted. - s64 u_w1 = s->mean + diff; - - if (!s->init) { - s->mean = x_w; - s->variance = 0; - } else { - s->mean = u_w1; - s->variance = ((var_w0 << w) - var_w0 + ((diff_w * (x_w - u_w1)) >> w)) >> w; - } - s->init = true; -} -EXPORT_SYMBOL_GPL(mean_and_variance_weighted_update); - -/** - * mean_and_variance_weighted_get_mean() - get mean from @s - * @s: mean and variance number of samples and their sums - */ -s64 mean_and_variance_weighted_get_mean(struct mean_and_variance_weighted s) -{ - return fast_divpow2(s.mean, s.weight); -} -EXPORT_SYMBOL_GPL(mean_and_variance_weighted_get_mean); - -/** - * mean_and_variance_weighted_get_variance() -- get variance from @s - * @s: mean and variance number of samples and their sums - */ -u64 mean_and_variance_weighted_get_variance(struct mean_and_variance_weighted s) -{ - // always positive don't need fast divpow2 - return s.variance >> s.weight; -} -EXPORT_SYMBOL_GPL(mean_and_variance_weighted_get_variance); - -/** - * mean_and_variance_weighted_get_stddev() - get standard deviation from @s - * @s: mean and variance number of samples and their sums - */ -u32 mean_and_variance_weighted_get_stddev(struct mean_and_variance_weighted s) -{ - return int_sqrt64(mean_and_variance_weighted_get_variance(s)); -} -EXPORT_SYMBOL_GPL(mean_and_variance_weighted_get_stddev); - -MODULE_AUTHOR("Daniel B. Hill"); -MODULE_LICENSE("GPL"); diff --git a/fs/bcachefs/mean_and_variance.h b/fs/bcachefs/mean_and_variance.h deleted file mode 100644 index 64df11ab422b..000000000000 --- a/fs/bcachefs/mean_and_variance.h +++ /dev/null @@ -1,201 +0,0 @@ -/* SPDX-License-Identifier: GPL-2.0 */ -#ifndef MEAN_AND_VARIANCE_H_ -#define MEAN_AND_VARIANCE_H_ - -#include -#include -#include -#include - -#define SQRT_U64_MAX 4294967295ULL - -/* - * u128_u: u128 user mode, because not all architectures support a real int128 - * type - * - * We don't use this version in userspace, because in userspace we link with - * Rust and rustc has issues with u128. - */ - -#if defined(__SIZEOF_INT128__) && defined(__KERNEL__) && !defined(CONFIG_PARISC) - -typedef struct { - unsigned __int128 v; -} __aligned(16) u128_u; - -static inline u128_u u64_to_u128(u64 a) -{ - return (u128_u) { .v = a }; -} - -static inline u64 u128_lo(u128_u a) -{ - return a.v; -} - -static inline u64 u128_hi(u128_u a) -{ - return a.v >> 64; -} - -static inline u128_u u128_add(u128_u a, u128_u b) -{ - a.v += b.v; - return a; -} - -static inline u128_u u128_sub(u128_u a, u128_u b) -{ - a.v -= b.v; - return a; -} - -static inline u128_u u128_shl(u128_u a, s8 shift) -{ - a.v <<= shift; - return a; -} - -static inline u128_u u128_square(u64 a) -{ - u128_u b = u64_to_u128(a); - - b.v *= b.v; - return b; -} - -#else - -typedef struct { - u64 hi, lo; -} __aligned(16) u128_u; - -/* conversions */ - -static inline u128_u u64_to_u128(u64 a) -{ - return (u128_u) { .lo = a }; -} - -static inline u64 u128_lo(u128_u a) -{ - return a.lo; -} - -static inline u64 u128_hi(u128_u a) -{ - return a.hi; -} - -/* arithmetic */ - -static inline u128_u u128_add(u128_u a, u128_u b) -{ - u128_u c; - - c.lo = a.lo + b.lo; - c.hi = a.hi + b.hi + (c.lo < a.lo); - return c; -} - -static inline u128_u u128_sub(u128_u a, u128_u b) -{ - u128_u c; - - c.lo = a.lo - b.lo; - c.hi = a.hi - b.hi - (c.lo > a.lo); - return c; -} - -static inline u128_u u128_shl(u128_u i, s8 shift) -{ - u128_u r; - - r.lo = i.lo << shift; - if (shift < 64) - r.hi = (i.hi << shift) | (i.lo >> (64 - shift)); - else { - r.hi = i.lo << (shift - 64); - r.lo = 0; - } - return r; -} - -static inline u128_u u128_square(u64 i) -{ - u128_u r; - u64 h = i >> 32, l = i & U32_MAX; - - r = u128_shl(u64_to_u128(h*h), 64); - r = u128_add(r, u128_shl(u64_to_u128(h*l), 32)); - r = u128_add(r, u128_shl(u64_to_u128(l*h), 32)); - r = u128_add(r, u64_to_u128(l*l)); - return r; -} - -#endif - -static inline u128_u u64s_to_u128(u64 hi, u64 lo) -{ - u128_u c = u64_to_u128(hi); - - c = u128_shl(c, 64); - c = u128_add(c, u64_to_u128(lo)); - return c; -} - -u128_u u128_div(u128_u n, u64 d); - -struct mean_and_variance { - s64 n; - s64 sum; - u128_u sum_squares; -}; - -/* expontentially weighted variant */ -struct mean_and_variance_weighted { - bool init; - u8 weight; /* base 2 logarithim */ - s64 mean; - u64 variance; -}; - -/** - * fast_divpow2() - fast approximation for n / (1 << d) - * @n: numerator - * @d: the power of 2 denominator. - * - * note: this rounds towards 0. - */ -static inline s64 fast_divpow2(s64 n, u8 d) -{ - return (n + ((n < 0) ? ((1 << d) - 1) : 0)) >> d; -} - -/** - * mean_and_variance_update() - update a mean_and_variance struct @s1 with a new sample @v1 - * and return it. - * @s1: the mean_and_variance to update. - * @v1: the new sample. - * - * see linked pdf equation 12. - */ -static inline void -mean_and_variance_update(struct mean_and_variance *s, s64 v) -{ - s->n++; - s->sum += v; - s->sum_squares = u128_add(s->sum_squares, u128_square(abs(v))); -} - -s64 mean_and_variance_get_mean(struct mean_and_variance s); -u64 mean_and_variance_get_variance(struct mean_and_variance s1); -u32 mean_and_variance_get_stddev(struct mean_and_variance s); - -void mean_and_variance_weighted_update(struct mean_and_variance_weighted *s, s64 v); - -s64 mean_and_variance_weighted_get_mean(struct mean_and_variance_weighted s); -u64 mean_and_variance_weighted_get_variance(struct mean_and_variance_weighted s); -u32 mean_and_variance_weighted_get_stddev(struct mean_and_variance_weighted s); - -#endif // MEAN_AND_VAIRANCE_H_ diff --git a/fs/bcachefs/mean_and_variance_test.c b/fs/bcachefs/mean_and_variance_test.c deleted file mode 100644 index 019583c3ca0e..000000000000 --- a/fs/bcachefs/mean_and_variance_test.c +++ /dev/null @@ -1,240 +0,0 @@ -// SPDX-License-Identifier: GPL-2.0 -#include - -#include "mean_and_variance.h" - -#define MAX_SQR (SQRT_U64_MAX*SQRT_U64_MAX) - -static void mean_and_variance_basic_test(struct kunit *test) -{ - struct mean_and_variance s = {}; - - mean_and_variance_update(&s, 2); - mean_and_variance_update(&s, 2); - - KUNIT_EXPECT_EQ(test, mean_and_variance_get_mean(s), 2); - KUNIT_EXPECT_EQ(test, mean_and_variance_get_variance(s), 0); - KUNIT_EXPECT_EQ(test, s.n, 2); - - mean_and_variance_update(&s, 4); - mean_and_variance_update(&s, 4); - - KUNIT_EXPECT_EQ(test, mean_and_variance_get_mean(s), 3); - KUNIT_EXPECT_EQ(test, mean_and_variance_get_variance(s), 1); - KUNIT_EXPECT_EQ(test, s.n, 4); -} - -/* - * Test values computed using a spreadsheet from the psuedocode at the bottom: - * https://fanf2.user.srcf.net/hermes/doc/antiforgery/stats.pdf - */ - -static void mean_and_variance_weighted_test(struct kunit *test) -{ - struct mean_and_variance_weighted s = { .weight = 2 }; - - mean_and_variance_weighted_update(&s, 10); - KUNIT_EXPECT_EQ(test, mean_and_variance_weighted_get_mean(s), 10); - KUNIT_EXPECT_EQ(test, mean_and_variance_weighted_get_variance(s), 0); - - mean_and_variance_weighted_update(&s, 20); - KUNIT_EXPECT_EQ(test, mean_and_variance_weighted_get_mean(s), 12); - KUNIT_EXPECT_EQ(test, mean_and_variance_weighted_get_variance(s), 18); - - mean_and_variance_weighted_update(&s, 30); - KUNIT_EXPECT_EQ(test, mean_and_variance_weighted_get_mean(s), 16); - KUNIT_EXPECT_EQ(test, mean_and_variance_weighted_get_variance(s), 72); - - s = (struct mean_and_variance_weighted) { .weight = 2 }; - - mean_and_variance_weighted_update(&s, -10); - KUNIT_EXPECT_EQ(test, mean_and_variance_weighted_get_mean(s), -10); - KUNIT_EXPECT_EQ(test, mean_and_variance_weighted_get_variance(s), 0); - - mean_and_variance_weighted_update(&s, -20); - KUNIT_EXPECT_EQ(test, mean_and_variance_weighted_get_mean(s), -12); - KUNIT_EXPECT_EQ(test, mean_and_variance_weighted_get_variance(s), 18); - - mean_and_variance_weighted_update(&s, -30); - KUNIT_EXPECT_EQ(test, mean_and_variance_weighted_get_mean(s), -16); - KUNIT_EXPECT_EQ(test, mean_and_variance_weighted_get_variance(s), 72); -} - -static void mean_and_variance_weighted_advanced_test(struct kunit *test) -{ - struct mean_and_variance_weighted s = { .weight = 8 }; - s64 i; - - for (i = 10; i <= 100; i += 10) - mean_and_variance_weighted_update(&s, i); - - KUNIT_EXPECT_EQ(test, mean_and_variance_weighted_get_mean(s), 11); - KUNIT_EXPECT_EQ(test, mean_and_variance_weighted_get_variance(s), 107); - - s = (struct mean_and_variance_weighted) { .weight = 8 }; - - for (i = -10; i >= -100; i -= 10) - mean_and_variance_weighted_update(&s, i); - - KUNIT_EXPECT_EQ(test, mean_and_variance_weighted_get_mean(s), -11); - KUNIT_EXPECT_EQ(test, mean_and_variance_weighted_get_variance(s), 107); -} - -static void do_mean_and_variance_test(struct kunit *test, - s64 initial_value, - s64 initial_n, - s64 n, - unsigned weight, - s64 *data, - s64 *mean, - s64 *stddev, - s64 *weighted_mean, - s64 *weighted_stddev) -{ - struct mean_and_variance mv = {}; - struct mean_and_variance_weighted vw = { .weight = weight }; - - for (unsigned i = 0; i < initial_n; i++) { - mean_and_variance_update(&mv, initial_value); - mean_and_variance_weighted_update(&vw, initial_value); - - KUNIT_EXPECT_EQ(test, mean_and_variance_get_mean(mv), initial_value); - KUNIT_EXPECT_EQ(test, mean_and_variance_get_stddev(mv), 0); - KUNIT_EXPECT_EQ(test, mean_and_variance_weighted_get_mean(vw), initial_value); - KUNIT_EXPECT_EQ(test, mean_and_variance_weighted_get_stddev(vw),0); - } - - for (unsigned i = 0; i < n; i++) { - mean_and_variance_update(&mv, data[i]); - mean_and_variance_weighted_update(&vw, data[i]); - - KUNIT_EXPECT_EQ(test, mean_and_variance_get_mean(mv), mean[i]); - KUNIT_EXPECT_EQ(test, mean_and_variance_get_stddev(mv), stddev[i]); - KUNIT_EXPECT_EQ(test, mean_and_variance_weighted_get_mean(vw), weighted_mean[i]); - KUNIT_EXPECT_EQ(test, mean_and_variance_weighted_get_stddev(vw),weighted_stddev[i]); - } - - KUNIT_EXPECT_EQ(test, mv.n, initial_n + n); -} - -/* Test behaviour with a single outlier, then back to steady state: */ -static void mean_and_variance_test_1(struct kunit *test) -{ - s64 d[] = { 100, 10, 10, 10, 10, 10, 10 }; - s64 mean[] = { 22, 21, 20, 19, 18, 17, 16 }; - s64 stddev[] = { 32, 29, 28, 27, 26, 25, 24 }; - s64 weighted_mean[] = { 32, 27, 22, 19, 17, 15, 14 }; - s64 weighted_stddev[] = { 38, 35, 31, 27, 24, 21, 18 }; - - do_mean_and_variance_test(test, 10, 6, ARRAY_SIZE(d), 2, - d, mean, stddev, weighted_mean, weighted_stddev); -} - -static void mean_and_variance_test_2(struct kunit *test) -{ - s64 d[] = { 100, 10, 10, 10, 10, 10, 10 }; - s64 mean[] = { 10, 10, 10, 10, 10, 10, 10 }; - s64 stddev[] = { 9, 9, 9, 9, 9, 9, 9 }; - s64 weighted_mean[] = { 32, 27, 22, 19, 17, 15, 14 }; - s64 weighted_stddev[] = { 38, 35, 31, 27, 24, 21, 18 }; - - do_mean_and_variance_test(test, 10, 6, ARRAY_SIZE(d), 2, - d, mean, stddev, weighted_mean, weighted_stddev); -} - -/* Test behaviour where we switch from one steady state to another: */ -static void mean_and_variance_test_3(struct kunit *test) -{ - s64 d[] = { 100, 100, 100, 100, 100 }; - s64 mean[] = { 22, 32, 40, 46, 50 }; - s64 stddev[] = { 32, 39, 42, 44, 45 }; - s64 weighted_mean[] = { 32, 49, 61, 71, 78 }; - s64 weighted_stddev[] = { 38, 44, 44, 41, 38 }; - - do_mean_and_variance_test(test, 10, 6, ARRAY_SIZE(d), 2, - d, mean, stddev, weighted_mean, weighted_stddev); -} - -static void mean_and_variance_test_4(struct kunit *test) -{ - s64 d[] = { 100, 100, 100, 100, 100 }; - s64 mean[] = { 10, 11, 12, 13, 14 }; - s64 stddev[] = { 9, 13, 15, 17, 19 }; - s64 weighted_mean[] = { 32, 49, 61, 71, 78 }; - s64 weighted_stddev[] = { 38, 44, 44, 41, 38 }; - - do_mean_and_variance_test(test, 10, 6, ARRAY_SIZE(d), 2, - d, mean, stddev, weighted_mean, weighted_stddev); -} - -static void mean_and_variance_fast_divpow2(struct kunit *test) -{ - s64 i; - u8 d; - - for (i = 0; i < 100; i++) { - d = 0; - KUNIT_EXPECT_EQ(test, fast_divpow2(i, d), div_u64(i, 1LLU << d)); - KUNIT_EXPECT_EQ(test, abs(fast_divpow2(-i, d)), div_u64(i, 1LLU << d)); - for (d = 1; d < 32; d++) { - KUNIT_EXPECT_EQ_MSG(test, abs(fast_divpow2(i, d)), - div_u64(i, 1 << d), "%lld %u", i, d); - KUNIT_EXPECT_EQ_MSG(test, abs(fast_divpow2(-i, d)), - div_u64(i, 1 << d), "%lld %u", -i, d); - } - } -} - -static void mean_and_variance_u128_basic_test(struct kunit *test) -{ - u128_u a = u64s_to_u128(0, U64_MAX); - u128_u a1 = u64s_to_u128(0, 1); - u128_u b = u64s_to_u128(1, 0); - u128_u c = u64s_to_u128(0, 1LLU << 63); - u128_u c2 = u64s_to_u128(U64_MAX, U64_MAX); - - KUNIT_EXPECT_EQ(test, u128_hi(u128_add(a, a1)), 1); - KUNIT_EXPECT_EQ(test, u128_lo(u128_add(a, a1)), 0); - KUNIT_EXPECT_EQ(test, u128_hi(u128_add(a1, a)), 1); - KUNIT_EXPECT_EQ(test, u128_lo(u128_add(a1, a)), 0); - - KUNIT_EXPECT_EQ(test, u128_lo(u128_sub(b, a1)), U64_MAX); - KUNIT_EXPECT_EQ(test, u128_hi(u128_sub(b, a1)), 0); - - KUNIT_EXPECT_EQ(test, u128_hi(u128_shl(c, 1)), 1); - KUNIT_EXPECT_EQ(test, u128_lo(u128_shl(c, 1)), 0); - - KUNIT_EXPECT_EQ(test, u128_hi(u128_square(U64_MAX)), U64_MAX - 1); - KUNIT_EXPECT_EQ(test, u128_lo(u128_square(U64_MAX)), 1); - - KUNIT_EXPECT_EQ(test, u128_lo(u128_div(b, 2)), 1LLU << 63); - - KUNIT_EXPECT_EQ(test, u128_hi(u128_div(c2, 2)), U64_MAX >> 1); - KUNIT_EXPECT_EQ(test, u128_lo(u128_div(c2, 2)), U64_MAX); - - KUNIT_EXPECT_EQ(test, u128_hi(u128_div(u128_shl(u64_to_u128(U64_MAX), 32), 2)), U32_MAX >> 1); - KUNIT_EXPECT_EQ(test, u128_lo(u128_div(u128_shl(u64_to_u128(U64_MAX), 32), 2)), U64_MAX << 31); -} - -static struct kunit_case mean_and_variance_test_cases[] = { - KUNIT_CASE(mean_and_variance_fast_divpow2), - KUNIT_CASE(mean_and_variance_u128_basic_test), - KUNIT_CASE(mean_and_variance_basic_test), - KUNIT_CASE(mean_and_variance_weighted_test), - KUNIT_CASE(mean_and_variance_weighted_advanced_test), - KUNIT_CASE(mean_and_variance_test_1), - KUNIT_CASE(mean_and_variance_test_2), - KUNIT_CASE(mean_and_variance_test_3), - KUNIT_CASE(mean_and_variance_test_4), - {} -}; - -static struct kunit_suite mean_and_variance_test_suite = { - .name = "mean and variance tests", - .test_cases = mean_and_variance_test_cases -}; - -kunit_test_suite(mean_and_variance_test_suite); - -MODULE_AUTHOR("Daniel B. Hill"); -MODULE_LICENSE("GPL"); diff --git a/fs/bcachefs/util.c b/fs/bcachefs/util.c index 56b815fd9fc6..d7ea95abb9df 100644 --- a/fs/bcachefs/util.c +++ b/fs/bcachefs/util.c @@ -22,9 +22,9 @@ #include #include #include +#include #include "eytzinger.h" -#include "mean_and_variance.h" #include "util.h" static const char si_units[] = "?kMGTPEZY"; diff --git a/fs/bcachefs/util.h b/fs/bcachefs/util.h index b414736d59a5..0059481995ef 100644 --- a/fs/bcachefs/util.h +++ b/fs/bcachefs/util.h @@ -17,8 +17,7 @@ #include #include #include - -#include "mean_and_variance.h" +#include #include "darray.h" diff --git a/include/linux/mean_and_variance.h b/include/linux/mean_and_variance.h new file mode 100644 index 000000000000..64df11ab422b --- /dev/null +++ b/include/linux/mean_and_variance.h @@ -0,0 +1,201 @@ +/* SPDX-License-Identifier: GPL-2.0 */ +#ifndef MEAN_AND_VARIANCE_H_ +#define MEAN_AND_VARIANCE_H_ + +#include +#include +#include +#include + +#define SQRT_U64_MAX 4294967295ULL + +/* + * u128_u: u128 user mode, because not all architectures support a real int128 + * type + * + * We don't use this version in userspace, because in userspace we link with + * Rust and rustc has issues with u128. + */ + +#if defined(__SIZEOF_INT128__) && defined(__KERNEL__) && !defined(CONFIG_PARISC) + +typedef struct { + unsigned __int128 v; +} __aligned(16) u128_u; + +static inline u128_u u64_to_u128(u64 a) +{ + return (u128_u) { .v = a }; +} + +static inline u64 u128_lo(u128_u a) +{ + return a.v; +} + +static inline u64 u128_hi(u128_u a) +{ + return a.v >> 64; +} + +static inline u128_u u128_add(u128_u a, u128_u b) +{ + a.v += b.v; + return a; +} + +static inline u128_u u128_sub(u128_u a, u128_u b) +{ + a.v -= b.v; + return a; +} + +static inline u128_u u128_shl(u128_u a, s8 shift) +{ + a.v <<= shift; + return a; +} + +static inline u128_u u128_square(u64 a) +{ + u128_u b = u64_to_u128(a); + + b.v *= b.v; + return b; +} + +#else + +typedef struct { + u64 hi, lo; +} __aligned(16) u128_u; + +/* conversions */ + +static inline u128_u u64_to_u128(u64 a) +{ + return (u128_u) { .lo = a }; +} + +static inline u64 u128_lo(u128_u a) +{ + return a.lo; +} + +static inline u64 u128_hi(u128_u a) +{ + return a.hi; +} + +/* arithmetic */ + +static inline u128_u u128_add(u128_u a, u128_u b) +{ + u128_u c; + + c.lo = a.lo + b.lo; + c.hi = a.hi + b.hi + (c.lo < a.lo); + return c; +} + +static inline u128_u u128_sub(u128_u a, u128_u b) +{ + u128_u c; + + c.lo = a.lo - b.lo; + c.hi = a.hi - b.hi - (c.lo > a.lo); + return c; +} + +static inline u128_u u128_shl(u128_u i, s8 shift) +{ + u128_u r; + + r.lo = i.lo << shift; + if (shift < 64) + r.hi = (i.hi << shift) | (i.lo >> (64 - shift)); + else { + r.hi = i.lo << (shift - 64); + r.lo = 0; + } + return r; +} + +static inline u128_u u128_square(u64 i) +{ + u128_u r; + u64 h = i >> 32, l = i & U32_MAX; + + r = u128_shl(u64_to_u128(h*h), 64); + r = u128_add(r, u128_shl(u64_to_u128(h*l), 32)); + r = u128_add(r, u128_shl(u64_to_u128(l*h), 32)); + r = u128_add(r, u64_to_u128(l*l)); + return r; +} + +#endif + +static inline u128_u u64s_to_u128(u64 hi, u64 lo) +{ + u128_u c = u64_to_u128(hi); + + c = u128_shl(c, 64); + c = u128_add(c, u64_to_u128(lo)); + return c; +} + +u128_u u128_div(u128_u n, u64 d); + +struct mean_and_variance { + s64 n; + s64 sum; + u128_u sum_squares; +}; + +/* expontentially weighted variant */ +struct mean_and_variance_weighted { + bool init; + u8 weight; /* base 2 logarithim */ + s64 mean; + u64 variance; +}; + +/** + * fast_divpow2() - fast approximation for n / (1 << d) + * @n: numerator + * @d: the power of 2 denominator. + * + * note: this rounds towards 0. + */ +static inline s64 fast_divpow2(s64 n, u8 d) +{ + return (n + ((n < 0) ? ((1 << d) - 1) : 0)) >> d; +} + +/** + * mean_and_variance_update() - update a mean_and_variance struct @s1 with a new sample @v1 + * and return it. + * @s1: the mean_and_variance to update. + * @v1: the new sample. + * + * see linked pdf equation 12. + */ +static inline void +mean_and_variance_update(struct mean_and_variance *s, s64 v) +{ + s->n++; + s->sum += v; + s->sum_squares = u128_add(s->sum_squares, u128_square(abs(v))); +} + +s64 mean_and_variance_get_mean(struct mean_and_variance s); +u64 mean_and_variance_get_variance(struct mean_and_variance s1); +u32 mean_and_variance_get_stddev(struct mean_and_variance s); + +void mean_and_variance_weighted_update(struct mean_and_variance_weighted *s, s64 v); + +s64 mean_and_variance_weighted_get_mean(struct mean_and_variance_weighted s); +u64 mean_and_variance_weighted_get_variance(struct mean_and_variance_weighted s); +u32 mean_and_variance_weighted_get_stddev(struct mean_and_variance_weighted s); + +#endif // MEAN_AND_VAIRANCE_H_ diff --git a/lib/Kconfig.debug b/lib/Kconfig.debug index 975a07f9f1cc..817ddfe132cd 100644 --- a/lib/Kconfig.debug +++ b/lib/Kconfig.debug @@ -2191,6 +2191,15 @@ config CPUMASK_KUNIT_TEST If unsure, say N. +config MEAN_AND_VARIANCE_UNIT_TEST + tristate "mean_and_variance unit tests" if !KUNIT_ALL_TESTS + depends on KUNIT + select MEAN_AND_VARIANCE + default KUNIT_ALL_TESTS + help + This option enables the kunit tests for mean_and_variance module. + If unsure, say N. + config TEST_LIST_SORT tristate "Linked list sorting test" if !KUNIT_ALL_TESTS depends on KUNIT diff --git a/lib/math/Kconfig b/lib/math/Kconfig index 0634b428d0cb..7530ae9a3584 100644 --- a/lib/math/Kconfig +++ b/lib/math/Kconfig @@ -15,3 +15,6 @@ config PRIME_NUMBERS config RATIONAL tristate + +config MEAN_AND_VARIANCE + tristate diff --git a/lib/math/Makefile b/lib/math/Makefile index 91fcdb0c9efe..8cdfa13a67ce 100644 --- a/lib/math/Makefile +++ b/lib/math/Makefile @@ -4,6 +4,8 @@ obj-y += div64.o gcd.o lcm.o int_log.o int_pow.o int_sqrt.o reciprocal_div.o obj-$(CONFIG_CORDIC) += cordic.o obj-$(CONFIG_PRIME_NUMBERS) += prime_numbers.o obj-$(CONFIG_RATIONAL) += rational.o +obj-$(CONFIG_MEAN_AND_VARIANCE) += mean_and_variance.o obj-$(CONFIG_TEST_DIV64) += test_div64.o obj-$(CONFIG_RATIONAL_KUNIT_TEST) += rational-test.o +obj-$(CONFIG_MEAN_AND_VARIANCE_UNIT_TEST) += mean_and_variance_test.o diff --git a/lib/math/mean_and_variance.c b/lib/math/mean_and_variance.c new file mode 100644 index 000000000000..ba90293204ba --- /dev/null +++ b/lib/math/mean_and_variance.c @@ -0,0 +1,164 @@ +// SPDX-License-Identifier: GPL-2.0 +/* + * Functions for incremental mean and variance. + * + * This program is free software; you can redistribute it and/or modify it + * under the terms of the GNU General Public License version 2 as published by + * the Free Software Foundation. + * + * This program is distributed in the hope that it will be useful, but WITHOUT + * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + * FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for + * more details. + * + * Copyright © 2022 Daniel B. Hill + * + * Author: Daniel B. Hill + * + * Description: + * + * This is includes some incremental algorithms for mean and variance calculation + * + * Derived from the paper: https://fanf2.user.srcf.net/hermes/doc/antiforgery/stats.pdf + * + * Create a struct and if it's the weighted variant set the w field (weight = 2^k). + * + * Use mean_and_variance[_weighted]_update() on the struct to update it's state. + * + * Use the mean_and_variance[_weighted]_get_* functions to calculate the mean and variance, some computation + * is deferred to these functions for performance reasons. + * + * see lib/math/mean_and_variance_test.c for examples of usage. + * + * DO NOT access the mean and variance fields of the weighted variants directly. + * DO NOT change the weight after calling update. + */ + +#include +#include +#include +#include +#include +#include +#include +#include + +u128_u u128_div(u128_u n, u64 d) +{ + u128_u r; + u64 rem; + u64 hi = u128_hi(n); + u64 lo = u128_lo(n); + u64 h = hi & ((u64) U32_MAX << 32); + u64 l = (hi & (u64) U32_MAX) << 32; + + r = u128_shl(u64_to_u128(div64_u64_rem(h, d, &rem)), 64); + r = u128_add(r, u128_shl(u64_to_u128(div64_u64_rem(l + (rem << 32), d, &rem)), 32)); + r = u128_add(r, u64_to_u128(div64_u64_rem(lo + (rem << 32), d, &rem))); + return r; +} +EXPORT_SYMBOL_GPL(u128_div); + +/** + * mean_and_variance_get_mean() - get mean from @s + * @s: mean and variance number of samples and their sums + */ +s64 mean_and_variance_get_mean(struct mean_and_variance s) +{ + return s.n ? div64_u64(s.sum, s.n) : 0; +} +EXPORT_SYMBOL_GPL(mean_and_variance_get_mean); + +/** + * mean_and_variance_get_variance() - get variance from @s1 + * @s1: mean and variance number of samples and sums + * + * see linked pdf equation 12. + */ +u64 mean_and_variance_get_variance(struct mean_and_variance s1) +{ + if (s1.n) { + u128_u s2 = u128_div(s1.sum_squares, s1.n); + u64 s3 = abs(mean_and_variance_get_mean(s1)); + + return u128_lo(u128_sub(s2, u128_square(s3))); + } else { + return 0; + } +} +EXPORT_SYMBOL_GPL(mean_and_variance_get_variance); + +/** + * mean_and_variance_get_stddev() - get standard deviation from @s + * @s: mean and variance number of samples and their sums + */ +u32 mean_and_variance_get_stddev(struct mean_and_variance s) +{ + return int_sqrt64(mean_and_variance_get_variance(s)); +} +EXPORT_SYMBOL_GPL(mean_and_variance_get_stddev); + +/** + * mean_and_variance_weighted_update() - exponentially weighted variant of mean_and_variance_update() + * @s: mean and variance number of samples and their sums + * @x: new value to include in the &mean_and_variance_weighted + * + * see linked pdf: function derived from equations 140-143 where alpha = 2^w. + * values are stored bitshifted for performance and added precision. + */ +void mean_and_variance_weighted_update(struct mean_and_variance_weighted *s, s64 x) +{ + // previous weighted variance. + u8 w = s->weight; + u64 var_w0 = s->variance; + // new value weighted. + s64 x_w = x << w; + s64 diff_w = x_w - s->mean; + s64 diff = fast_divpow2(diff_w, w); + // new mean weighted. + s64 u_w1 = s->mean + diff; + + if (!s->init) { + s->mean = x_w; + s->variance = 0; + } else { + s->mean = u_w1; + s->variance = ((var_w0 << w) - var_w0 + ((diff_w * (x_w - u_w1)) >> w)) >> w; + } + s->init = true; +} +EXPORT_SYMBOL_GPL(mean_and_variance_weighted_update); + +/** + * mean_and_variance_weighted_get_mean() - get mean from @s + * @s: mean and variance number of samples and their sums + */ +s64 mean_and_variance_weighted_get_mean(struct mean_and_variance_weighted s) +{ + return fast_divpow2(s.mean, s.weight); +} +EXPORT_SYMBOL_GPL(mean_and_variance_weighted_get_mean); + +/** + * mean_and_variance_weighted_get_variance() -- get variance from @s + * @s: mean and variance number of samples and their sums + */ +u64 mean_and_variance_weighted_get_variance(struct mean_and_variance_weighted s) +{ + // always positive don't need fast divpow2 + return s.variance >> s.weight; +} +EXPORT_SYMBOL_GPL(mean_and_variance_weighted_get_variance); + +/** + * mean_and_variance_weighted_get_stddev() - get standard deviation from @s + * @s: mean and variance number of samples and their sums + */ +u32 mean_and_variance_weighted_get_stddev(struct mean_and_variance_weighted s) +{ + return int_sqrt64(mean_and_variance_weighted_get_variance(s)); +} +EXPORT_SYMBOL_GPL(mean_and_variance_weighted_get_stddev); + +MODULE_AUTHOR("Daniel B. Hill"); +MODULE_LICENSE("GPL"); diff --git a/lib/math/mean_and_variance_test.c b/lib/math/mean_and_variance_test.c new file mode 100644 index 000000000000..f45591a169d8 --- /dev/null +++ b/lib/math/mean_and_variance_test.c @@ -0,0 +1,239 @@ +// SPDX-License-Identifier: GPL-2.0 +#include +#include + +#define MAX_SQR (SQRT_U64_MAX*SQRT_U64_MAX) + +static void mean_and_variance_basic_test(struct kunit *test) +{ + struct mean_and_variance s = {}; + + mean_and_variance_update(&s, 2); + mean_and_variance_update(&s, 2); + + KUNIT_EXPECT_EQ(test, mean_and_variance_get_mean(s), 2); + KUNIT_EXPECT_EQ(test, mean_and_variance_get_variance(s), 0); + KUNIT_EXPECT_EQ(test, s.n, 2); + + mean_and_variance_update(&s, 4); + mean_and_variance_update(&s, 4); + + KUNIT_EXPECT_EQ(test, mean_and_variance_get_mean(s), 3); + KUNIT_EXPECT_EQ(test, mean_and_variance_get_variance(s), 1); + KUNIT_EXPECT_EQ(test, s.n, 4); +} + +/* + * Test values computed using a spreadsheet from the psuedocode at the bottom: + * https://fanf2.user.srcf.net/hermes/doc/antiforgery/stats.pdf + */ + +static void mean_and_variance_weighted_test(struct kunit *test) +{ + struct mean_and_variance_weighted s = { .weight = 2 }; + + mean_and_variance_weighted_update(&s, 10); + KUNIT_EXPECT_EQ(test, mean_and_variance_weighted_get_mean(s), 10); + KUNIT_EXPECT_EQ(test, mean_and_variance_weighted_get_variance(s), 0); + + mean_and_variance_weighted_update(&s, 20); + KUNIT_EXPECT_EQ(test, mean_and_variance_weighted_get_mean(s), 12); + KUNIT_EXPECT_EQ(test, mean_and_variance_weighted_get_variance(s), 18); + + mean_and_variance_weighted_update(&s, 30); + KUNIT_EXPECT_EQ(test, mean_and_variance_weighted_get_mean(s), 16); + KUNIT_EXPECT_EQ(test, mean_and_variance_weighted_get_variance(s), 72); + + s = (struct mean_and_variance_weighted) { .weight = 2 }; + + mean_and_variance_weighted_update(&s, -10); + KUNIT_EXPECT_EQ(test, mean_and_variance_weighted_get_mean(s), -10); + KUNIT_EXPECT_EQ(test, mean_and_variance_weighted_get_variance(s), 0); + + mean_and_variance_weighted_update(&s, -20); + KUNIT_EXPECT_EQ(test, mean_and_variance_weighted_get_mean(s), -12); + KUNIT_EXPECT_EQ(test, mean_and_variance_weighted_get_variance(s), 18); + + mean_and_variance_weighted_update(&s, -30); + KUNIT_EXPECT_EQ(test, mean_and_variance_weighted_get_mean(s), -16); + KUNIT_EXPECT_EQ(test, mean_and_variance_weighted_get_variance(s), 72); +} + +static void mean_and_variance_weighted_advanced_test(struct kunit *test) +{ + struct mean_and_variance_weighted s = { .weight = 8 }; + s64 i; + + for (i = 10; i <= 100; i += 10) + mean_and_variance_weighted_update(&s, i); + + KUNIT_EXPECT_EQ(test, mean_and_variance_weighted_get_mean(s), 11); + KUNIT_EXPECT_EQ(test, mean_and_variance_weighted_get_variance(s), 107); + + s = (struct mean_and_variance_weighted) { .weight = 8 }; + + for (i = -10; i >= -100; i -= 10) + mean_and_variance_weighted_update(&s, i); + + KUNIT_EXPECT_EQ(test, mean_and_variance_weighted_get_mean(s), -11); + KUNIT_EXPECT_EQ(test, mean_and_variance_weighted_get_variance(s), 107); +} + +static void do_mean_and_variance_test(struct kunit *test, + s64 initial_value, + s64 initial_n, + s64 n, + unsigned weight, + s64 *data, + s64 *mean, + s64 *stddev, + s64 *weighted_mean, + s64 *weighted_stddev) +{ + struct mean_and_variance mv = {}; + struct mean_and_variance_weighted vw = { .weight = weight }; + + for (unsigned i = 0; i < initial_n; i++) { + mean_and_variance_update(&mv, initial_value); + mean_and_variance_weighted_update(&vw, initial_value); + + KUNIT_EXPECT_EQ(test, mean_and_variance_get_mean(mv), initial_value); + KUNIT_EXPECT_EQ(test, mean_and_variance_get_stddev(mv), 0); + KUNIT_EXPECT_EQ(test, mean_and_variance_weighted_get_mean(vw), initial_value); + KUNIT_EXPECT_EQ(test, mean_and_variance_weighted_get_stddev(vw),0); + } + + for (unsigned i = 0; i < n; i++) { + mean_and_variance_update(&mv, data[i]); + mean_and_variance_weighted_update(&vw, data[i]); + + KUNIT_EXPECT_EQ(test, mean_and_variance_get_mean(mv), mean[i]); + KUNIT_EXPECT_EQ(test, mean_and_variance_get_stddev(mv), stddev[i]); + KUNIT_EXPECT_EQ(test, mean_and_variance_weighted_get_mean(vw), weighted_mean[i]); + KUNIT_EXPECT_EQ(test, mean_and_variance_weighted_get_stddev(vw),weighted_stddev[i]); + } + + KUNIT_EXPECT_EQ(test, mv.n, initial_n + n); +} + +/* Test behaviour with a single outlier, then back to steady state: */ +static void mean_and_variance_test_1(struct kunit *test) +{ + s64 d[] = { 100, 10, 10, 10, 10, 10, 10 }; + s64 mean[] = { 22, 21, 20, 19, 18, 17, 16 }; + s64 stddev[] = { 32, 29, 28, 27, 26, 25, 24 }; + s64 weighted_mean[] = { 32, 27, 22, 19, 17, 15, 14 }; + s64 weighted_stddev[] = { 38, 35, 31, 27, 24, 21, 18 }; + + do_mean_and_variance_test(test, 10, 6, ARRAY_SIZE(d), 2, + d, mean, stddev, weighted_mean, weighted_stddev); +} + +static void mean_and_variance_test_2(struct kunit *test) +{ + s64 d[] = { 100, 10, 10, 10, 10, 10, 10 }; + s64 mean[] = { 10, 10, 10, 10, 10, 10, 10 }; + s64 stddev[] = { 9, 9, 9, 9, 9, 9, 9 }; + s64 weighted_mean[] = { 32, 27, 22, 19, 17, 15, 14 }; + s64 weighted_stddev[] = { 38, 35, 31, 27, 24, 21, 18 }; + + do_mean_and_variance_test(test, 10, 6, ARRAY_SIZE(d), 2, + d, mean, stddev, weighted_mean, weighted_stddev); +} + +/* Test behaviour where we switch from one steady state to another: */ +static void mean_and_variance_test_3(struct kunit *test) +{ + s64 d[] = { 100, 100, 100, 100, 100 }; + s64 mean[] = { 22, 32, 40, 46, 50 }; + s64 stddev[] = { 32, 39, 42, 44, 45 }; + s64 weighted_mean[] = { 32, 49, 61, 71, 78 }; + s64 weighted_stddev[] = { 38, 44, 44, 41, 38 }; + + do_mean_and_variance_test(test, 10, 6, ARRAY_SIZE(d), 2, + d, mean, stddev, weighted_mean, weighted_stddev); +} + +static void mean_and_variance_test_4(struct kunit *test) +{ + s64 d[] = { 100, 100, 100, 100, 100 }; + s64 mean[] = { 10, 11, 12, 13, 14 }; + s64 stddev[] = { 9, 13, 15, 17, 19 }; + s64 weighted_mean[] = { 32, 49, 61, 71, 78 }; + s64 weighted_stddev[] = { 38, 44, 44, 41, 38 }; + + do_mean_and_variance_test(test, 10, 6, ARRAY_SIZE(d), 2, + d, mean, stddev, weighted_mean, weighted_stddev); +} + +static void mean_and_variance_fast_divpow2(struct kunit *test) +{ + s64 i; + u8 d; + + for (i = 0; i < 100; i++) { + d = 0; + KUNIT_EXPECT_EQ(test, fast_divpow2(i, d), div_u64(i, 1LLU << d)); + KUNIT_EXPECT_EQ(test, abs(fast_divpow2(-i, d)), div_u64(i, 1LLU << d)); + for (d = 1; d < 32; d++) { + KUNIT_EXPECT_EQ_MSG(test, abs(fast_divpow2(i, d)), + div_u64(i, 1 << d), "%lld %u", i, d); + KUNIT_EXPECT_EQ_MSG(test, abs(fast_divpow2(-i, d)), + div_u64(i, 1 << d), "%lld %u", -i, d); + } + } +} + +static void mean_and_variance_u128_basic_test(struct kunit *test) +{ + u128_u a = u64s_to_u128(0, U64_MAX); + u128_u a1 = u64s_to_u128(0, 1); + u128_u b = u64s_to_u128(1, 0); + u128_u c = u64s_to_u128(0, 1LLU << 63); + u128_u c2 = u64s_to_u128(U64_MAX, U64_MAX); + + KUNIT_EXPECT_EQ(test, u128_hi(u128_add(a, a1)), 1); + KUNIT_EXPECT_EQ(test, u128_lo(u128_add(a, a1)), 0); + KUNIT_EXPECT_EQ(test, u128_hi(u128_add(a1, a)), 1); + KUNIT_EXPECT_EQ(test, u128_lo(u128_add(a1, a)), 0); + + KUNIT_EXPECT_EQ(test, u128_lo(u128_sub(b, a1)), U64_MAX); + KUNIT_EXPECT_EQ(test, u128_hi(u128_sub(b, a1)), 0); + + KUNIT_EXPECT_EQ(test, u128_hi(u128_shl(c, 1)), 1); + KUNIT_EXPECT_EQ(test, u128_lo(u128_shl(c, 1)), 0); + + KUNIT_EXPECT_EQ(test, u128_hi(u128_square(U64_MAX)), U64_MAX - 1); + KUNIT_EXPECT_EQ(test, u128_lo(u128_square(U64_MAX)), 1); + + KUNIT_EXPECT_EQ(test, u128_lo(u128_div(b, 2)), 1LLU << 63); + + KUNIT_EXPECT_EQ(test, u128_hi(u128_div(c2, 2)), U64_MAX >> 1); + KUNIT_EXPECT_EQ(test, u128_lo(u128_div(c2, 2)), U64_MAX); + + KUNIT_EXPECT_EQ(test, u128_hi(u128_div(u128_shl(u64_to_u128(U64_MAX), 32), 2)), U32_MAX >> 1); + KUNIT_EXPECT_EQ(test, u128_lo(u128_div(u128_shl(u64_to_u128(U64_MAX), 32), 2)), U64_MAX << 31); +} + +static struct kunit_case mean_and_variance_test_cases[] = { + KUNIT_CASE(mean_and_variance_fast_divpow2), + KUNIT_CASE(mean_and_variance_u128_basic_test), + KUNIT_CASE(mean_and_variance_basic_test), + KUNIT_CASE(mean_and_variance_weighted_test), + KUNIT_CASE(mean_and_variance_weighted_advanced_test), + KUNIT_CASE(mean_and_variance_test_1), + KUNIT_CASE(mean_and_variance_test_2), + KUNIT_CASE(mean_and_variance_test_3), + KUNIT_CASE(mean_and_variance_test_4), + {} +}; + +static struct kunit_suite mean_and_variance_test_suite = { + .name = "mean and variance tests", + .test_cases = mean_and_variance_test_cases +}; + +kunit_test_suite(mean_and_variance_test_suite); + +MODULE_AUTHOR("Daniel B. Hill"); +MODULE_LICENSE("GPL"); -- cgit v1.2.3