deletions | additions       diff --git a/Sternheimer.tex b/Sternheimer.tex  index d8b15b0..f9777e0 100644  --- a/Sternheimer.tex  +++ b/Sternheimer.tex 
...
resonance~\cite{Andrade_2007}. The method is suited for linear and  non-linear response; higher-order Sternheimer equations can be  obtained for higher-order variations. For second-order response,  however, we can apply the \(2n\,+\,1\) \(2N\,+\,1\)  theorem (also known as Wigner's \(2n\,+\,1\) \(2N\,+\,1\)  rule)~\cite{Gonze_1989,Corso_1996} to get the coefficients directly from first-order response variations.  In the Sternheimer formalism we consider the response to   monochromatic perturbative field \(\lambda  \delta{V}(\vec{r})\cos\left(\omega{t}\right)\). This perturbation induces a variation in the time-dependent Kohn-Sham orbitals, that we denote $\delta\varphi_{k}(\vec $\delta\varphi_{n}(\vec  r, \omega)$. These variations allows to calculate response quantities like the frequency-dependent polarization. In order to calculate the variations of the orbitals we need to solve a linear equation that only depends on the occupied orbitals (atomic units are used throughout)  %\begin{multline}  \begin{equation}  \label{eq:sternheimer}  \left\{H - \epsilon_k\pm\omega \epsilon_n\pm\omega  + \mathrm{i}\eta\right\}\delta\varphi_{k}(\vec \mathrm{i}\eta\right\}\delta\varphi_{n}(\vec  r, \pm\omega) =\\  -\mathrm{P}_c\,\delta{H}(\pm\omega) \varphi_k(\vec \varphi_n(\vec  r) \,,  \end{equation}  % 
...
\begin{equation}  \label{eq:varrho}  \delta{n}(\vec r, \omega) = \sum_k^{\rm occ.} \Big\{  \left[\varphi_k(\vec{r})\right]^*\delta\varphi_{k}(\vec \left[\varphi_n(\vec{r})\right]^*\delta\varphi_{n}(\vec  r, \omega)\\ + \left[\delta\varphi_{k}(\vec \left[\delta\varphi_{n}(\vec  r, -\omega)\right]^*\varphi_k(\vec{r}) -\omega)\right]^*\varphi_n(\vec{r})  \Big\}\ ,  \end{equation}  %\end{multline} 
...
correct position of the poles of the causal response function, and,  consequently obtain the imaginary part of the  po\-la\-ri\-za\-bi\-li\-ty and remove the divergences of  the equation for resonant frequencies. In the usual implementation of DFTP, $P_c = 1 - \sum^{\rm occ}_k \left| \varphi_k \varphi_n  \right> \left< \varphi_k \varphi_n  \right|$ which effectively removes the components of \(\delta\varphi_{k}(\vec \(\delta\varphi_{n}(\vec  r, \pm\omega)\) in the subspace of the occupied ground-state wave-functions. In linear  response, these components do not contribute to the variation of the  density.  We have found that it is not strictly necessary to project out the occupied subspace,   the crucial part is simply to remove the  projection of $\delta \varphi_k$ \varphi_n$  on $\varphi_k$ $\varphi_n$  (and any other states degenerate with it), which is not physically meaningful and arises from a phase convention. To fix the phase, it is  sufficient to apply a minimal projector  $P_k = 1 - \sum_l^{\epsilon_l = \epsilon_k} \epsilon_n}  \left| \varphi_l \right> \left< \varphi_l \right|$. We optionally use this approach to obtain the entire response wavefunction, not just the projection in the  unoccupied subspace, which is needed for obtaining effective masses. While the full projection can become time-consuming  for large systems, it saves time overall since it increases the condition number of the matrix for the linear solver, 
...
We can solve for linear response to various different perturbations. The most straight-forward case is the response to an electric field \(\mathcal{E}_i\) in finite systems, where the perturbation is \(\delta V = r_i\) . In this case the polarizability as  %  \begin{align}  \alpha_{ij} \left( \omega \right) = - \sum_k f_k \sum_n f_n  \left[ \left< \varphi_k \varphi_n  \left| r_i \right| \frac{\partial \varphi_k}{\partial \varphi_n}{\partial  \mathcal{E}_{j, \omega}} \right> + \left< \frac{\partial \varphi_k}{\partial \varphi_n}{\partial  \mathcal{E}_{j, -\omega}} \left| r_i \right| \varphi_k \varphi_n  \right> \right] \label{eq:sternheimer_polarizability}\ . \end{align}  %  The calculations of the polarizability yield optical response properties~\cite{Andrade_2007,Vila_2010} and, for imaginary frequencies, van der Waals coefficients~\cite{Botti_2008}.