/*
 * DF-5580 differential harness for sys/platform/pc64/x86_64/__udivti3.c
 *
 * Compiles the EXACT kernel source (unmodified, #include'd below) with the
 * same gcc that built the kernel (GCC 8.3 [DragonFly]) and checks every
 * vector against native unsigned __int128 hardware division.
 *
 * Batteries:
 *   E  edge battery        - cross product of ~230 interesting operands
 *   X  exhaustive small    - all (n,d) pairs in [0,2048)^2  (4.19M)
 *   S  sr==64 boundary     - vectors pinned to shift-register == 64 +/- 4
 *                            (the exact UB shift-count boundary)
 *   R  random              - biased-random 128-bit pairs (default 2^31)
 *
 * Checks both quotient (__udivti3) and remainder (static __udivmodti4),
 * reachable because we include the .c directly.
 *
 * Exit 0 = all vectors match. Exit 1 = mismatch (printed).
 */
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <stdint.h>
#include <inttypes.h>
#include <time.h>

typedef __int128 _int128_t;              /* match machine/stdint.h */
typedef unsigned __int128 _uint128_t;

#include "__udivti3.c"                   /* the kernel source, verbatim */

static uint64_t total, mism;
static int64_t badn_a, badd_a;           /* first mismatch, for the report */

static void
chk(_uint128_t n, _uint128_t d, const char *battery)
{
	_uint128_t q1, r1, qe, re;

	qe = n / d;                          /* hardware reference */
	re = n % d;
	q1 = __udivti3(n, d);
	r1 = 0;
	(void)__udivmodti4(n, d, &r1);
	total++;
	if (q1 != qe || r1 != re) {
		unsigned long long nh = (unsigned long long)(n >> 64),
		    nl = (unsigned long long)n,
		    dh = (unsigned long long)(d >> 64),
		    dl = (unsigned long long)d;
		if (mism < 8) {
			printf("MISMATCH[%s] n=%016llx%016llx d=%016llx%016llx "
			       "got q=%016llx%016llx r=%016llx%016llx "
			       "want q=%016llx%016llx r=%016llx%016llx\n",
			       battery, nh, nl, dh, dl,
			       (unsigned long long)(q1 >> 64), (unsigned long long)q1,
			       (unsigned long long)(r1 >> 64), (unsigned long long)r1,
			       (unsigned long long)(qe >> 64), (unsigned long long)qe,
			       (unsigned long long)(re >> 64), (unsigned long long)re);
		}
		if (mism == 0) {
			badn_a = (int64_t)nl; badd_a = (int64_t)dl;
		}
		mism++;
	}
}

/* ---- interesting-value set for the edge battery ---- */
static int
build_set(_uint128_t *s, int max)
{
	int k, i = 0;
	uint64_t v[64] = {0, 1, 2, 3, 4, 5, 6, 7, 15, 16, 17, 31, 32, 33,
	    62, 63, 64, 65, 127, 128, 129, 255, 256, 257, 1023, 1024, 1025,
	    0xffff, 0x10000, 0x10001, 0x7fffffff, 0x80000000, 0x80000001,
	    0xfffffffe, 0xffffffff};
	int nv = 38;

	for (k = 0; k < nv && i < max; k++) {
		s[i++] = v[k];                       /* small values */
		s[i++] = ((_uint128_t)v[k]) << 64;   /* in the high half */
		s[i++] = (((_uint128_t)v[k]) << 64) | v[k];
		s[i++] = ~(_uint128_t)0 - v[k];
	}
	for (k = 0; k < 64 && i < max; k++) {    /* every 2^k and 2^k +/- 1 */
		uint64_t b = (k >= 64) ? 0 : ((uint64_t)1 << k);
		s[i++] = b;
		s[i++] = b + 1;
		s[i++] = b - 1;
		s[i++] = ((_uint128_t)b) << 64;
		s[i++] = (((_uint128_t)b) << 64) + 1;
		s[i++] = (((_uint128_t)b) << 64) - 1;
	}
	s[i++] = ~(_uint128_t)0;                 /* 128-bit all ones */
	s[i++] = ~(_uint128_t)0 - 1;
	s[i++] = (_uint128_t)1 << 127;
	return i;
}

static void
battery_edge(void)
{
	static _uint128_t s[1024];
	int n = build_set(s, 1024);
	int a, b;

	for (a = 0; a < n; a++)
		for (b = 0; b < n; b++) {
			if (s[b] == 0)
				continue;       /* caller contract: divisor != 0 */
			chk(s[a], s[b], "E");
		}
	printf("E edge cross-product: %d x %d operands\n", n, n);
}

static void
battery_exhaustive(uint64_t lim)
{
	uint64_t a, b;

	for (a = 0; a < lim; a++)
		for (b = 1; b < lim; b++)
			chk(a, b, "X");
	printf("X exhaustive [0,%llu)^2\n", (unsigned long long)lim);
}

/*
 * Pin sr to exactly 64 (and +/- 4) for the d.s.high==0, d.s.low non-pow2
 * path: sr = 1 + 64 + clz(d.lo) - clz(n.hi)  =>  clz(n.hi) = clz(d.lo)+1-sr+64.
 * We brute-force small n.hi/d.lo pairs whose clz difference hits the target.
 */
static int clz64(uint64_t x) { return x ? __builtin_clzl(x) : 64; }

static void
battery_sr64(void)
{
	uint64_t nh, dl, nl;
	_uint128_t n, d;
	int c, hits = 0;

	for (c = 57; c <= 69; c++) {           /* target sr value */
		for (nh = 1; nh < (1ULL << 20); nh++) {
			for (dl = 3; dl < 4096; dl += 2) {   /* odd: never pow2 */
				if (dl && (dl & (dl - 1)) == 0)
					continue;
				if (nh == 0 || (nh & (nh - 1)) == 0)
					continue;               /* n.hi pow2 is fine, keep anyway */
				if ((int)(1 + 64 + clz64(dl) - clz64(nh)) != c)
					continue;
				nl = 0xdeadbeefcafebabeULL;
				n = ((_uint128_t)nh << 64) | nl;
				d = dl;                      /* dh == 0 */
				chk(n, d, "S");
				hits++;
				if (hits > 40000)
					return;
			}
		}
	}
	printf("S sr-boundary vectors: %d\n", hits);
}

/* xorshift128+ */
static uint64_t rs0, rs1;
static uint64_t rnd64(void)
{
	uint64_t x = rs0, y = rs1;
	rs0 = y;
	x ^= x << 23;
	rs1 = x ^ y ^ (x >> 17) ^ (y >> 26);
	return rs1 + y;
}

static _uint128_t
mask_to_len(_uint128_t v, int bits)
{
	if (bits <= 0)  return 0;
	if (bits >= 128) return v;
	return v & (((_uint128_t)1 << bits) - 1);
}

static void
battery_random(uint64_t count, uint64_t seed0, uint64_t seed1)
{
	uint64_t i;
	uint64_t t0 = time(NULL);

	rs0 = seed0 ^ 0x9E3779B97F4A7C15ULL;
	rs1 = seed1 ^ 0xBF58476D1CE4E5B9ULL;
	if (!rs0) rs0 = 1;
	for (i = 0; i < count; i++) {
		_uint128_t n, d;
		int nb, db;
		uint64_t u = rnd64();

		if ((u & 15) == 0) {             /* 1/16 fully uniform */
			n = ((_uint128_t)rnd64() << 64) | rnd64();
			d = ((_uint128_t)rnd64() << 64) | rnd64();
		} else {
			db = 1 + (int)(rnd64() % 128);
			if ((u >> 8) & 1)
				nb = db - 2 + (int)(rnd64() % 5);   /* near db */
			else
				nb = (int)(rnd64() % 129);
			if (nb < 0) nb = 0;
			if (nb > 128) nb = 128;
			n = mask_to_len(((_uint128_t)rnd64() << 64) | rnd64(), nb);
			d = mask_to_len(((_uint128_t)rnd64() << 64) | rnd64(), db);
			if (d == 0) d = 1;
		}
		if (d == 0)
			continue;                    /* caller contract */
		chk(n, d, "R");
		if ((i & 0xFFFFFFF) == 0xFFFFFFF)
			fprintf(stderr, "  ... %llu/%llu (%.0f vec/s)\n",
			    (unsigned long long)i, (unsigned long long)count,
			    (double)i / (double)(time(NULL) - t0 + 1));
	}
	printf("R random vectors: %llu (seeds %llx:%llx)\n",
	    (unsigned long long)count, (unsigned long long)seed0,
	    (unsigned long long)seed1);
}

int
main(int argc, char **argv)
{
	uint64_t rcnt = 1ULL << 31;
	uint64_t seed0 = 0x243F6A8885A308D3ULL, seed1 = 0x13198A2E03707344ULL;
	int only_edge = 0;

	if (argc > 1 && strcmp(argv[1], "-e") == 0)
		only_edge = 1;
	if (argc > 2) rcnt = strtoull(argv[2], NULL, 0);
	if (argc > 4) { seed0 = strtoull(argv[3], NULL, 0);
			seed1 = strtoull(argv[4], NULL, 0); }

	printf("DF-5580 differential: kernel __udivti3.c vs native __int128\n");
	battery_edge();
	if (!only_edge) {
		battery_exhaustive(2048);
		battery_sr64();
		battery_random(rcnt, seed0, seed1);
	}
	printf("TOTAL checked=%llu mismatches=%llu\n",
	    (unsigned long long)total, (unsigned long long)mism);
	printf("%s\n", mism == 0 ? "PASS: all vectors match hardware division"
	                         : "FAIL: mismatch against hardware division");
	return mism == 0 ? 0 : 1;
}