Tensor in QuantumOptics

Hello :slight_smile:
I have a program using QuantumOptics. It uses several bosonic modes. While writing the Hamiltonian, I wonder if I should explictly use the tensor product : I[1] tensor I[2] tensor a[3] for exemple, instead of only a[3] ? Here is my code below :slight_smile: :

    n = 3
    Nmax = 50
    var_t = 0.125
    var_U = 0.125
    var_mu = 0.125

    basis = [ FockBasis(Nmax) for i in 1:n]
    I = [ identityoperator(basis[i]) for i in 1:n]
    a = [ destroy(basis[i]) for i in 1:n]
    adag = [ create(basis[i]) for i in 1:n]
    _number = [ number(basis[i]) for i in 1:n]


    H = 0 * a[1]
    for i in 1:(n-1)
        #H += - var_t * ( b(i)' * b(i+1) + b(i) * b(i+1)' )
        H += -var_t * ( adag[i] * a[i+1] + adag[i+1] * a[i] )
    end

    for i in 1:n 
        #H += var_U/2.0 * (b(i)' * b(i)) * ( b(i)' * b(i)-1.0) - var_mu * b(i)'*b(i)
        H += var_U/2.0 * (_number[i] * (_number[i] - I[i] )) - var_mu * _number[i]
    end

is my code OK ? Thanks :slight_smile:

Finally I use the following :

    basis = [ FockBasis(Nmax) for i in 1:n]
    I = [ identityoperator(basis[i]) for i in 1:n]
    _I = tensor(I...)
    a = [ tensor( [ (i==j ? destroy(basis[i]) : I[j]) for j in 1:n ]... ) for i in 1:n] 
    adag = [ tensor( [ (i==j ? create(basis[i]) : I[j]) for j in 1:n ]... ) for i in 1:n] 
    _number = [ tensor( [ (i==j ? number(basis[i]) : I[j]) for j in 1:n ]... ) for i in 1:n]

    H = 0 * a[1]
    for i in 1:( n-1 )
        H += -var_t * ( adag[i] * a[i+1] + adag[i+1] * a[i] )
    end

    for i in 1:n 
        H += var_U/2.0 * ( _number[i] * (_number[i] - _I )) - var_mu * _number[i]
    end

    L = [ var_epsilon * a[1] , var_epsilon * adag[n]]

I think it is OK !

I think the embed function would make your life much easier. It is the main tool to use when you want to move from a subsystem Hamiltonian (e.g. for a given harmonic oscillator) to a Hamiltonian for the entire system (e.g. multiple modes).