def energy_tot(e, N, P, H, Vnn):
    """
    Compute the total energy.

    INPUT:
    e: MO energies
    N: num of electrons
    P: density matrix
    H: h_core matrix
    Vnn: nuclear nuclear repulsion energy
    """
    e_tot = 0
    for i in range(int(N/2)):
        e_tot += e[i].real
    e_tot = e_tot + 0.5 * (P * H).sum() + Vnn
    return e_tot