Sequence of primes by trial division: Difference between revisions

Added solution for EDSAC.
m (→‎{{header|C sharp}}: tweaked to work properly with limits other than 100.)
(Added solution for EDSAC.)
Line 1:
 
 
{{task|Prime Numbers}}
 
Line 737 ⟶ 739:
1000000000163
*stopped*</lang>
 
=={{header|EDSAC order code}}==
To save time, the program omits multiples of 2 and 3, and tests 5, 7, 11, 13, ..., the increments being alternately 2 and 4. Division is done implicitly by using wheels. The program stores possible prime divisors in a list, whose size can be edited. If a prime is found, and the list isn't yet full, the prime is added to the list.
 
Initially, no wheels are active. When the number under test reaches the "milestone" 5^2, the wheel for 5 is activated, and the milestone is updated to 7^2. When that is reached, the wheel for 7 is activated, and the milestone is updated to 11^2; and so on.
 
In the program as posted, the list of possible divisors holds 30 primes, from 5 to 131. The program finds primes less than 131^2 = 17161, the largest being 17159. Assuming 650 orders per second, this would have taken the original EDSAC about an hour.
<lang edsac>
[List of primes by trial division, for Rosetta Code website.]
[EDSAC program, Initial Orders 2.]
[Division is done implicitly by the use of wheels.
One wheel for each possible prime divisor, up to an editable limit.]
 
T51K [G parameter: print subroutine, 54 locations.]
P56F [must be even address]
T47K [M parameter: main routine.]
P110F [must be even address]
 
[============================= M parameter ===============================]
E25KTM
GK
[35-bit values. First clear them completely.
This is done to ensure that the middle bit ("sandwich digit") is zero.]
T#ZPF T2#ZPF T4#ZPF T6#ZPF T8#ZPF
[Back to normal loading]
TZ
[0] PDPF [number under test; initially to 1, pre-inc'd to 5]
[2] P1FPF [increment, alternately 2 and 4]
[4] P12DPF ['milestone', square of prime; initially 25]
[6] PDPF [constant 1]
[8] P3FPF [constant 6]
[17-bit values]
[10] P30F [*EDIT HERE* Number of primes to store (in address field)]
[11] PF [flag < 0 if number is prime; 0 if factor is found]
[12] #F [figure shift]
[13] @F [carriage return]
[14] &F [line feed]
[15] K4096F [null char]
[16] A112@ [A order for list{0}]
[17] T112@ [T order for list{0}]
[18] AF [limit of A order for testing primes]
[19] AF [limit of A order for loading wheels]
[20] TF [limit of T order for storing primes]
[21] O1F [subtract from T order to make A order for previous address]
[22] W1F [add to T order to make U order for next address]
[23] OF [add to A order to make T order for same address]
[24] P2F [to inc an address by 2]
 
[Enter with acc = 0]
[25] O12@ [set teleprinter to figures]
[Set limits for list of trial prime divisors.
The list contains wheels and negatives of primes, thus:
wheel for 5; -5; wheel for 7; -7; wheel for 11; -11; etc]
A10@ [number of items in prime list]
LD [times 2 words per item (wheel + prime)]
A17@ [add T order for list{0}]
U20@ [store T order for exclusive end of wheels]
S21@ [make A order for inclusive end of primes]
T19@ [store it]
A16@ [load A order for start of lise]
U18@ [store as exclusive end of active wheels]
A2F [inc address, exclusive end of active primes]
T100@ [plant in code]
A17@ [load T order to store first wheel]
T89@ [plant in code]
[Main loop: update increment, alternately 2 and 4]
[Assume acc = 0;]
[38] A8#@ [load 6]
S2#@ [subtract incremet]
T2#@ [store new increment]
[First priority: keep the wheels turning]
A16@ [load order that loads first wheel]
U11@ [set prime flag (any negative value will do)]
[43] U49@ [plant order in code]
S18@ [more wheels to test?]
E66@ [if not, jump with acc = 0]
A18@ [restore after test]
A23@ [make order to store wheel]
T62@ [plant in code]
[49] AF [load wheel]
A2@ [apply current inc as 17-bit 2 or 4]
G62@ [if wheel still < 0, just store updated value]
T1F [wheel >= 0, save in 1F]
S1F [wheel = 0?]
G56@ [no, skip next order]
T11@ [yes, so prime flag := 0]
[56] TF [clear acc]
A49@ [make A order for negative of prime]
A2F
T61@ [plant in code]
A1F [load wheel again]
[61] AF [add negative of prime to set wheel < 0]
[62] TF [store updated wheel]
A49@ [on to next wheel]
A24@
G43@ [always jump, since A < 0]
[Update the number under test. Assume acc = 0.]
[66] A#@ [add incrememnt to number under test]
A2#@
U#@ [store new number]
[Test whether we've reached the "milestone", i.e. number = p^2.]
S4#@ [subtract milestone]
E94@ [if reached milestone, jump with acc = 0]
TF [clear acc]
A11@ [acc < 0 if number is prime, 0 if composite]
E38@ [if composite, loop for next number with acc = 0]
[Here when number is found to be prime.]
TF [clear acc]
A#@ [load number]
TD [copy number 0D for printing]
[77] A77@
GG [call print routine, clears acc]
O13@O14@ [print CR, LF]
[If list of primes isn't yet full, store the prime just found.
It's slightly more convenient to store the negative of the prime.
Also, the wheel is initialized to the negative of the prime.]
A89@ [load T order to store ]
S20@ [compare with end of list]
E38@ [if list is full, loop with acc = 0]
A20@ [restore acc after test]
A22@ [make U order for wheel + 1, i.e. for prime]
T88@ [plant in code]
S@ [load negative of latest prime]
[88] UF [store in list]
[89] TF [initialize wheel for this prime]
A89@ [inc address by 2 for next time]
A24@
T89@
E38@ [loop with acc = 0]
[Here when number tested equals the "milestone" p^2 (p prime).
We need to activate the wheel for the prime p,
and update the milestone to the next prime after p.]
[Assume acc = 0.]
[94] A100@ [load A order below]
S19@ [test against A order for end of list]
E110@ [if reached end of list, exit]
A19@ [restore acc after test]
A24@ [inc address in A order]
T100@ [plant in next order]
[100] AF [load negative of prime from list]
TF [to 0F]
HF [to mult reg]
VF [acc := square of prime scaled by 2^(-32)]
R1F [scale by 2^(-34) for 35-bit value]
T4#@ [update]
A18@ [start testing next prime wheel]
A24@
T18@
E38@ [loop with acc = 0]
[Here on exit from program]
[110] O15@ [print null to flush printer buffer]
[111] ZF [stop]
[Array of wheels and primes, immediately after program code]
[112]
 
[============================ G parameter ==================================]
[Modified library subroutine P7.]
[Prints signed integer; up to 10 digits, left-justified.]
[Input: 0D = integer,]
[54 locations. Load at even address. Workspace 4D.]
E25KTG
GKA3FT42@A49@T31@ADE10@T31@A48@T31@SDTDH44#@NDYFLDT4DS43@
TFH17@S17@A43@G23@UFS43@T1FV4DAFG50@SFLDUFXFOFFFSFL4FT4D
A49@T31@A1FA43@G20@XFP1024FP610D@524D!FO46@O26@XFSFL8FT4DE39@
 
[========================= M parameter again ===============================]
E25KTM
GK
E25Z [define entry point]
PF [acc = 0 on entry]
</lang>
{{out}}
<pre>
5
7
11
13
[...]
17047
17053
17077
17093
17099
17107
17117
17123
17137
17159
</pre>
 
=={{header|Eiffel}}==
113

edits