Thanks for the question, apateonas. You’re on the right track, as be able to use your operator definition in functions like expect
or correlation_matrix
you need to define it as an op
or op!
function overload as you are trying to do.
Fortunately it’s pretty easy to do! Here is an example from our “Qubit” site type of a two-site operator (the CNOT gate in this case):
ITensors.op(::OpName"CNOT", ::SiteType"Qubit") = [
1 0 0 0
0 1 0 0
0 0 0 1
0 0 1 0
]
So you just have to have the op
method overload take the name of your operator and SiteType that you want it to work for, then have the function return the (appropriately sized) matrix defining your operator.
Please let me know if you need more information than that for your case. Even though your case involves fermions, you shouldn’t have to do anything special since your operator is “bosonic” overall, meaning it involves an even number of c operators.