How to write program for [A(drager)^n, A^n].

A(drager)= a(t)(drager) exp(-iwt) and A= a(t) exp(iwt)

[KALYAN ACHARJYA]

Here Drager is an operator or Just Symbol?

[Samrat Banerjee]

Drager is a operator.

Answer 1

[David Goodmanson]

Hi Samrat,

You can't represent Adag and A by finite-dimensional matrices. Let's assume adag and a are NxN matrices that have no explicit time dependence (which can be implemented in other ways as below). Then

a*adag - adag*a =  eye(N)

But from the cyclic property of the trace,

trace(a*adag) - trace(adag*a) =  0;    trace(eye(N)) = N.

So finite-dimensional matrices don’t work. However, you can do the next best thing. Assume you have a state that is a linear combination of harmonic oscillator states with coefficients c_j(t), j = 0 ... m, and you know from your problem that m is never larger than 100 or whatever. Then you can operate on a vector of c_j(t) coefficients of total length comfortably greater than 100, with the extra entries above that value set to zero. That works.

Here each c_j(t) is a scalar coefficient for a particular time t.

Since Matlab is one-based, the coefficient c_0 is entry 1, the coefficient c_1 is entry2, etc. Everything is off by one.

The functions a and adag are determined by their effect on the states that correspond to each coefficient:

a|n> = sqrt(n)|n-1>               adag|n> = sqrt(n+1)|n+1> 

So adag slides the coefficients in the c vector to the right by one and ‘a’ slides them to the left.

I don’t want to mess up the opportunity for you to create your own functions, so here is an outline: Assume c is a row vector.

    For a:  find the length of the c vector, = m.  Create the vector sqrt(1:m-1) and multiply it
    term-by term with c(2:end)  to create the vector d.  The result is the vector [d 0]. 
    For adag:  find the length of the c vector, = m.  Create the vector sqrt(1:m-1) and multiply it
    term-by term with c(1:end-1)  to create the vector d.  The result is the vector [0 d].

If everything is done correctly, then two checks are:

c = [1 1 1 1 1 1 1 1 1 etc, 0 0 0 0 0 etc]
adag(a( c)) =? [0 1 2 3 4 5 6 7 8 9 etc, 0 0 0 0  etc]                the number operator  

c = [1.4  2.6 7.7 0.4 4.6 (could be anything), 0 0 0 0 etc];
a(adag( c)) - adag(a( c))  =?  c                                      anticommutation identity

My versions of these functions also contain a warning if the last entry of c is nonzero, which means that there is not enough zerofilling.