Breaking truncated LCG with https://github.com/rkm0959/Inequality_Solving_with_CVP

solver.sage

from sage.modules.free_module_integer import IntegerLattice

# Directly taken from rbtree's LLL repository
# From https://oddcoder.com/LOL-34c3/, https://hackmd.io/@hakatashi/B1OM7HFVI
def Babai_CVP(mat, target):
	M = IntegerLattice(mat, lll_reduce=True).reduced_basis
	G = M.gram_schmidt()[0]
	diff = target
	for i in reversed(range(G.nrows())):
		diff -=  M[i] * ((diff * G[i]) / (G[i] * G[i])).round()
	return target - diff


def solve(mat, lb, ub, weight = None):
	num_var  = mat.nrows()
	num_ineq = mat.ncols()

	max_element = 0 
	for i in range(num_var):
		for j in range(num_ineq):
			max_element = max(max_element, abs(mat[i, j]))

	if weight == None:
		weight = num_ineq * max_element

    # sanity checker
	if len(lb) != num_ineq:
		print("Fail: len(lb) != num_ineq")
		return

	if len(ub) != num_ineq:
		print("Fail: len(ub) != num_ineq")
		return

	for i in range(num_ineq):
		if lb[i] > ub[i]:
			print("Fail: lb[i] > ub[i] at index", i)
			return

    	# heuristic for number of solutions
	DET = 0

	if num_var == num_ineq:
		DET = abs(mat.det())
		num_sol = 1
		for i in range(num_ineq):
			num_sol *= (ub[i] - lb[i])
		if DET == 0:
			print("Zero Determinant")
		else:
			num_sol //= DET
			# + 1 added in for the sake of not making it zero...
			print("Expected Number of Solutions : ", num_sol + 1)

	# scaling process begins
	max_diff = max([ub[i] - lb[i] for i in range(num_ineq)])
	applied_weights = []

	for i in range(num_ineq):
		ineq_weight = weight if lb[i] == ub[i] else max_diff // (ub[i] - lb[i])
		applied_weights.append(ineq_weight)
		for j in range(num_var):
			mat[j, i] *= ineq_weight
		lb[i] *= ineq_weight
		ub[i] *= ineq_weight

	# Solve CVP
	target = vector([(lb[i] + ub[i]) // 2 for i in range(num_ineq)])
	result = Babai_CVP(mat, target)

	for i in range(num_ineq):
		if (lb[i] <= result[i] <= ub[i]) == False:
			print("Fail : inequality does not hold after solving")
			break
    
    	# recover x
	fin = None

	if DET != 0:
		mat = mat.transpose()
		fin = mat.solve_right(result)
	
	## recover your result
	return result, applied_weights, fin

truncated_lcg.sage

a = 0x876387638763
b = 0xDEADBEEF
m = 2 ** 64
s = 314231423142
trunc = 2 ** 48
n = 4  # number of states


def lcg():
    global s
    s = (a * s + b) % m
    return s


xs = [lcg() for _ in range(n)]  # xs, actual states
zs = [x % trunc for x in xs]  # zs, truncated lower bits
ys = [x - z for x, z in zip(xs, zs)]  # ys, remaining upper bits, known to attacker
assert all(y + z == x for x, y, z in zip(xs, ys, zs))

ks = []

for z1, z2, y1, y2 in zip(zs, zs[1:], ys, ys[1:]):
    rhs = y2 - a * y1 - b
    lhs = z1 * a - z2
    assert (rhs - lhs) % m == 0
    ks.append((rhs - lhs) // m)

M = Matrix(ZZ, n * 2 - 1, n - 1 + n)
for i in range(n - 1):
    M[i, i] = a
    M[i + 1, i] = -1
    M[i, i + n - 1] = 1
    M[n + i, i] = m
M[n - 1, -1] = 1
print(M)
print(M.dimensions())

print("expected solution", vector(zs + ks) * M)

ub = []

for y1, y2 in zip(ys, ys[1:]):
    ub.append(y2 - a * y1 - b)

lb = list(ub)

lb += [0] * n
ub += [trunc] * n

print("lb", lb)
print("ub", ub)

load(
    "solver.sage"
)  # https://raw.githubusercontent.com/rkm0959/Inequality_Solving_with_CVP/main/solver.sage

result, applied_weights, fin = solve(M, lb, ub)
assert vector(fin[:n]) == vector(zs), "failed to find zs"
print(fin[:n])

Output example:

[     148861542958947                    0                    0                    1                    0                    0                    0]
[                  -1      148861542958947                    0                    0                    1                    0                    0]
[                   0                   -1      148861542958947                    0                    0                    1                    0]
[                   0                    0                   -1                    0                    0                    0                    1]
[18446744073709551616                    0                    0                    0                    0                    0                    0]
[                   0 18446744073709551616                    0                    0                    0                    0                    0]
[                   0                    0 18446744073709551616                    0                    0                    0                    0]
(7, 7)
expected solution (-2384867795891451280358037978463983, -1103876558545954614637563724283631, -995269686663205541008598889250543, 185138152245537, 192147573321138, 245302493464773, 233877055260190)
lb [-2384867795891451280358037978463983, -1103876558545954614637563724283631, -995269686663205541008598889250543, 0, 0, 0, 0]
ub [-2384867795891451280358037978463983, -1103876558545954614637563724283631, -995269686663205541008598889250543, 281474976710656, 281474976710656, 281474976710656, 281474976710656]
Expected Number of Solutions :  1
(185138152245537, 192147573321138, 245302493464773, 233877055260190)