Let $p = 10^9 + 7$.
Euler's theorem (the same one you cited) tells us $21^{p-1} \equiv 1 \mod p$. Thus when raising $21$ to some huge power, we can reduce that exponent mod $p-1$ without changing the final result. It turns out that $10^{18} \equiv 36 \mod (p-1)$, so $21^{10^{18}} \equiv 21^{36} \mod p$, and the modular inverse of $21^{10^{18}}$ will be the same as the modular inverse of $21^{36}$.
Next we should use the same method you were starting to try before. The inverse of $21^{36}$ mod $p$ will be $\left(21^{36}\right)^{p-2}$ mod $p$.
From here, we could brute force the answer using $p-2$ multiplications on a computer, but $10^9$ operations may take awhile depending which language you use. Instead, we can speed up the process using the exponentiation by squaring algorithm, which lets us use around $\log_2(p)$ operations instead of $p$ operations.
Here's a short Python script that puts all these steps together.
p = 10**9 + 7
e = 10**18 % (p-1) # aka 36
b = pow(21, e, p) # b is 21^(10^18) mod p
result = pow(b, p-2, p)
print(result)
The final answer returned is 753568836
. Note that in Python, pow(x, y, z)
means $x^y \mod z$, and it runs fast because it incorporates "exponentiation by squaring" (or something similar). The code runs instantly on my laptop.
pow()
function nowadays hides all the maths. It is a black box. Trypow(21, -1, 10**9 + 7)
– Gribouillis Dec 06 '22 at 09:26