# OK to post
# HW3, RDB, 2024-01-28

# Thursday homework.

# 2.
 
# EA and EAnew both compute the gcd, but EA does it via iteration rather than
# recursion. Calling a function results in a small amount of overhead. In a
# compiled language EA would likely be optimized to use iteration.

# EAnew is still not as good as Maple's builtin gcd. There are apparently lots
# of subtle algorithmic improvements that can be made to the Euclidean
# algorithm. I would guess that Maple spent a lot of time implementing these,
# and we just did the "naive" thing.

# 3.

EstimateProbGCD := proc(N, K)
    local relprimes, r, k, a, b:
    relprimes := 0:
    r := rand(1..N):
    for k from 1 to K do
        a := r():
        b := r():
        if igcd(a, b) = 1 then
            relprimes += 1:
        fi:
    od:

    relprimes / K:
end:

# I ran EstimateProbGCD(10^10, 10^6) a few times.
# The results were not too varied, somewhere around 0.60.

# 4.

# The actual limit is

#   1 / Zeta(2) = 6 / Pi^2.

# If you are happy with that, then stop reading now!

# We want to get the asymptotics of the sum

#   sum([gcd(a, b) = 1], a, b <= n).

# Maybe, if I was smarter, I would see some counting argument here about
# primes. Because I am not, I will try a number theory trick I know.

# If mu(d) is the Mobius function, then

#   [gcd(a, b) = 1] = sum(mu(d), d | gcd(a, b)).

# So our sum is

#   sum(mu(d), d | gcd(a, b), a, b <= n)

# Now, d | gcd(a, b) iff d | a and d | b. So we get

#   sum(mu(d) [d | a][d | b], d >= 1, a, b <= n).

# If we sum over a, then b, we get

#   sum(mu(d) floor(n / d)^2, d >= 1).

# If we write

#   floor(n / d)^2 = ((n / d) + {n / d})^2
#                  = n^2 / d^2 + (n / d){n / d} + {n/d}^2,

# then our sum is

#   n^2 sum(mu(d) / d^2, d >= 1) + sum(mu(d) ((n/d) {n/d} + {n/d}^2), d >= 1)
#   n^2 / Zeta(2) + (stuff).

# I'm almost certain that (stuff) is O(n).

# For terms with d >= n, we can use the triangle inequality to get a bound of

#   sum({n / d} (n / d + {n / d}), d >= n).
#   =
#   2 sum((n / d)^2, d >= n).
#   =
#   O(n).

# For terms with d < n, we have

#   sum(mu(d) {n / d} (n / d + {n / d}), d < n).

# I am pretty sure this is also O(n), but it is annoying to show.

# Putting everything together, our sum is

#   n^2 / Zeta(2) + O(n),

# and dividing by n^2 gives 1 / Zeta(2) + O(1/n).

# 6.

# To send a "signed" message M, I first apply my private key to M to produce
# PRIV(M). Everyone in the world can recover M from this, because my public key
# is available. But, assuming that RSA is secure, no one could have *produced*
# a PRIV(M) which decrypts to a reasonable message other than me.

MakeRSASignature := proc(M, sender_key, recipient_key)
    local sender_n, sender_priv, recip_n, recip_pub, signed_M, secret_signed_M:
    sender_n := sender_key[1]:
    sender_priv := sender_key[3]:

    recip_n := recipient_key[1]:
    recip_pub := recipient_key[2]:

    signed_M := M&^sender_priv mod sender_n:
    print(`signed M`, signed_M);
    secret_signed_M := signed_M&^recip_pub mod recip_n:
end:

read "C3.txt":

sender_key := MakeRSAkey(100):
recip_key := MakeRSAkey(150):

message := 1766:

secret := MakeRSASignature(message, sender_key, recip_key):

# Only the recipient can do this!
decoded := secret&^recip_key[3] mod recip_key[1];

# Only the sender could have this result in the correct message!
checked := decoded&^sender_key[2] mod sender_key[1]:

message;
checked;