Trivial examples:  Pythagorean triangles [[0, 0, 0], [2ab, 0, 0], [0, 0, a^2-b^2]] , where a and b - any unequal positive integers.

Nontrivial examples can be found using a simple code:

N:=10:

L:=[]:

for x1 to N do

for y1 from x1 to N do

for z1 from y1 to N do

for x2 to N do

for y2 to N do

for z2 to N do

a:=sqrt((x2-x1)^2+(y2-y1)^2+(z2-z1)^2):

if type(a, positive) and type(sqrt(x1^2+y1^2+z1^2)+ sqrt(x2^2+y2^2+z2^2)+a, integer) then L:=[op(L), {[0, 0, 0], [x1, y1, z1], [x2, y2, z2]}]: fi:

od: od: od: od: od: od:

nops(L);

L;

Yes, of course! See help for the  addcoords  command.

Your code can be significantly improved:

1) Since always  x2=x1+x3 and y2=y1+y3 and z2=z1+z3 , then variables x2, y2, z2  should not be written as a variables in a loop. This will speed up program execution.

2) You can exclude the  different trivial variants. If the points (x1, y1, z1) and (x3, y3, z3) lie in a coordinate plane, eg, xOy, then it will always be

(x1, y1, z1) = (x, y, 0) , (x3, y3, z3) = (-y, x, 0)  or  (x1, y1, z1) = (x, y, 0) , (x3, y3, z3) = (y, -x, 0)

(x and y - any numbers)

3) You can also foreclose variants, which are obtained from one another by simple transformations, such as reflections in the coordinate planes and axes.

4) Lists  [[0,0,0], [x1, y1, z1], [x2, y2, z2], [x3, y3, z3]]  and [[0,0,0], [x3, y3, z3], [x2, y2, z2], [x1, y1, z1]] give the same square, so from two lists, you can leave one list.

To implement  3) and 4) can be used  ListTools[Categorize]  command.

Write your 5 variables   (a,b,c,d,e)  in the list.

Example:

Case:= [ [1, 2, 3, 1.5, 6], [3, 3, 1, 3, 8], [5, 3, 3, 5, 3] ]:
Case[1], Case[2], Case[3];

       [1, 2, 3, 1.5, 6], [3, 3, 1, 3, 8], [5, 3, 3, 5, 3]

This problem can be solved similar to your problem of an equilateral triangle. Take advantage of the fact that if [[0,0,0], [x1, y1, z1], [x2, y2, z2], [x3, y3, z3]] are the vertices of a square, then x2=x1+x3, y2=y1+y3, z2=z1+z3. It is therefore sufficient to seek only two vertices [x1, y1, z1], [x3, y3, z3] .

Let us the solution to be non trivial if none of the vertices of other than [0, 0, 0] , does not lie in the coordinate plane, and the absolute values ​​of the coordinates of each vertex are unequal and integer numbers. The simplest nontrivial solution is  [[0, 0, 0], [6, 3, -2], [9, 1, 4], [3, -2, 6]] . Interestingly, the length of the side of the square is also an integer equal to 7.

If I understand your problem, it is simple. Should be in the loop to multiply all the elements of each row to the least common multiple of their denominators.

Code of the procedure:

DoIntegers:=proc(A)

local m, n, i, j, d, B, f;

uses LinearAlgebra;

m:=RowDimension(A);

n:=ColumnDimension(A);

for i to m do

d:=ilcm(seq(denom(A[i,j]), j=1..n));

for j to n do

B[i,j]:=A[i,j]*d;

od; od;

f:=(i,j)->B[i,j];

B:=Matrix(m, n, f);

end proc:

Example:

A:=Matrix([[1, 1/2, 1/3], [1/4, 1/5, 1/6], [7, 8, 9]]);

B=DoIntegers(A);

1) Your code is correct.

2) Since the triangle with the vertex at an arbitrary point can be obtained from the triangle with the vertex at the origin of the usual translation, there is no need to write new code.

3) Rewrote my previous code. Now the triangles, differing only in the numbering of the vertices are not distinct. In this code, the user himself determines the required number of triangles (parameter N). The search range is not fixed from the beginning, but gradually expands as new triangles. Therefore, if the number N is not too large, then the procedure is fast.

N_equil_triangles:=proc(N)

local x1, y1, z1, x2, y2, z2, k, n, a;

global L;

k:=0:

L:=[]:

for n while k<N do

for x1 from -n to n while k<N do

for y1 from -n to n while k<N do

for z1 from -n to n while k<N do

for x2 from -n to n while k<N do

for y2 from -n to n while k<N do

for z2 from -n to n while k<N do

a:=x1^2+y1^2+z1^2:

if a>0 and a=x2^2+y2^2+z2^2 and a=(x2-x1)^2+(y2-y1)^2+(z2-z1)^2 then

L:=[op(L),{[0,0,0], [x1,y1,z1], [x2,y2,z2]}]: k:=nops(L): fi:

od: od: od: od: od: od: od:

end proc:

Example:

t:=time():

N_equil_triangles(1000):

L[500..510];

Time=cat(convert(time()-t, symbol),`  sec`);

I'm not sure that all colors are identified correctly. Details can be improved:

R:=3: r:=1: h:=0.25:  # h specifies the width of the band

Colors:=["Red", "HotPink",  "Blue", "BlueViolet", "DeepSkyBlue", "LightBlue", "LightSteelBlue", "LimeGreen", "Lime", "Yellow", "Orange", "Gold"]:

E:=seq(seq(plottools[line]([(R+r*cos(11/3*t))*cos(t), (R+r*cos(11/3*t))*sin(t), r*sin(11/3*t)],[(R+r*cos(11/3*t))*cos(t+h), (R+r*cos(11/3*t))*sin(t+h), r*sin(11/3*t)], color=Colors[k+1], thickness=5), t=6*k/12*3.14..6*(k+1)/12*3.14, 0.01), k=0..11):

plots[display]( E, scaling=constrained);

Points A and B are taken arbitrarily, and point C is taken so that  vector C-A is 3 times longer than  B-A.

The set of points in space equidistant from two given planes is the union of two perpendicular planes, each of which bisects the corresponding dihedral angle between the given planes. In fact you have found this set. Its equation is E . To select from it each of the planes you can as follows:

with(geom3d):

plane(P,2*x+3*y-6*z+1=0, [x, y, z]):

plane(Q,x -4*y-8*z+1=0, [x,y,z]):

point(M, x, y, z):

E:=distance(M,P)= distance(M,Q);

simplify(E) assuming op([2,1], lhs(E))>0 and op([2,1], rhs(E))>0 ;  # equation of the first plane

simplify(E) assuming op([2,1], lhs(E))>0 and op([2,1], rhs(E))<0 ;  # equation of the second plane

If there is no  DirectSearch package, then the problem can be easily solved by exhaustive search . To reduce the amount of enumeration, we assume that one vertex locates at the origin:

restart;

L:=[]:

for x1 from -10 to 10 do

for y1 from -10 to 10 do

for z1 from -10 to 10 do

for x2 from -10 to 10 do

for y2 from -10 to 10 do

for z2 from -10 to 10 do

a:=x1^2+y1^2+z1^2:

if a>0 and a=x2^2+y2^2+z2^2 and a=(x2-x1)^2+(y2-y1)^2+(z2-z1)^2 then

L:=[op(L), [[0,0,0], [x1,y1,z1], [x2,y2,z2]]] fi;

od: od: od: od: od: od:

nops(L);  L[1..20];

We found 3936 solutions and displayed on the screen the first 20 decisions.

Let compare the solution obtained with the previous one by DirectSearch package. For correctness  we put one vertex at the origin (assume [х3, y3, z3]=[0, 0, 0]). Also we find the total number of solutions and display the first 20 ones:

restart;

interface(rtablesize=infinity):

L:=DirectSearch:-SolveEquations([x1^2+y1^2+z1^2 = (x2-x1)^2+(y2-y1)^2+(z2-z1)^2,

x1^2+y1^2+z1^2 = x2^2+y2^2+z2^2],  {x1^2+y1^2+z1^2 > 0, abs(x1) <= 10,

abs(y1) <= 10, abs(z1)<=10, abs(x2) <= 10, abs(y2) <= 10, abs(z2) <= 10}, assume = integer, AllSolutions):

LinearAlgebra[RowDimension](L);

seq(L[k,3], k=1..20);

We see that the program found only 98 solutions from 3936 solutions .

The procedure  N_points  finds N different random points in the plane eq, such that |x|<=m and |y|<=m. The number N should not be too large because number of cases is finite.

Npoints:=proc(N, eq, m)

local L, i, a, b;

L:=[];

for i while nops(L)<N do

a:=rand(-m..m)(); b:=rand(-m..m)();

L:={op(L),[a, b, solve(eval(eq, [x=a, y=b]), z)]};

od;

L;

end proc;

Example:

Npoints(2, 2*x + 2*y -z + 4=0, 5);

{[2,-1,6], [5,-5,4]}

The problem can easily be solved by a simple repetitive routine. To reduce the amount of exhaustive search we will search only the basic solutions. The integer solution will be called basic if the conditions 0 <=a<=b<=c<=d are hold. The program finds all the basic solutions (there are 726) and, for clarity, displays the first 100 solutions:

L:=[ ]:

for d from 0 to 100 do

for a from 0 to floor(d/sqrt(3)) do

for b from a to floor(sqrt(d^2-2*a^2)) do

for c from b to floor(sqrt(d^2-a^2-b^2)) do

if a^2+b^2+c^2=d^2 then L:=[op(L), [a, b, c, d]]; fi;

od: od: od: od:

nops(L);  L[1..100];

eq1:=4*x^n*y^m+(x+y)^n+x^3:

simplify(log[x+y](select(has, eq1, x+y)), symbolic);

                                        n

Try in Maple 16 interchange A and B. Maybe this will help:

plots[display](B, A, view=[2..7, 2..4], scaling=constrained, labels=[x, y]);