The ProjMPO-ITensor product error message inside tdvp

Zhen_H · May 4, 2025, 4:04am

Hi,

I am using the tdvp function and it gives me the following error:

ProjMPO-ITensor product P*v is not equal to the order of the ITensor v.

I’ve attached the full error message in the end. Note that I have not tweeked anything inside tdvp. In other words, I am using tdvp straightly.

Here is my julia code. I have a psi loaded from a h5 file. Then I do the following:

n = 9
sites = siteinds(psi)
cdag = op("Cdag", sites, n)
psi[n] = psi[n]*cdag

psi = tdvp(H, t, psi; time_step=dt)

I am a bit stuck on debugging this, and I am wondering if there are any suggestions on what could possibly lead to this error. Maybe it is worth mentioning that if I commented out psi[n] = psi[n]*cdag, then the whole thing runs without error.
The full error message is as follows:

ERROR: LoadError: The order of the ProjMPO-ITensor product P*v is not equal to the order of the ITensor v, this is probably due to an index mismatch.
Common reasons for this error: 
(1) You are trying to multiply the ProjMPO with the 2-site wave-function at the wrong position.
(2) `orthogonalize!` was called, changing the MPS without updating the ProjMPO.

P*v inds: ((dim=2|id=845|"Fermion,Site,n=9"), (dim=2|id=845|"Fermion,Site,n=9")'', (dim=4|id=535|"sum")', (dim=2|id=978|"Fermion,Site,n=2")', (dim=2|id=288|"Fermion,Site,n=1")') 

v inds: ((dim=2|id=288|"Fermion,Site,n=1"), (dim=2|id=978|"Fermion,Site,n=2"), (dim=4|id=535|"sum"))
Stacktrace:
  [1] error(s::String)
    @ Base ./error.jl:35
  [2] product(P::ProjMPO, v::ITensor)
    @ ITensorMPS ~/.julia/packages/ITensorMPS/hghpU/src/abstractprojmpo/abstractprojmpo.jl:73
  [3] AbstractProjMPO
    @ ~/.julia/packages/ITensorMPS/hghpU/src/abstractprojmpo/abstractprojmpo.jl:87 [inlined]
  [4] apply
    @ ~/.julia/packages/KrylovKit/jC5gU/src/apply.jl:2 [inlined]
  [5] expintegrator(A::ProjMPO, t::Float64, u::Tuple{ITensor, ITensor}, alg::KrylovKit.Arnoldi{KrylovKit.ModifiedGramSchmidt2, Float64})
    @ KrylovKit ~/.julia/packages/KrylovKit/jC5gU/src/matrixfun/expintegrator.jl:109
  [6] expintegrator
    @ ~/.julia/packages/KrylovKit/jC5gU/src/matrixfun/expintegrator.jl:102 [inlined]
  [7] expintegrator(::ProjMPO, ::Float64, ::ITensor; kwargs::@Kwargs{})
    @ KrylovKit ~/.julia/packages/KrylovKit/jC5gU/src/matrixfun/expintegrator.jl:98
  [8] expintegrator
    @ ~/.julia/packages/KrylovKit/jC5gU/src/matrixfun/expintegrator.jl:94 [inlined]
  [9] exponentiate
    @ ~/.julia/packages/KrylovKit/jC5gU/src/matrixfun/exponentiate.jl:83 [inlined]
 [10] #exponentiate_updater#765
    @ ~/.julia/packages/ITensorMPS/hghpU/src/solvers/tdvp.jl:5 [inlined]
 [11] exponentiate_updater
    @ ~/.julia/packages/ITensorMPS/hghpU/src/solvers/tdvp.jl:4 [inlined]
 [12] region_update!(nsite_val::Val{2}, reverse_step_val::Val{true}, reduced_operator::ProjMPO, state::MPS, b::Int64; updater::typeof(ITensorMPS.exponentiate_updater), updater_kwargs::@NamedTuple{}, current_time::Float64, time_step::Float64, outputlevel::Int64, normalize::Bool, direction::Base.Order.ForwardOrdering, noise::Bool, which_decomp::Nothing, svd_alg::Nothing, cutoff::Float64, maxdim::Int64, mindim::Int64, maxtruncerr::Float64)
    @ ITensorMPS ~/.julia/packages/ITensorMPS/hghpU/src/solvers/sweep_update.jl:408
 [13] region_update!(reduced_operator::ProjMPO, state::MPS, b::Int64; updater::Function, updater_kwargs::@NamedTuple{}, nsite::Int64, reverse_step::Bool, current_time::Float64, outputlevel::Int64, time_step::Float64, normalize::Bool, direction::Base.Order.ForwardOrdering, noise::Bool, which_decomp::Nothing, svd_alg::Nothing, cutoff::Float64, maxdim::Int64, mindim::Int64, maxtruncerr::Float64)
    @ ITensorMPS ~/.julia/packages/ITensorMPS/hghpU/src/solvers/sweep_update.jl:187
 [14] sub_sweep_update(direction::Base.Order.ForwardOrdering, reduced_operator::ProjMPO, state::MPS; updater::Function, updater_kwargs::@NamedTuple{}, which_decomp::Nothing, svd_alg::Nothing, sweep::Int64, current_time::Float64, time_step::Float64, nsite::Int64, reverse_step::Bool, normalize::Bool, observer!::ITensorMPS.EmptyObserver, outputlevel::Int64, maxdim::Int64, mindim::Int64, cutoff::Float64, noise::Bool)
    @ ITensorMPS ~/.julia/packages/ITensorMPS/hghpU/src/solvers/sweep_update.jl:107
 [15] sweep_update(order::ITensorMPS.TDVPOrder{2, Base.Order.ForwardOrdering()}, reduced_operator::ProjMPO, state::MPS; current_time::Float64, time_step::Float64, kwargs::@Kwargs{updater::typeof(ITensorMPS.exponentiate_updater), updater_kwargs::@NamedTuple{}, nsite::Int64, reverse_step::Bool, sweep::Int64, observer!::ITensorMPS.EmptyObserver, normalize::Bool, outputlevel::Int64, maxdim::Int64, mindim::Int64, cutoff::Float64, noise::Bool})
    @ ITensorMPS ~/.julia/packages/ITensorMPS/hghpU/src/solvers/sweep_update.jl:25
 [16] macro expansion
    @ ~/.julia/packages/ITensorMPS/hghpU/src/solvers/alternating_update.jl:69 [inlined]
 [17] macro expansion
    @ ./timing.jl:421 [inlined]
 [18] alternating_update(operator::MPO, init::MPS; updater::Function, updater_kwargs::@NamedTuple{}, nsweeps::Int64, checkdone::Returns{Bool}, write_when_maxdim_exceeds::Nothing, nsite::Int64, reverse_step::Bool, time_start::Float64, time_step::Float64, order::Int64, observer!::ITensorMPS.EmptyObserver, sweep_observer!::ITensorMPS.EmptyObserver, outputlevel::Int64, normalize::Bool, maxdim::Vector{Int64}, mindim::Int64, cutoff::Vector{Float64}, noise::Bool)
    @ ITensorMPS ~/.julia/packages/ITensorMPS/hghpU/src/solvers/alternating_update.jl:68
 [19] alternating_update
    @ ~/.julia/packages/ITensorMPS/hghpU/src/solvers/alternating_update.jl:23 [inlined]
 [20] #tdvp#767
    @ ~/.julia/packages/ITensorMPS/hghpU/src/solvers/tdvp.jl:86 [inlined]
 [21] top-level scope

Thanks in advance!

mtfishman · May 4, 2025, 11:10am

You should change psi[n] = psi[n] * cdag to psi[n] = apply(cdag, psi[n]).

cdag is an operator with a pair of primed and unprimed indices, so when you contract it with psi[n]the result is a tensor where the site index is primed, and therefore isn’t in the proper form when you put it into tdvp. apply makes sure that the site index gets unprimed.

Zhen_H · May 4, 2025, 5:15pm

Thanks Matthew! That seems to be the problem.

miles · May 4, 2025, 5:38pm

Also, just since I see you are working with fermions, make sure when you apply Cdag or C to a site of an MPS that unless it’s the first site, you also apply the Jordan-Wigner string “F” operator to all previous sites. Otherwise you will get incorrect results.

Zhen_H · May 4, 2025, 6:08pm

Thanks Miles, great point!

Just to make sure I understand correctly, when constructing MPOs using H = MPO(os, sites) where sites are fermions, this Jordan-Wigner business is taken care of automatically, right?

miles · May 4, 2025, 6:34pm

Yes, that is correct. But not when using individual operators outside an MPO.

system · May 14, 2025, 6:35pm

This topic was automatically closed 10 days after the last reply. New replies are no longer allowed.

Topic		Replies	Views
Generalization of the TDVP function DMRG and Numerical Methods julia , mpo	5	610	May 24, 2023
tdvp time developing ITensor Julia Questions julia , mps	4	739	December 28, 2022
Problem using tdvp in julia ITensor Julia Questions julia	2	199	December 3, 2023
Errors in TDVP function DMRG and Numerical Methods julia , mps	10	799	June 6, 2023
Regarding errors encountered when using the expand function ITensor Julia Questions tdvp	0	53	July 6, 2024

The ProjMPO-ITensor product error message inside tdvp

Related topics