GIER/Primtals jagt 2013

Fra DDHFwiki
Spring til navigation Spring til søgning

Lidt ligesom at køre veteranbiler, skal gamle maskiner have noget at rive i. Her beskriver jeg en primtalsjagt.

2013-Aug-08: GIER legestue: Primtal 1

Date: Thu, 8 Aug 2013 16:24:31 +0200
From: Thorkil Naur
To: torsdagsaktive

Kære TA,

Nu står vores GIER jo til almindelig skue og brug, så jeg tænkte, at så skulle den også have noget at rive i. Og eftersom intet mere sundt eller nyttigt kan tænkes end at jagte primtal, besluttede jeg, at vores GIER skulle have fornøjelsen af at finde et 100-cifret primtal. Det er planen, i hvert fald, vi er der ikke endnu, men ting og sager er sket undervejs, som måske kan interessere.

Metoden er af Pocklington (Pocklington, H. C. "The Determination of the Prime or Composite Nature of Large Numbers by Fermat's Theorem." Proc. Cambridge Phil. Soc. 18, 29-30, 1914/16) og kan benyttes til at indse, at et tal N er et primtal, hvis tilstrækkelig store primfaktorer af N-1 er kendte. Og beregningerne sker med Jørgen Kjærs DEMON-5 som udgangspunkt, passende justeret og udvidet.

Når man bruger Pocklington og tilsvarende metoder i faktorisering, er situationen typisk den, at man står med et givet N og vil undersøge, om N er et primtal. For endeligt at afgøre sagen med Pocklington, må man som nævnt have tilstrækkeligt med faktorer af N-1. Og det er jo ikke sikkert, at disse faktorer sådan lige lader sig skaffe til veje.

Anderledes hvis man bare vil finde et stort primtal: Lad os sige, at vi allerede har et primtal, B, i hånden. Så kan vi danne N = A*B+1 med håb om, at N er et primtal. Her giver Pocklington en effektiv måde at afgøre sagen, hvis blot A <= B, fordi B herved bliver en tilstrækkelig stor primfaktor i N-1 = A*B. Tanken er derfor, med et primtal, B, i hånden, at gennemsøge tallene B*B+1, (B-1)*B+1, (B-2)*B+1 og så fremdeles, og så bruge Pocklington til at finde det første primtal i denne serie af tal. Herved vil vi finde et primtal, der er kun lidt mindre end kvadratet af B, altså med omtrent dobbelt så mange cifre som B.

Eller vil vi? Tjah, der er jo ingen garantier, men matematikken fortæller jo, at primtalstætheden i nærheden af tallet N er 1/ln(N), så man må forvente, at ca. hvert ln(10^100) = 100*ln(10) ca. = 230'te tal er et primtal, og det er jo et absolut overkommeligt antal at gennemsøge.

Vi kan således forvente at finde et 100-cifret primtal, N = A*B+1, hvis vi allerede har et 50-cifret primtal, B, i hånden. Og B finder vi på tilsvarende måde, med udgangspunkt i et 25-cifret primtal. Og så fremdeles. Indtil vi når ned på jorden og må finde et primtal uden Pocklington, det kunne f.eks. være omtrent 6-cifret, så skal vi blot prøve divisorer til omkring 1000. (Alt dette forstået bagfra, naturligvis: Vi starter med et 6-cifret primtal, går videre med et 12-cifret, 25-cifret, 50-cifret og endelig et 100-cifret.)

Pocklington i teorien

I detaljer lyder det specielle tilfælde af Pocklington, som vi bruger, på denne måde: Vi har et tal N med N-1 = A*B, hvor B er et primtal. Vi skal nu finde et tal, d, f.eks. d = 2, hvor disse betingelser er opfyldt:

  1. d^(N-1) mod N = 1
  2. Største fælles divisor (sfd) af d^A-1 = d^((N-1)/B)-1 og N er lig med 1.

Og kan vi finde sådan et d, siger Pocklington, at alle primfaktorer i N er på formen k*B+1 for et eller andet k. Tilføjer vi betingelsen, at A <= B, må N være et primtal: For hvis N var sammensat, skulle N jo have mindst to primfaktorer på formen k*B+1, så N var mindst lig (B+1)*(B+1) = B^2+2*B+1. Men det duer jo ikke, for N = A*B+1, hvor A <= B, så den eneste mulighed er, at N kun har en enkelt primfaktor og derfor er et primtal.

Et eksempel: Lad B være primtallet 101, A = 96 og N = A*B+1 = 9697. Vi prøver d = 2 og beregner i første omgang: 

 d^(N-1) mod N = 2^9696 mod N = 1

Nu er det jo ikke en helt beskeden beregning: 2^9696 er et tal med 2919 cifre, som jo skal frembringes med, umiddelbart, 9695 multiplikationer. Men, som det måske er velkendt, kan udtryk med store eksponenter beregnes med et antal multiplikationer, der er proportionalt med logaritmen til eksponenten, f.eks. udtrykt rekursivt således: 

 d^0 = 1 for d <> 0, 0 for d = 0
 d^(2*m) = (d^m)^2, m > 0
 d^(2*m+1) = (d^m)^2 * d, m >= 0

Dertil kommer, at eftersom slutresultatet skal reduceres modulo N, kan alle mellemresultaterne også reduceres modulo N, så tallenes størrelse bliver behersket.

Tilbage til eksemplet: d^(N-1) mod N = 1 er jo den kendte Fermat-betingelse, der er opfyldt af alle primtal N, der ikke går op i d. Og selvom Fermat-betingelsen ikke i sig selv sikrer, at N er et primtal, er sådanne undtagelser sjældne. Og vores N = 9697 opfylder altså Fermat-betingelsen. Havde betingelsen ikke været opfyldt, havde vi kunnet konkludere, at 9697 var sammensat.

Næste skridt er Pocklington-betingelsen (2) sfd(d^A-1,N) = sfd(2^96-1,9697) = 1. Her skal vi igen beregne et udtryk med en stor eksponent (2^96), og fordi resultatet blot skal videre i en beregning af største fælles divisor med N = 9697, kan vi som tidligere reducere både slutresultatet og alle mellemresultaterne modulo N.

Vi finder: 

 d^A mod N = 2^96 mod 9697 = 4868
 sfd(d^A-1,N) = sfd(4868-1,9697) = 1 med Euclids algoritme.

De to betingelser (1) og (2) er således opfyldte og vi kan konkludere, at alle N's primfaktorer har formen k*B+1 = k*101+1. Og igen, hvis der var to sådanne faktorer, var den mindste mulighed 102*102 = 10404, hvilket er større end N. Så der kan kun være en enkelt sådan primfaktor og N = 9697 må derfor være et primtal.

I praksis arrangeres beregningen lidt anderledes: I forløbet skal vi jo både beregne d^A mod N og d^(N-1) mod N, men eftersom N-1 = A*B, er 

 d^(N-1) mod N = d^(A*B) mod N,

så hvis vi sætter 

 X = d^A mod N,

kan vi beregne Y = d^(N-1) mod N ved 

 Y = X^B mod N

og spare et betydeligt antal multiplikationer modulo N.

Pocklington i praksis

Pocklington testen for N = A*B+1, hvor B er primtal og A <= B, kommer således til at se sådan ud:

  1. Sæt d = 2.
  2. Beregn X = d^A mod N.
  3. Hvis X = 1, er også X^B mod N = d^(N-1) mod N = 1, så Fermat-betingelsen er opfyldt og vi kan ikke udelukke, at N er primtal. Men sfd(X-1,N) = sfd(0,N) = N, så den anden Pocklington-betingelse er ikke opfyldt. Vi må prøve med et nyt d, sætter d = d+1 og går tilbage til step 2.
  4. Hvis derimod X er forskellig fra 1, går vi videre og beregner Y = X^B mod N, som jo er d^(N-1) mod N.
  5. Hvis Y er forskellig fra 1 er vi færdige og kan konkludere, at N er sammensat.
  6. Hvis Y = 1, mangler vi blot at beregne Z = sfd(X-1,N). Hvis Z = 1, er N primtal i henhold til Pocklington og A <= B. Hvis Z er forskellig fra 1, har vi fundet en faktor i N og denne faktor kan ikke være N selv, eftersom 1 < X < N. Så N er sammensat.

Som det ses, er denne beregning ikke som udgangspunkt sikret nogensinde at slutte, eftersom vi kan risikere, at X = 1 for alle d = 2, 3, ... I praksis vil vi oftest finde, med d = 2, at N er primtal eller konkludere, at N er sammensat, fordi Y er forskellig fra 1.

DEMON-5 justeret

Beregningerne sker med Jørgen Kjær's DEMON-5 som udgangspunkt. De grundlæggende strukturer er bevaret, ikke mindst repræsentationen af et stort tal som et integer array med fast længde, hvor hver celle indeholder 10 decimale cifre af det store tal. Samt ideen med at holde meget store tal (mere end 400 decimale cifre) på tromlen.

Proceduren LONG MULT til multiplikation af store tal er erstattet af LONG MULT 2, der især for meget store tal sparer en del tromletransporter. Og den vigtige LONG DIVIDE til division af et stort tal med et andet stort tal, inklusive beregning af resten, er tilføjet.

Dernæst er de mere sammensatte procedurer som EXP, PI TO osv., der ikke benyttes i primtalsjagten, fjernet af pladshensyn.

Primtalsjagten omfatter disse nye procedurer:

POWER
Beregner A^B, eventuelt modulo C.
GCD 1
Beregner største fælles divisor ("Greatest Common Divisor") af to store tal med Euclids algoritme.
TRIAL DIVIDE
Testdividerer et stort tal med små divisorer, evt. til en grænse.
IS PRIME 1
Tester om et stort tal er primtal med testdivision.
PRIME NEXT 1
Finder næste primtal over en grænse, testet med IS PRIME 1 (testdivision).
IS PRIME 2
Tester om et tal A*B+1 er primtal med Pocklington, hvor B skal være et primtal og A <= B. Inden Pocklington sættes ind for alvor, afprøves små faktorer op til en givet grænse for om muligt at afgøre sagen i en fart.
PRIME NEXT 2
Finder næste primtal blandt (A-1)*B+1, (A-2)*B+1, ..., hvor B skal være et primtal og A <= B, testet med IS PRIME 2 (Pocklington).

Af mere administrativ karakter er tilføjet demon5LOAD og demon5READ, der kan benyttes til at indlæse DEMON-5 registre (A, B, C) via strimmellæseren i det format, der skrives af DEMON-5's PRINT. Samtidig er DEMON-5 ordre-dekoderen tilføjet ordrer til at udskrive alle tre registre. Tilsammen giver disse tilføjelser en rimelig måde at "dumpe" registrene med udskrift på strimmelhuller (aktiveret via krydspanelet) med henblik på senere genindlæsning og fortsættelse af en beregning.

Endelig er DEMON-5 ordre-dekoderen omskrevet, så ordrelisten kan udskrives direkte ud fra den dynamisk bestemte ordrerækkefølge. Det er således ikke længere nødvendigt at vedligeholde den redundante writetext af alle ordrerne, som blev brugt tidligere. Dette er bekvemt, når man arbejder med DEMON-5, ikke mindst når man ofte fjerner og tilføjer ordrer. Til gengæld er det typisk nødvendigt at have en udskrift af den aktuelt gældende ordreliste foran sig ved brug af programmet. Listen kan jo stadig udskrives på skrivemaskinen af programmet selv, men bliver det kun på opfordring.

Eksempel på kørsel

Her er et eksempel på en primtalsjagt, hvor først 101 findes frem med PRIME NEXT 1 og derefter 9697 med PRIME NEXT 2: 

 2013-Jul-03 09.41 / TN
 Select language: d: danish, e: english, f: french, g: german.: 
 Dansk
 PROGRAM DEMON-5. Beregning af store tal. Programmet simulerer en maskine med
 3 registre, A, B og C, som har D decimaler og E cifre før kommaet.
 Opgiv antal decimaler, D. -1 er stop:  0
 Og antallet af heltalscifre, E:  210
 No (-1 opremser ordrer):  -1
  0:  stop
  1:  A := r;
  2:  LOAD(A);
  3:  READ(A, B, or C);
  4:  write(A);
  5:  write(B);
  6:  write(C);
  7:  B := A;
  8:  C := A;
  9:  A := B;
  10:  C := B;
  11:  A := C;
  12:  B := C;
  13:  A := A + B;
  14:  A := A - B;
  15:  A := A*r;
  16:  A := A/r;
  17:  C := A<*2>B;
  18:  (C,A) := (A_:B,A _m_o_d B);
  19:  (C,A) := (A_:B,A _m_o_d B) with r decimals;
  20:  A := A |& B (_m_o_d C, if non-zero);
  21:  A := IS PRIME 1(B);
  22:  A := PRIME NEXT 1(A);
  23:  A := IS PRIME 2(B*C(prime)+1,dmax);
  24:  C := PRIME NEXT 2([A-1,A-2,..]*B(prime)+1,dmax);
  25:  A := GCD 1(A,B);
 No:  1 A := r;r := 1.000000  ' 2
 tt 103+18 -> 121
 No: 22 A := PRIME NEXT 1(A);
 TRIAL DIVIDE: A =
         101 
  Prime, divided to 11
 tt 142+84 -> 226
 No:  4 write(A);
         101 
 tt 236+9 -> 245
 No:  7 B := A;
 tt 256+8 -> 264
 No: 24 C := PRIME NEXT 2([A-1,A-2,..]*B(prime)+1,dmax);r :=        10
 IS PRIME 2: N =
        9899 
   has no factors _< 7
   Computing 2 |& A _m_o_d N for A =
          98 
   Result 1: X =
        6488 
   X |& B _m_o_d N for B =
         101 
   Result 2: Y =
          81 
   Y |= 1, so N composite.
 IS PRIME 2: N =
        9697 
   has no factors _< 7
   Computing 2 |& A _m_o_d N for A =
          96 
   Result 1: X =
        4868 
   X |& B _m_o_d N for B =
         101 
   Result 2: Y =
           1 
   Y = 1, so verify Z = GCD( X-1, N ) =
           1 
   Z = 1, so N is prime
 tt 288+3468 -> 3756
 No:  6 write(C);
        9697 
 tt 3768+10 -> 3778
 No:  0 stop
 Opgiv antal decimaler, D. -1 er stop: 
 image
 #029p2    7.7.10  17#062

Værdien af ga4 standard variablen tracks transferred før og efter udførslen af hver ordre udskrives (linierne med "tt ...+... -> ...") som et groft udtryk for køretiden. Informationer om beregningerne udskrives, som vist, hvis kbon. PRIME NEXT 2 kan afbrydes efter hver test med kaon og vil efterlade A og B klar til videre jagt med fornyet brug af ordre 24. Det er her, jeg forestiller mig, at man ville kunne "dumpe" registrene på strimmel og fortsætte jagten ved en senere lejlighed.

C preprocessor

Af bekvemmelighedsgrunde bruger jeg C preprocessoren på DEMON-5 i ASCII format til at producere den GIER-klare udgave: En del gentagelser af logikken til at administrere meget lange tal på tromlen skrev jeg i første omgang som Algol procedurer, men det stjal en betragtelig mængde køretid, så i stedet benytter jeg C preprocessor makroer til at generere denne Algol kode inline. I forbindelse med den aktuelle primtalsjagt er jeg også begyndt at bruge C preprocessor #include's.

At producere en hulstrimmel

For at køre programmer på den rigtige GIER, må vi have dem på hulstrimmel. Heldigvis har vi mulighed for at udskrive hulstrimler fra moderne apparatur, som beskrevet på http://datamuseum.dk/wiki/GIER/Testsoftware#Hulning_af_strimler_til_GIER. Inden jeg kom i gang, advarede Mogens om, at der måske var vanskeligheder med den beskrevne mekanisme, at strimlerne måske ikke blev helt som ønskede. Så han anbefalede at bruge sumkoder på strimlerne, for derved at gardere mod fejl i skrivninger og anden kommunikation af data på strimlerne.

GIER/Hjælp 3 sumkoder er en konvention for skrivning, fortolkning og kontrol af sumkoder på strimler, som i hvert fald ga4, slip og ikke mindst edit er enige om: Se A Manual of Gier Algol 4 (http://datamuseum.dk/w/images/f/fa/Ga4man.pdf) og A Manual of Help 3 (http://datamuseum.dk/w/images/1/15/Help3man.pdf) for flere detaljer om disse ting. Essensen er, at vi tidligt på vores strimmel skal introducere en CLEAR CODE (for at nulstille tegn-summen og aktivere summeringen) og dernæst til sidst (og måske også med jævne mellemrum) inkludere en SUM CODE, der så efterfølges af et sum-tegn, passende beregnet som summen af de foregående tegn. Både ga4 og slip kontrollerer sådanne summer, ligeså edit, men edit erstatter også et forkert sumtegn med det rigtige, så et gennemløb af data med edit, eventuelt helt uden at foretage egentlige rettelser, kan bruges til at generere en udgave af data med de rigtige summer.

Basis for generering af en rigtig strimmel til GIER fremkommer således ved:

  1. I tilfældet med DEMON-5 til primtalsjagt, kørsel af C preprocessoren til behandling af makroer og #include direktiver.
  2. Outputtet skal så gennem edit, kørt via GIER simulatoren, for at generere de rigtige sum-tegn. Her er det vigtigt, at slutresultatet er en .flx-fil og ikke en .asc-fil, så den er klar til at blive kopieret til strimmel-hulleren.

Med disse ting på plads, gik jeg så i gang med at hulle. Nu er DEMON-5, også med de angivne udeladelser, ikke et helt lille program, så hulningen tog tid og brugte en del strimmel. Og processen gik da også i stå på et tidspunkt, hvor der ikke var mere blank strimmel tilbage i hulle-maskinen. Min umiddelbare plan, for ikke at spilde den strimmel, der var blevet hullet, var at afslutte den partielle strimmel med noget blank strimmel, med henblik på at påklistre den manglende del. Det er sådan, at jeg aldrig i min tilværelse som hulstrimmel-jonglør er blevet god til at klistre strimler sammen, uden noget blank strimmel at løbe på. Mens sammenklistring af strimler med passende blanke passager er lykkedes i mange tilfælde.

Desværre lykkedes det ikke at inkludere blank strimmel i enden af den partielle strimmel, der i stedet blev hullet ud til sin yderste ende. Mogens foreslog så, at jeg i stedet for at satse på sammenklistring, simpelthen bibeholdt det samlede program på separate strimler, som blot skulle indlæses efter tur.

Som sagt så gjort: For let at identificere, hvor langt den partielle strimmel var nået, lagde jeg de sidste par meter i GIER strimmellæseren og lod edit udskrive på skrivemaskinen. Derefter var det en smal sag at klippe den over et passende sted og generere en ny strimmel med resten af programmet. Dvs. af årsager, som ikke er fundet, stoppede udskriften af anden del, inden udskriften var afsluttet, så det var nødvendigt at udskrive endnu en strimmel med resten af programmet. Denne tredie del viste sig at være komplet, så hele programmet foreligger således som tre hulstrimler, der blot skal indlæses efter hinanden.

Første afprøvning

Så går det løs med den rigtige GIER: Med dens begrænsede kapacitet, er vi overladt til transient compileren, dvs. den på hulstrimler, som beskrevet på http://datamuseum.dk/wiki/GIER/Testsoftware#GIER_ALGOL_4. Nu ligger mit program jo på tre hulstrimer, og under første forsøg med at oversætte, meldte oversætteren 2 gange "character" (ugyldigt tegn) et stykke tid efter, at første strimmel var læst og mens jeg var i færd med at lægge strimmel nummer 2 i læseren. Denne oversættelse gik derfor galt. I næste forsøg lykkedes det imidlertid og programmet kunne køre. Det oversatte program fylder omkring 190 tromlespor, free har lidt over 200 ledige spor, så det går lige.

Næste oplevelse var så, at i hvert fald tegnet '1' og vistnok også '-' ikke fungerer rigtigt: De kan tastes og skrives af skrivemaskinen, men tilsyneladende får GIER ikke fat i tegnet. Det betyder så i første omgang, at taster jeg "210" som det ønskede antal cifre før kommaet, får jeg kun 20, og det duer jo ikke. I stedet brugte jeg 220. (Jeg kunne vist også bruge f.eks. 209 og lade DEMON-5 runde op til nærmeste 10.)

Og så var dette første eksperiment jo allerede ved at kæntre inden søsætningen, for det er jo ordre 

  1:  A := r;

man bruger til at få data ind i DEMON-5's registre og kan jeg ikke taste 1, kan jeg jo ikke aktivere denne indlæsning. Heldigvis kunne jeg bruge en ordre 

  2:  LOAD(A);

der læser register A ind fra skrivemaskinen i DEMON-5 PRINT format. Denne test ordre skal egentlig ikke bruges, men jeg havde ikke fået fjernet den endnu.

Men på kort sigt, indtil '1'-tegnet måske kommer til at virke, må jeg jo finde en mulighed for at aktivere de forskellige ordrer, der indeholder tegnet '1'. Deriblandt er jo f.eks. begge de ordrer, der kopierer registeret C til andre registre og det havde jeg jo tænkt mig at bruge i min primtalsjagt. Min umiddelbare plan på dette felt er at justere programmet til at reducere den indtastede ordre modulo 222, så f.eks. 12: B := C; alternativt kunne indtastes som 12+222=234.

Tegnet '-' er mindre kritisk: "-1" angives ganske vist som den ordre, der udskriver ordrelisten, men i virkeligheden kan enhver ugyldig ordre bruges her.

Men videre, kom der faktisk en lille smule gang i en egentlig primtalsjagt: Men simple krumspring (fordi '1' ikke virker) fik jeg indlæst et passende 6-cifret tal i A, fundet det nærmest større primtal med ordre 22 (som blot dividerer), kopieret resultatet til register B og fundet et afledt 12-cifret (eller det var måske 13-cifret, jeg har ikke mine notater ved hånden i skrivende stund) primtal med ordre 24, der bruger Pocklington.

For at komme videre, skal jeg nu gøre det umulige, nemlig kopiere C registeret (der indeholder det 12-13-cifrede primtal) til register A og B. Derved kommer vi så til næste problemstilling, for i stedet for (som jeg selvfølgelig kunne, men nødigt ville) at indtaste tallet, forsøgte jeg at gå via strimmelhulleren, også for at få gang i andre aspekter af jagten. Planen var at udskrive C på strimmel, så rette strimlen til at ramme register A i stedet for register C (det ønskede register fremgår af data på strimlen) og så indlæse til register A ved hjælp af den modificerede strimmel.

(Et andet problem med manuel overførsel af tal er, at det ofte er vanskeligt at læse skrivemaskinens output. Undervejs greb jeg til at hulle data på strimmel, hvor de decimale cifre kunne aflæses med større sikkerhed end på skrivemaskinen.)

Desværre udviklede strimmelhulleren under forløbet en sygdom: Til at starte med virkede den i flere tilfælde, men efterhånden skete det ofte, at udskrift, der skulle komme til hulleren, fik GIER til at hænge. Frakobling af huller via krydspanelet samt tryk på Reset+Start fik programmet til at fortsætte, men også blot dette at aktivere generering af tom strimmel fik GIER til at hænge under omstændigheder, som jeg ikke fik udforsket fuldstændig. Den eneste afhjælpning af fejl, som jeg forsøgte, og den havde kun begrænset effekt, var at slukke og tænde for strømmen til hulleren, ved at hive ledningen ud og sætte den tilbage igen.

Men jeg fik dog til sidst tastet det 12-13-cifrede primtal ind (det indeholdt heldigvis ingen 1-taller) og brugt det til at generere et 25-26-cifret primtal).

Slut for denne gang

De seneste aktiviteter foregik torsdag den 25. juli 2013. Jeg efterlod papirsporet i skrivemaskinen, hvor man måske stadig kan se resultaterne. Spændende fortsættelser venter ved senere lejligheder.

Hilsner Thorkil

2013-Sep-26: 41371368561830689276952510478484206292283710462079

Denne torsdag kom der så gang i jagten for alvor: Programmet havde overlevet i GIER siden sidste torsdag, så det var bare at tænde, trykke Reset + HP og skrive run<, så var vi kørende.

For et par torsdage siden, den 15. august 2013, havde jeg fået produceret et alternativ til de to strimler, der udgør den sidste del af programmet, med den planlagte rettelse, hvor det indtastede ordrenummer reduceres modulo 222. Som nævnt har dette krumspring til hensigt at tillade alternative ordrenumre, sådan at man aldrig får behov for at indtaste 1-taller, der jo ikke for tiden modtages af GIER.

Torsdagen efter havde en gæst, Povl Kvols, foreslået en startværdi for jagten: Vi skal jo ende med et 100-cifret primtal, som jo fremkommer ved nogenlunde kvadratet på et 50-cifret primtal, der igen er omtrent kvadratet på et 25-cifret primtal, som opstår ved omtrent at kvadrere et 12-13-cifret primtal, som endelig beregnes som omtrent kvadratet på et 6-7-cifret primtal, som vi altså skal vælge en omtrentlig størrelse for. Og vi kan jo lige så godt vælge det 100-cifrede tal og så tage 4 kvadratrødder for at finde udgangspunktet for den første jagt. Og her var det så, at Povl forslog 69!*10. Det er jo velkendt blandt, nok især ældre, regnemaskinebrugere, at den største faktultet, der kan beregnes på en maskine med 2 decimale cifre i eksponenten er 69!, der er omtrent 1.71122452428142e+98, et tal på 99 cifre. 70! har jo 101 cifre, så Povl valgte at sige, jamen så starter vi med 69!*10, der har netop 100 cifre. Efter 4 kvadratrødder, kommer vi ned til 1592508.23661706, hvor jagten på det første primtal skal starte.

Allerede her er vi jo i vanskeligheder, for tallet indeholder et 1-ciffer, som vi jo stadig ikke kan taste, men 1592508 = 6370032/4, så jeg tastede bare 6370032 ind og lod DEMON-5 dividere med 4.

Efter kort tid, fandt programmet det næste større primtal, B1 = 1592533, hvor primtalsbeviset blot bestod i at dividere med faktorer op til kvadratroden. Med dette primtal i hånden, skal vi så gennemsøge A1*B1+1 for A1 = B1-1, B1-2, ... og teste, først ved division med små faktorer og dernæst med Pocklington, indtil vi finder et primtal.

Efter kort tid fandt programmet, at B2 = A1*B1 + 1 = 1592526*B1 + 1 = 2536150208359 var et primtal. Og endnu to primtal blev fundet: 

 B3 = A2*B2 + 1 = 2536150208310*B2 + 1 = 6432057879235127753263291
 B4 = A3*B3 + 1 = 6432057879235127753263258*B3 + 1 = 41371368561830689276952510478484206292283710462079

Testen for det 50-cifrede tal B4 tog omtrent en halv time, fordelt på denne måde: 

 Division med faktorer til 5000         5 minutter
 X = 2^A3 mod B4                       10 minutter
 Y = X^B3 mod B4                       10 minutter
 Største fælles divisor af X-1 og B4    5 minutter

Et par andre tal blev testet på samme måde og forkastet som sammensatte, men de fleste af de B3-A3 = 33 kandidater blev forkastet ved at finde små faktorer. ln(B4) = 114.2, så vi var lidt heldige, at kunne nå denne jagt på en enkelt torsdag.

For at komme videre ved en senere lejlighed, skal vi jo kunne indlæse B4 og strimmelhulleren har det ikke så godt. I stedet brugte jeg Poul-Henning og Mogens' mekanisme, der tillader styring af GIER skrivemaskinen fra en PC via USB porten. Som en sideeffekt, kan udskrevne tegn skrives på en fil på PC'en og dermed kan data fra GIER opfanges. Mekanismen er nu ikke ganske stabil og det krævede et par forsøg at få overført B4, her gengivet i den klassiske DEMON-5 stil: 

  6 write(C);
 41371 36856 18306 89276 95251 04784 84206 29228 37104 62079 
 tt 333064+10 -> 333075

Tanken er at hulle en strimmel med dette tal fra PC'en og så bruge den nye READ ordre til at indlæse tallet via RC2000.

Og vi skal jo videre og finde B5. Et gæt på, hvor lang tid det tager at teste sådan et 100-cifret tal kunne være dette:

  • Divisionerne kan jo reguleres ved at justere grænsen, men mon ikke det er fornuftigt at bruge en 10-15 minutter på det, i et forsøg på at spare brugen af Pocklington.
  • For potensopløftningerne, skal der udføres omtrent dobbelt så mange operationer for 100-cifrede tal som for 50-cifrede. Dertil kommer, at hver operation, en multiplikation efterfulgt af en division, vil tage noget, der ligner 4 gange så lang tid, med disse simple n^2 skole-metoder, som DEMON-5 bruger. Så i alt en faktor 8.
  • Største fælles divisor skal nok også have en faktor 8.

I alt 15 + 8*(10 + 10 + 5) minutter = 3 timer og 35 minutter. Så mere end en enkelt af de dyre kandidater kan det nok ikke blive til på en torsdag.

2013-Okt-03: 1711590136678832715194635872474013002927665545820488747130228034523451878255063681857184340779739581? Nej.

Vi skal i gang med de 100-cifrede primtalskandidater og jeg fik produceret en strimmel med sidste torsdags primtal B4 = 41371368561830689276952510478484206292283710462079 på hulstrimmel i det format, der kan læses af READ ordren. Igen var GIER klar til kamp efter Reset+HP+run<.

I løbet af aftenen blev i alt 59 tal fundet sammensatte og det drejer sig jo om de 59 tal (B4-1)*B4+1, (B4-2)*B4+1, ... (B4-59)*B4+1. De fleste af disse tal har små divisorer, der findes i den indledende fase med test-division (blandt andet er jo hvert andet tal lige) og jeg havde sat grænsen for divisorerne til 10000.

To af tallene, (B4-35)*B4+1 og (B4-59)*B4+1, det sidste er tallet 1711590136678832715194635872474013002927665545820488747130228034523451878255063681857184340779739581 fra overskriften, havde ikke divisorer under 10000. Med ca. 167 divisioner i minuttet, tog det omkring 20 minutter at dække intervallet til 10000: Alle divisorer på formen 6*k+1 og 6*k+5 afprøves, og dem er der omtrent 10000/3 = 3333 af.

(B4-35)*B4+1 og (B4-59)*B4+1 måtte således igennem den fulde Pocklington test, men begge blev fundet sammensatte. Testen bruger ca. 42 minutter om hver af de to potensopløftninger, der skal beregnes for en enkelt kandidat, så i alt ca. 104 minutter for den fulde test, når den melder "sammensat". Og det er jo noget mindre end vi gættede på tidligere. Formentlig skyldes dette, at der trods alt er væsentlige elementer i beregningerne, der kun tager tid proportionalt med ciffer-antallet og ikke med kvadratet på ciffer-antallet.

Separat har jeg fundet, at den mindste primfaktor i (B4-35)*B4+1 er 103613, mens det ikke er lykkedes (endnu) at finde en primfaktor i (B4-59)*B4+1. Vi kunne således være sluppet for den ene Pocklington ved at sætte divisorgrænsen til 104000, men derved ville divisionerne have taget 200 minutter i stedet for de 20, så det ville ikke have været en god forretning.

På den anden side har vi (B4-29)*B4+1, hvis mindste primfaktor, 2309, var størst blandt de 59 kandidater, hvor vi fandt faktorer. Ved at reducere divisorgrænsen til 2400, kunne vi således have nedsat tiden til at dividere med omkring 15 minutter for hver af de kandidater, der skulle igennem Pocklington.

Men overalt set synes jeg, at grænsen på 10000 for divisorer ser fornuftig ud.

Med udskrifterne 

  4 write(A);
 41371 36856 18306 89276 95251 04784 84206 29228 37104 62020
 tt 371192+12 -> 371205
 No:  5 write(B);
 41371 36856 18306 89276 95251 04784 84206 29228 37104 62079
 tt 371213+9 -> 371223

i hus er vi klar til næste runde ved en senere lejlighed.

Programkode

Programkoden ligger i Fil:Demon5sp1 20130815 1428.asc. Denne .asc udgave skal preprocesses med C preprocessoren og forudsætter, at visse #include'de filer er tilgængelige. For i første omgang at rode mig uden om at gøre disse #include-filer tilgængelige direkte, er i stedet den preprocessede, a2flx'ede og endelig, med edit opsummerede, fil lagt i Fil:Demon5sp1 20130815 1515 summed.flx. Denne .flx udgave af programmet er umiddelbart klar til at blive hældt gennem GIER Algol 4 oversætteren. Og hvis man flx2a'er den, kan jo de #include'de dele læses, mens, imidlertid, andre dele, der er resultatet af makro-ekspandering, bliver mindre læselige.

2013-Okt-10: Heller ikke (B4-60)*B4+1, (B4-61)*B4+1, ... (B4-71)*B4+1

Denne gang fik vi testet de 12 tal (B4-60)*B4+1, (B4-61)*B4+1, ... (B4-71)*B4+1, men alle var sammensatte. (B4-69)*B4+1 og (B4-71)*B4+1 måtte igennem Pocklington, de øvrige havde små primfaktorer, med 7 som den største.

Jeg har kontrolleret delresultaterne i Pocklington-beregningerne (X og Y) og de ser fine ud: GIER kan godt regne.

2013-Okt-17: (B4-72)*B4+1, (B4-73)*B4+1, ... (B4-95)*B4+1 er alle sammensatte

I samme stil som sidst, fik vi ryddet yderligere en bunke kandidater til side som sammensatte: (B4-89)*B4+1 og (B4-95)*B4+1 måtte igennem Pocklington.

Separate beregninger viser, at 

 (B4-89)*B4+1 = 443897087 * 23197238735767 * 166219200091901796074411640417851251206992863322487482773466962726386262909859
 (B5-95)*B4+1 = 492438485052764251 * 59948738022760755211 * 57978602540078325609253816919703795137834885219083284627247417

så det er jo lige til at regne efter.

2013-Okt-25: (B4-96)*B4+1, (B4-97)*B4+1, ... (B4-137)*B4+1 er sammensatte

Og ligesom tidligere, nåede vi to Pocklington tests på en aften: (B4-125)*B4+1 og (B4-137)*B4+1. Største faktor i de øvrige kandidater var 457, der gik op i (B4-117)*B4+1. Separat finder vi: 

 (B4-125)*B4+1 = 755634856579 * 23734168365012654583567657 * 95436339419294918172072387202547996523702671634455842824240589
 (B4-137)*B4+1 = 323027637853 * 5298587291337978476694055226779199373844933521671192135169779140166969601660791777855623

2013-Okt-31: (B4-138)*B4+1, (B4-139)*B4+1, ... (B4-179)*B4+1 er sammensatte

Pocklington for (B4-167)*B4+1 og (B4-179)*B4+1.

2013-Nov-07: (B4-180)*B4+1, (B4-181)*B4+1, ... (B4-225)*B4+1 er sammensatte

Pocklington for (B4-197)*B4+1 og (B4-225)*B4+1.

2013-Nov-21:(B4-237)*B4+1 = 1711590136678832715194635872474013002927665545813124643526222171832154331389893493137157840317489519 er primtal

I dag lykkedes det at finde et 100-cifret primtal. Indledningsvis måtte kandidaten 

 (B4-231)*B4+1 = 608831899142811973 * 2811268823280479704149933565714671436078989289258953607311981994872651492303830741

afvises med Pocklington, men så viste beregningerne, at (B4-237)*B4+1 = 1711590136678832715194635872474013002927665545813124643526222171832154331389893493137157840317489519 er et primtal. Den afsluttende beregning af en GCD (største fælles divisor), som vi ikke har prøvet før med tal af denne størrelse, tog omtrent 5 minutter.

2014-Jan-09: 200 cifre?

De seneste torsdage har GIER været skidt tilpas, men i dag fik Ole skiftet en sikring, der var røget, og det lykkedes at få min justerede udgave af DEMON-5 oversat. Med B4 indlæst fra strimmel, fik jeg beregnet vores 100-cifrede primtal B5 = (B4-237)*B4+1 og fik startet en jagt på et primtal på formen (B5-k)*B5+1. Disse tal har ikke helt 200 cifre, kun 199, men vi må lade os nøje i denne omgang. Men under alle omstændigheder: GIER fandt små faktorer i (B5-1)*B5+1, ... (B5-8)*B5+1 og en faktor 1213 i (B5-9)*B5+1.

På et tidspunkt blev jeg klar over, at jeg som vanligt havde givet DEMON-5 210 cifre at arbejde med og det er ikke nok til en Pocklington med 199-cifrede primtalskandidater, så jeg afbrød jagten.

2014-Jan-18: Stak overløb

Den første 199-cifrede kandidat til Pocklington kvaltes i stak overløb, så der skal fedtes nogle bits af de lokale variable for at komme videre.

2014-Jan-24: Fejl i LONG DIVIDE 2

Under arbejdet med at spare lidt på variabelforbruget, fandt jeg en fejl i LONG DIVIDE 2: Den virker ikke altid, når tallene ligger på tromlen. De hidtidige beregninger har jo været med tal, der var små nok til at holde i kernelageret, så denne fejl er vi ikke rendt ind i tidligere.

Hentet fra "https://datamuseum.dk/w/index.php?title=GIER/Primtals_jagt_2013&oldid=8177"