 ring r=32003,x,dp;
matrix A[3][3] = 1,3,2,5,0,3,2,4,5; // define a matrix
print(A); // nice printing of small matrices
==> 1,3,2,
==> 5,0,3,
==> 2,4,5
A[2,3]; // matrix entry
==> 3
A[2,3] = A[2,3] + 1; // change entry
A[2,1..3] = 1,2,3; // change 2nd row
print(A);
==> 1,3,2,
==> 1,2,3,
==> 2,4,5
matrix E[3][3]; E = E + 1; // the unit matrix
matrix B =x*E  A;
print(B);
==> x1,3, 2,
==> 1, x2,3,
==> 2, 4, x5
// the same (but xA does not work):
B = A+x;
print(B);
==> x1,3, 2,
==> 1, x2,3,
==> 2, 4, x5
det(B); // the characteristic polynomial of A
==> x38x22x1
A*A*A  8 * A*A  2*A == E; // CayleyHamilton
==> 1
vector v =[x,1,x2];
A*v; // multiplication of matrix and vector
==> _[1,1]=2x2+x3
==> _[2,1]=3x2+x2
==> _[3,1]=5x2+2x4
matrix m[2][2]=1,2,3;
print(mtranspose(m));
==> 0,1,
==> 1,0
