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]);

A:=plot(f(x), x=2..7, color=white, filled=true):

B:=plot(g(x), x=2..7, color=cyan, filled=true):

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

The following method works:

eq2:=subsop(seq([1,k]=0, k=2..4), [1,5]=q11, eq1);

If you want to build a pyramid of arbitrary, it is convenient to use a special procedure . Formal arguments are: X - a list of names for the vertices of the pyramid, Y - a list of the coordinates. After specifying the pyramid for it's building  plots[display] command should be used. You can add other objects.

Code of Procedure:

Pyramid := proc(X::(list(symbol)), Y::(list(list)))

local L, L1, T, T1, n;

n:=nops(X);

L:=plottools[polygon]([seq(Y[i],  i=1..n-1)], style=line, color=blue, thickness=2);

L1:=seq(plottools[line](Y[i],Y[n], color=blue, thickness=2), i=1..n-1);

T:=plots[textplot3d]([seq([op(Y[i]), X[i]], i=1..n-1)], align={left,below}, color=black, font=[TIMES,ROMAN,18]);

T1:=plots[textplot3d]([op(Y[n]),X[n]], align=above, color=black, font=[TIMES,ROMAN,18]);

L, L1, T, T1;

end proc:

Example of constructing a pyramid with a base:

plots[display](Pyramid([A,B,C,D,S],[[0,0,0],[2,0,0],[2,3,0],[0,3,0],[0,0,4]]),plot3d(0,x=-1..3,y=-1..4,style=surface,color=pink), scaling=constrained, orientation=[35, 65]);

Use subs command instead of algsubs command.

print(`Answer:`); print(sqrt(5));

I do not know how to write it in a single line.

L:=combinat[permute]([0$8, 1$8], 8):

V:=[seq(v||i, i=1..8)]:

M:=[seq(add(L[k,i]*V[i], i=1..8), k=1..2^8)];

M is the list of all required sums.