Refactor ode

[halvarsu]

Implement same changes as configuration-interaction/pull/#4 for the time-dependent solvers.

[halvarsu]

removed the template-function again

[halvarsu]

Any comments before I apply the same treatment to the other truncation levels and OA?:)

[Schoyen]

The only immediate comments I can think of is to get rid of self._amplitudes (the time-dependent solver should not store the amplitudes other than the template), this means that there is no more need for a set_initial_conditions-method. Instead this should be done in the integrator class. An example use-case (with potential bugs) is:

import numpy as np
import matplotlib.pyplot as plt

from scipy.integrate import complex_ode

from quantum_systems import GeneralOrbitalSystem, ODQD
from quantum_systems.time_evolution_operators import LaserField
from coupled_cluster import CCSD, TDCCSD
from gauss_integrator import GaussIntegrator

n = 2
l = 6
omega = 1

gos = GeneralOrbitalSystem(n, ODQD(l, 11, 201))
gos.set_time_evolution_operator(LaserField(lambda t: np.sin(omega * t)))

cc = CCSD(gos, verbose=True).compute_ground_state()
tdcc = TDCCSD(gos, verbose=True)

r = complex_ode(tdcc).set_integrator("GaussIntegrator")

# Ground state initial conditions
r.set_initial_value(cc.get_amplitudes(get_t_0=True))

t_final = 5
dt = 1e-2

num_steps = int(t_final / dt) + 1
energy = np.zeros(num_steps, dtype=np.complex128)

i = 0

while r.successful() and r.t <= t_final:
    if i % 100 == 0:
        print(f"{i} / {num_steps}")

    # All observables should be computed by passing in the time and the
    # amplitude vector from the integrator
    energy[i] = tdcc.compute_energy(r.t, r.y)
    r.integrate(r.t + dt)

    i += 1

plt.plot(np.linspace(0, t_final, num_steps), energy.real)
plt.show()

[halvarsu]

Pretty sure I did that for tdcc and tdccd in 3d58a43 :) I used something very similar to your script for TDCCD and it works. Can add it.

[halvarsu]

Ok, I think that should be it. Ready for review :)

[Schoyen]

I'll have a look at it in the morning!

[halvarsu]

This is beginning to look very good, but there are a few changes that needs to be addressed before merging. These are:

1. All `compute_energy` functions need to take in the current time to make sure that the Hamiltonian is updated to the correct time step. This might solve some of the failing tests.

2. With this new layout and a renaming for the matrix elements from `self.*_orig` to `self.*` there is the potential for more reuse of the parent functionality. Especially in the `update_hamiltonian` function and the constructor.

3. Good that you've added tests for `TDCCS`, however these tests needs to be fixed. For now I am content if you mark them with `@pytest.mark.skip` to skip these tests for now.

1. Done, it did not fix the failing tests (i.e. the numerical ones, right?). I am pretty sure they are failing due to difference in integrator, as the results are within the change I get from changing to yet other integrators.

2. update_hamiltonian is only defined twice, right? It should be fine to keep two versions of it for now, the large hurdle is all the compute_*-methods which all basically pass the arguments along to some amplitude-evaluated expression.

3. 👍