Items tagged with matrix

Dear All

I want to use subs command to replace products in expression with entries from matrix. When I use "subs" command it returns same result as original without any substitution. But when I use "algsubs", it works fine, but is complicated to apply if there are large number of substitutions.

Please see followings:


M := Matrix([[U[1], U[2], -epsilon*U[1]+U[3], U[4], U[5], U[6]], [U[1], U[2], -epsilon*U[2]+U[3], epsilon*U[1]+U[4], epsilon*U[2]+U[5], U[6]], [e^epsilon*U[1], e^epsilon*U[2], U[3], U[4], U[5], U[6]], [U[1], -epsilon*U[1]+U[2], U[3], U[4], -epsilon*U[4]+U[5], U[6]-epsilon*(U[3]+2*U[5])+epsilon^2*U[4], U[1], e^(-epsilon)*U[2], U[3], e^epsilon*U[4], U[5], e^(-epsilon)*U[6]], [U[1], e^(-epsilon)*U[2], U[3], e^epsilon*U[4], U[5], e^(-epsilon)*U[6]], [epsilon*U[2]+U[1], U[2], U[3], U[4]+epsilon*(U[3]+2*U[5])+epsilon^2*U[6], epsilon*U[6]+U[5], U[6]]])

M := Vector(4, {(1) = ` 6 x 12 `*Matrix, (2) = `Data Type: `*anything, (3) = `Storage: `*rectangular, (4) = `Order: `*Fortran_order})


M[5, 1]; 1; M[5, 3]








subs(T[5]*U[1] = M[5, 1], T[5]*U[2] = M[5, 2], T[5]*U[3] = M[5, 3], w*T[5]*U[1]*a[1]+w*T[5]*U[2]*a[2]+w*T[5]*U[3]*a[3]+w*T[5]*U[4]*a[4]+w*T[5]*U[5]*a[5]+w*T[5]*U[6]*a[6]+z*T[4]*U[1]*a[1]+z*T[4]*U[2]*a[2]+z*T[4]*U[3]*a[3]+z*T[4]*U[4]*a[4]+z*T[4]*U[5]*a[5]+z*T[4]*U[6]*a[6]+T[6]*U[1]*a[1]+T[6]*U[2]*a[2]+T[6]*U[3]*a[3]+T[6]*U[4]*a[4]+T[6]*U[5]*a[5]+T[6]*U[6]*a[6])



Why this command subs in not replacing  products  T[5]*U[1], T[5]*U[2], T[5]*U[3] from matrix (1)NULL?????

algsubs(T[5]*U[1] = M[5, 1], w*T[5]*U[1]*a[1]+w*T[5]*U[2]*a[2]+w*T[5]*U[3]*a[3]+w*T[5]*U[4]*a[4]+w*T[5]*U[5]*a[5]+w*T[5]*U[6]*a[6]+z*T[4]*U[1]*a[1]+z*T[4]*U[2]*a[2]+z*T[4]*U[3]*a[3]+z*T[4]*U[4]*a[4]+z*T[4]*U[5]*a[5]+z*T[4]*U[6]*a[6]+T[6]*U[1]*a[1]+T[6]*U[2]*a[2]+T[6]*U[3]*a[3]+T[6]*U[4]*a[4]+T[6]*U[5]*a[5]+T[6]*U[6]*a[6])



It worked !!!! But If I try "algsubs" like :

algsubs(T[5]*U[1] = M[5, 1], T[5]*U[2] = M[5, 2], T[5]*U[3] = M[5, 3], w*T[5]*U[1]*a[1]+w*T[5]*U[2]*a[2]+w*T[5]*U[3]*a[3]+w*T[5]*U[4]*a[4]+w*T[5]*U[5]*a[5]+w*T[5]*U[6]*a[6]+z*T[4]*U[1]*a[1]+z*T[4]*U[2]*a[2]+z*T[4]*U[3]*a[3]+z*T[4]*U[4]*a[4]+z*T[4]*U[5]*a[5]+z*T[4]*U[6]*a[6]+T[6]*U[1]*a[1]+T[6]*U[2]*a[2]+T[6]*U[3]*a[3]+T[6]*U[4]*a[4]+T[6]*U[5]*a[5]+T[6]*U[6]*a[6])

Error, invalid input: algsubs expects its 3rd argument, x, to be of type {name, function, list({name, function}), set({name, function})}, but received T[5]*U[3] = U[3]


See it is giving error !!! but still I can use "algsubs" as follow:


expand(algsubs(T[6]*U[6] = M[6, 6], expand(algsubs(T[6]*U[5] = M[6, 5], expand(algsubs(T[6]*U[4] = M[6, 4], expand(algsubs(T[6]*U[3] = M[6, 3], expand(algsubs(T[6]*U[2] = M[6, 2], expand(algsubs(T[6]*U[1] = M[6, 1], expand(algsubs(T[5]*U[6] = M[5, 6], expand(algsubs(T[5]*U[5] = M[5, 5], expand(algsubs(T[5]*U[4] = M[5, 4], expand(algsubs(T[5]*U[3] = M[5, 3], expand(algsubs(T[5]*U[2] = M[5, 2], expand(algsubs(T[5]*U[1] = M[5, 1], expand(algsubs(T[4]*U[6] = M[4, 6], expand(algsubs(T[4]*U[5] = M[4, 5], expand(algsubs(T[4]*U[4] = M[4, 4], expand(algsubs(T[4]*U[3] = M[4, 3], expand(algsubs(T[4]*U[2] = M[4, 2], expand(algsubs(T[4]*U[1] = M[4, 1], expand(algsubs(T[3]*U[6] = M[3, 6], expand(algsubs(T[3]*U[5] = M[3, 5], expand(algsubs(T[3]*U[4] = M[3, 4], expand(algsubs(T[3]*U[3] = M[3, 3], expand(algsubs(T[3]*U[2] = M[3, 2], expand(algsubs(T[3]*U[1] = M[3, 1], expand(algsubs(T[2]*U[6] = M[2, 6], expand(algsubs(T[2]*U[5] = M[2, 5], expand(algsubs(T[2]*U[4] = M[2, 4], expand(algsubs(T[2]*U[3] = M[2, 3], expand(algsubs(T[2]*U[3] = M[2, 3], expand(algsubs(T[2]*U[2] = M[2, 2], expand(algsubs(T[2]*U[1] = M[2, 1], expand(algsubs(T[1]*U[6] = M[1, 6], expand(algsubs(T[1]*U[5] = M[1, 5], expand(algsubs(T[1]*U[4] = M[1, 4], expand(algsubs(T[1]*U[3] = M[1, 3], expand(algsubs(T[1]*U[2] = M[1, 2], expand(algsubs(T[1]*U[1] = M[1, 1], w*T[5]*U[1]*a[1]+w*T[5]*U[2]*a[2]+w*T[5]*U[3]*a[3]+w*T[5]*U[4]*a[4]+w*T[5]*U[5]*a[5]+w*T[5]*U[6]*a[6]+z*T[4]*U[1]*a[1]+z*T[4]*U[2]*a[2]+z*T[4]*U[3]*a[3]+z*T[4]*U[4]*a[4]+z*T[4]*U[5]*a[5]+z*T[4]*U[6]*a[6]+T[6]*U[1]*a[1]+T[6]*U[2]*a[2]+T[6]*U[3]*a[3]+T[6]*U[4]*a[4]+T[6]*U[5]*a[5]+T[6]*U[6]*a[6]))))))))))))))))))))))))))))))))))))))))))))))))))))))))))))))))))))))))))



But this procedure is complex, I want to use "subs" only. I know it complexity of terms in (3) might me creating problem.

collect(a[5]*epsilon*U[6]+a[1]*epsilon*U[2]+U[6]*a[6]+U[2]*a[2]+U[3]*a[3]+U[4]*a[4]+U[5]*a[5]+z*a[1]*U[1]+epsilon^2*U[6]*a[4]+epsilon*U[3]*a[4]+2*epsilon*U[5]*a[4]+U[6]*w*a[6]/e^epsilon+U[4]*e^epsilon*w*a[4]+w*a[1]*U[1]+U[2]*w*a[2]/e^epsilon+epsilon^2*z*U[4]*a[6]-epsilon*z*U[3]*a[6]-2*epsilon*z*U[5]*a[6]+w*U[3]*a[3]+w*U[5]*a[5]-z*a[5]*epsilon*U[4]-z*a[2]*epsilon*U[1]+z*U[2]*a[2]+z*U[3]*a[3]+z*U[4]*a[4]+z*U[5]*a[5]+z*U[6]*a[6]+a[1]*U[1], [U[1], U[2], U[3], U[4], U[5], U[6]])







I thought I would share some code for computing sparse matrix products in Maple.  For floating point matrices this is done quickly, but for algebraic datatypes there is a performance problem with the builtin routines, as noted in

Download spm.txt

The code is fairly straightforward in that it uses op(1,A) to extract the dimensions and op(2,A) to extract the non-zero elements of a Matrix or Vector, and then loops over those elements.  I included a sparse map function for cases where you want to map a function (like expand) over non-zero elements only.

# sparse matrix vector product
spmv := proc(A::Matrix,V::Vector)
local m,n,Ae,Ve,Vi,R,e;
n, m := op(1,A);
if op(1,V) <> m then error "incompatible dimensions"; end if;
Ae := op(2,A);
Ve := op(2,V);
Vi := map2(op,1,Ve);
R := Vector(n, storage=sparse);
for e in Ae do
n, m := op(1,e);
if member(m, Vi) then R[n] := R[n] + A[n,m]*V[m]; end if;
end do;
return R;
end proc:
# sparse matrix product
spmm := proc(A::Matrix, B::Matrix)
local m,n,Ae,Be,Bi,R,l,e,i;
n, m := op(1,A);
i, l := op(1,B);
if i <> m then error "incompatible dimensions"; end if;
Ae := op(2,A);
Be := op(2,B);
R := Matrix(n,l,storage=sparse);
for i from 1 to l do
Bi, Be := selectremove(type, Be, (anything,i)=anything);
Bi := map2(op,[1,1],Bi);
for e in Ae do
n, m := op(1,e);
if member(m, Bi) then R[n,i] := R[n,i] + A[n,m]*B[m,i]; end if;
end do;
end do;
return R;
end proc:
# sparse map
smap := proc(f, A::{Matrix,Vector})
local B, Ae, e;
if A::Vector then
B := Vector(op(1,A),storage=sparse):
B := Matrix(op(1,A),storage=sparse):
end if;
Ae := op(2,A);
for e in Ae do
B[op(1,e)] := f(op(2,e),args[3..nargs]);
end do;
return B;
end proc:

As for how it performs, here is a demo inspired by the original post.

n := 674;
k := 6;
A := Matrix(n,n,storage=sparse):
for i to n do
  for j to k do
    A[i,irem(rand(),n)+1] := randpoly(x):
  end do:
end do:
V := Vector(n):
for i to k do
  V[irem(rand(),n)+1] := randpoly(x):
end do:
C := CodeTools:-Usage( spmv(A,V) ):  # 7ms, 25x faster
CodeTools:-Usage( A.V ):  # 174 ms
B := Matrix(n,n,storage=sparse):
for i to n do
  for j to k do
    B[i,irem(rand(),n)+1] := randpoly(x):
  end do:
end do:
C := CodeTools:-Usage( spmm(A,B) ):  # 2.74 sec, 50x faster
CodeTools:-Usage( A.B ):  # 2.44 min
# expand and collect like terms
C := CodeTools:-Usage( smap(expand, C) ):

I want to set the determinant of the coefficient matrix equal to zero and then solving for the roots. But I could not achieve it via Maple. Can you help me please? 

You can reach two examples in the following file.  Yeni_Microsoft_Word_Belgesi_(2).docx


Besides. how can i compute the following transcental equations via maple 










Sorry for disturbing you. I am wondering if there is an easier approach in Maple that could convert a system of second order differential equations into matrix form. Of course, we could do it by hand easily if the degrees of freedom is small. I would like to know if we could use Maple to do so. 

Here is an example with 6 degrees of freedom: the variables are u, v, w, alpha, beta and gamma. And, this is a uncoupled system.

Vector(6, {(1) = 2*R^2*(diff(w(t), t))*Pi*Omega*h*rho+R^2*(diff(u(t), t, t))*Pi*h*rho-R^2*u(t)*Pi*Omega^2*h*rho = 0, (2) = R^2*(diff(v(t), t, t))*Pi*h*rho = 0, (3) = -2*R^2*(diff(u(t), t))*Pi*Omega*h*rho+R^2*(diff(w(t), t, t))*Pi*h*rho-R^2*w(t)*Pi*Omega^2*h*rho = 0, (4) = (1/4)*R^4*Pi*(diff(alpha(t), t, t))*h*rho+(1/12)*R^2*Pi*(diff(alpha(t), t, t))*h^3*rho+(1/6)*R^2*Pi*(diff(gamma(t), t))*Omega*h^3*rho-(1/12)*R^2*Pi*alpha(t)*Omega^2*h^3*rho = 0, (5) = (1/2)*R^4*Pi*(diff(beta(t), t, t))*h*rho-(1/2)*R^4*Pi*beta(t)*Omega^2*h*rho = 0, (6) = (1/4)*R^4*Pi*(diff(gamma(t), t, t))*h*rho+(1/12)*R^2*Pi*(diff(gamma(t), t, t))*h^3*rho-(1/6)*R^2*Pi*(diff(alpha(t), t))*Omega*h^3*rho-(1/12)*R^2*Pi*gamma(t)*Omega^2*h^3*rho = 0});

The objective is to reform it into matrix form : M*diff(X(t), t, t)+C*diff(X(t), t)+K*X(t)=F.

Thank you in advance for taking a look. 

Eight matrices inside the list J:

J := [Matrix([[1, 0], [0, 1]]), Matrix([[0, -I], [-I, 0]]), Matrix([[0, -1], [1, 0]]), Matrix([[-I, 0], [0, I]]), Matrix([[-1, 0], [0, -1]]), Matrix([[0, I], [I, 0]]), Matrix([[0, 1], [-1, 0]]), Matrix([[I, 0], [0, -I]])];


Function member identifies J[2] as a member of J and returns its position in j:

member(J[2], J, 'j'); j;


Matrix multiplication inside a loop does not have a matrix type:

for i to 1 do for j to 2 do a := J[i].J[j]; member(a, J, 'k'); print(i, j, k, a, whattype(a)) end do end do;

Has anyone any ideia of what is going on?


Hello everyone.  Let n be a positive number; I seek how to write the matrix Z\in M_{2^n,n} that is defined as follows:

for every i,j, Z_{ij}\in \{-1,1\}.

When I consider the 2^n rows of Z, I want to find ALL the possible sequences of length n with the entries +-1. A toy example, when n=2, is Z=\begin{pmatrix}1&1\\1&-1\\-1&1\\-1&-1\end{pmatrix}.

Thanks in advance;



Anyone can solve this??



Many thanks!

Hi guys,
I want to import symbolic matrix from matlab to Maple, How I can do that ? 

I want to solve a system of equations using f-solve (two unknowns) and exporting the solutions to a matrix where the solutions are in seperate columns. How do I do this?

I have tried:

for i from 1 to 937 do


end if

end do


But this returns the solutions for X and Y in the same column. Also, for the values that are not possible to solve, it returns the entire expression instead of e.g. 0 or "undefined".

Thank you.

 from determinant's polynomial?                                                                                                       

I tried to make a procedure that would find the determinant of any 3x3 matrix but I keep getting unterminated procedure what should I change??

    local  i::integer,  j::integer,  x::integer,  y::integer,  A::integer,  B::integer,  indice::integer,  S::integer:
      A:=0:  B:=0:  S:=0: i:= 1:  j:=1:  x:=1:  y := 1:  indice:=1:

for x from 1 to 3 do
   indice := 1:  
   for i  from 1 to  3 do  
       if x = i then
           end do:
        else if x <>i then
           for j from 1 to 3 do
                  if y = j then  end do:  
                  else if indice = 1 then  A := matA(i,j):  indice := 2:  
                     else if indice = 2 then  B:=matA(i,j):  indice :=3:    
                       else if indice = 3 then  B :=B * matA(i,j):  indice :=4:    
                         else if indice = 4 then  A := A * matA(i,j):      
                             if x = 2 then  S := S + (-1 * matA(x,y)) * (A-B):
                               else  S := S + matA(x,y) * (A-B):   end if:   end if:
 end do:
 end if:
 end do:
 end do:
 end proc:

Hi, I have 2 Questions about programming in maple. I will be thankful if you help me as soon as possible.

First; How can I display a n*n HilbertMatrix?

Second: I wanna make a 2*2 matrix which its transpose is egual to its inverse, How can I do that by helping reflectionmatrix?

(I'm an amateur programmer, I use maple 13 on my pc)

I faced a very large eigenproblem during my research. The square matrix under consideration is of size more than 2^30 times 2^30. I have tried to deal with this problem by the QR algorithm with double implicit shift (more precisely, the Francis double step QR algorithm). I'm a very beginner of programming, but I tried as follows:


A := Matrix([[7, 3, 4, -11, -9, -2], [-6, 4, -5, 7, 1, 12], [-1, -9, 2, 2, 9, 1], [-8, 0, -1, 5, 0, 8], [-4, 3, -5, 7, 2, 10], [6, 1, 4, -11, -7, -1]]):
H := HessenbergForm(A):
for p while p>2 do: 
for k from 0 to p-3 do:  
if k<3 then z:=H(k+4,k+1):   
end if: 
if abs(H(p,q))<10^(-20)*(abs(H(q,q))+abs(H(p,p))) then    H(p,q):=0: p:=p-1:q=p-1:  
elif abs(H(p-1,q-1))<10^(-20)*(abs(H(q-1,q-1))+abs(H(q,q))) then    H(p-1,q-1):=0: p:=p-2:q:=p-1:  
end if:  od:

It seemed that replacing 0 in a Hessenberg matrix by a non-zero element is not allowed. How can I remedy this?

Plus, can anyone tell me the problem of the above thing(it's not really a programming...;( ), please?

I would also appreciate it if someone let me know a better idea for a huge eigenproblem.

Thanks in advance.

in LinearAlgebra Eigenvectors calculation.

Maple 2015 Error



So the above output startled me.  I have used the Maple Linear Algebra Eigenvalues, Eigenvectors commands many times with no problem.  Can any one explain to me what is going on.  The program correctly calculates the eigenvalues for the matrix which are all distinct for a real symmetric matrix, and thus should have three distinct non-zero eigenvectors, yet the eigenvectore command returns only zeros for the eigenvectors.  I calculated an eigenvector by hand corresponding to the eigenvalue of 1 and obtained (1, -sqrt(2)/sqrt(3), -1/sqrt(3).


So this is either a serious bug or I am going completely insane. 

Sorry for the uninformative title. I've never used Maple, but I'm willing to buy a student license and learn it. But before spending too much effort and money I need to know if it suits my needs.

Basically what I need to do is:

1) I have a positive definite symmetric matrix of size nxn, where n can range from 2 to inf. I don't know the elements, except the fact that the diagonal has ones everywhere. All I know is that the elements out of the diagonal are in the range [0,1)

2) I have to compute the lower triangular cholesky decomposition of this matrix, lets call it L.

3) I need to subtract from each element of L the mean of the elements in the respective column. Lets call this matrix L*

4) Then I need to evaluate another nxn matrix computed from the elements of L* following a simple pattern.

5) Finally I need to find the eigenvalues of this last matrix.

What I would ideally want is to get a symbolic representation of the n eigenvalues as symbolic functions of the (unknown) elements of the matrix at point 1.

I can drop the assumption of n being unknown, i.e. fix n=3 and get the 3 functions that, after replacing the right values, give me the eigenvalues, then fix n=4 and get 4 functions, etc.

Is this possible to do in maple?

Thank you

5 6 7 8 9 10 11 Last Page 7 of 46