@mkonto 

The vector allSolutions depends on the variables t[1], t[2], ...
You must substutute each of them.
eval(allSolutions,[t[1]=7,t[2]=5]);

eval(allSolutions, t=0) works because 0[x] evaluates to 0; for t=1 it does not work because  1[x] remains unevaluated.

If you want a random vector from allSolutions, use e.g.

eval(allSolutions,[seq(u='rand()', u=indets(allSolutions))]);

@Carl Love 

In wikipedia b=0. But after a translation by a vector u_0 (where Au_0 = b), one obtains Rouben's formula.

@mkonto 

I have changed the procedure to work with Gradient without changing the Gradient vector to free.

with(VectorCalculus);
SetCoordinates(cartesian[x,y,z]):
n1 := Gradient(a*x^2+y^2+z^2):
nf:=ConvertVector(n1,free):
vf:=Proj(<<1|2|3>,<2|4|6>>, <1,2>, nf);
v1:=Proj(<<1|2|3>,<2|4|6>>, <1,2>, n1);


V1 and vf will be displayed as

But if Proj is defined after loading VectorCalculus, they will will be displayed:

@Markiyan Hirnyk 

Thank you.

@Dave L 

I did not noticed there the interdiction to post (and discuss) the result of a Maple command such as

showstat(int)

which BTW does not contain any copyright mention. After all, Maplesoft offered this possibilty, it is not a "reverese engineering" made by a user. Otherwise it would be problematic to post the result of int(x^2,x).

@MariaL 

Just inspect the Cayley table. Foget about 90, it was probably a typo.

@krismalo 

fsolve without location finds only U__C <0 and the problem needs >0. That is why {U__C=0..1000, X=0..360} was added, but fsolve failed (not unusual!).

_Z1 and _B1 are arbitrary integer and binary constants contained in the general solution (see ?solve,details) .  They must be specified if a numeric solution is wanted.

@zia9206314 

p is your polynomial. I converted it in 1D math. In the original 2D worksheet the name Nu was somehow distorted and I replaced it by N [but you may revert to Nu if you want]:

p:=-(1/416179814400)*(N^4*k1^4*k2^7+7*N^4*k1^4*k2^6+21*N^4*k1^4*k2^5+35*N^4*k1^4*k2^4+35*N^4*k1^4*k2^3+364*N^3*k1^3*k2^5+21*N^4*k1^4*k2^2+1820*N^3*k1^3*k2^4+7*N^4*k1^4*k2+3640*N^3*k1^3*k2^3-23744*N^2*k1^2*k2^5+N^4*k1^4+3640*N^3*k1^3*k2^2-94976*N^2*k1^2*k2^4+1820*N^3*k1^3*k2-118160*N^2*k1^2*k2^3+364*N^3*k1^3-22064*N^2*k1^2*k2^2+49168*N^2*k1^2*k2-1435648*N*k1*k2^3+24304*N^2*k1^2-2871296*N*k1*k2^2-1216192*N*k1*k2+9625600*k2^3+219456*N*k1+9625600*k2^2-3990528*k2)*N^3*k1^3+(1/147456)*(N^2*k1^2*k2^4+4*N^2*k1^2*k2^3+6*N^2*k1^2*k2^2+4*N^2*k1^2*k2+N^2*k1^2+56*N*k1*k2^2+112*N*k1*k2+56*N*k1-640*k2^2-640*k2+144)*N^2*k1^2-(1/14745600)*(N^3*k1^3*k2^5+5*N^3*k1^3*k2^4+10*N^3*k1^3*k2^3+10*N^3*k1^3*k2^2+5*N^3*k1^3*k2+120*N^2*k1^2*k2^3+N^3*k1^3+360*N^2*k1^2*k2^2+360*N^2*k1^2*k2-2944*N*k1*k2^3+120*N^2*k1^2-5888*N*k1*k2^2-1520*N*k1*k2+1424*N*k1-13312*k2)*N^2*k1^2+(1/2123366400)*(N^4*k1^4*k2^6+6*N^4*k1^4*k2^5+15*N^4*k1^4*k2^4+20*N^4*k1^4*k2^3+15*N^4*k1^4*k2^2+220*N^3*k1^3*k2^4+6*N^4*k1^4*k2+880*N^3*k1^3*k2^3+N^4*k1^4+1320*N^3*k1^3*k2^2-9344*N^2*k1^2*k2^4+880*N^3*k1^3*k2-28032*N^2*k1^2*k2^3+220*N^3*k1^3-21008*N^2*k1^2*k2^2+4704*N^2*k1^2*k2+7024*N^2*k1^2-205312*N*k1*k2^2-205312*N*k1*k2+14400*N*k1+409600*k2^2)*N^2*k1^2-(1/4)*(k2+1)*N*k1+(1/64)*(N*k1*k2^2+2*N*k1*k2+N*k1+4)*N*k1-(1/2304)*(N^2*k1^2*k2^3+3*N^2*k1^2*k2^2+3*N^2*k1^2*k2+N^2*k1^2+20*N*k1*k2+20*N*k1-64*k2)*N*k1+1;

Then the command works.
(It is much better to use 1D math because you can see exactly the content)

collect-expand-Nu.mw

@Markiyan Hirnyk 

I inserted some comments. Note that the code is short and does not use tricks, so I think that more comments are not cecessary.

@dorna01 

I have spotted an error in M2[2,1]; you have omega(1-25000000000000/...)
instead of (probably)  omega*(1-25000000000000/...)
but this does not change essentially the problem (i.e. the equation is almost 0=0).

@Axel Vogt 

Ok, but what about symmetry?

If you want the coeffs expanded use:

sort(collect(p,N,expand));

I suggest to use the D1-math. Note that I used N instead of Nu (during the conversion D2-Math tio D1-Math).

You may of course write

% = 0

if you want.

@sunflower 

It is not that hard.
a)
k is the length of the largest number in the Pascal triangle. For example, if n=6, the largest number is binomial(6,3)=20, so k=2.

b) Each number in the triangle will be printed in a field of length 2k.

c) You shoud look at ?printf
E.g.
printf("   %4d%4d%4d\n", 5,10,200);


prints
      5  10 200
(notice the number of spaces; \n is for the newline).

Unfortunately both solutions may not work if det<=0; in John's, the symmetry is distroyed by RowOperation and in Carl's the determinant could be -2 (if n is even). So, RandomMatrix should be repeted until det>0. [the case det=0 is anyway "improbable"].
Another solution would be to start with an arbitrary nxn matrix B and compute const.B.Transpose(B) ; but in this case the result is not quite random, the matrix being nonnegative.

@Markiyan Hirnyk 

What general case? I have already said that all the cases are isomorphic, so one is enough.
Ok, I will stop here.