Regarding errors encountered when using the expand function

Hello everyone. Recently, I encountered some unresolved errors while using the expand function for tdvp evolution.
As shown in the following example, I am using a two-component SSH model Hamiltonian. When I use 2-TDVP, it can get the correct results.

function hamiltonian(L,t1,t2,U,sites)
    os = OpSum()
    for i in 1:2:2L-2
        if isodd((i+1)/2)
            os += -t1, "Adag", i, "A", i+2
            os += -t1, "Adag", i+2, "A", i
        elseif iseven((i+1)/2)
            os += -t2, "Adag", i, "A", i+2
            os += -t2, "Adag", i+2, "A", i
        end
    end
    for i in 1:2:2L
        os += U/2, "N",i,"N",i
        os += -U/2, "N",i
    end
    Heff = MPO(os,sites)
    return Heff
end;
function increase_bond_dimension!(psi,dimension)
    for bond in 1:(length(psi)-1)
        bond_index = commonind(psi[bond], psi[bond + 1])
        aux = Index(dimension; tags=tags(bond_index))
        psi[bond] = psi[bond] * delta(bond_index, aux)
        psi[bond + 1] = psi[bond + 1] * delta(bond_index, aux)
    end
    return psi
end;

L = 4
t1 = 0.2
t2 = 1
U = 8
ttotal = 10
tau = 0.1
dimension = 50

sites = siteinds("Boson",2*L)
states = ["1","0","1","0"]
psi0 = MPS(sites,states)
psi0 = increase_bond_dimension!(psi0,dimension)
Ham = hamiltonian(L,t1,t2,U,sites)

for j in 1:tsteps
    psi = tdvp(Ham,-im*tau,psi;time_step=-im*tau,normalize=true,maxdim=500,cutoff=0,nsite=2,outputlevel=1)
end


But when I added the expand function and changed it to 1-TDVP, it was obvious that there was an error in his result.

for j in 1:tsteps
        psi=expand(psi, kry_mpo;alg="global_krylov", krylovdim=4, cutoff=0)
        psi = tdvp(Ham,-im*tau,psi;time_step=-im*tau,normalize=true,maxdim=500,cutoff=0,nsite=1,outputlevel=1)
end

I am curious if this is because I did not use the expand function correctly, or if there were other errors. I’m glad everyone is willing to teach.