LÖVE

I created a small helper function for mapping vectors to whole numbers.
https://github.com/IAmSegfault/N3N This is useful for creating a deterministic rng for each chunk in your game world, for example the procedurally generated universe in the 1984 version of Elite. I haven't tested it thoroughly, so I do not know yet if there are duplicate numbers due to rounding errors, etc. All you need to do to use it is modify the magic number inside and require the module. Let me know what you think of it!

Nice idea. A few comments:

• It would be nice to mention that it maps vectors of integers to integers.
• In double precision, there are about 2^54 consecutive integers; having three degrees of freedom gives about 2^(54*3) = 2^162 possible variations, each of which must map to a different number. You can't possibly have a reversible mapping from every number between 0 and 2^162 to a number from -2^53 to 2^53, for the same reason that you can't make a reversible mapping from integers between 0 and 3 to integers between 0 and 5. There MUST be duplicate numbers. Perhaps suggest shorter magic numbers, so that they don't grow too much. Float precision decreases as the numbers get bigger, which increases the chance of collision.
• In the comment, you use / to indicate combinations. That confused the hell out of me until I read the code. Maybe specify it as e.g. binom(x,1)+binom(x+y+1, 2)+...? (or combine(), though that's not a very used name to my knowledge)
• x1 and x2 are global. I don't think you intended that?
• Why : instead of . in the function definition? It seems like an isolated function, kinda like e.g. math.sin, and it doesn't use self.

Thanks for posting info on the magic number affecting collision chances. I didn't consider it too important to mention that it is only for mapping integers since I wrote it mainly to work with a chunking system. Also thanks for pointing out that I accidentally added a few globals.

If i'm understanding right, then the point still stands about your code not being truly bijective, since as pgimeno said, and due to lua (and Löve) using double precision floats for all number representations (unless you use luaJIT FFI and C integer types specifically), so you will never be able to force 2^216 numbers onto 2^54 bijectively;

you could do 2^18, or 262144, since 2^(18*3) would be 2^(54) and that should be completely bijective, if my math is right.

Thanks Zorg, this clarifies even further. I guess the stackexchange post I found was wrong about the formula being bijective. Math isn't exactly my strong point.

I started adding constraints on the coordinates, I should have this actually working in a day or two.
*edit* I added a lazy constraint on the coordinates and I calculated the maximum absolute value of the magic number.

There is also another bug in combine(), you should divide by i, not by k. And the x variable defined there is global too.

I found another problem. You map negatives to positives in combine(), not on the x, y, z. This results in a map that is not bijective for negative input, and indeed there are many collisions, even after fixing the previous bug.

The fix is to apply the mapping from negative to positive on the x, y, z inputs, not on the numbers that you send to combine().

With all that fixed, it's harder to get collisions with low numbers. One I've found is (-125000, -125000, 0) and (125000, 125001, 0).

Ok, so with the constraint of r never exceeding 262144 the limits on the coordinates are (-17, -17, -17) to (18,18,18) for 3d and (-27, -27) to (28, 28) for 2d. I comapred every mapping in the 2d area coordinates and the 3d volume coordinates and found no collisions. This function is now truly bijective. The small limits on the coordinates now mean that the function is only useful for sectors or possibly superchunks, not chunks.
I also added an update so you can change the magic number at runtime within limits of 262144 * magic <= 2^54.
If LOVE and LuaJIT ever update to Lua 5.3 I will update this library to take advantage of 64 bit integers.

You can take advantage of 64 bit integers by using FFI. It's not too easy to handle, and I'm not sure if you can carry on all of the operations necessary to succeed in this case.

Warning: wall of text ahead.

I've now understood the formula. It's basically a method of numbering spheres arranged as a tetrahedral pyramid, which when adequately transformed and using 1x1x1 cubes instead of spheres, fills N³, that's why it's bijective for the naturals. The pre-mapping of negative to positive numbers makes it count in the shape of an octahedron centred at the origin.

To be picky, your implementation is actually not bijective, but it's injective, which is good enough, as that guarantees that you get a different number for each coordinate, even if there are numbers that don't correspond to any coordinate. The reason it's not bijective is that your function 'constraint' does not produce zero as an output.

Since zero isn't a possible input to the next phase, it also doesn't produce any of the numbers that the tetrahedral function would return for an input value of 0 in any coordinate, for example for the numbers from 0 to 39 it can only generate four numbers: 14, 24, 27 and 28 (10%), but fortunately the percentage grows as the numbers do. It's basically discarding all numbers that lay in the X, Y and Z planes, i.e. three layers of balls, if you prefer to see it that way. That makes it an injection from Z³ to N, not a bijection.

Here's my implementation of the tetrahedral pyramid method: Note how my z_to_n function can return 0 and therefore in theory biject() is a bijection between the N and Z³.

The tetrahedral function is written in a way that will never suffer 53-bit overflow in the intermediate results if the final theoretical result fits in 53 bits (something that very few binomial calculation algorithms care about, thus reducing their range of proper operation). Division by 2 is always exact for floats, because they work in binary, but division by 3 isn't, so that one needs a case-by-case check in order to only divide by 3 the number that is already a multiple of 3 (actually by 6 which is 3*2, but again, division by 2 is exact).

But tetrahedral numbering is not the only approach to a bijective map, and possibly not the best approach either.

The tetrahedral method is cool, in that you can have one different seed per supported square that lies within an octahedron, and you're fine as long as you don't hit a coordinate that doesn't fit within the octahedron for which there's no overflow (I'll call this the "maximal lossless octahedron" or MLO). But most of the time, in games, you want to work with an NxNxN cube or even NxNxM (if height range is different than area range). And that admits a simpler bijective formula, for which you can easily know the range of values it can take, and can therefore know exactly how far you can push it. Working with the tetrahedron formula, to do the same you'd have to find the biggest cube that can fit within the MLO, and that would mark the limits of your playing area. It's much better to define the limits in advance and play with them until you're satisfied, than to have to stick to limits that are intrinsic to the algorithm being employed, don't you think?

The method is basically to take the coordinates as a mixed radix representation of the seed. It works as follows (in pseudocode):

RangeX = MaxX - MinX + 1
RangeY = MaxY - MinY + 1
RangeZ = MaxZ - MinZ + 1 -- not used in the seed formula, but used to calculate the range of the seed
Seed = (z - MinZ) * RangeY * RangeX + (y - MinY) * RangeX + (x - MinX)

or, to save a multiplication:
Seed = ((z - MinZ) * RangeY + (y - MinY)) * RangeX + (x - MinX)

Now, the maximum value of Seed is exactly RangeX * RangeY * RangeZ - 1, that is, you have exactly RangeX * RangeY * RangeZ possible seeds when x, y, z are within their respective ranges. When RangeX, RangeY and RangeZ are powers of 2, you can count bits, which even simplifies the calculations more.

For example, love.math.newRandomGenerator can use two 32-bit input seeds. If you restrict the ranges to, say, 10, 10 and 10 bits, you can have a 30-bit number in the low part of the seed and any arbitrary 32-bit magic number in the high part. That gives you a range of 2^10 = 1024 numbers in each direction, which you can distribute for example like this: from -512 to +511 inclusive, plus a magic number from 0 to 4294967295. Or you can have a cuboid instead of a cube, using 11 bits for X and Y and 10 bits for Z. Or 16 bit for x, 16 bits for y, 16 bits for z and 16 for the magic number; the latter two would go into the high part. That would give you a range of -32768 to 32767 (or from 0 to 65535 or however you want to distribute it) in every direction, and a range of 0 to 65535 for the magic number.

If you've read all the way up to here, congratulations!

Wow, this implementation outdoes what I came up with, great work! I'll rename my function to injection with the new info you provided.