Find all integer pairs that produce a given Loeschian number

Inspired by and drawns from Is this number Loeschian?

A positive integer $$\k\$$ is a Loeschian number if

• $$\k\$$ can be expressed as $$\i^2 + j^2 + i\times j\$$ for $$\i\$$, $$\j\$$ integers.

For example, the first positive Loeschian numbers are: $$\1\$$ ($$\i=1, j=0\$$); $$\3\$$ ($$\i=j=1\$$); $$\4\$$ ($$\i=2, j=0\$$); $$\7\$$ ($$\i=2, j=1\$$); $$\9\$$ ($$\i=-3, j=3\$$)1; ... Note that $$\i, j\$$ for a given $$\k\$$ are not unique. For example, $$\9\$$ can also be generated with $$\i=3, j=0\$$.

Other equivalent characterizations of these numbers are:

• $$\k\$$ can be expressed as $$\i^2 + j^2 + i\times j\$$ for $$\i, j\$$ non-negative integers. (For each pair of integers $$\i, j\$$ there's a pair of non-negative integers that gives the same $$\k\$$)

• There is a set of $$\k\$$ contiguous hexagons that forms a tesselation on a hexagonal grid (see illustrations for $$\k = 4\$$ and for $$\k = 7\$$). (Because of this property, these numbers find application in mobile cellular communication networks.)

• See more characterizations in the OEIS page of the sequence.

The first few Loeschian numbers are

0, 1, 3, 4, 7, 9, 12, 13, 16, 19, 21, 25, 27, 28, 31, 36, 37, 39, 43, 48, 49, 52, 57, 61, 63, 64, 67, 73, 75, 76, 79, 81, 84, 91, 93, 97, 100, 103, 108, 109, 111, 112, 117, 121, 124, 127, 129, 133, 139, 144, 147, 148, 151, 156, 157, 163, 169, 171, 172, 175, 181, 183, 189, 192...

1while ($$\i=-3, j=3\$$) produces 9, stick to non-negative integers, so ($$\i=0, j=3\$$).

Loeschian numbers also appear in determining if a coincident point in a pair of rotated hexagonal lattices is closest to the origin?

The challenge

Given a non-negative integer $$\k\$$, output all pairs of non-negative integers $$\i, j\$$ such that $$\i^2 + j^2 + i\times j=k\$$. If none are found (i.e. $$\k\$$ is not Loeschian) then return nothing or some suitable flag other than $$\(0, 0)\$$ since that produces the first Loeschian number, $$\0\$$.

For reversed order pairs like $$\(0, 4)\$$ and $$\(4, 0)\$$ either include both, or one member of the pair, but it should be the same for all cases (i.e. not sometimes one and other times both).

The program or function should handle (say in less than a minute) inputs up to $$\100,000\$$, or up to data type limitations.

This is code golf so shortest code wins.

Test cases

in       out
0      (0, 0)
1      (0, 1), (1, 0)
3      (1, 1)
4      (0, 2), (2, 0)
9      (0, 3), (3, 0)
12     (2, 2)
16     (0, 4), (4, 0)
27     (3, 3)
49     (0, 7), (3, 5), (5, 3), (7, 0)
147    (2, 11), (7, 7), (11, 2)
169    (0, 13), (7, 8), (8, 7), (13, 0)
196    (0, 14), (6, 10), (10, 6), (14, 0)
361    (0, 19), (5, 16), (16, 5), (19, 0)
507    (1, 22), (13, 13), (22, 1)
2028   (2, 44), (26, 26), (44, 2)
8281   (0, 91), (11, 85), (19, 80), (39, 65), (49, 56), (56, 49), (65, 39), (80, 19), (85, 11), (91, 0)
12103  (2, 109), (21, 98), (27, 94), (34, 89), (49, 77), (61, 66), (66, 61), (77, 49), (89, 34), (94, 27), (98, 21), (109, 2)
• Nice first challenge! For output, it's fine how it is, although allowing both is probably better choice. Oct 11 '20 at 13:35
• I think the best option is to just leave it undefined (as in answers pick how they handle it). Most challenges here tend not to be too strict with I/O stuff. Oct 11 '20 at 13:41
• @RedwolfPrograms done, thanks!
– uhoh
Oct 11 '20 at 13:44
• Shouldn't (2,109) be included for 12103?
– kjs
Oct 11 '20 at 14:12
• No problem, and welcome to the site! Oct 11 '20 at 15:08

05AB1E, 11 9 bytes

Thanks to @ovs for -2!

ÝãʒãÀ¦POQ

Oh my, I am getting in the hang of code golfing!

Prints a list of all the valid pairs (e.g [[1, 0], [0, 1]]). If there are none, the list is empty ([]). Also outputs both of any reverse integer pairs.

Try it online!

How?

You may count this as a port of the other answers, but I took a look only at the Husk answer before writing the program!

Ý             # Push a list of all numbers from 0 to the input.
ã            # Push the cartesian power of lists. (Basically, finding all possible pairs)
ʒ           # For each pair...
ãÀ¦        # Find all other permutations of the pair.
P       # Multiply each permutation.
Q     # If the result is not equal to the input, yeet (throw) them from the list.
# Automatically print the pairs not yeeted.
• ã (Cartesian power) can replace . Using the same builtin, ãÀ¦PO is one byte shorter than DPsnO+.
– ovs
Oct 11 '20 at 19:30
• @ovs Nice one! (How did I not see ã in info.txt right beside â) Oct 11 '20 at 19:49
• Nice answer! Just a small FYI, the DPsnO+ could also have been PynO+. Although I like the current ãÀ¦PO of @ovs more. :) Oct 12 '20 at 9:21

f x|l<-[0..x]=[(i,j)|i<-l,j<-l,i*i+j*j+i*j==x]

Try it online!

Husk, 15 13 bytes

fo=¹§+Πṁ□π2…0

Try it online!

-2 bytes from Zgarb.

Outputs [] for non-Loeschians.

Explanation

fo=¹§+Πṁ□π2…0
…0 range from 0..n
π2   create all possible pairs using 0..n
fo            filter by the following two functions:
§         f: fork: § f g h x = f (g x) (h x)
ṁ□        sum of squares
Π          and fold by multiplication
=¹          g: is that equal to 1?
• F* can be Π and ´×e can be π2 Oct 11 '20 at 17:51

Jelly, 1110 9 bytes

ŻpḋÄ$=¥Ƈ Try it online! Outputs [] for non-Loeschian numbers -1 byte thanks to Sisyphus Not particularly efficient, but that can be fixed for an additional 2 bytes. Uses the fact that a Loeschian number can be expressed as $$\i\times i + j\times(i+j)\$$ by using Jelly’s vectorisation and cumulative sum. How it works ŻpḋÄ$=¥Ƈ - Main link. Takes n on the left
Ż         - Yield [0, 1, ..., n]
p       - Cartesian product with itself, yielding [[0, 0], [0, 1], ..., [n, n]]
¥Ƈ - Filter the pairs, keeping those where the following is true:
$= - The pair equals n after the following is done: Ä - Cumulative sum. Yield [i, i+j] ḋ - Dot product with [i, j]; Yields i×i + j×(i+j) • Would the downvoter mind explaining why they cast their vote? This is a perfectly valid answer Oct 11 '20 at 15:32 • Nice answer. You can do ŻpḋÄ$=¥Ƈ for 9. Oct 18 '20 at 2:54
• @Sisyphus Very nice, I knew Jelly had a shorter way to multiply-then-sum! That should also let me save a byte on my answer to the other Loeschian numbers question Oct 18 '20 at 11:27

Scala, 64 62 bytes

k=>0.to(k)flatMap(i=>0.to(k)filter(j=>i*i+j*j+i*j==k)map(i->))

Try it online!

Thanks to user for -2

Japt, 13 bytes

ô ï f@¶Xx²+X×

Try it

ô ï f@¶Xx²+X×     :Implicit input of integer U
ô                 :Range [0,U]
ï               :Cartesian product
f             :Filter by
@            :Passing each X through the following function
¶           :  Is U equal to
Xx         :  X reduced by addition
²        :  After squaring each
+X×     :  Plus X reduced by multiplication

Wolfram Language (Mathematica), 46 bytes

Solve[i^2+j^2+i*j==#&&i>=j>=0,{i,j},Integers]&

Try it online!

Perl 5, 8177 74 bytes

sub{map{$i=$_;grep{$k==$i**2+$i*$_+$_**2&&($_=[$i,$_])}$i..$k}0..($k=pop)} Try it online! A bit ungolfed: sub f { my$k=pop;                      #gangnam style, k=pop from input
grep { $k==pop@$_ }             #pop last of three elems
#...in the candidate array
#...and return as result
#...if last = i*i+i*j+j*j = k
map  {                          #two loops from 0 to sqrt $k my$i=$_; #outer loop var map { my$j=$_; #inner loop var [$i, $j,$i*$i+$i*$j+$j*$j] #result candidate } 0..sqrt$k                     #or  $i..sqrt$k  to return only i<=j
}
0..sqrt$k } Note: Saving bytes by removing the two sqrt makes it run A LOT slower, but it will still return the correct result. • Can you make the inner loop go from 0 to i (inclusive) instead? That might also make your code faster. – Neil Oct 11 '20 at 17:19 • Yes, that halves the result sets and the runtime. Still ok though. – kjs Oct 11 '20 at 18:50 • Well, I was thinking it also makes your code shorter, which is good, given that this is code-golf... – Neil Oct 11 '20 at 18:52 Scala, 53 bytes | =>for(i<-0 to|;j<-0 to|if| ==i*i+j*j+i*j)yield(i,j) Try it online! C (gcc), 75 bytes i,j;f(x){for(i=j=x;~j;i-=!i?j--,-x:1)i*i+j*j+i*j-x||printf("(%d,%d)",i,j);} Try it online! JavaScript (V8), 57 bytes Prints the pairs $$\(x,y),\:x\le y\$$. n=>{for(y=n+1;x=y--;)for(;x--;)x*x+y*y+x*y-n||print(x,y)} Try it online! Charcoal, 25 bytes ＮθＦ⊕₂θＦ⊕ι¿⁼θ⁻Ｘ⁺ικ²×ικＩ⟦ικ Try it online! Link is to verbose version of code. Only outputs those pairs where i>=j. speeds the code up so that the larger test cases complete within a minute, but it is not needed for smaller test cases. Explanation: Ｎθ Input k. Ｆ⊕₂θ Loop i from 0 to √k inclusive. Ｆ⊕ι Loop j from 0 to i inclusive. ¿⁼θ⁻Ｘ⁺ικ²×ικ If k=(i+j)²-ij, then... Ｉ⟦ικ Output i and j on separate lines. Just for fun, here's a 73-byte Retina 1.0 answer that only finds nontrivial solutions (i.e. neither i nor j is zero): .+ * L$w`^((_)+)(?=(?<-2>\1)+(?(2)$.)(_(_)*)(?<-4>\1\3)*$(?(4).))
$.1$.3

Try it online! Very slow, so don't try anything over about 500.

• Wow, running Try it seems incredibly fast compared to several other answers here, certainly part of that is due to the efficiency of this approach but is that all of it or is Charcoal also inherently faster somehow? I'm new here and I don't really understand how TIO works or implements the languages.
– uhoh
Oct 11 '20 at 18:07
• @uhoh Charcoal is written in Python 3 so any Python answer should beat it. It's not trying to be an efficient approach but it doesn't have to allocate memory like a Cartesian product approach would.
– Neil
Oct 11 '20 at 18:56
• @uhoh Also, it's much much slower without the of course.
– Neil
Oct 11 '20 at 18:57
• I'm intrigued, will read up on this. Thanks!
– uhoh
Oct 11 '20 at 19:07

Wolfram Language (Mathematica), 42 bytes

Array[(+##)^2-##&,{#,#}+1,0]~Position~#-1&

Try it online!

Gets slow on larger inputs.

Array[          (* Create a table of *)
(+##)^2-##&,    (* (i+j)^2-i j *)
{#,#}+1,0]      (* for i,j = 0...k *)
~Position~#-1   (* and find where that expression equals k *)

Python 3 - 79 70 Bytes

Just the trivial solution to get started with Python 3

k=9  # not counted
def f(k):
a=range(k);return [(i,j) for i in a for j in a if k==i*i+j*j+i*j]
print(f(k))   # not counted

Try it online!

Thanks to @Jakque

lambda, testing framework in TIO, spaces, k+1, factorization

lambda k:[(i,j)for i in range(k+1)for j in range(k+1)if(i+j)*i+j*j==k]

Try it online!

• Yon need to use range(k+1) otherwise you don't pass the test for k=0 and k=1. You can also gain some bytes by removing some useless space : def f(k):a=range(k+1);return[(i,j)for i in a for j in a if k==i*i+j*j+i*j] works perfectly and is cheaper ;-) May 19 '21 at 17:02
• You can also do better with an anonymous lambda function : Try it online! May 19 '21 at 17:04
• by foctorizing i*i+j*j+i*j you can have (i+j)*i+j*j. It doesn't gain any byte in itself but now you can remove the space after the if like that : if(i+j)*i+j*j==k to save one byte May 20 '21 at 8:52