Discover MakerZone

MATLAB and Simulink resources for Arduino, LEGO, and Raspberry Pi

Learn more

Discover what MATLAB® can do for your career.

Opportunities for recent engineering grads.

Apply Today

Thread Subject:
Difficult matrix construction

Subject: Difficult matrix construction

From: Dennis Belleter

Date: 7 Oct, 2010 16:20:06

Message: 1 of 11

I'm having some trouble designing a MatLab algorithm for the construction of a rather complicated matrix. The marix structure is as follows:

| B 0 0 .... 0 |
| A*B B 0 .... 0 |
| A^2*B A*B B .... 0 |
| : : : : : |
| A^(N-1)*B A^(N-2)*B A^(N-3)*B .... B |

where A and B are matrices as well.

I want to construct and use this matrix in an m-file for any arbitrary N.
I hope that someone has an algorithm for the construction of this matrix and can help me.
Thanks in advance.

Kind Regards

Dennis

Subject: Difficult matrix construction

From: Sean

Date: 7 Oct, 2010 16:35:04

Message: 2 of 11

"Dennis Belleter" <d.j.w.belleter@student.tue.nl> wrote in message <i8krvm$dm3$1@fred.mathworks.com>...
> I'm having some trouble designing a MatLab algorithm for the construction of a rather complicated matrix. The marix structure is as follows:
>
> | B 0 0 .... 0 |
> | A*B B 0 .... 0 |
> | A^2*B A*B B .... 0 |
> | : : : : : |
> | A^(N-1)*B A^(N-2)*B A^(N-3)*B .... B |
>
> where A and B are matrices as well.
>
> I want to construct and use this matrix in an m-file for any arbitrary N.
> I hope that someone has an algorithm for the construction of this matrix and can help me.

Like this?

n = 4; %Must equal size of both dimensions of A and B
fancy_mat = B.*tril(A.^toeplitz(0:n-1));

Subject: Difficult matrix construction

From: Dennis Belleter

Date: 7 Oct, 2010 17:14:05

Message: 3 of 11

"Sean " <sean.dewolski@nospamplease.umit.maine.edu> wrote in message <i8ksro$b5v$1@fred.mathworks.com>...
> "Dennis Belleter" <d.j.w.belleter@student.tue.nl> wrote in message <i8krvm$dm3$1@fred.mathworks.com>...
> > I'm having some trouble designing a MatLab algorithm for the construction of a rather complicated matrix. The marix structure is as follows:
> >
> > | B 0 0 .... 0 |
> > | A*B B 0 .... 0 |
> > | A^2*B A*B B .... 0 |
> > | : : : : : |
> > | A^(N-1)*B A^(N-2)*B A^(N-3)*B .... B |
> >
> > where A and B are matrices as well.
> >
> > I want to construct and use this matrix in an m-file for any arbitrary N.
> > I hope that someone has an algorithm for the construction of this matrix and can help me.
>
> Like this?
>
> n = 4; %Must equal size of both dimensions of A and B
> fancy_mat = B.*tril(A.^toeplitz(0:n-1));

Thanks Sean,

This works brilliant for scalar values, the only problem I have is dat A is a square matrix (2x2) and B is a collum (2x1).
Is it possible to do it in this manner for such a system?

Kind Regards,

Dennis

Subject: Difficult matrix construction

From: Sean

Date: 7 Oct, 2010 17:19:04

Message: 4 of 11

"Dennis Belleter" <d.j.w.belleter@student.tue.nl> wrote in message <i8kv4t$cr1$1@fred.mathworks.com>...
> "Sean " <sean.dewolski@nospamplease.umit.maine.edu> wrote in message <i8ksro$b5v$1@fred.mathworks.com>...
> > "Dennis Belleter" <d.j.w.belleter@student.tue.nl> wrote in message <i8krvm$dm3$1@fred.mathworks.com>...
> > > I'm having some trouble designing a MatLab algorithm for the construction of a rather complicated matrix. The marix structure is as follows:
> > >
> > > | B 0 0 .... 0 |
> > > | A*B B 0 .... 0 |
> > > | A^2*B A*B B .... 0 |
> > > | : : : : : |
> > > | A^(N-1)*B A^(N-2)*B A^(N-3)*B .... B |
> > >
> > > where A and B are matrices as well.
> > >
> > > I want to construct and use this matrix in an m-file for any arbitrary N.
> > > I hope that someone has an algorithm for the construction of this matrix and can help me.
> >
> > Like this?
> >
> > n = 4; %Must equal size of both dimensions of A and B
> > fancy_mat = B.*tril(A.^toeplitz(0:n-1));
>
> Thanks Sean,
>
> This works brilliant for scalar values, the only problem I have is dat A is a square matrix (2x2) and B is a collum (2x1).
> Is it possible to do it in this manner for such a system?
>
> Kind Regards,
>
> Dennis

It's definitely possible but you didn't illustrate what elements of A go where and same for B.

Subject: Difficult matrix construction

From: Jonathan

Date: 7 Oct, 2010 17:29:06

Message: 5 of 11

Dennis,

Perhaps the following will work. The A and B I use as a test make the result easy to check by hand.

********************************************
B = ones(4,2);
A = 10*eye(size(B,1));
n = 4;

result = zeros(size(B)*n);
for k = 0:n-1
    result = result + kron(diag(ones(1,n-k),-k),A^k*B);
end
disp(result);
********************************************

This will work with matrix B of any dimensions, any square matrix A compatible with B and any non-negative integer n.

Regards,
Jonathan


"Dennis Belleter" <d.j.w.belleter@student.tue.nl> wrote in message <i8krvm$dm3$1@fred.mathworks.com>...
> I'm having some trouble designing a MatLab algorithm for the construction of a rather complicated matrix. The marix structure is as follows:
>
> | B 0 0 .... 0 |
> | A*B B 0 .... 0 |
> | A^2*B A*B B .... 0 |
> | : : : : : |
> | A^(N-1)*B A^(N-2)*B A^(N-3)*B .... B |
>
> where A and B are matrices as well.
>
> I want to construct and use this matrix in an m-file for any arbitrary N.
> I hope that someone has an algorithm for the construction of this matrix and can help me.
> Thanks in advance.
>
> Kind Regards
>
> Dennis

Subject: Difficult matrix construction

From: Bruno Luong

Date: 7 Oct, 2010 17:29:07

Message: 6 of 11

% Data
A = rand(2);
B = rand(2,1);
N = 3;

% Engine
Ap = allmatpow(A,N-1);
C = arrayfun(@(p) Ap(:,:,p+1)*B, tril(toeplitz(0:N-1)),'Unif', false);
C(triu(true(N),1))={zeros(size(B))};
cell2mat(C)

function powers = allmatpow(A, pmax)
% function powers = allmatpow(A, pmax)
%
% compute all powers of exponent P = (0,1,...,PMAX) of matrix A.
% return the result in 2d array POWERS having size [length(X) x PMAX+1)]

powers= zeros([size(A) pmax+1]);
powers(:,:,1) = eye(size(A));
if pmax>0
    powers(:,:,2) = A;
    p2 = 1;
    r = 1;
    for k=2:pmax
        powers(:,:,k+1) = powers(:,:,p2+1)*powers(:,:,r+1);
        if (p2==r)
            p2 = 2*p2;
            r = 1;
        else
            r = r+1;
        end
    end % for-loop
end

end %% allmatpow

% Bruno

Subject: Difficult matrix construction

From: Roger Stafford

Date: 7 Oct, 2010 17:37:03

Message: 7 of 11

"Dennis Belleter" <d.j.w.belleter@student.tue.nl> wrote in message <i8krvm$dm3$1@fred.mathworks.com>...
> I'm having some trouble designing a MatLab algorithm for the construction of a rather complicated matrix. The marix structure is as follows:
>
> | B 0 0 .... 0 |
> | A*B B 0 .... 0 |
> | A^2*B A*B B .... 0 |
> | : : : : : |
> | A^(N-1)*B A^(N-2)*B A^(N-3)*B .... B |
>
> where A and B are matrices as well.
>
> I want to construct and use this matrix in an m-file for any arbitrary N.
> I hope that someone has an algorithm for the construction of this matrix and can help me.
> Thanks in advance.
>
> Kind Regards
>
> Dennis
- - - - - - -
  Probably the most efficient way is to use a for-loop. That way you can avoid needless repetition of matrix multiplication by A. Let the sizes of A and B be n x n and n x m, respectively. The size of the result, C, will be n*N by m*N.

 C = zeros(n*N,m*N);
 C(1:n,1:m) = B;
 for k = 2:N
  C((k-1)*n+1:k*n,1:m) = A * C((k-2)*n+1:(k-1)*n,1:m);
  C((k-1)*n+1:k*n,m+1:k*m) = C((k-2)*n+1:(k-1)*n,1:(k-1)*m);
 end

Roger Stafford

Subject: Difficult matrix construction

From: Dennis Belleter

Date: 7 Oct, 2010 17:38:04

Message: 8 of 11

"Sean " <sean.dewolski@nospamplease.umit.maine.edu> wrote in message <i8kve8$2f7$1@fred.mathworks.com>...
> "Dennis Belleter" <d.j.w.belleter@student.tue.nl> wrote in message <i8kv4t$cr1$1@fred.mathworks.com>...
> > "Sean " <sean.dewolski@nospamplease.umit.maine.edu> wrote in message <i8ksro$b5v$1@fred.mathworks.com>...
> > > "Dennis Belleter" <d.j.w.belleter@student.tue.nl> wrote in message <i8krvm$dm3$1@fred.mathworks.com>...
> > > > I'm having some trouble designing a MatLab algorithm for the construction of a rather complicated matrix. The marix structure is as follows:
> > > >
> > > > | B 0 0 .... 0 |
> > > > | A*B B 0 .... 0 |
> > > > | A^2*B A*B B .... 0 |
> > > > | : : : : : |
> > > > | A^(N-1)*B A^(N-2)*B A^(N-3)*B .... B |
> > > >
> > > > where A and B are matrices as well.
> > > >
> > > > I want to construct and use this matrix in an m-file for any arbitrary N.
> > > > I hope that someone has an algorithm for the construction of this matrix and can help me.
> > >
> > > Like this?
> > >
> > > n = 4; %Must equal size of both dimensions of A and B
> > > fancy_mat = B.*tril(A.^toeplitz(0:n-1));
> >
> > Thanks Sean,
> >
> > This works brilliant for scalar values, the only problem I have is dat A is a square matrix (2x2) and B is a collum (2x1).
> > Is it possible to do it in this manner for such a system?
> >
> > Kind Regards,
> >
> > Dennis
>
> It's definitely possible but you didn't illustrate what elements of A go where and same for B.

Oké I'll try to explain it more clearly

the first row of the Matrix so be as such:
[ B 0 0 ... 0] (where the row is n-elements long)

the second row is this:
[A*B B 0 .... 0] which is again n-elements long

then the third row is like this:
[A^2*B A*B B .... 0] (n-elements long)

and the n-th row is like this:
[A^(n-1)*B A^(n-2)*B A^(n-3)*B .... B] (n-elements long)

where the 0's are (2x1) zero collums and every A and B is the full matrix or collum so the rows are actually rows of (2x1) collums.

Subject: Difficult matrix construction

From: Dennis Belleter

Date: 7 Oct, 2010 17:45:04

Message: 9 of 11

"Jonathan " <bubbathekid@gmail.com> wrote in message <i8l012$b8g$1@fred.mathworks.com>...
> Dennis,
>
> Perhaps the following will work. The A and B I use as a test make the result easy to check by hand.
>
> ********************************************
> B = ones(4,2);
> A = 10*eye(size(B,1));
> n = 4;
>
> result = zeros(size(B)*n);
> for k = 0:n-1
> result = result + kron(diag(ones(1,n-k),-k),A^k*B);
> end
> disp(result);
> ********************************************
>
> This will work with matrix B of any dimensions, any square matrix A compatible with B and any non-negative integer n.
>
> Regards,
> Jonathan
>
>
> "Dennis Belleter" <d.j.w.belleter@student.tue.nl> wrote in message <i8krvm$dm3$1@fred.mathworks.com>...
> > I'm having some trouble designing a MatLab algorithm for the construction of a rather complicated matrix. The marix structure is as follows:
> >
> > | B 0 0 .... 0 |
> > | A*B B 0 .... 0 |
> > | A^2*B A*B B .... 0 |
> > | : : : : : |
> > | A^(N-1)*B A^(N-2)*B A^(N-3)*B .... B |
> >
> > where A and B are matrices as well.
> >
> > I want to construct and use this matrix in an m-file for any arbitrary N.
> > I hope that someone has an algorithm for the construction of this matrix and can help me.
> > Thanks in advance.
> >
> > Kind Regards
> >
> > Dennis

Thanks Jonathan that really did the job in an easy way!!

Subject: Difficult matrix construction

From: Matt J

Date: 7 Oct, 2010 21:44:03

Message: 10 of 11

"Dennis Belleter" <d.j.w.belleter@student.tue.nl> wrote in message <i8krvm$dm3$1@fred.mathworks.com>...
> I'm having some trouble designing a MatLab algorithm for the construction of a rather complicated matrix. The marix structure is as follows:
>
> | B 0 0 .... 0 |
> | A*B B 0 .... 0 |
> | A^2*B A*B B .... 0 |
> | : : : : : |
> | A^(N-1)*B A^(N-2)*B A^(N-3)*B .... B |
>
> where A and B are matrices as well.
>
========

I don't know if this is worth pointing out, since you haven't said how large N typically is, or how this matrix will be used.

However, the matrix is basically two Toeplitz matrices interleaved.
You could construct it as shown below, but if you will be using it for matrix-vector multiplications, it is basically equivalent to 2 convolutions, and could be more efficiently implemented both memory-wise and speed-wise using FFTs.

K=zeros(2,N);
tmp=B;
for ii=1:N,
    K(:,ii)=tmp ;
    tmp=A*tmp;
end


finalMatrix(1:2:2*N,1:N)=toeplitz(K(1,:),[K(1,1),zeros(1,N-1)]);
finalMatrix(2:2:2*N,1:N)=toeplitz(K(2,:),[K(2,1), zeros(1,N-1)]);

Subject: Difficult matrix construction

From: Matt J

Date: 8 Oct, 2010 16:33:03

Message: 11 of 11

"Matt J " <mattjacREMOVE@THISieee.spam> wrote in message <i8lev3$fdg$1@fred.mathworks.com>...
>
> However, the matrix is basically two Toeplitz matrices interleaved.
> You could construct it as shown below, but if you will be using it for matrix-vector multiplications, it is basically equivalent to 2 convolutions, and could be more efficiently implemented both memory-wise and speed-wise using FFTs.
=======

For kicks, I decided to implement the FFT approach using using this tool

http://www.mathworks.com/matlabcentral/fileexchange/26611-on-the-fly-definition-of-custom-matrix-objects

It allows a matrix multiplication to be done using normal matrix algebra syntax, but using FFTs internally. Notice the improvement in CPU time and the orders of magnitude saying in memory:


%fake data
N=2000;
A=rand(2);
B=rand(2,1);
x=rand(N);

%%convolution data
K=zeros(2,N);
tmp=B;
for ii=1:N,
    K(:,ii)=tmp ;
    tmp=A*tmp;
end



tic; %Approach 1 - compute the matrix and multiply normally

 M(2:2:2*N,1:N)=toeplitz(K(2,:),[K(2,1), zeros(1,N-1)]);
 M(1:2:2*N,1:N)=toeplitz(K(1,:),[K(1,1),zeros(1,N-1)]);
 
 z1=M*x;

toc;
%Elapsed time is 2.241007 seconds.


tic; %Approach 2 -- using FFTs

 Mobj=MatrixObj; %object equivalent to M
  Mobj.Ops.mtimes=@SpecialMult;

 z2=Mobj*x;

toc;
%Elapsed time is 1.252811 seconds.


>>whos M Mobj

  Name Size Bytes Class Attributes

  M 4000x2000 64000000 double
  Mobj 1x1 4412 MatrixObj


%%%%%%%%%%%%%
    function Y=SpecialMult(Q,X) %Nested function

        XX=fft(X,2*N);
        Y(2:2:2*N,1:N)=fftconv(K(2,:));
        Y(1:2:2*N,1:N)=fftconv(K(1,:));
        
       
        
            function CC=fftconv(b) %Nested function
        
                b=fft(b(:),2*N);
                CC = ifft( bsxfun(@times,b,XX));
                CC = real(CC(1:N,:));
        
            end
        
    end

Tags for this Thread

What are tags?

A tag is like a keyword or category label associated with each thread. Tags make it easier for you to find threads of interest.

Anyone can tag a thread. Tags are public and visible to everyone.

Contact us