# Trithagorean Triples

A Pythagorean Triple is a positive integer solution to the equation:

A Trithagorean triple is a positive integer solution to the equation:

Where Δn finds the nth triangular number. All Trithagorean triples are also solutions to the equation:

Given a positive integer c, output all the pairs of positive integers a,b such that the sum of the ath and bth triangular numbers is the cth triangular number. You may output the pairs in whatever way is most convenient. You should only output each pair once.

This is

## Test Cases

2: []
3: [(2, 2)]
21: [(17, 12), (20, 6)]
23: [(18, 14), (20, 11), (21, 9)]
78: [(56, 54), (62, 47), (69, 36), (75, 21), (77, 12)]
153: [(111, 105), (122, 92), (132, 77), (141, 59), (143, 54), (147, 42), (152, 17)]
496: [(377, 322), (397, 297), (405, 286), (427, 252), (458, 190), (469, 161), (472, 152), (476, 139), (484, 108), (493, 54), (495, 31)]
1081: [(783, 745), (814, 711), (836, 685), (865, 648), (931, 549), (954, 508), (979, 458), (989, 436), (998, 415), (1025, 343), (1026, 340), (1053, 244), (1066, 179), (1078, 80), (1080, 46)]
1978: [(1404, 1393), (1462, 1332), (1540, 1241), (1582, 1187), (1651, 1089), (1738, 944), (1745, 931), (1792, 837), (1826, 760), (1862, 667), (1890, 583), (1899, 553), (1917, 487), (1936, 405), (1943, 370), (1957, 287), (1969, 188)]
2628: [(1880, 1836), (1991, 1715), (2033, 1665), (2046, 1649), (2058, 1634), (2102, 1577), (2145, 1518), (2204, 1431), (2300, 1271), (2319, 1236), (2349, 1178), (2352, 1172), (2397, 1077), (2418, 1029), (2426, 1010), (2523, 735), (2547, 647), (2552, 627), (2564, 576), (2585, 473), (2597, 402), (2622, 177), (2627, 72)]
9271: [(6631, 6479), (6713, 6394), (6939, 6148), (7003, 6075), (7137, 5917), (7380, 5611), (7417, 5562), (7612, 5292), (7667, 5212), (7912, 4832), (7987, 4707), (8018, 4654), (8180, 4363), (8207, 4312), (8374, 3978), (8383, 3959), (8424, 3871), (8558, 3565), (8613, 3430), (8656, 3320), (8770, 3006), (8801, 2914), (8900, 2596), (8917, 2537), (9016, 2159), (9062, 1957), (9082, 1862), (9153, 1474), (9162, 1417), (9207, 1087), (9214, 1026), (9229, 881), (9260, 451), (9261, 430), (9265, 333)]

• Can we output repeated pairs? Example, for 21 output [(17, 12), (20, 6), (12, 17), (6, 20)] Commented May 24, 2017 at 15:52
• I thought you were asking us to find a^3+ b^3 = c^3. :D Commented May 24, 2017 at 15:54
• @LuisMendo No. I'll include this in the question. Commented May 24, 2017 at 15:55
• @BetaDecay MATL, 0 bytes Commented May 24, 2017 at 16:05
• @EriktheOutgolfer a^3+ b^3 = c^3 is known to have no integer solutions; see Fermat's last theorem Commented May 24, 2017 at 16:10

# Mathematica, 5349 48 bytes

Solve[{x.(x+1)==#^2+#,a>=b>0},x={a,b},Integers]&


Example:

In[1]:= Solve[{x.(x+1)==#^2+#,a>=b>0},x={a,b},Integers]&[21]

Out[1]= {{a -> 17, b -> 12}, {a -> 20, b -> 6}}

• ooohh, nice vectorization, way better than what I would have done Commented May 24, 2017 at 17:07

# MATL, 17 13 bytes

:Ys&+G:s=R&fh


Each pair is output with the smaller number first.

Try it online!

### Explanation

Consider input 3.

:      % Implicitly input n. Push [1 2 ... n]
% STACK: [1 2 3]
Ys     % Comulative sum
% STACK: [1 3 6]
&+     % All pairwise sums
% STACK: [2 4 7; 4 6 9; 7 9 12]
G:s    % Push 1+2+...+n
% STACK: [2 4 7; 4 6 9; 7 9 12], 6
=      % Is equal?
% STACK: [0 0 0; 0 1 0; 0 0 0]
R      % Upper triangular part of matrix. This removes duplicate pairs
% STACK: [0 0 0; 0 1 0; 0 0 0]
&f     % Push row and column indices (1-based) of non-zero entries
% STACK: 2, 2
h      % Concatenate horizontally. Implicitly display
% STACK: [2, 2]

• Explanation please? Commented May 24, 2017 at 16:01
• @EriktheOutgolfer Sure, I'll add it later in the day Commented May 24, 2017 at 16:03
• Just make sure you output the unique pairs, that mainly why I asked. Commented May 24, 2017 at 16:05
• @EriktheOutgolfer Yes, they are unique (R takes care of that) Commented May 24, 2017 at 16:07
• @EriktheOutgolfer Explanation added Commented May 24, 2017 at 17:30

# Jelly, 12 bytes

j‘c2ḅ-
ŒċçÐḟ


Try it online!

### How it works

ŒċçÐḟ   Main link. Argument: c

Œċ      Yield all 2-combinations w/repetition of elements of [1, ..., c].
çÐḟ   Filterfalse; keep only 2-combinations for which the helper link returns 0.

j‘c2ḅ-  Helper link. Left argument: [a, b]. Right argument: c

j       Join [a, b] with separator c, yielding [a, c, b].
‘      Increment; yield [a+1, c+1, b+1].
c2    Combination count; compute [C(a+1,c), C(c+1,c), C(b+1,c)], yielding
[½a(a+1), ½c(c+1), ½b(b+1)].
ḅ-  Convert from base -1 to integer, yielding
½(-1)²a(a+1) + ½(-1)¹c(c+1) + ½(-1)⁰b(b+1) = ½(a(a+1) - c(c+1) + b(b+1)),
which is 0 if and only if a(a+1) + b(b+1) = c(c+1).


## Python 2, 69 bytes

Try it online

lambda c:[(a,b)for a in range(c)for b in range(a+1)if~a*a==c*~c-~b*b]


-9 bytes, thanks to @WheatWizard

• And ~a*a==c*~c-~b*b is a byte shorter than that. Try it online! Commented May 24, 2017 at 16:12
• @WheatWizard Good stuff :D Commented May 24, 2017 at 16:28

# Jelly, 16 14 bytes

RS
ŒċÇ€S$⁼¥ÐfÇ  Try it online! This is too long for sure... Explanation: ŒċÇ€S$⁼¥ÐfÇ (main) Arguments: z
Œċ             Return [[1, 1], [1, 2], ..., [1, z], [2, 2], ..., [z, z]]
Ç    Return T(z)
Ç€S$⁼¥Ðf Only keep the pairs such as ΣT(a, b)=T(z) RS (helper 1) Arguments: z R [1, 2, ..., z] S Take the sum  # AWK, 72 bytes {for(n=$1;++i<=n;)for(j=i;j<=n;++j)if(i^2+j^2+i+j==n^2+n)$0=$0" "i":"j}1


Try it online!

Output is c a1:b1 a2:b2 .... The TIO link has 4 extra bytes i=0; to allow for multiline input.

This isn't efficient at all, but it works. :)

# PHP, 94 Bytes

for($a=$c=$argn;$a--;)for($b=$a;$b;$b--)$a**2+$a+$b**2+$b!=$c**2+$c?:$e[]=[$a,$b];print_r($e);


Try it online!

f c=[(a,b)|a<-[1..c],b<-[1..a],a^2+a+b^2+b==c^2+c]


Usage example: f 21 -> [(17,12),(20,6)]. Try it online!

Uses the 2nd equation.

# J, 35 bytes

1+[:~.,~/:~@#:[:I.@,@({:=+/~)2!2+i.


Try it online!

# Axiom, 281204196 191 bytes

q(b,m)==(r:=1+4*m;v:=4.*b*(b+1);r<v=>0;(sqrt(r-v)-1)/2);g(c:NNI):Any==(r:List List INT:=[];i:=0;repeat(i:=i+1;i>=c=>break;w:=q(i,c^2+c);w>=i and fractionPart(w)=0=>(r:=cons([w::INT,i],r)));r)


test and ungolf

-- if m=c^2+c than a^2+a+b^2+b-m=0 has the solutions [a,b] with a>0,b>0
-- if it is used a=(-1+sqrt(1+4*m-4*(b)*(b-1)))/2 because the other return a<0
-- o(b,m) return that solution if 1+4*m-4*(b)*(b-1)>0 [so exist in R sqrt] else return 0
o(b,m)==(r:=1+4*m;v:=4.*b*(b+1);r<v=>0;(sqrt(r-v)-1)/2)

--it Gets one not negative integer c; return one list of list(ordered) of 2 integers
--[a,b] with  a^2+a+b^2+b=c^2+c
gg(c:NNI):List List INT==
r:List List INT:=[]  -- initialize the type make program more fast at last it seems 10x
i:=0
repeat
i:=i+1
i>=c=>break
w:=o(i,c^2+c)
w>=i and fractionPart(w)=0=>(r:=cons([w::INT,i],r))
r

(6) -> [[i,g(i)]  for i in [2,3,21,23,78,153,496,1081,1978,2628,9271]]
(6)
[[2,[]], [3,[[2,2]]], [21,[[17,12],[20,6]]], [23,[[18,14],[20,11],[21,9]]],
[78,[[56,54],[62,47],[69,36],[75,21],[77,12]]],
[153,[[111,105],[122,92],[132,77],[141,59],[143,54],[147,42],[152,17]]],

[496,
[[377,322], [397,297], [405,286], [427,252], [458,190], [469,161],
[472,152], [476,139], [484,108], [493,54], [495,31]]
]
,

[1081,
[[783,745], [814,711], [836,685], [865,648], [931,549], [954,508],
[979,458], [989,436], [998,415], [1025,343], [1026,340], [1053,244],
[1066,179], [1078,80], [1080,46]]
]
,

[1978,
[[1404,1393], [1462,1332], [1540,1241], [1582,1187], [1651,1089],
[1738,944], [1745,931], [1792,837], [1826,760], [1862,667], [1890,583],
[1899,553], [1917,487], [1936,405], [1943,370], [1957,287], [1969,188]]
]
,

[2628,
[[1880,1836], [1991,1715], [2033,1665], [2046,1649], [2058,1634],
[2102,1577], [2145,1518], [2204,1431], [2300,1271], [2319,1236],
[2349,1178], [2352,1172], [2397,1077], [2418,1029], [2426,1010],
[2523,735], [2547,647], [2552,627], [2564,576], [2585,473], [2597,402],
[2622,177], [2627,72]]
]
,

[9271,
[[6631,6479], [6713,6394], [6939,6148], [7003,6075], [7137,5917],
[7380,5611], [7417,5562], [7612,5292], [7667,5212], [7912,4832],
[7987,4707], [8018,4654], [8180,4363], [8207,4312], [8374,3978],
[8383,3959], [8424,3871], [8558,3565], [8613,3430], [8656,3320],
[8770,3006], [8801,2914], [8900,2596], [8917,2537], [9016,2159],
[9062,1957], [9082,1862], [9153,1474], [9162,1417], [9207,1087],
[9214,1026], [9229,881], [9260,451], [9261,430], [9265,333]]
]
]
Type: List List Any


# CJam, 30 28 bytes

{_,2m*:$_&f+{{_)*}%~+=},1f>}  Anonymous block expecting its argument on the stack and leaving the result on the stack. Try it online! ### Explanation I will refer to the input as n _, e# Copy n, and get the range from 0 to n-1. 2m* e# Get the 2nd Cartesian power of this range. :$_&   e# Sort the pairs and deduplicate, to get all unique pairs.
f+     e# Prepend n to each pair.
{      e# Filter these triplets; keep only those that give a truthy result:
{     e#  Map this block over the triplet:
_)*  e#   Multiply x by x+1. (i.e. x^2 + x)
}%    e#  (end map)
~+=   e#  Check if the sum of the second and third is equal to the first.
},     e# (end filter)
1f>    e# Remove the first element from all remaining triplets.


# Pyth - 23 21 bytes

L*bhbfqyQ+yhTyeT.CUQ2


Try it

L*bhbfqyQ+yhTyeT.CUQ2
L*bhb                     Define y(b)=b*(b+1)
.CUQ2     All pairs of numbers less than the input
fqyQ+yhTyeT          Filter based on whether y(input) == y(1st elem. of pair) + y(2nd elem. of pair)


## JavaScript (ES6), 83 bytes

c=>[...Array(c*c)].map((_,x)=>[x%c,x/c|0]).filter(([a,b])=>a>=b&a++*a+b++*b==c*c+c)


### Test cases

Omitting here the largest inputs that take too much time for the snippet.

let f =

c=>[...Array(c*c)].map((_,x)=>[x%c,x/c|0]).filter(([a,b])=>a>=b&a++*a+b++*b==c*c+c)

console.log(JSON.stringify(f(2)))
console.log(JSON.stringify(f(3)))
console.log(JSON.stringify(f(21)))
console.log(JSON.stringify(f(23)))
console.log(JSON.stringify(f(78)))
console.log(JSON.stringify(f(153)))
console.log(JSON.stringify(f(496)))

# Husk, 12 11 bytes

fo=Σ¹ṁΣüOπ2


Try it online!

Σ has a mode where it returns the nth triangular number, which is very useful here.

-1 byte from Dominic van Essen.

## Explanation

fo=Σ¹ṁΣüOπ2
π2  cartesian power 2 for range 1..n(all possible pairs)
üO    uniquify, ignoring ordering
fo           filter by the following 2 functions:
ṁΣ      triangular numbers of the pair summed
=          equals
Σ¹        nth triangular number?

• 11 bytes without the ḣ (using cpowN function of π) Commented Oct 18, 2020 at 12:49
• Obviously that's a thing aaaaaa Commented Oct 18, 2020 at 13:29