deletions | additions       diff --git a/Sternheimer2.tex b/Sternheimer2.tex  new file mode 100644  index 0000000..a8ca614  --- /dev/null  +++ b/Sternheimer2.tex 
...
We can solve for linear response to various different perturbations. The most straight-forward case is the response of a finite system to an electric field \(\mathcal{E}_{i,\omega}\) with frequency \(\omega\) in the direction \(i\) , where the perturbation operator is \(\delta \hat v = \hat r_i\) . In this case the polarizability can be calculated as  %  \begin{align}  \alpha_{ij} \left( \omega \right) = - \sum_n f_n \left[ \left< \varphi_n \left| r_i \right| \frac{\partial \varphi_n}{\partial \mathcal{E}_{j, \omega}} \right> + \left< \frac{\partial \varphi_n}{\partial \mathcal{E}_{j, -\omega}} \left| r_i \right| \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}.  It is also possible to compute vibrational frequencies~\cite{Baroni_2001,Kadantsev_2005},  in this case we calculate the response to a ionic displacement perturbation  $\partial H/\partial R_{i \alpha} = \partial V_{\alpha}/\partial x_i$, for each direction $i$ and atom $\alpha$.  % FIXME: make clear this is a non-local operator in general  %currently only implemented for finite systems.  The elements of the dynamical matrix or Hessian are  \begin{align}  D_{i \alpha, j \beta} = \frac{\partial^2 E}{\partial R_{i \alpha} \partial R_{j \beta}} = D_{i \alpha, j \beta}^{\rm ion-ion} + \sum_n \left[ f_n \left< \varphi_k \left| \frac{\partial V_{\alpha}}{\partial x_i} \right| \frac{\partial \varphi_n}{\partial R_{j \beta}} \right> + {\rm cc.}  + \delta_{\alpha \beta} \left< \varphi_n \left| \frac{\partial^2 V_{\alpha}}{\partial x_i \partial x_j} \right| \varphi_n \right> \right]  \end{align}  The ion-ion term is  \begin{align}  D_{i \alpha, j \beta}^{\rm ion-ion} =  \begin{cases}  Z_\alpha \sum_\gamma Z_\gamma \left[ \frac{\delta_{ij}}{\left|R_\alpha - R_\gamma \right|^3} - 3 \frac{\left( R_{i \alpha} - R_{i \gamma} \right)}{\left|R_\alpha - R_\gamma \right|^4} \right] & \textrm{if}\ \alpha = \beta \\  -Z_\alpha Z_\beta \left[ \frac{\delta_{ij}}{\left|R_\alpha - R_\beta \right|^3} - 3 \frac{\left( R_{i \alpha} - R_{i \beta} \right)}{\left|R_\alpha - R_\gamma \right|^4} \right] & \textrm{if}\ \alpha \ne \beta  \end{cases}  \end{align}  where $Z_\alpha$ is the ionic charge of atom $\alpha$.  The matrix elements of the perturbation are best calculated by differentiating the wavefunctions, as discussed for the forces in section \ref{sec:forces}. % FIXME  \begin{align}  \left< \varphi_n \left| \frac{\partial V_{\alpha}}{\partial x_i} \right| \frac{\partial \varphi_n}{\partial R_{j \beta}} \right> =  \left< \varphi_n \left| V_{\alpha} \right| \frac{\partial^2 \varphi_n}{\partial R_{j \beta} \partial x_i} \right> \\  \left< \varphi_n \left| \frac{\partial^2 V_{\alpha}}{\partial x_i \partial x_j} \right| \varphi_n \right> =  \left< \frac{\partial^2 \varphi_n}{\partial x_i \partial x_j} \left| V_{\alpha} \right| \varphi_n \right> + {\rm cc.} +  \left< \frac{\partial \varphi_n}{\partial x_i} \left| V_{\alpha} \right| \frac{\partial \varphi_n}{\partial x_j} \right> + {\rm cc.}  \end{align}  The dynamical matrix is symmetrized to enforce  \begin{align}  D_{i \alpha, j \beta} = D_{j \beta, i \alpha}  \end{align}  and then vibrational frequencies $\omega$ are obtained by solving the eigenvalue equation  \begin{align}  \frac{1}{\sqrt{m_\alpha m_\beta}} D_{i \alpha, j \beta} y_{j \beta} = -\omega^2 y_{j \beta}  \end{align}  For a finite system of $N$ atoms, there should be 3 zero-frequency translational modes and 3  zero-frequency rotational modes. However, they may appear at positive or imaginary frequencies,  due to the finite size of the simulation domain, the discretization of the grid, and finite precision  in solution of the ground state and Sternheimer equation. Improving convergence brings them closer to zero.  The Born effective charges can be computed from the response of the dipole moment to ionic displacement:  \begin{align}  Z^{*}_{i j \alpha} = -\frac{\partial^2 E}{\partial \mathcal{E}_i \partial R_{j \alpha}}  = \frac{\partial \mu_i}{\partial R_{j \alpha}} = Z_{\alpha} \delta_{ij} + \sum_k \left< \varphi \left| x_i \right| \frac{\partial \varphi_k}{\partial R_{j \alpha}} \right>  \end{align}  The intensity of each mode for absorption of radiation polarized in direction $i$ is calculated by multiplying by the normal mode eigenvector $y$:  \begin{align}  I_i = \sum_{j \alpha} Z^{*}_{ij \alpha} y_{j \alpha}  \end{align}  These intensities can be used to predict infrared spectra.  The Born charges must obey the acoustic sum rule, from translational invariance:  \begin{align}  \sum_{\alpha} Z^{*}_{i j \alpha} = Q_{\rm tot} \delta_{ij}  \end{align}  For each $ij$, we enforce this sum rule by distributing the discrepancy equally among the atoms, and thus obtaining corrected Born charges: % cite Chang, Gonze  \begin{align}  \tilde{Z}^{*}_{i j \alpha} = Z^{*}_{i j \alpha} + \left( Q_{\rm tot} \delta_{ij} - \sum_{\alpha} Z^{*}_{i j \alpha} \right) / N  \end{align}  The discrepancy arises from the same causes as the non-zero translational and rotational modes.  The Sternheimer equation can be used in conjunction with $\vec{k} \cdot \vec{p}$ perturbation theory~\cite{Cardona_1966} to obtain  band velocities and effective masses, as well as to apply electric fields via the quantum theory of polarization.  We consider the effective Hamiltonian for the \(\vec{k}\)-point.  % FIXME: use of k for bands interferes with k-points here!!  We compute band group velocities in a periodic system from  \begin{align}  v_{n k} = \frac{1}{\hbar} \frac{\partial \epsilon_k}{\partial k}   = \frac{1}{\hbar} \left< \varphi_k \left| \frac{\partial H_k}{\partial k} \right| \varphi_k \right> \\  \frac{\partial H_k}{\partial k} = -i \nabla + \vec{k} + \left[ V_{\alpha}, \vec{r} \right]  \end{align}  where $H_k$ is the effective Hamiltonian for the k-point, \textit{i.e.} with $\nabla \rightarrow \nabla + i \vec{k}$.  non-SCF Sternheimer, one-shot  we use dH/d(ik) so real-> real.  adding in occupied contributions  Inverse effective mass tensors can be calculated by solving the Sternheimer equation with the perturbation  \begin{align}  m^{-1}_{ijnk} = \frac{1}{\hbar^2} \frac{\partial^2 \epsilon_{nk}}{\partial k_i \partial k_j}  = \delta_{ij} + \left< \varphi_{nk} \left| dH/dk \right| d \varphi_{nk}/dk \right> + {\rm cc.}  + \left< \varphi_{nk} \left| \left[ \vec{r}, \left[ \vec{r}, V_{\alpha} \right] \right] \right| \varphi_{nk} \right>\ .  \end{align}  The $\vec{k} \cdot \vec{p}$ wavefunctions can be used to compute response to electric fields in periodic systems.  In finite systems, a homogeneous electric field can be represented simply via the position operator $\vec{r}$.  However, this operator is not well defined in a periodic system and cannot be used. According to the quantum theory of polarization,  the solution is to replace $\vec{r} \varphi$ with $-i d/dk \varphi$. % cite Baroni, also Gonze papers  The response to this perturbation can be used to compute the polarizability using from eq.~\ref{eq:sternheimer_polarizability}  \begin{align}  \alpha_{ij} \left( \omega \right) = i \sum_k f_k \left[ \left. \left< \frac{\partial \varphi_k}{\partial k_i} \right| \frac{\partial \varphi_k}{\partial \mathcal{E}_{j, \omega}} \right> + \left. \left< \frac{\partial \varphi_k}{\partial \mathcal{E}_{j, -\omega}} \right| \frac{\partial \varphi_k}{\partial k_j} \right> \right]  \end{align}  The polarizability is most usefully represented in a periodic system via the dielectric constant  \begin{align}  \epsilon_{ij} = \delta_{ij} + \frac{4 \pi}{V} \alpha_{ij}  \end{align}  where $V$ is the volume of the unit cell.  We can compute the Born charges from the electric-field response in either finite or periodic systems (as a complementary approach to using the vibrational response):  \begin{align}  Z^{*}_{i j \alpha} = -\frac{\partial^2 E}{\partial \mathcal{E}_i \partial R_{j \alpha}} = \frac{\partial F_{j \alpha}}{\partial \mathcal{E}_i}  =   \end{align}  single-point Berry phase for polarization.  'symmetrization' of $\epsilon$ and Born charge tensors  solvers: cg, bicgstab, qmr. when applicable.  % ought to discuss dichroism?  % Arto's stuff  % magnetic is non-self-consistent if wfns are real.  cite my chapter on linear response