Sweeps and bond dimension in TDVP

Hi,
It is generally recommended that for DMRG calculation to use larger sweeps start with smaller bond dimensions and then increase slowly for better convergence. But I am confused about TDVP. Because of the GitHub example, only one sweep has been used with 2000 bond dimensions. Is it just an example to show, or does it have some specific reason to choose only one seep with a large bond dimension? What is the optimum choice to start with in TDVP?

In the following code snippit

void mean_current_evolve(int L, int N, int max_occ, double V, double U_ONSITE, double U_NN, double pump_delay, double probe_delay, double A0_pump, double A0_probe)
{
    Real dt = Time[1] - Time[0];
    auto sweeps = Sweeps(50);
    sweeps.maxdim() = 10, 10, 10, 10, 10, 20, 20, 20, 20, 20, 20, 20, 20, 20, 20, 40, 40, 40, 40, 40, 40, 40, 40, 40, 40, 80, 80, 80, 80, 80, 80, 80, 80, 80, 80, 160, 160, 160, 160, 160, 300, 300, 600, 600, 1200, 1200, 2400, 2400, 4800, 4800;
    sweeps.niter() = 6;
    sweeps.noise() = 0.0;
    sweeps.cutoff() = 1E-8;
    for (int time_point = 0; time_point < Time.size(); ++time_point)
    {
        System Sp(L, N, max_occ, V, U_ONSITE, U_NN, pump_delay, probe_delay, A0_pump, A0_probe, time_point, "p");
        System Sp_pr(L, N, max_occ, V, U_ONSITE, U_NN, pump_delay, probe_delay, A0_pump, A0_probe, Time[time_point], "p_pr");
        Sp.current();
        Sp_pr.current();
        auto Hp = Sp.H;
        auto Hp_pr = Sp_pr.H;
        auto Jp = Sp.J;
        auto Jp_pr = Sp_pr.J;
        MPS psip;
        MPS psip_pr;
        if (time_point == 0)
        {
            Sp.ground_state();
            Sp_pr.ground_state();
            psip = Sp.psi0;
            psip_pr = Sp_pr.psi0;
        }

        if (time_point < 3)
        {
            // Global subspace expansion
            vector<Real> epsilonK = {1E-12, 1E-12};
            addBasis(psip, Hp, epsilonK, {"Cutoff", 1E-8, "Method", "DensityMatrix", "KrylovOrd", 3, "DoNormalize", true, "Quiet", true});
            addBasis(psip_pr, Hp_pr, epsilonK, {"Cutoff", 1E-8, "Method", "DensityMatrix", "KrylovOrd", 3, "DoNormalize", true, "Quiet", true});
        }

        tdvp(psip, Hp, -dt, sweeps, {"Truncate", true, "DoNormalize", true, "Quiet", true, "NumCenter", 1});
        tdvp(psip_pr, Hp_pr, -dt, sweeps, {"Truncate", true, "DoNormalize", true, "Quiet", true, "NumCenter", 1});

        auto M_Jpsip = real(innerC(psip, Hp, psip));
        auto M_Jpsip_pr = real(innerC(psip_pr, Hp_pr, psip_pr));
        // auto M_Jpsip = real(innerC(psip, Jp, psip));
        // auto M_Jpsip_pr = real(innerC(psip_pr, Jp_pr, psip_pr));

        // auto Jpsip = applyMPO(Jp, psip);
        // auto M_Jpsip = inner(psip, Jpsip);
        // auto Jpsip_pr = applyMPO(Jp_pr, psip_pr);
        // auto M_Jpsip_pr = inner(psip_pr, Jpsip_pr);
        // Current.push_back(M_Jpsip_pr - M_Jpsip);
    }
}

The error comes as

From line 132, file mps/mpoalgs.cc

Error in applyMPO, MPS is uninitialized.

Error in applyMPO, MPS is uninitialized.
Aborted (core dumped)

If I use

auto M_Jpsip = real(innerC(psip, Jp, psip));
auto M_Jpsip_pr = real(innerC(psip_pr, Jp_pr, psip_pr));

error is

From line 568, file indexset.cc

In findIndex: more than one Index found, consider using findInds instead

In findIndex: more than one Index found, consider using findInds instead
Aborted (core dumped)

and If I use

        auto Jpsip = applyMPO(Jp, psip);
        auto M_Jpsip = inner(psip, Jpsip);
        auto Jpsip_pr = applyMPO(Jp_pr, psip_pr);
        auto M_Jpsip_pr = inner(psip_pr, Jpsip_pr);

The error is

From line 568, file indexset.cc

In findIndex: more than one Index found, consider using findInds instead

In findIndex: more than one Index found, consider using findInds instead
Aborted (core dumped)

I have no idea how to resolve this issue. I am Using TDVP to evolve \psi of time-dependent Hamiltonian H(t). In simple words, I am trying to solve the Schrodinger equation.

I found one mistake as

MPS psip;
MPS psip_pr;

Inside the for loop, which causes problem, but If I modify it like

void mean_current_evolve(int L, int N, int max_occ, double V, double U_ONSITE, double U_NN, double pump_delay, double probe_delay, double A0_pump, double A0_probe)
{
    Real dt = Time[1] - Time[0];
    auto sweeps = Sweeps(1);
    sweeps.maxdim() = 1000;
    sweeps.niter() = 6;
    sweeps.noise() = 0.0;
    sweeps.cutoff() = 1E-8;
    MPS psip;
    MPS psip_pr;
    for (int time_point = 0; time_point < Time.size(); ++time_point)
    {
        System Sp(L, N, max_occ, V, U_ONSITE, U_NN, pump_delay, probe_delay, A0_pump, A0_probe, time_point, "p");
        System Sp_pr(L, N, max_occ, V, U_ONSITE, U_NN, pump_delay, probe_delay, A0_pump, A0_probe, Time[time_point], "p_pr");
        Sp.current();
        Sp_pr.current();
        auto Hp = Sp.H;
        auto Hp_pr = Sp_pr.H;
        auto Jp = Sp.J;
        auto Jp_pr = Sp_pr.J;
        if (time_point == 0)
        {
            Sp.ground_state();
            Sp_pr.ground_state();
            psip = Sp.psi0;
            psip_pr = Sp_pr.psi0;
        }

        if (time_point < 3)
        {
            // Global subspace expansion
            vector<Real> epsilonK = {1E-12, 1E-12};
            addBasis(psip, Hp, epsilonK, {"Cutoff", 1E-8, "Method", "DensityMatrix", "KrylovOrd", 3, "DoNormalize", true, "Quiet", true});
            addBasis(psip_pr, Hp_pr, epsilonK, {"Cutoff", 1E-8, "Method", "DensityMatrix", "KrylovOrd", 3, "DoNormalize", true, "Quiet", true});
        }

        tdvp(psip, Hp, -dt, sweeps, {"Truncate", true, "DoNormalize", true, "Quiet", true, "NumCenter", 1});
        tdvp(psip_pr, Hp_pr, -dt, sweeps, {"Truncate", true, "DoNormalize", true, "Quiet", true, "NumCenter", 1});

        // auto M_Jpsip = real(innerC(psip, Hp, psip));
        // auto M_Jpsip_pr = real(innerC(psip_pr, Hp_pr, psip_pr));
        // auto M_Jpsip = real(innerC(psip, Jp, psip));
        // auto M_Jpsip_pr = real(innerC(psip_pr, Jp_pr, psip_pr));

        // auto Jpsip = applyMPO(Jp, psip);
        // auto M_Jpsip = inner(psip, Jpsip);
        // auto Jpsip_pr = applyMPO(Jp_pr, psip_pr);
        // auto M_Jpsip_pr = inner(psip_pr, Jpsip_pr);
        // Current.push_back(M_Jpsip_pr - M_Jpsip);
    }
}

Now the error appearing as

From line 226, file mps/mpoalgs.cc

MPS and MPO have different site indices in applyMPO method ‘DensityMatrix’

MPS and MPO have different site indices in applyMPO method ‘DensityMatrix’
Aborted (core dumped)

Hi, sorry for the slow reply on this. For an error of this type (saying the MPS and MPO have different indices), a good thing to do is to print out the MPS and MPO you are using or print out a few of their tensors. Are the site indices different? I.e. do they have different id numbers? Or are their prime levels for some reason not the conventional ones (usually two copies of each site on the MPO, one with prime level 0 and the other with prime level 1, and MPS site indices unprimed).