getE()

function E=getE(W)

Input:

• W: Expansion coefficients for $N_s$ unconstrained wave functions, stored as an $\prod S_k\,\times\,N_s$ matrix

Global variables:

• gbl_Vdual: Dual potential coefficients stored as a $S_k\,\times\,1$ column vector.
• gbl_f: State fillings stored as an $N_s\,\times\,1$ column vector.

Output:

• E: Energies summed over $N_s$ states

Modify getE() to accept variable fillings as above. Note that, because getE() now simply passes gbl_f to getn(), there is nothing to change for the calculation of the electron density and associated terms. The calculation of the kinetic energy, however, must change. In our matrix language, the new expression is

$$T=-\frac{1}{2} \mbox{Tr\,}\left( (\mbox{\bf Diag\,} \vec f) \left(Y^\dagger L Y \right)\right),$$

where $Y$ is the matrix of expansion coefficients for the orthonormal wave functions.
Hint: The extra parentheses in the above expression limit the size of and number computations required in evaluating intermediate expressions.