@sand15 I am sorry, I don't know how he did it. 

@sand15 My friend uses  asymptote https://asymptote.sourceforge.io/ and get the result

@mmcdara I don't say "with n = 14, The answer of Maple and your result are different." but I said There are differences between the two programs at the end of the questions when I compare two cases n = 14 and n = 15. 

@mmcdara with n = 14, The answer of Mathamatica and your result are different.

@nm How can I write solution to this mylist? 
We have, vector AB = (),  vector AC = ( ). Normal vector of the plane ABC is Crossproduct(AB, AC). The equation of equation ABC is 

@nm Thank you very much.

@vv Thank you very much. 

@Kitonum My code is slower than your. How can I reduce the time compile?

restart;
with(Student[MultivariateCalculus]);
A := [0, 0, 4];
B := [6, -2, 6];
M := [x, y, z];
d := Line(<1, -1, 2>, [4, -8, 4]);
f := Distance(M, d);
p := (Student[MultivariateCalculus]):-`.`(<x, y, z> - <0, 0, 4>, <x, y, z> - <6, -2, 6>);
extrema(f, {p = 0, x^2 + (y - 2)^2 + (z + 1)^2 = 29}, {x, y, z}, s);
s;

@acer Thank you very much

@acer I tried this code and get massage "Warning, solve may be ignoring assumptions on the input variables." 

with(Student:-MultivariateCalculus);
A := [0, 0, 0];
B := [1, 0, 0];
C := [1, 1, 0];
DD := [0, 2, 0];
assume(0 < h);
S := [0, 0, h];
d1 := Line(S, C);
d2 := Line(B, DD);
cos(Angle(d1, d2));
P := Plane(A, B, C);
temp := Angle(d1, P);
solve(temp = Pi/3);

I see this question https://mathematica.stackexchange.com/questions/296690/how-can-i-tell-mathematica-create-heronian-triangles-in-2d-like-maplepost.

@acer Thank you very much.

@acer I have mylist := [(x - 7)*(x^2 + (-2*m - 30)*x + m), (x - 8)*(x^2 + (-m + 8)*x + 2*m)]

I want to write this mylist in the form
L:=[[(x - 7)*(x^2 -2 (m 15)*x + m), (x - 8)*(x^2 - (m - 8)*x + 2*m)]]
How can I writr that?
I tried
sort~(map~(normal, mylist), x)

and get 

[(x - 7)*(x^2 - 2*m*x - 30*x + m), (x - 8)*(x^2 - m*x + 8*x + 2*m)]

@acer  I use Mathematica to forming those entries. Can you find them by Maple? 

With Maple, I tried
restart;
f := x -> x^4 + (a*m + b)*x^3 + (c*m + d)*x^2 + k*m + t;
n := 0;
for a from -4 to 1 do
for b from -7 to -1 do

for c to 3 do

for d to 5 do

for k from 3 to 7 do

for t from -10 to -5 do

for m from -5 to -1 do

if eval(diff(f(x), x), x = 2) = 0 and 0 < eval(diff(f(x), x, x), x = 2) and igcd(a, b) = 1 and igcd(c, d) = 1 and igcd(k, m) = 1 then n := n + 1; L[n] := [a, b, c, d, k, t, m]; end if; end do; end do; end do; end do; end do; end do;
end do;
L := convert(L, list);
nops(L);

@acer I like this two methods and I want this form. Thanks.