TSGCTF 2021 Writeups

5 Crypto + 2 PPC = 7 solves. Favorite Challenges = This is DSA, Lumberjack Against Nature.

 

Beginner's Crypto

from secret import e
from Crypto.Util.number import getStrongPrime, isPrime
 
p = getStrongPrime(1024)
q = getStrongPrime(1024)
N = p * q
phi = (p - 1) * (q - 1)
 
with open('flag.txt', 'rb') as f:
    flag = int.from_bytes(f.read(), 'big')
 
assert(isPrime(e))
assert(isPrime(e + 2))
assert(isPrime(e + 4))
 
e1 = pow(e, 0x10001, phi)
e2 = pow(e + 2, 0x10001, phi)
e3 = pow(e + 4, 0x10001, phi)
 
c1 = pow(flag, e1, N)
c2 = pow(flag, e2, N)
c3 = pow(flag, e3, N)
 
print(f'p = {p}')
print(f'q = {q}')
print(f'c1 = {c1}')
print(f'c2 = {c2}')
print(f'c3 = {c3}')
 

 

As it will be soon mentioned in the flag of this challenge, we will solve this without $p, q$. 

Since $e, e+2, e+4$ is all prime and at least one of them has to be a multiple of $3$, we see $e=3$. 

Now we can see that $c_1, c_2, c_3$ can be computed as $$c_1 \equiv m^{3^{65537}} \pmod{n}, \quad c_2 \equiv m^{5^{65537}} \pmod{n}, \quad c_3 \equiv m^{7^{65537}} \pmod{n}$$ Since $\gcd(3^{65537}, 5^{65537}) = 1$, we can use extended Euclidean algorithm on the exponents to find $m$. 

 

Since $3^{65537}$ and $5^{65537}$ are large numbers, it is recommended to use GMPY or Sagemath's integers.

 

p = 167710954518007348037383082265231465648795974011761905177264545864288011527333715495850532989338171489309608848431113452814709692343039027970312735521415071265608660628968391884287240987858607818275329135585153511665148279408708087727501421558738163577629329044315775019460018956186674179846621352371150072281
q = 130354329753344570838569091064852072757046774566775609047544069941246798511317343102715733555464772099991834579660053860799207243561908291522943696711982657846373844514551117658179060004064010647453939332217996817580433587341521331941287365948919907797478197717562721233289937471168288241937022054501586986443
c1 = 2560344169447809042170685026483682125499025654554670516499742981486615082413150123244985585751880264831112089324011804397189638172356179296987581738515619297036118472798499254785110885662931526277474101787493114656242031264678448394380651657330967744585361662315313462698221954777506355498445242300193032704972074020068699180111637362566860530694807230108024167631423062629721393506643291591971626450262144814424411172618188943774725105690851574922374544865628890948773274109561622040022136970632948166009941425683576381155722191980954262373394704682297682490061906408535261437100820855976015526295573831744458528440
c2 = 9041231631916227099296501948589424780380702196870972231114747229225732542137483840187783630590878594711315671224997985975031038623195921968945234067183003568830416719957054703139219879265482072634572699299971785171441858501409377942183918216246312330291820452436486171483461790388518159980027140392750222843449604265528929311978655519463562520038992870162220913137870017065557254099767583925177889051326144499369420594398043223307161794788085369471538477803421726790780799629276012701406231535048423554314287152404245482928538931953627397633165453319078105028671410039195670727134471011040601278722143504641171853743
c3 = 3193069356811106774640161554961405075257002069448498144279061282023129342916422283816661697787316681475161942522570615456264481238277711114193792510286127129056376618422336477707825009085263623755329815306483253646072909132096678270667136193038337386976289222105363398033633185639402128949635525665502328717781718263894690234837016959581149138917064108193064639981137359869717065147934752707676203651598070046066514316196771853484143158367616177332902152347890310640338106015356361617700741042461419248117687350565094928451141103632305400493998164788411031832078388030194992306440474662871408938796429927990102583837
n = p * q
 
# using n only
 
e1 = pow(Integer(3), 0x10001)
e2 = pow(Integer(5), 0x10001)
 
g1 = inverse_mod(e1, e2)
g2 = Integer((e1 * g1 - 1) // e2)
 
flag = (pow(c1, g1, n) * inverse_mod(Integer(pow(c2, g2, n)), n)) % n 
 
print(long_to_bytes(int(flag)))

 

Minimalist's Private

from Crypto.Util.number import isPrime
from random import randrange
from secret import p, q, L, e, d
 
class RSA:
    def __init__(self, p, q, L, e, d):
        assert(isPrime(p) and isPrime(q))
        self.N = p * q
        self.L = L
        self.e = e
        self.d = d
 
        # these are the normal RSA conditions
        for _ in range(100):
            assert(pow(randrange(1, self.N), self.L, self.N) == 1)
        assert(self.e * self.d % self.L == 1)
 
        # minimal is the best
        assert(self.L * self.L <= 10000 * self.N)
 
    def gen_private_key(self):
        return (self.N, self.d)
 
    def gen_public_key(self):
        return (self.N, self.e)
 
    def encrypt(self, msg):
        return pow(msg, self.e, self.N)
 
    def decrypt(self, c):
        return pow(c, self.d, self.N)
 
flag = open('flag.txt', 'rb').read()
msg = int.from_bytes(flag, byteorder='big')
assert(msg < p * q)
 
rsa = RSA(p, q, L, e, d)
encrypted = rsa.encrypt(msg)
assert(rsa.decrypt(encrypted) == msg)
 
print(f'N, e = {rsa.gen_public_key()}')
print(f'c = {encrypted}')
 

 

We see that $L \ge \text{lcm}(p-1, q-1)$. If we let $G = \gcd(p-1, q-1)$ and $p-1 = Ga$, $q-1 = Gb$, we have $$(Gab)^2 \le L^2 \le 10^4 n = 10^4(Ga + 1)(Gb + 1)$$ which shows us that $$ab \le 10^4 \left(1 + \frac{1}{Ga} \right) \left(1 + \frac{1}{Gb} \right) \le 4 \cdot 10^4$$ There are not a lot of pairs $(a, b)$ with $ab \le 4 \cdot 10^4$, in fact, the number of pairs $(a, b)$ with $ab \le n$ is around $\mathcal{O}(n \log n)$, so we can brute force all such $(a, b)$ and try to solve for $G$ with the quadratic equation $$n = (Ga+1)(Gb+1)$$

 

def inthroot(a, nn):
    return a.nth_root(nn, truncate_mode=True)[0]
 
n, e = (1108103848370322618250236235096737547381026108763302516499816051432801216813681568375319595638932562835292256776016949573972732881586209527824393027428125964599378845347154409633878436868422905300799413838645686430352484534761305185938956589612889463246508935994301443576781452904666072122465831585156151, 65537)
c = 254705401581808316199469430068831357413481187288921393400711004895837418302514065107811330660948313420965140464021505716810909691650540609799307500282957438243553742714371028405100267860418626513481187170770328765524251710154676478766892336610743824131087888798846367363259860051983889314134196889300426
 
for i in tqdm(range(1, 5000)):
    for j in range(1, 5000 // i + 5):
        aa = i * j 
        bb = i + j 
        cc = 1 - n 
        try:
            tt = (-bb + inthroot(Integer(bb * bb - 4 * aa * cc), 2)) // (2 * aa)
            p = i * tt + 1
            q = j * tt + 1 
            if p * q == n:
                print("HEY")
                print(p, q)
                phi = LCM(p - 1, q - 1)
                d = inverse(e, phi)
                print(d)
                print(long_to_bytes(pow(c, d, n)))
                exit()
        except:
            pass
 
 

 

Baba is Flag

require 'openssl'
require 'digest'
 
STDOUT.sync = true
 
class OpenSSL::PKey::EC::Point
  def xy
    n = to_bn(:uncompressed).to_i
    mask = (1 << group.degree) - 1
    return (n >> group.degree) & mask, n & mask
  end
  alias_method :+, :add
  alias_method :*, :mul
end
 
class ECDSA
  def initialize
    @curve = OpenSSL::PKey::EC::Group.new('secp256k1')
    @G = @curve.generator
    @n = @curve.order.to_i
    @d = OpenSSL::BN.rand(@curve.degree).to_i
    @Q = @G * @d
  end
 
  def inv(x)
    x.pow(@n - 2, @n)
  end
 
  def sign(msg)
    z = Digest::SHA256.hexdigest(msg).hex
    k = OpenSSL::BN.rand(@curve.degree / 3).to_s.unpack1('H*').hex
    x, y = (@G * k).xy
 
    # We all like hacks, ain't we?
    # s = (z + x * @d) * inv(k) % @n
    s = (z + @d) * inv(k) % @n
 
    return x, s
  end
 
  def verify(msg, x, s)
    return false if x % @n == 0 || s % @n == 0
    z = Digest::SHA256.hexdigest(msg).hex
 
    # ditto
    # x2, y2 = (@G * (z * inv(s)) + @Q * (x * inv(s))).xy
    x2, y2 = (@G * (z * inv(s)) + @Q * inv(s)).xy
 
    return x == x2
  end
end
 
ecdsa = ECDSA.new
 
5.times do
  puts <<~EOS
    1. Sign
    2. Find rule
    3. Exit
  EOS
 
  print 'choice? '
 
  case gets.chomp
  when '1'
    x, s = ecdsa.sign('Baba')
    puts 'Baba is:'
    puts "x = #{x}"
    puts "s = #{s}"
  when '2'
    print 'Which rule do you want to know? '; msg = gets.chomp
    print 'x? '; x = gets.to_i
    print 's? '; s = gets.to_i
 
    if ecdsa.verify(msg, x, s)
      if msg == 'Baba'
        puts 'Baba is you'
      elsif msg == 'Flag'
        puts "Flag is #{ENV['FLAG']}"
      else
        puts 'Not Found :('
      end
    else
      puts 'Invalid :('
    end
  else
    exit
  end
end
 
puts 'You is defeat.'
 

 

Here, we want to forge $(x, s)$ such that $$s \cdot \text{lift_x}(x) = Q + z \cdot G$$ where $z$ is the hash of the message and $\text{lift_x}(x)$ is the point with $x$-coordinate equal to $x$. By asking for the signature of 'Baba', we get a pair $(x, s)$ that corresponds to the hash of 'Baba'. Since $x, s, z, G$ are all known, we can recover the value of $Q$. 

 

Now we can simply take $s = 1$ and $x$ to be the $x$-coordinate of $Q + z \cdot G$, where $z$ is the hash of 'Flag' to make a valid signature.

This solves the problem, and this vuln is of course from the missing $x$ in the signature formula.

 

r = remote('34.146.212.53', 65434)
 
p = (1 << 256) - (1 << 32) - (1 << 9) - (1 << 8) - (1 << 7) - (1 << 6) - (1 << 4) - 1
 
E = EllipticCurve(GF(p), [0, 7])
Gx = 0x79BE667EF9DCBBAC55A06295CE870B07029BFCDB2DCE28D959F2815B16F81798
Gy = 0x483ADA7726A3C4655DA4FBFC0E1108A8FD17B448A68554199C47D08FFB10D4B8
 
G = E(Gx, Gy)
n = E.order()
print(isPrime(n))
 
h1 = bytes_to_long(hashlib.sha256(b'Baba').digest())
h2 = bytes_to_long(hashlib.sha256(b'Flag').digest())
 
for i in range(3):
    r.recvline()
r.sendline(b"1")
r.recvline()
X1 = int(r.recvline().split()[-1])
S1 = int(r.recvline().split()[-1])
 
print(X1)
print(S1)
 
 
target1 = S1 * E.lift_x(GF(p)(X1))
 
target2 = target1 + (h2 - h1) * G
for i in range(3):
    r.recvline()
r.sendline(b"2")
r.sendline("Flag")
r.sendline(str(int(target2.xy()[0])))
r.sendline(b"1")
print(r.recvline())
print(r.recvline())

 

This is DSA

# See also https://github.com/tsg-ut/pycryptodome
from Crypto.PublicKey import DSA
from Crypto.Signature import DSS
from Crypto.Hash import SHA256
from Crypto.Util.number import getPrime
from Crypto.Random.random import randrange
from base64 import b64decode
from signal import alarm
import os
 
alarm(15)
 
q = getPrime(256)
print(f'q = {q}')
 
p = int(input('p? '))
h = int(input('h? '))
 
g = pow(h, (p - 1) // q, p)
x = randrange(q)
y = pow(g, x, p)
 
print(f'g = {g}')
print(f'y = {y}')
 
dsa = DSA.construct((y, g, p, q, x))
dss = DSS.new(dsa, 'fips-186-3')
 
print('Thank you for helping me with DSA! Now give me the base64-encoded signature of sha256("flag")')
sign = b64decode(input('sign? '))
 
dss.verify(SHA256.new(b'flag'), sign)
print(f"Awesome! {os.environ.get('FLAG')}")

 

I took ridiculously long to solve this challenge, for several reasons. Here's my story. 

 

Step 1 : removing all "irrational" ideas

In standard DSA, we are forced to solve a discrete logarithm in a subgroup of $\mathbb{F}_p^{\star}$ with order $q$. 

Since $q$ is 256 bits, this is very hard to do, and there are no tricks that I know of to make a nice $p$ that makes it possible. 

Therefore, directly attacking this problem is not possible. We have to go around it. 

 

The first thing that caught my eye was that there was no check that $(y, g, p, q, x)$ was valid DSA tuple. 

If I could do send some nasty values on $p, h$, then maybe this problem would be very easy to solve. 

 

At this point, about 10 minutes have passed. I decided to look at pycryptodome's code for DSA construction.

Then https://github.com/Legrandin/pycryptodome/blob/master/lib/Crypto/PublicKey/DSA.py#L489 happened. 

It turns out that pycryptodome does check everything automatically, even when not specified. 20 minutes have passed.

 

Step 2 : actually knowing what the challenge is

After wasting an additional 40 minutes, I found that the library was patched. 

https://github.com/tsg-ut/pycryptodome/commit/22388c5fec4607e8e255926c3e95724a2f070e76  

 

So it doesn't check $p \equiv 1 \pmod{q}$ anymore! This is good stuff. However, one thing still bugged me. 

After sending $h$, the server computes $g \equiv h^{\lfloor (p-1)/q \rfloor} \pmod{p}$. It then checks

• $g \not\equiv 1 \pmod{p}$
• $g^q \equiv 1 \pmod{p}$

If $p$ was a prime, then such $g$ can only exist if $p \equiv 1 \pmod{q}$. This forces $p$ to not be prime.

 

However, the pycryptodome library mentions that it checks for the primality of $p$, and there were no patches on that.

 

So I looked at the primality test function used in the repository. It consisted of a few Miller-Rabin tests and one Lucas test. 

If there was no Lucas test, it was not very hard to bypass this with some very large semiprimes. Because the number of Miller-Rabin iterations in the repository did not consider adversarial attacks, if we send a very large "carmichael like" semiprime, then we would only do one round of Miller-Rabin, and have good probability of passing the Miller-Rabin part of the primality test. Of course, the downside is 

• We actually need $p$ to be exactly 2048 bits, but I didn't know that at the time 
• We also need to pass Lucas test, and adversarial examples passing both Miller-Rabin and Lucas is not known, I think?

 

After deciding that finding a composite $p$ that passes the primality check is as hard as writing a conference paper, I looked back.

 

fmt_error = test_probable_prime(p) == COMPOSITE
fmt_error = test_probable_prime(q) == COMPOSITE

 

That second line had to be OR'ed, not substituted, yet it was substituted. This was also in the original repository.

This meant that the primality check of $p$ is completely ignored, which means I didn't need to think about Miller-Rabin and stuff.

 

Step 3 : the finish

Since $p$ doesn't have to be prime, we will use the variable $n$ to avoid confusion.

OK, so I want to solve a discrete logarithm in a subgroup of $\mathbb{Z}_{n}^\star$ with a group order of $q$. 

This is a classical one - we can always take $n$ to be a power of $q$ and use Binomial Theorem. To be exact, we take $n = q^8$ and $h = 1 + q^7$.

 

Now $$g \equiv h^{q^7 - 1} \equiv 1 - q^7 \pmod{q^8}$$ and we see that $$y \equiv g^x \equiv 1 - x q^7 \pmod{q^8}$$ which means we can recover $0 \le x < q$ easily from $y$. Check out Paillier Cryptosystem.

  

Since $n$ needs to be exactly $2048$ bits, $n = q^8$ may fail. You can either reconnect until successful, or try to multiply $2$ until $n$ is exactly $2048$ bits. In this case, you also need to patch $h$ a bit as well. This is left as an exercise.

 

This was an astounding problem (thanks to the author!), as one of the main vulns was in the original repository as well. 

I briefly thought about whether this unpatched vuln is enough to create some issues, but I don't think so. 

 

 
 
r = remote('34.146.212.53', 61234)
 
s = r.recvline()
q = int(s.split()[-1])
 
p = q ** 8
while p.bit_length() < 2048:
    p = 2 * p 
 
h = 1 + 16 * q ** 7
r.sendline(str(p))
r.sendline(str(h))
 
g = int(r.recvline().split()[-1])
y = int(r.recvline().split()[-1])
 
print(2 <= g < p)
print(pow(g, q, p) == 1)
 
gs = ((g - 1) // (q ** 7)) % q
ys = ((y - 1) // (q ** 7)) % q
 
x = (ys * inverse(gs, q)) % q 
 
res = bytes_to_long(hashlib.sha256(b'flag').digest())
 
k = 1
rr = g % q
ss = (res + x * rr) % q
 
print(r.recvline())
 
 
res = long_to_bytes(rr, 32) + long_to_bytes(ss, 32)
 
r.sendline(b64encode(res))
 
print(r.recvline())
print(r.recvline())

 

Flag is Win

require 'openssl'
require 'digest'
 
STDOUT.sync = true
 
class OpenSSL::PKey::EC::Point
  def xy
    n = to_bn(:uncompressed).to_i
    mask = (1 << group.degree) - 1
    return (n >> group.degree) & mask, n & mask
  end
  alias_method :+, :add
  alias_method :*, :mul
end
 
class ECDSA
  def initialize
    @curve = OpenSSL::PKey::EC::Group.new('secp256k1')
    @G = @curve.generator
    @n = @curve.order.to_i
    @d = OpenSSL::BN.rand(@curve.degree).to_i
    @Q = @G * @d
  end
 
  def inv(x)
    x.pow(@n - 2, @n)
  end
 
  def sign(msg)
    z = Digest::SHA256.hexdigest(msg).hex
    k = OpenSSL::BN.rand(@curve.degree / 3).to_s.unpack1('H*').hex
    x, y = (@G * k).xy
 
    # We should discourage every evil hacks
    s = (z + x * @d) * inv(k) % @n
 
    return x, s
  end
 
  def verify(msg, x, s)
    return false if x % @n == 0 || s % @n == 0
    z = Digest::SHA256.hexdigest(msg).hex
 
    # ditto
    x2, y2 = (@G * (z * inv(s)) + @Q * (x * inv(s))).xy
 
    return x == x2
  end
end
 
ecdsa = ECDSA.new
 
5.times do
  puts <<~EOS
    1. Sign
    2. Find rule
    3. Exit
  EOS
 
  print 'choice? '
 
  case gets.chomp
  when '1'
    x, s = ecdsa.sign('Baba')
    puts 'Baba is:'
    puts "x = #{x}"
    puts "s = #{s}"
  when '2'
    print 'Which rule do you want to know? '; msg = gets.chomp
    print 'x? '; x = gets.to_i
    print 's? '; s = gets.to_i
 
    if ecdsa.verify(msg, x, s)
      if msg == 'Baba'
        puts 'Baba is you'
      elsif msg == 'Flag'
        puts "Flag is #{ENV['FLAG']}"
      else
        puts 'Not Found :('
      end
    else
      puts 'Invalid :('
    end
  else
    exit
  end
end
 
puts 'You is defeat.'
 
 

 

This challenge also took me ridiculously long because I made many mistakes and my intuition on lattices is not solid. 

 

It took me very long to realize that I have ruby installed on WSL and I could test what that whole unpack hex thing is. 

Of course, experienced CTF players may notice that unpack hex thing from last year's SECCON, but I didn't solve that challenge :P

 

Anyways, if you test out that unpack hex thing, we can see that $k$ has the form of $$ 48 \cdot \sum_{m=0}^{26} 256^m + \sum_{m=0}^{26} v_m \cdot 256^m$$ where $0 \le v_m \le 9$. We also know that $$k_1 s_1 \equiv z + x_1 d \pmod{n}, \quad k_2 s_2 \equiv z + x_2 d \pmod{n}$$ which, after canceling $d$ out, gives $$k_1(s_1x_2) - k_2(s_2x_1) \equiv z(x_2-x_1) \pmod{n}$$ This can be written as a linear equation of $26 \times 2$ variables between $0$ and $9$ inclusive, and we can shove it into CVP repository.

It seems like you need BKZ instead of LLL to find the correct values, which is understandable since BKZ is very strong.

 

I took a lot of time trying to use as many signatures as possible, leading to very large matrix size and longer runtime. 

I also tried a lot of various hacks which worked very well for ACSC Share the Flag, but they didn't work here. Lattices are hard...

 

# Directly taken from rbtree's LLL repository
# From https://oddcoder.com/LOL-34c3/, https://hackmd.io/@hakatashi/B1OM7HFVI
def Babai_CVP(mat, target):
    M = mat.BKZ(block_size = 35)
    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)
 
    print(result[num_ineq - 1] - target[num_ineq-1])
 
    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
 
 
 
 
r = remote('34.146.212.53', 35719)
 
p = (1 << 256) - (1 << 32) - (1 << 9) - (1 << 8) - (1 << 7) - (1 << 6) - (1 << 4) - 1
 
E = EllipticCurve(GF(p), [0, 7])
Gx = 0x79BE667EF9DCBBAC55A06295CE870B07029BFCDB2DCE28D959F2815B16F81798
Gy = 0x483ADA7726A3C4655DA4FBFC0E1108A8FD17B448A68554199C47D08FFB10D4B8
 
G = E(Gx, Gy)
n = E.order()
print(isPrime(n))
 
h1 = bytes_to_long(hashlib.sha256(b'Baba').digest())
h2 = bytes_to_long(hashlib.sha256(b'Flag').digest())
 
X = []
S = []
for _ in range(4):
    for i in range(3):
        r.recvline()
    r.sendline(b"1")
    r.recvline()
    X.append(int(r.recvline().split()[-1]))
    S.append(int(r.recvline().split()[-1]))
 
NUM_EQ = 4
test = False
 
D = 26
 
supp = []
if test:
    d = rand.randint(1, n)
    for i in range(NUM_EQ):
        cc = []
        k = 0
        for j in range(2 * D):
            if j % 2 == 0:
                u = rand.randint(0, 9)
                supp.append(u)
                k += u * (16 ** j)
                cc.append(u)
            else:
                k += 3 * (16 ** j)
        x = int((k * G).xy()[0])
        s = ((h1 + x * d) * inverse(k, n)) % n 
        X[i] = x
        S[i] = s 
    supp.append(d)
 
print(supp)
M = Matrix(ZZ, 2 * D + 1, 2 * D + 1)
lb = [0] * (2 * D + 1)
ub = [0] * (2 * D + 1) 
 
base_k = 0
for i in range(D):
    base_k += 3 * 16 * (256 ** i)
 
for i in range(2 * D):
    M[i, i] = 1
    lb[i] = 0
    ub[i] = 16 
 
for i in range(D):
    M[i, 2 * D] = int(((256 ** i) * (S[0] * X[1])) % n)
    M[i + D, 2 * D] = int(n - ((256 ** i) * (S[1] * X[0])) % n) 
    M[2 * D, 2 * D] = int(n)
    lb[2 * D] = int((h1 * (X[1] - X[0]) - base_k * S[0] * X[1] + base_k * S[1] * X[0]) % n)
    ub[2 * D] = int((h1 * (X[1] - X[0]) - base_k * S[0] * X[1] + base_k * S[1] * X[0]) % n)
 
 
result, applied_weights, fin = solve(M, lb, ub)
print(fin)
 
k0 = base_k 
for i in range(26):
    k0 += (256 ** i) * int(fin[i]) 
 
d = (inverse(X[0], n) * (k0 * S[0] - h1)) % n 
 
x = Gx 
s = (h2 + x * d) % n 
 
for i in range(3):
    print(r.recvline())
r.sendline(b"2")
r.sendline(b"Flag")
r.sendline(str(x))
r.sendline(str(s))
print(r.recvline())
print(r.recvline())
print(r.recvline())

 

Lumberjack in Nature

from mpmath import mp, power, ln
import json
 
mp.dps = 1000000000
 
def decode(enc):
    return int(power(2, enc * ln(2)))
 
s, e = json.load(open('encoded.json'))
flag = decode(s << e)
 
print(flag.to_bytes((flag.bit_length() + 7) // 8, 'big')[:74])

 

To solve this problem, we need to know the higher bits of $$2^{s \cdot 2^e \cdot \ln 2}$$ where $e = 13371337$ is very large. This is equivalent to finding the decimal part of $s \cdot 2^e \cdot \ln 2$ to a very high precision. 

 

Since you need the decimal part, and $2^e$ is very large, if we want to do direct computation we would need at least $10^7$ binary digits of precision, which seems like too much to handle, even for SageMath. We would like the computation to be easier to do.

 

UPDATE : Never underestimate SageMath! Using $2 \cdot 10^7$ binary digits of precision works very well and fast.

 

The key idea is to approximate $\ln 2$ using the Taylor series $$\ln 2 = \sum_{n=1}^\infty \frac{1}{2^n n}$$ This implies that $$s \cdot 2^e \cdot \ln 2 = \sum_{n=1}^\infty \frac{s 2^{e-n}}{n}$$ and we can compute the decimal part of this as follows. We will sum from $n=1$ to $n=14000000$ as it is enough for precision.

• If $e > n$, then compute $r = s \cdot 2^{e-n} \pmod{n}$ and add $r/n$ to the sum
• If $n \le e \le n+600$, then compute $r = s \pmod{n \cdot 2^{n-e}}$ and add $r / (n \cdot 2^{n-e})$ to the sum
• If $e > n+600$, then add $s / (n \cdot 2^{n-e})$ to the sum 
• After each addition, if the value is larger than $1$, subtract $1$ from it

This is enough to compute the decimal part of $s \cdot 2^e \cdot \ln(2)$ with $10^4$ bits of precision in a few minutes. 

 

Now we can compute the higher bits of $2^{s \cdot 2^e \cdot \ln(2)}$ as well, and shift it and make it into bytes to get our flag.

 

R = RealField(10000)
s, e = 1644076701048410800736598044521957621165075009220047353598695908154798545574669213628822983013478131645136107829991107147966592023064802684205984935604556580976432082255148549763, 13371337
print(s.bit_length())
 
res = R(0)
 
for i in tqdm(range(1, 14000000)):
    # s / i* 2^(e-i)
    if i <= e:
        cc = int(  (s * int(pow(2, e - i, i)) ) % i )
        res += R(cc) / R(i)
    elif i <= e + 600:
        cc = s % (i * pow(2, i-e))
        res += R(cc) / R(i * (R(2) ** (i - e)))
    else:
        res += R(s) / R(i * (R(2) ** (i - e)))
    if res >= R(1):
        res -= R(1)
    
print(res)
res = R(2) ** res
 
for i in range(70 * 8, 80 * 8):
    cc = int(res * R(2 ** i))
    print(cc.to_bytes((cc.bit_length() + 7) // 8, 'big'))

 

Lumberjack Against Nature

from mpmath import power, ln
from random import SystemRandom
from string import ascii_letters
from signal import alarm
 
from secret import decode_fast, flag
 
alarm(10)
 
def to_string(number):
    return number.to_bytes((number.bit_length() + 7) // 8, 'big')[:74]
def decode(enc):
    return to_string(int(power(2, enc * ln(2))))
 
assert(decode(1337 << 5) == decode_fast(1337, 5))
 
 
plaintext = ''.join(SystemRandom().choice(ascii_letters) for _ in range(74)).encode()
e = 13371337
 
print(f'decode(s << {e}) == {plaintext}')
s = int(input('s = ? '))
 
if 0 < s < 2 ** (75 * 8) and decode_fast(s, e) == plaintext:
    print(f'Congrats! {flag}')
else:
    print(":P")
 

 

Now we have to go around. Denote the decimal term of $2^{13371337} \cdot \ln(2)$ as $t$, and the target plaintext viewed as a integer as $v$.

 

We want to find an $0 \le s < 2^{600}$ such that $$2^{s \cdot t - z} \approx v$$ for some integer $z$. To solve this, we take the logarithm again and multiply $2^{5000}$, giving us $$s \cdot \lfloor 2^{5000} t \rfloor - 2^{5000} z \approx \lfloor \log_2(v) \cdot 2^{5000} \rfloor$$ It's clear that we can compute the two values $$T = \lfloor 2^{5000} t \rfloor , \quad V = \lfloor \log_2(v) \cdot 2^{5000} \rfloor$$ using the methods described in the challenge above and arbitrary precision logarithms from SageMath. Now we want something like $$ sT \pmod{2^{5000}} \approx V \pmod{2^{5000}}$$ If we plug in the values of the challenge above, we see that $$ V - 2^{4409} \pmod{2^{5000}} \le sT \pmod{2^{5000}} \le V + 2^{4409} \pmod{2^{5000}}$$ Here, note that $$L \pmod{M} \le x \pmod{M} \le R \pmod{M}$$ should be regarded as $x$ lies in the clockwise strip from $L$ to $R$ in a clock divided into $M$ pieces. Check the link below. 

 

Anyways, it's now clear that we want to solve the system $$0 \le s < 2^{600}, \quad V - 2^{4409} \pmod{2^{5000}} \le sT \pmod{2^{5000}} \le V + 2^{4409} \pmod{2^{5000}}$$ which is possible with the "special case variation" of the CVP repository. This will give around $2^{10}$ candidates for $s$. 

 

Since we have already precomputed $t$ and $T$, we can check the validity of each $s$ very easily.

While this solution works with a relatively low probability, it still works well enough to first blood this challenge. :)

 

def ceil(n, m): # returns ceil(n/m)
    return (n + m - 1) // m
 
def is_inside(L, R, M, val): # is L <= val <= R in mod M context?
    if L <= R:
        return L <= val <= R
    else:
        R += M
        if L <= val <= R:
            return True
        if L <= val + M <= R:
            return True 
        return False
 
## some notes : it's good idea to check for gcd(A, M) = 1
## in CTF context, if gcd(A, M) != 1, we can factorize M and sometimes we can solve the challenge
## in competitive programming context, we need to check gcd(A, M) = 1 and decide whether solution even exists..
def optf(A, M, L, R): # minimum nonnegative x s.t. L <= Ax mod M <= R
    if L == 0:
        return 0
    if 2 * A > M:
        L, R = R, L
        A, L, R = M - A, M - L, M - R
    cc_1 = ceil(L, A)
    if A * cc_1 <= R:
        return cc_1
    cc_2 = optf(A - M % A, A, L % A, R % A)
    return ceil(L + M * cc_2, A)
 
# check if L <= Ax (mod M) <= R has a solution
def sol_ex(A, M, L, R):
    if L == 0 or L > R:
        return True
    g = GCD(A, M)
    if (L - 1) // g == R // g:
        return False
    return True
 
## find all solutions for L <= Ax mod M <= R, S <= x <= E:
def solve(A, M, L, R, S, E):
    # this is for estimate only : if very large, might be a bad idea to run this
    # print("Expected Number of Solutions : ", ((E - S + 1) * (R - L + 1)) // M + 1)
    if sol_ex(A, M, L, R) == False:
        return []
    cur = S - 1
    ans = []
    num_sol = 0
    while cur <= E:
        NL = (L - A * (cur + 1)) % M
        NR = (R - A * (cur + 1)) % M
        if NL > NR:
            cur += 1
        else:
            val = optf(A, M, NL, NR)
            cur += 1 + val
        if cur <= E:
            ans.append(cur)
            # remove assert for performance if needed
            assert is_inside(L, R, M, (A * cur) % M)
            num_sol += 1
    print("Actual Number of Solutions : ", num_sol)
    return ans
 
R = RealField(10000)
s, e = 1, 13371337
res = R(0)
 
for i in tqdm(range(1, 14000000)):
    # s / i* 2^(e-i)
    if i <= e:
        cc = int(  (s * int(pow(2, e - i, i)) ) % i )
        res += R(cc) / R(i)
    elif i <= e + 600:
        cc = s % (i * pow(2, i-e))
        res += R(cc) / R(i * (R(2) ** (i - e)))
    else:
        res += R(s) / R(i * (R(2) ** (i - e)))
    if res >= R(1):
        res -= R(1)
 
v = int(res * R(2 ** 5000))
print(v)
 
sys.setrecursionlimit(10 ** 6)
 
while True:
    r = remote('34.146.212.53', 53928)
    s = r.recvline()
    print(s)
    s = s[-76:-2]
    print(s)
 
    cc = bytes_to_long(s)
    res = R(cc).log() / R(2).log()
    res = int(res * R(2 ** 5000))
 
    # enc * v - integer * 2^5000 = ln_2(val) * 2^5000
    # enc * v - integer * 2^5000 = res 
    fin = solve(v, 1 << 5000, (res - (1 << 4409)) % (1 << 5000), (res + (1 << 4409)) % (1 << 5000), 0, 1 << 600)
    dec = R(v) / R(2 ** 5000)
 
    finished = False
    for cand in fin:
        if finished:
            break
        val = dec * R(cand)
        val = val - val.floor()
        val = R(2) ** val
        for i in range(70 * 8, 80 * 8):
            flag = int(val * R(2 ** i))
            flag = flag.to_bytes((flag.bit_length() + 7) // 8, 'big')
            if s == flag[:74]:
                print(s)
                print(cand.bit_length())
                print(flag)
                print(cand)
                r.sendline(str(cand))
                ff = r.recvline()
                if b"? :P" in ff:
                    finished = True
                    break
                else:
                    print(ff)
    r.close()