| Groth16 + Powers of Tau + Tornado Cash (0) | 2022.03.14 |
|---|---|
| Lattice Attacks Introductory Survey (0) | 2021.07.22 |
| Sigma Protocol (0) | 2021.02.15 |
| Schnorr's Protocol & Signature (0) | 2021.02.13 |
| Diophantine Argument of Knowledge (4) | 2021.01.20 |
We played as "rbtree fan club" and placed 3rd overall.
babycrypto 1
|
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
|
#!/usr/bin/env python
from base64 import b64decode
from base64 import b64encode
import socket
import multiprocessing
from Crypto.Cipher import AES
from Crypto.Random import get_random_bytes
from Crypto.Util.Padding import pad, unpad
import hashlib
import sys
class AESCipher:
def __init__(self, key):
self.key = key
def encrypt(self, data):
iv = get_random_bytes(AES.block_size)
self.cipher = AES.new(self.key, AES.MODE_CBC, iv)
return b64encode(iv + self.cipher.encrypt(pad(data,
AES.block_size)))
def encrypt_iv(self, data, iv):
self.cipher = AES.new(self.key, AES.MODE_CBC, iv)
return b64encode(iv + self.cipher.encrypt(pad(data,
AES.block_size)))
def decrypt(self, data):
raw = b64decode(data)
self.cipher = AES.new(self.key, AES.MODE_CBC, raw[:AES.block_size])
return unpad(self.cipher.decrypt(raw[AES.block_size:]), AES.block_size)
flag = open("flag", "rb").read().strip()
COMMAND = [b'test',b'show']
def run_server(client, aes_key, token):
client.send(b'test Command: ' + AESCipher(aes_key).encrypt(token+COMMAND[0]) + b'\n')
client.send(b'**Cipher oracle**\n')
client.send(b'IV...: ')
iv = b64decode(client.recv(1024).decode().strip())
client.send(b'Message...: ')
msg = b64decode(client.recv(1024).decode().strip())
client.send(b'Ciphertext:' + AESCipher(aes_key).encrypt_iv(msg,iv) + b'\n\n')
while(True):
client.send(b'Enter your command: ')
tt = client.recv(1024).strip()
tt2 = AESCipher(aes_key).decrypt(tt)
client.send(tt2 + b'\n')
if tt2 == token+COMMAND[1]:
client.send(b'The flag is: ' + flag)
client.close()
break
if __name__ == '__main__':
server = socket.socket(socket.AF_INET, socket.SOCK_STREAM)
server.setsockopt(socket.SOL_SOCKET, socket.SO_REUSEADDR, 1)
server.bind(('0.0.0.0', 16001))
server.listen(1)
while True:
client, address = server.accept()
aes_key = get_random_bytes(AES.block_size)
token = b64encode(get_random_bytes(AES.block_size*10))[:AES.block_size*10]
process = multiprocessing.Process(target=run_server, args=(client, aes_key, token))
process.daemon = True
process.start()
|
cs |
Here, we are given
and we are supposed to find the AES-CBC encrypted result of the token appended with "show".
We can easily see that we just need to change the last block. Therefore, we only need to change the last block of the ciphertext as well.
We can simply XOR the next to last block of the ciphertext and the padded "show" string, and send it to the encryption oracle with 0 as IV.
The result of this encryption can be used as the last block of our final ciphertext, and this solves the problem. Here's my dirty code.
|
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
|
def bytexor(a, b):
assert len(a) == len(b)
return bytes(x ^ y for x, y in zip(a, b))
r = remote('35.200.115.41', '16001')
def get_data():
t = r.recvline()
print(t)
t = t.rstrip(b"\n")
t = t.split()[-1]
if t[:11] == b"Ciphertext:":
t = t[11:]
print(t)
return t
inp = b64decode(get_data())
IV = inp[:16]
CTXT = inp[16:]
print(CTXT)
print(r.recvline())
print(r.recv(1024))
ss = b64encode(b"\x00" * 16)
r.sendline(ss)
print(r.recv(1024))
val = bytexor(pad(b'show', 16), CTXT[9*16:10*16])
val = b64encode(val)
r.sendline(val)
cc = b64decode(get_data())
cc = cc[16:32]
print(len(cc))
vv = IV + CTXT[:160] + cc
vv = b64encode(vv)
print(r.recv(1024))
r.sendline(vv)
print(r.recv(1024))
|
cs |
babycrypto 2
|
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
|
#!/usr/bin/env python
from base64 import b64decode
from base64 import b64encode
import socket
import multiprocessing
from Crypto.Cipher import AES
from Crypto.Random import get_random_bytes
from Crypto.Util.Padding import pad, unpad
import hashlib
import sys
class AESCipher:
def __init__(self, key):
self.key = key
def encrypt(self, data):
iv = get_random_bytes(AES.block_size)
self.cipher = AES.new(self.key, AES.MODE_CBC, iv)
return b64encode(iv + self.cipher.encrypt(pad(data,
AES.block_size)))
def encrypt_iv(self, data, iv):
self.cipher = AES.new(self.key, AES.MODE_CBC, iv)
return b64encode(iv + self.cipher.encrypt(pad(data,
AES.block_size)))
def decrypt(self, data):
raw = b64decode(data)
self.cipher = AES.new(self.key, AES.MODE_CBC, raw[:AES.block_size])
return unpad(self.cipher.decrypt(raw[AES.block_size:]), AES.block_size)
flag = open("flag", "rb").read().strip()
AES_KEY = get_random_bytes(AES.block_size)
TOKEN = b64encode(get_random_bytes(AES.block_size*10-1))
COMMAND = [b'test',b'show']
PREFIX = b'Command: '
def run_server(client):
client.send(b'test Command: ' + AESCipher(AES_KEY).encrypt(PREFIX+COMMAND[0]+TOKEN) + b'\n')
while(True):
client.send(b'Enter your command: ')
tt = client.recv(1024).strip()
tt2 = AESCipher(AES_KEY).decrypt(tt)
client.send(tt2 + b'\n')
if tt2 == PREFIX+COMMAND[1]+TOKEN:
client.send(b'The flag is: ' + flag)
client.close()
break
if __name__ == '__main__':
server = socket.socket(socket.AF_INET, socket.SOCK_STREAM)
server.setsockopt(socket.SOL_SOCKET, socket.SO_REUSEADDR, 1)
server.bind(('0.0.0.0', 16002))
server.listen(1)
while True:
client, address = server.accept()
process = multiprocessing.Process(target=run_server, args=(client, ))
process.daemon = True
process.start()
|
cs |
This time, we don't have access to any oracle, but we now have to change the first block of the decrypted text.
This can be easily done by simply changing the IV as needed. Here's another dirty code that does the job.
|
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
|
def bytexor(a, b):
assert len(a) == len(b)
return bytes(x ^ y for x, y in zip(a, b))
r = remote('35.200.39.68', '16002')
def get_data():
t = r.recvline()
t = t.rstrip(b"\n")
t = t.split()[-1]
return t
inp = b64decode(get_data())
IV = inp[:16]
CTXT = inp[16:]
print(r.recv(1024))
IV = bytexor(IV, pad(b"Command: test", 16))
IV = bytexor(IV, pad(b"Command: show", 16))
vv = IV + CTXT
vv = b64encode(vv)
r.sendline(vv)
print(r.recv(1024))
print(r.recv(1024))
|
cs |
babycrypto 3
The value of $n, e, c$ are given, and we have to decrypt it.
Using RsaCTFTool, we can factorize $n$ (cm-factor is used) and RSA-decrypt the value of $c$.
The resulting value of $m$, when converted to bytes, looks like a base64 encoded value of a string.
We take some suffix of the results and base64 decode it, to get "CLOSING THE DISTANCE.".
|
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
|
n = 31864103015143373750025799158312253992115354944560440908105912458749205531455987590931871433911971516176954193675507337
e = 65537
c = 10879776433900426660843740332190892429769159527203392037251077478777616065501519198653853699716123394455804888854401574
# python3 RsaCtfTool.py -n 31864103015143373750025799158312253992115354944560440908105912458749205531455987590931871433911971516176954193675507337 -e 65537 --uncipher 10879776433900426660843740332190892429769159527203392037251077478777616065501519198653853699716123394455804888854401574 --private
val = 0x00026067ff851ecdcb61e50b83a515e3005130785055306c4f527942555345556752456c545645464f5130557543673d3d0a
p = 291664785919250248097148750343149685985101
q = 109249057662947381148470526527596255527988598887891132224092529799478353198637
val = long_to_bytes(val)
print(val)
tt = b'Q0xPU0lORyBUSEUgRElTVEFOQ0UuCg==\n' # remove prefix
print(b64decode(tt))
|
cs |
babycrypto 4
We are given a secp160r1 ECDSA signature data, with 32-bit nonces. Also, the first 16 bits of nonces are leaked.
It's trivial from the ECDSA definition that given one "full" data of $r, s, k, h$ we can retrieve the value of the secret key.
Therefore, by brute forcing the remaining unleaked 16 bits, we can create a set of 65536 candidates of the secret key for each dataset.
By taking the intersection of such datasets, we can find the secret key, solving the problem. Here's another code.
|
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
|
r = [0] * 20
s = [0] * 20
k = [0] * 20
h = [0] * 20
f = open('output.txt', 'r')
for i in range(0, 20):
v = f.readline()
cc = v.split()
r[i] = int(cc[0][2:], 16)
s[i] = int(cc[1][2:], 16)
k[i] = int(cc[2][2:], 16)
h[i] = int(cc[3][2:], 16)
n = 0x0100000000000000000001F4C8F927AED3CA752257
s1 = set()
s2 = set()
for i in range(0, 1 << 16):
# sk = h + r * pvk mod n
cc = ((s[0] * (k[0] + i) - h[0]) * inverse(r[0], n)) % n
s1.add(cc)
for i in range(0, 1 << 16):
cc = ((s[1] * (k[1] + i) - h[1]) * inverse(r[1], n)) % n
s2.add(cc)
print(s1 & s2)
|
cs |
| zer0lfsr 2021 Full Writeup (0CTF/TCTF 2021) (2) | 2021.07.05 |
|---|---|
| zer0lfsr 2021 Quick Solutions/Thoughts (0) | 2021.07.05 |
| zer0pts CTF 2021 Crypto Writeups (1) | 2021.03.07 |
| AeroCTF 2021 Crypto Writeups (0) | 2021.02.28 |
| UnionCTF Crypto Writeups (2) | 2021.02.21 |
빠르게 Nesterov's Estimate Sequence를 정리한다. 출처는 "Introductory Lectures on Convex Optimization"
정의 : $\{\phi_k(x)\}$와 $\{\lambda_k\}$, $\lambda_k \ge 0$이 함수 $f(x)$의 estimate sequence라는 것은, $\lambda_k \rightarrow 0$이고 임의의 $x \in \mathbb{R}^n$과 $k \ge 0$에 대하여 $$\phi_k(x) \le (1-\lambda_k) f(x) + \lambda_k \phi_0(x)$$라는 것이다. 이 경우, 만약 적당한 수열 $\{x_k\}$가 있어 $$f(x_k) \le \phi_k^{*} \equiv \min_{x \in \mathbb{R}^n} \phi_k(x)$$가 성립한다면, $$f(x_k) \le \phi_k(x_{*}) \le (1-\lambda_k) f_{*} + \lambda_k \phi_0(x_{*})$$가 되어 이를 정리하면 결과적으로 $$f(x_k) - f_{*} \le \lambda_k (\phi_0(x_{*}) - f_{*})$$를 얻는다. 그러니 $x_k$를 잘 잡고, $\lambda_k$가 빠르게 $0$으로 수렴하기만 하면 된다.
이제 estimate sequence를 하나 만들어보자. $f \in \mathcal{F}_{\mu, L}$이라 가정한다. 즉, $f$는 $\mu$-strongly convex, $L$-smooth function.
$\phi_0$를 임의의 함수, $\{y_k\}$를 임의의 $\mathbb{R}^n$ 위 수열, $\{\alpha_k\}$를 $\alpha_k \in (0, 1)$과 $\sum_{k} \alpha_k \rightarrow \infty$를 만족하는 수열이라 하자. $\lambda_0 = 1$이라 하자.
이제 $\{\phi_k\}$와 $\{\lambda_k\}$를 다음과 같이 재귀적으로 정의한다. $$\begin{equation*} \begin{aligned} \lambda_{k+1} &= (1 - \alpha_k) \lambda_k \\ \phi_{k+1}(x) & = (1-\alpha_k)\phi_k(x) + \alpha_k(f(y_k) + \langle \nabla f(y_k), x - y_k \rangle + \frac{\mu}{2} ||x-y_k||^2) \end{aligned} \end{equation*}$$
이는 estimate sequence가 된다. 증명은 $\phi_{k+1}(x) \le (1-\alpha_k) \phi_k(x) + \alpha_k f(x)$임을 활용하여 귀납법. 특히, $\sum_k \alpha_k \rightarrow \infty$ 조건에서 $\lambda_k \rightarrow 0$이 나온다. 이렇게 되면 estimate sequence의 큰 family를 하나 얻었으니, 여기서 입맛에 맞게 $\lambda_k$가 빠르게 감소하도록 $\{\alpha_k\}$와 $\{y_k\}$를 잘 잡아주면 되겠다.
$\phi_0$도 잡아줘야 하는데, 뒤에 계속해서 quadratic function과 섞어먹고 있으니, $\phi_0$도 quadratic하게 잡자.
$\phi_0 = \phi_0^{*} + \frac{\gamma_0}{2} ||x-v_0||^2$이라고 정의하면, $\phi_k$가 전부 quadratic function이 된다. 귀찮은 계산을 하면, $$ \begin{equation*} \begin{aligned} \gamma_{k+1} &= (1-\alpha_k) \gamma_k + \alpha_k \mu \\ v_{k+1} &= \frac{1}{\gamma_{k+1}} \left( (1-\alpha_k) \gamma_k v_k + \alpha_k \mu y_k - \alpha_k \nabla f(y_k) \right) \\ \phi_{k+1}^{*} &= (1-\alpha_k) \phi_k^{*} + \alpha_k f(y_k) - \frac{\alpha^2_k}{2 \gamma_{k+1}} ||\nabla f(y_k)||^2 \\ & + \frac{\alpha_k(1-\alpha_k) \gamma_k}{\gamma_{k+1}} \left( \frac{\mu}{2} ||y_k - v_k||^2 + \langle \nabla f(y_k), v_k - y_k \rangle \right) \end{aligned} \end{equation*}$$라 하면 $\phi_k(x) = \phi_k^{*} + \frac{\gamma_k}{2} ||x - v_k||^2$이 된다.
이 framework에 속하는 알고리즘에서 계산하기 가장 어렵게 생긴 것은 $x_k$를 설정하는 단계이다. 첫 $x_0$야 $\lambda_0 = 1$ 조건에서 아무렇게나 잡아도 충분하니, $x_k$에서 $x_{k+1}$로 쉽게 넘어가는 방법이 있으면 좋을 것 같다. $\phi_k^{*} \ge f(x_k)$를 가정하자. 위 식과 $f$의 convexity를 활용하면, $$\begin{equation*} \begin{aligned} \phi_{k+1}^{*} & \ge (1-\alpha_k) f(x_k) + \alpha_k f(y_k) - \frac{\alpha_k^2}{2\gamma_{k+1}} ||\nabla f(y_k)||^2 + \frac{\alpha_k (1-\alpha_k) \gamma_k}{\gamma_{k+1}} \langle \nabla f(y_k), v_k - y_k \rangle \\ & \ge f(y_k) - \frac{\alpha_k^2}{2\gamma_{k+1}} ||\nabla f(y_k)||^2 + (1-\alpha_k) \langle \nabla f(y_k), x_k - y_k + \frac{\alpha_k \gamma_k}{\gamma_{k+1}} (v_k - y_k) \rangle \end{aligned} \end{equation*}$$ 그러니 우리는 $$f(y_k) - \frac{\alpha_k^2}{2\gamma_{k+1}} ||\nabla f(y_k)||^2 + (1-\alpha_k) \langle \nabla f(y_k), x_k - y_k + \frac{\alpha_k \gamma_k}{\gamma_{k+1}} (v_k - y_k) \ge f(x_{k+1})$$인 $x_{k+1}$이 잡기 쉬운 값이었으면 좋겠다. 그런데 우리는 이미 강력한 Descent Lemma, 즉 $$f \left(x - \frac{1}{L} \nabla f(x) \right) \le f(x) - \frac{1}{2L} || \nabla f(x)||^2$$이라는 식을 알고 있으니, 이를 위 부등식에 끼워맞추면 좋겠다는 생각을 할 수 있다. 이러면 필요한 식이 $$\frac{\alpha_k^2}{2 \gamma_{k+1}} = \frac{1}{2L}, \quad x_k - y_k + \frac{\alpha_k \gamma_k}{\gamma_{k+1}} (v_k - y_k) = 0, \quad x_{k+1} = y_k - \frac{1}{L} \nabla f(y_k)$$
이를 모두 순서에 맞게 잘 합치면, 다음과 같은 iteration을 얻는다. $\gamma_0 > 0$과 $x_0$를 임의로 잡고, $v_0 = x_0$라 하자.
$$ \begin{equation*} \begin{aligned} \alpha_k \in (0, 1) &: L\alpha_k^2 = (1-\alpha_k) \gamma_k + \alpha_k \mu \\ \gamma_{k+1} &= (1-\alpha_k)\gamma_k + \alpha_k \mu \\ y_k &= \frac{\alpha_k \gamma_k v_k + \gamma_{k+1} x_k}{\gamma_k + \alpha_k \mu} \\ x_{k+1} &= y_k - \frac{1}{L} \nabla f(y_k) \\ v_{k+1} &= \frac{1}{\gamma_{k+1}} \left( (1-\alpha_k)y_k v_k + \alpha_k \mu y_k - \alpha_k \nabla f(y_k) \right) \end{aligned} \end{equation*}$$를 반복한다. 여기서 $\gamma_0 = L$을 주고 $\lambda_k$의 값에 대한 bound를 열심히 계산하면, 이 방법이 acceleration을 보임을 확인할 수 있다.
이제 복잡한 식을 간단하게 하기 위해서 열심히 식 조작을 하면
$$ \begin{equation*} \begin{aligned} v_{k+1} &= x_k + \frac{1}{\alpha_k} (x_{k+1} - x_k) \\ \alpha_{k+1}^2 &= (1-\alpha_{k+1}) \alpha_k^2 + \alpha_{k+1} \cdot \mu /L \\ \beta_k &= \frac{\alpha_k (1-\alpha_k)}{\alpha_k^2 + \alpha_{k+1}} \\ y_{k+1} &= x_{k+1} + \beta_k (x_{k+1} - x_k) \end{aligned} \end{equation*} $$을 얻는다. 결국 모든 식이 다 사라지고 남는 것은 $\alpha_k, x_k, y_k$ 뿐이다. 결과는 $$\begin{equation*} \begin{aligned} x_{k+1} &= y_k - \frac{1}{L} \nabla f(y_k) \\ \alpha_{k+1}^2 &= (1-\alpha_{k+1})\alpha_k^2 + \alpha_{k+1} \cdot \mu / L \\ \beta_k &= \frac{\alpha_k (1-\alpha_k)}{\alpha_k^2 + \alpha_{k+1}} \\ y_{k+1} &= x_{k+1} + \beta_k (x_{k+1} - x_k) \end{aligned} \end{equation*} $$ 특히 $\alpha_0 \ge \sqrt{\mu/L}$이면 acceleration이 일어난다. 또한, $\alpha_0 = \sqrt{\mu/L}$이면 $\alpha, \beta$의 값이 고정되어 fixed momentum인 경우가 나타난다.
물론 아예 $\alpha_0 = \sqrt{\mu / L}$로 잡아버리면 $\mu = 0$인 경우에 위험하므로 주의하자. 이렇게 해서 유도 과정이 끝난다.
다음에 기회가 되면 여기서 아이디어를 일부 얻어간 Catalyst에 대해 알아보고자 한다.
| 4-5월 정리 (1) | 2021.05.31 |
|---|---|
| 3월 정리 (0) | 2021.04.02 |
| 2월 정리 (0) | 2021.02.28 |
| 1월 정산 + 2월 계획 + 1학기 수강신청 (1) | 2021.02.01 |
| 2020 : A Retrospective (9) | 2020.12.27 |
I participated in zer0pts CTF in Super Guesser, and we reached 1st place.
Since there are many challenges to cover, this writeup will be very concise.
The solution codes may be a bit messy, since I wrote this in a hurry :P
Also, r4j solved signme, so I'll leave that writeup to r4j :D
|
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
|
from Crypto.Util.number import getStrongPrime, GCD
from random import randint
from flag import flag
import os
def pad(m: int, n: int):
# PKCS#1 v1.5 maybe
ms = m.to_bytes((m.bit_length() + 7) // 8, "big")
ns = n.to_bytes((n.bit_length() + 7) // 8, "big")
assert len(ms) <= len(ns) - 11
ps = b""
while len(ps) < len(ns) - len(ms) - 3:
p = os.urandom(1)
if p != b"\x00":
ps += p
return int.from_bytes(b"\x00\x02" + ps + b"\x00" + ms, "big")
while True:
p = getStrongPrime(512)
q = getStrongPrime(512)
n = p * q
phi = (p-1)*(q-1)
e = 1337
if GCD(phi, e) == 1:
break
m = pad(int.from_bytes(flag, "big"), n)
c1 = pow(m, e, n)
c2 = pow(m // 2, e, n)
print("n =", n)
print("e =", e)
print("c1=", c1)
print("c2=", c2)
|
cs |
Writing $\lfloor m/2 \rfloor = x$, we see that $m$ is either $2x$ or $2x+1$. If $m = 2x$, we can perform the related message attack with polynomials $$p_1(x) = (2x)^e - c_1, \quad p_2(x) = x^e - c_2$$ Of course, in this case we would have $c_1 \equiv 2^e c_2 \pmod{n}$ and polynomial GCD would be useless. If $m = 2x+1$, we can do it with $$p_1(x) = (2x+1)^e - c_1, \quad p_2(x) = x^e - c_2$$ which works and gives our desired value of $x$. $m$ can be recovered as $2x+1 \pmod{n}$.
|
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
|
# sage
n = 113135121314210337963205879392132245927891839184264376753001919135175107917692925687745642532400388405294058068119159052072165971868084999879938794441059047830758789602416617241611903275905693635535414333219575299357763227902178212895661490423647330568988131820052060534245914478223222846644042189866538583089
e = 1337
c1 = 89077537464844217317838714274752275745737299140754457809311043026310485657525465380612019060271624958745477080123105341040804682893638929826256518881725504468857309066477953222053834586118046524148078925441309323863670353080908506037906892365564379678072687516738199061826782744188465569562164042809701387515
c2 = 18316499600532548540200088385321489533551929653850367414045951501351666430044325649693237350325761799191454032916563398349042002392547617043109953849020374952672554986583214658990393359680155263435896743098100256476711085394564818470798155739552647869415576747325109152123993105242982918456613831667423815762
P.<x> = PolynomialRing(Zmod(n))
def GCD(f, g, n):
g = g % f
if g == 0:
return f
t = g.lc()
if gcd(t, n) != 1:
print(t)
exit()
tt = inverse_mod(Integer(t), n)
g = g * tt
return GCD(g, f, n)
f = (2 * x + 1)^e - c1
g = x^e - c2
print(GCD(f, g, n))
# back to python
m = -113131683468284119607964196541575421728552328779560800910707074969146190674728631210678514366174908879152073378298687396133007061413453069629241334690374397179449596372765780570923850759894857943624531609494119853246060991087070835093327088258517606688574536921269027653158032526646833175962482267255713310835
n = 113135121314210337963205879392132245927891839184264376753001919135175107917692925687745642532400388405294058068119159052072165971868084999879938794441059047830758789602416617241611903275905693635535414333219575299357763227902178212895661490423647330568988131820052060534245914478223222846644042189866538583089
e = 1337
c1= 89077537464844217317838714274752275745737299140754457809311043026310485657525465380612019060271624958745477080123105341040804682893638929826256518881725504468857309066477953222053834586118046524148078925441309323863670353080908506037906892365564379678072687516738199061826782744188465569562164042809701387515
c2= 18316499600532548540200088385321489533551929653850367414045951501351666430044325649693237350325761799191454032916563398349042002392547617043109953849020374952672554986583214658990393359680155263435896743098100256476711085394564818470798155739552647869415576747325109152123993105242982918456613831667423815762
print(n+m)
fin = (2 * (n+m) + 1) % n
print(long_to_bytes(fin))
|
cs |
|
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
|
import os
import signal
import random
from base64 import b64encode
from Crypto.Util.number import getStrongPrime, bytes_to_long
from Crypto.Util.Padding import pad
from Crypto.Cipher import AES
from flag import flag
p = getStrongPrime(1024)
key = os.urandom(32)
iv = os.urandom(AES.block_size)
aes = AES.new(key=key, mode=AES.MODE_CBC, iv=iv)
c = aes.encrypt(pad(flag, AES.block_size))
key = bytes_to_long(key)
print("Encrypted flag: {}".format(b64encode(iv + c).decode()))
print("p = {}".format(p))
print("key.bit_length() = {}".format(key.bit_length()))
signal.alarm(600)
while key > 0:
r = random.randint(2, p-1)
s = random.randint(2, p-1)
t = random.randint(2, p-1)
print("t = {}".format(t))
a = int(input("a = ")) % p
b = int(input("b = ")) % p
c = int(input("c = ")) % p
d = int(input("d = ")) % p
assert all([a > 1 , b > 1 , c > 1 , d > 1])
assert len(set([a,b,c,d])) == 4
u = pow(a, r, p) * pow(c, s, p) % p
v = pow(b, r, p) * pow(c, s, p) % p
x = u ^ (key & 1)
y = v ^ ((key >> 1) & 1)
z = pow(d, r, p) * pow(t, s, p) % p
key = key >> 2
print("x = {}".format(x))
print("y = {}".format(y))
print("z = {}".format(z))
|
cs |
The key idea is to make a small subgroup which the values of $x, y \in \mathbb{F}_p$ must belong to.
To do this, we wait until $p \equiv 1 \pmod{4}$, then take any $t \in \mathbb{F}_p$ with order $4$.
This can be done by taking some random $g \in \mathbb{F}_p$ and choosing $t = g^{(p-1)/4}$.
Checking if order of $t$ is $4$ can now be done simply by verifying $t^2 \neq 1$ in $\mathbb{F}_p$.
Now we simply send $a = t$, $b = t^2$, $c = t^3$. Since $a^4 = b^4 = c^4 = 1$, we also have $u^4 = v^4 = 1$.
Therefore, we can see if $x = u$ and if $y= v$ by simply checking if $x^4=1$ and if $y^4=1$.
Of course, it should be noted that getStrongPrime does not generate a prime of the form $p=2q+1$ with $q$ prime. :wink:
|
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
|
r = remote('crypto.ctf.zer0pts.com', '10130')
def recvint():
s = r.recvline()[:-1]
s = s.split()[-1]
return int(s)
S = r.recvline()[:-1]
S = S.split()[-1]
V = base64.b64decode(S)
iv = V[:16]
ctxt = V[16:]
p = recvint()
keylen = recvint()
if p % 4 != 1:
print("BAD PRIME")
exit()
t = 0
for g in range(2, 2000):
t = pow(g, (p-1)//4, p)
if (t * t) % p != 1:
break
print("Begin Key Finding")
keyv = 0
add = 1
for i in tqdm(range(0, 128)):
t = recvint()
r.sendline(str(t))
r.sendline(str((t * t) % p))
r.sendline(str((t * t * t) % p))
r.sendline(str(5))
x = recvint()
y = recvint()
z = recvint()
if pow(x, 4, p) != 1:
keyv += add
add = add * 2
if pow(y, 4, p) != 1:
keyv += add
add = add * 2
if i >= 126:
keyv = long_to_bytes(keyv)
aes = AES.new(key=keyv, mode=AES.MODE_CBC, iv=iv)
ptxt = aes.decrypt(ctxt)
print(ptxt)
|
cs |
|
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
|
import random
import signal
from flag import flag
from Crypto.Util.number import getStrongPrime, inverse
HANDNAMES = {
1: "Rock",
2: "Scissors",
3: "Paper"
}
def commit(m, key):
(g, p), (x, _) = key
r = random.randint(2, p-1)
c1 = pow(g, r, p)
c2 = m * pow(g, r*x, p) % p
return (c1, c2)
def decrypt(c, key):
c1, c2 = c
_, (x, p)= key
m = c2 * inverse(pow(c1, x, p), p) % p
return m
def keygen(size):
p = getStrongPrime(size)
g = random.randint(2, p-1)
x = random.randint(2, p-1)
return (g, p), (x, p)
signal.alarm(3600)
key = keygen(1024)
(g, p), _ = key
print("[yoshiking]: Hello! Let's play Janken(RPS)")
print("[yoshiking]: Here is g: {}, and p: {}".format(g, p))
round = 0
wins = 0
while True:
round += 1
print("[system]: ROUND {}".format(round))
yoshiking_hand = random.randint(1, 3)
c = commit(yoshiking_hand, key)
print("[yoshiking]: my commitment is={}".format(c))
hand = input("[system]: your hand(1-3): ")
print("")
try:
hand = int(hand)
if not (1 <= hand <= 3):
raise ValueError()
except ValueError:
print("[yoshiking]: Ohhhhhhhhhhhhhhhh no! :(")
exit()
yoshiking_hand = decrypt(c, key)
print("[yoshiking]: My hand is ... {}".format(HANDNAMES[yoshiking_hand]))
print("[yoshiking]: Your hand is ... {}".format(HANDNAMES[hand]))
result = (yoshiking_hand - hand + 3) % 3
if result == 0:
print("[yoshiking]: Draw, draw, draw!!!")
elif result == 1:
print("[yoshiking]: Yo! You win!!! Ho!")
wins += 1
print("[system]: wins: {}".format(wins))
if wins >= 100:
break
elif result == 2:
print("[yoshiking]: Ahahahaha! I'm the winnnnnnner!!!!")
print("[yoshiking]: You, good loser!")
print("[system]: you can check that yoshiking doesn't cheat")
print("[system]: here's the private key: {}".format(key[1][0]))
exit()
print("[yoshiking]: Wow! You are the king of roshambo!")
print("[yoshiking]: suge- flag ageru")
print(flag)
|
cs |
The first observation is that if $g$ is a quadratic residue (QR), the commitment leaks the legendre symbol of $m$.
Indeed, if $g$ is QR, then the legendre symbol of $m$ is simply the legendre symbol of $c_2$.
We generate instances until we find a prime $p$ that not all $1, 2, 3$ are QR modulo $p$.
In this case, by leaking the legendre symbol of $m$, we can reduce the number of candidates for yoshiking's hand from 3 to 2.
Therefore, we can at least force a tie, and occasionally win. This is enough to solve the challenge.
|
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
|
r = remote('crypto.ctf.zer0pts.com', '10463')
r.recvline()
s = r.recvline()
T = s.split()
g = int(T[4][:-1])
p = int(T[-1].rstrip())
QR = []
NQR = []
if pow(g, (p-1)//2 , p) != 1:
print("BAD PRIME 1")
exit()
QR.append(1)
if pow(2, (p-1)//2, p-1) == 1:
QR.append(2)
else:
NQR.append(2)
if pow(3, (p-1)//2, p-1) == 1:
QR.append(3)
else:
NQR.append(3)
if len(QR) == 3:
print("lol")
exit()
def get_commit():
s = r.recvline()
t = s.split(b'(')[-1]
c1 = int(t.split(b',')[0])
c2 = int(t.split(b',')[1][1:-2])
return c1, c2
for i in range(0, 500):
print(r.recvline())
c1, c2 = get_commit()
if pow(c2, (p-1) // 2, p) == 1:
if len(QR) == 1:
r.sendline(b"3")
else:
if 2 in QR:
r.sendline(b"1")
else:
r.sendline(b"3")
else:
if len(NQR) == 1:
if 2 in NQR:
r.sendline(b"1")
if 3 in NQR:
r.sendline(b"2")
else:
r.sendline(b"2")
r.recvline()
r.recvline()
r.recvline()
s = r.recvline()
if s.split()[-1].rstrip() == b"draw!!!":
continue
else:
s = r.recvline()
print(s)
if s.split()[-1].rstrip() == b"100":
for t in range(0, 4):
print(r.recvline())
|
cs |
|
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
|
from Crypto.Util.number import*
from flag import flag
nbits = 256
p = random_prime(1 << nbits)
Fp = Zmod(p)
P.<v> = PolynomialRing(Fp)
b = randrange(p)
d = 2
F = v^2 + b
v0 = randrange(p)
v1 = F(v0)
k = ceil(nbits * (d / (d + 1)))
w0 = (v0 >> (nbits - k))
w1 = (v1 >> (nbits - k))
# encrypt
m = bytes_to_long(flag)
v = v1
for i in range(5):
v = F(v)
m ^^= int(v)
print(f"p = {p}")
print(f"b = {b}")
print(f"m = {m}")
print(f"w0 = {w0}")
print(f"w1 = {w1}")
|
cs |
It's clear that we just need to find the value of $v_0$. Denote $l = 256$. Using the given info, we can write $$\begin{equation*} \begin{aligned} v_0 &= w_0 \cdot 2^{l - k} + x \\ v_1 &= w_1 \cdot 2^{l-k} + y \\ v_1 & \equiv v_0^2 + b \pmod{n} \end{aligned} \end{equation*}$$ Combining, this gives us the key relation $$ w_1 \cdot 2^{l-k} + y \equiv w_0^2 \cdot 2^{2l-2k} + b + w_0 \cdot 2^{l-k+1} \cdot x + x^2 \pmod{n}$$ The key idea is to remove the quadratic term $x^2$. Denote $$\begin{equation*} \begin{aligned} A &= w_0 \cdot 2^{l-k+1} \\ B &= w_1 \cdot 2^{l-k} - w_0^2 \cdot 2^{2l-2k} - b \end{aligned} \end{equation*}$$ Note that $A, B$ are constants that we know. Our relation is now $$ Ax \equiv B + y - x^2 \pmod{n}$$ Since $x, y$ are all below $2^{l-k}$, we can see that $$ B - 2^{2l-2k} \le Ax \pmod{n} \le B + 2^{l-k}$$ which is a type of relation that can be solved. Check out my writeups for PBCTF Special Gift and SECCON crypto01.
Note that there are $2^{l-k}$ possible values of $x$, and there are approximately $2^{2l-2k}$ possible values of $Ax \pmod{n}$.
Since $l-k + 2l-2k = 3l-3k \approx l$, we can decide that there must be fairly small number of $x$ that satisfy the desired inequality.
This type of analysis was done in both of writeups above, and recently highlighted by Mystiz in AeroCTF writeup.
Now that there are three challenges I solved with this method, I think I'll add this to my CVP repository as well...
|
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
|
# notation a bit different from solution, but same idea
# simplify polynomial in sage
p = 86160765871200393116432211865381287556448879131923154695356172713106176601077
b = 71198163834256441900788553646474983932569411761091772746766420811695841423780
m = 88219145192729480056743197897921789558305761774733086829638493717397473234815
w0 = 401052873479535541023317092941219339820731562526505
w1 = 994046339364774179650447057905749575131331863844814
nbits = 256
P.<x, y> = PolynomialRing(GF(p))
k = (nbits * 2 + 2) // 3
f = w0 * 2^(nbits - k) + x
g = w1 * 2^(nbits - k) + y
h = g - f^2 - b
print(h)
# python
def ceil(n, m):
return (n + m - 1) // m
def optf(A, M, L, R):
if L == 0:
return 0
if 2 * A > M:
L, R = R, L
A = M - A
L = M - L
R = 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)
x = 18867904637006146022735447
y = 19342813113834066795298816
p = 86160765871200393116432211865381287556448879131923154695356172713106176601077
b = 71198163834256441900788553646474983932569411761091772746766420811695841423780
m = 88219145192729480056743197897921789558305761774733086829638493717397473234815
w0 = 401052873479535541023317092941219339820731562526505
w1 = 994046339364774179650447057905749575131331863844814
C1 = 55130802749277213576496911760053178817655787149958046010477129311148596128757
C2 = 78083221913223461198494116323396529665894773452683783127339675579334647310194
nbits = 256
k = (nbits * 2 + 2) // 3
delt = nbits - k
# C1 x + y + C2 - x^2 == 0 mod p
# C1 x == x^2 - y - C2 mod p
L = (0 - (1 << (delt)) - C2 + p) % p
R = ((1 << (2 * delt)) - C2 + p) % p
lst = 0
while lst <= (1 << delt):
NL = (L - C1 * (lst + 1)) % p
NR = (R - C1 * (lst + 1)) % p
if NL > NR:
lst = lst + 1
else:
cc = optf(C1, p, NL, NR)
lst = lst + 1 + cc
mm = m
x = lst
v0 = w0 * (1 << (nbits - k)) + x
v1 = (v0 * v0 + b) % p
v = v1
# try out!
for i in range(5):
v = (v * v + b) % p
mm ^= v
print(long_to_bytes(mm))
|
cs |
|
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
|
from Crypto.Cipher import AES
from Crypto.Random import get_random_bytes
from binascii import hexlify, unhexlify
from hashlib import md5
import os
import signal
from flag import flag
keys = [md5(os.urandom(3)).digest() for _ in range(3)]
def get_ciphers(iv1, iv2):
return [
AES.new(keys[0], mode=AES.MODE_ECB),
AES.new(keys[1], mode=AES.MODE_CBC, iv=iv1),
AES.new(keys[2], mode=AES.MODE_CFB, iv=iv2, segment_size=8*16),
]
def encrypt(m: bytes, iv1: bytes, iv2: bytes) -> bytes:
assert len(m) % 16 == 0
ciphers = get_ciphers(iv1, iv2)
c = m
for cipher in ciphers:
c = cipher.encrypt(c)
return c
def decrypt(c: bytes, iv1: bytes, iv2: bytes) -> bytes:
assert len(c) % 16 == 0
ciphers = get_ciphers(iv1, iv2)
m = c
for cipher in ciphers[::-1]:
m = cipher.decrypt(m)
return m
signal.alarm(3600)
while True:
print("==== MENU ====")
print("1. Encrypt your plaintext")
print("2. Decrypt your ciphertext")
print("3. Get encrypted flag")
choice = int(input("> "))
if choice == 1:
plaintext = unhexlify(input("your plaintext(hex): "))
iv1, iv2 = get_random_bytes(16), get_random_bytes(16)
ciphertext = encrypt(plaintext, iv1, iv2)
ciphertext = b":".join([hexlify(x) for x in [iv1, iv2, ciphertext]]).decode()
print("here's the ciphertext: {}".format(ciphertext))
elif choice == 2:
ciphertext = input("your ciphertext: ")
iv1, iv2, ciphertext = [unhexlify(x) for x in ciphertext.strip().split(":")]
plaintext = decrypt(ciphertext, iv1, iv2)
print("here's the plaintext(hex): {}".format(hexlify(plaintext).decode()))
elif choice == 3:
plaintext = flag
iv1, iv2 = get_random_bytes(16), get_random_bytes(16)
ciphertext = encrypt(plaintext, iv1, iv2)
ciphertext = b":".join([hexlify(x) for x in [iv1, iv2, ciphertext]]).decode()
print("here's the encrypted flag: {}".format(ciphertext))
exit()
else:
exit()
|
cs |
Let's consider plaintexts/ciphertexts of one block first. Denote $E_k$, $D_k$ as ECB encryption/decryption with key $k$.
The encryption process can be written as $$ P \rightarrow E_{k_1}(P) \rightarrow E_{k_2}(E_{k_1}(P) \oplus IV_1) \rightarrow E_{k_2}(E_{k_1}(P) \oplus IV_1) \oplus E_{k_3}(IV_2) = C$$ The key relation is $$E_{k_1}(P) \oplus IV_1 = D_{k_2}(C \oplus E_{k_3}(IV_2))$$ For a fixed value of $C, IV_2$, the right hand side is fixed. Therefore, the left hand side is also fixed.
We start by using the decryption oracle with $C, IV_1, IV_2$ being all 16 bytes of \x00. Denote $P_0$ as the resulting plaintext.
Then, we use the decryption oracle with $C, IV_1, IV_2$ all \x00, but $IV_1$ has the last byte \x01. Denote $P_1$ as the resulting plaintext.
Our observations lead to $$\text{bytes_to_long}(E_{k_1}(P_0) \oplus E_{k_1}(P_1)) = 1$$ This can be tested for all $2^{24}$ possible keys, and with that we find $k_1$.
Now that we know $k_1$, we also know the value of $D_{k_2} \circ E_{k_3}$ applied to 16 bytes of \x00.
This is a perfect scenario for a meet-in-the-middle attack, and with this we can find $k_2$, $k_3$.
retrieve data from server
|
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
|
# get data from server
r = remote('crypto.ctf.zer0pts.com', '10929')
def get_enc(ptxt):
r.recvuntil(b">")
r.sendline(b"1")
ptxt = bytes.hex(ptxt)
r.sendline(ptxt)
s = r.recvline().split()[-1].rstrip()
iv1, iv2, ctxt = s.split(b":")
iv1 = bytes.fromhex(iv1.decode())
iv2 = bytes.fromhex(iv2.decode())
ctxt = bytes.fromhex(ctxt.decode())
return iv1, iv2, ctxt
def get_dec(ctxt, iv1, iv2):
r.recvuntil(b">")
r.sendline(b"2")
ctxt = bytes.hex(ctxt)
iv1 = bytes.hex(iv1)
iv2 = bytes.hex(iv2)
goal = iv1 + ":" + iv2 + ":" + ctxt
r.sendline(goal)
s = r.recvline()
print(s)
s = s.split()[-1].rstrip()
ptxt = bytes.fromhex(s.decode())
return ptxt
def get_flag():
r.recvuntil(b">")
r.sendline(b"3")
s = r.recvline().split()[-1].rstrip()
iv1, iv2, ctxt = s.split(b":")
iv1 = bytes.fromhex(iv1.decode())
iv2 = bytes.fromhex(iv2.decode())
ctxt = bytes.fromhex(ctxt.decode())
assert len(iv1) == 16
assert len(iv2) == 16
assert len(ctxt) == 48
return iv1, iv2, ctxt
f = open("lol.txt", "w")
ptxt_0 = get_dec(b"\x00" * 16, b"\x00" * 16, b"\x00" * 16)
ptxt_1 = get_dec(b"\x00" * 16, b"\x00" * 15 + b"\x01", b"\x00" * 16)
f.write(str(bytes_to_long(ptxt_0)) + "\n")
f.write(str(bytes_to_long(ptxt_1)) + "\n")
iv1, iv2, ctxt = get_flag()
f.write(str(bytes_to_long(iv1)) + "\n")
f.write(str(bytes_to_long(iv2)) + "\n")
f.write(str(bytes_to_long(ctxt)) + "\n")
|
cs |
find the first key
|
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
|
def bytexor(a, b):
assert len(a) == len(b)
return bytes(x ^ y for x, y in zip(a, b))
for i in tqdm(range(0, 1 << 24)):
cc = long_to_bytes(i)
if len(cc) < 3:
cc = b"\x00" * (3 - len(cc)) + cc
k1 = hashlib.md5(cc).digest()
cipher = AES.new(k1, mode = AES.MODE_ECB)
val_1 = cipher.encrypt(ptxt_0)
val_2 = cipher.encrypt(ptxt_1)
det = bytexor(val_1, val_2)
if bytes_to_long(det) == 1:
print(i)
|
cs |
meet in the middle preparation
|
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
|
# we now know the first key
cipher = AES.new(key = key_1, mode = AES.MODE_ECB)
vv = cipher.encrypt(ptxt_0)
f = open("Data1.txt", "w")
for i in tqdm(range(0, 1 << 24)):
cc = long_to_bytes(i)
if len(cc) < 3:
cc = b"\x00" * (3 - len(cc)) + cc
assert bytes_to_long(cc) == i
k2 = hashlib.md5(cc).digest()
cipher = AES.new(key = k2, mode = AES.MODE_ECB)
res = cipher.encrypt(vv)
f.write(str(bytes_to_long(res)) + "\n")
for i in tqdm(range(0, 1 << 24)):
cc = long_to_bytes(i)
if len(cc) < 3:
cc = b"\x00" * (3 - len(cc)) + cc
k3 = hashlib.md5(cc).digest()
cipher = AES.new(key = k3, mode = AES.MODE_ECB)
res = cipher.encrypt(b"\x00" * 16)
f.write(str(bytes_to_long(res)) + "\n")
f.close()
|
cs |
meet in the middle in C++
|
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
|
map<string, int> M;
int main(void)
{
fio; int i, j;
freopen("Data1.txt", "r", stdin);
for(i=0 ; i < (1<<24) ; i++)
{
if(i % 1000000 == 0) cout << i / 1000000 << endl;
string s; cin >> s;
M[s] = i;
}
for(i=0 ; i<(1<<24) ; i++)
{
if(i % 1000000 == 0) cout << i / 1000000 << endl;
string s; cin >> s;
if(M[s] != 0)
{
cout << M[s] << " " << i << endl;
return 0;
}
}
return 0;
}
|
cs |
|
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
|
from Crypto.Util.number import bytes_to_long
from secrets import k, FLAG
p = 13046889097521646369087469608188552207167764240347195472002158820809408567610092324592843361428437763328630003678802379234688335664907752858268976392979073
a = 10043619664651911066883029686766120169131919507076163314397915307085965058341170072938120477911396027902856306859830431800181085603701181775623189478719241
b = 12964455266041997431902182249246681423017590093048617091076729201020090112909200442573801636087298080179764338147888667898243288442212586190171993932442177
E = EllipticCurve(GF(p),[a,b])
P = E(11283606203023552880751516189906896934892241360923251780689387054183187410315259518723242477593131979010442607035913952477781391707487688691661703618439980, 12748862750577419812619234165922125135009793011470953429653398381275403229335519006908182956425430354120606424111151410237675942385465833703061487938776991)
Q = k*P
R = (k+1)*P
p = int(Q[0])
q = int(R[0])
assert is_prime(p)
assert is_prime(q)
e = 0x10001
N = p*q
m = bytes_to_long(FLAG)
c = pow(m,e,N)
print(f'N = {N}')
print(f'c = {c}')
|
cs |
Denote $P = (u, v)$, $Q = (x, y)$. We can utilize the sum of points formula to get the x-coordinate of $R$ as $$ \left( \frac{y-v}{x-u} \right)^2 - (x + u)$$ Therefore, we have the relation $$ g(x, y) = x(y-v)^2 - (x+u)(x-u)^2 - N(x-u)^2 \equiv 0 \pmod{p}$$ We now have to combine this relation with $$y^2 \equiv x^3 + ax + b \pmod{p}$$ Initially I tried to calculate the resultant, but it bugged out (maybe due to incorrect parameters given) so I changed my approach.
If we print out the monomials of $g(x, y)$, we see that $y$ appears in two monomials, $xy^2$ and $xy$.
We can change $y^2$ as $x^3 + ax + b$ in the first one. Now we can drag the $xy$ term into the right hand side, square the equation, and change $y^2$ to $x^3 + ax + b$ to obtain a equation with only one variable $x$. For more details, please refer to the solution code.
This polynomial equation can now be solved, which gives us the candidates for the prime divisor of $N$.
|
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
|
pytp = 13046889097521646369087469608188552207167764240347195472002158820809408567610092324592843361428437763328630003678802379234688335664907752858268976392979073
a = 10043619664651911066883029686766120169131919507076163314397915307085965058341170072938120477911396027902856306859830431800181085603701181775623189478719241
b = 12964455266041997431902182249246681423017590093048617091076729201020090112909200442573801636087298080179764338147888667898243288442212586190171993932442177
E = EllipticCurve(GF(p),[a,b])
G = E(11283606203023552880751516189906896934892241360923251780689387054183187410315259518723242477593131979010442607035913952477781391707487688691661703618439980, 12748862750577419812619234165922125135009793011470953429653398381275403229335519006908182956425430354120606424111151410237675942385465833703061487938776991)
u = G[0]
v = G[1]
N = 22607234899418506929126001268361871457071114354768385952661316782742548112938224795906631400222949082488044126564531809419277303594848211922000498018284382244900831520857366772119155202621331079644609558409672584261968029536525583401488106146231216232578818115404806474812984250682928141729397248414221861387
c = 15850849981973267982600456876579257471708532525108633915715902825196241000151529259632177065183069032967782114646012018721535909022877307131272587379284451827627191021621449090672315265556221217089055578013603281682705976215360078119427612168005716370941190233189775697324558168779779919848728188151630185987
P.<x, y> = PolynomialRing(GF(p))
f = y ** 2 - x ** 3 - a * x - b
g = x * ((y - v) * (y - v) - (x + u) * (x-u) * (x-u)) - N * (x-u) * (x-u)
print(g.monomials())
# [x^4, x^3, xy^2, x^2, xy, x, 1]
T = g.coefficients()
P.<x> = PolynomialRing(Zmod(p))
f = 0
g = 0
f += T[0] * x^4
f += T[1] * x^3
f += T[2] * x * (x^3 + a* x + b)
f += T[3] * x^2
f += T[5] * x
f += T[6]
g += T[4] * T[4] * x^2 * (x^3 + a*x + b)
h = f * f - g
print(h.roots())
N = 22607234899418506929126001268361871457071114354768385952661316782742548112938224795906631400222949082488044126564531809419277303594848211922000498018284382244900831520857366772119155202621331079644609558409672584261968029536525583401488106146231216232578818115404806474812984250682928141729397248414221861387
c = 15850849981973267982600456876579257471708532525108633915715902825196241000151529259632177065183069032967782114646012018721535909022877307131272587379284451827627191021621449090672315265556221217089055578013603281682705976215360078119427612168005716370941190233189775697324558168779779919848728188151630185987
# cand = h.roots()
cand = [(5266647903652352665309561331835186152327627163271331811555419978564191000470060566535428497675116887002541568904535904345037425011015457585262022604897451, 1), (4292528248136861387890911319917455946841411872473250675409509735620572311636407361858881556677385609500178629430025710517411214702704597103005396234440737, 1), (11283606203023552880751516189906896934892241360923251780689387054183187410315259518723242477593131979010442607035913952477781391707487688691661703618439980, 2)]
for p, t in cand:
q = N // p
phi = (p-1) * (q-1)
print(p * q - N)
d = inverse(65537, phi)
print(long_to_bytes(pow(c, d, N)))
|
cs |
|
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
|
from Crypto.Util.number import *
from flag import flag
assert len(flag) == 40
p = 74894047922780452080480621188147614680859459381887703650502711169525598419741
a1 = 22457563127094032648529052905270083323161530718333104214029365341184039143821
a2 = 82792468191695528560800352263039950790995753333968972067250646020461455719312
Fp = Zmod(p^3)
E = EllipticCurve(Fp, [a1, a2])
m = bytes_to_long(flag)
S = E(201395103510950985196528886887600944697931024970644444173327129750000389064102542826357168547230875812115987973230106228243893553395960867041978131850021580112077013996963515239128729448812815223970675917812499157323530103467271226, 217465854493032911836659600850860977113580889059985393999460199722148747745817726547235063418161407320876958474804964632767671151534736727858801825385939645586103320316229199221863893919847277366752070948157424716070737997662741835)
T = m * S
with open("output.txt", "w") as f:
f.write("T: {}".format(T))
|
cs |
The required paper for this challenge is https://arxiv.org/pdf/2010.15543.pdf.
Basically, (like a few challenges this year) we need a homomorphism to a group we can efficiently calculate discrete logarithm on.
It turns out that the given elliptic curve over $\mathbb{Z}_{p^3}$ has a homomorphism to $\mathbb{Z}_{p^2}$.
UPD : It should also be noted that the binary operator in $\mathbb{Z}_{p^2}$ is sum, not product.
To calculate this, we need to know how to perform point multiplication in a projective curve/division free setting.
The formula for this can be found online, for example, here. I wonder if this computation can be done in sage...
|
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
|
N = 74894047922780452080480621188147614680853399736985793708596679454987247185378
# N is the order of the same elliptic curve, except the underlying field is GF(p) instead of ring Z/p^3Z
p = 74894047922780452080480621188147614680859459381887703650502711169525598419741
a = 22457563127094032648529052905270083323161530718333104214029365341184039143821
b = 82792468191695528560800352263039950790995753333968972067250646020461455719312
mod = p * p * p
def point_double(x, y, z):
t = (3 * x * x + a * z * z) % mod
u = (2 * y * z) % mod
v = (2 * u * x * y) % mod
w = (t * t - 2 * v + mod * 3) % mod
x2 = (u * w) % mod
y2 = (t * (v - w) - 2 * (u * u * y * y)) % mod
y2 = (y2 + mod) % mod
z2 = (u * u * u) % mod
return x2, y2, z2
def point_add(x0, y0, z0, x1, y1, z1):
if x0 == 0 and y0 == 0:
return x1, y1, z1
if x1 == 0 and y1 == 0:
return x0, y0, z0
t0 = (y0 * z1) % mod
t1 = (y1 * z0) % mod
t = (t0 - t1 + mod) % mod
u0 = (x0 * z1) % mod
u1 = (x1 * z0) % mod
u = (u0 - u1 + mod) % mod
u2 = (u * u) % mod
v = (z0 * z1) % mod
w = (t * t * v - u2 * (u0 + u1)) % mod
w = (w + mod) % mod
u3 = (u * u2) % mod
x2 = (u * w) % mod
y2 = (t * (u0 * u2 - w) - t0 * u3) % mod
y2 = (y2 + mod) % mod
z2 = (u3 * v) % mod
return x2, y2, z2
def point_mul(x, y, z, n):
xr, yr, zr = 0, 0, 1
while n:
if n % 2 == 1:
xr, yr, zr = point_add(xr, yr, zr, x, y, z)
n -= 1
n = n // 2
x, y, z = point_double(x, y, z)
return xr, yr, zr
x1 = 201395103510950985196528886887600944697931024970644444173327129750000389064102542826357168547230875812115987973230106228243893553395960867041978131850021580112077013996963515239128729448812815223970675917812499157323530103467271226
y1 = 217465854493032911836659600850860977113580889059985393999460199722148747745817726547235063418161407320876958474804964632767671151534736727858801825385939645586103320316229199221863893919847277366752070948157424716070737997662741835
x2 = 49376632602749543055345783411902198690599351794957124343389298933965298693663616388441379424236401744560279599744281594405742089477317921152802669021421009909184865835968088427615238677007575776072993333868804852765473010336459028
y2 = 342987792080103175522504176026047089398654876852013925736156942540831035248585067987750805770826115548602899760190394686399806864247192767745458016875262322097116857255158318478943025083293316585095725393024663165264177646858125759
print((y1 ** 2 - x1 ** 3 - a * x1 - b) % mod)
print((y2 ** 2 - x2 ** 3 - a * x2 - b) % mod)
xf, yf, zf = point_mul(x1, y1, 1, N)
xs, ys, zs = point_mul(x2, y2, 1, N)
A = xs * yf
B = ys * xf
g = GCD(A, B)
A //= g
B //= g
gg = (A * inverse(B, p * p)) % p
for i in range(0, 2):
cc = bytes_to_long(b"zer0pts{")
val = (gg - cc * (1 << 256)) % p + i * p
flag = b"zer0pts{" + long_to_bytes(val)
print(flag)
|
cs |
|
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
|
import base64
import random
from Crypto.Cipher import AES
from Crypto.Random import get_random_bytes
from Crypto.Util.Padding import pad
from qulacs import QuantumState, QuantumCircuit
from qulacs.gate import *
from secret import flag
GATES = {
'CNOT': (CNOT, 2),
'H': (H, 1), 'X': (X, 1), 'Y': (Y, 1), 'Z': (Z, 1),
'S': (S, 1), 'SDAG': (Sdag, 1), 'T': (T, 1), 'TDAG': (Tdag, 1)
}
def parse_circuit(asm, qbits=9):
""" Convert assembly into quantum circuit
i.e. q0 ---+--[Z]--
| <= CNOT 0,1; Z 0; H 1;
q1 --[X]-[H]--
"""
def apply(gate, args):
return gate(*args)
circuit = QuantumCircuit(qbits)
cnt = 0
for instruction in asm.replace('\n', '').split(';'):
t = instruction.strip().split()
if t == []:
continue
if len(t) < 2:
print("[-] Invalid instruction")
exit(0)
opecode, operand = t[0].upper(), t[1:]
if opecode not in GATES:
print("[-] Invalid gate")
exit(0)
operand = list(map(lambda x: int(x), ''.join(t[1:]).split(',')))
if not all(map(lambda x: 0 <= x < qbits, operand)):
print("[-] Invalid quantum bit specified")
exit(0)
if GATES[opecode][1] != len(operand):
print("[-] Invalid number of operands")
exit(0)
gate = apply(GATES[opecode][0], operand)
circuit.add_gate(gate)
cnt += 1
if cnt > 100:
print("[-] Too large circuit")
exit(0)
return circuit
def transfer_and_measure(q):
# Oops, noise might happen :(
noise = QuantumCircuit(9)
idx = random.randrange(0, 9)
noise.add_gate(DephasingNoise(idx, 0.31337))
noise.add_gate(DepolarizingNoise(idx, 0.31337))
noise.add_gate(BitFlipNoise(idx, 0.31337))
noise.update_quantum_state(q)
# Quantum arrived! You can update (or keep) the state :)
circuit = parse_circuit(input('[?] Your Circuit: '))
circuit.update_quantum_state(q)
# Now you measure the quantum state :P
return random.choices(range(2**9),
map(lambda x: abs(x), q.get_vector()))[0] & 1
if __name__ == '__main__':
N = 128
xi, xip = 0.98, 0.98
p = (xi * (1 + xi))**0.5 - xi
Np = int(N * (1 + 2*xi + 2*(xi*(1+xi))**0.5 + xip))
print("[+] Np = " + str(Np))
encoder = parse_circuit("""
CNOT 0,3; CNOT 0,6;
H 0; H 3; H 6;
CNOT 0,1; CNOT 0,2;
CNOT 3,4; CNOT 3,5;
CNOT 6,7; CNOT 6,8;
""")
measured = 0
ra, ba = 0, 0
for i in range(Np):
ra = (ra << 1) | random.choice([0, 1])
ba = (ba << 1) | random.choices([1, 0], [p, 1-p])[0]
q = QuantumState(9)
q.set_zero_state()
if ra & 1:
X(0).update_quantum_state(q)
if ba & 1:
H(0).update_quantum_state(q)
encoder.update_quantum_state(q)
measured = (measured << 1) | transfer_and_measure(q)
del q
print("[+] Measured state: " + bin(measured))
bb = int(input('[?] bb = '), 2)
print("[+] ba = " + bin(ba))
Nx, Nz = 0, 0
for i in range(Np):
if (ba >> i) & 1 == (bb >> i) & 1 == 1:
Nx += 1
elif (ba >> i) & 1 == (bb >> i) & 1 == 0:
Nz += 1
if Nx < N * xi or Nz - N < N * xi:
print("[-] Something went wrong :(")
exit(0)
xa = 0
for i in range(Np):
if (ba >> i) & 1 == (bb >> i) & 1 == 1:
xa = (xa << 1) | ((ra >> i) & 1)
print("[+] xa = " + bin(xa))
l = []
for i in range(Np):
if (ba >> i) & 1 == (bb >> i) & 1 == 0:
l.append(i)
chosen = random.sample(l, Nz - N)
m = 0
for i in chosen:
m |= 1 << i
print("[+] m = " + bin(m))
k = 0
for i in sorted(list(set(l) - set(chosen))):
k = (k << 1) | ((ra >> i) & 1)
key = int.to_bytes(k, N // 8, 'big')
iv = get_random_bytes(AES.block_size)
cipher = AES.new(key, AES.MODE_CBC, iv)
ct = cipher.encrypt(pad(flag, AES.block_size))
print("Result: " + base64.b64encode(iv + ct).decode())
|
cs |
The context of this challenge is quite nice - it's quantum error correction with Shor and QKD with BB84.
I realized the BB84 part of the program after I read the flag, but if you look at it knowing it's BB84, it's quite clear.
Working backwards, we observe the following and eventually find the solution
|
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
|
from qulacs import QuantumCircuit, QuantumState
from qulacs.gate import *
import random
from pwn import *
from Crypto.Cipher import AES
import base64
# CNOT gate
def ADD_CNOT(a, b):
return "CNOT " + str(a) + "," + str(b) + "; "
# single qubit gates
def ADD_single(s, a):
return s + " " + str(a) + "; "
# toffoli
def ADD_toffoli(a, b, c):
ret = ""
ret += ADD_single("H", c)
ret += ADD_CNOT(b, c)
ret += ADD_single("TDAG", c)
ret += ADD_CNOT(a, c)
ret += ADD_single("T", c)
ret += ADD_CNOT(b, c)
ret += ADD_single("TDAG", c)
ret += ADD_CNOT(a, c)
ret += ADD_single("T", b)
ret += ADD_single("T", c)
ret += ADD_CNOT(a, b)
ret += ADD_single("H", c)
ret += ADD_single("T", a)
ret += ADD_single("TDAG", b)
ret += ADD_CNOT(a, b)
return ret
N = 128
xi, xip = 0.98, 0.98
p = (xi * (1 + xi))**0.5 - xi
Np = int(N * (1 + 2*xi + 2*(xi*(1+xi))**0.5 + xip))
# shor error correction
decoder = ""
decoder += ADD_CNOT(0, 1)
decoder += ADD_CNOT(3, 4)
decoder += ADD_CNOT(6, 7)
decoder += ADD_CNOT(0, 2)
decoder += ADD_CNOT(3, 5)
decoder += ADD_CNOT(6, 8)
decoder += ADD_toffoli(1, 2, 0)
decoder += ADD_toffoli(4, 5, 3)
decoder += ADD_toffoli(7, 8, 6)
decoder += ADD_single("H", 0)
decoder += ADD_single("H", 3)
decoder += ADD_single("H", 6)
decoder += ADD_CNOT(0, 3)
decoder += ADD_CNOT(0, 6)
decoder += ADD_toffoli(3, 6, 0)
r = remote('others.ctf.zer0pts.com', '11099')
def get_bin():
s = r.recvline()
s = s.split()[-1].decode()
return int(s, 2)
r.recvline()
for i in range(0, 860):
r.sendline(decoder)
measure = get_bin()
bb = (1 << 400) - 1
r.sendline(bin(bb))
ba = get_bin()
xa = get_bin()
m = get_bin()
l = []
for i in range(Np):
if (ba >> i) & 1 == (bb >> i) & 1 == 0:
l.append(i)
chosen = []
for i in range(Np):
if (m & (1 << i)) > 0:
chosen.append(i)
k = 0
for i in sorted(list(set(l) - set(chosen))):
k = (k << 1) | ((measure >> i) & 1)
key = int.to_bytes(k, N // 8, 'big')
s = r.recvline()
vv = s.split()[-1]
fin = base64.b64decode(vv)
iv = fin[:16]
ct = fin[16:]
cipher = AES.new(key, AES.MODE_CBC, iv)
print(cipher.decrypt(ct))
|
cs |
| zer0lfsr 2021 Quick Solutions/Thoughts (0) | 2021.07.05 |
|---|---|
| LINE CTF Crypto Writeups (0) | 2021.03.22 |
| AeroCTF 2021 Crypto Writeups (0) | 2021.02.28 |
| UnionCTF Crypto Writeups (2) | 2021.02.21 |
| CryptoHack All Solve, Again (10) | 2021.01.31 |
2월에는 나름 잘 쓴 것 같은 좋은 블로그 글을 많이 쓴 것 같아서 기분이 좋다 :)
건강하고 재밌는 한 학기를 보낼 수 있었으면 좋겠다. 부담없이 과목을 수강하니, 재밌는 공부를 더 하고 싶다.
| 3월 정리 (0) | 2021.04.02 |
|---|---|
| 노트북 받고 한 일들 (8) | 2021.03.09 |
| 1월 정산 + 2월 계획 + 1학기 수강신청 (1) | 2021.02.01 |
| 2020 : A Retrospective (9) | 2020.12.27 |
| 1학기 후기 || 최근에 한/할 것들 || 2학기 계획 (0) | 2020.06.23 |
I participated in AeroCTF as a member of Super Guesser. We reached 1st place :)
We were able to grab first bloods on all three cryptography challenges.
Phoenix was solved by rbtree, and I solved horcrux and boggart.
|
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
|
#!/usr/bin/env sage
p = 2
F = GF(p)
P.<x> = PolynomialRing(F)
class Cipher:
def __init__(self, size, params):
self.size = size
self.params = params
def sequence(self, key):
while True:
key = key * self.params[0]
yield key + self.params[1]
def encrypt(self, key, data, strength):
for value, pbit in zip(self.sequence(key), data):
xbit = sum(value[i] for i in range(0, strength, 2))
ybit = mul(value[i] for i in range(1, strength, 2))
yield int(pbit) ^^ int(xbit) ^^ int(ybit)
def main():
size = 256
length = 1024
strength = 10
q = P.irreducible_element(size, 'minimal_weight')
R.<x> = P.quo(q)
key, a, b = [R.random_element() for _ in range(3)]
with open('flag.txt', 'rb') as file:
flag = file.read().strip()
message = int.from_bytes(flag, 'big')
assert message.bit_length() < size
plaintext = list(map(int, bin(message)[2:]))
padding = [0] * (length - len(plaintext))
cipher = Cipher(size, [a, b])
ciphertext = list(cipher.encrypt(key, padding + plaintext, strength))
result = int(''.join(map(str, ciphertext)), 2)
print(a)
print(b)
print(result)
if __name__ == '__main__':
main()
|
cs |
Let's look at the encryption scheme.
We have a key $k \in GF(2^{256})$ and known parameters $a, b \in GF(2^{256})$.
To encrypt a plaintext $m$, it is first converted to binary then padded at the front with zeros to have length 1024.
Then, we take the first 10 bits of the polynomial representation of $a^i k +b$ for each $i$, and calculate $$xbit_i = \sum \text{(even index bits)}$$ $$ybit_i = \prod \text{(odd index bits)}$$ where calculations are done in $GF(2)$.
Finally, we encrypt our text with a stream cipher, i.e. XOR our plaintext with $\{xbit_i \oplus ybit_i\}_{i \ge 1}$.
We now make some observations.
Now we have our solution idea.
We need to have $ybit = 0$ for all 256 bits/equations we selected. This has a probability of $$(31/32)^{256} \approx 0.0003$$ of succeeding, and we also need our system to be uniquely solvable. Thankfully, our algorithm terminates in about 2000 trials, and each trial doesn't take too long. It seems like SageMath is very efficient when it comes to matrices over $GF(2)$ :)
Here's the solution script, by rbtree. It runs in a few minutes.
|
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
|
import random
from Crypto.Util.number import *
p = 2
F = GF(p)
P.<x> = PolynomialRing(F)
size = 256
length = 1024
strength = 10
q = P.irreducible_element(size, 'minimal_weight')
R.<x> = P.quo(q)
class Cipher:
def __init__(self, size, params):
self.size = size
self.params = params
def sequence(self, key):
while True:
key = key * self.params[0]
yield key + self.params[1]
def encrypt(self, key, data, strength):
for value, pbit in zip(self.sequence(key), data):
xbit = sum(value[i] for i in range(0, strength, 2))
ybit = mul(value[i] for i in range(1, strength, 2))
yield int(pbit) ^^ int(xbit) ^^ int(ybit)
enc = 69824286833704501471834043923417254326103912707315595840737453739249974863266259092449058810542265536810346421685955365128856715192808287450464619418781355923155781710833586631897182535937891456025282049302526058466298304955387306232279075295308862156912873485647349272079984781574084434511227361370780842056
enc = list(map(int, bin(enc)[2:]))
enc = [0] * (length - len(enc)) + enc
a = x^255 + x^252 + x^246 + x^245 + x^240 + x^239 + x^236 + x^234 + x^233 + x^232 + x^231 + x^229 + x^228 + x^227 + x^224 + x^223 + x^218 + x^217 + x^216 + x^215 + x^210 + x^209 + x^204 + x^202 + x^200 + x^198 + x^197 + x^196 + x^192 + x^188 + x^187 + x^182 + x^181 + x^180 + x^178 + x^177 + x^176 + x^174 + x^173 + x^172 + x^167 + x^166 + x^161 + x^160 + x^157 + x^155 + x^154 + x^151 + x^150 + x^149 + x^148 + x^147 + x^146 + x^144 + x^140 + x^137 + x^135 + x^133 + x^132 + x^130 + x^129 + x^126 + x^122 + x^119 + x^118 + x^115 + x^112 + x^111 + x^109 + x^107 + x^106 + x^105 + x^104 + x^101 + x^100 + x^99 + x^97 + x^96 + x^94 + x^92 + x^87 + x^86 + x^84 + x^83 + x^81 + x^79 + x^75 + x^71 + x^69 + x^68 + x^67 + x^66 + x^65 + x^63 + x^62 + x^61 + x^56 + x^55 + x^53 + x^52 + x^50 + x^46 + x^44 + x^43 + x^41 + x^39 + x^38 + x^37 + x^36 + x^35 + x^34 + x^33 + x^32 + x^30 + x^29 + x^27 + x^24 + x^21 + x^17 + x^16 + x^14 + x^13 + x^12 + x^11 + x^10 + x^9 + x^5 + x^4 + x^3 + x + 1
b = x^255 + x^254 + x^250 + x^247 + x^243 + x^242 + x^241 + x^238 + x^235 + x^232 + x^229 + x^227 + x^222 + x^221 + x^219 + x^218 + x^217 + x^216 + x^215 + x^211 + x^207 + x^206 + x^204 + x^202 + x^201 + x^197 + x^195 + x^193 + x^192 + x^190 + x^189 + x^188 + x^186 + x^184 + x^181 + x^180 + x^179 + x^178 + x^176 + x^173 + x^172 + x^169 + x^167 + x^165 + x^161 + x^160 + x^158 + x^149 + x^147 + x^146 + x^145 + x^140 + x^138 + x^137 + x^134 + x^133 + x^132 + x^130 + x^129 + x^128 + x^126 + x^125 + x^124 + x^121 + x^120 + x^118 + x^117 + x^114 + x^112 + x^111 + x^110 + x^109 + x^108 + x^107 + x^106 + x^105 + x^101 + x^96 + x^95 + x^94 + x^93 + x^92 + x^90 + x^89 + x^88 + x^86 + x^85 + x^84 + x^83 + x^81 + x^80 + x^79 + x^78 + x^77 + x^76 + x^71 + x^70 + x^69 + x^68 + x^67 + x^64 + x^63 + x^59 + x^56 + x^55 + x^53 + x^50 + x^46 + x^43 + x^42 + x^40 + x^38 + x^37 + x^35 + x^34 + x^33 + x^25 + x^23 + x^22 + x^21 + x^18 + x^16 + x^14 + x^13 + x^12 + x^11 + x^10 + x^8 + x^3 + x^2 + x
t = sum(b[i] for i in range(0, 10, 2))
vs = [ x^i for i in range(size) ]
vecs = []
for idx in range(length - size):
for i in range(size):
vs[i] = vs[i] * a
arr = [int( sum(vs[j][i] for i in range(0, 10, 2)) ) for j in range(size)]
val = int(t) ^^ enc[idx]
vecs.append( (arr, val) )
cipher = Cipher(size, [a, b])
cnt = 0
while True:
cnt += 1
random.shuffle(vecs)
mat = Matrix(F, [vecs[i][0] for i in range(size)])
vec = Matrix(F, [vecs[i][1] for i in range(size)]).transpose()
try:
vec = mat.inverse() * vec
key = 0
for i in range(size):
if int(vec[i][0]):
key += x^i
plaintext = list(cipher.encrypt(key, enc, strength))
plaintext = int(''.join(map(str, plaintext)), 2)
plaintext = long_to_bytes(plaintext)
if b'Aero{' in plaintext:
while True:
print(plaintext)
else:
print("FAIL2", cnt)
except:
print("FAIL", cnt)
pass
|
cs |
|
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
|
#!/usr/bin/env python3.8
from os import urandom
from gmpy2 import next_prime
from random import randrange, getrandbits
from Crypto.Cipher import AES
from fastecdsa.curve import Curve
def bytes_to_long(data):
return int.from_bytes(data, 'big')
def generate_random_point(p):
while True:
a, x, y = (randrange(0, p) for _ in range(3))
b = (pow(y, 2, p) - pow(x, 3, p) - a * x) % p
if (4 * pow(a, 3, p) + 27 * pow(b, 2, p)) % p != 0:
break
return Curve(None, p, a, b, None, x, y).G
class DarkWizard:
def __init__(self, age):
self.power = int(next_prime(getrandbits(age)))
self.magic = generate_random_point(self.power)
self.soul = randrange(0, self.power)
def create_horcrux(self, location, weakness):
# committing a murder
murder = self.cast_spell(b'AVADA KEDAVRA')
# splitting the soul in half
self.soul = self.soul * pow(2, -1, self.power) % self.power
# making a horcrux
horcrux = (self.soul + murder) * self.magic
# nobody should know location and weakness of the horcrux
horcrux.x ^= location
horcrux.y ^= weakness
return horcrux
def cast_spell(self, spell_name):
spell = bytes_to_long(spell_name)
return spell %~ spell
def encrypt(key, plaintext):
cipher = AES.new(key=key, mode=AES.MODE_ECB)
padding = b'\x00' * (AES.block_size - len(plaintext) % AES.block_size)
return cipher.encrypt(plaintext + padding)
def main():
wizard_age = 3000
horcruxes_count = 2
wizard = DarkWizard(wizard_age)
print(f'Wizard\'s power:\n{hex(wizard.power)}\n')
print(f'Wizard\'s magic:\n{wizard.magic}\n')
key = urandom(AES.key_size[0])
horcrux_length = len(key) // horcruxes_count
for i in range(horcruxes_count):
key_part = key[i * horcrux_length:(i + 1) * horcrux_length]
horcrux_location = bytes_to_long(key_part[:horcrux_length // 2])
horcrux_weakness = bytes_to_long(key_part[horcrux_length // 2:])
horcrux = wizard.create_horcrux(horcrux_location, horcrux_weakness)
print(f'Horcrux #{i + 1}:\n{horcrux}\n')
with open('flag.txt', 'rb') as file:
flag = file.read().strip()
ciphertext = encrypt(key, flag)
print(f'Ciphertext:\n{ciphertext.hex()}')
if __name__ == '__main__':
main()
|
cs |
We first take a look at what's going on.
Writing this in mathematical equations, we have
Since XOR is hard to deal with, we change it to normal addition - of course, it's clear that
Therefore, we can change each XOR to addition, by changing our constraints from $0 \le \alpha_i, \beta_i < 2^{32}$ to $|\alpha_i|, |\beta_i| < 2^{32}$.
Now we have to solve, given $p, x_0, y_0, x_1, y_1, x_2, y_2$, the system of equations $$\begin{equation*} \begin{aligned} y_0^2 &= x_0^3 + ax_0 + b \\ (y_1 + \beta_1)^2 &= (x_1 + \alpha_1)^3 + a(x_1 + \alpha_1) + b \\ (y_2 + \beta_2)^2 &= (x_2 + \alpha_2)^3 + a(x_2 + \alpha_2) + b \end{aligned} \end{equation*}$$ over $\mathbb{F}_p$ and $|\alpha_i|, |\beta_i| < 2^{32}$. Now what do we do? There were a few ideas at this point...
Idea 1 : Lattices
Of course, "small solutions" imply that we have to use lattices. Maybe we can directly use them here?
For example, using github.com/rkm0959/Inequality_Solving_with_CVP we can solve a system of linear inequalities.
Turns out, it's quite difficult to do so, since our system has terms like $a \alpha_1$ and $a \alpha_2$.
Because we don't know the value of $a$ and its size can be large, we can't deal with this term using this method.
To use the tool in the link, we need all terms of our equations to be either
In this case, $a \alpha_1$ and $a \alpha_2$ term does not satisfy this requirement.
Idea 2 : Determinant
The failure of our Idea 1 leads to the key idea : removing $a$ in our system completely.
This can be done by realizing that $[1, -a, -b]^{T}$ is in the kernel of the matrix $$ \left[ \begin{matrix} y_0^2 - x_0^3 & x_0 & 1 \\ (y_1 + \beta_1)^2 - (x_1 + \alpha_1)^3 & (x_1 + \alpha_1) & 1 \\ (y_2 + \beta_2)^2 - (x_2 + \alpha_2)^3 & (x_2 + \alpha_2) & 1 \end{matrix} \right]$$ which means its determinant must be $0$ in $\mathbb{F}_p$. This gives us a polynomial equation in $\alpha_1, \alpha_2, \beta_1, \beta_2$. Now what?
The first idea was multivariable coppersmith - I tried to use defund's implementation, but the matrix size got way too large.
The second idea is, of course, lattices from Idea 1. If we expand our determinant, we will get a equation of the form $$\sum \text{(constant coefficient)} \cdot \text{(some monomial)} = 0$$ we can simply bound each monomials with $|\alpha_i|, |\beta_i| < 2^{32}$, then solve using the CVP toolkit above.
To be more precise : our equation will be in a form of $$\sum_{i=1}^k c_i m_i \equiv -c_{k+1} \pmod{p}$$ where $c_i$ is a known constant and $m_i$ is some non-constant monomial over $\alpha_i$, $\beta_i$. (Consider $m_{k+1}$ as $1$)
We solve this system by utilizing the following inequalities and plugging them into the CVP toolkit :
$$-2^{32 \cdot \deg(m_i)} < m_i < 2^{32 \cdot \deg(m_i)}$$ $$ \sum_{i=1}^k c_i m_i+ pX = -c_{k+1}$$ where $\deg(m_i)$ is the degree of the monomial. Note that abuse of notation was used - $m_i$ is used both as a monomial, and as the actual value of the monomial with the true values of $\alpha_i$, $\beta_i$ plugged in. For further details, please take a look in the github repo and the solution.
UPD : Mystiz gives a good insight on entropy to explain why this line of thinking works. Link is here.
This is a very useful way to think about whether or not the above CVP toolkit will work, when it comes to the uniqueness of the solution.
|
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
|
from sage.modules.free_module_integer import IntegerLattice
# https://github.com/rkm0959/Inequality_Solving_with_CVP
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
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
# 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 your result
return result, applied_weights
p = 0x27266a1284b761f793a529b9664693a6b1db36864a8664898d98d8010ec51afffaef2d1c79ab2078c1e0b289a1719d34b0a081ca325ba2b367017ee8e0824aaa9488409e76e923d0f0fc917ddfc1b0534e93a74246405dadfd1683f0dc31682eed0fa6b95fc235c845e16d2ef40463b7e746668dad82981fc4e05933aca410c65b36f89738f7d97502f6626c38f595338f3864638d8613fb74c16b63f3969a49ebd103ef354ed756c3cd5e67f1d2dbe5acdbc088bd6c1d503acef4ec59e4a09efac4729ca796ad25217fe74e7a0c7ef5a3e1fcd9eb9288fb89e842ef0b16642f7e84e27df4bcb623726e2c44ef46be07f9b5a5f92fe2c77d0de79fa6d46193b064207125d8935c2ff04f63e2f858e98d2518077dc58e13307f01d65ae953efd70980f3aeed320b7a6b66eb0c578dc3f05d426f412c4e3c7a9bcc68f27fe236fde41400371a39f53828824f5de3d5902cd3e7dcaee58b89c1a234188e391d412e7bc598d4f10b2bcb26aab7cd09194e80be046022ee8471
x0 = 0x4f69c16e493693a9342b7b9bb0ae42d0e7425996151c631a5620ad4d7ed372d285f04df82975a379080af08322fd927cf8ea9702f4533b5981351dd12358ca61c6e34af3f902b13d2981cae7cd2416e910bfe2c90f3b64c02bf933b1d11ff8ee384f5b507d9a0b1be38b15e7db3ba700ed32122a0163ec66530911fe142e098187771f5a248627a1709800069d402c1a61a17bd053ff6541f981898df6670420b0f3ea6ebd6d01d24f6fe9e3091a7c36df0c1ad7596e8ae090c54528385422911ff8e1dc201ab56844bc92acabc495c442104de5c34a5fe6661b7213500f4ccef6e76b94c34d60772fcf0c1250e66812d85efcda4b323250f740d0d57d40697fedae0bdbaceb0582cd7c82a27610e7f9fba5ac8f84e006d034dcf481f2dc9338acabd51c28afc1b4fc1d0cd7c1cd5ad21109e8708e9458e5301868c68920ab1aefb1e9184383002e6c893f1793443f8388f3f1b6e1f4ecaecb4cf000f57f677156efdacb20d60a35b8337ee416957aadffdc889dce487
y0 = 0x2393c955e5f44ab3ac052dafd83554222728393b8fd20630f3c4f122d8c86d872a7b692b3782f12a6e57352c46328f85fad2d73eaccfcbf7beb47f97da2148c4f3fc7d6751bf66569c97a64f732e7cc71767ba11b1419732035c0285ff6973fba230a3d315fba7820855208e07e6fc5bc4b2cc3868aacde3b8e9b2075bfe861b2ebccd9b6c836d85dc319290263961344d348fe8faa0aa3ebca76fe514a7b7ba313a8200727b22a8714f8bba9e5ec2c549e5bfb5857c050d1ff0471d4b01426c2ee583a7d4b6c30df5ad2d9a6902f574416f8d55ed192b29d521e5d23a54e5062400b539468dca9aa3dc558feac63c88fc696d42434ad9a83551e6860aeb6a4d84ca80713387fa8c56a1e473a82af63ec03a71202aba0ae46fb9a97132f4c92d332327e2c11b79008586b22d92d60c3155e88b5e1f9c193b363ca28f0990400afb6e8148458708b89c6023c0d5ebb746bcd754fe37a84ee4dea6baa273b2ef31a864e9586f01bae855cc0d6f6055b2546b2a918664bdb6
x1 = 0x2700ff7abe81679c770d3171c993c55da4a47cda3360e33d2b763e65517f307a39d401a256ee0634f3f1c6c244aad34422dccfe9ebd1803339972095eb2b0a57c82ee0db9ba94365cfb5270c4924c5c1ec00db2aff529aa923b113d2ca8ad6a6774fb7c655cd101ea63bc5a6ea0261f8da82d455219c7584d7de0757b2fda627c4684d3bac8f899c24178c7e0ecef6c226892b86043d3853cdb777889ce2901d8496bf0232dd000208eb2ce77e953c478551b1112ebf4b02f0086726210a50dc20ac08eb15d846084f9324b4f1f5ec73e3ef7e4a5207e04ebb866f731201e5f626084d1b61c158cfd0fdbaae8b8dd23ba599689c74a790933a3e77daa1c95fbde63b74381aff2b98f41cfebb7f1b220a4f2e8c3361734d7bd6648c720efc0a5f978917d8b84b4764e416762884f00104981e62d876d460bd8c1095cd755d8f31c5377e5ef935da77ff82e823b49d817a1f91bd4d155306173eb07efa4567d362e1cddbf0483873d5efc9286c36a58e5b3d995f3d01ca60
y1 = 0x1bf9a696ef976644dd903fb148892044fe65dd16bf12966a6b6e43be4b3b52aaa348a76ed5a53bbb400a366c59c96cf88a9c6a713273a722c0c8cd6f42b7c6f5e1911d2d323780bce65be80eef4d375dd3425a9ff832e4166765c7687aaf3c6f3c1ae7c561c46c49f0075bb68d48b95be5fd2c94996a42ba255eb638ae5ce449064cc81831f2239dbb93f4944693ee42a507dc2878988928dc9ed5bb5cb3b3eba9ac414b4171e505d285edaf02d13c0748e47b17a498cb0c15883977a11cde98995a022999a555fb3f1a3b9226dc4122f3db6c86036b1b0e6ada10b23b8c89ddb590f2186191c0f1148a04a8e35c905d4053943554c58e8f0e2ccac42c2307532e2b8f55b74fccb8ab24a977c557d10bd7c8621bef3e7f326d00c3ff28a52ee85ec61bb8aef14f67ac80da5885a93f2840a5e44d9800b0352163667ba851bb15c4083955e1fba5cd9d6e7bf103a2bbb0fbc868136cf2871815762514c691e6352af18324d1ccdc21da3e1c4c6aa771aa9fccf343ad64b8
x2 = 0xa439486ae1dbd96ca0796947609757b9dc068d8ed8287f98261c0c77e0940d29e25d27e16b97bc6983be505b8295886a5e880664ecb2d9161036f5beb1dbdac08da0cc2797d574e0feac67bbcb43275ab9e9d936aa17fd9a0a13a2f95e111b5e0893c4d4c3afae3bfa4a60fc5d0673215dd876b8bc960fc765401135dc708562cbeb572fcefaef1fef51782f6e0b495c0799cf89a9f47ffec11b0e780fb41e80d75681f7c9dd5fbf5e93e1dcda8091ca6a84de5b82b3abc9194913119e429781c12c11ebc6b6e80c01683344319879dec8fe59dd2f9c5b7b6e5e211153ae3b5fd161edae59777f3e76ecdeaef518495e512119451c6a97a8f471c4349ba7df8c3e9c1ba67f7a4fef57551d58b8c87607dc830c8074eec50841a1e365be90b528e4f39b8e19e7617191ac5118bb1abc44739d65384d915cb72e226122965ebae1581f55ecb12e523b0e904b3c0a1b26cf506ca030a68fe40d34c0912272277cd0d605ab057dbb5423cca28629661ad163232e766e80949
y2 = 0x29cfc760529c48d68bc086a5b4403f1d4446db04abe243c99baf659b5e67cd6cdac9a658f273f4c682b9e13dda72aaa1ede42c69c98640f2fc58eabeef143a65334e4a236a6e72a157d92ab1a541c9bcf7d953b386a40d68880312ed1e900f6ba481bf134515f5a4c245a0d924db0e5105fc44fae71fc991df85219403ba5d48f68d5d91151f8411b4cbb6971bd63030b523989fa7ec94c14136ac712d4e91eca4c930d79caab7c328043c762917c4b3f868ce037cae9f19a579272c6b13ca8c19d349f5777bf6ac9c078d8472c92582daf96b30d4f7b8fdc004a36b792c133e2b6511956480892636ea91a3361afd8a3afbe3a5bf889feb5a5dc143f5a917347591634218066f2f71b36afd5257f6637152ac9a0965daa881ddc86ca8a8545b255cbb14e7738297a55c7428d9a0e79dee37f801fc1d49d205be809e0cd1b8a1b9d1d5b1b1a9ae98e7c564718cb17d10a3dd01c2914d7bb96c45a4fc942c2ed628354882464407c3208938282471aa2b8016e46c998e4
R = Integers(p)
# determinant calculation
# a, b = alpha_1, beta_1
# c, d = alpha_2, beta_2
P.<a, b, c, d> = PolynomialRing(R)
f = 0
f += (y0 ** 2 - x0 ** 3) * (x1 + a) * 1
f += (x0) * 1 * ((y2 + d) ** 2 - (x2 + c) ** 3)
f += 1 * ((y1 + b) ** 2 - (x1 + a) ** 3) * (x2 + c)
f -= (y0 ** 2 - x0 ** 3) * 1 * (x2 + c)
f -= x0 * ((y1 + b) ** 2 - (x1 + a) ** 3) * (1)
f -= 1 * (x1 + a) * ((y2 + d) ** 2 - (x2 + c) ** 3)
# coefficient
C = f.coefficients()
# degree
print(f.monomials())
# [a^3*c, a*c^3, a^3, a^2*c, b^2*c, a*c^2, c^3, a*d^2, a^2, b^2, a*c, b*c, c^2, a*d, d^2, a, b, c, d, 1]
D = [poly.degree() for poly in f.monomials()]
print(D)
# [4, 4, 3, 3, 3, 3, 3, 3, 2, 2, 2, 2, 2, 2, 2, 1, 1, 1, 1, 0]
# CVP setup
print(len(C)) # 20
M = Matrix(ZZ, 20, 20)
lb = [0] * 20
ub = [0] * 20
# monomial bounds
for i in range(0, 19):
M[i, i] = 1
lb[i] = -(2 ** (32 * D[i]))
ub[i] = (2 ** (32 * D[i]))
# the polynomial
for i in range(0, 19):
M[i, 19] = C[i]
M[19, 19] = p
lb[19] = int((-C[19]) % p)
ub[19] = int((-C[19]) % p)
# solve CVP
result, applied_weights = solve(M, lb, ub)
# recover variables
alpha_1 = result[15] // applied_weights[15]
beta_1 = result[16] // applied_weights[16]
alpha_2 = result[17] // applied_weights[17]
beta_2 = result[18] // applied_weights[18]
print(alpha_1, beta_1, alpha_2, beta_2)
# recover actual key
alpha_1 = x1 ^^ (x1 + alpha_1)
beta_1 = y1 ^^ (y1 + beta_1)
alpha_2 = x2 ^^ (x2 + alpha_2)
beta_2 = y2 ^^ (y2 + beta_2)
print(alpha_1, beta_1, alpha_2, beta_2)
# the rest is trivial decryption
|
cs |
|
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
|
#!/usr/bin/env python3.8
from gmpy import next_prime
from random import getrandbits
def bytes_to_long(data):
return int.from_bytes(data, 'big')
class Wardrobe:
@staticmethod
def create_boggarts(fear, danger):
# for each level of danger we're increasing fear
while danger > 0:
fear = next_prime(fear)
if getrandbits(1):
yield fear
danger -= 1
class Wizard:
def __init__(self, magic, year, experience):
self.magic = magic
self.knowledge = year - 1 # the wizard is currently studying the current year
self.experience = experience
def cast_riddikulus(self, boggart):
# increasing the wizard's experience by casting the riddikulus charm
knowledge, experience = self.knowledge, self.experience
while boggart > 1:
knowledge, experience = experience, (experience * self.experience - knowledge) % self.magic
boggart -= 1
self.experience = experience
def main():
year = 3
bits = 512
boggart_fear = 31337
boggart_danger = 16
neutral_magic, light_magic, dark_magic = [getrandbits(bits) for _ in range(3)]
magic = next_prime(neutral_magic | light_magic) * next_prime(neutral_magic | dark_magic)
print('Hello. I am Professor Remus Lupin. Today I\'m going to show you how to deal with the boggart.')
print(neutral_magic)
print(magic)
with open('flag.txt', 'rb') as file:
flag = file.read().strip()
# some young wizards without knowledge of the riddikulus charm
harry_potter = Wizard(magic, year, bytes_to_long(b'the boy who lived'))
you = Wizard(magic, year, bytes_to_long(flag))
for boggart in Wardrobe.create_boggarts(boggart_fear, boggart_danger):
# wizards should train to learn the riddikulus charm
harry_potter.cast_riddikulus(boggart)
you.cast_riddikulus(boggart)
# wizard's experience should be increased
print(harry_potter.experience)
print(you.experience)
if __name__ == '__main__':
main()
|
cs |
There is a lot going on in this code, so let's take a look at it.
Now, the text "the boy who lived" and the flag goes through some encryption process, and we know their ciphertexts.
cast_riddikulus algorithm
Naturally, we have to look at the cast_riddikulus algorithm and how to reverse/decrypt it.
We note that the value of self.knowledge is fixed at $2$ at all times. Denote the self.experience as $x$, and the boggart as $n$.
The algorithm calculates a linear recurrence - it starts with $a_0=2$, $a_1 = x$, and moves up by $$a_{i} = xa_{i-1} - a_{i-2}$$ Let's forget about the modulos and solve this recurrence. Using standard techniques, we get $$a_n = \left( \frac{x - \sqrt{x^2 - 4}}{2} \right)^n + \left( \frac{x + \sqrt{x^2- 4}}{2} \right)^n $$ This looks kinda scary, but then we can see that this is just $$x = t + t^{-1} \implies a_n = t^n + t^{-n}$$ Note that this implication holds even when we consider the fact that we are working modulo $M$.
This has a lot of meaning - for example, the encryption function is commutative in a sense that encrypting a text (experience) with (boggart) $m$ then $n$ is equivalent to encrypting with $n$ then $m$, which is equivalent to encryption once with $mn=nm$.
It also implies that under some "mild" conditions, we can reverse the encryption function as follows -
In the above algorithm, we needed
Let's review our requirements and how difficult they are
Therefore, it's time to go back and try to factorize $M$ with our information of $M_N$.
factorization of $M$
The key idea is that the value of $M_N$ gives us known bits of each prime $p, q$ of $M$.
Consider $M'_N = M_N - (M_N \pmod{2^{16}})$, i.e. $M_N$ but changed the $16$ LSB to 0.
Therefore, we write $M''_N = M'_N \& (M'_N + 2^{16})$ and claim that $p, q$ all have the bits of $M''_N$.
It is known that with around 57% of the bits of $p, q$ known, we can factorize $N$ using the following algorithm -
We require that $M''_N$ has a lot of ones in its binary representation. We simply restart the connection until we get good parameters.
For the details of the factorization algorithm, please refer to the solution script below.
how did we encrypt this?
Now we know the factorization of $M$. We can now check if $ctxt^2 - 4$ is a quadratic residue modulo $M$.
We test this for the ciphertext of the "the boy who lived" and the ciphertext of the flag. If this check fails, we run our algorithm again with a new connection for a good set of parameters. We now assume $ctxt^2 - 4$ is a quadratic residue.
We reasonably guess that our encryption algorithm will have used 16 primes out of the 35 consecutive primes larger than 31337.
We also verify that $p-1$ and $q-1$ are all coprime to these 35 candidate primes. If not, we look for a good set of parameters.
Now we can encrypt and decrypt as we like. Let's find out which 16 primes were used for encryption.
To do this, we use meet-in-the-middle algorithm. Basically what we do is
Now we know which primes were used for encryption. All we need to do is decrypt the ciphertext of the flag.
additional notes
Note 1 : speeding up cast_riddikulus algorithm
Note 2 : requirements relaxation
Note 3 : meet-in-the-middle, but not quite middle
The solution code is displayed in parts for reader's understanding. Feel free to comment and ask questions...
reading parameters
|
1
2
3
4
5
6
7
8
9
10
11
12
13
|
def recvint():
T = r.recvline()
return int(T[:-1].decode(), 10)
def get_vals():
r.recvline()
C = recvint()
N = recvint()
H = bytes_to_long(b'the boy who lived')
H_AFTER = recvint()
F_AFTER = recvint()
return C, N, H_AFTER, F_AFTER
|
cs |
factorization of $M$
|
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
|
def bsum(x):
ret = 0
while x > 0:
ret += x % 2
x //= 2
return ret
def factor(N, C, cand, k):
if len(cand) == 0:
return
print(k, len(cand))
# fix pq mod 2^{k+1} by adding 2^k or not
ret = []
if (C >> k) % 2 == 1: # we know p, q has 2^k
for p, q in cand:
if p * q == N:
print(p)
print(q)
return
pp = p + (1 << k)
qq = q + (1 << k)
if (pp * qq - N) % (1 << (k + 1)) == 0 and pp * qq <= N:
ret.append((pp, qq))
factor(N, C, ret, k + 1)
else:
for p, q in cand:
if p * q == N:
print(p)
print(q)
return
for i in range(0, 2):
for j in range(0, 2):
pp = p + i * (1 << k)
qq = q + j * (1 << k)
if (pp * qq - N) % (1 << (k + 1)) == 0 and pp * qq <= N:
ret.append((pp, qq))
factor(N, C, ret, k + 1)
# M_N to M''_N
C -= C % (1 << 16)
C &= (C + (1 << 16))
# check number of 1s, if bsum(C) is small (<256), stop and rerun
print(C.bit_length(), bsum(C))
# factorize
factor(N, C, [(1, 1)], 1)
|
cs |
encryption and decryption
|
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
|
def mlucas(v, a, n): # encryption fast, from Wikipedia
v1, v2 = v, (v**2 - 2) % n
for bit in bin(a)[3:]:
if bit == "0":
v1, v2 = ((v1**2 - 2) % n, (v1*v2 - v) % n)
else:
v1, v2 = ((v1*v2 - v) % n, (v2**2 - 2) % n)
return v1
def REV(val, ex, pr): # decryption modulo a prime
# t^n + 1/t^n = val
if val == 0: # if invalid
return 0 # return invalid
cc = tonelli((val * val - 4) % pr, pr) # modular sqrt
if cc == -1: # failure
return 0 # return failure
tn = ((val + cc) * inverse(2, pr)) % pr
t = pow(tn, inverse(ex, pr-1), pr)
return ((t + inverse(t, pr))) % pr
|
cs |
meet-in-the-middle
|
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
|
primes = [31337, 31357, 31379, 31387, 31391, 31393, 31397, 31469, 31477, 31481, 31489, 31511, 31513, 31517, 31531, 31541, 31543, 31547, 31567, 31573, 31583, 31601, 31607, 31627, 31643, 31649, 31657, 31663, 31667, 31687, 31699, 31721, 31723, 31727, 31729, 31741, 31751, 31769, 31771, 31793]
L = primes[15:35]
R = primes[:15]
LV = []
RV = []
# generally, encryption is faster than decryption
# therefore, we leave the "heavy" work for the encryption part...
for combo_L in tqdm(itertools.combinations(L, 8)):
CC = list(combo_L)
cur = H
# you can also just multiply all primes and encrypt at once
for pr in combo_L:
cur = mlucas(cur, pr, p)
CC.append(cur)
LV.append(CC)
print("LEFT DONE")
for combo_R in tqdm(itertools.combinations(R, 8)):
CC = list(combo_R)
cur = HAFTER
# you can also just multiply all primes and decrypt at once
for pr in combo_R:
cur = REV(cur, pr, p)
CC.append(cur)
RV.append(CC)
print("RIGHT DONE")
for i in tqdm(range(len(LV))):
for j in range(len(RV)):
if LV[i][8] == RV[j][8]:
print("FOUND!!")
print(LV[i])
print(RV[j])
|
cs |
the finish, and the parameters I got
|
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
|
H = bytes_to_long(b'the boy who lived')
sys.setrecursionlimit(10 ** 6)
C = 12124756848659098434025945489515506912896022954145117746118560512007665385702760439414990812257455576297349156226093149988609289245714223348281989890389750
N = 157597985389833012337654133040126048344064622845161536236706459270343905778002470548499258715513540516431526518690532156245280894778788692043941237295606686168037171464988128988463706375526180496632421973522548093894845498612792150825707672843107252573999144787226703076358545319417530365329436368718460943493
HAFTER = 66161881147822169408519540711422900962287264738494143175834051626001922954586728648835878096124744364195826536091510407493007528877139856387261499433277826944946254511824024047480941829026088269865298686453128715170657018128276813244425143986311708022950785583195028647859774987948632731985531259912781472862
FAFTER = 149186530719822614329126547638374064715014252925601014676661223009475822460330945440469384214084001910035138025738722725987466200681944900264994344927428683388976167111544750466576538413516454786176229441173029050647235653998791477157269246962955063391947778663841551982999293815571149539542758304215156142104
C -= C % (1 << 16)
C &= (C + (1 << 16))
# factor(N, C, [(1, 1)], 1)
p = 12558711464274427739528720572494472142909592647353129013838950445222814801805965383430302364628487022743397586481672449715551542652546057434522020868473011
q = 12548897698474048380978452887676419841595083766206501465313606366388795637681128899285066184154566275021196792336453455837893284460576050862858626214885863
# recovered set of primes
pr = [31543, 31567, 31573, 31607, 31649, 31667, 31687, 31699, 31357, 31379, 31397, 31469, 31477, 31513, 31517, 31531]
hp = FAFTER % p
hq = FAFTER % q
# reverse
for i in pr:
hp = REV(hp, i, p)
hq = REV(hq, i, q)
flag, tt = CRT(hp, p, hq, q) # CRT
print(long_to_bytes(flag))
|
cs |
| LINE CTF Crypto Writeups (0) | 2021.03.22 |
|---|---|
| zer0pts CTF 2021 Crypto Writeups (1) | 2021.03.07 |
| UnionCTF Crypto Writeups (2) | 2021.02.21 |
| CryptoHack All Solve, Again (10) | 2021.01.31 |
| 0x41 CTF SPN Writeup (0) | 2021.01.31 |
이어서 간다. 이번 글에서 다룰 주제는
결과부터 소개한다. 우리가 앞에서 준비했던, 논문의 정리 4.1.만 소개한다. 그래도 논문 링크를 열고 정리 4.1.~4.4.를 한 번 확인하는 것을 추천.
결과도 매우 중요하지만, 아무래도 일단 관심이 가는 것은 대체 이렇게 복잡한 것을 어떻게 증명했냐는 것이다. 논문의 말을 옮기면,
Appendix를 구경하면 이 결과의 엄밀한 증명을 직접 볼 수 있다. ㄹㅇ 가슴이 웅장해진다...
여기서는 operator class $\mathcal{Q}_1$, $\mathcal{Q}_2$를 알고 있을 때 parameter $\alpha, \theta$를 optimal 하게 설정하는 방법을 다룬다.
이는 결국 $\rho^2_{*} (\alpha, \theta)$를 해당 $\alpha, \theta$에 대한 optimal contraction factor라고 할 때, $$\rho^2_{*} = \inf_{\alpha, \theta > 0} \rho^2_{*}(\alpha, \theta)$$와 대응되는 $\alpha_{*}$, $\theta_{*}$를 찾는 문제와 같다. 원래 $\theta \in (0, 2)$라는 조건이 있으나, 이 조건은 relax 해도 무방하다.
$\rho^2_{*}(\alpha, \theta)$의 값 자체는 SDP를 풀어서 계산할 수 있으니, 이차원 위에서 미분 불가능한 함수를 최소화하는 일반적인 알고리즘이 있다면 이를 사용하여 문제를 해결할 수도 있을 것이다. 그러나, 조금 더 작업을 거치면 더욱 편하게 문제를 풀 수 있다. $\theta$를 최적화 문제 안에 집어넣자.
Dual Problem $$\begin{equation*} \begin{aligned} \text{minimize} \quad & \rho^2 \\ \text{subject to} \quad & \lambda_\mu^A, \lambda_\beta^B \ge 0 \\ \quad & S = -M_O - \lambda_\mu^A M_\mu^A - \lambda_\beta^B M_{\beta}^B + \rho^2 M_I \succeq 0 \end{aligned} \end{equation*}$$을 다시 생각해보자. 우리가 $\alpha, \theta$를 고정된 값이 아니라 최적화 문제의 변수로 두지 못하는 이유는, 행렬 $S$의 entry가 $\alpha, \theta$의 nonlinear term을 포함하고 있기 때문이다. $\alpha$는 아예 분모에 있고, $\theta$는 $\theta^2$ 항이 $S$에 있다. 여기서 $\theta$에 대한 quadratic term을 제거하기 위해, Schur Complement를 사용한다.
$$\begin{equation*} \begin{aligned} S &= -M_O - \lambda_\mu^A M_\mu^A - \lambda_\beta^B M_{\beta}^B + \rho^2 M_I \\ S' &= - \lambda_\mu^A M_\mu^A - \lambda_\beta^B M_{\beta}^B + \rho^2 M_I \end{aligned} \end{equation*} $$라고 정의하면, $$M_O = (1, \theta, -\theta)^T (1, \theta, -\theta)$$이므로, Schur Complement의 성질에서 $S \succeq 0$는 $v = (1, \theta, -\theta)^T$라 했을 때 $$\overline{S} = \left[ \begin{matrix} S' & v \\ v^T & 1 \end{matrix} \right] \succeq 0$$과 동치이다. 또한, $S'$은 각 entry에 $\theta$에 대한 non-linear term이 없으므로, $\theta$까지 SDP에 넣을 수 있다.
즉, 새로운 문제를 $$\begin{equation*} \begin{aligned} \text{minimize} \quad & \rho^2 \\ \text{subject to} \quad & \lambda_\mu^A, \lambda_\beta^B, \theta \ge 0 \\ \quad & S' = - \lambda_\mu^A M_\mu^A - \lambda_\beta^B M_{\beta}^B + \rho^2 M_I \\ \quad & v = (1, \theta, -\theta)^T \\ \quad & \overline{S} = \left[ \begin{matrix} S' & v \\ v^T & 1 \end{matrix} \right] \succeq 0 \end{aligned} \end{equation*}$$로 바꾸면, $\theta$에 대한 최적화까지 자동으로 된다. 이제 남은 것은 $\alpha$에 대한 최적화다.
이렇게 해서 얻은 주어진 $\alpha$에 대한 optimal contraction factor를 $\rho^2_{*}(\alpha)$라고 하자. 이를 최소화하면 된다.
변수 하나에 대한 (미분 불가능한) 함수의 minimization은 다양한 방법이 있다. 그런데, 실험 결과에 따르면
그러므로, 이를 최소화하기 위해서 (PS를 했던 사람의 입장에서 쉽게 생각하면) 삼분탐색을 사용해도 좋다.
사실 scipy나 matlab 등에서 이미 derivative-free optimization이 구현이 되어있기 때문에, 이걸 가져다 쓰면 된다.
이건 논문에는 없고 (정확하게는 CVXPY로 하지 않고 MATLAB 등으로 한 것 같다) 그냥 내가 심심해서 하는 것이다.
논문의 저자들이 사용한 코드를 보고 싶다면, github.com/AdrienTaylor/OperatorSplittingPerformanceEstimation로 가자.
내가 여기서 짠 코드들은 github.com/rkm0959/rkm0959_implements에 가면 있다.
CVXPY를 이용하여 Theorem 4.1.을 검증해보자. 방법은 다음과 같다.
이를 구현한 코드는 다음과 같다.
|
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
|
import cvxpy as cp
import matplotlib.pyplot as plt
import numpy as np
import math
def THM_4_1_scaled(theta, mu, beta):
assert 0 < theta < 2 and mu > 0 and beta > 0
if mu * beta - mu + beta < 0 and theta <= 2 * (beta + 1) * (mu - beta - mu * beta) / (mu + mu * beta - beta - beta * beta - 2 * mu * beta * beta):
return 1, abs(1 - theta * beta / (beta + 1))
if mu * beta - mu - beta > 0 and theta <= 2 * (mu * mu + beta * beta + mu * beta + mu + beta - mu * mu * beta * beta) / (mu * mu + beta * beta + mu * mu * beta + mu * beta * beta + mu + beta - 2 * mu * mu * beta * beta):
return 2, abs(1 - theta * (1 + mu * beta) / (1 + mu) / (1 + beta))
if theta >= 2 * (mu * beta + mu + beta) / (2 * mu * beta + mu + beta):
return 3, abs(1 - theta)
if mu * beta + mu - beta < 0 and theta <= 2 * (mu + 1) * (beta - mu - mu * beta) / (beta + mu * beta - mu - mu * mu - 2 * mu * mu * beta):
return 4, abs(1 - theta * mu / (mu + 1))
ret = (2 - theta) / 4
ret = ret * ((2 - theta) * mu * (beta + 1) + theta * beta * (1 - mu))
ret = ret * ((2 - theta) * beta * (mu + 1) + theta * mu * (1 - beta))
ret = ret / mu / beta
ret = ret / (2 * mu * beta * (1 - theta) + (2 - theta) * (mu + beta + 1))
return 5, math.sqrt(ret)
def THM_4_1_unscaled(alpha, theta, mu, beta):
return THM_4_1_scaled(theta, mu * alpha, beta / alpha)
def gen_M_O(theta):
ret = [[1, theta, -theta],
[theta, theta * theta, -theta * theta],
[-theta, -theta * theta, theta * theta]]
return np.array(ret)
def gen_M_I():
ret = np.zeros((3, 3))
ret[0, 0] = 1
return ret
def gen_M_mu_A(alpha, mu):
ret = [[0, -1/2, 0],
[-1/2, -(1+alpha*mu), 1],
[0, 1, 0]]
return np.array(ret)
def gen_M_beta_B(alpha, beta):
ret = [[-beta/alpha, 0, beta/alpha + 1/2],
[0, 0, 0],
[beta/alpha + 1/2, 0, -beta/alpha-1]]
return np.array(ret)
def Primal(alpha, theta, mu, beta):
G = cp.Variable((3, 3), symmetric = True)
M_O = gen_M_O(theta)
M_I = gen_M_I()
M_mu_A = gen_M_mu_A(alpha, mu)
M_beta_B = gen_M_beta_B(alpha, beta)
constraints = []
constraints.append(G >> 0)
constraints.append(cp.trace(M_I @ G) == 1)
constraints.append(cp.trace(M_mu_A @ G) >= 0)
constraints.append(cp.trace(M_beta_B @ G) >= 0)
objective = cp.Maximize(cp.trace(M_O @ G))
problem = cp.Problem(objective, constraints)
problem.solve()
return math.sqrt(problem.value)
def Dual(alpha, theta, mu, beta):
lambda_mu_A = cp.Variable()
lambda_beta_B = cp.Variable()
rho_sq = cp.Variable()
M_O = gen_M_O(theta)
M_I = gen_M_I()
M_mu_A = gen_M_mu_A(alpha, mu)
M_beta_B = gen_M_beta_B(alpha, beta)
constraints = []
constraints.append(lambda_mu_A >= 0)
constraints.append(lambda_beta_B >= 0)
constraints.append(-M_O - lambda_mu_A * M_mu_A - lambda_beta_B * M_beta_B + rho_sq * M_I >> 0)
objective = cp.Minimize(rho_sq)
problem = cp.Problem(objective, constraints)
problem.solve()
return math.sqrt(problem.value)
eps = 0.0005
checked = [0, 0, 0, 0, 0]
for _ in range(300):
alpha = np.random.uniform(low = 0.5, high = 2.0)
theta = np.random.uniform(low = 0.05, high = 1.95)
mu = np.random.uniform(low = 0.1, high = 3.9)
beta = np.random.uniform(low = 0.1, high = 3.9)
assert alpha > 0 and theta > 0
assert mu > 0 and beta > 0
val_1 = Primal(alpha, theta, mu, beta)
val_2 = Dual(alpha, theta, mu, beta)
whi, val_3 = THM_4_1_unscaled(alpha, theta, mu, beta)
checked[whi - 1] += 1
if abs(val_2 - val_1) > eps or abs(val_3 - val_1) > eps or abs(val_3 - val_2) > eps:
print(val_1, val_2, val_3)
print("Failed at Case", whi)
print("Parameters:", alpha, theta, mu, beta)
break
print("Checked Cases")
for i in range(5):
print("Case #{} Verified {} Times!".format(i + 1, checked[i]))
print("Finished!")
|
cs |
이렇게 하여 Theorem 4.1.의 결과를 검수할 수 있다. SDP solver가 문제를 잘 해결할 수 있도록 극단적인 값들은 넣지 않았다.
각 코드의 의미는 거의 자명하니 (식을 그대로 옮긴 것에 불과하다) 자세한 설명은 생략하도록 하겠다. Task 1은 이렇게 끝.
이번에는 Part 4에서 다루었던 내용을 구현해보자. 우리가 할 일은
Part 4의 내용을 그대로 CVXPY로 옮기면, 다음과 같은 코드를 얻을 수 있다.
|
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
|
import cvxpy as cp
import matplotlib.pyplot as plt
import numpy as np
import math
def THM_4_1_scaled(theta, mu, beta):
assert 0 < theta < 2 and mu > 0 and beta > 0
if mu * beta - mu + beta < 0 and theta <= 2 * (beta + 1) * (mu - beta - mu * beta) / (mu + mu * beta - beta - beta * beta - 2 * mu * beta * beta):
return 1, abs(1 - theta * beta / (beta + 1))
if mu * beta - mu - beta > 0 and theta <= 2 * (mu * mu + beta * beta + mu * beta + mu + beta - mu * mu * beta * beta) / (mu * mu + beta * beta + mu * mu * beta + mu * beta * beta + mu + beta - 2 * mu * mu * beta * beta):
return 2, abs(1 - theta * (1 + mu * beta) / (1 + mu) / (1 + beta))
if theta >= 2 * (mu * beta + mu + beta) / (2 * mu * beta + mu + beta):
return 3, abs(1 - theta)
if mu * beta + mu - beta < 0 and theta <= 2 * (mu + 1) * (beta - mu - mu * beta) / (beta + mu * beta - mu - mu * mu - 2 * mu * mu * beta):
return 4, abs(1 - theta * mu / (mu + 1))
ret = (2 - theta) / 4
ret = ret * ((2 - theta) * mu * (beta + 1) + theta * beta * (1 - mu))
ret = ret * ((2 - theta) * beta * (mu + 1) + theta * mu * (1 - beta))
ret = ret / mu / beta
ret = ret / (2 * mu * beta * (1 - theta) + (2 - theta) * (mu + beta + 1))
return 5, math.sqrt(ret)
def THM_4_1_unscaled(alpha, theta, mu, beta):
return THM_4_1_scaled(theta, mu * alpha, beta / alpha)
def gen_Base():
ret = np.zeros((4, 4))
ret[0, 3] = 1
ret[3, 0] = 1
ret[3, 3] = 1
return ret
def gen_M_theta():
ret = np.zeros((4, 4))
ret[1, 3] = 1
ret[2, 3] = -1
ret[3, 1] = 1
ret[3, 2] = -1
return ret
def gen_M_I():
ret = np.zeros((4, 4))
ret[0, 0] = 1
return ret
def gen_M_mu_A(alpha, mu):
ret = [[0, -1/2, 0, 0],
[-1/2, -(1+alpha*mu), 1, 0],
[0, 1, 0, 0],
[0, 0, 0, 0]]
return np.array(ret)
def gen_M_beta_B(alpha, beta):
ret = [[-beta/alpha, 0, beta/alpha + 1/2, 0],
[0, 0, 0, 0],
[beta/alpha + 1/2, 0, -beta/alpha-1, 0],
[0, 0, 0, 0]]
return np.array(ret)
def opt_val(alpha, mu, beta, retrieve = False):
lambda_mu_A = cp.Variable()
lambda_beta_B = cp.Variable()
theta = cp.Variable()
rho_sq = cp.Variable()
Base = gen_Base()
M_theta = gen_M_theta()
M_I = gen_M_I()
M_mu_A = gen_M_mu_A(alpha, mu)
M_beta_B = gen_M_beta_B(alpha, beta)
constraints = []
constraints.append(lambda_mu_A >= 0)
constraints.append(lambda_beta_B >= 0)
constraints.append(theta >= 0)
constraints.append(Base + theta * M_theta - lambda_mu_A * M_mu_A - lambda_beta_B * M_beta_B + rho_sq * M_I >> 0)
objective = cp.Minimize(rho_sq)
problem = cp.Problem(objective, constraints)
problem.solve()
if retrieve == True:
return math.sqrt(problem.value), theta.value
else:
return math.sqrt(problem.value)
mu = 0.53
beta = 1.35
alpha_L = 0.05
alpha_R = 3.95
print("Begin Ternary Search")
while alpha_R - alpha_L >= 0.0001:
left = alpha_L + (alpha_R - alpha_L) / 3
right = alpha_R - (alpha_R - alpha_L) / 3
if opt_val(left, mu, beta) < opt_val(right, mu, beta):
alpha_R = right
else:
alpha_L = left
opt_alpha = (alpha_L + alpha_R) / 2
opt_rho, opt_theta = opt_val(opt_alpha, mu, beta, True)
print("Optimal alpha, theta", opt_alpha, opt_theta)
print("Optimal contraction factor", opt_rho)
print("Check with Theorem 4.1.", THM_4_1_unscaled(opt_alpha, opt_theta, mu, beta))
print("Begin Graphing")
X = np.linspace(0.25, 3.75, 50)
Y = np.array([opt_val(alpha, mu, beta) for alpha in X])
plt.plot(X, Y)
plt.show()
plt.savefig("Graph.PNG")
print("Finished!")
|
cs |
여기서 얻는 결론은 $\mu = 0.53$, $\beta = 1.35$일 때, optimal parameter는 $$\alpha_{*} = 1.5949751, \quad \theta_{*} = 1.4244099$$이고, 이때 얻어지는 optimal contraction factor는 $$\rho_{*} = 0.5010598$$이라는 것이다. 또한, 이 경우에 대한 $\rho^2_{*} (\alpha)$의 그래프를 그려서 unimodal function임을 확인할 수 있다. 그래프는 아래와 같다.
당연하지만, $x$축은 $\alpha$의 값, $y$축은 $\rho^2_{*}(\alpha)$의 값이다. 이렇게 Task 2와 논문 정리도 끝 :)
| Acceleration Overview (2) | 2021.03.28 |
|---|---|
| Nesterov's Estimate Sequence (0) | 2021.03.15 |
| [Computer-Assisted Proofs in Optimization] Operator Splitting Performance Estimation Problem - 2 (0) | 2021.02.23 |
| [Computer-Assisted Proofs in Optimization] Operator Splitting Performance Estimation Problem - 1 (0) | 2021.02.23 |
| [Computer-Assisted Proofs in Optimization] Performance Estimation Problem (0) | 2021.02.23 |