do local function t(t,o,n,P)local b,d,M,h=table.create(P),0,0,1 for P=1,P,1 do for b=math.max(1,(P+1)-#o),math.min(P,#t),1 do d=d+(n*t[b])*o[(P+1)-b]end local I=d%TWO_POW_24 b[P]=math.floor(I)d=(d-I)/TWO_POW_24 M=M+I*h h=h*TWO_POW_24 end return b,M end local o,n,P,b,d,M=0,{4,1,2;-2;2},4,{1},sha2_H_hi,sha2_H_lo repeat P=P+n[P%6]local h=1 repeat h=h+n[h%6]if h*h>P then local n=P^.33333333333333 local h=n*TWO_POW_40 h=t(table.create(1,math.floor(h)),b,1,2)local I,E=t(h,t(h,h,1,4),-1,4)local u=(h[2]%65536)*65536+math.floor(h[1]/256)local H=(h[1]%256)*16777216+math.floor(((E*(TWO_POW_NEG_56/3))*n)/P)if o<16 then n=math.sqrt(P)h=n*TWO_POW_40 h=t(table.create(1,math.floor(h)),b,1,2)I,E=t(h,h,-1,2)local u=(h[2]%65536)*65536+math.floor(h[1]/256)local H=(h[1]%256)*16777216+math.floor((E*TWO_POW_NEG_17)/n)local x=o%8+1 sha2_H_ext256[224][x]=H d[x],M[x]=u,H+u*hi_factor if x>7 then d,M=sha2_H_ext512_hi[384],sha2_H_ext512_lo[384]end end o=o+1 sha2_K_hi[o],sha2_K_lo[o]=u,H%K_lo_modulo+u*hi_factor break end until P%h==0 until o>79 end
--------------------------------------------------------------------------------
-- MAGIC NUMBERS CALCULATOR
--------------------------------------------------------------------------------
-- Q:
--    Is 53-bit "double" math enough to calculate square roots and cube roots of primes with 64 correct bits after decimal point?
-- A:
--    Yes, 53-bit "double" arithmetic is enough.
--    We could obtain first 40 bits by direct calculation of p^(1/3) and next 40 bits by one step of Newton's method.
do
	local function mul(src1, src2, factor, result_length)
		-- src1, src2 - long integers (arrays of digits in base TWO_POW_24)
		-- factor - small integer
		-- returns long integer result (src1 * src2 * factor) and its floating point approximation
		local result, carry, value, weight = table.create(result_length), 0, 0, 1
		for j = 1, result_length do
			for k = math.max(1, j + 1 - #src2), math.min(j, #src1) do
				carry = carry + factor * src1[k] * src2[j + 1 - k] -- "int32" is not enough for multiplication result, that's why "factor" must be of type "double"
			end

			local digit = carry % TWO_POW_24
			result[j] = math.floor(digit)
			carry = (carry - digit) / TWO_POW_24
			value = value + digit * weight
			weight = weight * TWO_POW_24
		end

		return result, value
	end

	local idx, step, p, one, sqrt_hi, sqrt_lo = 0, { 4, 1, 2, -2, 2 }, 4, { 1 }, sha2_H_hi, sha2_H_lo
	repeat
		p = p + step[p % 6]
		local d = 1
		repeat
			d = d + step[d % 6]
			if d * d > p then
				-- next prime number is found
				local root = p ^ (1 / 3)
				local R = root * TWO_POW_40
				R = mul(table.create(1, math.floor(R)), one, 1, 2)
				local _, delta = mul(R, mul(R, R, 1, 4), -1, 4)
				local hi = R[2] % 65536 * 65536 + math.floor(R[1] / 256)
				local lo = R[1] % 256 * 16777216 + math.floor(delta * (TWO_POW_NEG_56 / 3) * root / p)

				if idx < 16 then
					root = math.sqrt(p)
					R = root * TWO_POW_40
					R = mul(table.create(1, math.floor(R)), one, 1, 2)
					_, delta = mul(R, R, -1, 2)
					local hi = R[2] % 65536 * 65536 + math.floor(R[1] / 256)
					local lo = R[1] % 256 * 16777216 + math.floor(delta * TWO_POW_NEG_17 / root)
					local idx = idx % 8 + 1
					sha2_H_ext256[224][idx] = lo
					sqrt_hi[idx], sqrt_lo[idx] = hi, lo + hi * hi_factor
					if idx > 7 then
						sqrt_hi, sqrt_lo = sha2_H_ext512_hi[384], sha2_H_ext512_lo[384]
					end
				end

				idx = idx + 1
				sha2_K_hi[idx], sha2_K_lo[idx] = hi, lo % K_lo_modulo + hi * hi_factor
				break
			end
		until p % d == 0
	until idx > 79
end
