6 views (last 30 days)

Show older comments

I want to integrate of an exponential of a matrix ( E= exp (Mmat.*(t0-t))). First of all ‘expm’ built-in function in MATLAB has complexity of order N^3 ; while my matrix is a circulant matrix, so I can use a theorem instead that makes complexity of the order of NlogN;

The algorithm is :

( in this algorithm, E=exp(A) which A is a matrix )

By this algorithm I can easily calculate the exp of a matrix but my problem now is the calculation of the integral of this function I mean

load filem.mat

t0=1e-4;

Mmatt=@(t)Mmat.*(t0-t);

col1=Mmatt(:,1);

FFtcol1=fft(col1);

expFF=exp(FFtcol1);

expMt=ifft(expFF);% the first coulmn of exp (mt)

expMmat(:,1)=expMt;

for x=2:322

expMmat(:,x) = circshift(expMDelt,x-1,1);

end

fun=expMmat*Dmat;

f=integral(fun,0,1e-4,'ArrayValued',true);

I don’t know how can I define it as a function. I can use handle function or function defining in a separate file code but how can I do it. Also I can use integral built in finction in matlab.

Also I think I have an another option and calculate the amount of the function in each time and use the trapz function for integration but I think Integral is more suitable for this situation. I attach my code and the data. I really appreciate any help

thanks in advance

David Goodmanson
on 20 Jul 2021

Edited: David Goodmanson
on 20 Jul 2021

Hi raha,

I think you can go with the following technique, which does integration, differentiation etc. at the eigenvalue level. I was casual about it and created the circulant matrices by concatentation which of course would have to be improved for large matrices.

t = 4;

m1 = 2*rand(5,1)-1

M = circ(m1)

expM = expm(M*t)

% calculate N = expm(M*t)

c1 = M(:,1);

lam = fft(c1);

lamfun = exp(lam*t);

i1 = ifft(lamfun);

N = circ(i1)

max(abs(N-expM),[],'all') % check should be small

% calculate D = d/dt expm(M*t)

c1 = M(:,1);

lam = fft(c1);

lamfun = (exp(lam*t)).*lam;

i1 = ifft(lamfun);

D = circ(i1)

max(abs(D - M*expM),[],'all')

max(abs(D - expM*M),[],'all')

% calculate I = Int expm(M*t) dt (indefinite integral)

c1 = M(:,1);

lam = fft(c1);

lamfun = (exp(lam*t))./lam;

i1 = ifft(lamfun);

I = circ(i1)

max(abs(I*M -expM),[],'all')

max(abs(M*I -expM),[],'all')

function M = circ(m1)

% create circulant matrix by circular shift of columns

n = size(m1,1);

M = m1;

for k = 1:n-1

M = [M circshift(m1,k,1)];

end

end

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!