ProGramMA bij artiKel vaN maRtIjn LeIsInk Over EindCIJfers

Waarom is het getal $N= 823\ 252\ 653\ 353\ 834\ 464\ 227\ 313\ 169\ 615\ 015\ 138\ 644\ 385\ 446\ 477$ zo bijzonder?

Het antwoord is dat $N^3$ eindigt op een groot aantal drieën. Het programma hieronder berekent als je uitgaat van een reeks cijfers welk getal je moet kiezen om na tot de macht drie verheffen deze reeks cijfers te verkrijgen. Het programma is algemener. Je hoeft niet perse de macht $3$ gebruiken, je kunt een willekeurige macht $M$ kiezen. Belangrijk is wel dat $M$ niet $2$ of $5$ deelt. In dat geval kan het gebeuren dat er geen oplossing is, of dat er juist meerdere oplossingen zijn. Ga dat zelf na. Bedenk hoe je het programma kunt aanpassen dat het programma alle oplossingen genereert voor willekeurige machten. Het getal dat eindigt met de bijzondere cijferreeks noemen we $R$. In principe kan het programma ook worden uitgevoerd als $R$ even is of een vijfvoud. Echter de kans is groot (groter) dat er geen oplossing is. Dat laat ik over aan de lezer! Het programma bestaat uit twee stappen:

1. Je bepaalt het laatste cijfer (Dit gaat gewoon met trial and error)
2. Iteratief bepaal je steeds een extra cijfer. Hierna gaan we er van uit dat inmiddels $n$ cijfers al bekend zijn (in $A$, eigenlijk $A_n$). We schrijven $N=A+10^n\cdot B$. We maken hier gebruik van de formule: $R=(A+10^n\cdot B)^M = A^M + M\cdot A^{M-1}\cdot 10^n \cdot B + \mbox{rest}$ waarbij $\mbox{rest}$ een veelvoud is van $10^{2n}$. Hieruit volgt dat $B=(R-A^M)/(M\cdot A^{M-1}\cdot10^n)$. Van $B$ hebben we alleen het laatste cijfer nodig, en dat voegen we toe aan $A$.

We schrijven twee versies van hetzelfde programma. Bij de eerste versie proberen we bij elke stap de cijfers $0$ t/m $9$ uit. Bij de tweede versie berekenen we expliciet $B=(R-A^M)/(M\cdot A^{M-1}\cdot10^n)$. Het delen is niet simpelweg delen. $R-A^M$ bevat een factor $10^n$. Die factor kan uitgedeeld worden. Het restant moet worden gedeeld door $M\cdot A^{M-1}$. We gaan er van uit dat zowel $M$ als $A$ geen delers $2$ en $5$ bevatten. In dat geval kunnen we mod_inverse in sympy gebruiken.

We berekenen $B=\mbox{teller}/\mbox{noemer} \mod 10$, waarbij $\mbox{teller}=(R-A^M)/10^n$ en $\mbox{noemer}=M\cdot A^{M-1}$. Het laatste cijfer van $B$ voegen we (iteratief) toe aan $A$.

N = 823252653353834464227313169615015138644385446477
print(N**3)

 

## Eerste Programma
import math
Eindgetal = 123456789
Eindcijfer = Eindgetal % 10
Macht = 3
if math.gcd(Macht, 10) > 1:
    print("Het is raadzaam een andere macht te kiezen zonder factoren 2 en 5.")
if math.gcd(Eindcijfer, 10) > 1:
    print("Het is raadzaam een ander eindgetal te kiezen zonder factoren 2 en 5.")
Lengte = math.ceil(math.log(Eindgetal)/math.log(10))
print("Lengte Eindgetal:", Lengte)
print("Eindgetal:", Eindgetal)
print("Laatste cijfer:", Eindcijfer)
N = 0
# Trial and error: nagaan of 1 t/m 9 als laatste cijfer kan voorkomen
for i in range(9):
    NN = i+1
    NM = NN**Macht
    if NM % 10 == Eindcijfer:
        N = NN
if N == 0:
    print("Er is iets fout gegaan.")
else:
    print("Stap 1: N =", N, "=> N^", Macht, "=", N**Macht)
    for j in range(Lengte-1):
        for i in range(9):
            NN0 = i+1
            NN = N + 10**(j+1)*NN0
            NM = NN**Macht
            if NM % (10**(j+2)) == Eindgetal % (10**(j+2)):
                N = NN
        print("Stap", j+2, ": N =", N, "=> N^", Macht, "=", N**Macht)

 

## Tweede Programma
import math
import sympy
# Let op de namen van de variabelen zijn aangepast, 
# zodat ze in overeenstemming zijn met de inleiding
R = 123456789
R0 = R % 10
M = 3
if math.gcd(M, 10) > 1:
    print("Het is raadzaam een andere macht te kiezen zonder factoren 2 en 5.")
if math.gcd(R0, 10) > 1:
    print("Het is raadzaam een ander eindgetal te kiezen zonder factoren 2 en 5.")
Lengte = math.ceil(math.log(R)/math.log(10))
print("Lengte Eindgetal:", Lengte)
print("Eindgetal:", R)
print("Laatste cijfer:", R0)
A = 0
# Trial and error: nagaan of 1 t/m 9 als laatste cijfer kan voorkomen
for i in range(9):
    AA = i+1
    AM = AA**M
    if AM % 10 == R0:
        A = AA
if A == 0:
    print("Er is iets fout gegaan.")
else:
    print("Stap 1: N =", A, "=> N^", M, "=", A**M)
    for j in range(Lengte-1):
        teller = (R - A**M) // 10**(j+1)
        noemer = M * A**(M-1)
        inv_noemer = sympy.mod_inverse(noemer,10)
        B = teller*inv_noemer % 10
        A += B * 10**(j+1)
        print("Stap", j+2, ": N =", A, "=> N^", M, "=", A**M)