# vrandom - 硬件随机数生成库
# 基于 x86 RDRAND / RDSEED 指令，使用内嵌汇编实现
import t, c
import vipermath


# 重试次数限制（RDRAND 通常不需要重试，RDSEED 可能需要）
RDRAND_RETRY_LIMIT: t.CDefine = 10
RDSEED_RETRY_LIMIT: t.CDefine = 100


# ============================================================
# RDRAND - 硬件随机数（伪随机数发生器输出）
# ============================================================

def rdrand16() -> t.CUInt16T:
    val: t.CUInt32T = 0
    c.Asm(f"rdrand {c.AsmOut(val, t.ASM_DESCR.OUTPUT_REG)}",
        out=[c.AsmOut(val, t.ASM_DESCR.OUTPUT_REG)],
        op=[t.ASM_DESCR.CLOBBER_CC])
    return t.CUInt16T(val)


def rdrand32() -> t.CUInt32T:
    val: t.CUInt32T = 0
    c.Asm(f"rdrand {c.AsmOut(val, t.ASM_DESCR.OUTPUT_REG)}",
        out=[c.AsmOut(val, t.ASM_DESCR.OUTPUT_REG)],
        op=[t.ASM_DESCR.CLOBBER_CC])
    return val


def rdrand64() -> t.CUInt64T:
    val: t.CUInt64T = 0
    c.Asm(f"rdrand {c.AsmOut(val, t.ASM_DESCR.OUTPUT_REG)}",
        out=[c.AsmOut(val, t.ASM_DESCR.OUTPUT_REG)],
        op=[t.ASM_DESCR.CLOBBER_CC])
    return val


# ============================================================
# RDSEED - 硬件随机种子（熵源直接输出）
# ============================================================

def rdseed16() -> t.CUInt16T:
    val: t.CUInt32T = 0
    c.Asm(f"rdseed {c.AsmOut(val, t.ASM_DESCR.OUTPUT_REG)}",
        out=[c.AsmOut(val, t.ASM_DESCR.OUTPUT_REG)],
        op=[t.ASM_DESCR.CLOBBER_CC])
    return t.CUInt16T(val)


def rdseed32() -> t.CUInt32T:
    val: t.CUInt32T = 0
    c.Asm(f"rdseed {c.AsmOut(val, t.ASM_DESCR.OUTPUT_REG)}",
        out=[c.AsmOut(val, t.ASM_DESCR.OUTPUT_REG)],
        op=[t.ASM_DESCR.CLOBBER_CC])
    return val


def rdseed64() -> t.CUInt64T:
    val: t.CUInt64T = 0
    c.Asm(f"rdseed {c.AsmOut(val, t.ASM_DESCR.OUTPUT_REG)}",
        out=[c.AsmOut(val, t.ASM_DESCR.OUTPUT_REG)],
        op=[t.ASM_DESCR.CLOBBER_CC])
    return val


# ============================================================
# Step 版本（带成功标志，匹配 GCC builtin 接口）
# 返回 1 表示成功，0 表示失败；随机值写入 *p
# ============================================================

def rdrand16_step(p: t.CUInt16T | t.CPtr) -> t.CInt:
    val: t.CUInt32T = 0
    ok: t.CUInt32T = 0
    c.Asm(f"""rdrand {c.AsmOut(val, t.ASM_DESCR.OUTPUT_REG)}
              setc al
              movzx {c.AsmOut(ok, t.ASM_DESCR.OUTPUT_REG)}, al""",
        out=[c.AsmOut(val, t.ASM_DESCR.OUTPUT_REG), c.AsmOut(ok, t.ASM_DESCR.OUTPUT_REG)],
        op=[t.ASM_DESCR.CLOBBER_CC, t.ASM_DESCR.CLOBBER_AL])
    p[0] = t.CUInt16T(val)
    return t.CInt(ok)


def rdrand32_step(p: t.CUInt32T | t.CPtr) -> t.CInt:
    val: t.CUInt32T = 0
    ok: t.CUInt32T = 0
    c.Asm(f"""rdrand {c.AsmOut(val, t.ASM_DESCR.OUTPUT_REG)}
              setc al
              movzx {c.AsmOut(ok, t.ASM_DESCR.OUTPUT_REG)}, al""",
        out=[c.AsmOut(val, t.ASM_DESCR.OUTPUT_REG), c.AsmOut(ok, t.ASM_DESCR.OUTPUT_REG)],
        op=[t.ASM_DESCR.CLOBBER_CC, t.ASM_DESCR.CLOBBER_AL])
    p[0] = val
    return t.CInt(ok)


def rdrand64_step(p: t.CUInt64T | t.CPtr) -> t.CInt:
    val: t.CUInt64T = 0
    ok: t.CUInt32T = 0
    c.Asm(f"""rdrand {c.AsmOut(val, t.ASM_DESCR.OUTPUT_REG)}
              setc al
              movzx {c.AsmOut(ok, t.ASM_DESCR.OUTPUT_REG)}, al""",
        out=[c.AsmOut(val, t.ASM_DESCR.OUTPUT_REG), c.AsmOut(ok, t.ASM_DESCR.OUTPUT_REG)],
        op=[t.ASM_DESCR.CLOBBER_CC, t.ASM_DESCR.CLOBBER_AL])
    p[0] = val
    return t.CInt(ok)


def rdseed16_step(p: t.CUInt16T | t.CPtr) -> t.CInt:
    val: t.CUInt32T = 0
    ok: t.CUInt32T = 0
    c.Asm(f"""rdseed {c.AsmOut(val, t.ASM_DESCR.OUTPUT_REG)}
              setc al
              movzx {c.AsmOut(ok, t.ASM_DESCR.OUTPUT_REG)}, al""",
        out=[c.AsmOut(val, t.ASM_DESCR.OUTPUT_REG), c.AsmOut(ok, t.ASM_DESCR.OUTPUT_REG)],
        op=[t.ASM_DESCR.CLOBBER_CC, t.ASM_DESCR.CLOBBER_AL])
    p[0] = t.CUInt16T(val)
    return t.CInt(ok)


def rdseed32_step(p: t.CUInt32T | t.CPtr) -> t.CInt:
    val: t.CUInt32T = 0
    ok: t.CUInt32T = 0
    c.Asm(f"""rdseed {c.AsmOut(val, t.ASM_DESCR.OUTPUT_REG)}
              setc al
              movzx {c.AsmOut(ok, t.ASM_DESCR.OUTPUT_REG)}, al""",
        out=[c.AsmOut(val, t.ASM_DESCR.OUTPUT_REG), c.AsmOut(ok, t.ASM_DESCR.OUTPUT_REG)],
        op=[t.ASM_DESCR.CLOBBER_CC, t.ASM_DESCR.CLOBBER_AL])
    p[0] = val
    return t.CInt(ok)


def rdseed64_step(p: t.CUInt64T | t.CPtr) -> t.CInt:
    val: t.CUInt64T = 0
    ok: t.CUInt32T = 0
    c.Asm(f"""rdseed {c.AsmOut(val, t.ASM_DESCR.OUTPUT_REG)}
              setc al
              movzx {c.AsmOut(ok, t.ASM_DESCR.OUTPUT_REG)}, al""",
        out=[c.AsmOut(val, t.ASM_DESCR.OUTPUT_REG), c.AsmOut(ok, t.ASM_DESCR.OUTPUT_REG)],
        op=[t.ASM_DESCR.CLOBBER_CC, t.ASM_DESCR.CLOBBER_AL])
    p[0] = val
    return t.CInt(ok)


# ============================================================
# 带重试的便捷接口（自动重试直到成功或达到上限）
# ============================================================

def random_u16() -> t.CUInt16T:
    r: t.CUInt16T = 0
    i: t.CInt = 0
    while i < RDRAND_RETRY_LIMIT:
        if rdrand16_step(c.Addr(r)) == 1:
            return r
        i += 1
    return 0


def random_u32() -> t.CUInt32T:
    r: t.CUInt32T = 0
    i: t.CInt = 0
    while i < RDRAND_RETRY_LIMIT:
        if rdrand32_step(c.Addr(r)) == 1:
            return r
        i += 1
    return 0


def random_u64() -> t.CUInt64T:
    r: t.CUInt64T = 0
    i: t.CInt = 0
    while i < RDRAND_RETRY_LIMIT:
        if rdrand64_step(c.Addr(r)) == 1:
            return r
        i += 1
    return 0


def seed_u16() -> t.CUInt16T:
    r: t.CUInt16T = 0
    i: t.CInt = 0
    while i < RDSEED_RETRY_LIMIT:
        if rdseed16_step(c.Addr(r)) == 1:
            return r
        i += 1
    return 0


def seed_u32() -> t.CUInt32T:
    r: t.CUInt32T = 0
    i: t.CInt = 0
    while i < RDSEED_RETRY_LIMIT:
        if rdseed32_step(c.Addr(r)) == 1:
            return r
        i += 1
    return 0


def seed_u64() -> t.CUInt64T:
    r: t.CUInt64T = 0
    i: t.CInt = 0
    while i < RDSEED_RETRY_LIMIT:
        if rdseed64_step(c.Addr(r)) == 1:
            return r
        i += 1
    return 0


# ============================================================
# Python random 标准库函数实现
# 默认类型：t.CUInt64T (u64), t.CDouble (f64), t.CInt (i32)
# ============================================================

# ---------- 一、随机整数 ----------

def randint(a: t.CInt, b: t.CInt) -> t.CInt:
    return t.CInt(a + t.CInt(random_u64() % t.CUInt64T(b - a + 1)))


def randrange(start: t.CInt, stop: t.CInt, step: t.CInt = 1) -> t.CInt:
    width: t.CInt = 0
    if step > 0:
        width = (stop - start + step - 1) // step
    else:
        width = (start - stop - step - 1) // (-step)
    return t.CInt(start + t.CInt(random_u64() % t.CUInt64T(width)) * step)


def getrandbits(k: t.CInt) -> t.CUInt64T:
    r: t.CUInt64T = rdrand64()
    if k >= 64:
        return r
    mask: t.CUInt64T = (t.CUInt64T(1) << k) - t.CUInt64T(1)
    return r & mask


# ---------- 二、随机浮点数 ----------

def random() -> t.CDouble:
    return t.CDouble(rdrand64() >> 11) * (t.CDouble(1.0) / t.CDouble(9007199254740992.0))


def uniform(a: t.CDouble, b: t.CDouble) -> t.CDouble:
    return a + (b - a) * random()


def triangular(low: t.CDouble, high: t.CDouble, mode: t.CDouble) -> t.CDouble:
    u: t.CDouble = random()
    c: t.CDouble = (mode - low) / (high - low)
    if u > c:
        u = t.CDouble(1.0) - u
        c = t.CDouble(1.0) - c
        tmp: t.CDouble = low
        low = high
        high = tmp
    return low + (high - low) * vipermath.sqrt(u * c)


def betavariate(alpha: t.CDouble, beta: t.CDouble) -> t.CDouble:
    y: t.CDouble = gammavariate(alpha, t.CDouble(1.0))
    if y != t.CDouble(0.0):
        return y / (y + gammavariate(beta, t.CDouble(1.0)))
    return t.CDouble(0.0)


def expovariate(lambd: t.CDouble) -> t.CDouble:
    u: t.CDouble = random()
    if u < t.CDouble(1e-300):
        u = t.CDouble(1e-300)
    return -vipermath.log(t.CDouble(1.0) - u) / lambd


def gammavariate(alpha: t.CDouble, beta: t.CDouble) -> t.CDouble:
    u: t.CDouble = t.CDouble(0.0)
    d: t.CDouble = t.CDouble(0.0)
    c: t.CDouble = t.CDouble(0.0)
    x: t.CDouble = t.CDouble(0.0)
    v: t.CDouble = t.CDouble(0.0)
    if alpha <= t.CDouble(0.0) or beta <= t.CDouble(0.0):
        return t.CDouble(0.0)
    if alpha < t.CDouble(1.0):
        u = random()
        if u < t.CDouble(1e-300):
            u = t.CDouble(1e-300)
        return gammavariate(alpha + t.CDouble(1.0), beta) * vipermath.pow(u, t.CDouble(1.0) / alpha)
    d = alpha - t.CDouble(1.0) / t.CDouble(3.0)
    c = t.CDouble(1.0) / vipermath.sqrt(t.CDouble(9.0) * d)
    while 1:
        x = gauss(t.CDouble(0.0), t.CDouble(1.0))
        v = t.CDouble(1.0) + c * x
        if v > t.CDouble(0.0):
            v = v * v * v
            u = random()
            if u < t.CDouble(1.0) - t.CDouble(0.0331) * x * x * x * x:
                return d * v * beta
            if vipermath.log(u) < t.CDouble(0.5) * x * x + d * (t.CDouble(1.0) - v + vipermath.log(v)):
                return d * v * beta


def gauss(mu: t.CDouble, sigma: t.CDouble) -> t.CDouble:
    u1: t.CDouble = random()
    u2: t.CDouble = random()
    if u1 < t.CDouble(1e-300):
        u1 = t.CDouble(1e-300)
    z0: t.CDouble = vipermath.sqrt(t.CDouble(-2.0) * vipermath.log(u1)) * vipermath.cos(t.CDouble(2.0) * vipermath.U_M_PI * u2)
    return mu + sigma * z0


def normalvariate(mu: t.CDouble, sigma: t.CDouble) -> t.CDouble:
    NV_MAGICCONST: t.CDouble = t.CDouble(1.7155277699214135)
    u1: t.CDouble = t.CDouble(0.0)
    u2: t.CDouble = t.CDouble(0.0)
    z: t.CDouble = t.CDouble(0.0)
    zz: t.CDouble = t.CDouble(0.0)
    while 1:
        u1 = random()
        u2 = random()
        if u2 < t.CDouble(1e-300):
            continue
        z = NV_MAGICCONST * (u1 - t.CDouble(0.5)) / u2
        zz = z * z / t.CDouble(4.0)
        if zz <= t.CDouble(0.0) - vipermath.log(u2):
            break
    return mu + z * sigma


def lognormvariate(mu: t.CDouble, sigma: t.CDouble) -> t.CDouble:
    return vipermath.exp(gauss(mu, sigma))


def vonmisesvariate(mu: t.CDouble, kappa: t.CDouble) -> t.CDouble:
    if kappa <= t.CDouble(1e-6):
        return vipermath.U_M_PI * t.CDouble(2.0) * random()
    s: t.CDouble = t.CDouble(0.5) / kappa
    r: t.CDouble = s + vipermath.sqrt(t.CDouble(1.0) + s * s)
    u1: t.CDouble = t.CDouble(0.0)
    u2: t.CDouble = t.CDouble(0.0)
    u3: t.CDouble = t.CDouble(0.0)
    z: t.CDouble = t.CDouble(0.0)
    d: t.CDouble = t.CDouble(0.0)
    q: t.CDouble = t.CDouble(0.0)
    f: t.CDouble = t.CDouble(0.0)
    theta: t.CDouble = t.CDouble(0.0)
    two_pi: t.CDouble = vipermath.U_M_PI * t.CDouble(2.0)
    while 1:
        u1 = random()
        z = vipermath.cos(vipermath.U_M_PI * u1)
        d = t.CDouble(1.0) / (r + z)
        u2 = random()
        if u2 < t.CDouble(1.0) - d * d:
            break
        if u2 < (t.CDouble(1.0) - z) * vipermath.exp(kappa * (z - t.CDouble(1.0)) - d * kappa * r):
            break
    q = t.CDouble(1.0) / r
    f = (q + z) / (t.CDouble(1.0) + q * z)
    u3 = random()
    if u3 > t.CDouble(0.5):
        theta = mu + vipermath.acos(f)
    else:
        theta = mu - vipermath.acos(f)
    while theta < t.CDouble(0.0):
        theta += two_pi
    while theta >= two_pi:
        theta -= two_pi
    return theta


def paretovariate(alpha: t.CDouble) -> t.CDouble:
    u: t.CDouble = random()
    if u < t.CDouble(1e-300):
        u = t.CDouble(1e-300)
    return t.CDouble(1.0) / vipermath.pow(t.CDouble(1.0) - u, t.CDouble(1.0) / alpha)


def weibullvariate(alpha: t.CDouble, beta: t.CDouble) -> t.CDouble:
    u: t.CDouble = random()
    if u < t.CDouble(1e-300):
        u = t.CDouble(1e-300)
    return alpha * vipermath.pow(t.CDouble(-1.0) * vipermath.log(t.CDouble(1.0) - u), t.CDouble(1.0) / beta)