Pythagorean triples: Difference between revisions
m
→{{header|Wren}}: Changed to Wren S/H
(Add Nimrod) |
m (→{{header|Wren}}: Changed to Wren S/H) |
||
| (152 intermediate revisions by 63 users not shown) | |||
Line 1:
{{task}}
A [[wp:Pythagorean_triple|Pythagorean triple]] is defined as three positive integers <math>(a, b, c)</math> where <math>a < b < c</math>, and <math>a^2+b^2=c^2.</math>
They are called primitive triples if <math>a, b, c</math> are co-prime, that is, if their pairwise greatest common divisors <math>{\rm gcd}(a, b) = {\rm gcd}(a, c) = {\rm gcd}(b, c) = 1</math>.
Because of their relationship through the Pythagorean theorem, a, b, and c are co-prime if a and b are co-prime (<math>{\rm gcd}(a, b) = 1</math>).
Each triple forms the length of the sides of a right triangle, whose perimeter is <math>P=a+b+c</math>.
;Task:
The task is to determine how many Pythagorean triples there are with a perimeter no larger than 100 and the number of these that are primitive.
;Extra credit:
Deal with large values. Can your program handle a maximum perimeter of 1,000,000? What about 10,000,000? 100,000,000?
Note: the extra credit is not for you to demonstrate how fast your language is compared to others; you need a proper algorithm to solve them in a timely manner.
;Related tasks:
* [[Euler's sum of powers conjecture]]
* [[List comprehensions]]
* [[Pythagorean quadruples]]
<br><br>
=={{header|11l}}==
{{trans|D}}
<syntaxhighlight lang="11l">Int64 nTriples, nPrimitives, limit
F countTriples(Int64 =x, =y, =z)
L
V p = x + y + z
I p > :limit
R
:nPrimitives++
:nTriples += :limit I/ p
V t0 = x - 2 * y + 2 * z
V t1 = 2 * x - y + 2 * z
V t2 = t1 - y + z
countTriples(t0, t1, t2)
t0 += 4 * y
t1 += 2 * y
t2 += 4 * y
countTriples(t0, t1, t2)
z = t2 - 4 * x
y = t1 - 4 * x
x = t0 - 2 * x
L(p) 1..8
limit = Int64(10) ^ p
nTriples = nPrimitives = 0
countTriples(3, 4, 5)
print(‘Up to #11: #11 triples, #9 primitives.’.format(limit, nTriples, nPrimitives))</syntaxhighlight>
{{out}}
<pre>
Up to 10: 0 triples, 0 primitives.
Up to 100: 17 triples, 7 primitives.
Up to 1000: 325 triples, 70 primitives.
Up to 10000: 4858 triples, 703 primitives.
Up to 100000: 64741 triples, 7026 primitives.
Up to 1000000: 808950 triples, 70229 primitives.
Up to 10000000: 9706567 triples, 702309 primitives.
Up to 100000000: 113236940 triples, 7023027 primitives.
</pre>
=={{header|360 Assembly}}==
{{trans|Perl}}
<syntaxhighlight lang="360asm">* Pythagorean triples - 12/06/2018
PYTHTRI CSECT
USING PYTHTRI,R13 base register
B 72(R15) skip savearea
DC 17F'0' savearea
SAVE (14,12) save previous context
ST R13,4(R15) link backward
ST R15,8(R13) link forward
LR R13,R15 set addressability
MVC PMAX,=F'1' pmax=1
LA R6,1 i=1
DO WHILE=(C,R6,LE,=F'6') do i=1 to 6
L R5,PMAX pmax
MH R5,=H'10' *10
ST R5,PMAX pmax=pmax*10
MVC PRIM,=F'0' prim=0
MVC COUNT,=F'0' count=0
L R1,PMAX pmax
BAL R14,ISQRT isqrt(pmax)
SRA R0,1 /2
ST R0,NMAX nmax=isqrt(pmax)/2
LA R7,1 n=1
DO WHILE=(C,R7,LE,NMAX) do n=1 to nmax
LA R9,1(R7) m=n+1
LR R5,R9 m
AR R5,R7 +n
MR R4,R9 *m
SLA R5,1 *2
LR R8,R5 p=2*m*(m+n)
DO WHILE=(C,R8,LE,PMAX) do while p<=pmax
LR R1,R9 m
LR R2,R7 n
BAL R14,GCD gcd(m,n)
IF C,R0,EQ,=F'1' THEN if gcd(m,n)=1 then
L R2,PRIM prim
LA R2,1(R2) +1
ST R2,PRIM prim=prim+1
L R4,PMAX pmax
SRDA R4,32 ~
DR R4,R8 /p
A R5,COUNT +count
ST R5,COUNT count=count+pmax/p
ENDIF , endif
LA R9,2(R9) m=m+2
LR R5,R9 m
AR R5,R7 +n
MR R4,R9 *m
SLA R5,1 *2
LR R8,R5 p=2*m*(m+n)
ENDDO , enddo n
LA R7,1(R7) n++
ENDDO , enddo n
L R1,PMAX pmax
XDECO R1,XDEC edit pmax
MVC PG+15(9),XDEC+3 output pmax
L R1,COUNT count
XDECO R1,XDEC edit count
MVC PG+33(9),XDEC+3 output count
L R1,PRIM prim
XDECO R1,XDEC edit prim
MVC PG+55(9),XDEC+3 output prim
XPRNT PG,L'PG print
LA R6,1(R6) i++
ENDDO , enddo i
L R13,4(0,R13) restore previous savearea pointer
RETURN (14,12),RC=0 restore registers from calling sav
NMAX DS F nmax
PMAX DS F pmax
COUNT DS F count
PRIM DS F prim
PG DC CL80'Max Perimeter: ........., Total: ........., Primitive:'
XDEC DS CL12
GCD EQU * --------------- function gcd(a,b)
STM R2,R7,GCDSA save context
LR R3,R1 c=a
LR R4,R2 d=b
GCDLOOP LR R6,R3 c
SRDA R6,32 ~
DR R6,R4 /d
LTR R6,R6 if c mod d=0
BZ GCDELOOP then leave loop
LR R5,R6 e=c mod d
LR R3,R4 c=d
LR R4,R5 d=e
B GCDLOOP loop
GCDELOOP LR R0,R4 return(d)
LM R2,R7,GCDSA restore context
BR R14 return
GCDSA DS 6A context store
ISQRT EQU * --------------- function isqrt(n)
STM R3,R10,ISQRTSA save context
LR R6,R1 n=r1
LR R10,R6 sqrtn=n
SRA R10,1 sqrtn=n/2
IF LTR,R10,Z,R10 THEN if sqrtn=0 then
LA R10,1 sqrtn=1
ELSE , else
LA R9,0 snm2=0
LA R8,0 snm1=0
LA R7,0 sn=0
LA R3,0 okexit=0
DO UNTIL=(C,R3,EQ,=A(1)) do until okexit=1
AR R10,R7 sqrtn=sqrtn+sn
LR R9,R8 snm2=snm1
LR R8,R7 snm1=sn
LR R4,R6 n
SRDA R4,32 ~
DR R4,R10 /sqrtn
SR R5,R10 -sqrtn
SRA R5,1 /2
LR R7,R5 sn=(n/sqrtn-sqrtn)/2
IF C,R7,EQ,=F'0',OR,CR,R7,EQ,R9 THEN if sn=0 or sn=snm2 then
LA R3,1 okexit=1
ENDIF , endif
ENDDO , enddo until
ENDIF , endif
LR R5,R10 sqrtn
MR R4,R10 *sqrtn
IF CR,R5,GT,R6 THEN if sqrtn*sqrtn>n then
BCTR R10,0 sqrtn=sqrtn-1
ENDIF , endif
LR R0,R10 return(sqrtn)
LM R3,R10,ISQRTSA restore context
BR R14 return
ISQRTSA DS 8A context store
YREGS
END PYTHTRI</syntaxhighlight>
{{out}}
<pre>
Max Perimeter: 10, Total: 0, Primitive: 0
Max Perimeter: 100, Total: 17, Primitive: 7
Max Perimeter: 1000, Total: 325, Primitive: 70
Max Perimeter: 10000, Total: 4858, Primitive: 703
Max Perimeter: 100000, Total: 64741, Primitive: 7026
Max Perimeter: 1000000, Total: 808950, Primitive: 70229
</pre>
=={{header|Action!}}==
<syntaxhighlight lang="action!">DEFINE PTR="CARD"
DEFINE ENTRY_SIZE="3"
TYPE TRIPLE=[BYTE a,b,c]
TYPE TRIPLES=[
PTR buf ;BYTE ARRAY
BYTE count]
PTR FUNC GetItemAddr(TRIPLES POINTER arr BYTE index)
PTR addr
addr=arr.buf+index*ENTRY_SIZE
RETURN (addr)
PROC PrintTriples(TRIPLES POINTER arr)
INT i
TRIPLE POINTER t
FOR i=0 TO arr.count-1
DO
t=GetItemAddr(arr,i)
PrintF("(%B %B %B) ",t.a,t.b,t.c)
OD
RETURN
PROC Init(TRIPLES POINTER arr BYTE ARRAY b)
arr.buf=b
arr.count=0
RETURN
PROC AddItem(TRIPLES POINTER arr TRIPLE POINTER t)
TRIPLE POINTER p
p=GetItemAddr(arr,arr.count)
p.a=t.a
p.b=t.b
p.c=t.c
arr.count==+1
RETURN
PROC FindTriples(TRIPLES POINTER res BYTE limit)
BYTE ARRAY data(100)
BYTE half,i,j,k
TRIPLE t
Init(res,data)
half=limit/2
FOR i=1 TO half
DO
FOR j=i TO half
DO
FOR k=j TO limit
DO
IF i+j+k<limit AND i*i+j*j=k*k THEN
t.a=i t.b=j t.c=k
AddItem(res,t)
FI
OD
OD
OD
RETURN
BYTE FUNC Gcd(BYTE a,b)
BYTE tmp
IF a<b THEN
tmp=a a=b b=tmp
FI
WHILE b#0
DO
tmp=a MOD b
a=b b=tmp
OD
RETURN (a)
BYTE FUNC IsPrimitive(TRIPLE POINTER t)
IF Gcd(t.a,t.b)>1 THEN RETURN (0) FI
IF Gcd(t.b,t.c)>1 THEN RETURN (0) FI
IF Gcd(t.a,t.c)>1 THEN RETURN (0) FI
RETURN (1)
PROC FindPrimitives(TRIPLES POINTER arr,res)
BYTE ARRAY data(100)
INT i
TRIPLE POINTER t
Init(res,data)
FOR i=0 TO arr.count-1
DO
t=GetItemAddr(arr,i)
IF IsPrimitive(t) THEN
AddItem(res,t)
FI
OD
RETURN
PROC Main()
DEFINE LIMIT="100"
TRIPLES res,res2
FindTriples(res,LIMIT)
PrintF("There are %B pythagorean triples with a perimeter less than %B:%E%E",res.count,LIMIT)
PrintTriples(res)
FindPrimitives(res,res2)
PrintF("%E%E%E%B of them are primitive:%E%E",res2.count)
PrintTriples(res2)
RETURN</syntaxhighlight>
{{out}}
[https://gitlab.com/amarok8bit/action-rosetta-code/-/raw/master/images/Pythagorean_triples.png Screenshot from Atari 8-bit computer]
<pre>
There are 17 pythagorean triples with a perimeter less than 100:
(3 4 5) (5 12 13) (6 8 10) (7 24 25) (8 15 17) (9 12 15) (9 40 41) (10 24 26) (12 16 20)
(12 35 37) (15 20 25) (15 36 39) (16 30 34) (18 24 30) (20 21 29) (21 28 35) (24 32 40)
7 of them are primitive:
(3 4 5) (5 12 13) (7 24 25) (8 15 17) (9 40 41) (12 35 37) (20 21 29)
</pre>
=={{header|Ada}}==
Line 16 ⟶ 336:
Translation of efficient method from C, see [[wp:Pythagorean_triple#Parent.2Fchild_relationships|the WP article]]. Compiles on gnat/gcc.
<
procedure Pythagorean_Triples is
Line 48 ⟶ 368:
Large_Natural'Image(P_Cnt) & " Primitives");
end loop;
end Pythagorean_Triples;</
Output:
Line 61 ⟶ 381:
Up to 10 ** 8 : 113236940 Triples, 7023027 Primitives
Up to 10 ** 9 : 1294080089 Triples, 70230484 Primitives</pre>
=={{header|ALGOL 68}}==
Uses a table of square roots so OK for perimeters up to 1000 (or possibly 10 000).
<syntaxhighlight lang="algol68">
BEGIN # find some Pythagorean triples ( a, b, c ) #
# where a < b < c and a^2 + b^2 = c^2 #
INT max perimeter = 100; # maximum a + b + c we will consider #
INT max square = max perimeter * max perimeter;
# form a table of square roots of numbers to max perimeter ^ 2 #
[ 1 : max square ]INT sr;
FOR i TO UPB sr DO sr[ i ] := 0 OD;
FOR i TO max perimeter DO sr[ i * i ] := i OD;
PROC gcd = ( INT x, y )INT: # iterative gcd #
BEGIN
INT a := ABS x, b := ABS y;
WHILE b /= 0 DO
INT next a = b;
b := a MOD b;
a := next a
OD;
a
END # gcd # ;
# count the Pythagorean triples #
INT t count := 0, p count := 0;
FOR a TO max perimeter DO
INT a2 = a * a;
FOR b FROM a + 1 TO max perimeter - a
WHILE INT c = sr[ a2 + ( b * b ) ];
a + b + c <= max perimeter
DO IF c > b THEN # have a triple #
t count +:= 1;
IF gcd( a, b ) = 1 THEN # have a primitive triple #
p count +:= 1
FI
FI
OD
OD;
print( ( "Pythagorean triples with perimeters up to ", whole( max perimeter, 0 ), ":", newline ) );
print( ( " Primitive: ", whole( p count, 0 ), newline ) );
print( ( " Total: ", whole( t count, 0 ), newline ) )
END
</syntaxhighlight>
{{out}}
<pre>
Pythagorean triples with perimeters up to 100:
Primitive: 7
Total: 17
</pre>
=={{header|Arturo}}==
<syntaxhighlight lang="rebol">triples: new []
loop 1..50 'x [
loop 1..50 'y [
loop (max @[x y])..100 'z [
if 100 > sum @[x y z] [
if (z^2) = add x^2 y^2 ->
'triples ++ @[sort @[x y z]]
]
]
]
]
unique 'triples
print ["Found" size triples "pythagorean triples with a perimeter no larger than 100:"]
print triples
primitive: select triples => [1 = gcd]
print ""
print [size primitive "of them are primitive:"]
print primitive</syntaxhighlight>
{{out}}
<pre>Found 17 pythagorean triples with a perimeter no larger than 100:
[3 4 5] [5 12 13] [6 8 10] [7 24 25] [8 15 17] [9 12 15] [9 40 41] [10 24 26] [12 16 20] [12 35 37] [15 20 25] [15 36 39] [16 30 34] [18 24 30] [20 21 29] [21 28 35] [24 32 40]
7 of them are primitive:
[3 4 5] [5 12 13] [7 24 25] [8 15 17] [9 40 41] [12 35 37] [20 21 29]</pre>
=={{header|APL}}==
<syntaxhighlight lang="apl">
⍝ Determine whether given list of integers has GCD = 1
primitive←∧/1=2∨/⊢
⍝ Filter list given as right operand by applying predicate given as left operand
filter←{⍵⌿⍨⍺⍺ ⍵}
⍝ Function pytriples finds all triples given a maximum perimeter
∇res←pytriples maxperimeter;sos;sqrt;cartprod;ascending;ab_max;c_max;a_b_pairs;sos_is_sq;add_c;perimeter_rule
⍝ Input parameter maxperimeter is the maximum perimeter
⍝ Sum of squares of given list of nrs
sos←+/(×⍨⊢)
⍝ Square root
sqrt←(÷2)*⍨⊢
⍝ (cartesian product) all possible pairs of integers
⍝ from 1 to ⍵
cartprod←{,{⍺∘.,⍵}⍨⍳⍵}
⍝ Predicate: are values in given list ascending
⍝ Given e.g. pair a, b, c: is a ≤ b ≤ c?
ascending←∧/2≤/⊢
ab_max←⌊maxperimeter÷2
c_max←⌈maxperimeter×sqrt 2
⍝ Selects from all a,b combinations (a<abmax, b<abmax)
⍝ only those pairs where a ≤ b.
a_b_pairs←ascending filter¨cartprod(ab_max)
⍝ Predicate: is the sum of squares of a and b
⍝ itself a square? (does it occur in the squares list)
sos_is_sq←{{⍵≠1+c_max}(×⍨⍳c_max)⍳sos¨⍵}
⍝ Given a pair a,b add corresponding c to form a triple
add_c←{⍵,sqrt sos ⍵}
⍝ Predicate: sum of items less than or equal to max
perimeter_rule←{maxperimeter≥+/⍵}
res←perimeter_rule¨filter add_c¨sos_is_sq filter a_b_pairs
∇</syntaxhighlight>
{{out}}
<pre>
⍝ Get the number of triples
≢pytriples 100
17
⍝ Get the number of primitive triples
≢primitive¨filter pytriples 100
7
</pre>
=={{header|AutoHotkey}}==
<
SetBatchLines, -1
#SingleInstance, Force
Line 92 ⟶ 540:
Loop, 8
Msgbox % 10**A_Index ": " count_triples(10**A_Index)</
{{out}}
Line 104 ⟶ 552:
100000000: 7023027 primitives out of 113236940 triples</pre>
=={{header|
<syntaxhighlight lang="awk">
# syntax: GAWK -f PYTHAGOREAN_TRIPLES.AWK
# converted from Go
BEGIN {
printf("%5s %11s %11s %11s %s\n","limit","limit","triples","primitives","seconds")
for (max_peri=10; max_peri<=1E9; max_peri*=10) {
t = systime()
prim = 0
total = 0
new_tri(3,4,5)
printf("10^%-2d %11d %11d %11d %d\n",++n,max_peri,total,prim,systime()-t)
}
exit(0)
}
function new_tri(s0,s1,s2, p) {
p = s0 + s1 + s2
if (p <= max_peri) {
prim++
total += int(max_peri / p)
new_tri(+1*s0-2*s1+2*s2,+2*s0-1*s1+2*s2,+2*s0-2*s1+3*s2)
new_tri(+1*s0+2*s1+2*s2,+2*s0+1*s1+2*s2,+2*s0+2*s1+3*s2)
new_tri(-1*s0+2*s1+2*s2,-2*s0+1*s1+2*s2,-2*s0+2*s1+3*s2)
}
}
</syntaxhighlight>
{{out}}
<pre>
limit limit triples primitives seconds
10^1 10 0 0 0
10^2 100 17 7 0
10^3 1000 325 70 0
10^4 10000 4858 703 0
10^5 100000 64741 7026 0
10^6 1000000 808950 70229 0
10^7 10000000 9706567 702309 2
10^8 100000000 113236940 7023027 12
10^9 1000000000 1294080089 70230484 116
</pre>
=={{header|BASIC}}==
==={{header|ANSI BASIC}}===
{{trans|BBC BASIC}}
{{works with|Decimal BASIC}}
<syntaxhighlight lang="basic">100 DECLARE EXTERNAL SUB tri
110 !
120 PUBLIC NUMERIC U0(3,3), U1(3,3), U2(3,3), all, prim
130 DIM seed(3)
140 MAT READ U0, U1, U2
150 DATA 1, -2, 2, 2, -1, 2, 2, -2, 3
160 DATA 1, 2, 2, 2, 1, 2, 2, 2, 3
170 DATA -1, 2, 2, -2, 1, 2, -2, 2, 3
180 !
190 MAT READ seed
200 DATA 3, 4, 5
210 FOR power = 1 TO 7
220 LET all = 0
230 LET prim = 0
240 CALL tri(seed, 10^power , all , prim)
250 PRINT "Up to 10^";power,
260 PRINT USING "######### triples ######### primitives":all,prim
270 NEXT power
280 END
290 !
300 EXTERNAL SUB tri(i(), mp, all, prim)
310 DECLARE EXTERNAL FUNCTION SUM
320 DECLARE NUMERIC t(3)
330 !
340 IF SUM(i) > mp THEN EXIT SUB
350 LET prim = prim + 1
360 LET all = all + INT(mp / SUM(i))
370 !
380 MAT t = U0 * i
390 CALL tri(t, mp , all , prim)
400 MAT t = U1 * i
410 CALL tri(t, mp , all , prim)
420 MAT t = U2 * i
430 CALL tri(t, mp , all , prim)
440 END SUB
450 !
460 EXTERNAL FUNCTION SUM(a())
470 LET temp = 0
480 FOR i=LBOUND(a) TO UBOUND(a)
490 LET temp = temp + a(i)
500 NEXT i
510 LET SUM = temp
520 END FUNCTION</syntaxhighlight>
{{out}}
<pre>
Up to 10^ 1 0 triples 0 primitives
Up to 10^ 2 17 triples 7 primitives
Up to 10^ 3 325 triples 70 primitives
Up to 10^ 4 4858 triples 703 primitives
Up to 10^ 5 64741 triples 7026 primitives
Up to 10^ 6 808950 triples 70229 primitives
Up to 10^ 7 9706567 triples 702309 primitives
</pre>
==={{header|BBC BASIC}}===
The built-in array arithmetic is very well suited to this task!
<
U0%() = 1, -2, 2, 2, -1, 2, 2, -2, 3
U1%() = 1, 2, 2, 2, 1, 2, 2, 2, 3
Line 132 ⟶ 678:
t%() = U2%() . i%()
PROCtri(t%(), mp%, all%, prim%)
ENDPROC</
'''Output:'''
<pre>
Line 147 ⟶ 693:
=={{header|Bracmat}}==
{{trans|C}}
<
total prim max-peri U
. (.(1,-2,2) (2,-1,2) (2,-2,3))
Line 200 ⟶ 746:
pythagoreanTriples$;
</syntaxhighlight>
Output (under Linux):
Line 215 ⟶ 761:
With very few changes we can get rid of the stack exhausting recursion. Instead of calling <code>new-tri</code> recursively, be push the triples to test onto a stack and return to the <code>Main</code> function. In the innermost loop we pop a triple from the stack and call <code>new-tri</code>. The memory overhead is only a few megabytes for a max perimeter of 100,000,000. On my Windows XP box the whole computation takes at least 15 minutes! Given enough time (and memory), the program can compute results for larger perimeters.
<
total prim max-peri U stack
. (.(1,-2,2) (2,-1,2) (2,-2,3))
Line 268 ⟶ 814:
);
pythagoreanTriples$;</
=={{header|C}}==
Sample implemention; naive method, patentedly won't scale to larger numbers, despite the attempt to optimize it. Calculating up to 10000 is already a test of patience.
<
#include <stdlib.h>
Line 318 ⟶ 864:
return 0;
}</
Efficient method, generating primitive triples only as described in [[wp:Pythagorean_triple#Parent.2Fchild_relationships|the same WP article]]:<
#include <stdlib.h>
#include <stdint.h>
Line 366 ⟶ 912:
}
return 0;
}</
Up to 100: 17 triples, 7 primitives.
Up to 1000: 325 triples, 70 primitives.
Line 373 ⟶ 919:
Up to 1000000: 808950 triples, 70229 primitives.
Up to 10000000: 9706567 triples, 702309 primitives.
Up to 100000000: 113236940 triples, 7023027 primitives.</
Same as above, but with loop unwound and third recursion eliminated:
<
#include <stdlib.h>
#include <stdint.h>
Line 426 ⟶ 972:
}
return 0;
}</
=={{header|C++}}==
<syntaxhighlight lang="cpp">#include <cmath>
#include <iostream>
#include <numeric>
#include <tuple>
#include <vector>
using namespace std;
auto CountTriplets(unsigned long long maxPerimeter)
{
unsigned long long totalCount = 0;
unsigned long long primitveCount = 0;
auto max_M = (unsigned long long)sqrt(maxPerimeter/2) + 1;
for(unsigned long long m = 2; m < max_M; ++m)
{
for(unsigned long long n = 1 + m % 2; n < m; n+=2)
{
if(gcd(m,n) != 1)
{
continue;
}
// The formulas below will generate primitive triples if:
// 0 < n < m
// m and n are relatively prime (gcd == 1)
// m + n is odd
auto a = m * m - n * n;
auto b = 2 * m * n;
auto c = m * m + n * n;
auto perimeter = a + b + c;
if(perimeter <= maxPerimeter)
{
primitveCount++;
totalCount+= maxPerimeter / perimeter;
}
}
}
return tuple(totalCount, primitveCount);
}
int main()
{
vector<unsigned long long> inputs{100, 1000, 10'000, 100'000,
1000'000, 10'000'000, 100'000'000, 1000'000'000,
10'000'000'000}; // This last one takes almost a minute
for(auto maxPerimeter : inputs)
{
auto [total, primitive] = CountTriplets(maxPerimeter);
cout << "\nMax Perimeter: " << maxPerimeter << ", Total: " << total << ", Primitive: " << primitive ;
}
}
</syntaxhighlight>
{{out}}
<pre>
Max Perimeter: 100, Total: 17, Primitive: 7
Max Perimeter: 1000, Total: 325, Primitive: 70
Max Perimeter: 10000, Total: 4858, Primitive: 703
Max Perimeter: 100000, Total: 64741, Primitive: 7026
Max Perimeter: 1000000, Total: 808950, Primitive: 70229
Max Perimeter: 10000000, Total: 9706567, Primitive: 702309
Max Perimeter: 100000000, Total: 113236940, Primitive: 7023027
Max Perimeter: 1000000000, Total: 1294080089, Primitive: 70230484
Max Perimeter: 10000000000, Total: 14557915466, Primitive: 702304875
</pre>
=={{header|C sharp|C#}}==
Based on Ada example, which is a translation of efficient method from C, see [[wp:Pythagorean_triple#Parent.2Fchild_relationships|the WP article]].
<syntaxhighlight lang="c sharp">using System;
namespace RosettaCode.CSharp
{
class Program
{
static void Count_New_Triangle(ulong A, ulong B, ulong C, ulong Max_Perimeter, ref ulong Total_Cnt, ref ulong Primitive_Cnt)
{
ulong Perimeter = A + B + C;
if (Perimeter <= Max_Perimeter)
{
Primitive_Cnt = Primitive_Cnt + 1;
Total_Cnt = Total_Cnt + Max_Perimeter / Perimeter;
Count_New_Triangle(A + 2 * C - 2 * B, 2 * A + 2 * C - B, 2 * A + 3 * C - 2 * B, Max_Perimeter, ref Total_Cnt, ref Primitive_Cnt);
Count_New_Triangle(A + 2 * B + 2 * C, 2 * A + B + 2 * C, 2 * A + 2 * B + 3 * C, Max_Perimeter, ref Total_Cnt, ref Primitive_Cnt);
Count_New_Triangle(2 * B + 2 * C - A, B + 2 * C - 2 * A, 2 * B + 3 * C - 2 * A, Max_Perimeter, ref Total_Cnt, ref Primitive_Cnt);
}
}
static void Count_Pythagorean_Triples()
{
ulong T_Cnt, P_Cnt;
for (int I = 1; I <= 8; I++)
{
T_Cnt = 0;
P_Cnt = 0;
ulong ExponentNumberValue = (ulong)Math.Pow(10, I);
Count_New_Triangle(3, 4, 5, ExponentNumberValue, ref T_Cnt, ref P_Cnt);
Console.WriteLine("Perimeter up to 10E" + I + " : " + T_Cnt + " Triples, " + P_Cnt + " Primitives");
}
}
static void Main(string[] args)
{
Count_Pythagorean_Triples();
}
}
}</syntaxhighlight>
Output:
<pre>Perimeter up to 10E1 : 0 Triples, 0 Primitives
Perimeter up to 10E2 : 17 Triples, 7 Primitives
Perimeter up to 10E3 : 325 Triples, 70 Primitives
Perimeter up to 10E4 : 4858 Triples, 703 Primitives
Perimeter up to 10E5 : 64741 Triples, 7026 Primitives
Perimeter up to 10E6 : 808950 Triples, 70229 Primitives
Perimeter up to 10E7 : 9706567 Triples, 702309 Primitives
Perimeter up to 10E8 : 113236940 Triples, 7023027 Primitives</pre>
=={{header|Clojure}}==
Line 432 ⟶ 1,104:
for each pair ''(m,n)'' such that ''m>n>0'', ''m'' and ''n'' coprime and of opposite polarity (even/odd),
there is a primitive Pythagorean triple. It can be proven that the converse is true as well.
<
(defn pyth [peri]
Line 448 ⟶ 1,120:
(reduce (fn [[total prims] t] [(inc total), (if (first t) (inc prims) prims)])
[0 0]
ts))</
To handle really large perimeters, we can dispense with actually generating the triples and just calculate the counts:
<
(reduce (fn [[total prims] k] [(+ total k), (inc prims)]) [0 0]
(for [m (range 2 (Math/sqrt (/ peri 2)))
Line 457 ⟶ 1,129:
:while (<= p peri)
:when (= 1 (gcd m n))]
(quot peri p))))</
=={{header|CoffeeScript}}==
This algorithm scales linearly with the max perimeter. It uses two loops that are capped by the square root of the half-perimeter to examine/count provisional values of m and n, where m and n generate a, b, c, and p using simple number theory.
<
gcd = (x, y) ->
return x if y == 0
Line 496 ⟶ 1,168:
max_perim = Math.pow 10, 9 # takes under a minute
count_triples(max_perim)
</syntaxhighlight>
output
<pre>
Line 505 ⟶ 1,177:
=={{header|Common Lisp}}==
<
(loop for x in a collect
(loop for y in x
for z in b sum (* y z))))
(defun count-tri (lim &aux (prim 0) (cnt 0))
(
(count1 (mmul '(( 1 -2 2) ( 2 -1 2) ( 2 -2 3)) tr))
(
(count1 (mmul '((
(loop for p from 2 do (count-tri (expt 10 p)))</
1000: 70 prim, 325 all
10000: 703 prim, 4858 all
Line 529 ⟶ 1,199:
1000000: 70229 prim, 808950 all
10000000: 702309 prim, 9706567 all
...</
=={{header|Crystal}}==
{{trans|Ruby}}
<syntaxhighlight lang="ruby">class PythagoranTriplesCounter
def initialize(limit = 0)
@limit = limit
@total = 0
@primitives = 0
generate_triples(3, 4, 5)
end
def total; @total end
def primitives; @primitives end
private def generate_triples(a, b, c)
perim = a + b + c
return if perim > @limit
@primitives += 1
@total += @limit // perim
generate_triples( a-2*b+2*c, 2*a-b+2*c, 2*a-2*b+3*c )
generate_triples( a+2*b+2*c, 2*a+b+2*c, 2*a+2*b+3*c )
generate_triples(-a+2*b+2*c,-2*a+b+2*c,-2*a+2*b+3*c )
end
end
perim = 10
while perim <= 100_000_000
c = PythagoranTriplesCounter.new perim
p [perim, c.total, c.primitives]
perim *= 10
end</syntaxhighlight>
output
<pre>[10, 0, 0]
[100, 17, 7]
[1000, 325, 70]
[10000, 4858, 703]
[100000, 64741, 7026]
[1000000, 808950, 70229]
[10000000, 9706567, 702309]
[100000000, 113236940, 7023027]</pre>
=={{header|D}}==
===Lazy Functional Version===
With hints from the Haskell solution.
<
import std.stdio, std.range, std.algorithm, std.typecons, std.numeric;
enum triples = (in
iota(1, n + 1)
.map!(
.joiner.joiner
.filter!(t => t[0] ^^ 2 + t[1] ^^ 2 == t[2] ^^ 2 && t[].only.sum <= n)
.map!(t => tuple(t[0 .. 2].gcd == 1, t[]));
auto xs = triples(100);
writeln("Up to 100 there are ", xs.
xs.filter!q{ a[0] }.
}</
{{out}}
<pre>Up to 100 there are 17 triples, 7 are primitive.</pre>
===Shorter Version===
<syntaxhighlight lang="d">ulong[2] tri(ulong lim, ulong a=3, ulong b=4, ulong c=5)
pure nothrow @safe @nogc {
immutable l = a + b + c;
return [0, 0];
typeof(return) r = [
r[] += tri(lim, a - 2*b + 2*c, 2*a - b + 2*c, 2*a - 2*b + 3*c)[];
r[] += tri(lim, a
r[] += tri(lim,
return r;
}
void main() /*@safe*/ {
import std.stdio;
writeln(10 ^^ p, ' ', tri(10 ^^ p));
}</syntaxhighlight>
{{out}}
<pre>10 [0, 0]
Line 580 ⟶ 1,290:
100000 [7026, 64741]
1000000 [70229, 808950]
10000000 [702309, 9706567]
100000000 [7023027, 113236940]</pre>
Run-time (32 bit system): about 0.80 seconds with ldc2.
===Short SIMD Version===
With LDC compiler this is a little faster than the precedent version (remove @nogc to compile it with the current version of LDC compiler).
<
ulong2 tri(in ulong lim, in ulong a=3, in ulong b=4, in ulong c=5)
pure nothrow @safe @nogc {
immutable l = a + b + c;
if (l > lim)
return [0, 0];
typeof(return) r = [1, lim / l];
r += tri(lim, a - 2*b + 2*c, 2*a - b + 2*c, 2*a - 2*b + 3*c);
r += tri(lim, a + 2*b + 2*c, 2*a + b + 2*c, 2*a + 2*b + 3*c);
r += tri(lim, -a + 2*b + 2*c, -2*a + b + 2*c, -2*a + 2*b + 3*c);
return r;
}
void main() /*@safe*/ {
foreach (immutable p; 1 ..
writeln(10 ^^ p,
}</
The output is the same. Run-time (32 bit system): about 0.67 seconds with ldc2.
===Faster Version===
{{trans|C}}
<
alias Xuint = uint; // ulong if going over 1 billion.
Line 645 ⟶ 1,357:
limit, nTriples, nPrimitives);
}
}</
{{out}}
<pre>Up to 10: 0 triples, 0 primitives.
Line 655 ⟶ 1,367:
Up to 10000000: 9706567 triples, 702309 primitives.
Up to 100000000: 113236940 triples, 7023027 primitives.</pre>
Run-time: about 0.27 seconds with ldc2.
Using the power p up to 11, using ulong for xuint, and compiling with the dmd <code>-L/STACK:10000000</code> switch to increase the stack size to about 10MB:
Line 681 ⟶ 1,394:
Up to 10000000000: 14557915466 triples, 702304875 primitives.
Up to 100000000000: 161750315680 triples, 7023049293 primitives.</pre>
=={{header|Delphi}}==
See [[#Pascal|Pascal]].
=={{header|EasyLang}}==
{{trans|C}}
<syntaxhighlight>
global total prim maxperi .
proc newtri s0 s1 s2 . .
p = s0 + s1 + s2
if p <= maxperi
prim += 1
total += maxperi div p
newtri s0 - 2 * s1 + 2 * s2 2 * s0 - s1 + 2 * s2 2 * s0 - 2 * s1 + 3 * s2
newtri s0 + 2 * s1 + 2 * s2 2 * s0 + s1 + 2 * s2 2 * s0 + 2 * s1 + 3 * s2
newtri -s0 + 2 * s1 + 2 * s2 -2 * s0 + s1 + 2 * s2 -2 * s0 + 2 * s1 + 3 * s2
.
.
for maxperi in [ 100 10000000 ]
prim = 0
total = 0
newtri 3 4 5
print "Up to " & maxperi & ": " & total & " triples, " & prim & " primitives"
.
</syntaxhighlight>
=={{header|EDSAC order code}}==
Not much optimization is done in this code, which is best run on a fast simulator
if the maximum perimeter is more than 10^5 or so.
A maximum perimeter of 10^7 takes 340 million EDSAC orders
(about 6 days on the original EDSAC).
The number of primitive triples divided by the maximum perimeter seems to tend to a limit,
which looks very much like (ln 2)/pi^2 = 0.07023049277 (see especially the FreeBASIC output below).
<syntaxhighlight lang="edsac">
[Pythagorean triples for Rosetta code.
Counts (1) all Pythagorean triples (2) primitive Pythagorean triples,
with perimeter not greater than a given value.
Library subroutine M3, Prints header and is then overwritten.
Here, the last character sets the teleprinter to figures.]
..PZ [simulate blank tape]
PFGKIFAFRDLFUFOFE@A6FG@E8FEZPF
@&*!MAX!PERIM!!!!!TOTAL!!!!!!PRIM@&#.
..PZ
[Library subroutine P7, prints long strictly positive integer;
10 characters, right justified, padded left with spaces.
Closed, even; 35 storage locations; working position 4D.]
T 56 K
GKA3FT26@H28#@NDYFLDT4DS27@TFH8@S8@T1FV4DAFG31@SFLDUFOFFFSFL4F
T4DA1FA27@G11@XFT28#ZPFT27ZP1024FP610D@524D!FO30@SFL8FE22@
[Subroutine for positive integer division.
Input: 4D = dividend, 6D = divisor.
Output: 4D = remainder, 6D = quotient.
37 locations; working locations 0D, 8D.]
T 100 K
GKA3FT35@A6DU8DTDA4DRDSDG13@T36@ADLDE4@T36@T6DA4DSDG23@
T4DA6DYFYFT6DT36@A8DSDE35@T36@ADRDTDA6DLDT6DE15@EFPF
[Subroutine to return GCD of two non-negative 35-bit integers.
Input: Integers at 4D, 6D.
Output: GCD at 4D; changes 6D.
41 locations; working location 0D.]
T 200 K
GKA3FT39@S4DE37@T40@A4DTDA6DRDSDG15@T40@ADLDE6@T40@A6DSDG20@T6D
T40@A4DSDE29@T40@ADRDTDE16@S6DE39@TDA4DT6DSDT4DE5@A6DT4DEFPF
[************************ ROSETTA CODE TASK *************************
Subroutine to count Pythagorean triples with given maximum perimeter.
Input: 0D = maximum perimeter.
Output: 4D = number of triples, 6D = number of primitive.
0D is changed.
Must be loaded at an even address.
Uses the well-known fact that a primitive Pythagorean triple is of the form
(m^2 - n^2, 2*m*n, m^2 + n^2) where m, n are coprime and of opposite parity.]
T 300 K
G K
A 3 F [make link]
E 16 @ [jump over variables and constants]
[Double values are put here to ensure even address]
[Variables]
[2] P F P F [maximum perimeter]
[4] P F P F [total number of Pythagorean triples]
[6] P F P F [number of primitive Pythagorean triples]
[8] P F P F [m]
[10] P F P F [n]
[Constants]
T12#Z PF T12Z [clears sandwich digit between 12 and 13]
[12] P D P F [double-value 1]
T14#Z PF T14Z [clears sandwich digit between 14 and 15]
[14] P1F P F [double-value 2]
[Continue with code]
[16] T 69 @ [plant link for return]
A D [load maximum perimeter]
T 2#@ [store locally]
T 4#@ [initialize counts of triangles to 0]
T 6#@
A 12#@ [load 1]
T 8#@ [m := 1]
[Next m, inc by 1]
[23] T F [clear acc]
A 8#@ [load m]
A 12#@ [add 1]
T 8#@ [update m]
H 8#@ [mult reg := m]
C 12#@ [acc := m AND 1]
A 12#@ [add 1]
T 10#@ [n := 1 if m even, 2 if m odd]
[Here to count triangles arising from m, n.
It's assumed m and n are known coprime.]
[31] A 31 @ [call the count subroutine,]
G 70 @ [result is in 6D]
S 6 D [load negative count]
G 40 @ [jump if count > 0]
[No triangles found for this n.
If n = 1 or 2 then whole thing is finished.
Else move on to next m.]
T F [clear acc]
A 14#@ [load 2]
S 10#@ [2 - n]
G 23 @ [if n > 2, go to next m]
E 64 @ [if n <= 2, exit]
[Found triangles, count is in 6D]
[40] T F [clear acc]
A 4#@ [load total count]
A 6 D [add count just found]
T 4#@ [update total count]
A 6#@ [load primitive count]
A 12#@ [add 1]
T 6#@ [update primitive count]
[47] T F [clear acc]
A 10#@ [load n]
A 14#@ [add 2]
U 10#@ [update n]
S 8#@ [is n > m?]
E 23 @ [if so, loop back for next m]
[Test whether m and n are coprime.]
T F [clear acc]
A 8#@ [load m]
T 4 D [to 4D for GCD routine]
A 10#@ [load n]
T 6 D [to 6D for GCD routine]
A 58 @ [call GCD routine,]
G 200 F [GCD is returned in 4D]
A 4 D [load GCD]
S 14#@ [is GCD = 1? (test by subtracting 2)]
E 47 @ [no, go straight to next n]
G 31 @ [yes, count triangles, then next n]
[64] T F [exit, clear acc]
A 4#@ [load total number of triples]
T 4 D [return in 4D]
A 6#@ [load number of primitive triples]
T 6 D [return in 6D]
[69] E F
[2nd-level subroutine to count triangles arising from m, n.
Assumes m, n are coprime and of opposite parity,
and m is in the multiplier register.
Result is returned in 6D.]
[70] A 3 F [make and plant link for return]
T 91 @
A 2#@ [acc := maximum perimeter]
T 4 D [to 4D for division routine]
A 8#@ [load m]
A 10#@ [add n]
T D [m + n to 0D]
V D [acc := m*(m + n)]
[Need to shift product 34 left to restore integer scaling.
Since we want 2*m*(m+n), shift 35 left.]
L F [13 left (maximum possible)]
L F [13 more]
L 128 F [9 more]
T 6 D [perimeter to 6D for division routine]
A 4 D [load maximum perimeter]
S 6 D [is perimeter > maximum?]
G 89 @ [quick exit if so]
T F [clear acc]
A 86 @ [call division routine,]
G 100 F [leaves count in 6D]
E 91 @ [jump to exit]
[89] T F [acc := 0]
T 6 D [return count = 0]
[91] E F
[Main routine. Load at an even address.]
T 500 K
G K
[The initial maximum perimeter is repeatedly multiplied by 10]
T#Z PF TZ [clears sandwich digit between 0 and 1]
[0] P50F PF [initial maximum perimeter <---------- EDIT HERE]
[2] P 3 F [number of values to calculate <---------- EDIT HERE]
[3] P D [1]
[4] P F P F [maximum perimeter]
[6] P F P F [total number of triples]
[8] P F P F [number of primitive triples]
[10] P F [negative count of values]
[11] # F [figures shift]
[12] @ F [carriage return]
[13] & F [line feed]
[14] K 4096 F [null char]
[Enter with acc = 0]
[15] S 2 @ [initialize a negative counter]
T 10 @ [(standard EDSAC practice)]
A #@ [initialize maximum perimeter]
T 4#@
[19] T F [clear acc]
A 4#@ [load maximum perimeter]
T D [to 0D for subroutine]
A 22 @ [call subroutine to count triples]
G 300 F
A 4 D [returns total number in 4D]
T 6#@ [save locally]
A 6 D [returns number of primitive in 6D]
T 8#@ [save locally]
[Print the result]
A 4#@ [load maximum perimeter]
T D [to 0D for print subroutine]
A 30 @ [call print subroutine]
G 56 F
A 6#@ [repeat for total number of triples]
T D
A 34 @
G 56 F
A 8#@ [repeat for number of primitive triples]
T D
A 38 @
G 56 F
O 12 @
O 13 @
A 10 @ [load negative count]
A 3 @ [add 1]
E 53 @ [out if reached 0]
T 10 @ [else update count]
A 4#@ [load max perimeter]
U D [temp store]
L 1 F [times 4]
A D [times 5]
L D [times 10]
T 4#@ [update]
E 19 @ [loop back]
[53] O 14 @ [done; print null to flush printer buffer]
Z F [stop]
E 15 Z [define entry point]
P F [acc = 0 on entry]
</syntaxhighlight>
{{out}}
<pre>
MAX PERIM TOTAL PRIM
100 17 7
1000 325 70
10000 4858 703
100000 64741 7026
1000000 808950 70229
10000000 9706567 702309
</pre>
=={{header|Eiffel}}==
<syntaxhighlight lang="eiffel">
class
APPLICATION
create
make
feature
make
local
perimeter: INTEGER
do
perimeter := 100
from
until
perimeter > 1000000
loop
total := 0
primitive_triples := 0
count_pythagorean_triples (3, 4, 5, perimeter)
io.put_string ("There are " + total.out + " triples, below " + perimeter.out + ". Of which " + primitive_triples.out + " are primitives.%N")
perimeter := perimeter * 10
end
end
count_pythagorean_triples (a, b, c, perimeter: INTEGER)
-- Total count of pythagorean triples and total count of primitve triples below perimeter.
local
p: INTEGER
do
p := a + b + c
if p <= perimeter then
primitive_triples := primitive_triples + 1
total := total + perimeter // p
count_pythagorean_triples (a + 2 * (- b + c), 2 * (a + c) - b, 2 * (a - b + c) + c, perimeter)
count_pythagorean_triples (a + 2 * (b + c), 2 * (a + c) + b, 2 * (a + b + c) + c, perimeter)
count_pythagorean_triples (- a + 2 * (b + c), 2 * (- a + c) + b, 2 * (- a + b + c) + c, perimeter)
end
end
feature {NONE}
primitive_triples: INTEGER
total: INTEGER
end
</syntaxhighlight>
{{out}}
<pre>
There are 17 triples, below 100. Of which 7 are primitives.
There are 325 triples, below 1000. Of which 70 are primitives.
There are 4858 triples, below 10000. Of which 703 are primitives.
There are 64741 triples, below 100000. Of which 7026 are primitives.
There are 808950 triples, below 1000000. Of which 70229 are primitives.
</pre>
=={{header|Elixir}}==
{{trans|Ruby}}
<syntaxhighlight lang="elixir">defmodule RC do
def count_triples(limit), do: count_triples(limit,3,4,5)
defp count_triples(limit, a, b, c) when limit<(a+b+c), do: {0,0}
defp count_triples(limit, a, b, c) do
{p1, t1} = count_triples(limit, a-2*b+2*c, 2*a-b+2*c, 2*a-2*b+3*c)
{p2, t2} = count_triples(limit, a+2*b+2*c, 2*a+b+2*c, 2*a+2*b+3*c)
{p3, t3} = count_triples(limit,-a+2*b+2*c,-2*a+b+2*c,-2*a+2*b+3*c)
{1+p1+p2+p3, div(limit, a+b+c)+t1+t2+t3}
end
end
list = for n <- 1..8, do: Enum.reduce(1..n, 1, fn(_,acc)->10*acc end)
Enum.each(list, fn n -> IO.inspect {n, RC.count_triples(n)} end)</syntaxhighlight>
{{out}}
<pre>
{10, {0, 0}}
{100, {7, 17}}
{1000, {70, 325}}
{10000, {703, 4858}}
{100000, {7026, 64741}}
{1000000, {70229, 808950}}
{10000000, {702309, 9706567}}
{100000000, {7023027, 113236940}}
</pre>
=={{header|Erlang}}==
<syntaxhighlight lang="erlang">%%
%% Pythagorian triples in Erlang, J.W. Luiten
%%
Line 728 ⟶ 1,788:
L = lists:seq(1, Max),
Answer = lists:map(fun(X) -> count_triples(X) end, L),
lists:foreach(fun(Result) -> display_result(Result) end, Answer).</
Output:
Line 744 ⟶ 1,804:
Up to 10 ** 11 : 161750315680 triples, 7023049293 primitives
</pre>
=={{header|ERRE}}==
<syntaxhighlight lang="erre">PROGRAM PIT
BEGIN
PRINT(CHR$(12);) !CLS
PRINT(TIME$)
FOR POWER=1 TO 7 DO
PLIMIT=10#^POWER
UPPERBOUND=INT(1+PLIMIT^0.5)
PRIMITIVES=0
TRIPLES=0
EXTRAS=0 ! will count the in-range multiples of any primitive
FOR M=2 TO UPPERBOUND DO
FOR N=1+(M MOD 2=1) TO M-1 STEP 2 DO
TERM1=2*M*N
TERM2=M*M-N*N
TERM3=M*M+N*N
PERIMETER=TERM1+TERM2+TERM3
IF PERIMETER<=PLIMIT THEN TRIPLES=TRIPLES+1
A=TERM1
B=TERM2
REPEAT
R=A-B*INT(A/B)
A=B
B=R
UNTIL R<=0
! we've found a primitive triple if a = 1, since hcf =1.
! and it is inside perimeter range. Save it in an array
IF (A=1) AND (PERIMETER<=PLIMIT) THEN
PRIMITIVES=PRIMITIVES+1
!-----------------------------------------------
!swap so in increasing order of side length
!-----------------------------------------------
IF TERM1>TERM2 THEN SWAP(TERM1,TERM2)
!-----------------------------------------------
!we have the primitive & removed any multiples.
!Now calculate ALL the multiples in range.
!-----------------------------------------------
NEX=INT(PLIMIT/PERIMETER)
EXTRAS=EXTRAS+NEX
END IF
!scan
END FOR
END FOR
PRINT("Primit. with perimeter <=";10#^power;"is";primitives;"&";extras;"non-prim.triples.")
PRINT(TIME$)
END FOR
PRINT PRINT("** End **")
END PROGRAM</syntaxhighlight>
{{out}}
<pre>16:08:39
Primit. with perimeter <= 10 is 0 & 0 non-prim.triples.
16:08:39
Primit. with perimeter <= 100 is 7 & 17 non-prim.triples.
16:08:39
Primit. with perimeter <= 1000 is 70 & 325 non-prim.triples.
16:08:39
Primit. with perimeter <= 10000 is 703 & 4858 non-prim.triples.
16:08:39
Primit. with perimeter <= 100000 is 7026 & 64741 non-prim.triples.
16:08:41
Primit. with perimeter <= 1000000 is 70229 & 808950 non-prim.triples.
16:09:07
Primit. with perimeter <= 10000000 is 702309 & 9706567 non-prim.triples.
16:13:10
** End **</pre>
=={{header|Euphoria}}==
{{trans|D}}
<
sequence r
atom p
Line 767 ⟶ 1,906:
? tri(max_peri, {3, 4, 5})
max_peri *= 10
end while</
Output:
Line 779 ⟶ 1,918:
100000000: {7023027,113236940}
</pre>
=={{header|F_Sharp|F#}}==
{{trans|OCaml}}
<
let rec iter t =
let d = n - t*t
Line 814 ⟶ 1,952:
printfn "For perimeters up to %d there are %d total and %d primitive" i s p;;
List.iter show [ 100; 1000; 10000; 100000; 1000000; 10000000; 100000000 ]</
{{out}}
<pre>For perimeters up to 100 there are 17 total and 7 primitive
Line 827 ⟶ 1,965:
Pretty slow (100 times slower than C)...
<
math.functions math.matrices math.ranges sequences ;
IN: rosettacode.pyth
Line 869 ⟶ 2,007:
"Up to %d: %d triples, %d primitives.\n" printf ;
: pyth ( -- )
8 [1,b] [ 10^ dup count-triplets pprint-triplet-count ] each ;</
<pre>Up to 10: 0 triples, 0 primitives.
Line 880 ⟶ 2,018:
Up to 100000000: 113236940 triples, 7023027 primitives.
Running time: 57.968821207 seconds</pre>
=={{header|Forth}}==
<syntaxhighlight lang="forth ">
Line 1,059 ⟶ 2,195:
ok
</syntaxhighlight>
=={{header|Fortran}}==
{{works with|Fortran|90 and later}}
{{trans|C efficient method}}
<
implicit none
Line 1,109 ⟶ 2,245:
max_peri = max_peri * 10
end do
end program Pythagorean</
Output:<pre>Up to 10 0 triples 0 primitives
Up to 100 17 triples 7 primitives
Line 1,119 ⟶ 2,255:
Up to 100000000 113236940 triples 7023027 primitives</pre>
=={{header|
The upper limit is set to 12(10^12), this will take about 3-4 hours.
If you can't wait that long better lower it to 11(10^11).
===Version 1===
Normal version
<syntaxhighlight lang="freebasic">' version 30-05-2016
' compile with: fbc -s console
' primitive pythagoras triples
' a = m^2 - n^2, b = 2mn, c = m^2 + n^2
' m, n are positive integers and m > n
' m -
' p = a
' max m for give perimeter
' p = m^2 - n^2 + 2mn + m^2 + n^2
' p = 2mn + m^2 + m^2 + n^2 - n^2 = 2mn + 2m^2
' m >> n and n = 1 ==> p = 2m + 2m^2 = 2m(1 + m)
' m >> 1 ==> p = 2m(m) = 2m^2
' max m for given perimeter = sqr(p / 2)
Function gcd(x As UInteger, y As UInteger) As UInteger
Dim As UInteger t
x = t
Wend
Return x
End Function
Sub pyth_trip(limit As ULongInt, ByRef trip As ULongInt, ByRef prim As ULongInt)
Dim As ULongInt perimeter, lby2 = limit Shr 1
Dim As UInteger m, n
Dim As ULongInt a, b, c
For m = 2 To Sqr(limit / 2)
For n = 1 + (m And 1) To (m - 1) Step 2
' common divisor, try next n
If (gcd(m, n) > 1) Then Continue For
a = CULngInt(m) * m - n * n
b = CULngInt(m) * n * 2
c = CULngInt(m) * m + n * n
perimeter = a + b + c
' perimeter > limit, since n goes up try next m
If perimeter >= limit Then Continue For, For
prim += 1
If perimeter < lby2 Then
trip += limit \ perimeter
Else
trip += 1
End If
Next n
Next m
End Sub
' ------=< MAIN >=------
Dim As String str1, buffer = Space(14)
Dim As ULongInt limit, trip, prim
Dim As Double t, t1 = Timer
Print "below triples primitive time"
Print
For x As UInteger = 1 To 12
t = Timer
limit = 10 ^ x : trip = 0 : prim = 0
pyth_trip(limit, trip, prim)
LSet buffer, Str(prim) : str1 = buffer
Print Using "10^## ################ "; x; trip;
If x > 7 Then
Print str1;
Print Using " ######.## sec."; Timer - t
Else
Print str1
End If
Next x
Print : Print
Print Using "Total time needed #######.## sec."; Timer - t1
' empty keyboard buffer
While InKey <> "" : Wend
Print : Print "hit any key to end program"
Sleep
End</syntaxhighlight>
{{out}}
<pre>below triples primitive time
10^ 1 0 0
10^ 2 17 7
10^ 3 325 70
10^ 4 4858 703
10^ 5 64741 7026
10^ 6 808950 70229
10^ 7 9706567 702309
10^ 8 113236940 7023027 0.94 sec.
10^ 9 1294080089 70230484 10.13 sec.
10^10 14557915466 702304875 109.75 sec.
10^11 161750315680 7023049293 1204.56 sec.
10^12 1779214833461 70230492763 13031.84 sec.
Total time needed 14357.31 sec.</pre>
===Version 2===
Attempt to make a faster version (about 20% faster)
<syntaxhighlight lang="freebasic">' version 30-05-2016
' compile with: fbc -s console
' max m for give perimeter
' p = m^2 - n^2 + 2mn + m^2 + n^2
' p = 2mn + m^2 + m^2 + n^2 - n^2 = 2mn + 2m^2
' m >> n and n = 1 ==> p = 2m + 2m^2 = 2m(1 + m)
' m >> 1 ==> p = 2m(m) = 2m^2
' max m for given perimeter = sqr(p / 2)
Function gcd(x As UInteger, y As UInteger) As UInteger
Dim As UInteger t
While y
t = y
y = x Mod y
x = t
Wend
Return x
End Function
Sub pyth_trip_fast(limit As ULongInt, ByRef trip As ULongInt, ByRef prim As ULongInt)
Dim As ULongInt perimeter, lby2 = limit Shr 1
Dim As UInteger mx2 = 4
For m As UInteger = 2 To Sqr(limit / 2)
perimeter = (CULngInt(m) * m * 2) - IIf(m And 1, 0, m * 2)
mx2 = mx2 + 4
For n As UInteger = 1 + (m And 1) To (m - 1) Step 2
perimeter += mx2
' common divisor, try next n
If (gcd(m, n) > 1) Then Continue For
' perimeter > limit, since n goes up try next m
If perimeter >= limit Then Continue For, For
prim += 1
If perimeter < lby2 Then
trip += limit \ perimeter
Else
trip += 1
End If
Next n
Next m
End Sub
' ------=< MAIN >=------
Dim As String str1, buffer = Space(14)
Dim As ULongInt limit, trip, prim
Dim As Double t, t1 = Timer
Print "below triples primitive time"
Print
For x As UInteger = 1 To 12
t = Timer
limit = 10 ^ x : trip = 0 : prim = 0
pyth_trip_fast(limit, trip, prim)
LSet buffer, Str(prim) : str1 = buffer
Print Using "10^## ################ "; x; trip;
If x > 7 Then
Print str1;
Print Using " ######.## sec."; Timer - t
Else
Print str1
End If
Next x
Print : Print
Print Using "Total time needed #######.## sec."; Timer - t1
' empty keyboard buffer
While InKey <> "" : Wend
Print : Print "hit any key to end program"
Sleep
End</syntaxhighlight>
{{out}}
<pre>below triples primitive time
10^ 1 0 0
10^ 2 17 7
10^ 3 325 70
10^ 4 4858 703
10^ 5 64741 7026
10^ 6 808950 70229
10^ 7 9706567 702309
10^ 8 113236940 7023027 0.66 sec.
10^ 9 1294080089 70230484 7.48 sec.
10^10 14557915466 702304875 83.92 sec.
10^11 161750315680 7023049293 945.95 sec.
10^12 1779214833461 70230492763 10467.94 sec.
Total time needed 11506.01 sec.</pre>
The time needed is about 11 times the time needed for the previous limit.
To calculate 10^12 with the fast version that would take about 32 hours.
The variable that holds the number of triples will eventual overflow
at 10^18 - 10^19. To reach that stage you need the program to run
for a few years.
=={{header|Go}}==
<
import "fmt"
Line 1,202 ⟶ 2,503:
maxPeri, total, prim)
}
}</
Output:
<pre>
Line 1,216 ⟶ 2,517:
Up to 100000000000: 161750315680 triples, 7023049293 primitives
</pre>
=={{header|Groovy}}==
===Parent/Child Algorithm===
Solution:
<syntaxhighlight lang="groovy">class Triple {
BigInteger a, b, c
def getPerimeter() { this.with { a + b + c } }
boolean isValid() { this.with { a*a + b*b == c*c } }
}
def initCounts (def n = 10) {
(n..1).collect { 10g**it }.inject ([:]) { Map map, BigInteger perimeterLimit ->
map << [(perimeterLimit): [primative: 0g, total: 0g]]
}
}
def findPythagTriples, findChildTriples
findPythagTriples = {Triple t = new Triple(a:3, b:4, c:5), Map counts = initCounts() ->
def p = t.perimeter
def currentCounts = counts.findAll { pLimit, tripleCounts -> p <= pLimit }
if (! currentCounts || ! t.valid) { return }
currentCounts.each { pLimit, tripleCounts ->
tripleCounts.with { primative ++; total += pLimit.intdiv(p) }
}
findChildTriples(t, currentCounts)
counts
}
findChildTriples = { Triple t, Map counts ->
t.with {
[
[ a - 2*b + 2*c, 2*a - b + 2*c, 2*a - 2*b + 3*c],
[ a + 2*b + 2*c, 2*a + b + 2*c, 2*a + 2*b + 3*c],
[-a + 2*b + 2*c, -2*a + b + 2*c, -2*a + 2*b + 3*c]
]*.sort().each { aa, bb, cc ->
findPythagTriples(new Triple(a:aa, b:bb, c:cc), counts)
}
}
}</syntaxhighlight>
Test:
<syntaxhighlight lang="groovy">printf (' LIMIT PRIMATIVE ALL\n')
findPythagTriples().sort().each { perimeterLimit, result ->
def exponent = perimeterLimit.toString().size() - 1
printf ('a+b+c <= 10E%2d %9d %12d\n', exponent, result.primative, result.total)
}</syntaxhighlight>
Output:
<pre> LIMIT PRIMATIVE ALL
a+b+c <= 10E 1 0 0
a+b+c <= 10E 2 7 17
a+b+c <= 10E 3 70 325
a+b+c <= 10E 4 703 4858
a+b+c <= 10E 5 7026 64741
a+b+c <= 10E 6 70229 808950
a+b+c <= 10E 7 702309 9706567
a+b+c <= 10E 8 7023027 113236940
a+b+c <= 10E 9 70230484 1294080089
a+b+c <= 10E10 702304875 14557915466</pre>
=={{header|Haskell}}==
<syntaxhighlight lang="haskell">pytr :: Int -> [(Bool, Int, Int, Int)]
pytr n =
filter
(\(_, a, b, c) -> a + b + c <= n)
[ (prim a b c, a, b, c)
| a <- xs,
b <- drop a xs,
c <- drop b xs,
a ^ 2 + b ^ 2 == c ^ 2
]
where
xs = [1 .. n]
prim a b _ = gcd a b == 1
main :: IO ()
main =
putStrLn $
"Up to 100 there are "
<> show (length xs)
<> " triples, of which "
<> show (length $ filter (\(x, _, _, _) -> x) xs)
<> " are primitive."
where
xs = pytr 100</syntaxhighlight>
{{Out}}
<pre>Up to 100 there are 17 triples, of which 7 are primitive.</pre>
Or equivalently (desugaring the list comprehension down to nested concatMaps, and pruning back the search space a little):
<syntaxhighlight lang="haskell">------------------- PYTHAGOREAN TRIPLES ------------------
pythagoreanTriplesBelow :: Int -> [[Int]]
pythagoreanTriplesBelow n =
concatMap
( \x ->
concatMap
(\y -> concatMap (go x y) [y + 1 .. m])
[x + 1 .. m]
)
[1 .. m]
where
m = quot n 2
go x y z
| x + y + z <= n && x ^ 2 + y ^ 2 == z ^ 2 =
[[x, y, z]]
| otherwise = []
--------------------------- TEST -------------------------
main :: IO ()
main =
mapM_
(print . length)
( [id, filter (\[x, y, _] -> gcd x y == 1)]
<*> [pythagoreanTriplesBelow 100]
)</syntaxhighlight>
{{Out}}
<pre>17
7</pre>
Recursive primitive generation:
<syntaxhighlight lang="haskell">triangles :: Int -> [[Int]]
triangles max_peri
| max_peri < 12 = []
| otherwise = concat tiers
where
tiers = takeWhile (not . null) $ iterate tier [[3, 4, 5]]
tier = concatMap (filter ((<= max_peri) . sum) . tmul)
tmul t =
map
(map (sum . zipWith (*) t))
[ [[1, -2, 2], [2, -1, 2], [2, -2, 3]],
[[1, 2, 2], [2, 1, 2], [2, 2, 3]],
[[-1, 2, 2], [-2, 1, 2], [-2, 2, 3]]
]
triangleCount max_p = (length t, sum $ map ((max_p `div`) . sum) t)
where
t = triangles max_p
main :: IO ()
main =
mapM_
((putStrLn . (\n -> show n <> " " <> show (triangleCount n))) . (10 ^))
[1 .. 7]</syntaxhighlight>
{{out}}
<pre>10 (0,0)
100 (7,17)
1000 (70,325)
10000 (703,4858)
100000 (7026,64741)
1000000 (70229,808950)
10000000 (702309,9706567)</pre>
=={{header|Icon}} and {{header|Unicon}}==
This uses the elegant formula (#IV) from [[wp:Formulas_for_generating_Pythagorean_triples|Formulas for generating Pythagorean triples]]
<syntaxhighlight lang="icon">
link numbers
link printf
Line 1,258 ⟶ 2,714:
every (s := "") ||:= !sort(x) do s ||:= ","
return s[1:-1]
end</
{{libheader|Icon Programming Library}}
Line 1,287 ⟶ 2,743:
Under perimiter=10000000: Pythagorean Triples=9706567 including primitives=702309
Time=560625, Collections: total=16 string=8 block=8</pre>
=={{header|J}}==
Line 1,331 ⟶ 2,748:
Brute force approach:
<
r=. i. 0 3
for_a. 1 + i. <.(y-1)%3 do.
Line 1,344 ⟶ 2,761:
)
prim=: 1 = 2 +./@{. |:</
Example use:
Line 1,350 ⟶ 2,767:
First column indicates whether the triple is primitive, and the remaining three columns are a, b and c.
<
1 3 4 5
1 5 12 13
Line 1,375 ⟶ 2,792:
325 70
(# , [: {. +/) pytr 10000
4858 703</
pytr 10000 takes 4 seconds on this laptop, and time to complete grows with square of perimeter, so pytr 1e6 should take something like 11 hours using this algorithm on this machine.
Line 1,381 ⟶ 2,798:
A slightly smarter approach:
<
'm n'=. |:(#~ 1 = 2 | +/"1)(#~ >/"1) ,/ ,"0/~ }. i. <. %: y
prim=. (#~ 1 = 2 +./@{. |:) (#~ y >: +/"1)m (-&*: ,. +:@* ,. +&*:) n
/:~ ; <@(,.~ # {. 1:)@(*/~ 1 + y i.@<.@% +/)"1 prim
)</
usage for trips is the same as for pytr. Thus:
<
0 0
(# , 1 {. +/) trips 100
Line 1,402 ⟶ 2,819:
808950 70229
(# , 1 {. +/) trips 10000000
9706567 702309</
The last line took about 16 seconds.
Line 1,408 ⟶ 2,825:
That said, we do not actually have to generate all the triples, we just need to count them. Thus:
<
'm n'=. |:(#~ 1 = 2 | +/"1)(#~ >/"1) ,/ ,"0/~ }. i. <. %: y
<.y%+/"1 (#~ 1 = 2 +./@{. |:) (#~ y >: +/"1)m (-&*: ,. +:@* ,. +&*:) n
)</
The result is a list of positive integers, one number for each primitive triple which fits within the limit, giving the number of triples which are multiples of that primitive triple whose perimeter is no greater than the limiting perimeter.
<syntaxhighlight lang="text"> (#,+/)trc 1e8
7023027 113236940</
But note that J's memory footprint reached 6.7GB during the computation, so to compute larger values the computation would have to be broken up into reasonable sized blocks.
===Traversal of the Tree of Primitive Pythagorean Triples===
On my laptop this code takes 1.35 seconds for a perimeter of up to 1 million, though it takes 2 minutes for 10 million, so performance is between the previous code "brute force" and what is called "slightly smarter approach" It could probably be sped up slightly by not sorting the triples (there is no need with this problem.)
<syntaxhighlight lang="j">
mp =: +/ . * "2 1
T =: 3 3 3$ 1 _2 2 2 _1 2 2 _2 3 1 2 2 2 1 2 2 2 3 _1 2 2 _2 1 2 _2 2 3
branch =: dyad define NB. Go down one branch of the tree, usage: <perimeter> branch <triple>
(x >: +/"1 next) # next =. T (/:~ @ mp) y
)
pythag =: monad define NB. pythagorean triples with max perimeter
t1 =. 0 3$ 0
if. y >: 12 do.
t0 =. 1 3$ 3 4 5
while. #t0 > 0 do.
t =. {. t0
t1 =. t1, t
t0 =. (}. t0), y branch t
end.
end.
/:~ t1
)
count =: monad define "0 NB. count triples with max perimeter
y, (#t), +/ <. y % +/"1 t =. pythag y
)
(9!:11) 7 NB. change output precision
echo 'Counts of primitive and total number of Pythagorean triples with perimeter ≤ 10^n.'
echo count 10 ^ >: i.6
exit ''
</syntaxhighlight>
{{Out}}
<pre>
Counts of primitive and total number of Pythagorean triples with perimeter ≤ 10^n.
10 0 0
100 7 17
1000 70 325
10000 703 4858
100000 7026 64741
1000000 70229 808950
</pre>
=={{header|Java}}==
===Brute force===
[[Category:Arbitrary precision]]Theoretically, this can go "forever", but it takes a while, so only the minimum is shown. Luckily, <code>BigInteger</code> has a GCD method built in.
<
import java.math.BigInteger;
import static java.math.BigInteger.ONE;
Line 1,474 ⟶ 2,935:
+ tripCount + " triples, of which " + primCount + " are primitive.");
}
}</
Output:
<pre>3, 4, 5 primitive
Line 1,497 ⟶ 2,958:
[[Pythagorean triples/Java/Brute force primitives]]
===Parent/child===
{{trans|
{{works with|Java|1.5+}}
This can also go "forever" theoretically. Letting it go to another order of magnitude overflowed the stack on the computer this was tested on. This version also does not show the triples as it goes, it only counts them.
<
public class Triples{
Line 1,537 ⟶ 2,998:
}
}
}</
Output:
<pre>100: 17 triples, 7 primitive.
Line 1,546 ⟶ 3,007:
10000000: 9706567 triples, 702309 primitive.</pre>
=={{header|JavaScript}}==
===ES6===
Exhaustive search of a full cartesian product. Not scalable.
<syntaxhighlight lang="javascript">(() => {
"use strict";
// Arguments: predicate, maximum perimeter
// pythTripleCount :: ((Int, Int, Int) -> Bool) -> Int -> Int
const pythTripleCount = p =>
maxPerim => {
const
xs = enumFromTo(1)(
Math.floor(maxPerim / 2)
);
return xs.flatMap(
x => xs.slice(x).flatMap(
y => xs.slice(y).flatMap(
z => ((x + y + z <= maxPerim) &&
((x * x) + (y * y) === z * z) &&
p(x, y, z)) ? [
[x, y, z]
] : []
)
)
).length;
};
// ---------------------- TEST -----------------------
const main = () => [10, 100, 1000]
.map(n => ({
maxPerimeter: n,
triples: pythTripleCount(() => true)(n),
primitives: pythTripleCount(
(x, y) => gcd(x)(y) === 1
)(n)
}));
// ---------------- GENERIC FUNCTIONS ----------------
// abs :: Num -> Num
const abs =
// Absolute value of a given number
// without the sign.
x => 0 > x ? (
-x
) : x;
// enumFromTo :: Int -> Int -> [Int]
const enumFromTo = m =>
n => Array.from({
length: 1 + n - m
}, (_, i) => m + i);
// gcd :: Integral a => a -> a -> a
const gcd = x =>
y => {
const zero = x.constructor(0);
const go = (a, b) =>
zero === b ? (
a
) : go(b, a % b);
return go(abs(x), abs(y));
};
// MAIN ---
return main();
})();</syntaxhighlight>
{{Out}}
<syntaxhighlight lang="javascript">[{"maxPerimeter":10, "triples":0, "primitives":0},
{"maxPerimeter":100, "triples":17, "primitives":7},
{"maxPerimeter":1000, "triples":325, "primitives":70}]</syntaxhighlight>
=={{header|jq}}==
{{Works with |jq |1.4}}
The jq program presented here is based on Euclid's formula, and uses the same algorithm and notation as in the AutoHotKey section.
The implementation illustrates how an inner function with arity 0 can
attain a high level of efficiency with both jq 1.4 and later. A simpler implementation is possible with versions of jq greater than 1.4.
<syntaxhighlight lang="jq">def gcd(a; b):
def _gcd:
if .[1] == 0 then .[0]
else [.[1], .[0] % .[1]] | _gcd
end;
[a,b] | _gcd ;
# Return: [total, primitives] for pythagorean triangles having
# perimeter no larger than peri.
# The following uses Euclid's formula with the convention: m > n.
def count(peri):
# The inner function can be used to count for a given value of m:
def _count:
# state [n,m,p, [total, primitives]]
.[0] as $n | .[1] as $m | .[2] as $p
| if $n < $m and $p <= peri then
if (gcd($m;$n) == 1)
then .[3] | [ (.[0] + ((peri/$p)|floor) ), (.[1] + 1)]
else .[3]
end
| [$n+2, $m, $p+4*$m, .] | _count
else .
end;
# m^2 < m*(m+1) <= m*(m+n) = perimeter/2
reduce range(2; (peri/2) | sqrt + 1) as $m
( [1, 2, 12, [0,0]];
(($m % 2) + 1) as $n
| (2 * $m * ($m + $n) ) as $p # a+b+c for this (m,n)
| [$n, $m, $p, .[3]] | _count
) | .[3] ;
# '''Example''':
def pow(i): . as $in | reduce range(0; i) as $j (1; . * $in);
range(1; 9) | . as $i | 10|pow($i) as $i | "\($i): \(count($i) )"
</syntaxhighlight>
{{Out}}
<syntaxhighlight lang="sh">$ jq -M -c -r -n -f Pythagorean_triples.jq
10: [0,0]
100: [17,7]
1000: [325,70]
10000: [4858,703]
100000: [64741,7026]
1000000: [808950,70229]
10000000: [9706567,702309]
100000000: [113236940,7023027]
</syntaxhighlight>
=={{header|Julia}}==
This solution uses the the Euclidian concept of m and n as generators of Pythagorean triplets. When m and n are coprime and have opposite parity, the generated triplets are primitive. It works reasonably well up to a limit of 10^10.
<syntaxhighlight lang="julia">
function primitiven{T<:Integer}(m::T)
1 < m || return T[]
m != 2 || return T[1]
!isprime(m) || return T[2:2:m-1]
rp = trues(m-1)
if isodd(m)
rp[1:2:m-1] = false
end
for p in keys(factor(m))
rp[p:p:m-1] = false
end
T[1:m-1][rp]
end
function pythagoreantripcount{T<:Integer}(plim::T)
primcnt = 0
fullcnt = 0
11 < plim || return (primcnt, fullcnt)
for m in 2:plim
p = 2m^2
p+2m <= plim || break
for n in primitiven(m)
q = p + 2m*n
q <= plim || break
primcnt += 1
fullcnt += div(plim, q)
end
end
return (primcnt, fullcnt)
end
println("Counting Pythagorian Triplets within perimeter limits:")
println(" Limit All Primitive")
for om in 1:10
(pcnt, fcnt) = pythagoreantripcount(10^om)
println(@sprintf " 10^%02d %11d %9d" om fcnt pcnt)
end
</syntaxhighlight>
{{out}}
<pre>
Counting Pythagorian Triplets within perimeter limits:
Limit All Primitive
10^01 0 0
10^02 17 7
10^03 325 70
10^04 4858 703
10^05 64741 7026
10^06 808950 70229
10^07 9706567 702309
10^08 113236940 7023027
10^09 1294080089 70230484
10^10 14557915466 702304875
</pre>
=={{header|Kotlin}}==
{{trans|Go}}
Due to deep recursion, I needed to increase the stack size to 4MB to get up to a maximum perimeter of 10 billion. Expect a run time of around 30 seconds on a typical laptop.
<syntaxhighlight lang="scala">// version 1.1.2
var total = 0L
var prim = 0L
var maxPeri = 0L
fun newTri(s0: Long, s1: Long, s2: Long) {
val p = s0 + s1 + s2
if (p <= maxPeri) {
prim++
total += maxPeri / p
newTri( s0 - 2 * s1 + 2 * s2, 2 * s0 - s1 + 2 * s2, 2 * s0 - 2 * s1 + 3 * s2)
newTri( s0 + 2 * s1 + 2 * s2, 2 * s0 + s1 + 2 * s2, 2 * s0 + 2 * s1 + 3 * s2)
newTri(-s0 + 2 * s1 + 2 * s2, -2 * s0 + s1 + 2 * s2, -2 * s0 + 2 * s1 + 3 * s2)
}
}
fun main(args: Array<String>) {
maxPeri = 100
while (maxPeri <= 10_000_000_000L) {
prim = 0
total = 0
newTri(3, 4, 5)
println("Up to $maxPeri: $total triples, $prim primatives")
maxPeri *= 10
}
}</syntaxhighlight>
{{out}}
<pre>
Up to 100: 17 triples, 7 primatives
Up to 1000: 325 triples, 70 primatives
Up to 10000: 4858 triples, 703 primatives
Up to 100000: 64741 triples, 7026 primatives
Up to 1000000: 808950 triples, 70229 primatives
Up to 10000000: 9706567 triples, 702309 primatives
Up to 100000000: 113236940 triples, 7023027 primatives
Up to 1000000000: 1294080089 triples, 70230484 primatives
Up to 10000000000: 14557915466 triples, 702304875 primatives
</pre>
=={{header|Lasso}}==
<
define num_pythagorean_triples(max_perimeter::integer) => {
local(max_b) = (#max_perimeter / 3)*2
Line 1,566 ⟶ 3,261:
stdout(`Number of Pythagorean Triples in a Perimeter of 100: `)
stdoutnl(num_pythagorean_triples(100))
</syntaxhighlight>
Output:
<pre>Number of Pythagorean Triples in a Perimeter of 100: 17
</pre>
=={{header|Liberty BASIC}}==
<syntaxhighlight lang="lb">
print time$()
Line 1,618 ⟶ 3,312:
print "End"
end
</syntaxhighlight>
<pre>
17:59:34
Line 1,635 ⟶ 3,329:
End
</pre>
=={{header|Mathematica}}/{{header|Wolfram Language}}==
=== Brute force ===
<syntaxhighlight lang="mathematica">pythag[n_]:=Block[{soln=Solve[{a^2+b^2==c^2,a+b+c<=n,0<a<b<c},{a,b,c},Integers]},{Length[soln],Count[GCD[a,b]/.soln,1]}]</syntaxhighlight>
Now prepare timings
<syntaxhighlight lang="mathematica">
pTiming[n_] := With[{comp = Timing@pythag@(10^n)},
{HoldForm[10^n], comp[[2, 1]], comp[[2, 2]], Round@comp[[1]]}];
{{"n", "Triples", "Primitives", "Timing(s)"}}~Join~(pTiming /@ Range@5) // Grid
</syntaxhighlight>
{{out}}
<pre>
n Triples Primitives Time(s)
10^1 0 0 3
10^2 17 7 5
10^3 325 70 7
10^4 4858 703 12
10^5 64741 7026 175
</pre>
=== Faster Primitives ===
The following uses generating formulae and is adapted from [[:https://mathematica.stackexchange.com/a/15904/2249]]
<syntaxhighlight lang="mathematica">primitivePythag[p_] := Join @@ Table[If[CoprimeQ[m, n], {2 m n, m^2 - n^2, m^2 + n^2}, ## &[]],{m, 2, Floor @ Sqrt @ p},{n, 1 + m ~Mod~ 2, m, 2}] // Select[Total[#] <= p &] // Length</syntaxhighlight>
Now prepare timings
<syntaxhighlight lang="mathematica">ppTiming[n_] := With[{comp = Timing@primitivePythag@(10^n)},{HoldForm[10^n], comp[[2]], Round@comp[[1]]}];
{{"n", "Primitives", "Timing(s)"}}~Join~(ppTiming /@ Range@9) // Grid</syntaxhighlight>
{{out}}
<pre>
n Primitives Time(s)
10^1 0 0
10^2 7 0
10^3 70 0
10^4 703 0
10^5 7026 0
10^6 70229 1
10^7 702309 10
10^8 7023027 111
10^9 70230484 1111
</pre>
The primitive counts are where the computational grunt-work lies (multiples under the limit can be readily computed) and this meets task challenge up to 10^8. "Faster Primitive" generates all triples and further, doesn't take advantage of built-in compilation which would need to be exploited if further order-of-magnitude improvements were required.
=={{header|MATLAB}} / {{header|Octave}}==
<
a = 1:N;
b = a(ones(N,1),:).^2;
Line 1,668 ⟶ 3,395:
printf('There are %i Pythagorean Triples and %i primitive triples with a perimeter smaller than %i.\n',...
sum(ix), length(p), N); </
Output:
<pre> There are 17 Pythagorean Triples and 7 primitive triples with a perimeter smaller than 100.</pre>
=={{header|Maxima}}==
<syntaxhighlight lang="maxima">
/* Function that returns a pythagorean triple from two parameters */
pythag(u,v):=[u^2-v^2,2*u*v,u^2+v^2]$
/* Predicate function to check for primitivity */
primitivep(lst):=if lreduce('gcd,lst)=1 then true$
/* Function that returns perimeter */
perim(lst):=apply("+",lst)$
/* Function to return a list of triples by parameter u */
/* Parameter v is controlled to be lesser or equal than u, and when equal are deleted */
param_pythag(n):=block(
create_list(lambda([x,y],pythag(x,y) and x#y)(i,j),i,1,n,j,1,i),
delete(false,%%))$
/* Test case */
/* With the function param_pythag as it is some non primitive triples are missing, but not the primitives */
sublist(param_pythag(6),lambda([x],primitivep(x) and perim(x)<=100));
/* The number of triples, primitive or not, can be recovered from the primitives */
block(
apply(append,makelist(%*i,i,1,8)),
sublist(%%,lambda([x],perim(x)<=100)));
</syntaxhighlight>
{{out}}
<pre>
[[3,4,5],[5,12,13],[15,8,17],[7,24,25],[21,20,29],[9,40,41],[35,12,37]]
[[3,4,5],[5,12,13],[15,8,17],[7,24,25],[21,20,29],[9,40,41],[35,12,37],[6,8,10],[10,24,26],[30,16,34],[9,12,15],[15,36,39],[12,16,20],[15,20,25],[18,24,30],[21,28,35],[24,32,40]]
</pre>
=={{header|Mercury}}==
From [[List comprehensions]]:
<syntaxhighlight lang="mercury">
:- module comprehension.
:- interface.
:- import_module io.
:- import_module int.
:- type triple ---> triple(int, int, int).
:- pred pythTrip(int::in,triple::out) is nondet.
:- pred main(io::di, io::uo) is det.
:- implementation.
:- import_module solutions.
pythTrip(Limit,triple(X,Y,Z)) :-
nondet_int_in_range(1,Limit,X),
nondet_int_in_range(X,Limit,Y),
nondet_int_in_range(Y,Limit,Z),
pow(Z,2) = pow(X,2) + pow(Y,2).
main(!IO) :-
solutions((pred(Triple::out) is nondet :- pythTrip(20,Triple)),Result),
write(Result,!IO).
</syntaxhighlight>
=={{header|Modula-3}}==
Note that this code only works on 64bit machines (where <tt>INTEGER</tt> is 64 bits). Modula-3 provides a <tt>LONGINT</tt> type, which is 64 bits on 32 bit systems, but there is a bug in the implementation apparently.
<
IMPORT IO, Fmt;
Line 1,705 ⟶ 3,495:
i := i * 10;
UNTIL i = 10000000;
END PyTriple64.</
Output:
Line 1,716 ⟶ 3,506:
</pre>
=={{header|
===Brute force method===
{{trans|Java}}
<syntaxhighlight lang="nanoquery">import math
// a function to check if three numbers are a valid triple
def is_triple(a, b, c)
if not (a < b) and (b < c)
return false
end
return (a^2 + b^2) = c^2
end
// a function to check if the numbers are coprime
def is_coprime(a, b, c)
global math
return (math.gcd(a, b)=1) && (math.gcd(a, c)=1) && (math.gcd(b, c)=1)
end
// the maximum perimeter to check
perimeter = 100
perimeter2 = int(perimeter / 2) - 1
perimeter3 = int(perimeter / 3) - 1
// loop though and look for pythagorean triples
ts = 0
ps = 0
for a in range(1, perimeter3)
for b in range(a + 1, perimeter2)
for c in range(b + 1, perimeter2)
if (a + b + c) <= perimeter
if is_triple(a,b,c)
ts += 1
print a + ", " + b + ", " + c
if is_coprime(a,b,c)
ps += 1
print " primitive"
end
println
end
end
end
end
end
print "Up to a perimeter of " + perimeter + ", there are " + ts
println " triples, of which " + ps + " are primitive."</syntaxhighlight>
{{out}}
<pre>3, 4, 5 primitive
5, 12, 13 primitive
6, 8, 10
7, 24, 25 primitive
8, 15, 17 primitive
9, 12, 15
9, 40, 41 primitive
10, 24, 26
12, 16, 20
12, 35, 37 primitive
15, 20, 25
15, 36, 39
16, 30, 34
18, 24, 30
20, 21, 29 primitive
21, 28, 35
24, 32, 40
Up to a perimeter of 100, there are 17 triples, of which 7 are primitive.</pre>
=={{header|Nim}}==
Compile with option <code>-d:release</code>. Without release option (i.e. in debug mode), the programs ends prematurely by reaching the recursion depth limit.
{{trans|C}}
<
[ 1, 2, 2, 2, 1, 2, 2, 2, 3],
[-1, 2, 2, -2, 1, 2, -2, 2, 3]]
Line 1,726 ⟶ 3,585:
maxPeri = 10
proc newTri(ins: array[0..2, int]) =
var p = ins[0] + ins[1] + ins[2]
if p > maxPeri: return
Line 1,742 ⟶ 3,601:
newTri([3, 4, 5])
echo "Up to ", maxPeri, ": ", total, " triples, ", prim, " primitives"
maxPeri *= 10</
Output:
<pre>Up to 10: 0 triples, 0 primitives
Line 1,754 ⟶ 3,613:
=={{header|OCaml}}==
<
let rec iter t =
let d = n - t*t in
Line 1,786 ⟶ 3,645:
Printf.printf "For perimeters up to %d there are %d total and %d primitive\n%!" i s p;;
List.iter show [ 100; 1000; 10000; 100000; 1000000; 10000000; 100000000 ]</
Output:
<pre>For perimeters up to 100 there are 17 total and 7 primitive
Line 1,795 ⟶ 3,654:
For perimeters up to 10000000 there are 9706567 total and 702309 primitive
For perimeters up to 100000000 there are 113236940 total and 7023027 primitive</pre>
=={{header|Ol}}==
<syntaxhighlight lang="scheme">
; triples generator based on Euclid's formula, creates lazy list
(define (euclid-formula max)
(let loop ((a 3) (b 4) (c 5) (tail #null))
(if (<= (+ a b c) max)
(cons (tuple a b c) (lambda ()
(let ((d (- b)) (z (- a)))
(loop (+ a d d c c) (+ a a d c c) (+ a a d d c c c) (lambda ()
(loop (+ a b b c c) (+ a a b c c) (+ a a b b c c c) (lambda ()
(loop (+ z b b c c) (+ z z b c c) (+ z z b b c c c) tail))))))))
tail)))
; let's do calculations
(define (calculate max)
(let loop ((p 0) (t 0) (ll (euclid-formula max)))
(cond
((null? ll)
(cons p t))
((function? ll)
(loop p t (ll)))
(else
(let ((triple (car ll)))
(loop (+ p 1) (+ t (div max (apply + triple)))
(cdr ll)))))))
; print values for 10..100000
(for-each (lambda (max)
(print max ": " (calculate max)))
(map (lambda (n) (expt 10 n)) (iota 6 1)))
</syntaxhighlight>
{{out}}
<pre>
10: (0 . 0)
100: (7 . 17)
1000: (70 . 325)
10000: (703 . 4858)
100000: (7026 . 64741)
1000000: (70229 . 808950)
</pre>
=={{header|PARI/GP}}==
This version is reasonably efficient and can handle inputs like a million quickly.
<
my(prim,total,P);
lim\=1;
Line 1,812 ⟶ 3,713:
[prim,total]
};
do(100)</
=={{header|Pascal}}==
<
var
Line 1,845 ⟶ 3,746:
maxPeri := maxPeri * 10;
end;
end.</
Output (on Core2Duo 2GHz laptop):
<pre>time ./PythagoreanTriples
Line 1,860 ⟶ 3,761:
=={{header|Perl}}==
<
my ($n, $m) = @_;
while($n){
Line 1,886 ⟶ 3,787:
tripel 10**$_ for 1..8;
</syntaxhighlight>
{{out}}
<pre>Max. perimeter: 10, Total: 0, Primitive: 0
Line 1,898 ⟶ 3,799:
</pre>
=={{header|Phix}}==
{{Trans|Pascal}}
<!--<syntaxhighlight lang="phix">(phixonline)-->
<span style="color: #008080;">with</span> <span style="color: #008080;">javascript_semantics</span>
<span style="color: #004080;">atom</span> <span style="color: #000000;">total</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">prim</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">maxPeri</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">10</span>
<span style="color: #008080;">procedure</span> <span style="color: #000000;">tri</span><span style="color: #0000FF;">(</span><span style="color: #004080;">atom</span> <span style="color: #000000;">s0</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">s1</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">s2</span><span style="color: #0000FF;">)</span>
<span style="color: #004080;">atom</span> <span style="color: #000000;">p</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">s0</span> <span style="color: #0000FF;">+</span> <span style="color: #000000;">s1</span> <span style="color: #0000FF;">+</span> <span style="color: #000000;">s2</span>
<span style="color: #008080;">if</span> <span style="color: #000000;">p</span><span style="color: #0000FF;"><=</span><span style="color: #000000;">maxPeri</span> <span style="color: #008080;">then</span>
<span style="color: #000000;">prim</span> <span style="color: #0000FF;">+=</span> <span style="color: #000000;">1</span>
<span style="color: #000000;">total</span> <span style="color: #0000FF;">+=</span> <span style="color: #7060A8;">floor</span><span style="color: #0000FF;">(</span><span style="color: #000000;">maxPeri</span><span style="color: #0000FF;">/</span><span style="color: #000000;">p</span><span style="color: #0000FF;">)</span>
<span style="color: #000000;">tri</span><span style="color: #0000FF;">(</span> <span style="color: #000000;">s0</span><span style="color: #0000FF;">+</span><span style="color: #000000;">2</span><span style="color: #0000FF;">*(-</span><span style="color: #000000;">s1</span><span style="color: #0000FF;">+</span><span style="color: #000000;">s2</span><span style="color: #0000FF;">),</span> <span style="color: #000000;">2</span><span style="color: #0000FF;">*(</span> <span style="color: #000000;">s0</span><span style="color: #0000FF;">+</span><span style="color: #000000;">s2</span><span style="color: #0000FF;">)-</span><span style="color: #000000;">s1</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">2</span><span style="color: #0000FF;">*(</span> <span style="color: #000000;">s0</span><span style="color: #0000FF;">-</span><span style="color: #000000;">s1</span><span style="color: #0000FF;">+</span><span style="color: #000000;">s2</span><span style="color: #0000FF;">)+</span><span style="color: #000000;">s2</span><span style="color: #0000FF;">);</span>
<span style="color: #000000;">tri</span><span style="color: #0000FF;">(</span> <span style="color: #000000;">s0</span><span style="color: #0000FF;">+</span><span style="color: #000000;">2</span><span style="color: #0000FF;">*(</span> <span style="color: #000000;">s1</span><span style="color: #0000FF;">+</span><span style="color: #000000;">s2</span><span style="color: #0000FF;">),</span> <span style="color: #000000;">2</span><span style="color: #0000FF;">*(</span> <span style="color: #000000;">s0</span><span style="color: #0000FF;">+</span><span style="color: #000000;">s2</span><span style="color: #0000FF;">)+</span><span style="color: #000000;">s1</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">2</span><span style="color: #0000FF;">*(</span> <span style="color: #000000;">s0</span><span style="color: #0000FF;">+</span><span style="color: #000000;">s1</span><span style="color: #0000FF;">+</span><span style="color: #000000;">s2</span><span style="color: #0000FF;">)+</span><span style="color: #000000;">s2</span><span style="color: #0000FF;">);</span>
<span style="color: #000000;">tri</span><span style="color: #0000FF;">(-</span><span style="color: #000000;">s0</span><span style="color: #0000FF;">+</span><span style="color: #000000;">2</span><span style="color: #0000FF;">*(</span> <span style="color: #000000;">s1</span><span style="color: #0000FF;">+</span><span style="color: #000000;">s2</span><span style="color: #0000FF;">),</span> <span style="color: #000000;">2</span><span style="color: #0000FF;">*(-</span><span style="color: #000000;">s0</span><span style="color: #0000FF;">+</span><span style="color: #000000;">s2</span><span style="color: #0000FF;">)+</span><span style="color: #000000;">s1</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">2</span><span style="color: #0000FF;">*(-</span><span style="color: #000000;">s0</span><span style="color: #0000FF;">+</span><span style="color: #000000;">s1</span><span style="color: #0000FF;">+</span><span style="color: #000000;">s2</span><span style="color: #0000FF;">)+</span><span style="color: #000000;">s2</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;">procedure</span>
<span style="color: #008080;">while</span> <span style="color: #000000;">maxPeri</span><span style="color: #0000FF;"><=</span><span style="color: #000000;">1e8</span> <span style="color: #008080;">do</span>
<span style="color: #000000;">prim</span> <span style="color: #0000FF;">:=</span> <span style="color: #000000;">0</span><span style="color: #0000FF;">;</span>
<span style="color: #000000;">total</span> <span style="color: #0000FF;">:=</span> <span style="color: #000000;">0</span><span style="color: #0000FF;">;</span>
<span style="color: #000000;">tri</span><span style="color: #0000FF;">(</span><span style="color: #000000;">3</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">4</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">5</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;">"Up to %d: %d triples, %d primitives.\n"</span><span style="color: #0000FF;">,</span> <span style="color: #0000FF;">{</span><span style="color: #000000;">maxPeri</span><span style="color: #0000FF;">,</span><span style="color: #000000;">total</span><span style="color: #0000FF;">,</span><span style="color: #000000;">prim</span><span style="color: #0000FF;">})</span>
<span style="color: #000000;">maxPeri</span> <span style="color: #0000FF;">*=</span> <span style="color: #000000;">10</span><span style="color: #0000FF;">;</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">while</span>
<!--</syntaxhighlight>-->
{{out}}
<pre>
Up to 10: 0 triples, 0 primitives.
Up to 100: 17 triples, 7 primitives.
Up to 1000: 325 triples, 70 primitives.
Up to 10000: 4858 triples, 703 primitives.
Up to 100000: 64741 triples, 7026 primitives.
Up to 1000000: 808950 triples, 70229 primitives.
Up to 10000000: 9706567 triples, 702309 primitives.
Up to 100000000: 113236940 triples, 7023027 primitives.
</pre>
=={{header|
<syntaxhighlight lang="php"><?php
function gcd($a, $b)
{
if ($a == 0)
return $b;
if ($b == 0)
return $a;
if($a == $b)
return $a;
if($a > $b)
return gcd($a-$b, $b);
return gcd($a, $b-$a);
}
$pytha = 0;
$prim = 0;
$max_p = 100;
for
$aa = $a**2;
for ($b = $a + 1; $b < $max_p/2; $b++) {
$bb = $b**2;
for ($c = $b + 1; $c < $max_p/2; $c++) {
$cc = $c**2;
if ($aa + $bb < $cc) break;
if ($a + $b + $c > $max_p) break;
if ($aa + $bb == $cc) {
$pytha++;
if (gcd($a, $b) == 1) $prim++;
}
}
}
}
echo 'Up to ' . $max_p . ', there are ' . $pytha . ' triples, of which ' . $prim . ' are primitive.';</syntaxhighlight>
=={{header|Picat}}==
===Prolog-style===
{{trans|Prolog}}
<syntaxhighlight lang="picat">main :-
garbage_collect(300_000_000),
Data = [100, 1_000, 10_000, 100_000, 1_000_000, 10_000_000, 100_000_000],
member(Max, Data),
count_triples(Max, Total, Prim),
printf("upto %d, there are %d Pythagorean triples (%d primitive.)%n", Max, Total, Prim),
fail,
count_triples(Max, Total, Prims) :-
Ps = findall(S, (triple(Max, A, B,
Prims = Ps.len,
Total = sum([Max div P : P in Ps]).
% - between_by/4
between_by(A, B, N, K) :-
C = (B - A) div N,
between(0, C, J),
K = N*J + A.
% - Pythagorean triple generator
triple(P, A, B, C) :-
Max = floor(sqrt(P/2)) - 1,
between(0, Max, M),
between_by(Start, Pm, 2, N),
gcd(M, N) == 1,
X = M*M - N*N,
Y = 2*M*N,
C = M*M + N*N,
order2(X, Y, A, B),
(A + B + C) =< P.
order2(A, B, A, B) :- A < B, !.
order2(A, B, B, A).</syntaxhighlight>
{{out}}
<pre>upto 100, there are 17 Pythagorean triples (7 primitive.)
upto 1000, there are 325 Pythagorean triples (70 primitive.)
upto 10000, there are 4857 Pythagorean triples (702 primitive.)
upto 100000, there are 64741 Pythagorean triples (7026 primitive.)
upto 1000000, there are 808950 Pythagorean triples (70229 primitive.)
upto 10000000, there are 9706567 Pythagorean triples (702309 primitive.)
upto 100000000, there are 113236940 Pythagorean triples (7023027 primitive.)
CPU time 4.612 seconds.</pre>
===Another approach===
{{trans|Go}}
Picat doesn't have global variables, so all parameters are placed in the call to <code>newTri/6</code>. This is slightly faster than the Prolog port.
<syntaxhighlight lang="picat">main =>
foreach(MaxPeri in [10**I : I in 2..8])
[Total, Prim] = newTri(MaxPeri,0,0,3,4,5),
printf("Up to %d: %d triples, %d primitives\n", MaxPeri, Total, Prim)
end.
newTri(MaxPeri,Prim,Total,S0, S1, S2) = [PrimRet,TotalRet] =>
P =
if P <= MaxPeri then
Prim2 = Prim + 1,
Total2 = Total + MaxPeri div P,
[Prim3,Total3] = newTri(MaxPeri,Prim2,Total2, +1*S0-2*S1+2*S2, +2*S0-1*S1+2*S2, +2*S0-2*S1+3*S2),
[Prim4,Total4] = newTri(MaxPeri,Prim3,Total3, +1*S0+2*S1+2*S2, +2*S0+1*S1+2*S2, +2*S0+2*S1+3*S2),
[Prim5,Total5] = newTri(MaxPeri,Prim4,Total4, -1*S0+2*S1+2*S2, -2*S0+1*S1+2*S2, -2*S0+2*S1+3*S2),
PrimRet = Prim5,
TotalRet = Total5
else
PrimRet = Prim,
TotalRet = Total
end.</syntaxhighlight>
{{out}}
<pre>Up to 100: 7 triples, 17 primitives
Up to 1000: 70 triples, 325 primitives
Up to 10000: 703 triples, 4858 primitives
Up to 100000: 7026 triples, 64741 primitives
Up to 1000000: 70229 triples, 808950 primitives
Up to 10000000: 702309 triples, 9706567 primitives
Up to 100000000: 7023027 triples, 113236940 primitives
CPU time 4.373 seconds.</pre>
===With "global variables"===
Actually, Picat has some support for global variables, using the global available map (<code>get_global_map</code>).
<syntaxhighlight lang="picat">main =>
foreach(MaxPeri in [10**I : I in 2..8])
Map = get_global_map(),
Map.put(max_peri,MaxPeri),
Map.put(prim,0),
Map.put(total,0),
newTri3(3,4,5),
printf("Up to %d: %d triples, %d primitives\n", MaxPeri, Map.get(total), Map.get(prim))
end.
newTri2(S0, S1, S2) =>
P = S0 + S1 + S2,
Map = get_global_map(),
if P <= Map.get(max_peri) then
Map.put(prim, Map.get(prim)+1),
Map.put(total,Map.get(total) + Map.get(max_peri) div P),
newTri2(+1*S0-2*S1+2*S2, +2*S0-1*S1+2*S2, +2*S0-2*S1+3*S2),
newTri2(+1*S0+2*S1+2*S2, +2*S0+1*S1+2*S2, +2*S0+2*S1+3*S2),
newTri2(-1*S0+2*S1+2*S2, -2*S0+1*S1+2*S2, -2*S0+2*S1+3*S2)
end.</syntaxhighlight>
This version is - however - slower: 7.401s.
=={{header|PicoLisp}}==
{{trans|C}}
<
(let (Total 0 Prim 0 In (3 4 5))
(recur (In)
Line 2,030 ⟶ 4,004:
(recurse
(mapcar '((U) (sum * U In)) Row) ) ) ) ) )
(prinl "Up to " Max ": " Total " triples, " Prim " primitives.") ) )</
Output:
<pre>Up to 10: 0 triples, 0 primitives.
Line 2,043 ⟶ 4,017:
=={{header|PL/I}}==
Version 1
<
/*********************************************************************
* REXX pgm counts number of Pythagorean triples
Line 2,132 ⟶ 4,106:
End;
End;</
Version 2
<syntaxhighlight lang="pl/i">
pythagorean: procedure options (main, reorder); /* 23 January 2014 */
declare (a, b, c) fixed (3);
Line 2,166 ⟶ 4,140:
end GCD;
end pythagorean;</
Output:
<pre>
Line 2,174 ⟶ 4,148:
<pre>
TRIPLES= 325 PRIMITIVES= 70;
</pre>
=={{header|PowerShell}}==
<syntaxhighlight lang="powershell">
function triples($p) {
if($p -gt 4) {
# ai + bi + ci = pi <= p
# ai < bi < ci --> 3ai < pi <= p and ai + 2bi < pi <= p
$pa = [Math]::Floor($p/3)
1..$pa | foreach {
$ai = $_
$pb = [Math]::Floor(($p-$ai)/2)
($ai+1)..$pb | foreach {
$bi = $_
$pc = $p-$ai-$bi
($bi+1)..$pc | where {
$ci = $_
$pi = $ai + $bi + $ci
$ci*$ci -eq $ai*$ai + $bi*$bi
} |
foreach {
[pscustomobject]@{
a = "$ai"
b = "$bi"
c = "$ci"
p = "$pi"
}
}
}
}
}
else {
Write-Error "$p is not greater than 4"
}
}
function gcd ($a, $b) {
function pgcd ($n, $m) {
if($n -le $m) {
if($n -eq 0) {$m}
else{pgcd $n ($m%$n)}
}
else {pgcd $m $n}
}
$n = [Math]::Abs($a)
$m = [Math]::Abs($b)
(pgcd $n $m)
}
$triples = (triples 100)
$coprime = $triples |
where {((gcd $_.a $_.b) -eq 1) -and ((gcd $_.a $_.c) -eq 1) -and ((gcd $_.b $_.c) -eq 1)}
"There are $(($triples).Count) Pythagorean triples with perimeter no larger than 100
and $(($coprime).Count) of them are coprime."
</syntaxhighlight>
<b>Output:</b>
<pre>
There are 17 Pythagorean triples with perimeter no larger than 100 and 7 of them are coprime.
</pre>
=={{header|Prolog}}==
<syntaxhighlight lang="prolog">
show :-
Data = [100, 1_000, 10_000, 100_000, 1_000_000, 10_000_000, 100_000_000],
forall(
member(Max, Data),
(count_triples(Max, Total, Prim),
format("upto ~D, there are ~D Pythagorean triples (~D primitive.)~n", [Max, Total, Prim]))).
div(A, B, C) :- C is A div B.
count_triples(Max, Total, Prims) :-
findall(S, (triple(Max, A, B, C), S is A + B + C), Ps),
length(Ps, Prims),
maplist(div(Max), Ps, Counts), sumlist(Counts, Total).
% - between_by/4
between_by(A, B, N, K) :-
C is (B - A) div N,
between(0, C, J),
K is N*J + A.
% - Pythagorean triple generator
triple(P, A, B, C) :-
Max is floor(sqrt(P/2)) - 1,
between(0, Max, M),
Start is (M /\ 1) + 1, succ(Pm, M),
between_by(Start, Pm, 2, N),
gcd(M, N) =:= 1,
X is M*M - N*N,
Y is 2*M*N,
C is M*M + N*N,
order2(X, Y, A, B),
(A + B + C) =< P.
order2(A, B, A, B) :- A < B, !.
order2(A, B, B, A).
</syntaxhighlight>
{{Out}}
<pre>
?- show.
upto 100, there are 17 Pythagorean triples (7 primitive.)
upto 1,000, there are 325 Pythagorean triples (70 primitive.)
upto 10,000, there are 4,857 Pythagorean triples (702 primitive.)
upto 100,000, there are 64,741 Pythagorean triples (7,026 primitive.)
upto 1,000,000, there are 808,950 Pythagorean triples (70,229 primitive.)
upto 10,000,000, there are 9,706,567 Pythagorean triples (702,309 primitive.)
upto 100,000,000, there are 113,236,940 Pythagorean triples (7,023,027 primitive.)
true.
</pre>
Line 2,185 ⟶ 4,271:
<math>n \le 10,000</math><br /><br />
<
Procedure.i ConsoleWrite(t.s) ; compile using /CONSOLE option
Line 2,254 ⟶ 4,340:
ConsoleWrite("")
et=ElapsedMilliseconds()-st:ConsoleWrite("Elapsed time = "+str(et)+" milliseconds")
</syntaxhighlight>
Line 2,278 ⟶ 4,364:
=={{header|Python}}==
Two methods, the second of which is much faster
<
Line 2,334 ⟶ 4,420:
for maxperimeter in range(mn, mx+1, mn):
printit(maxperimeter, algo)
</syntaxhighlight>
;Output:
Line 2,394 ⟶ 4,480:
Up to a perimeter of 19500 there are 10388 triples, of which 1373 are primitive
Up to a perimeter of 20000 there are 10689 triples, of which 1408 are primitive</pre>
Barebone minimum for this task:<
setrecursionlimit(2000) # 2000 ought to be big enough for everybody
Line 2,407 ⟶ 4,493:
for peri in [10 ** e for e in range(1, 8)]:
print peri, triples(peri)</
100 (7, 17)
1000 (70, 325)
Line 2,413 ⟶ 4,499:
100000 (7026, 64741)
1000000 (70229, 808950)
10000000 (702309, 9706567)</
=={{header|Quackery}}==
Based on the method shown in the Mathologer video "Fibonacci = Pythagoras: Help save a stunning discovery from oblivion!". https://youtu.be/94mV7Fmbx88
From the first solution "(1, 1, 2, 3)" (in quackery notation <code>[ 1 1 2 3 ]</code>) use functions <code>f1</code>, <code>f2</code>, and <code>f3</code> for recursive depth-first tree traversal. Terminating condition is that perimeter > limit (i.e. 100). At each node add one to the count of primitive pythagorean triples, and n to the count of all pythagorean triples, where n = floor(limit/perimeter).
<pre> f1 corresponds to the right hand branch of the tree in the mathologer video
f2 corresponds to the left hand branch...
f3 corresponds to the central branch...
a b c d A B C D A B C D
f1 [ 1 1 2 3 ] --> [ 1 2 3 5 ] [ a c A+B B+C ]
f2 [ 1 1 2 3 ] --> [ 3 1 4 5 ] [ d b A+B B+C ]
f3 [ 1 1 2 3 ] --> [ 3 2 5 7 ] [ d c A+B B+C ]
a b c d
[ 1 1 2 3 ] --> [ 3 4 5 ] [ a*d 2b*c (b*d)+(a*c) ] pythagorean triple
a*d + 2b*c + b*d + a*c perimeter
= (a*(c+d))+(b*(2c+d))</pre>
The recursive part of the word <code>task</code> is
<pre> [ dup perimeter
limit share over < iff 2drop done ( end on terminating condition )
1 primitives tally
limit share swap / triples tally
dup f1 recurse
dup f2 recurse
f3 again ]</pre>
Note that <code>f3 again</code> is equivalent to <code>f3 recurse</code> but a smidgeon faster by optimising tail-end recursion.
The ancillary stacks <code>limit</code>,<code>primitives</code>, and <code>triples</code> can be integer variables in languages that have variables.
<syntaxhighlight lang="Quackery"> [ dup 0 peek
swap 2 peek
2dup + 2dup +
join join join ] is f1 ( [ --> [ )
[ dup 3 peek
swap 1 peek
2dup + 2dup +
join join join ] is f2 ( [ --> [ )
[ dup 3 peek
swap 2 peek
2dup + 2dup +
join join join ] is f3 ( [ --> [ )
[ do over + tuck + rot * unrot * + ] is perimeter ( [ --> n )
[ stack ] is limit ( --> s )
[ stack ] is primitives ( --> s )
[ stack ] is triples ( --> s )
[ limit put
0 primitives put
0 triples put
' [ 1 1 2 3 ]
[ dup perimeter
limit share over < iff 2drop done
1 primitives tally
limit share swap / triples tally
dup f1 recurse
dup f2 recurse
f3 again ]
say "Pythagorean triples, perimeter < "
limit take echo
say ": "
triples take echo
say ", of which "
primitives take echo
say " are primitive." cr ] is task ( n --> )
7 times [ 10 i^ 2 + ** task ]</syntaxhighlight>
{{out}}
<pre>Pythagorean triples, perimeter < 100: 17, of which 7 are primitive.
Pythagorean triples, perimeter < 1000: 325, of which 70 are primitive.
Pythagorean triples, perimeter < 10000: 4858, of which 703 are primitive.
Pythagorean triples, perimeter < 100000: 64741, of which 7026 are primitive.
Pythagorean triples, perimeter < 1000000: 808950, of which 70229 are primitive.
Pythagorean triples, perimeter < 10000000: 9706567, of which 702309 are primitive.
Pythagorean triples, perimeter < 100000000: 113236940, of which 7023027 are primitive.
</pre>
=={{header|Racket}}==
<
#| Euclid's enumeration formula and counting is fast enough for extra credit.
Line 2,460 ⟶ 4,633:
113236940, 7023027.
cpu time: 11976 real time: 12215 gc time: 2381
|#</
=={{header|
(formerly Perl 6)
{{works with|Rakudo|2018.09}}
Here is a straight-forward, naïve brute force implementation:
<syntaxhighlight lang="raku" line>constant limit = 100;
for [X] [^limit] xx 3 -> (\a, \b, \c) {
say [a, b, c] if a < b < c and a + b + c <= limit and a*b + b*b == c*c
}</syntaxhighlight>
{{out}}
<pre style="height:25ex">[3 4 5]
[5 12 13]
[6 8 10]
[7 24 25]
[8 15 17]
[9 12 15]
[9 40 41]
[10 24 26]
[12 16 20]
[12 35 37]
[15 20 25]
[15 36 39]
[16 30 34]
[18 24 30]
[20 21 29]
[21 28 35]
[24 32 40]
</pre>
Here is a slightly less naive brute force implementation, but still not practical for large perimeter limits.
<syntaxhighlight lang="raku" line>my $limit = 10000;
my atomicint $i = 0;
my @triples[$limit/2];
(3 .. $limit/2).race.map: -> $c {
for 1 .. $c -> $a {
my $b = ($c * $c - $a * $a).sqrt;
last if $c + $a + $b > $limit;
last if $a > $b;
@triples[$i⚛++] = ([gcd] $c, $a, $b) > 1 ?? 0 !! 1 if $b == $b.Int;
}
}
say my $result = "There are {+@triples.grep:{$_ !eqv Any}} Pythagorean Triples with a perimeter <= $limit,"
~"\nof which {[+] @triples.grep: so *} are primitive.";</syntaxhighlight>
{{out}}
<pre>There are 4858 Pythagorean Triples with a perimeter <= 10000,
of which 703 are primitive.</pre>
Here's a much faster version. Hint, "oyako" is Japanese for "parent/child". <tt>:-)</tt>
<syntaxhighlight lang="raku" line>sub triples($limit) {
my $primitive = 0;
my $civilized = 0;
sub oyako($a, $b, $c) {
my $perim = $a + $b + $c;
return if $perim > $limit;
++$primitive; $civilized += $limit div $perim;
oyako( $a - 2*$b + 2*$c, 2*$a - $b + 2*$c, 2*$a - 2*$b + 3*$c);
oyako( $a + 2*$b + 2*$c, 2*$a + $b + 2*$c, 2*$a + 2*$b + 3*$c);
oyako(-$a + 2*$b + 2*$c, -2*$a + $b + 2*$c, -2*$a + 2*$b + 3*$c);
}
oyako(3,4,5);
"$limit => ($primitive $civilized)";
}
for 10,100,1000 ... * -> $limit {
say triples $limit;
}</syntaxhighlight>
Output:
<pre>10 => (0 0)
100 => (7 17)
1000 => (70 325)
10000 => (703 4858)
100000 => (7026 64741)
1000000 => (70229 808950)
10000000 => (702309 9706567)
100000000 => (7023027 113236940)
1000000000 => (70230484 1294080089)
^C</pre>
The geometric sequence of limits will continue on forever, so eventually when you get tired of waiting (about a billion on my computer), you can just stop it. Another efficiency trick of note: we avoid passing the limit in as a parameter to the inner helper routine, but instead supply the limit via the lexical scope. Likewise, the accumulators are referenced lexically, so only the triples themselves need to be passed downward, and nothing needs to be returned.
Here is an alternate version that avoids naming any scalars that can be handled by vector processing instead. Using vectorized ops allows a bit more potential for parallelization in theory, but techniques like the use of complex numbers to add two numbers in parallel, and the use of <tt>gather</tt>/<tt>take</tt> generate so much overhead that this version runs 70-fold slower than the previous one.
<syntaxhighlight lang="raku" line>constant @coeff = [[+1, -2, +2], [+2, -1, +2], [+2, -2, +3]],
sub triples($limit) {
sub oyako(@trippy) {
my $perim = [+] @trippy;
return if $perim > $limit;
take (1 + ($limit div $perim)i);
for @coeff -> @nine {
oyako (map -> @three { [+] @three »*« @trippy }, @nine);
}
return;
}
my $complex = 0i + [+] gather oyako([3,4,5]);
"$limit => ({$complex.re, $complex.im})";
}
for 10, 100, 1000, 10000 -> $limit {
say triples $limit;
}</syntaxhighlight>
{{out}}
<pre>10 => (0 0)
100 => (7 17)
1000 => (70 325)
10000 => (703 4858)</pre>
=={{header|REXX}}==
===using GCD for determinacy===
<syntaxhighlight lang="rexx">/*REXX program counts the number of Pythagorean triples that exist given a maximum */
/*──────────────────── perimeter of N, and also counts how many of them are primitives.*/
parse arg N . /*obtain optional argument from the CL.*/
if N=='' | N=="," then N= 100 /*Not specified? Then use the default.*/
do j=1 for N; @.j= j*j; end /*pre-compute some squares. */
N66= N * 2%3 /*calculate 2/3 of N (for a+b). */
T= 0; P= 0 /*set the number of Triples, Primitives*/
do a=3 to N%3 /*limit side to 1/3 of the perimeter.*/
do b= a+1 /*the triangle can't be isosceles. */
ab= a + b /*compute a partial perimeter (2 sides)*/
if ab>=N66 then iterate a /*is a+b≥66% perimeter? Try different A*/
aabb= @.a + @.b /*compute the sum of a²+b² (shortcut)*/
do c=b+1 /*compute the value of the third side. */
if ab+c > N then iterate a /*is a+b+c>perimeter ? Try different A.*/
if @.c >aabb then iterate b /*is c² > a²+b² ? Try " B.*/
if @.c\==aabb then iterate /*is c² ¬= a²+b² ? Try " C.*/
T= T + 1 /*eureka. We found a Pythagorean triple*/
P= P + (gcd(a, b)==1) /*is this triple a primitive triple? */
end /*c*/
end /*b*/
end /*a*/
_= left('', 7) /*for padding the output with 7 blanks.*/
say 'max perimeter =' N _ "Pythagorean triples =" T _ 'primitives =' P
exit /*stick a fork in it, we're all done. */
/*──────────────────────────────────────────────────────────────────────────────────────*/
gcd: procedure; parse arg x,y; do until y==0; parse value x//y y with y x; end; return x</syntaxhighlight>
{{out|output|text= when using the default input of: <tt> 100 </tt>}}
<pre>
max perimeter = 100 Pythagorean triples = 17 primitives = 7
</pre>
<pre>
max perimeter = 1000 Pythagorean triples = 325 primitives = 70
</pre>
===using single evenness for determinacy===
This REXX version takes advantage that primitive Pythagorean triples must have one and only one ''even'' number.
This REXX version is about '''10%''' faster than the 1<sup>st</sup> REXX version.
Non-primitive Pythagorean triples are generated after a primitive triple is found.
<syntaxhighlight lang="rexx">/*REXX program counts the number of Pythagorean triples that exist given a maximum */
/*──────────────────── perimeter of N, and also counts how many of them are primitives.*/
parse arg N .
if N=='' | N=="," then N= 100 /*Not specified? Then use the default.*/
@.= 0; do j=1 for N; @.j= j*j; end /*pre-compute some squares. */
N66= N * 2%3 /*calculate 2/3 of N (for a+b). */
P= 0; T= 0; do a=3 to N%3 /*limit side to 1/3 of the perimeter.*/
aEven= a//2==0 /*set variable to 1 if A is even. */
do b=a+1 by 1+aEven; ab= a + b /*the triangle can't be isosceles. */
if ab>=N66 then iterate a /*is a+b≥66% perimeter? Try different A*/
aabb= @.a + @.b /*compute the sum of a²+b² (shortcut)*/
do c=b + 1 /*compute the value of the third side. */
if aEven then if c//2==0 then iterate /*both A&C even? Skip it*/
if ab+c>n then iterate a /*a+b+c > perimeter? Try different A. */
if @.c > aabb then iterate b /*is c² > a²+b² ? " " B. */
if @.c\==aabb then iterate /*is c² ¬= a²+b² ? " " C. */
if @.a.b.c then iterate /*Is this a duplicate? Then try again.*/
T= T + 1 /*Eureka! We found a Pythagorean triple*/
P= P + 1 /*count this also as a primitive triple*/
do m=2 while a*m+b*m+c*m<=N /*generate non-primitives Pythagoreans.*/
T= T + 1 /*Eureka! We found a Pythagorean triple*/
am= a*m; bm= b*m; cm= c*m /*create some short-cut variable names.*/
@.am.bm.cm= 1 /*mark Pythagorean triangle as a triple*/
end /*m*/
end /*c*/
end /*b*/
end /*a*/ /*stick a fork in it, we're all done. */
_= left('', 7) /*for padding the output with 7 blanks.*/
say 'max perimeter =' N _ "Pythagorean triples =" T _ 'primitives =' P</syntaxhighlight>
{{out|output|text= is identical to the 1<sup>st</sup> REXX version.}}<br><br>
{{out|output|text= when using the input of: <tt> 10000 </tt>}}
<pre>
max perimeter = 10000
</pre>
=={{header|Ring}}==
<syntaxhighlight lang="ring">
size = 100
sum = 0
prime = 0
for i = 1 to size
for j =
if pow(i,2) + pow(j,2) = pow(k,2) and (i+j+k) < 101
if gcd(i,j) = 1 prime += 1 ok
sum += 1
see "" + i + " " + j + " " + k + nl ok
next
next
next
see "Total : " + sum + nl
see "Primitives : " + prime + nl
func gcd gcd, b
while b
gcd = b
b = c % b
end
return gcd
</syntaxhighlight>
Output:
<pre>
3 4 5
5 12 13
6 8 10
7 24 25
8 15 17
9 12 15
9 40 41
10 24 26
12 16 20
12 35 37
15 20 25
15 36 39
16 30 34
18 24 30
20 21 29
21 28 35
24 32 40
Total : 17
Primitives : 7
</pre>
=={{header|Ruby}}==
{{trans|Java}}
<
def initialize(limit)
@limit = limit
@total = 0
@
generate_triples(3, 4, 5)
end
attr_reader :total, :
private
Line 2,560 ⟶ 4,889:
return if perim > @limit
@
@total += @limit / perim
Line 2,572 ⟶ 4,901:
while perim <= 100_000_000
c = PythagoranTriplesCounter.new perim
p [perim, c.total, c.
perim *= 10
end</
output
Line 2,585 ⟶ 4,914:
[10000000, 9706567, 702309]
[100000000, 113236940, 7023027]</pre>
=={{header|Rust}}==
<syntaxhighlight lang="rust">use std::thread;
fn f1 (a : u64, b : u64, c : u64, d : u64) -> u64 {
let mut primitive_count = 0;
for triangle in [[a - 2*b + 2*c, 2*a - b + 2*c, 2*a - 2*b + 3*c],
[a + 2*b + 2*c, 2*a + b + 2*c, 2*a + 2*b + 3*c],
[2*b + 2*c - a, b + 2*c - 2*a, 2*b + 3*c - 2*a]] .iter() {
let l = triangle[0] + triangle[1] + triangle[2];
if l > d { continue; }
primitive_count += 1 + f1(triangle[0], triangle[1], triangle[2], d);
}
primitive_count
}
fn f2 (a : u64, b : u64, c : u64, d : u64) -> u64 {
let mut triplet_count = 0;
for triangle in [[a - 2*b + 2*c, 2*a - b + 2*c, 2*a - 2*b + 3*c],
[a + 2*b + 2*c, 2*a + b + 2*c, 2*a + 2*b + 3*c],
[2*b + 2*c - a, b + 2*c - 2*a, 2*b + 3*c - 2*a]] .iter() {
let l = triangle[0] + triangle[1] + triangle[2];
if l > d { continue; }
triplet_count += (d/l) + f2(triangle[0], triangle[1], triangle[2], d);
}
triplet_count
}
fn main () {
let new_th_1 = thread::Builder::new().stack_size(32 * 1024 * 1024).spawn (move || {
let mut i = 100;
while i <= 100_000_000_000 {
println!(" Primitive triples below {} : {}", i, f1(3, 4, 5, i) + 1);
i *= 10;
}
}).unwrap();
let new_th_2 =thread::Builder::new().stack_size(32 * 1024 * 1024).spawn (move || {
let mut i = 100;
while i <= 100_000_000_000 {
println!(" Triples below {} : {}", i, f2(3, 4, 5, i) + i/12);
i *= 10;
}
}).unwrap();
new_th_1.join().unwrap();
new_th_2.join().unwrap();
}</syntaxhighlight>
{{out}}
<pre> Primitive triples below 100 : 7
Triples below 100 : 17
Primitive triples below 1000 : 70
Triples below 1000 : 325
Primitive triples below 10000 : 703
Triples below 10000 : 4858
Primitive triples below 100000 : 7026
Triples below 100000 : 64741
Primitive triples below 1000000 : 70229
Triples below 1000000 : 808950
Primitive triples below 10000000 : 702309
Triples below 10000000 : 9706567
Primitive triples below 100000000 : 7023027
Triples below 100000000 : 113236940
Primitive triples below 1000000000 : 70230484
Triples below 1000000000 : 1294080089
Primitive triples below 10000000000 : 702304875
Triples below 10000000000 : 14557915466
Primitive triples below 100000000000 : 7023049293
Triples below 100000000000 : 161750315680
real 2m22.676s
user 3m39.239s
sys 0m0.024s</pre>
=={{header|Scala}}==
{{Out}}Best seen running in your browser either by [https://scalafiddle.io/sf/CAz60TW/0 ScalaFiddle (ES aka JavaScript, non JVM)] or [https://scastie.scala-lang.org/soOLJ673Q82l78OCgIx4oQ Scastie (remote JVM)].
<syntaxhighlight lang="scala">object PythagoreanTriples extends App {
println(" Limit Primatives All")
for {e <- 2 to 7
limit = math.pow(10, e).longValue()
} {
var primCount, tripCount = 0
def parChild(a: BigInt, b: BigInt, c: BigInt): Unit = {
val perim = a + b + c
val (a2, b2, c2, c3) = (2 * a, 2 * b, 2 * c, 3 * c)
if (limit >= perim) {
primCount += 1
tripCount += (limit / perim).toInt
parChild(a - b2 + c2, a2 - b + c2, a2 - b2 + c3)
parChild(a + b2 + c2, a2 + b + c2, a2 + b2 + c3)
parChild(-a + b2 + c2, -a2 + b + c2, -a2 + b2 + c3)
}
}
parChild(BigInt(3), BigInt(4), BigInt(5))
println(f"a + b + c <= ${limit.toFloat}%3.1e $primCount%9d $tripCount%12d")
}
}</syntaxhighlight>
=={{header|Scheme}}==
{{works with|Gauche Scheme}}
<syntaxhighlight lang="scheme">(use srfi-42)
(define (py perim)
(define prim 0)
(values
(sum-ec
(: c perim) (: b c) (: a b)
(if (and (<= (+ a b c) perim)
(= (square c) (+ (square b) (square a)))))
(begin (when (= 1 (gcd a b)) (inc! prim)))
1)
prim))</syntaxhighlight>
<b>Testing:</b>
<pre>
gosh> (py 100)
17
7
</pre>
=={{header|Scratch}}==
Scratch is a visual programming language. Click the link, then "see inside" to see the code.
https://scratch.mit.edu/projects/79066598/
Output: 17 Pythagorean triples with a perimeter less than 100, 7 of which are primitive.
=={{header|Seed7}}==
Line 2,590 ⟶ 5,048:
The example below uses [http://seed7.sourceforge.net/libraries/bigint.htm bigInteger] numbers:
<
include "bigint.s7i";
Line 2,620 ⟶ 5,078:
max_peri *:= 10_;
end while;
end func;</
Output:
Line 2,632 ⟶ 5,090:
Up to 10000000: 9706567 triples, 702309 primitives.
Up to 100000000: 113236940 triples, 7023027 primitives.
</pre>
=={{header|Sidef}}==
{{trans|Raku}}
<syntaxhighlight lang="ruby">func triples(limit) {
var primitive = 0
var civilized = 0
func oyako(a, b, c) {
(var perim = a+b+c) > limit || (
primitive++
civilized += int(limit / perim)
oyako( a - 2*b + 2*c, 2*a - b + 2*c, 2*a - 2*b + 3*c)
oyako( a + 2*b + 2*c, 2*a + b + 2*c, 2*a + 2*b + 3*c)
oyako(-a + 2*b + 2*c, -2*a + b + 2*c, -2*a + 2*b + 3*c)
)
}
oyako(3,4,5)
"#{limit} => (#{primitive} #{civilized})"
}
for n (1..Inf) {
say triples(10**n)
}</syntaxhighlight>
{{out}}
<pre>
10 => (0 0)
100 => (7 17)
1000 => (70 325)
10000 => (703 4858)
100000 => (7026 64741)
1000000 => (70229 808950)
^C
</pre>
=={{header|Swift}}==
{{trans|Pascal}}
<syntaxhighlight lang="swift">var total = 0
var prim = 0
var maxPeri = 100
func newTri(s0:Int, _ s1:Int, _ s2: Int) -> () {
let p = s0 + s1 + s2
if p <= maxPeri {
prim += 1
total += maxPeri / p
newTri( s0 + 2*(-s1+s2), 2*( s0+s2) - s1, 2*( s0-s1+s2) + s2)
newTri( s0 + 2*( s1+s2), 2*( s0+s2) + s1, 2*( s0+s1+s2) + s2)
newTri(-s0 + 2*( s1+s2), 2*(-s0+s2) + s1, 2*(-s0+s1+s2) + s2)
}
}
while maxPeri <= 100_000_000 {
prim = 0
total = 0
newTri(3, 4, 5)
print("Up to \(maxPeri) : \(total) triples \( prim) primitives.")
maxPeri *= 10
}</syntaxhighlight>
{{out}}
<pre>
Up to 100 : 17 triples 7 primitives.
Up to 1000 : 325 triples 70 primitives.
Up to 10000 : 4858 triples 703 primitives.
Up to 100000 : 64741 triples 7026 primitives.
Up to 1000000 : 808950 triples 70229 primitives.
Up to 10000000 : 9706567 triples 702309 primitives.
Up to 100000000 : 113236940 triples 7023027 primitives.
</pre>
Line 2,637 ⟶ 5,167:
Using the efficient method based off the Wikipedia article:
<!--There's no technical reason to limit the code to just these values, but generation does get progressively slower with larger maximum perimiters. 10M is about as much as I have patience for; I'm generally impatient! -->
<
lappend q 3 4 5
set idx [set count [set prim 0]]
Line 2,661 ⟶ 5,191:
lassign [countPythagoreanTriples $i] count primitive
puts "perimeter limit $i => $count triples, $primitive primitive"
}</
Output:
<pre>
Line 2,671 ⟶ 5,201:
perimeter limit 1000000 => 808950 triples, 70229 primitive
perimeter limit 10000000 => 9706567 triples, 702309 primitive
</pre>
=={{header|VBA}}==
{{trans|Pascal}}<syntaxhighlight lang="vb">Dim total As Variant, prim As Variant, maxPeri As Variant
Private Sub newTri(s0 As Variant, s1 As Variant, s2 As Variant)
Dim p As Variant
p = CDec(s0) + CDec(s1) + CDec(s2)
If p <= maxPeri Then
prim = prim + 1
total = total + maxPeri \ p
newTri s0 + 2 * (-s1 + s2), 2 * (s0 + s2) - s1, 2 * (s0 - s1 + s2) + s2
newTri s0 + 2 * (s1 + s2), 2 * (s0 + s2) + s1, 2 * (s0 + s1 + s2) + s2
newTri -s0 + 2 * (s1 + s2), 2 * (-s0 + s2) + s1, 2 * (-s0 + s1 + s2) + s2
End If
End Sub
Public Sub Program_PythagoreanTriples()
maxPeri = CDec(100)
Do While maxPeri <= 10000000#
prim = CDec(0)
total = CDec(0)
newTri 3, 4, 5
Debug.Print "Up to "; maxPeri; ": "; total; " triples, "; prim; " primitives."
maxPeri = maxPeri * 10
Loop
End Sub</syntaxhighlight>{{out}}
<pre>Up to 100 : 17 triples, 7 primitives.
Up to 1000 : 325 triples, 70 primitives.
Up to 10000 : 4858 triples, 703 primitives.
Up to 100000 : 64741 triples, 7026 primitives.
Up to 1000000 : 808950 triples, 70229 primitives.
Up to 10000000 : 9706567 triples, 702309 primitives.
</pre>
=={{header|VBScript}}==
{{trans|Perl}}
<syntaxhighlight lang="vb">
For i=1 To 8
WScript.StdOut.WriteLine triples(10^i)
Next
Function triples(pmax)
prim=0 : count=0 : nmax=Sqr(pmax)/2 : n=1
Do While n <= nmax
m=n+1 : p=2*m*(m+n)
Do While p <= pmax
If gcd(m,n)=1 Then
prim=prim+1
count=count+Int(pmax/p)
End If
m=m+2
p=2*m*(m+n)
Loop
n=n+1
Loop
triples = "Max Perimeter: " & pmax &_
", Total: " & count &_
", Primitive: " & prim
End Function
Function gcd(a,b)
c = a : d = b
Do
If c Mod d > 0 Then
e = c Mod d
c = d
d = e
Else
gcd = d
Exit Do
End If
Loop
End Function
</syntaxhighlight>
=={{header|Visual Basic}}==
{{Trans|VBA}}
{{works with|Visual Basic|5}}
{{works with|Visual Basic|6}}
{{works with|VBA|Access 97}}
{{works with|VBA|6.5}}
{{works with|VBA|7.1}}
<syntaxhighlight lang="vb">Option Explicit
Dim total As Long, prim As Long, maxPeri As Long
Public Sub NewTri(ByVal s0 As Long, ByVal s1 As Long, ByVal s2 As Long)
Dim p As Long, x1 As Long, x2 As Long
p = s0 + s1 + s2
If p <= maxPeri Then
prim = prim + 1
total = total + maxPeri \ p
x1 = s0 + s2
x2 = s1 + s2
NewTri s0 + 2 * (-s1 + s2), 2 * x1 - s1, 2 * (x1 - s1) + s2
NewTri s0 + 2 * x2, 2 * x1 + s1, 2 * (x1 + s1) + s2
NewTri -s0 + 2 * x2, 2 * (-s0 + s2) + s1, 2 * (-s0 + x2) + s2
End If
End Sub
Public Sub Main()
maxPeri = 100
Do While maxPeri <= 10& ^ 8
prim = 0
total = 0
NewTri 3, 4, 5
Debug.Print "Up to "; maxPeri; ": "; total; " triples, "; prim; " primitives."
maxPeri = maxPeri * 10
Loop
End Sub</syntaxhighlight>
{{out}}
<pre>Up to 100 : 17 triples, 7 primitives.
Up to 1000 : 325 triples, 70 primitives.
Up to 10000 : 4858 triples, 703 primitives.
Up to 100000 : 64741 triples, 7026 primitives.
Up to 1000000 : 808950 triples, 70229 primitives.
Up to 10000000 : 9706567 triples, 702309 primitives.
Up to 100000000 : 113236940 triples, 7023027 primitives.</pre>
=={{header|Wren}}==
{{trans|Go}}
Limited to a maximum perimeter of 10 billion in order to finish in a reasonable time.
<syntaxhighlight lang="wren">var sc = System.clock
var total = 0
var prim = 0
var maxPeri = 0
var newTri // recursive function so needs to be declared before it can be called
newTri = Fn.new { |s0, s1, s2|
var p = s0 + s1 + s2
if (p <= maxPeri) {
prim = prim + 1
total = total + (maxPeri/p).floor
newTri.call( 1*s0-2*s1+2*s2, 2*s0-1*s1+2*s2, 2*s0-2*s1+3*s2)
newTri.call( 1*s0+2*s1+2*s2, 2*s0+1*s1+2*s2, 2*s0+2*s1+3*s2)
newTri.call(-1*s0+2*s1+2*s2, -2*s0+1*s1+2*s2, -2*s0+2*s1+3*s2)
}
}
maxPeri = 100
while (maxPeri <= 1e10) {
prim = 0
total = 0
newTri.call(3, 4, 5)
var secs = (System.clock - sc).round
System.print("Up to %(maxPeri): %(total) triples, %(prim) primitives, %(secs) seconds")
maxPeri = 10 * maxPeri
}</syntaxhighlight>
{{out}}
Timings are for an Intel Core i7-8565U machine running Wren 0.4.0 on Ubuntu 20.04.
<pre>
Up to 100: 17 triples, 7 primitives, 0 seconds
Up to 1000: 325 triples, 70 primitives, 0 seconds
Up to 10000: 4858 triples, 703 primitives, 0 seconds
Up to 100000: 64741 triples, 7026 primitives, 0 seconds
Up to 1000000: 808950 triples, 70229 primitives, 0 seconds
Up to 10000000: 9706567 triples, 702309 primitives, 0 seconds
Up to 100000000: 113236940 triples, 7023027 primitives, 4 seconds
Up to 1000000000: 1294080089 triples, 70230484 primitives, 45 seconds
Up to 10000000000: 14557915466 triples, 702304875 primitives, 463 seconds
</pre>
=={{header|XPL0}}==
Simple minded algorithm:
<syntaxhighlight lang="xpl0">func GCD(N, D); \Return the greatest common divisor of N and D
int N, D, R; \numerator and denominator
[if D > N then
[R:=D; D:=N; N:=R];
while D > 0 do
[R:= rem(N/D);
N:= D;
D:= R;
];
return N;
];
int Max, PrimCnt, TripCnt, M, N, A, B, C, K, Prim;
[Max:= 10;
repeat PrimCnt:= 0; TripCnt:= 0;
for M:= 2 to Max do
for N:= 1 to M do
[if GCD(M,N) = 1 \coprime\ and
((M&1) = 0 xor (N&1) = 0) \one even\ then
[A:= M*M - N*N;
B:= 2*M*N;
C:= M*M + N*N;
Prim:= A+B+C;
if Prim <= Max then PrimCnt:= PrimCnt+1;
for K:= Max/Prim downto 1 do
if K*Prim <= Max then TripCnt:= TripCnt+1;
];
];
Format(6, 0);
Text(0, "Up to"); RlOut(0, float(Max));
RlOut(0, float(TripCnt)); Text(0, " triples,");
RlOut(0, float(PrimCnt)); Text(0, " primitives.^m^j");
Max:= Max*10;
until Max > 10_000;
]</syntaxhighlight>
{{out}}
<pre>
Up to 10 0 triples, 0 primitives.
Up to 100 17 triples, 7 primitives.
Up to 1000 325 triples, 70 primitives.
Up to 10000 4858 triples, 703 primitives.
</pre>
{{trans|Go}}
<syntaxhighlight lang="xpl0">int Total, Prim, MaxPeri;
proc NewTri(S0, S1, S2);
int S0, S1, S2, P;
[P:= S0 + S1 + S2;
if P <= MaxPeri then
[Prim:= Prim+1;
Total:= Total + MaxPeri/P;
NewTri(+1*S0-2*S1+2*S2, +2*S0-1*S1+2*S2, +2*S0-2*S1+3*S2);
NewTri(+1*S0+2*S1+2*S2, +2*S0+1*S1+2*S2, +2*S0+2*S1+3*S2);
NewTri(-1*S0+2*S1+2*S2, -2*S0+1*S1+2*S2, -2*S0+2*S1+3*S2);
];
];
[MaxPeri:= 10;
while MaxPeri <= 100_000_000 do
[Prim:= 0;
Total:= 0;
NewTri(3, 4, 5);
Format(10, 0);
Text(0, "Up to"); RlOut(0, float(MaxPeri)); Text(0, ":");
RlOut(0, float(Total)); Text(0, " triples,");
RlOut(0, float(Prim)); Text(0, " primitives.^m^j");
MaxPeri:= MaxPeri*10;
];
]</syntaxhighlight>
{{out}}
<pre>
Up to 10: 0 triples, 0 primitives.
Up to 100: 17 triples, 7 primitives.
Up to 1000: 325 triples, 70 primitives.
Up to 10000: 4858 triples, 703 primitives.
Up to 100000: 64741 triples, 7026 primitives.
Up to 1000000: 808950 triples, 70229 primitives.
Up to 10000000: 9706567 triples, 702309 primitives.
Up to 100000000: 113236940 triples, 7023027 primitives.
</pre>
=={{header|zkl}}==
{{trans|D}}
<syntaxhighlight lang
p:=a + b + c;
if(p>lim)
tri(lim, a - 2*b + 2*c, 2*a - b + 2*c, 2*a - 2*b + 3*c),
tri(lim, a + 2*b + 2*c, 2*a + b + 2*c, 2*a + 2*b + 3*c),
tri(lim, -a + 2*b + 2*c, -2*a + b + 2*c, -2*a + 2*b + 3*c)
);
}</
<
{{out}}
<pre>10: L(0,0)
Line 2,704 ⟶ 5,478:
</pre>
Max stack size is arbitrary but not adjustable.
=={{header|ZX Spectrum Basic}}==
ZX Spectrum: 8 bit microprocessor 3.5 Mhz doing all the work.
In an effort to get some decent speed the program is made to be as fast as it can.
It takes about 90 seconds for limit = 10 000 and 17 minutes for limit=100 000 and 3.5 hours for limit = 1000 000.
To set the limits.
Set in line nr: 1 L to the starting limit.
Set in line nr: 11 IF L<=(last limit to calculate)
Ex. start as limit 100 and end on limit 1000.
Set in line nr: 1 LET L=100.
Set in line nr: 11 IF L<=1000 THEN GO TO 2
<syntaxhighlight lang="zxbasic"> 1 LET Y=0: LET X=0: LET Z=0: LET V=0: LET U=0: LET L=10: LET T=0: LET P=0: LET N=4: LET M=0: PRINT "limit trip. prim."
2 FOR U=2 TO INT (SQR (L/2)): LET Y=U-INT (U/2)*2: LET N=N+4: LET M=U*U*2: IF Y=0 THEN LET M=M-U-U
3 FOR V=1+Y TO U-1 STEP 2: LET M=M+N: LET X=U: LET Y=V
4 LET Z=Y: LET Y=X-INT (X/Y)*Y: LET X=Z: IF Y<>0 THEN GO TO 4
5 IF X>1 THEN GO TO 8
6 IF M>L THEN GO TO 9
7 LET P=P+1: LET T=T+INT (L/M)
8 NEXT V
9 NEXT U
10 PRINT L;TAB 8;T;TAB 16;P
11 LET N=4: LET T=0: LET P=0: LET L=L*10: IF L<=100000 THEN GO TO 2</syntaxhighlight>
{{out}}
<pre>limit trip. prim.
10 0 0
100 17 7
1000 325 70
10000 4858 703
100000 64741 7026
1000000 808950 70229</pre>
[[Category:Geometry]]
| |||