AKS test for primes: Difference between revisions
Content deleted Content added
Ada version |
|||
(42 intermediate revisions by 23 users not shown) | |||
Line 38:
* [http://www.youtube.com/watch?v=HvMSRWTE2mI Fool-Proof Test for Primes] - Numberphile (Video). The accuracy of this video is disputed -- at best it is an oversimplification.
<br><br>
=={{header|11l}}==
{{trans|Python}}
<syntaxhighlight lang="11l">F expand_x_1(p)
V ex = [BigInt(1)]
L(i) 0 .< p
ex.append(ex.last * -(p - i) I/ (i + 1))
R reversed(ex)
F aks_test(p)
I p < 2
R 0B
V ex = expand_x_1(p)
ex[0]++
R !any(ex[0 .< (len)-1].map(mult -> mult % @p != 0))
print(‘# p: (x-1)^p for small p’)
L(p) 12
print(‘#3: #.’.format(p, enumerate(expand_x_1(p)).map((n, e) -> ‘#.#.#.’.format(‘+’ * (e >= 0), e, I n {(‘x^#.’.format(n))} E ‘’)).join(‘ ’)))
print("\n# small primes using the aks test")
print((0..100).filter(p -> aks_test(p)))</syntaxhighlight>
{{out}}
<pre>
# p: (x-1)^p for small p
0: +1
1: -1 +1x^1
2: +1 -2x^1 +1x^2
3: -1 +3x^1 -3x^2 +1x^3
4: +1 -4x^1 +6x^2 -4x^3 +1x^4
5: -1 +5x^1 -10x^2 +10x^3 -5x^4 +1x^5
6: +1 -6x^1 +15x^2 -20x^3 +15x^4 -6x^5 +1x^6
7: -1 +7x^1 -21x^2 +35x^3 -35x^4 +21x^5 -7x^6 +1x^7
8: +1 -8x^1 +28x^2 -56x^3 +70x^4 -56x^5 +28x^6 -8x^7 +1x^8
9: -1 +9x^1 -36x^2 +84x^3 -126x^4 +126x^5 -84x^6 +36x^7 -9x^8 +1x^9
10: +1 -10x^1 +45x^2 -120x^3 +210x^4 -252x^5 +210x^6 -120x^7 +45x^8 -10x^9 +1x^10
11: -1 +11x^1 -55x^2 +165x^3 -330x^4 +462x^5 -462x^6 +330x^7 -165x^8 +55x^9 -11x^10 +1x^11
# small primes using the aks test
[2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43, 47, 53, 59, 61, 67, 71, 73, 79, 83, 89, 97]
</pre>
=={{header|8th}}==
<
with: a
Line 96 ⟶ 139:
"The primes upto 50 are (via AKS): " . .primes-via-aks cr
bye</
{{out}}
<pre>
Line 109 ⟶ 152:
The primes upto 50 are (via AKS): 2 3 5 7 11 13 17 19 23 29 31 37 41 43 47
</pre>
=={{header|AArch64 Assembly}}==
{{works with|as|Raspberry Pi 3B version Buster 64 bits <br> or android 64 bits with application Termux }}
<syntaxhighlight lang="aarch64 assembly">
/* ARM assembly AARCH64 Raspberry PI 3B or android 64 bits */
/* program AKS64.s */
/*******************************************/
/* Constantes file */
/*******************************************/
/* for this file see task include a file in language AArch64 assembly*/
.include "../includeConstantesARM64.inc"
.equ MAXI, 64
.equ NUMBERLOOP, 10
/*********************************/
/* Initialized data */
/*********************************/
.data
szMessResult: .asciz " (x-1)^@ = "
szMessResult1: .asciz " @ x^@ "
szMessResPrime: .asciz "Number @ is prime. \n"
szCarriageReturn: .asciz "\n"
/*********************************/
/* UnInitialized data */
/*********************************/
.bss
sZoneConv: .skip 24
qTabCoef: .skip 8 * MAXI
/*********************************/
/* code section */
/*********************************/
.text
.global main
main: // entry of program
mov x4,#1
1: // loop
mov x0,x4
bl computeCoef // compute coefficient
ldr x0,qAdrqTabCoef
mov x0,x4
bl displayCoef // display coefficient
add x4,x4,1
cmp x4,NUMBERLOOP
blt 1b
mov x4,1
2:
mov x0,x4
bl isPrime // is prime ?
cmp x0,1
bne 3f
mov x0,x4
ldr x1,qAdrsZoneConv
bl conversion10 // call decimal conversion
add x1,x1,x0
strb wzr,[x1]
ldr x0,qAdrszMessResPrime
ldr x1,qAdrsZoneConv // insert value conversion in message
bl strInsertAtCharInc
bl affichageMess
3:
add x4,x4,1
cmp x4,MAXI
blt 2b
100: // standard end of the program
mov x0,0 // return code
mov x8,EXIT // request to exit program
svc 0 // perform the system call
qAdrszCarriageReturn: .quad szCarriageReturn
qAdrsZoneConv: .quad sZoneConv
qAdrqTabCoef: .quad qTabCoef
qAdrszMessResPrime: .quad szMessResPrime
/***************************************************/
/* display coefficients */
/***************************************************/
// x0 contains a number
displayCoef:
stp x1,lr,[sp,-16]! // save registres
stp x2,x3,[sp,-16]! // save registres
stp x4,x5,[sp,-16]! // save registres
stp x6,x7,[sp,-16]! // save registres
mov x2,x0
ldr x1,qAdrsZoneConv //
bl conversion10 // call decimal conversion
add x1,x1,x0
strb wzr,[x1]
ldr x0,qAdrszMessResult
ldr x1,qAdrsZoneConv // insert value conversion in message
bl strInsertAtCharInc
bl affichageMess
ldr x3,qAdrqTabCoef
1:
ldr x0,[x3,x2,lsl #3]
ldr x1,qAdrsZoneConv //
bl conversion10S // call decimal conversion
2: // removing spaces
ldrb w6,[x1]
cmp x6,' '
cinc x1,x1,eq
beq 2b
ldr x0,qAdrszMessResult1
bl strInsertAtCharInc
mov x4,x0
mov x0,x2
ldr x1,qAdrsZoneConv // else display odd message
bl conversion10 // call decimal conversion
add x1,x1,x0
strb wzr,[x1]
mov x0,x4
ldr x1,qAdrsZoneConv // insert value conversion in message
bl strInsertAtCharInc
bl affichageMess
subs x2,x2,#1
bge 1b
ldr x0,qAdrszCarriageReturn
bl affichageMess
100:
ldp x6,x7,[sp],16 // restaur des 2 registres
ldp x4,x5,[sp],16 // restaur des 2 registres
ldp x2,x3,[sp],16 // restaur des 2 registres
ldp x1,lr,[sp],16 // restaur des 2 registres
ret
qAdrszMessResult: .quad szMessResult
qAdrszMessResult1: .quad szMessResult1
/***************************************************/
/* compute coefficient */
/***************************************************/
// x0 contains a number
computeCoef:
stp x1,lr,[sp,-16]! // save registres
stp x2,x3,[sp,-16]! // save registres
stp x4,x5,[sp,-16]! // save registres
stp x6,x7,[sp,-16]! // save registres
ldr x1,qAdrqTabCoef // address coefficient array
mov x2,1
str x2,[x1] // store 1 to coeff [0]
mov x3,0 // indice 1
1:
add x4,x3,1
mov x5,1
str x5,[x1,x4,lsl #3]
mov x6,x3 // indice 2 = indice 1
2:
cmp x6,0 // zero ? -> end loop
ble 3f
sub x4,x6,1
ldr x5,[x1,x4,lsl 3]
ldr x4,[x1,x6,lsl 3]
sub x5,x5,x4
str x5,[x1,x6,lsl 3]
sub x6,x6,1
b 2b
3:
ldr x2,[x1] // inversion coeff [0]
neg x2,x2
str x2,[x1]
add x3,x3,1
cmp x3,x0
blt 1b
100:
ldp x6,x7,[sp],16 // restaur des 2 registres
ldp x4,x5,[sp],16 // restaur des 2 registres
ldp x2,x3,[sp],16 // restaur des 2 registres
ldp x1,lr,[sp],16 // restaur des 2 registres
ret
/***************************************************/
/* verify number is prime */
/***************************************************/
// x0 contains a number
isPrime:
stp x1,lr,[sp,-16]! // save registres
stp x2,x3,[sp,-16]! // save registres
stp x4,x5,[sp,-16]! // save registres
bl computeCoef
ldr x4,qAdrqTabCoef // address coefficient array
ldr x2,[x4]
add x2,x2,1
str x2,[x4]
ldr x2,[x4,x0,lsl 3]
sub x2,x2,#1
str x2,[x4,x0,lsl 3]
mov x5,x0 // number start
1:
ldr x1,[x4,x5,lsl 3] // load one coeff
sdiv x2,x1,x0
msub x3,x2,x0,x1 // compute remainder
cmp x3,#0 // remainder = zéro ?
bne 99f // if <> no prime
subs x5,x5,#1 // next coef
bgt 1b // and loop
mov x0,#1 // prime
b 100f
99:
mov x0,0 // no prime
100:
ldp x4,x5,[sp],16 // restaur des 2 registres
ldp x2,x3,[sp],16 // restaur des 2 registres
ldp x1,lr,[sp],16 // restaur des 2 registres
ret
/********************************************************/
/* File Include fonctions */
/********************************************************/
/* for this file see task include a file in language AArch64 assembly */
.include "../includeARM64.inc"
</syntaxhighlight>
{{Output}}
<pre>
(x-1)^1 = +1 x^1 -1 x^0
(x-1)^2 = +1 x^2 -2 x^1 +1 x^0
(x-1)^3 = +1 x^3 -3 x^2 +3 x^1 -1 x^0
(x-1)^4 = +1 x^4 -4 x^3 +6 x^2 -4 x^1 +1 x^0
(x-1)^5 = +1 x^5 -5 x^4 +10 x^3 -10 x^2 +5 x^1 -1 x^0
(x-1)^6 = +1 x^6 -6 x^5 +15 x^4 -20 x^3 +15 x^2 -6 x^1 +1 x^0
(x-1)^7 = +1 x^7 -7 x^6 +21 x^5 -35 x^4 +35 x^3 -21 x^2 +7 x^1 -1 x^0
(x-1)^8 = +1 x^8 -8 x^7 +28 x^6 -56 x^5 +70 x^4 -56 x^3 +28 x^2 -8 x^1 +1 x^0
(x-1)^9 = +1 x^9 -9 x^8 +36 x^7 -84 x^6 +126 x^5 -126 x^4 +84 x^3 -36 x^2 +9 x^1 -1 x^0
Number 1 is prime.
Number 2 is prime.
Number 3 is prime.
Number 5 is prime.
Number 7 is prime.
Number 11 is prime.
Number 13 is prime.
Number 17 is prime.
Number 19 is prime.
Number 23 is prime.
Number 29 is prime.
Number 31 is prime.
Number 37 is prime.
Number 41 is prime.
Number 43 is prime.
Number 47 is prime.
Number 53 is prime.
Number 59 is prime.
Number 61 is prime.
</pre>
=={{header|Ada}}==
<
procedure Test_For_Primes is
Line 199 ⟶ 487:
Show_Pascal_Triangles;
Show_Primes;
end Test_For_Primes;</
{{out}}
Line 218 ⟶ 506:
The code below uses Algol 68 Genie which provides arbitrary precision arithmetic for LONG LONG modes.
<
BEGIN
COMMENT
Line 305 ⟶ 593:
print (newline)
END
</syntaxhighlight>
{{out}}
<pre>
Line 319 ⟶ 607:
And just to show off, the primes between 900 and 1000 are 907 911 919 929 937 941 947 953 967 971 977 983 991 997
</pre>
=={{header|ARM Assembly}}==
{{works with|as|Raspberry Pi <br> or android 32 bits with application Termux}}
<syntaxhighlight lang="arm assembly">
/* ARM assembly Raspberry PI or android 32 bits */
/* program AKS.s */
/* REMARK 1 : this program use routines in a include file
see task Include a file language arm assembly
for the routine affichageMess conversion10
see at end of this program the instruction include */
/* for constantes see task include a file in arm assembly */
/************************************/
/* Constantes */
/************************************/
.include "../constantes.inc"
.equ MAXI, 32
.equ NUMBERLOOP, 10
/*********************************/
/* Initialized data */
/*********************************/
.data
szMessResult: .asciz " (x-1)^@ = "
szMessResult1: .asciz " @ x^@ "
szMessResPrime: .asciz "Number @ is prime. \n"
szCarriageReturn: .asciz "\n"
/*********************************/
/* UnInitialized data */
/*********************************/
.bss
sZoneConv: .skip 24
iTabCoef: .skip 4 * MAXI
/*********************************/
/* code section */
/*********************************/
.text
.global main
main: @ entry of program
mov r4,#1
1: @ loop
mov r0,r4
bl computeCoef @ compute coefficient
ldr r0,iAdriTabCoef
mov r0,r4
bl displayCoef @ display coefficient
add r4,r4,#1
cmp r4,#NUMBERLOOP
blt 1b
mov r4,#1
2:
mov r0,r4
bl isPrime @ is prime ?
cmp r0,#1
bne 3f
mov r0,r4
ldr r1,iAdrsZoneConv
bl conversion10 @ call decimal conversion
add r1,r0
mov r5,#0
strb r5,[r1]
ldr r0,iAdrszMessResPrime
ldr r1,iAdrsZoneConv @ insert value conversion in message
bl strInsertAtCharInc
bl affichageMess
3:
add r4,r4,#1
cmp r4,#MAXI
blt 2b
100: @ standard end of the program
mov r0, #0 @ return code
mov r7, #EXIT @ request to exit program
svc #0 @ perform the system call
iAdrszCarriageReturn: .int szCarriageReturn
iAdrsZoneConv: .int sZoneConv
iAdriTabCoef: .int iTabCoef
iAdrszMessResPrime: .int szMessResPrime
/***************************************************/
/* display coefficients */
/***************************************************/
// r0 contains a number
displayCoef:
push {r1-r6,lr} @ save registers
mov r2,r0
ldr r1,iAdrsZoneConv @
bl conversion10 @ call decimal conversion
add r1,r0
mov r5,#0
strb r5,[r1]
ldr r0,iAdrszMessResult
ldr r1,iAdrsZoneConv @ insert value conversion in message
bl strInsertAtCharInc
bl affichageMess
ldr r3,iAdriTabCoef
1:
ldr r0,[r3,r2,lsl #2]
ldr r1,iAdrsZoneConv @
bl conversion10S @ call decimal conversion
2: @ removing spaces
ldrb r6,[r1]
cmp r6,#' '
addeq r1,#1
beq 2b
ldr r0,iAdrszMessResult1
bl strInsertAtCharInc
mov r4,r0
mov r0,r2
ldr r1,iAdrsZoneConv @ else display odd message
bl conversion10 @ call decimal conversion
add r1,r0
mov r5,#0
strb r5,[r1]
mov r0,r4
ldr r1,iAdrsZoneConv @ insert value conversion in message
bl strInsertAtCharInc
bl affichageMess
subs r2,r2,#1
bge 1b
ldr r0,iAdrszCarriageReturn
bl affichageMess
100:
pop {r1-r6,lr} @ restaur registers
bx lr @ return
iAdrszMessResult: .int szMessResult
iAdrszMessResult1: .int szMessResult1
/***************************************************/
/* compute coefficient */
/***************************************************/
// r0 contains a number
computeCoef:
push {r1-r6,lr} @ save registers
ldr r1,iAdriTabCoef @ address coefficient array
mov r2,#1
str r2,[r1] @ store 1 to coeff [0]
mov r3,#0 @ indice 1
1:
add r4,r3,#1
mov r5,#1
str r5,[r1,r4,lsl #2]
mov r6,r3 @ indice 2 = indice 1
2:
cmp r6,#0 @ zero ? -> end loop
ble 3f
sub r4,r6,#1
ldr r5,[r1,r4,lsl #2]
ldr r4,[r1,r6,lsl #2]
sub r5,r5,r4
str r5,[r1,r6,lsl #2]
sub r6,r6,#1
b 2b
3:
ldr r2,[r1] @ inversion coeff [0]
neg r2,r2
str r2,[r1]
add r3,r3,#1
cmp r3,r0
blt 1b
100:
pop {r1-r6,lr} @ restaur registers
bx lr @ return
/***************************************************/
/* verify number is prime */
/***************************************************/
// r0 contains a number
isPrime:
push {r1-r5,lr} @ save registers
bl computeCoef
ldr r4,iAdriTabCoef @ address coefficient array
ldr r2,[r4]
add r2,r2,#1
str r2,[r4]
ldr r2,[r4,r0,lsl #2]
sub r2,r2,#1
str r2,[r4,r0,lsl #2]
mov r5,r0 @ number start
mov r1,r0 @ divisor
1:
ldr r0,[r4,r5,lsl #2] @ load one coeff
cmp r0,#0 @ if negative inversion
neglt r0,r0
bl division @ because this routine is number positive only
cmp r3,#0 @ remainder = zéro ?
movne r0,#0 @ if <> no prime
bne 100f
subs r5,r5,#1 @ next coef
bgt 1b
mov r0,#1 @ prime
100:
pop {r1-r5,lr} @ restaur registers
bx lr @ return
/***************************************************/
/* ROUTINES INCLUDE */
/***************************************************/
.include "../affichage.inc"
</syntaxhighlight>
{{Output}}
<pre>
(x-1)^1 = +1 x^1 -1 x^0
(x-1)^2 = +1 x^2 -2 x^1 +1 x^0
(x-1)^3 = +1 x^3 -3 x^2 +3 x^1 -1 x^0
(x-1)^4 = +1 x^4 -4 x^3 +6 x^2 -4 x^1 +1 x^0
(x-1)^5 = +1 x^5 -5 x^4 +10 x^3 -10 x^2 +5 x^1 -1 x^0
(x-1)^6 = +1 x^6 -6 x^5 +15 x^4 -20 x^3 +15 x^2 -6 x^1 +1 x^0
(x-1)^7 = +1 x^7 -7 x^6 +21 x^5 -35 x^4 +35 x^3 -21 x^2 +7 x^1 -1 x^0
(x-1)^8 = +1 x^8 -8 x^7 +28 x^6 -56 x^5 +70 x^4 -56 x^3 +28 x^2 -8 x^1 +1 x^0
(x-1)^9 = +1 x^9 -9 x^8 +36 x^7 -84 x^6 +126 x^5 -126 x^4 +84 x^3 -36 x^2 +9 x^1 -1 x^0
Number 1 is prime.
Number 2 is prime.
Number 3 is prime.
Number 5 is prime.
Number 7 is prime.
Number 11 is prime.
Number 13 is prime.
Number 17 is prime.
Number 19 is prime.
Number 23 is prime.
Number 29 is prime.
Number 31 is prime.
</pre>
=={{header|AutoHotkey}}==
{{works with|AutoHotkey L}}
<
; Function modified from http://rosettacode.org/wiki/Pascal%27s_triangle#AutoHotkey
pascalstriangle(n=8) ; n rows of Pascal's triangle
Line 374 ⟶ 889:
}
Msgbox % t
Return</
{{out}}
<pre>(x-1)^0=+1
Line 393 ⟶ 908:
The primality test uses a pattern that looks for a fractional factor. If such a factor is found, the test fails. Otherwise it succeeds.
<
& (expandx-1P=.forceExpansion$((x+-1)^!arg))
& ( isPrime
Line 419 ⟶ 934:
& out$"2 <= Primes <= 50:"
& out$!primes
);</
Output:
<pre>Polynomial representations of (x-1)^p for p <= 7 :
Line 441 ⟶ 956:
2 3 5 7 11 13 17 19 23 29 31 37 41 43 47</pre>
The AKS test kan be written more concisely than the task describes. This prints the primes between 980 and 1000:
<
& 980:?n
& whl
Line 449 ⟶ 964:
)
)
);</
Output:
<pre>Primes between 980 and 1000, short version:
Line 457 ⟶ 972:
=={{header|C}}==
<
#include <stdlib.h>
Line 507 ⟶ 1,022:
putchar('\n');
return 0;
}</
The ugly output:
Line 524 ⟶ 1,039:
primes (never mind the 1): 1 2 3 5 7 11 13 17 19 23 29 31 37 41 43 47 53 59 61
</pre>
=={{header|C sharp|C#}}==
{{trans|C}}
<syntaxhighlight lang="csharp">
using System;
public class AksTest
{
static long[] c = new long[100];
static void Main(string[] args)
{
for (int n = 0; n < 10; n++) {
coef(n);
Console.Write("(x-1)^" + n + " = ");
show(n);
Console.WriteLine("");
}
Console.Write("Primes:");
for (int n = 1; n <= 63; n++)
if (is_prime(n))
Console.Write(n + " ");
Console.WriteLine('\n');
Console.ReadLine();
}
static void coef(int n)
{
int i, j;
if (n < 0 || n > 63) System.Environment.Exit(0);// gracefully deal with range issue
for (c[i = 0] = 1L; i < n; c[0] = -c[0], i++)
for (c[1 + (j = i)] = 1L; j > 0; j--)
c[j] = c[j - 1] - c[j];
}
static bool is_prime(int n)
{
int i;
coef(n);
c[0] += 1;
c[i = n] -= 1;
while (i-- != 0 && (c[i] % n) == 0) ;
return i < 0;
}
static void show(int n)
{
do {
Console.Write("+" + c[n] + "x^" + n);
}while (n-- != 0);
}
}
</syntaxhighlight>
=={{header|C++}}==
{{trans|Pascal}}
<
#include <iomanip>
#include <iostream>
Line 616 ⟶ 1,189:
cout << endl;
}
</syntaxhighlight>
{{out}}
<pre>
Line 631 ⟶ 1,204:
2 3 5 7 11 13 17 19 23 29 31 37 41 43 47 53 59 61
</pre>
=={{header|Clojure}}==
{{incorrect|Clojure|(x-1)(x-1) is x^2-2x+1 so should the result 2 (-1 2 1) rather be 2 (1 2 1) and similarly for several others}}
The *' function is an arbitrary precision multiplication.
<
"kth coefficient of (x - 1)^n"
[n k]
Line 709 ⟶ 1,225:
(doseq [n (range 10)] (println n (cs n)))
(println)
(println "primes < 50 per AKS:" (filter aks? (range 2 50)))</
{{output}}
<pre>coefficient series n (k[0] .. k[n])
Line 726 ⟶ 1,242:
=={{header|CoffeeScript}}==
<
a = []
return () ->
Line 764 ⟶ 1,280:
console.log ""
console.log "The primes upto 50 are: #{primes}"</
{{out}}
<pre>(x - 1)^0 = + 1
Line 778 ⟶ 1,294:
=={{header|Common Lisp}}==
<
(cond
((= p 0) #(1))
Line 813 ⟶ 1,329:
(loop for i from 0 to 50
do (when (primep i) (format t "~D " i)))
(format t "~%"))</
{{out}}
<pre># p: (x-1)^p for small p:
Line 828 ⟶ 1,344:
=={{header|Crystal}}==
{{trans|Ruby}}
<
p.times.reduce([1]) do |ex, _|
([0_i64] + ex).zip(ex + [0]).map { |x, y| x - y }
end
end
Line 836 ⟶ 1,352:
def prime?(p)
return false if p < 2
coeff = x_minus_1_to_the(p)[1..p//2] # only need half of coeff terms
coeff.all?{ |n| n%p == 0 }
end
Line 847 ⟶ 1,363:
end
puts "\nPrimes below 50:", 50.times.select { |n| prime? n }.join(',')
</syntaxhighlight>
{{out}}
Line 867 ⟶ 1,383:
=={{header|D}}==
{{trans|Python}}
<
BigInt[] expandX1(in uint p) pure /*nothrow*/ {
Line 895 ⟶ 1,411:
"\nSmall primes using the AKS test:".writeln;
101.iota.filter!aksTest.writeln;
}</
{{out}}
<pre># p: (x-1)^p for small p:
Line 913 ⟶ 1,429:
Small primes using the AKS test:
[2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43, 47, 53, 59, 61, 67, 71, 73, 79, 83, 89, 97]</pre>
=={{header|Delphi}}==
See [[#Pascal]].
=={{header|EchoLisp}}==
We use the math.lib library and the poly functions to compute and display the required polynomials. A polynomial P(x) = a0 +a1*x + .. an*x^n is a list of coefficients (a0 a1 .... an).
<
(lib 'math.lib)
;; 1 - x^p : P = (1 0 0 0 ... 0 -1)
Line 936 ⟶ 1,453:
(test (lambda(a) (zero? (modulo a p))))) ;; p divides a[i] ?
(apply and (map test P)))) ;; returns #t if true for all a[i]
</syntaxhighlight>
{{Output}}
<
(show-them 13) →
p 1 x -1
Line 960 ⟶ 1,477:
2 3 5 7 11 13 17 19 23 29 31 37 41 43 47 53 59 61 67 71 73 79 83 89 97
</syntaxhighlight>
=={{header|Elena}}==
{{trans|C#}}
ELENA
<
singleton AksTest
Line 978 ⟶ 1,496:
c[i] := 1l;
for (int i := 0
c[1 + i] := 1l;
for (int j := i
c[j] := c[j - 1] - c[j]
};
Line 1,018 ⟶ 1,536:
public program()
{
for (int n := 0
AksTest.coef(n);
Line 1,027 ⟶ 1,545:
console.print("Primes:");
for (int n := 1
if (AksTest.is_prime(n))
{
Line 1,035 ⟶ 1,553:
console.printLine().readChar()
}</
{{out}}
<pre>
Line 1,053 ⟶ 1,571:
=={{header|Elixir}}==
{{trans|Erlang}}
<
def iterate(f, x), do: fn -> [x | iterate(f, f.(x))] end
Line 1,091 ⟶ 1,609:
end
AKS.main</
{{out}}
Line 1,112 ⟶ 1,630:
The Erlang io module can print out lists of characters with any level of nesting as a flat string. (e.g. ["Er", ["la", ["n"]], "g"] prints as "Erlang") which is useful when constructing the strings to print out for the binomial expansions. The program also shows how lazy lists can be implemented in Erlang.
<
-import(lists, [all/2, seq/2, zip/2]).
Line 1,151 ⟶ 1,669:
[[N || {Row, N} <- zip(tl(tl(take(51, pascal()))), seq(2, 50)),
primerow(Row, N)]]).
</syntaxhighlight>
{{out}}
<pre>(x - 1)^0 = + 1
Line 1,163 ⟶ 1,681:
The primes upto 50: [2,3,5,7,11,13,17,19,23,29,31,37,41,43,47]
</pre>
=={{header|F_Sharp|F#}}==
<syntaxhighlight lang="fsharp">
// N-grams. Nigel Galloway: April 11th., 2024
let fN g=let n,g=g@[0I],0I::g in List.map2(fun n g->n-g) n g
Seq.unfold(fun g->Some(g,fN g))[1I]|>Seq.take 9|>Seq.iteri(fun n g->printfn "%d -> %A" n g); printfn ""
let fG (n::g) l=(n+1I)%l=0I && g|>List.forall(fun n->n%l=0I)
Seq.unfold(fun(n,g)->Some((n,g),(n+1I,fN g)))(0I,[1I])|>Seq.skip 2|>Seq.filter(fun(n,_::g)->fG (List.rev g) n)|>Seq.takeWhile(fun(n,_)->n<100I)|>Seq.iter(fun(n,_)->printf "%A " n); printfn ""
</syntaxhighlight>
{{out}}
<pre>
0 -> [1]
1 -> [1; -1]
2 -> [1; -2; 1]
3 -> [1; -3; 3; -1]
4 -> [1; -4; 6; -4; 1]
5 -> [1; -5; 10; -10; 5; -1]
6 -> [1; -6; 15; -20; 15; -6; 1]
7 -> [1; -7; 21; -35; 35; -21; 7; -1]
8 -> [1; -8; 28; -56; 70; -56; 28; -8; 1]
2 3 5 7 11 13 17 19 23 29 31 37 41 43 47 53 59 61 67 71 73 79 83 89 97
</pre>
=={{header|Factor}}==
<
math.polynomials prettyprint sequences ;
IN: rosetta-code.aks-test
Line 1,204 ⟶ 1,745:
: aks-test ( -- ) part1 nl part2 ;
MAIN: aks-test</
{{out}}
<pre>
Line 1,221 ⟶ 1,762:
=={{header|Forth}}==
<
1 swap 1+ dup 1 ?do over over i - i */ negate swap loop drop ;
Line 1,245 ⟶ 1,786:
11 0 ?do i . ." : " i .poly cr loop cr
50 1 ?do i prime? if i . then loop
cr ;</
{{out}}
<pre>
Line 1,261 ⟶ 1,802:
2 3 5 7 11 13 17 19 23 29 31 37 41 43 47</pre>
=={{header|Fortran}}==
<
program aks
implicit none
Line 1,401 ⟶ 1,943:
end function
end program aks
</syntaxhighlight>
{{out}}
<pre>
Line 1,414 ⟶ 1,956:
2 3 5 7 11 13 17 19 23 29 31
2 3 5 7 11 13 17 19 23 29 31 37 41 43 47 53 59 61 67 71 73 79 83 89 97 101 103 107 109 113
</pre>
=={{header|FreeBASIC}}==
<
'UPPER LIMIT OF FREEBASIC ULONGINT GETS PRIMES UP TO 70
Sub string_split(s_in As String,char As String,result() As String)
Line 1,571 ⟶ 2,033:
Next n
Sleep</
{{out}}
<pre>(x-1)^1 = -1 +x
Line 1,585 ⟶ 2,047:
Primes up to 70:
2 3 5 7 11 13 17 19 23 29 31 37 41 43 47 53 59 61 67</pre>
=={{header|Go}}==
<syntaxhighlight lang="go">package main
import "fmt"
func bc(p int) []int64 {
c := make([]int64, p+1)
r := int64(1)
for i, half := 0, p/2; i <= half; i++ {
c[i] = r
c[p-i] = r
r = r * int64(p-i) / int64(i+1)
}
for i := p - 1; i >= 0; i -= 2 {
c[i] = -c[i]
}
return c
}
func main() {
for p := 0; p <= 7; p++ {
fmt.Printf("%d: %s\n", p, pp(bc(p)))
}
for p := 2; p < 50; p++ {
if aks(p) {
fmt.Print(p, " ")
}
}
fmt.Println()
}
var e = []rune("²³⁴⁵⁶⁷")
func pp(c []int64) (s string) {
if len(c) == 1 {
return fmt.Sprint(c[0])
}
p := len(c) - 1
if c[p] != 1 {
s = fmt.Sprint(c[p])
}
for i := p; i > 0; i-- {
s += "x"
if i != 1 {
s += string(e[i-2])
}
if d := c[i-1]; d < 0 {
s += fmt.Sprintf(" - %d", -d)
} else {
s += fmt.Sprintf(" + %d", d)
}
}
return
}
func aks(p int) bool {
c := bc(p)
c[p]--
c[0]++
for _, d := range c {
if d%int64(p) != 0 {
return false
}
}
return true
}</syntaxhighlight>
{{out}}
<pre>
0: 1
1: x - 1
2: x² - 2x + 1
3: x³ - 3x² + 3x - 1
4: x⁴ - 4x³ + 6x² - 4x + 1
5: x⁵ - 5x⁴ + 10x³ - 10x² + 5x - 1
6: x⁶ - 6x⁵ + 15x⁴ - 20x³ + 15x² - 6x + 1
7: x⁷ - 7x⁶ + 21x⁵ - 35x⁴ + 35x³ - 21x² + 7x - 1
2 3 5 7 11 13 17 19 23 29 31 37 41 43 47
</pre>
=={{header|Haskell}}==
<
Line 1,606 ⟶ 2,147:
putStrLn $ unlines [show i ++ ": " ++ printPoly (expand i) | i <- [0..10]]
putStrLn "-- Primes up to 100:"
print (filter test [1..100])</
{{out}}
<pre>-- p: (x-1)^p for small p
Line 1,625 ⟶ 2,166:
=={{header|Idris}}==
<
-- Computes Binomial Coefficients
Line 1,694 ⟶ 2,235:
printExpansions 10
putStrLn "\n-- Primes Up To 100:"
putStrLn $ show $ filter isPrime [0..100]</
{{out}}
<pre>-- p: (x-1)^p for small p
Line 1,714 ⟶ 2,255:
=={{header|J}}==
'''Solution''':<
testAKS =: 0 *./ .= ] | binomialExpansion NB. 3) Use that function to create another which determines whether p is prime using AKS.</
'''Examples''':<
+-++--+----+-------+-----------+---------------+------------------+
|0||_2|_3 3|_4 6 _4|_5 10 _10 5|_6 15 _20 15 _6|_7 21 _35 35 _21 7|
Line 1,726 ⟶ 2,267:
2 3 5 7 11 13 17 19 23 29 31 37 41 43 47
i.&.:(_1&p:) 50 NB. Double-check our results using built-in prime filter.
2 3 5 7 11 13 17 19 23 29 31 37 41 43 47</
=={{header|Java}}==
{{trans|C}}
'''Solution''':<
private static final long[] c = new long[64];
Line 1,775 ⟶ 2,316:
System.out.println();
}
}</
Output:
<pre>
Line 1,793 ⟶ 2,334:
=={{header|JavaScript}}==
{{trans|CoffeeScript}}
<
pascal = function() {
Line 1,871 ⟶ 2,412:
console.log("");
console.log("The primes upto 50 are: " + primes);</
{{out}}
<pre>(x - 1)^0 = + 1
Line 1,884 ⟶ 2,425:
The primes upto 50 are: 2,3,5,7,11,13,17,19,23,29,31,37,41,43,47</pre>
Reviewed (ES6):
<
var cs = []; if (n) while (n--) coef(); return coef
function coef() {
Line 1,906 ⟶ 2,447:
document.write('<br>Primes: ');
for (var coef=pascal(2), n=2; n<=50; n+=1) if (isPrime(coef())) document.write(' ', n)</
{{output}}
(x-1)<sup>0</sup> = +1
Line 1,919 ⟶ 2,460:
Primes: 2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43, 47
{{trans|C}}
<
for (var c=[1], i=0; i<n; c[0]=-c[0], i+=1) {
c[i+1]=1; for (var j=i; j; j-=1) c[j] = c[j-1]-c[j]
Line 1,940 ⟶ 2,481:
document.write('<br>Primes: ');
for (var n=2; n<=50; n++) if (isPrime(n)) document.write(' ', n)</
{{output}}
(x-1)<sup>0</sup> = +1
Line 1,968 ⟶ 2,509:
easily be adapted to use a BigInt library such as the one at
https://github.com/joelpurra/jq-bigint
<
# Input: an OptPascal array
# Output: the next OptPascal array (obtained by adding adjacent items,
Line 2,002 ⟶ 2,543:
def pascal: nth(.; pascals);
def optpascal: nth(.; optpascals);</
'''Task 1:''' "A method to generate the coefficients of (x-1)^p"
<
def alternate_signs: . as $in
| reduce range(0; length) as $i ([]; . + [$in[$i] * (if $i % 2 == 0 then 1 else -1 end )]);
(.+1) | pascal | alternate_signs;</
'''Task 2:''' "Show here the polynomial expansions of (x − 1)^p for p in the range 0 to at least 7, inclusive."
<
{{out}}
<
Coefficients for (x - 1)^1: [1,-1]
Coefficients for (x - 1)^2: [1,-2,1]
Line 2,020 ⟶ 2,561:
Coefficients for (x - 1)^5: [1,-5,10,-10,5,-1]
Coefficients for (x - 1)^6: [1,-6,15,-20,15,-6,1]
Coefficients for (x - 1)^7: [1,-7,21,-35,35,-21,7,-1]</
'''Task 3:''' Prime Number Test
For brevity, we show here only the relatively efficient solution based on optpascal/0:
<
. as $N
| if . < 2 then false
else (1+.) | optpascal
| all( .[2:][]; . % $N == 0 )
end;</
'''Task 4:''' "Use your AKS test to generate a list of all primes under 35."
<
{{out}}
<
3
5
Line 2,045 ⟶ 2,586:
23
29
31</
'''Task 5:''' "As a stretch goal, generate all primes under 50."
<
{{out}}
<
=={{header|Julia}}==
'''Task 1'''
<syntaxhighlight lang="julia">
function polycoefs(n::Int64)
pc = typeof(n)[]
Line 2,068 ⟶ 2,609:
return pc
end
</syntaxhighlight>
Perhaps this should be done with a comprehension, but properly accounting for the sign is tricky in that case.
Line 2,074 ⟶ 2,615:
'''Task 2'''
<
function stringpoly(n::Int64)
Line 2,103 ⟶ 2,644:
return st
end
</syntaxhighlight>
Of course this could be simpler, but this produces a nice payoff in typeset equations that do on include extraneous characters (leading pluses and coefficients of 1).
Line 2,109 ⟶ 2,650:
'''Task 3'''
<syntaxhighlight lang="julia">
function isaksprime(n::Int64)
if n < 2
Line 2,121 ⟶ 2,662:
return true
end
</syntaxhighlight>
'''Task 4'''
<syntaxhighlight lang="julia">
println("<math>")
println("\\begin{array}{lcl}")
Line 2,144 ⟶ 2,685:
end
println()
</syntaxhighlight>
{{out}}
Line 2,167 ⟶ 2,708:
=={{header|Kotlin}}==
<
fun binomial(n: Int, k: Int): Long = when {
Line 2,222 ⟶ 2,763:
println("\nThe prime numbers under 62 are:")
println(primes)
}</
{{out}}
Line 2,240 ⟶ 2,781:
[2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43, 47, 53, 59, 61]
</pre>
=={{header|LFE}}==
Loosely translated from Erlang.
<syntaxhighlight lang="lisp">
(defun next-row (row)
(cons 1
(cl:maplist
(match-lambda
(((list a)) a)
(((cons a (cons b _))) (+ a b)))
row)))
(defun pascal (n)
(pascal n '(())))
(defun pascal
((0 rows) (cdr (lists:reverse rows)))
((n (= (cons row _) rows))
(pascal (- n 1) (cons (next-row row) rows))))
(defun show-x
((0) "")
((1) "x")
((n) (++ "x^" (integer_to_list n))))
(defun rhs
(('() _ _) "")
(((cons coef coefs) sgn exp)
(++
(if (< sgn 0) " - " " + ")
(integer_to_list coef)
(show-x exp)
(rhs coefs (- sgn) (- exp 1)))))
(defun binomial-text (row)
(let ((degree (- (length row) 1)))
(++ "(x - 1)^" (integer_to_list degree) " =" (rhs row 1 degree))))
(defun primerow
(('() _)
'true)
((`(1 . ,rest) n)
(primerow rest n))
(((cons a (cons a _)) n) ; stop when we've checked half the list
(=:= 0 (rem a n)))
(((cons a rest) n)
(andalso
(=:= 0 (rem a n))
(primerow rest n))))
(defun main (_)
(list-comp
((<- row (pascal 8)))
(lfe_io:format "~s~n" (list (binomial-text row))))
(lfe_io:format "~nThe primes upto 50: ~p~n"
(list
(list-comp
((<- (tuple row n) (lists:zip (cddr (pascal 51)) (lists:seq 2 50)))
(primerow row n))
n))))
</syntaxhighlight>
{{Out}}
<pre>
(x - 1)^0 = + 1
(x - 1)^1 = + 1x - 1
(x - 1)^2 = + 1x^2 - 2x + 1
(x - 1)^3 = + 1x^3 - 3x^2 + 3x - 1
(x - 1)^4 = + 1x^4 - 4x^3 + 6x^2 - 4x + 1
(x - 1)^5 = + 1x^5 - 5x^4 + 10x^3 - 10x^2 + 5x - 1
(x - 1)^6 = + 1x^6 - 6x^5 + 15x^4 - 20x^3 + 15x^2 - 6x + 1
(x - 1)^7 = + 1x^7 - 7x^6 + 21x^5 - 35x^4 + 35x^3 - 21x^2 + 7x - 1
The primes upto 50: (2 3 5 7 11 13 17 19 23 29 31 37 41 43 47)
</pre>
=={{header|Lua}}==
<syntaxhighlight lang="lua">-- AKS test for primes, in Lua, 6/23/2020 db
local function coefs(n)
local list = {[0]=1}
for k = 0, n do list[k+1] = math.floor(list[k] * (n-k) / (k+1)) end
for k = 1, n, 2 do list[k] = -list[k] end
return list
end
local function isprimeaks(n)
local c = coefs(n)
c[0], c[n] = c[0]-1, c[n]+1
for i = 0, n do
if (c[i] % n ~= 0) then return false end
end
return true
end
local function pprintcoefs(n, list)
local result = ""
for i = 0, n do
local s = i==0 and "" or list[i]>=0 and " + " or " - "
local c, e = math.abs(list[i]), n-i
if (c==1 and e > 0) then c = "" end
local x = e==0 and "" or e==1 and "x" or "x^"..e
result = result .. s .. c .. x
end
print("(x-1)^" .. n .." : " .. result)
end
for i = 0, 9 do
pprintcoefs(i, coefs(i))
end
local primes = {}
for i = 2, 53 do
if (isprimeaks(i)) then primes[#primes+1] = i end
end
print(table.concat(primes, ", "))</syntaxhighlight>
{{out}}
<pre>(x-1)^0 : 1
(x-1)^1 : x - 1
(x-1)^2 : x^2 - 2x + 1
(x-1)^3 : x^3 - 3x^2 + 3x - 1
(x-1)^4 : x^4 - 4x^3 + 6x^2 - 4x + 1
(x-1)^5 : x^5 - 5x^4 + 10x^3 - 10x^2 + 5x - 1
(x-1)^6 : x^6 - 6x^5 + 15x^4 - 20x^3 + 15x^2 - 6x + 1
(x-1)^7 : x^7 - 7x^6 + 21x^5 - 35x^4 + 35x^3 - 21x^2 + 7x - 1
(x-1)^8 : x^8 - 8x^7 + 28x^6 - 56x^5 + 70x^4 - 56x^3 + 28x^2 - 8x + 1
(x-1)^9 : x^9 - 9x^8 + 36x^7 - 84x^6 + 126x^5 - 126x^4 + 84x^3 - 36x^2 + 9x - 1
2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43, 47, 53</pre>
=={{header|Lambdatalk}}==
<syntaxhighlight lang="scheme">
{require lib_BN} // for big numbers
1) pascalian binomial coefficient C(n,p) = n!/(p!(n-p)!) = (n*(n-1)...(n-p+1))/(p*(p-1)...2*1)
{def coeff
{lambda {:n :p}
{BN.intPart
{BN./ {S.reduce BN.* {S.serie :n {- :n :p -1} -1}}
{S.reduce BN.* {S.serie :p 1 -1}}}}}}
-> coeff
2) polynomial expansions of (x − 1)^p
{def sign
{lambda {:n}
{if {= {% :n 2} 0} then + else -}}}
-> sign
{def coeffs
{lambda {:n}
{br}(x - 1)^:n =
{if {= :n 0}
then + 1x^0
else {if {= :n 1}
then + 1x^1 - 1x^0
else {sign 0} 1x^:n
{S.map {{lambda {:p :n} {sign {- :p :n}} {coeff :p :n}x^{- :p :n}} :n}
{S.serie 1 {- :n 1}}}
{sign :n} 1x^0}}}}
-> coeffs
{S.map coeffs {S.serie 0 7}}
->
(x - 1)^0 = + 1x^0
(x - 1)^1 = + 1x^1 - 1x^0
(x - 1)^2 = + 1x^2 - 2x^1 + 1x^0
(x - 1)^3 = + 1x^3 + 3x^2 - 3x^1 - 1x^0
(x - 1)^4 = + 1x^4 - 4x^3 + 6x^2 - 4x^1 + 1x^0
(x - 1)^5 = + 1x^5 + 5x^4 - 10x^3 + 10x^2 - 5x^1 - 1x^0
(x - 1)^6 = + 1x^6 - 6x^5 + 15x^4 - 20x^3 + 15x^2 - 6x^1 + 1x^0
(x - 1)^7 = + 1x^7 + 7x^6 - 21x^5 + 35x^4 - 35x^3 + 21x^2 - 7x^1 - 1x^0
3) primality test
Taking into account the symmetry of the list of coefficients and the uselessness of the sign
in the calculation of the divisibility, one can limit the tests to half of the list,
and define a simplified function, aks_coeffs:
{def aks_coeffs
{lambda {:n}
{S.map {coeff :n} {S.serie 1 {+ {/ {- :n 1} 2} 1}}}}}
-> aks_coeffs
{def divide
{lambda {:a :b}
{= {BN.compare {BN.% :b :a} 0} 0}}}
-> divide
{def isprime
{lambda {:n}
{if {and {S.map {divide :n} {aks_coeffs :n}}} then :n else .}}}
-> isprime
{S.map isprime {S.serie 2 100}}
-> 2 3 . 5 . 7 . . . 11 . 13 . . . 17 . 19 . . . 23 . . . . . 29 . 31 . . . . . 37 . . . 41 . 43 . . . 47
. . . . . 53 . . . . . 59 . 61 . . . . . 67 . . . 71 . 73 . . . . . 79 . . . 83 . . . . . 89 . . . . . . . 97 . . .
</syntaxhighlight>
=={{header|Liberty BASIC}}==
{{trans|Pascal}}
{{works with|Just BASIC|any}}
<syntaxhighlight lang="lb">
global pasTriMax
pasTriMax = 61
Line 2,324 ⟶ 3,060:
wend
end sub
</syntaxhighlight>
{{out}}
<pre>
Line 2,342 ⟶ 3,078:
=={{header|Maple}}==
Maple handles algebraic manipulation of polynomials natively.
<
1
Line 2,364 ⟶ 3,100:
7 6 5 4 3 2
x - 7 x + 21 x - 35 x + 35 x - 21 x + 7 x - 1
</syntaxhighlight>
To implement the primality test, we write the following procedure that uses the (built-in) polynomial expansion to generate a list of coefficients of the expanded polynomial.
<
Use <code>polc</code> to implement <code>prime?</code> which does the primality test.
<
Of course, rather than calling <code>polc</code>, we can inline it, just for the sake of making the whole thing a one-liner (while adding argument type-checking for good measure):
<
This agrees with the built-in primality test <code>isprime</code>:
<
true
</syntaxhighlight>
Use <code>prime?</code> with the built-in Maple <code>select</code> procedure to pick off the primes up to 50:
<
[2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43, 47]
</syntaxhighlight>
=={{header|Mathematica}} / {{header|Wolfram Language}}==
Algebraic manipulation is built into Mathematica, so there's no need to create a function to do (x-1)^p
<
(x - 1)^( Range[0, 7]) // Expand // TableForm
Print["primes under 50"]
Line 2,389 ⟶ 3,124:
coefflist[p_Integer] := Coefficient[poly[p], x, #] & /@ Range[0, p - 1];
AKSPrimeQ[p_Integer] := (Mod[coefflist[p] , p] // Union) == {0};
Select[Range[1, 50], AKSPrimeQ]</
{{out}}
<pre>powers of (x-1)
Line 2,403 ⟶ 3,138:
primes under 50
{1, 2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43, 47}</pre>
=={{header|Maxima}}==
<syntaxhighlight lang="maxima">
/* Function that lists coefficients given the exponent */
pol_to_list(n):=if n=0 then [1] else block(pol:expand((x-1)^n),
makelist(numfactor(part(pol,i)),i,1,length(pol)))$
/* Function to expand the polynomial (x-1)^n */
expansion(n):=block(
coeflist:pol_to_list(n),
makelist(x^(length(%%)-i),i,1,length(%%)),
%%*coeflist,
apply("+",%%))$
/* AKS based predicate function */
aks_primep(n):=if n=1 then false else block(sample:expansion(n)-(x^n-1),
makelist(numfactor(part(sample,i)),i,1,length(sample)),
if not member(false,makelist(integerp(%%[i]/n),i,1,length(%%))) then true)$
/* Test cases */
makelist(expansion(i),i,0,7);
block(
makelist(aks_primep(i),i,1,50),
sublist_indices(%%,lambda([x],x=true)));
</syntaxhighlight>
{{out}}
<pre>
[1,x-1,x^2-2*x+1,x^3-3*x^2+3*x-1,x^4-4*x^3+6*x^2-4*x+1,x^5-5*x^4+10*x^3-10*x^2+5*x-1,x^6-6*x^5+15*x^4-20*x^3+15*x^2-6*x+1,x^7-7*x^6+21*x^5-35*x^4+35*x^3-21*x^2+7*x-1]
[2,3,5,7,11,13,17,19,23,29,31,37,41,43,47]
</pre>
=={{header|Nim}}==
<syntaxhighlight lang="nim">
from math import binom
import strutils
# Table of unicode superscript characters.
const Exponents: array[0..9, string] = ["⁰", "¹", "²", "³", "⁴", "⁵", "⁶", "⁷", "⁸", "⁹"]
iterator coeffs(n: int): int =
## Yield the coefficients of the expansion of (x - 1)ⁿ.
var sign = 1
for k in 0..n:
yield binom(n, k) * sign
sign = -sign
iterator polyExpansion(n: int): tuple[c, e: int] =
## Yield the coefficients and the exponents of the expansion of (x - 1)ⁿ.
var e = n
for c in coeffs(n):
yield(c, e)
dec e
proc termString(c, e: int): string =
## Return the string for the term c * e^n.
if e == 0:
result.addInt(c)
else:
if c != 1:
result.addInt(c)
result.add('x')
if e != 1:
result.add(Exponents[e])
proc polyString(n: int): string =
## Return the string for the expansion of (x - 1)ⁿ.
for (c, e) in polyExpansion(n):
if c < 0:
result.add(" - ")
elif e != n:
result.add(" + ")
result.add(termString(abs(c), e))
proc isPrime(n: int): bool =
## Check if a number is prime using the polynome expansion.
result = true
for (c, e) in polyExpansion(n):
if e in 1..(n-1): # xⁿ and 1 are eliminated by the subtraction.
if c mod n != 0:
return false
#---------------------------------------------------------------------------------------------------
echo "Polynome expansions:"
for p in 0..9:
echo "(x - 1)$1 = $2".format(Exponents[p], polyString(p))
var primes: string
for p in 2..34:
if p.isPrime():
primes.addSep(", ", 0)
primes.addInt(p)
echo "\nPrimes under 35: ", primes
for p in 35..50:
if p.isPrime():
primes.add(", ")
primes.addInt(p)
echo "\nPrimes under 50: ", primes
</syntaxhighlight>
{{out}}
<pre>
Polynome expansions:
(x - 1)⁰ = 1
(x - 1)¹ = x - 1
(x - 1)² = x² - 2x + 1
(x - 1)³ = x³ - 3x² + 3x - 1
(x - 1)⁴ = x⁴ - 4x³ + 6x² - 4x + 1
(x - 1)⁵ = x⁵ - 5x⁴ + 10x³ - 10x² + 5x - 1
(x - 1)⁶ = x⁶ - 6x⁵ + 15x⁴ - 20x³ + 15x² - 6x + 1
(x - 1)⁷ = x⁷ - 7x⁶ + 21x⁵ - 35x⁴ + 35x³ - 21x² + 7x - 1
(x - 1)⁸ = x⁸ - 8x⁷ + 28x⁶ - 56x⁵ + 70x⁴ - 56x³ + 28x² - 8x + 1
(x - 1)⁹ = x⁹ - 9x⁸ + 36x⁷ - 84x⁶ + 126x⁵ - 126x⁴ + 84x³ - 36x² + 9x - 1
Primes under 35: 2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31
Primes under 50: 2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43, 47
</pre>
=={{header|Objeck}}==
{{trans|Java}}
<
@c : static : Int[];
Line 2,462 ⟶ 3,318:
} while (n-- <> 0);
}
}</
Output:
Line 2,483 ⟶ 3,339:
{{trans|Clojure}}
Uses [http://github.com/c-cube/gen gen] library for lazy streams and [https://forge.ocamlcore.org/projects/zarith zarith] for arbitrarily sized integers. Runs as is through the [https://github.com/diml/utop utop] REPL.
<
#require "zarith"
open Z
Line 2,524 ⟶ 3,380:
print_endline ("primes < 50 per AKS: " ^
(Gen.filter aks (range (of_int 2) (of_int 50)) |>
Gen.map to_string |> Gen.to_list |> String.concat " "))</
{{output}}
Line 2,540 ⟶ 3,396:
primes < 50 per AKS: (2 3 5 7 11 13 17 19 23 29 31 37 41 43 47)</pre>
=={{header|Oforth}}==
<
: nextCoef( prev -- [] )
Line 2,560 ⟶ 3,417:
0 10 for: i [ System.Out "(x-1)^" << i << " = " << coefs( i ) << cr ]
50 seq filter( #prime? ) apply(#.) printcr
;</
{{out}}
Line 2,579 ⟶ 3,436:
=={{header|PARI/GP}}==
<
vector(8,n,getPoly(n-1))
AKS_slow(n)=my(P=getPoly(n));for(i=1,n-1,if(polcoeff(P,i)%n,return(0))); 1;
AKS(n)=my(X=('x-1)*Mod(1,n));X^n=='x^n-1;
select(AKS, [1..50])</
{{out}}
<pre> [1, x - 1, x^2 - 2*x + 1, x^3 - 3*x^2 + 3*x - 1, x^4 - 4*x^3 + 6*x^2 - 4*x + 1, x^5 - 5*x^4 + 10*x^3 - 10*x^2 + 5*x - 1, x^6 - 6*x^5 + 15*x^4 - 20*x^3 + 15*x^2 - 6*x + 1, x^7 - 7*x^6 + 21*x^5 - 35*x^4 + 35*x^3 - 21*x^2 + 7*x - 1]
Line 2,590 ⟶ 3,447:
=={{header|Pascal}}==
{{works with|Free Pascal}}
<
const
pasTriMax = 61;
Line 2,683 ⟶ 3,540:
Write(n:3);
WriteLn;
end.</
;output:
<pre>
Line 2,699 ⟶ 3,556:
=={{header|Perl}}==
<
use warnings;
# Select one of these lines. Math::BigInt is in core, but quite slow.
Line 2,727 ⟶ 3,584:
print "expansions of (x-1)^p:\n";
print binpoly($_),"\n" for 0..9;
print "Primes to 80: [", join(",", grep { binprime($_) } 2..80), "]\n";</
{{out}}
<pre>expansions of (x-1)^p:
Line 2,746 ⟶ 3,603:
{{libheader|ntheory}}
<
# Uncomment next line to see the r and s values used. Set to 2 for more detail.
# prime_set_config(verbose => 1);
say join(" ", grep { is_aks_prime($_) } 1_000_000_000 .. 1_000_000_100);</
{{out}}
<pre>1000000007 1000000009 1000000021 1000000033 1000000087 1000000093 1000000097</pre>
=={{header|Phix}}==
<!--<syntaxhighlight lang="phix">-->
<span style="color: #000080;font-style:italic;">-- demo/rosetta/AKSprimes.exw
-- Does not work for primes above 53, which is actually beyond the original task anyway.
-- Translated from the C version, just about everything is (working) out-by-1, what fun.</span>
<span style="color: #004080;">sequence</span> <span style="color: #000000;">c</span> <span style="color: #0000FF;">=</span> <span style="color: #7060A8;">repeat</span><span style="color: #0000FF;">(</span><span style="color: #000000;">0</span><span style="color: #0000FF;">,</span><span style="color: #000000;">100</span><span style="color: #0000FF;">)</span>
<span style="color: #008080;">procedure</span> <span style="color: #000000;">coef</span><span style="color: #0000FF;">(</span><span style="color: #004080;">integer</span> <span style="color: #000000;">n</span><span style="color: #0000FF;">)</span>
<span style="color: #000080;font-style:italic;">-- out-by-1, ie coef(1)==^0, coef(2)==^1, coef(3)==^2 etc.</span>
<span style="color: #000000;">c</span><span style="color: #0000FF;">[</span><span style="color: #000000;">n</span><span style="color: #0000FF;">]</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">1</span>
<span style="color: #008080;">for</span> <span style="color: #000000;">i</span><span style="color: #0000FF;">=</span><span style="color: #000000;">n</span><span style="color: #0000FF;">-</span><span style="color: #000000;">1</span> <span style="color: #008080;">to</span> <span style="color: #000000;">2</span> <span style="color: #008080;">by</span> <span style="color: #0000FF;">-</span><span style="color: #000000;">1</span> <span style="color: #008080;">do</span>
<span style="color: #000000;">c</span><span style="color: #0000FF;">[</span><span style="color: #000000;">i</span><span style="color: #0000FF;">]</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">c</span><span style="color: #0000FF;">[</span><span style="color: #000000;">i</span><span style="color: #0000FF;">]+</span><span style="color: #000000;">c</span><span style="color: #0000FF;">[</span><span style="color: #000000;">i</span><span style="color: #0000FF;">-</span><span style="color: #000000;">1</span><span style="color: #0000FF;">]</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">for</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">procedure</span>
<span style="color: #008080;">function</span> <span style="color: #000000;">is_aks_prime</span><span style="color: #0000FF;">(</span><span style="color: #004080;">integer</span> <span style="color: #000000;">n</span><span style="color: #0000FF;">)</span>
<span style="color: #000000;">coef</span><span style="color: #0000FF;">(</span><span style="color: #000000;">n</span><span style="color: #0000FF;">+</span><span style="color: #000000;">1</span><span style="color: #0000FF;">);</span> <span style="color: #000080;font-style:italic;">-- (I said it was out-by-1)</span>
<span style="color: #008080;">for</span> <span style="color: #000000;">i</span><span style="color: #0000FF;">=</span><span style="color: #000000;">2</span> <span style="color: #008080;">to</span> <span style="color: #000000;">n</span><span style="color: #0000FF;">-</span><span style="color: #000000;">1</span> <span style="color: #008080;">do</span> <span style="color: #000080;font-style:italic;">-- (technically "to n" is more correct)</span>
<span style="color: #008080;">if</span> <span style="color: #7060A8;">remainder</span><span style="color: #0000FF;">(</span><span style="color: #000000;">c</span><span style="color: #0000FF;">[</span><span style="color: #000000;">i</span><span style="color: #0000FF;">],</span><span style="color: #000000;">n</span><span style="color: #0000FF;">)!=</span><span style="color: #000000;">0</span> <span style="color: #008080;">then</span>
<span style="color: #008080;">return</span> <span style="color: #000000;">0</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">if</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">for</span>
<span style="color: #008080;">return</span> <span style="color: #000000;">1</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">function</span>
<span style="color: #008080;">procedure</span> <span style="color: #000000;">show</span><span style="color: #0000FF;">(</span><span style="color: #004080;">integer</span> <span style="color: #000000;">n</span><span style="color: #0000FF;">)</span>
<span style="color: #000080;font-style:italic;">-- (As per coef, this is (working) out-by-1)</span>
<span style="color: #004080;">object</span> <span style="color: #000000;">ci</span>
<span style="color: #008080;">for</span> <span style="color: #000000;">i</span><span style="color: #0000FF;">=</span><span style="color: #000000;">n</span> <span style="color: #008080;">to</span> <span style="color: #000000;">1</span> <span style="color: #008080;">by</span> <span style="color: #0000FF;">-</span><span style="color: #000000;">1</span> <span style="color: #008080;">do</span>
<span style="color: #000000;">ci</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">c</span><span style="color: #0000FF;">[</span><span style="color: #000000;">i</span><span style="color: #0000FF;">]</span>
<span style="color: #008080;">if</span> <span style="color: #000000;">ci</span><span style="color: #0000FF;">=</span><span style="color: #000000;">1</span> <span style="color: #008080;">then</span>
<span style="color: #008080;">if</span> <span style="color: #7060A8;">remainder</span><span style="color: #0000FF;">(</span><span style="color: #000000;">n</span><span style="color: #0000FF;">-</span><span style="color: #000000;">i</span><span style="color: #0000FF;">,</span><span style="color: #000000;">2</span><span style="color: #0000FF;">)=</span><span style="color: #000000;">0</span> <span style="color: #008080;">then</span>
<span style="color: #008080;">if</span> <span style="color: #000000;">i</span><span style="color: #0000FF;">=</span><span style="color: #000000;">1</span> <span style="color: #008080;">then</span>
<span style="color: #008080;">if</span> <span style="color: #000000;">n</span><span style="color: #0000FF;">=</span><span style="color: #000000;">1</span> <span style="color: #008080;">then</span>
<span style="color: #000000;">ci</span> <span style="color: #0000FF;">=</span> <span style="color: #008000;">"1"</span>
<span
<span style="color: #000000;">ci</span> <span style="color: #0000FF;">=</span> <span style="color: #008000;">"+1"</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">if</span>
<span
<span style="color: #000000;">ci</span> <span style="color: #0000FF;">=</span> <span style="color: #008000;">""</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">if</span>
<span
<span style="color: #000000;">ci</span> <span style="color: #0000FF;">=</span> <span style="color: #008000;">"-1"</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">if</span>
<span style="color: #008080;">else</span>
<span style="color: #008080;">if</span> <span style="color: #7060A8;">remainder</span><span style="color: #0000FF;">(</span><span style="color: #000000;">n</span><span style="color: #0000FF;">-</span><span style="color: #000000;">i</span><span style="color: #0000FF;">,</span><span style="color: #000000;">2</span><span style="color: #0000FF;">)=</span><span style="color: #000000;">0</span> <span style="color: #008080;">then</span>
<span style="color: #000000;">ci</span> <span style="color: #0000FF;">=</span> <span style="color: #7060A8;">sprintf</span><span style="color: #0000FF;">(</span><span style="color: #008000;">"+%d"</span><span style="color: #0000FF;">,</span><span style="color: #000000;">ci</span><span style="color: #0000FF;">)</span>
<span
<span style="color: #000000;">ci</span> <span style="color: #0000FF;">=</span> <span style="color: #7060A8;">sprintf</span><span style="color: #0000FF;">(</span><span style="color: #008000;">"-%d"</span><span style="color: #0000FF;">,</span><span style="color: #000000;">ci</span><span style="color: #0000FF;">)</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">if</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">if</span>
<span style="color: #008080;">if</span> <span style="color: #000000;">i</span><span style="color: #0000FF;">=</span><span style="color: #000000;">1</span> <span style="color: #008080;">then</span> <span style="color: #000080;font-style:italic;">-- ie ^0</span>
<span style="color: #7060A8;">printf</span><span style="color: #0000FF;">(</span><span style="color: #000000;">1</span><span style="color: #0000FF;">,</span><span style="color: #008000;">"%s"</span><span style="color: #0000FF;">,{</span><span style="color: #000000;">ci</span><span style="color: #0000FF;">})</span>
<span style="color: #008080;">elsif</span> <span style="color: #000000;">i</span><span style="color: #0000FF;">=</span><span style="color: #000000;">2</span> <span style="color: #008080;">then</span> <span style="color: #000080;font-style:italic;">-- ie ^1</span>
<span style="color: #7060A8;">printf</span><span style="color: #0000FF;">(</span><span style="color: #000000;">1</span><span style="color: #0000FF;">,</span><span style="color: #008000;">"%sx"</span><span style="color: #0000FF;">,{</span><span style="color: #000000;">ci</span><span style="color: #0000FF;">})</span>
<span style="color:
<span style="color: #7060A8;">printf</span><span style="color: #0000FF;">(</span><span style="color: #000000;">1</span><span style="color: #0000FF;">,</span><span style="color: #008000;">"%sx^%d"</span><span style="color: #0000FF;">,{</span><span style="color: #000000;">ci</span><span style="color: #0000FF;">,</span><span style="color: #000000;">i</span><span style="color: #0000FF;">-</span><span style="color: #000000;">1</span><span style="color: #0000FF;">})</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">if</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">for</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">procedure</span>
<span style="color: #008080;">procedure</span> <span style="color: #000000;">main</span><span style="color: #0000FF;">()</span>
<span style="color: #008080;">for</span> <span style="color: #000000;">n</span><span style="color: #0000FF;">=</span><span style="color: #000000;">1</span> <span style="color: #008080;">to</span> <span style="color: #000000;">10</span> <span style="color: #008080;">do</span> <span style="color: #000080;font-style:italic;">-- (0 to 9 really)</span>
<span style="color: #000000;">coef</span><span style="color: #0000FF;">(</span><span style="color: #000000;">n</span><span style="color: #0000FF;">);</span>
<span style="color: #7060A8;">printf</span><span style="color: #0000FF;">(</span><span style="color: #000000;">1</span><span style="color: #0000FF;">,</span><span style="color: #008000;">"(x-1)^%d = "</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">n</span><span style="color: #0000FF;">-</span><span style="color: #000000;">1</span><span style="color: #0000FF;">);</span>
<span style="color: #000000;">show</span><span style="color: #0000FF;">(</span><span style="color: #000000;">n</span><span style="color: #0000FF;">);</span>
<span style="color: #7060A8;">puts</span><span style="color: #0000FF;">(</span><span style="color: #000000;">1</span><span style="color: #0000FF;">,</span><span style="color: #008000;">'\n'</span><span style="color: #0000FF;">);</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">for</span>
<span style="color: #7060A8;">puts</span><span style="color: #0000FF;">(</span><span style="color: #000000;">1</span><span style="color: #0000FF;">,</span><span style="color: #008000;">"\nprimes (<=53):"</span><span style="color: #0000FF;">);</span>
<span style="color: #000080;font-style:italic;">-- coef(2); -- (needed to reset c, if we want to avoid saying 1 is prime...)</span>
<span style="color: #000000;">c</span><span style="color: #0000FF;">[</span><span style="color: #000000;">2</span><span style="color: #0000FF;">]</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">1</span> <span style="color: #000080;font-style:italic;">-- (this manages "", which is all that call did anyway...)</span>
<span style="color: #008080;">for</span> <span style="color: #000000;">n</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">2</span> <span style="color: #008080;">to</span> <span style="color: #000000;">53</span> <span style="color: #008080;">do</span>
<span style="color: #008080;">if</span> <span style="color: #000000;">is_aks_prime</span><span style="color: #0000FF;">(</span><span style="color: #000000;">n</span><span style="color: #0000FF;">)</span> <span style="color: #008080;">then</span>
<span style="color: #7060A8;">printf</span><span style="color: #0000FF;">(</span><span style="color: #000000;">1</span><span style="color: #0000FF;">,</span><span style="color: #008000;">" %d"</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">n</span><span style="color: #0000FF;">);</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">if</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">for</span>
<span style="color: #7060A8;">puts</span><span style="color: #0000FF;">(</span><span style="color: #000000;">1</span><span style="color: #0000FF;">,</span><span style="color: #008000;">'\n'</span><span style="color: #0000FF;">);</span>
<span style="color: #008080;">if</span> <span style="color: #7060A8;">getc</span><span style="color: #0000FF;">(</span><span style="color: #000000;">0</span><span style="color: #0000FF;">)</span> <span style="color: #008080;">then</span> <span style="color: #008080;">end</span> <span style="color: #008080;">if</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">procedure</span>
<span style="color: #000000;">main</span><span style="color: #0000FF;">()</span>
<!--</syntaxhighlight>-->
{{out}}
<pre>
Line 2,910 ⟶ 3,707:
primes (<=53): 2 3 5 7 11 13 17 19 23 29 31 37 41 43 47 53
</pre>
=={{header|Picat}}==
<syntaxhighlight lang="picat">
pascal([]) = [1].
pascal(L) = [1|sum_adj(L)].
sum_adj(Row) = Next =>
Next = L,
while (Row = [A,B|_])
L = [A+B|Rest],
L := Rest,
Row := tail(Row)
end,
L = Row.
show_x(0) = "".
show_x(1) = "x".
show_x(N) = S, N > 1 => S = [x, '^' | to_string(N)].
show_term(Coef, Exp) = cond((Coef != 1; Exp == 0), Coef.to_string, "") ++ show_x(Exp).
expansions(N) =>
Row = [],
foreach (I in 0..N-1)
Row := pascal(Row),
writef("(x - 1)^%d = ", I),
Exp = I,
Sgn = '+',
foreach (Coef in Row)
if Exp != I then
writef(" %w ", Sgn)
end,
writef("%s", show_term(Coef, Exp)),
Exp := Exp - 1,
Sgn := cond(Sgn == '+', '-', '+')
end,
nl
end.
primerow([], _).
primerow([A,A|_], N) :- A mod N == 0. % end when we've seen half the list.
primerow([A|As], N) :- (A mod N == 0; A == 1), primerow(As, N).
primes_upto(N) = Primes =>
Primes = L,
Row = [1, 1],
foreach (K in 2..N)
Row := pascal(Row),
if primerow(Row, K) then
L = [K|Rest],
L := Rest
end
end,
L = [].
main =>
expansions(8),
writef("%nThe primes upto 50 (via AKS) are: %w%n", primes_upto(50)).
</syntaxhighlight>
{{Out}}
<pre>
(x - 1)^0 = 1
(x - 1)^1 = x - 1
(x - 1)^2 = x^2 - 2x + 1
(x - 1)^3 = x^3 - 3x^2 + 3x - 1
(x - 1)^4 = x^4 - 4x^3 + 6x^2 - 4x + 1
(x - 1)^5 = x^5 - 5x^4 + 10x^3 - 10x^2 + 5x - 1
(x - 1)^6 = x^6 - 6x^5 + 15x^4 - 20x^3 + 15x^2 - 6x + 1
(x - 1)^7 = x^7 - 7x^6 + 21x^5 - 35x^4 + 35x^3 - 21x^2 + 7x - 1
The primes upto 50 (via AKS) are: [2,3,5,7,11,13,17,19,23,29,31,37,41,43,47]
</pre>
=={{header|PicoLisp}}==
<syntaxhighlight lang="text">(de pascal (N)
(let D 1
(make
Line 2,932 ⟶ 3,801:
(range 2 50) ) )
(bye)</
{{out}}
Line 2,950 ⟶ 3,819:
=={{header|PL/I}}==
<syntaxhighlight lang="pl/i">
AKS: procedure options (main, reorder); /* 16 September 2015, derived from Fortran */
Line 3,086 ⟶ 3,955:
end AKS;
</syntaxhighlight>
Results obtained:
<pre>
Line 3,108 ⟶ 3,977:
connection can be expressed directly in Prolog by the following prime number
generator:
<
prime(P) :-
pascal([1,P|Xs]),
append(Xs, [1], Rest),
forall( member(X,Xs), 0 is X mod P).
</syntaxhighlight>
where pascal/1 is a generator of rows of the Pascal triangle, for
example as defined below; the other predicates used above are
Line 3,147 ⟶ 4,016:
===Pascal Triangle Generator===
<
% To generate the n-th row of a Pascal triangle
% pascal(+N, Row)
Line 3,191 ⟶ 4,060:
S is X1 + X2,
add_pairs( [X2, X3|Xs], Ys).
</syntaxhighlight>
===Solutions===
Solutions with output from SWI-Prolog:
<
%%% Task 1: "A method to generate the coefficients of (1-X)^p"
Line 3,219 ⟶ 4,088:
% As required by the problem statement, but necessarily very inefficient:
:- between(0, 7, N), coefficients(N, Coefficients), writeln(Coefficients), fail ; true.
</syntaxhighlight>
<pre>
[1]
Line 3,238 ⟶ 4,107:
</pre>
<
%%% Task 3. Use the previous function in creating [sic]
%%% another function that when given p returns whether p is prime
Line 3,250 ⟶ 4,119:
append(Cs, [_], Coefficients),
forall( member(C, Cs), 0 is C mod N).
</syntaxhighlight>
The following is more efficient (because it relies on optpascal
Line 3,257 ⟶ 4,126:
recomputation):
<
prime(N) :- optpascal([1,N|Xs]), forall( member(X,Xs), 0 is X mod N).
</syntaxhighlight>
<
%%% Task 4. Use your AKS test to generate a list of all primes under 35.
Line 3,273 ⟶ 4,142:
% Output: 1 2 3 5 7 11 13 17 19 23 29 31 37 41 43 47
</syntaxhighlight>
=== Alternate version using symbolic computation ===
This version uses Prolog's symbolic pattern matching to generate the polynomial equations that are the subject of the first part of the task. The rows of Pascal's triangle are generated, then combined with the "a" and "b" terms of the expression (a + b)^n. The expression is then algebraically simplified to remove things like x^0, x^1, 1*x, a + 0, etc.
It would be possible to go further and expand (a + b)^n and collect terms, but that turns out to be rather difficult because terms have to be reordered into a canonical form first -- for example, an expansion can produce both j*a*b and k*b*a, and these terms would have to be reordered to get (j+k)*a*b.
<syntaxhighlight lang="prolog">
main :- task1(8), nl, task2(50), halt.
task1(N) :-
pascal(Z),
length(Rows, N),
prefix(Rows, Z),
forall(member(Row, Rows),
(length(Row, K), succ(DecK, K),
binomial(x, -1, Row, Expr),
format("(x-1)**~w = ~w~n", [DecK, Expr]))).
task2(Upto) :-
primes_upto(Upto, Ps),
format("The primes upto ~w (via AKS) are: ~p~n", [Upto, Ps]).
pascal(Lz) :-
lazy_list(pascal_row, [], Lz).
pascal_row([], R1, R1) :- R1 = [1], !.
pascal_row(R0, R1, R1) :-
sum_adj(R0, Next), R1 = [1|Next].
sum_adj(L, L) :- L = [_], !.
sum_adj([A|As], [C|Cs]) :-
As = [B|_], C is A + B,
sum_adj(As, Cs).
% First part of task -- create textual representation of (x-1)^n
% here we generate expression trees
%
binomial(A, B, Coefs, Expr) :-
length(Coefs, N), succ(DecN, N),
binomial(B, DecN, A, 0, Coefs, Exp0),
reduce(Exp0, Exp1),
addition_to_subtraction(Exp1, Expr).
binomial(_, _, _, _, [], 0) :- !.
binomial(A, PowA, B, PowB, [N|Ns], Ts + T) :-
T = N * A**PowA * B**PowB,
IncPow is PowB + 1,
DecPow is PowA - 1,
binomial(A, DecPow, B, IncPow, Ns, Ts).
addition_to_subtraction(A + B, X) :-
addition_to_subtraction(A, C),
(make_positive(B, D) -> X = C - D; X = C + B), !.
addition_to_subtraction(X, X).
make_positive(N, Term) :- integer(N), N < 0, !, Term is -N.
make_positive(A*B, Term) :-
make_positive(A, PosA),
(PosA = 1 -> Term = B, !; Term = PosA*B).
reduce(A, C) :-
simplify(A, B),
(B = A -> C = A; reduce(B, C)).
simplify(_**0, 1) :- !.
simplify(1**_, 1) :- !.
simplify(-1**N, Z) :- integer(N), (0 is N /\ 1 -> Z = 1; Z = -1), !.
simplify(X**1, X) :- !.
simplify(0 + A, A) :- !.
simplify(A + 0, A) :- !.
simplify(A + B, C) :-
integer(A),
integer(B), !,
C is A + B.
simplify(A + B, C + D) :- !,
simplify(A, C),
simplify(B, D).
simplify(0 * _, 0) :- !.
simplify(_ * 0, 0) :- !.
simplify(1 * A, A) :- !.
simplify(A * 1, A) :- !.
simplify(A * B, C) :-
integer(A),
integer(B), !,
C is A * B.
simplify(A * B, C * D) :- !,
simplify(A, C),
simplify(B, D).
simplify(X, X).
% Second part of task -- Use the coefficients of Pascal's Triangle to check primality.
%
primerow([1, N| Rest]) :- primerow(N, Rest).
primerow(_, End) :- (End = []; End = [1]), !.
primerow(_, [A,A|_]) :- !. % end when we've seen half the list.
primerow(N, [A|As]) :- A mod N =:= 0, primerow(N, As).
second([_,N|_], N).
primes_upto(N, Ps) :-
pascal(Z),
Z = [_, _ | Rows], % we only care about 2nd row on up. ([1,2,1])
succ(DecN, N), length(CheckRows, DecN), prefix(CheckRows, Rows),
include(primerow, CheckRows, PrimeRows),
maplist(second, PrimeRows, Ps).
?- main.
</syntaxhighlight>
{{Out}}
<pre>
$ swipl aks.pl
(x-1)**0 = 1
(x-1)**1 = x-1
(x-1)**2 = x**2-2*x+1
(x-1)**3 = x**3-3*x**2+3*x-1
(x-1)**4 = x**4-4*x**3+6*x**2-4*x+1
(x-1)**5 = x**5-5*x**4+10*x**3-10*x**2+5*x-1
(x-1)**6 = x**6-6*x**5+15*x**4-20*x**3+15*x**2-6*x+1
(x-1)**7 = x**7-7*x**6+21*x**5-35*x**4+35*x**3-21*x**2+7*x-1
The primes upto 50 (via AKS) are: [2,3,5,7,11,13,17,19,23,29,31,37,41,43,47]
</pre>
=={{header|PureBasic}}==
<
Define vzr.b = -1, vzc.b = ~vzr, nMAX.i = 10, n.i , k.i
Line 3,319 ⟶ 4,312:
If isPrime(n,pd()) : Print(Str(n)+Space(2)) : EndIf
Next
Input()</
{{out}}
<pre> 0: + 1
Line 3,338 ⟶ 4,331:
=={{header|Python}}==
<
# This version uses a generator and thus less computations
c =1
Line 3,353 ⟶ 4,346:
# we stop without computing all possible solutions
return False
return True</
or equivalently:
<
if p==2:return True
c=1
Line 3,361 ⟶ 4,354:
c=c*(p-i)//(i+1)
if c%p:return False
return True</
alternatively:
<
ex = [1]
for i in range(p):
Line 3,382 ⟶ 4,375:
print('\n# small primes using the aks test')
print([p for p in range(101) if aks_test(p)])</
{{out}}
Line 3,405 ⟶ 4,398:
Using a wikitable and math features with the following additional code produces better formatted polynomial output:
<
{| class="wikitable" style="text-align:left;"
|+ Polynomial Expansions and AKS prime test
Line 3,419 ⟶ 4,412:
for n,e in enumerate(expand_x_1(p))),
aks_test(p)))
print('|}')</
{{out}}
Line 3,479 ⟶ 4,472:
|-
|}
=={{header|Quackery}}==
<syntaxhighlight lang="Quackery"> [ dup 0 peek swap
[] 0 rot 0 join
witheach
[ abs tuck +
rot join swap ]
drop
[] swap witheach
[ rot negate
dup dip unrot
* join ]
nip ] is nextcoeffs ( [ --> [ )
[ ' [ 1 ] swap
times nextcoeffs ] is coeffs ( n --> [ )
[ dup 2 < iff
[ drop false ]
done
true swap
dup coeffs
behead drop
-1 split drop
witheach
[ over mod
0 != if
[ dip not
conclude ] ]
drop ] is prime ( n --> b )
[ behead
0 < iff
[ say "-" ]
else sp
dup size
dup 0 = iff
[ 1 echo 2drop ]
done
say "x"
dup 1 = iff
[ 2drop
say " + 1" ]
done
say "^"
echo
witheach
[ dup 0 < iff
[ say " - " ]
else
[ say " + " ]
abs echo
i 0 = if done
say "x"
i 1 = if done
say "^"
i echo ] ] is echocoeffs ( [ --> )
8 times
[ i^ echo
say ": "
i^ coeffs
echocoeffs
cr ]
cr
50 times
[ i^ prime if
[ i^ echo sp ] ]</syntaxhighlight>
{{out}}
<pre>0: 1
1: -x + 1
2: x^2 - 2x + 1
3: -x^3 + 3x^2 - 3x + 1
4: x^4 - 4x^3 + 6x^2 - 4x + 1
5: -x^5 + 5x^4 - 10x^3 + 10x^2 - 5x + 1
6: x^6 - 6x^5 + 15x^4 - 20x^3 + 15x^2 - 6x + 1
7: -x^7 + 7x^6 - 21x^5 + 35x^4 - 35x^3 + 21x^2 - 7x + 1
2 3 5 7 11 13 17 19 23 29 31 37 41 43 47 </pre>
=={{header|R}}==
Line 3,485 ⟶ 4,560:
<
i<-2:p-1
l<-unique(factorial(p) / (factorial(p-i) * factorial(i)))
Line 3,493 ⟶ 4,568:
print(noquote("It isn't prime."))
}
}</
=={{header|Racket}}==
With copious use of the math/number-theory library...
<
(require math/number-theory)
Line 3,559 ⟶ 4,634:
(displayln (for/list ((i (in-range lowest-tested-number 50)) #:when (prime?/AKS i)) i))
(displayln (for/list ((i (in-range lowest-tested-number 100)) #:when (prime?/AKS i)) i))
</syntaxhighlight>
{{out}}
Line 3,577 ⟶ 4,652:
(2 3 5 7 11 13 17 19 23 29 31 37 41 43 47)
(2 3 5 7 11 13 17 19 23 29 31 37 41 43 47 53 59 61 67 71 73 79 83 89 97)</pre>
=={{header|Raku}}==
(formerly Perl 6)
{{works with|Rakudo|2015-12}}
The expansions are generated similarly to how most FP languages generate sequences that resemble Pascal's triangle, using a zipwith meta-operator (Z) with subtraction, applied between two lists that add a 0 on either end to the prior list. Here we define a constant infinite sequence using the <tt>...</tt> sequence operator with a "whatever" endpoint. In fact, the second term <tt>[1,-1]</tt> could have been generated from the first term, but we put it in there for documentation so the reader can see what direction things are going.
The <tt>polyprime</tt> function pretty much reads like the original description. Is it "so" that the p'th expansion's coefficients are all divisible by p? The <tt>.[1 ..^ */2]</tt> slice is done simply to weed out divisions by 1 or by factors we've already tested (since the coefficients are symmetrical in terms of divisibility). If we wanted to write <tt>polyprime</tt> even more idiomatically, we could have made it another infinite constant list that is just a mapping of the first list, but we decided that would just be showing off. <tt>:-)</tt>
<syntaxhighlight lang="raku" line>constant expansions = [1], [1,-1], -> @prior { [|@prior,0 Z- 0,|@prior] } ... *;
sub polyprime($p where 2..*) { so expansions[$p].[1 ..^ */2].all %% $p }
# Showing the expansions:
say ' p: (x-1)ᵖ';
say '-----------';
sub super ($n) {
$n.trans: '0123456789'
=> '⁰¹²³⁴⁵⁶⁷⁸⁹';
}
for ^13 -> $d {
say $d.fmt('%2i: '), (
expansions[$d].kv.map: -> $i, $n {
my $p = $d - $i;
[~] gather {
take < + - >[$n < 0] ~ ' ' unless $p == $d;
take $n.abs unless $p == $d > 0;
take 'x' if $p > 0;
take super $p - $i if $p > 1;
}
}
)
}
# And testing the function:
print "\nPrimes up to 100:\n { grep &polyprime, 2..100 }\n";</syntaxhighlight>
{{out}}
<pre> p: (x-1)ᵖ
-----------
0: 1
1: x - 1
2: x² - 2x + 1
3: x³ - 3x² + 3x - 1
4: x⁴ - 4x³ + 6x² - 4x + 1
5: x⁵ - 5x⁴ + 10x³ - 10x² + 5x - 1
6: x⁶ - 6x⁵ + 15x⁴ - 20x³ + 15x² - 6x + 1
7: x⁷ - 7x⁶ + 21x⁵ - 35x⁴ + 35x³ - 21x² + 7x - 1
8: x⁸ - 8x⁷ + 28x⁶ - 56x⁵ + 70x⁴ - 56x³ + 28x² - 8x + 1
9: x⁹ - 9x⁸ + 36x⁷ - 84x⁶ + 126x⁵ - 126x⁴ + 84x³ - 36x² + 9x - 1
10: x¹⁰ - 10x⁹ + 45x⁸ - 120x⁷ + 210x⁶ - 252x⁵ + 210x⁴ - 120x³ + 45x² - 10x + 1
11: x¹¹ - 11x¹⁰ + 55x⁹ - 165x⁸ + 330x⁷ - 462x⁶ + 462x⁵ - 330x⁴ + 165x³ - 55x² + 11x - 1
12: x¹² - 12x¹¹ + 66x¹⁰ - 220x⁹ + 495x⁸ - 792x⁷ + 924x⁶ - 792x⁵ + 495x⁴ - 220x³ + 66x² - 12x + 1
Primes up to 100:
2 3 5 7 11 13 17 19 23 29 31 37 41 43 47 53 59 61 67 71 73 79 83 89 97</pre>
=={{header|REXX}}==
===version 1===
<
* 09.02.2014 Walter Pachl
* 22.02.2014 WP fix 'accounting' problem (courtesy GS)
Line 3,628 ⟶ 4,763:
Say ' ' subword(pl,29)
Say 'Largest coefficient:' mmm
Say 'This has' length(mmm) 'digits' </
{{out}}
<pre>(x-1)**0 = 1
Line 3,649 ⟶ 4,784:
The program determines programmatically the required number of decimal digits (precision) for the large coefficients.
<
parse arg Z . /*obtain optional argument from the CL.*/
if Z=='' | Z=="," then Z= 200 /*Not specified? Then use the default.*/
Line 3,691 ⟶ 4,826:
if #=='' then #= "none"; return # /*if null, return "none"; else return #*/
/*──────────────────────────────────────────────────────────────────────────────────────*/
isAKSp: if z==word(#,words(#)) then return ' is a prime.'; else return " isn't a prime."</
{{out|output|text= for task requirement #2, showing thirty-one expressions using as input: <tt> -31 </tt>}}
<br>(Shown at five-sixth size.)
Line 3,757 ⟶ 4,892:
Found 95 primes and the largest coefficient has 150 decimal digits.
</pre>
=={{header|RPL}}==
{| class="wikitable"
! RPL code
! Comment
|-
|
≪ → p
≪ { } 1
0 p '''FOR''' j
p j COMB OVER * ROT SWAP +
SWAP NEG
'''NEXT''' DROP
≫ ≫ ‘'''XPM1'''’ STO
≪ DUP SIZE → coeffs size
≪ 0 1 size '''FOR''' j
coeffs j GET 'X' size j - ^ * + '''NEXT'''
≫ ≫ ‘'''SHOWA'''’ STO
≪ DUP '''XPM1''' → p coeffs
≪ 1
2 p '''FOR''' j
coeffs j GET ABS p MOD NOT AND '''NEXT'''
≫ ≫ ‘'''PRIM?'''’ STO
|
'''XPM1''' ''( p -- { (x-1)^p_coefficients } ) ''
coeffs = { } ; sign = 1
loop for j=0 to p
append C(p,j)*sign to coeffs
sign = -sign
end loop
return coeffs
'''SHOWA''' ''( { coeffs } -- 'polynom' )''
loop for each coeff
append jth polynomial term
return polynom
'''PRIM?''' ''( p -- boolean ) ''
res = true
loop for j=2 to p
res = res AND (mod(coeffs[j],p)=0)
return res
|}
{{in}}
<pre>
≪ { } 1 7 FOR n n XPM1 SHOWA + NEXT ≫
≪ { } 2 35 FOR n IF n PRIM? THEN n + END NEXT ≫
</pre>
{{out}}
<pre>
8: '-1+X'
7: '1-2*X+X^2'
6: '-1+3*X-3*X^2+X^3'
5: '1-4*X-4*X^3+6*X^2+X^4'
4: '-1+5*X-5*X^4-10*X^2+10*X^3+X^5'
3: '1-6*X-6*X^5+15*X^2-20*X^3+15*X^4+X^6'
2: '-1+7*X-7*X^6-21*X^2+35*X^3-35*X^4+21*X^5+X^7'
1: { 2 3 5 7 11 13 17 19 23 29 31 }
</pre>
Line 3,762 ⟶ 4,958:
Using the `polynomial` Rubygem, this can be written directly from the definition in the description:
<
def x_minus_1_to_the(p)
Line 3,779 ⟶ 4,975:
end
puts "\nPrimes below 50:", 50.times.select {|n| prime? n}.join(',')</
Or without the dependency:
<
p.times.inject([1]) do |ex, _|
([0] + ex).zip(ex + [0]).map { |x,y| x - y }
Line 3,803 ⟶ 4,999:
end
puts "\nPrimes below 50:", 50.times.select {|n| prime? n}.join(',')</
{{out}}
Line 3,820 ⟶ 5,016:
=={{header|Rust}}==
<syntaxhighlight lang="rust">fn aks_coefficients(k: usize) -> Vec<i64> {
let mut coefficients = vec![0i64; k + 1];
coefficients[0] = 1;
for i in 1..(k + 1) {
Line 3,840 ⟶ 5,034:
} else {
let c = aks_coefficients(p);
(1
}
}
Line 3,848 ⟶ 5,042:
println!("{}: {:?}", i, aks_coefficients(i));
}
for i in (1..
print!("{} ", i);
}
}</
{{out}}
<pre>0: [1]
Line 3,865 ⟶ 5,059:
An alternative version which computes the coefficients in a more functional but less efficient way.
<
fn aks_coefficients(k: usize) -> Vec<i64> {
if k == 0 {
Line 3,878 ⟶ 5,072:
}
}
</syntaxhighlight>
=={{header|Scala}}==
<
val pascal = (( Vector(Vector(BigInt(1))) /: (1 to 50)) { (rows, i) =>
Line 3,910 ⟶ 5,104:
val primes = (2 to 50).filter(isPrime)
println
println(primes mkString " ")</
{{out}}
Line 3,926 ⟶ 5,120:
=={{header|Scheme}}==
<syntaxhighlight lang="scheme">
;; implement mod m arithmetic with polnomials in x
;; as lists of coefficients, x^0 first.
Line 3,983 ⟶ 5,177:
;; > (filter rosetta-aks-test (iota 50))
;; '(2 3 5 7 11 13 17 19 23 29 31 37 41 43 47)
</syntaxhighlight>
=={{header|Scilab}}==
<syntaxhighlight lang="text">
clear
xdel(winsid())
Line 4,055 ⟶ 5,249:
end
disp(R)
</syntaxhighlight>
{{out}}
Line 4,079 ⟶ 5,273:
=={{header|Seed7}}==
<
const func array integer: expand_x_1 (in integer: p) is func
Line 4,102 ⟶ 5,296:
ex := expand_x_1(p);
ex[1] +:= 1;
for key idx range
noop;
end for;
Line 4,138 ⟶ 5,332:
end for;
writeln;
end func;</
{{out}}
Line 4,162 ⟶ 5,356:
=={{header|Sidef}}==
{{trans|Perl}}
<
p >= 2 || return false
for i in (1 .. p>>1) {
Line 4,184 ⟶ 5,378:
say "expansions of (x-1)^p:"
for i in ^10 { say binpoly(i) }
say "Primes to 80: [#{2..80 -> grep { binprime(_) }.join(' ')}]"</
{{out}}
<pre>
Line 4,204 ⟶ 5,398:
Using the [https://ideas.repec.org/c/boc/bocode/s455001.html moremata] library to print the polynomial coefficients. They are in decreasing degree order. To install moremata, type '''ssc install moremata''' in Stata. Since Stata is using double precision floating-point instead of 32 bit integers, the polynomials are exact up to p=54.
<
function pol(n) {
a=J(1,n+1,1)
Line 4,266 ⟶ 5,460:
2 3 5 7 11 13 17 19 23 29 31 37 41 43 47
end</
=={{header|Swift}}==
<
var result = [Int](count : n+1, repeatedValue : 0)
Line 4,349 ⟶ 5,543:
}
}
</syntaxhighlight>
{{out}}
Line 4,370 ⟶ 5,564:
=={{header|Tcl}}==
A recursive method with memorization would be more efficient, but this is sufficient for small-scale work.
<
set clist 1
for {set i 0} {$i < $p} {incr i} {
Line 4,420 ⟶ 5,614:
}
}
puts "primes under 50: [join $sub50primes ,]"</
{{out}}
<pre>
Line 4,433 ⟶ 5,627:
primes under 35: 2,3,5,7,11,13,17,19,23,29,31
primes under 50: 2,3,5,7,11,13,17,19,23,29,31,37,41,43,47
</pre>
=={{header|Transd}}==
{{trans|Python}}
<syntaxhighlight lang="scheme">
#lang transd
MainModule: {
poly: (λ n Long()
(with v Vector<Long>([1])
(for i in Range(n) do
(append v (/ (* (get v -1) (- (- n i))) (to-Long (+ i 1))))
)
(reverse v)
(ret v)
)
),
aks_test: (λ n Long() -> Bool()
(if (< n 2) (ret false))
(with v (poly n)
(set-el v 0 (+ (get v 0) 1))
(ret (not (any Range(in: v 0 -1) (λ (mod @it n)))))
)
),
_start: (λ (with v Vector<Long>()
(for p in Seq( 12 ) do
(set v (poly p))
(textout "(x - 1)^" p " = ")
(for i in v do (textout :sign i "x^" @idx " "))
(textout "\n")
)
(textout "\nList of primes in 2-62 interval:\n")
(for p in Seq( 2 62 ) do
(if (aks_test p) (textout p " "))
)
))
}</syntaxhighlight>{{out}}
<pre>
(x - 1)^0 = +1x^0
(x - 1)^1 = -1x^0 +1x^1
(x - 1)^2 = +1x^0 -2x^1 +1x^2
(x - 1)^3 = -1x^0 +3x^1 -3x^2 +1x^3
(x - 1)^4 = +1x^0 -4x^1 +6x^2 -4x^3 +1x^4
(x - 1)^5 = -1x^0 +5x^1 -10x^2 +10x^3 -5x^4 +1x^5
(x - 1)^6 = +1x^0 -6x^1 +15x^2 -20x^3 +15x^4 -6x^5 +1x^6
(x - 1)^7 = -1x^0 +7x^1 -21x^2 +35x^3 -35x^4 +21x^5 -7x^6 +1x^7
(x - 1)^8 = +1x^0 -8x^1 +28x^2 -56x^3 +70x^4 -56x^5 +28x^6 -8x^7 +1x^8
(x - 1)^9 = -1x^0 +9x^1 -36x^2 +84x^3 -126x^4 +126x^5 -84x^6 +36x^7 -9x^8 +1x^9
(x - 1)^10 = +1x^0 -10x^1 +45x^2 -120x^3 +210x^4 -252x^5 +210x^6 -120x^7 +45x^8 -10x^9 +1x^10
(x - 1)^11 = -1x^0 +11x^1 -55x^2 +165x^3 -330x^4 +462x^5 -462x^6 +330x^7 -165x^8 +55x^9 -11x^10 +1x^11
List of primes in 2-62 interval:
2 3 5 7 11 13 17 19 23 29 31 37 41 43 47 53 59 61
</pre>
=={{header|uBasic/4tH}}==
<syntaxhighlight lang="text">For n = 0 To 9
Push n : Gosub _coef : Gosub _drop
Print "(x-1)^";n;" = ";
Line 4,507 ⟶ 5,755:
_drop ' ( n --)
If Pop() Endif
Return</
{{out}}
<pre>
Line 4,526 ⟶ 5,774:
=={{header|VBA}}==
{{trans|Phix}}
<syntaxhighlight lang="vb">
'-- Does not work for primes above 97, which is actually beyond the original task anyway.
'-- Translated from the C version, just about everything is (working) out-by-1, what fun.
Line 4,602 ⟶ 5,850:
Next n
End Sub
</
(x-1)^0 = 1
(x-1)^1 = x-1
Line 4,615 ⟶ 5,863:
primes (<= 99 ):
2 3 5 7 11 13 17 19 23 29 31 37 41 43 47 53 59 61 67 71 73 79 83 89 97
</pre>
=={{header|V (Vlang)}}==
{{trans|go}}
<syntaxhighlight lang="v (vlang)">fn bc(p int) []i64 {
mut c := []i64{len: p+1}
mut r := i64(1)
for i, half := 0, p/2; i <= half; i++ {
c[i] = r
c[p-i] = r
r = r * i64(p-i) / i64(i+1)
}
for i := p - 1; i >= 0; i -= 2 {
c[i] = -c[i]
}
return c
}
fn main() {
for p := 0; p <= 7; p++ {
println("$p: ${pp(bc(p))}")
}
for p := 2; p < 50; p++ {
if aks(p) {
print("$p ")
}
}
println('')
}
const e = [`²`,`³`,`⁴`,`⁵`,`⁶`,`⁷`]
fn pp(c []i64) string {
mut s := ''
if c.len == 1 {
return c[0].str()
}
p := c.len - 1
if c[p] != 1 {
s = c[p].str()
}
for i := p; i > 0; i-- {
s += "x"
if i != 1 {
s += e[i-2].str()
}
d := c[i-1]
if d < 0 {
s += " - ${-d}"
} else {
s += " + $d"
}
}
return s
}
fn aks(p int) bool {
mut c := bc(p)
c[p]--
c[0]++
for d in c {
if d%i64(p) != 0 {
return false
}
}
return true
}</syntaxhighlight>
{{out}}
<pre>
0: 1
1: x - 1
2: x² - 2x + 1
3: x³ - 3x² + 3x - 1
4: x⁴ - 4x³ + 6x² - 4x + 1
5: x⁵ - 5x⁴ + 10x³ - 10x² + 5x - 1
6: x⁶ - 6x⁵ + 15x⁴ - 20x³ + 15x² - 6x + 1
7: x⁷ - 7x⁶ + 21x⁵ - 35x⁴ + 35x³ - 21x² + 7x - 1
2 3 5 7 11 13 17 19 23 29 31 37 41 43 47
</pre>
=={{header|Wren}}==
{{trans|Go}}
<syntaxhighlight lang="wren">var bc = Fn.new { |p|
var c = List.filled(p+1, 0)
var r = 1
var half = (p/2).floor
for (i in 0..half) {
c[i] = r
c[p-i] = r
r = (r * (p-i) / (i+1)).floor
}
var j = p - 1
while (j >= 0) {
c[j] = -c[j]
j = j - 2
}
return c
}
var e = "²³⁴⁵⁶⁷".codePoints.toList
var pp = Fn.new { |c|
if (c.count == 1) return "%(c[0])"
var p = c.count - 1
var s = ""
if (c[p] != 1) s = "%(c[p])"
if (p == 0) return s
for (i in p..1) {
s = s + "x"
if (i != 1) s = s + String.fromCodePoint(e[i-2])
var d = c[i-1]
s = s + ((d < 0) ? " - %(-d)" : " + %(d)")
}
return s
}
var aks = Fn.new { |p|
var c = bc.call(p)
c[p] = c[p] - 1
c[0] = c[0] + 1
for (d in c) {
if (d%p != 0) return false
}
return true
}
for (p in 0..7) System.print("%(p): %(pp.call(bc.call(p)))")
System.print("\nAll primes under 50:")
for (p in 2..49) if (aks.call(p)) System.write("%(p) ")
System.print()</syntaxhighlight>
{{out}}
<pre>
0: 1
1: x - 1
2: x² - 2x + 1
3: x³ - 3x² + 3x - 1
4: x⁴ - 4x³ + 6x² - 4x + 1
5: x⁵ - 5x⁴ + 10x³ - 10x² + 5x - 1
6: x⁶ - 6x⁵ + 15x⁴ - 20x³ + 15x² - 6x + 1
7: x⁷ - 7x⁶ + 21x⁵ - 35x⁴ + 35x³ - 21x² + 7x - 1
All primes under 50:
2 3 5 7 11 13 17 19 23 29 31 37 41 43 47
</pre>
=={{header|Yabasic}}==
{{trans|Phix}}
<
// Translated from the C version, just about everything is (working) out-by-1, what fun.
Line 4,703 ⟶ 6,095:
end sub
AKS_test_for_primes()</
=={{header|Zig}}==
Uses Zig's compile-time interpreter to create Pascal's triangle at compile-time.
<syntaxhighlight lang="zig">
const std = @import("std");
const assert = std.debug.assert;
pub fn main() !void {
const stdout = std.io.getStdOut().writer();
var i: u6 = 0;
while (i < 8) : (i += 1)
try showBinomial(stdout, i);
try stdout.print("\nThe primes upto 50 (via AKS) are: ", .{});
i = 2;
while (i <= 50) : (i += 1) if (aksPrime(i))
try stdout.print("{} ", .{i});
try stdout.print("\n", .{});
}
fn showBinomial(writer: anytype, n: u6) !void {
const row = binomial(n).?;
var sign: u8 = '+';
var exp = row.len;
try writer.print("(x - 1)^{} =", .{n});
for (row) |coef| {
try writer.print(" ", .{});
if (exp != row.len)
try writer.print("{c} ", .{sign});
exp -= 1;
if (coef != 1 or exp == 0)
try writer.print("{}", .{coef});
if (exp >= 1) {
try writer.print("x", .{});
if (exp > 1)
try writer.print("^{}", .{exp});
}
sign = if (sign == '+') '-' else '+';
}
try writer.print("\n", .{});
}
fn aksPrime(n: u6) bool {
return for (binomial(n).?) |coef| {
if (coef > 1 and coef % n != 0)
break false;
} else true;
}
pub fn binomial(n: u32) ?[]const u64 {
if (n >= rmax)
return null
else {
const k = n * (n + 1) / 2;
return pascal[k .. k + n + 1];
}
}
const rmax = 68;
// evaluated and created at compile-time
const pascal = build: {
@setEvalBranchQuota(100_000);
var coefficients: [(rmax * (rmax + 1)) / 2]u64 = undefined;
coefficients[0] = 1;
var j: u32 = 0;
var k: u32 = 1;
var n: u32 = 1;
while (n < rmax) : (n += 1) {
var prev = coefficients[j .. j + n];
var next = coefficients[k .. k + n + 1];
next[0] = 1;
var i: u32 = 1;
while (i < n) : (i += 1)
next[i] = prev[i] + prev[i - 1];
next[i] = 1;
j = k;
k += n + 1;
}
break :build coefficients;
};
</syntaxhighlight>
{{Out}}
<pre>
$ zig run aks.zig
(x - 1)^0 = 1
(x - 1)^1 = x - 1
(x - 1)^2 = x^2 - 2x + 1
(x - 1)^3 = x^3 - 3x^2 + 3x - 1
(x - 1)^4 = x^4 - 4x^3 + 6x^2 - 4x + 1
(x - 1)^5 = x^5 - 5x^4 + 10x^3 - 10x^2 + 5x - 1
(x - 1)^6 = x^6 - 6x^5 + 15x^4 - 20x^3 + 15x^2 - 6x + 1
(x - 1)^7 = x^7 - 7x^6 + 21x^5 - 35x^4 + 35x^3 - 21x^2 + 7x - 1
The primes upto 50 (via AKS) are: 2 3 5 7 11 13 17 19 23 29 31 37 41 43 47
</pre>
=={{header|zkl}}==
{{trans|Python}}
<
fcn expand_x_1(p){
ex := L(BN(1));
Line 4,726 ⟶ 6,215:
println("\n# small primes using the aks test");
println([0..110].filter(aks_test).toString(*));</
{{out}}
<pre>
|