Thank you Miles, that clarifies a lot.

But somehow the productMPS is still not working and gives me the error "LoadError: ArgumentError: Overload of "state" function not found for Index tags "CFermion,Site,n=1", even though I can print out the states and sites in productMPS(sites, states).

I copy and paste my code below for your reference: (I followed the "Electron" site type closely)

using ITensors

ITensors.space(::SiteType"CFermion") = 4

function ITensors.space(

::SiteType"CFermion";

conserve_qns=false,

conserve_Q1=conserve_qns,

conserve_nfparity=conserve_qns,

qnname_Q1="Q1",

qnname_nfparity="NfParity",

)

if conserve_nfparity && conserve_Q1

return [

QN((qnname_nfparity, 0, -2), (qnname_Q1, 0)) => 1,

QN((qnname_nfparity, 1, -2), (qnname_Q1, +1)) => 1,

QN((qnname_nfparity, 1, -2), (qnname_Q1, +2)) => 1,

QN((qnname_nfparity, 0, -2), (qnname_Q1, +3)) => 1

]

end

return 4

end

#val(::ValName"Emp", ::SiteType"CFermion") = 1

#val(::ValName"f1", ::SiteType"CFermion") = 2

#val(::ValName"f2", ::SiteType"CFermion") = 3

#val(::ValName"f1f2", ::SiteType"CFermion") = 4

ITensors.state(::StateName"Emp", ::SiteType"CFermion") = [1.0,0,0,0]

ITensors.state(::StateName"f1", ::SiteType"CFermion") = [0.0,1,0,0]

ITensors.state(::StateName"f2", ::SiteType"CFermion") = [0.0,0,1,0]

ITensors.state(::StateName"f1f2", ::SiteType"CFermion") = [0.0,0,0,1]

function ITensors.op!(Op::ITensor, ::OpName"C1", ::SiteType"CFermion", s::Index)

Op[s' => 1, s => 2] = 1.0

return Op[s' => 3, s => 4] = 1.0

end

function ITensors.op!(Op::ITensor, ::OpName"C1d", ::SiteType"CFermion", s::Index)

Op[s' => 2, s => 1] = 1.0

return Op[s' => 4, s => 3] = 1.0

end

function ITensors.op!(Op::ITensor, ::OpName"C2", ::SiteType"CFermion", s::Index)

Op[s' => 1, s => 3] = 1.0

return Op[s' => 2, s => 4] = 1.0

end

function ITensors.op!(Op::ITensor, ::OpName"C2d", ::SiteType"CFermion", s::Index)

Op[s' => 3, s => 1] = 1.0

return Op[s' => 4, s => 2] = 1.0

end

function ITensors.op!(Op::ITensor, ::OpName"N1", ::SiteType"CFermion", s::Index)

Op[s' => 2, s => 2] = 1.0

return Op[s' => 4, s => 4] = 1.0

end

function ITensors.op!(Op::ITensor, ::OpName"N2", ::SiteType"CFermion", s::Index)

Op[s' => 3, s => 3] = 1.0

return Op[s' => 4, s => 4] = 1.0

end

has_fermion_string(::OpName"C1", ::SiteType"CFermion") = true

has_fermion_string(::OpName"C1d", ::SiteType"CFermion") = true

has_fermion_string(::OpName"C2", ::SiteType"CFermion") = true

has_fermion_string(::OpName"C2d", ::SiteType"CFermion") = true

let

N = 6

sites = siteinds("CFermion",N; conserve_qns = true)

println(sites)

states = [mod(i,2)==0 ? "f1" : "f2" for i=1:N]

println(states)

psi0 = productMPS(sites,states)

println(psi0)

end