/*
 * merkle.c -- RFC 9162 Merkle tree hash, inclusion and consistency proofs,
 * checked against published test vectors, the RFC 9162 Section 2.1.5
 * example, and exhaustive small-tree tests; plus two construction flaws
 * (missing domain separation, Bitcoin-style last-node duplication).
 *
 * Build: gcc -O2 -Wall -Wextra -std=c11 -o merkle merkle.c -lcrypto
 * Run:   ./merkle
 */
#include <openssl/evp.h>
#include <stdint.h>
#include <stdio.h>
#include <stdlib.h>
#include <string.h>

typedef struct { unsigned char b[32]; } hash_t;
typedef struct { const unsigned char *p; size_t len; } entry_t;

static EVP_MD_CTX *g_ctx;
static int g_fail;

#define CHECK(cond, ...)                                        \
    do {                                                        \
        if (!(cond)) {                                          \
            g_fail++;                                           \
            fprintf(stderr, "FAIL %s:%d: ", __FILE__, __LINE__); \
            fprintf(stderr, __VA_ARGS__);                       \
            fputc('\n', stderr);                                \
        }                                                       \
    } while (0)

static void sha256_parts(const void *const *parts, const size_t *lens, int k,
                         hash_t *out)
{
    unsigned int n = 0;
    if (EVP_DigestInit_ex(g_ctx, EVP_sha256(), NULL) != 1) abort();
    for (int i = 0; i < k; i++)
        if (EVP_DigestUpdate(g_ctx, parts[i], lens[i]) != 1) abort();
    if (EVP_DigestFinal_ex(g_ctx, out->b, &n) != 1 || n != 32) abort();
}

static int heq(const hash_t *a, const hash_t *b) { return memcmp(a->b, b->b, 32) == 0; }

/* ---------------- RFC 9162 Section 2.1.1 ---------------- */

static const unsigned char LEAF_PREFIX = 0x00, NODE_PREFIX = 0x01;

static hash_t leaf_hash(const entry_t *e)
{
    const void *p[2] = { &LEAF_PREFIX, e->p };
    size_t l[2] = { 1, e->len };
    hash_t h;
    sha256_parts(p, l, 2, &h);
    return h;
}

static hash_t node_hash(const hash_t *left, const hash_t *right)
{
    const void *p[3] = { &NODE_PREFIX, left->b, right->b };
    size_t l[3] = { 1, 32, 32 };
    hash_t h;
    sha256_parts(p, l, 3, &h);
    return h;
}

/* largest power of two strictly smaller than n (n >= 2) */
static size_t split_point(size_t n)
{
    size_t k = 1;
    while (k << 1 < n) k <<= 1;
    return k;
}

static hash_t mth(const entry_t *d, size_t n)
{
    if (n == 0) {
        hash_t h;
        sha256_parts(NULL, NULL, 0, &h);
        return h;
    }
    if (n == 1) return leaf_hash(&d[0]);
    size_t k = split_point(n);
    hash_t l = mth(d, k), r = mth(d + k, n - k);
    return node_hash(&l, &r);
}

/* PATH(m, D_n), Section 2.1.3.1; appends to out[*len] */
static void path(size_t m, const entry_t *d, size_t n, hash_t *out, size_t *len)
{
    if (n <= 1) return;
    size_t k = split_point(n);
    if (m < k) {
        path(m, d, k, out, len);
        out[(*len)++] = mth(d + k, n - k);
    } else {
        path(m - k, d + k, n - k, out, len);
        out[(*len)++] = mth(d, k);
    }
}

/* SUBPROOF(m, D_n, b), Section 2.1.4.1 */
static void subproof(size_t m, const entry_t *d, size_t n, int b, hash_t *out,
                     size_t *len)
{
    if (m == n) {
        if (!b) out[(*len)++] = mth(d, n);
        return;
    }
    size_t k = split_point(n);
    if (m <= k) {
        subproof(m, d, k, b, out, len);
        out[(*len)++] = mth(d + k, n - k);
    } else {
        subproof(m - k, d + k, n - k, 0, out, len);
        out[(*len)++] = mth(d, k);
    }
}

static void consistency_proof(size_t m, const entry_t *d, size_t n, hash_t *out,
                              size_t *len)
{
    *len = 0;
    subproof(m, d, n, 1, out, len);
}

/* Section 2.1.3.2 */
static int verify_inclusion(uint64_t leaf_index, uint64_t tree_size,
                            const hash_t *leaf, const hash_t *proof, size_t plen,
                            const hash_t *root)
{
    if (leaf_index >= tree_size) return 0;
    uint64_t fn = leaf_index, sn = tree_size - 1;
    hash_t r = *leaf;
    for (size_t i = 0; i < plen; i++) {
        if (sn == 0) return 0;
        if ((fn & 1) || fn == sn) {
            r = node_hash(&proof[i], &r);
            while (!(fn & 1) && fn != 0) { fn >>= 1; sn >>= 1; }
        } else {
            r = node_hash(&r, &proof[i]);
        }
        fn >>= 1; sn >>= 1;
    }
    return sn == 0 && heq(&r, root);
}

/* Section 2.1.4.2, for 0 < first < second */
static int verify_consistency(uint64_t first, uint64_t second,
                              const hash_t *first_hash, const hash_t *second_hash,
                              const hash_t *proof, size_t plen)
{
    if (first == 0 || first >= second || plen == 0) return 0;
    hash_t p[130];
    size_t n = 0;
    if ((first & (first - 1)) == 0) p[n++] = *first_hash;
    if (plen + n > sizeof p / sizeof p[0]) return 0;
    memcpy(p + n, proof, plen * sizeof(hash_t));
    n += plen;

    uint64_t fn = first - 1, sn = second - 1;
    while (fn & 1) { fn >>= 1; sn >>= 1; }
    hash_t fr = p[0], sr = p[0];
    for (size_t i = 1; i < n; i++) {
        if (sn == 0) return 0;
        if ((fn & 1) || fn == sn) {
            fr = node_hash(&p[i], &fr);
            sr = node_hash(&p[i], &sr);
            while (!(fn & 1) && fn != 0) { fn >>= 1; sn >>= 1; }
        } else {
            sr = node_hash(&sr, &p[i]);
        }
        fn >>= 1; sn >>= 1;
    }
    return heq(&fr, first_hash) && heq(&sr, second_hash) && sn == 0;
}

/* ---------------- proof sizes without hashing ---------------- */

static size_t path_len(size_t m, size_t n)
{
    if (n <= 1) return 0;
    size_t k = split_point(n);
    return 1 + (m < k ? path_len(m, k) : path_len(m - k, n - k));
}

static size_t subproof_len(size_t m, size_t n, int b)
{
    if (m == n) return b ? 0 : 1;
    size_t k = split_point(n);
    return 1 + (m <= k ? subproof_len(m, k, b) : subproof_len(m - k, n - k, 0));
}

static int bitlen64(uint64_t x) { return x ? 64 - __builtin_clzll(x) : 0; }

/* closed form used by transparency-dev/merkle proof.decompInclProof */
static size_t path_len_closed(uint64_t m, uint64_t n)
{
    int inner = bitlen64(m ^ (n - 1));
    return (size_t)inner + (size_t)__builtin_popcountll(m >> inner);
}

static int ceil_log2(uint64_t n) { return n <= 1 ? 0 : bitlen64(n - 1); }

/* ---------------- test vectors ---------------- */

static int hex2bin(const char *s, unsigned char *out, size_t cap)
{
    size_t n = strlen(s) / 2;
    if (n > cap) return -1;
    for (size_t i = 0; i < n; i++) {
        unsigned v;
        if (sscanf(s + 2 * i, "%2x", &v) != 1) return -1;
        out[i] = (unsigned char)v;
    }
    return (int)n;
}

static hash_t H(const char *hex)
{
    hash_t h;
    if (hex2bin(hex, h.b, 32) != 32) abort();
    return h;
}

static void print_hash(const hash_t *h)
{
    for (int i = 0; i < 32; i++) printf("%02x", h->b[i]);
}

/* transparency-dev/merkle v0.0.2, testonly/constants.go (Apache-2.0) */
static const char *LEAF_HEX[8] = {
    "", "00", "10", "2021", "3031", "40414243", "5051525354555657",
    "606162636465666768696a6b6c6d6e6f",
};
static const char *ROOT_HEX[9] = {
    "e3b0c44298fc1c149afbf4c8996fb92427ae41e4649b934ca495991b7852b855",
    "6e340b9cffb37a989ca544e6bb780a2c78901d3fb33738768511a30617afa01d",
    "fac54203e7cc696cf0dfcb42c92a1d9dbaf70ad9e621f4bd8d98662f00e3c125",
    "aeb6bcfe274b70a14fb067a5e5578264db0fa9b51af5e0ba159158f329e06e77",
    "d37ee418976dd95753c1c73862b9398fa2a2cf9b4ff0fdfe8b30cd95209614b7",
    "4e3bbb1f7b478dcfe71fb631631519a3bca12c9aefca1612bfce4c13a86264d4",
    "76e67dadbcdf1e10e1b74ddc608abd2f98dfb16fbce75277b5232a127f2087ef",
    "ddb89be403809e325750d3d263cd78929c2942b7942a34b77e122c9594a74c8c",
    "5dc9da79a70659a9ad559cb701ded9a2ab9d823aad2f4960cfe370eff4604328",
};

/* proof/verify_test.go; index is 0-based here (the Go file stores index+1) */
struct pv { size_t a, n, len; const char *h[3]; };
static const struct pv INCL[] = {
    { 0, 8, 3, { "96a296d224f285c67bee93c30f8a309157f0daa35dc5b87e410b78630a09cfc7",
                 "5f083f0a1a33ca076a95279832580db3e0ef4584bdff1f54c8a360f50de3031e",
                 "6b47aaf29ee3c2af9af889bc1fb9254dabd31177f16232dd6aab035ca39bf6e4" } },
    { 5, 8, 3, { "bc1a0643b12e4d2d7c77918f44e0f4f79a838b6cf9ec5b5c283e1f4d88599e6b",
                 "ca854ea128ed050b41b35ffc1b87b8eb2bde461e9e3b5596ece6b9d5975a0ae0",
                 "d37ee418976dd95753c1c73862b9398fa2a2cf9b4ff0fdfe8b30cd95209614b7" } },
    { 2, 3, 1, { "fac54203e7cc696cf0dfcb42c92a1d9dbaf70ad9e621f4bd8d98662f00e3c125" } },
    { 1, 5, 3, { "6e340b9cffb37a989ca544e6bb780a2c78901d3fb33738768511a30617afa01d",
                 "5f083f0a1a33ca076a95279832580db3e0ef4584bdff1f54c8a360f50de3031e",
                 "bc1a0643b12e4d2d7c77918f44e0f4f79a838b6cf9ec5b5c283e1f4d88599e6b" } },
};
/* here a = first tree size */
static const struct pv CONS[] = {
    { 1, 8, 3, { "96a296d224f285c67bee93c30f8a309157f0daa35dc5b87e410b78630a09cfc7",
                 "5f083f0a1a33ca076a95279832580db3e0ef4584bdff1f54c8a360f50de3031e",
                 "6b47aaf29ee3c2af9af889bc1fb9254dabd31177f16232dd6aab035ca39bf6e4" } },
    { 6, 8, 3, { "0ebc5d3437fbe2db158b9f126a1d118e308181031d0a949f8dededebc558ef6a",
                 "ca854ea128ed050b41b35ffc1b87b8eb2bde461e9e3b5596ece6b9d5975a0ae0",
                 "d37ee418976dd95753c1c73862b9398fa2a2cf9b4ff0fdfe8b30cd95209614b7" } },
    { 2, 5, 2, { "5f083f0a1a33ca076a95279832580db3e0ef4584bdff1f54c8a360f50de3031e",
                 "bc1a0643b12e4d2d7c77918f44e0f4f79a838b6cf9ec5b5c283e1f4d88599e6b" } },
    { 6, 7, 3, { "0ebc5d3437fbe2db158b9f126a1d118e308181031d0a949f8dededebc558ef6a",
                 "b08693ec2e721597130641e8211e7eedccb4c26413963eee6c1e2ed16ffb1a5f",
                 "d37ee418976dd95753c1c73862b9398fa2a2cf9b4ff0fdfe8b30cd95209614b7" } },
};

static unsigned char g_leafbuf[8][16];
static entry_t g_vec[8];

static void load_vectors(void)
{
    for (int i = 0; i < 8; i++) {
        int n = hex2bin(LEAF_HEX[i], g_leafbuf[i], sizeof g_leafbuf[i]);
        if (n < 0) abort();
        g_vec[i].p = g_leafbuf[i];
        g_vec[i].len = (size_t)n;
    }
}

static void test_vectors(void)
{
    for (size_t n = 0; n <= 8; n++) {
        hash_t r = mth(g_vec, n), want = H(ROOT_HEX[n]);
        CHECK(heq(&r, &want), "root of size %zu", n);
    }
    hash_t pr[64];
    for (size_t t = 0; t < sizeof INCL / sizeof INCL[0]; t++) {
        const struct pv *v = &INCL[t];
        size_t len = 0;
        path(v->a, g_vec, v->n, pr, &len);
        CHECK(len == v->len, "inclusion (%zu,%zu) length", v->a, v->n);
        for (size_t i = 0; i < len && i < v->len; i++) {
            hash_t w = H(v->h[i]);
            CHECK(heq(&pr[i], &w), "inclusion (%zu,%zu) element %zu", v->a, v->n, i);
        }
        hash_t lh = leaf_hash(&g_vec[v->a]), root = H(ROOT_HEX[v->n]);
        CHECK(verify_inclusion(v->a, v->n, &lh, pr, len, &root), "verify incl vector %zu", t);
    }
    for (size_t t = 0; t < sizeof CONS / sizeof CONS[0]; t++) {
        const struct pv *v = &CONS[t];
        size_t len = 0;
        consistency_proof(v->a, g_vec, v->n, pr, &len);
        CHECK(len == v->len, "consistency (%zu,%zu) length", v->a, v->n);
        for (size_t i = 0; i < len && i < v->len; i++) {
            hash_t w = H(v->h[i]);
            CHECK(heq(&pr[i], &w), "consistency (%zu,%zu) element %zu", v->a, v->n, i);
        }
        hash_t r1 = H(ROOT_HEX[v->a]), r2 = H(ROOT_HEX[v->n]);
        CHECK(verify_consistency(v->a, v->n, &r1, &r2, pr, len), "verify cons vector %zu", t);
    }
    printf("[vectors] roots for sizes 0..8, %zu inclusion and %zu consistency proofs "
           "match transparency-dev/merkle v0.0.2\n",
           sizeof INCL / sizeof INCL[0], sizeof CONS / sizeof CONS[0]);
}

/* RFC 9162 Section 2.1.5: node names a..l for the 7-leaf tree */
static void test_rfc_example(void)
{
    hash_t lf[7];
    for (int i = 0; i < 7; i++) lf[i] = leaf_hash(&g_vec[i]);
    hash_t a = lf[0], b = lf[1], c = lf[2], d = lf[3], e = lf[4], f = lf[5], j = lf[6];
    hash_t g = node_hash(&a, &b), h = node_hash(&c, &d), i = node_hash(&e, &f);
    hash_t k = node_hash(&g, &h), l = node_hash(&i, &j), root = node_hash(&k, &l);
    hash_t hash0 = node_hash(&g, &c), hash2 = node_hash(&k, &i);

    hash_t t = mth(g_vec, 7);
    CHECK(heq(&t, &root), "rfc example root");
    t = mth(g_vec, 3);
    CHECK(heq(&t, &hash0), "rfc example hash0");
    t = mth(g_vec, 4);
    CHECK(heq(&t, &k), "rfc example hash1");
    t = mth(g_vec, 6);
    CHECK(heq(&t, &hash2), "rfc example hash2");

    struct { size_t m; int len; const hash_t *want[4]; const char *s; } incl[] = {
        { 0, 3, { &b, &h, &l }, "[b, h, l]" },
        { 3, 3, { &c, &g, &l }, "[c, g, l]" },
        { 4, 3, { &f, &j, &k }, "[f, j, k]" },
        { 6, 2, { &i, &k }, "[i, k]" },
    };
    hash_t pr[16];
    for (size_t x = 0; x < sizeof incl / sizeof incl[0]; x++) {
        size_t len = 0;
        path(incl[x].m, g_vec, 7, pr, &len);
        int ok = (int)len == incl[x].len;
        for (size_t y = 0; ok && y < len; y++) ok = heq(&pr[y], incl[x].want[y]);
        CHECK(ok, "rfc inclusion d%zu", incl[x].m);
        CHECK(verify_inclusion(incl[x].m, 7, &lf[incl[x].m], pr, len, &root),
              "rfc inclusion verify d%zu", incl[x].m);
        printf("[rfc9162 2.1.5] inclusion proof for d%zu = %s: %s\n", incl[x].m, incl[x].s,
               ok ? "ok" : "MISMATCH");
    }
    struct { size_t m; int len; const hash_t *want[4]; const hash_t *old; const char *s; } cons[] = {
        { 3, 4, { &c, &d, &g, &l }, &hash0, "[c, d, g, l]" },
        { 4, 1, { &l }, &k, "[l]" },
        { 6, 3, { &i, &j, &k }, &hash2, "[i, j, k]" },
    };
    for (size_t x = 0; x < sizeof cons / sizeof cons[0]; x++) {
        size_t len = 0;
        consistency_proof(cons[x].m, g_vec, 7, pr, &len);
        int ok = (int)len == cons[x].len;
        for (size_t y = 0; ok && y < len; y++) ok = heq(&pr[y], cons[x].want[y]);
        CHECK(ok, "rfc consistency %zu", cons[x].m);
        CHECK(verify_consistency(cons[x].m, 7, cons[x].old, &root, pr, len),
              "rfc consistency verify %zu", cons[x].m);
        printf("[rfc9162 2.1.5] consistency proof PROOF(%zu, D7) = %s: %s\n", cons[x].m,
               cons[x].s, ok ? "ok" : "MISMATCH");
    }
    printf("[rfc9162 2.1.5] root of D7 = ");
    print_hash(&root);
    printf("\n");
}

/* ---------------- exhaustive small trees ---------------- */

#define NMAX 128

static void test_exhaustive(void)
{
    static unsigned char buf[NMAX][8];
    entry_t d[NMAX];
    for (size_t i = 0; i < NMAX; i++) {
        memcpy(buf[i], &i, sizeof i < 8 ? sizeof i : 8);
        d[i].p = buf[i];
        d[i].len = 8;
    }
    hash_t roots[NMAX + 1], lf[NMAX];
    for (size_t n = 1; n <= NMAX; n++) roots[n] = mth(d, n);
    for (size_t i = 0; i < NMAX; i++) lf[i] = leaf_hash(&d[i]);

    long n_incl = 0, n_cons = 0, n_forged = 0, n_rejected = 0;
    hash_t pr[64];
    for (size_t n = 1; n <= NMAX; n++) {
        for (size_t m = 0; m < n; m++) {
            size_t len = 0;
            path(m, d, n, pr, &len);
            n_incl++;
            CHECK(len == path_len_closed(m, n), "closed form n=%zu m=%zu", n, m);
            CHECK((int)len <= ceil_log2(n), "incl bound n=%zu m=%zu", n, m);
            CHECK(verify_inclusion(m, n, &lf[m], pr, len, &roots[n]), "incl n=%zu m=%zu", n, m);
            /* every single-bit flip in the proof must be rejected */
            for (size_t e = 0; e < len; e++) {
                pr[e].b[e % 32] ^= 1;
                n_forged++;
                n_rejected += !verify_inclusion(m, n, &lf[m], pr, len, &roots[n]);
                pr[e].b[e % 32] ^= 1;
            }
            /* same proof, wrong index or wrong size */
            if (m + 1 < n) {
                n_forged++;
                n_rejected += !verify_inclusion(m + 1, n, &lf[m], pr, len, &roots[n]);
            }
            if (n + 1 <= NMAX) {
                n_forged++;
                n_rejected += !verify_inclusion(m, n + 1, &lf[m], pr, len, &roots[n + 1]);
            }
        }
        for (size_t m = 1; m < n; m++) {
            size_t len = 0;
            consistency_proof(m, d, n, pr, &len);
            n_cons++;
            CHECK(len == subproof_len(m, n, 1), "subproof_len n=%zu m=%zu", n, m);
            CHECK((int)len <= ceil_log2(n) + 1, "cons bound n=%zu m=%zu", n, m);
            CHECK(verify_consistency(m, n, &roots[m], &roots[n], pr, len), "cons n=%zu m=%zu", n, m);
            for (size_t e = 0; e < len; e++) {
                pr[e].b[e % 32] ^= 1;
                n_forged++;
                n_rejected += !verify_consistency(m, n, &roots[m], &roots[n], pr, len);
                pr[e].b[e % 32] ^= 1;
            }
            /* a proof for (m, n) must not pass as a proof for another old size */
            if (m + 1 < n) {
                n_forged++;
                n_rejected += !verify_consistency(m + 1, n, &roots[m + 1], &roots[n], pr, len);
            }
        }
    }
    CHECK(n_forged == n_rejected, "forged proofs accepted: %ld", n_forged - n_rejected);
    printf("[exhaustive n<=%d] %ld inclusion + %ld consistency proofs verified; "
           "%ld tampered proofs, %ld rejected\n",
           NMAX, n_incl, n_cons, n_forged, n_rejected);
}

/* Section 2.1.2: stack algorithm must agree with the recursive MTH */
static void test_stack_algorithm(void)
{
    enum { N = 1000 };
    static unsigned char buf[N][4];
    entry_t d[N];
    for (int i = 0; i < N; i++) {
        memcpy(buf[i], &i, 4);
        d[i].p = buf[i];
        d[i].len = 4;
    }
    hash_t stack[64];
    int bad = 0;
    for (size_t n = 1; n <= N; n++) {
        int sp = 0;
        for (size_t i = 0; i < n; i++) {
            stack[sp++] = leaf_hash(&d[i]);
            for (size_t x = i; x & 1; x >>= 1) {
                hash_t r = stack[--sp], l = stack[--sp];
                stack[sp++] = node_hash(&l, &r);
            }
        }
        while (sp > 1) {
            hash_t r = stack[--sp], l = stack[--sp];
            stack[sp++] = node_hash(&l, &r);
        }
        hash_t m = mth(d, n);
        bad += !heq(&stack[0], &m);
    }
    CHECK(bad == 0, "stack algorithm mismatches: %d", bad);
    printf("[stack 2.1.2] incremental stack root == recursive MTH for n=1..%d\n", N);
}

/* ---------------- proof sizes ---------------- */

static void proof_sizes(void)
{
    enum { LIM = 4096 };
    long checked = 0;
    int max_c = 0, tight = 0;
    for (size_t n = 1; n <= LIM; n++) {
        for (size_t m = 0; m < n; m++) {
            CHECK(path_len(m, n) == path_len_closed(m, n), "closed form n=%zu m=%zu", n, m);
            checked++;
        }
        for (size_t m = 1; m < n; m++) {
            int c = (int)subproof_len(m, n, 1);
            CHECK(c <= ceil_log2(n) + 1, "cons bound n=%zu m=%zu", n, m);
            if (c > max_c) max_c = c;
            if (c == ceil_log2(n) + 1) tight++;
        }
    }
    printf("[sizes] closed-form inclusion length == recursive PATH length for all "
           "%ld (m, n), n<=%d\n", checked, LIM);
    printf("[sizes] consistency length <= ceil(log2 n)+1 for all 0<m<n<=%d; "
           "max %d; bound attained by %d pairs\n", LIM, max_c, tight);

    const uint64_t ns[] = { 7, 1000, 1024, 1025, 1000000, 1000000000ULL };
    printf("[sizes] inclusion proof length over all leaf indices (SHA-256, 32 B/hash)\n");
    printf("  %12s %9s %5s %5s %8s %10s\n", "n", "ceil_lg", "min", "max", "mean", "max_bytes");
    for (size_t t = 0; t < sizeof ns / sizeof ns[0]; t++) {
        uint64_t n = ns[t];
        size_t mn = 99, mx = 0;
        double sum = 0;
        for (uint64_t m = 0; m < n; m++) {
            size_t l = path_len_closed(m, n);
            if (l < mn) mn = l;
            if (l > mx) mx = l;
            sum += (double)l;
        }
        printf("  %12llu %9d %5zu %5zu %8.3f %10zu\n", (unsigned long long)n, ceil_log2(n),
               mn, mx, sum / (double)n, mx * 32);
    }
    const uint64_t cn[] = { 1000, 1000000 };
    printf("[sizes] consistency proof length over all old sizes 0<m<n\n");
    printf("  %12s %9s %5s %5s %8s\n", "n", "ceil_lg+1", "min", "max", "mean");
    for (size_t t = 0; t < sizeof cn / sizeof cn[0]; t++) {
        uint64_t n = cn[t];
        size_t mn = 99, mx = 0;
        double sum = 0;
        for (uint64_t m = 1; m < n; m++) {
            size_t l = subproof_len(m, n, 1);
            if (l < mn) mn = l;
            if (l > mx) mx = l;
            sum += (double)l;
        }
        printf("  %12llu %9d %5zu %5zu %8.3f\n", (unsigned long long)n, ceil_log2(n) + 1, mn, mx,
               sum / (double)(n - 1));
    }
}

/* ---------------- flaw 1: no domain separation ---------------- */

static hash_t plain_leaf(const unsigned char *p, size_t len)
{
    const void *pp[1] = { p };
    size_t l[1] = { len };
    hash_t h;
    sha256_parts(pp, l, 1, &h);
    return h;
}

static hash_t plain_node(const hash_t *l, const hash_t *r)
{
    const void *pp[2] = { l->b, r->b };
    size_t ll[2] = { 32, 32 };
    hash_t h;
    sha256_parts(pp, ll, 2, &h);
    return h;
}

static void flaw_second_preimage(void)
{
    const char *txt[4] = { "d0", "d1", "d2", "d3" };
    hash_t lf[4];
    for (int i = 0; i < 4; i++) lf[i] = plain_leaf((const unsigned char *)txt[i], 2);
    hash_t g = plain_node(&lf[0], &lf[1]), h = plain_node(&lf[2], &lf[3]);
    hash_t root4 = plain_node(&g, &h);

    /* two forged 64-byte "leaves": the concatenated child hashes */
    unsigned char x0[64], x1[64];
    memcpy(x0, lf[0].b, 32); memcpy(x0 + 32, lf[1].b, 32);
    memcpy(x1, lf[2].b, 32); memcpy(x1 + 32, lf[3].b, 32);
    hash_t fx0 = plain_leaf(x0, 64), fx1 = plain_leaf(x1, 64);
    hash_t root2 = plain_node(&fx0, &fx1);
    int same_plain = heq(&root4, &root2);
    /* the one-element "proof" [h] for leaf x0 at index 0 of a size-2 tree */
    hash_t via_proof = plain_node(&fx0, &h);
    int proof_ok = heq(&via_proof, &root4);

    entry_t e4[4], e2[2];
    for (int i = 0; i < 4; i++) { e4[i].p = (const unsigned char *)txt[i]; e4[i].len = 2; }
    e2[0].p = x0; e2[0].len = 64;
    e2[1].p = x1; e2[1].len = 64;
    hash_t r4 = mth(e4, 4), r2 = mth(e2, 2);
    int same_rfc = heq(&r4, &r2);

    CHECK(same_plain && proof_ok && !same_rfc, "second preimage demo");
    printf("[flaw: no domain separation] root([d0,d1,d2,d3]) == root([H(d0)||H(d1), "
           "H(d2)||H(d3)]): %s; forged inclusion proof accepted: %s\n",
           same_plain ? "yes" : "no", proof_ok ? "yes" : "no");
    printf("[flaw: no domain separation] same two lists under RFC 9162 prefixes: roots %s\n",
           same_rfc ? "EQUAL" : "differ");
}

/* ---------------- flaw 2: Bitcoin-style duplication ---------------- */

static hash_t sha256d(const void *p, size_t len)
{
    const void *pp[1] = { p };
    size_t l[1] = { len };
    hash_t h1, h2;
    sha256_parts(pp, l, 1, &h1);
    pp[0] = h1.b; l[0] = 32;
    sha256_parts(pp, l, 1, &h2);
    return h2;
}

/* mirrors ComputeMerkleRoot() in Bitcoin Core v28.0 src/consensus/merkle.cpp */
static hash_t btc_root(const hash_t *in, size_t n, int *mutated)
{
    hash_t *v = malloc((n + 1) * sizeof *v);
    if (!v) abort();
    memcpy(v, in, n * sizeof *v);
    int mut = 0;
    while (n > 1) {
        for (size_t pos = 0; pos + 1 < n; pos += 2)
            if (heq(&v[pos], &v[pos + 1])) mut = 1;
        if (n & 1) { v[n] = v[n - 1]; n++; }
        for (size_t i = 0; i < n / 2; i++) {
            unsigned char cat[64];
            memcpy(cat, v[2 * i].b, 32);
            memcpy(cat + 32, v[2 * i + 1].b, 32);
            v[i] = sha256d(cat, 64);
        }
        n /= 2;
    }
    hash_t r = v[0];
    free(v);
    *mutated = mut;
    return r;
}

static void flaw_duplication(void)
{
    char txt[8][8];
    hash_t txid[8];
    const int order[8] = { 1, 2, 3, 4, 5, 6, 5, 6 };
    for (int i = 0; i < 8; i++) {
        snprintf(txt[i], sizeof txt[i], "tx%d", order[i]);
        txid[i] = sha256d(txt[i], strlen(txt[i]));
    }
    int m6, m8;
    hash_t r6 = btc_root(txid, 6, &m6), r8 = btc_root(txid, 8, &m8);
    int same_btc = heq(&r6, &r8);

    entry_t e[8];
    for (int i = 0; i < 8; i++) { e[i].p = txid[i].b; e[i].len = 32; }
    hash_t f6 = mth(e, 6), f8 = mth(e, 8);
    int same_rfc = heq(&f6, &f8);

    CHECK(same_btc && !m6 && m8 && !same_rfc, "duplication demo");
    printf("[flaw: duplicate last node] Bitcoin-style root([1..6]) == root([1..6,5,6]): %s; "
           "mutation flag: %d vs %d\n", same_btc ? "yes" : "no", m6, m8);
    printf("[flaw: duplicate last node] same txid lists under RFC 9162: roots %s\n",
           same_rfc ? "EQUAL" : "differ");
}

int main(void)
{
    g_ctx = EVP_MD_CTX_new();
    if (!g_ctx) return 2;
    load_vectors();
    test_vectors();
    test_rfc_example();
    test_stack_algorithm();
    test_exhaustive();
    proof_sizes();
    flaw_second_preimage();
    flaw_duplication();
    EVP_MD_CTX_free(g_ctx);
    printf(g_fail ? "FAILED: %d checks\n" : "all checks passed\n", g_fail);
    return g_fail ? 1 : 0;
}