Convert co-variance to Hamiltonian matrix

hey,
I was trying to reconstruct the Hamiltonian from a covariance matrix but I always get nan values. this function implements the conversion. I checked and I realized that it’s the “atanh” that causes the problem. is there a way to overcome this issue??
I’m using the latest version of strawberryfields from the master

here is a simple script to produce the same error:

prog = sf.Program(1)
amp = .7
np.random.seed(5)
x = np.random.rand(1)

with prog.context as q:

    sf.ops.DisplacedSqueezed(r=amp,p=x) | q[0]

eng = sf.Engine('gaussian')

result = eng.run(prog)    

cov = result.state.cov()
print(sf.decompositions.covmat_to_hamil(cov))

Thanks @kareem_essafty for pointing this out. I’m not sure if this is a bug or not. I’ll consult with the team and we’ll get back to you when we understand more what’s causing this issue.

hey, i have read the paper that is pointed out in the docs and I believe by changing only this part of this line

H = (1j * omega @ (v @ np.diag(np.arctanh(1.0/l.real)) @ np.linalg.inv(v))).real

to

H = (1j * omega @ (v @ np.diag(np.arctanh(1.0/l)) @ np.linalg.inv(v))).real

it will lead to overcoming the issue of nan values

is that right?

Hi @kareem_essafty — Thanks for the question!
I believe sf is giving you the correct output, in that a Gibbs state, i.e., a state of the form \exp(-\beta H) is only pure at infinite inverse temperature \beta = 1/(k_B T) where k_B is the Boltzmann constant and T is the temperature.
A quick solution to this problem is to start with a thermal state with a small mean number of photons and then extrapolate to the zero temperature limit.

1 Like

hey @Nicolas_Quesada, I tried the following:

prog = sf.Program(2)
np.random.seed(5)
x = np.random.rand(1)
print(x[0])
with prog.context as q:

    sf.ops.Thermal(n=1e-015) | q[0]
    sf.ops.Thermal(n=1e-015) | q[1]
    #sf.ops.S2gate(r=0.7,phi=0) | q
    #sf.ops.DisplacedSqueezed(alpha=1.,r=.7,p=x[0]) | q[0]
    #sf.ops.DisplacedSqueezed(alpha=2,r=.7,p=x[0]+np.pi) |q[1]
    sf.ops.Sgate(r=0.7,phi=x[0]) | q[0]
    sf.ops.Sgate(r=0.7,phi=x[0]+np.pi/4) | q[1]
    sf.ops.Dgate(7.2) | q[0]
    sf.ops.Dgate(-0.6-0.3j) | q[1]

eng = sf.Engine('gaussian')

result = eng.run(prog)    

cov = result.state.cov()
print(cov)
print(" ")
print(sf.decompositions.covmat_to_hamil(cov)) 

the result didn’t include any nan and was perfect except for the fact that the hamiltonian was not symmetric.
besides that, if i used a displacedsqueezed operation instead of the S and D gates i get nan values

Hi — I think you should use bigger values of the mean photon number in the thermal state. Maybe try 10^-3

I still get nan when I use the displacedsqueezed preparetion after the thermal operation

another question but please bare with me I’m sorry,

mean photon number for S gate is sinh(squezzing amp)**2 and the phase doesn’t affect right?
and for the D gate is abs(alpha)**2 and alpha is a complex number right?