daniel_w_

0 Reputation

One Badge

12 years, 40 days

MaplePrimes Activity


These are replies submitted by daniel_w_

@acer

That would explain why it does not work...

I found an 'Inverse' command for inert matrix inverses but I have not even been calling it correctly. The calling sequence is: Inverse(A) mod n. A being a matrix, n being an integer - the modulus.

I was hoping to get a qualitative insight into the form of the eigenvalues, or at least the real parts of them as these are of bigger importance to me. But from what I can see so far, I might have to settle for doing this numerically...

Thanks!

 

@Carl Love 

Hi, the only package I am using here is the 'Linear Algebra' package.

I have attached my Maple script, it is quite short.

My main aim is to obtain the analytical form of the eigenvalues of the matrix exponential of 'Lss', howver, when I first executed the 'Matrix Exponential' command, it took over 70 000 seconds at which stage I got fed up and decided to do it step by step to see what is causing the problem.

I am trying to diagonalize the 4 by 4 matrix 'Lss' first, for which purpose the matrix of eigenvectors must be obtained as well as its inverse. What I find is that the command 'Inverse' works almost instantaneously. 'Matrix Inverse' takes much longer and I really don't know why. The other problem is that right at the end of the script I test whether the inverse is a true inverse and it seems doubtful as I expect an identity matrix as a result. That is why the thought that 'Inverse' is an approximation popped into my head. I think the main problem is that the eigenvector matrix terms are very long but the main problem is probably that I am not very experienced with Maple and most probably am doing something stupid.

 

Thank you,

Daniel

restart

with(LinearAlgebra):

unprotect(GAMMA):

assume(GAMMA[1] >= 0, GAMMA[2] >= 0, GAMMA[3] >= 0, omega >= 0, Delta >= 0):

interface(showassumed = 0):

tran := proc (x) options operator, arrow; HermitianTranspose(x) end proc:

kron := proc (x, y) options operator, arrow; KroneckerProduct(x, y) end proc:

mult := proc (x, y) options operator, arrow; Multiply(x, y) end proc:

Id := Matrix(2, [1, 0, 0, 1]):#Identity``

Lm := Matrix(2, [0, 0, sqrt(GAMMA[2]), 0]):

Lp := Matrix(2, [0, sqrt(GAMMA[1]), 0, 0]):

Lz := Matrix(2, [(1/2)*sqrt(GAMMA[3]), 0, 0, -(1/2)*sqrt(GAMMA[3])]):

K := -1/2*(mult(tran(Lp), Lp)+mult(tran(Lm), Lm)+mult(tran(Lz), Lz))

L1 := kron(Id, K)+kron(K, Id)

L1 := Matrix(4, 4, {(1, 1) = -GAMMA[2]-(1/4)*GAMMA[3], (1, 2) = 0, (1, 3) = 0, (1, 4) = 0, (2, 1) = 0, (2, 2) = -(1/2)*GAMMA[1]-(1/4)*GAMMA[3]-(1/2)*GAMMA[2], (2, 3) = 0, (2, 4) = 0, (3, 1) = 0, (3, 2) = 0, (3, 3) = -(1/2)*GAMMA[1]-(1/4)*GAMMA[3]-(1/2)*GAMMA[2], (3, 4) = 0, (4, 1) = 0, (4, 2) = 0, (4, 3) = 0, (4, 4) = -GAMMA[1]-(1/4)*GAMMA[3]})

(1)

H := Matrix(2, [(1/2)*Delta, (1/2)*omega, (1/2)*omega, -(1/2)*Delta]):

# For testing purposesNULL

L2 := -i*(kron(Id, H)-kron(H, Id))

L2 := Matrix(4, 4, {(1, 1) = 0, (1, 2) = -(1/2)*i*omega, (1, 3) = (1/2)*i*omega, (1, 4) = 0, (2, 1) = -(1/2)*i*omega, (2, 2) = i*Delta, (2, 3) = 0, (2, 4) = (1/2)*i*omega, (3, 1) = (1/2)*i*omega, (3, 2) = 0, (3, 3) = -i*Delta, (3, 4) = -(1/2)*i*omega, (4, 1) = 0, (4, 2) = (1/2)*i*omega, (4, 3) = -(1/2)*i*omega, (4, 4) = 0})

(2)

L3 := kron(Lz, Lz)+kron(Lp, Lp)+kron(Lm, Lm)

L3 := Matrix(4, 4, {(1, 1) = (1/4)*GAMMA[3], (1, 2) = 0, (1, 3) = 0, (1, 4) = GAMMA[1], (2, 1) = 0, (2, 2) = -(1/4)*GAMMA[3], (2, 3) = 0, (2, 4) = 0, (3, 1) = 0, (3, 2) = 0, (3, 3) = -(1/4)*GAMMA[3], (3, 4) = 0, (4, 1) = GAMMA[2], (4, 2) = 0, (4, 3) = 0, (4, 4) = (1/4)*GAMMA[3]})

(3)

Lss := `<,>`(L1+L2+L3)#Full Lindblad SuperoperatorNULL

Lss := Matrix(4, 4, {(1, 1) = -GAMMA[2], (1, 2) = -(1/2)*i*omega, (1, 3) = (1/2)*i*omega, (1, 4) = GAMMA[1], (2, 1) = -(1/2)*i*omega, (2, 2) = -(1/2)*GAMMA[1]-(1/2)*GAMMA[3]-(1/2)*GAMMA[2]+i*Delta, (2, 3) = 0, (2, 4) = (1/2)*i*omega, (3, 1) = (1/2)*i*omega, (3, 2) = 0, (3, 3) = -(1/2)*GAMMA[1]-(1/2)*GAMMA[3]-(1/2)*GAMMA[2]-i*Delta, (3, 4) = -(1/2)*i*omega, (4, 1) = GAMMA[2], (4, 2) = (1/2)*i*omega, (4, 3) = -(1/2)*i*omega, (4, 4) = -GAMMA[1]})

(4)

eigval, eigvec := Eigenvectors(Lss):

eigvec[4]#NULLaccess specific eigenvectorsNULL

Vector[row](4, {(1) = 1, (2) = 1, (3) = 1, (4) = 1})

(5)

eigvmatrix := `<,>`(eigvec[1], eigvec[2], eigvec[3], eigvec[4])

P := Transpose(eigvmatrix)

Pinv := Inverse(P)#eigenvector matrix in the appropriate form

simplify(mult(Pinv, P))#check that the inverse 'works'

`[Length of output exceeds limit of 1000000]`

(6)

``


Download Electron_Lindblad_su.mwElectron_Lindblad_su.mw

Page 1 of 1