Atomic orbital normalization done within unit cell instead of full real-space domain

[nichousha233]

here
In my opinion, the normalization of the wavefunction should be performed using the global real-space wavefunction before folding, rather than the one confined to a single unit cell after folding.

[1040403998]

Hi @nichousha233, thanks a lot for taking the time to carefully look into the PDOS code — really appreciate your attention to detail!

Regarding the normalization, the current implementation performs normalization after folding the atomic orbitals into the unit cell with the correct Bloch phase. This ensures that the projected states are consistent with the periodic boundary conditions and that the projection coefficients remain physically meaningful (i.e., ≤ 1).

If normalization were instead done before folding, the normalization factors would be independent of the k-point, but the resulting projections could exceed unity when the cutoff radius is large. To be more precise, in our current workflow the Kohn-Sham (KS) orbitals are normalized within the unit cell, i.e.

$$ \int_{\text {cell }}\left|\psi_{n \mathbf{k}}(\mathbf{r})\right|^2 d \mathbf{r}=1 $$

To make the projector consistent with that inner product, we first fold the localized atomic orbital with the Bloch phase

$$ \tilde{\phi}_\mu^{(\mathbf{k})}(\mathbf{r})=\sum_{\mathbf{R}} \phi_\mu(\mathbf{r}-\mathbf{R}) e^{\mathbf{i} \mathbf{k} \cdot \mathbf{R}} $$

and then normalize $\tilde{\phi}_\mu^{(\mathbf{k})}$ over the unit cell:

$$ \int_{\text {cell }}\left|\tilde{\phi}_\mu^{(\mathbf{k})}(\mathbf{r})\right|^2 d \mathbf{r}=1 $$

If instead one normalizes the unfolded $\phi_\mu$ in real space (over all space) before folding, the two states being overlapped are normalized over different domains/measures (global $\mathbb{R}^3$ vs. unit cell). In that case, the standard Cauchy–Schwarz bound

$$ \lvert \langle \psi_{n,\mathbf{k}} \mid \tilde{\phi}_\mu^{(\mathbf{k})} \rangle \rvert \le 1 $$

no longer applies, and for large cutoff radii the overlap $\lvert \langle \psi_{n,\mathbf{k}} \mid \tilde{\phi}_\mu^{(\mathbf{k})} \rangle \rvert^2$ can exceed 1.

When the cutoff is small, both procedures are numerically close, which is why this can be easy to miss. But for robustness across systems, we keep the post-fold, unit-cell normalization as the default. In our previous comparisons with Quantum ESPRESSO (QE), the PDOS results agree very well, and I believe QE adopts the same normalization convention.

[nichousha233]

Hi @1040403998 , thanks for your kind words! I'm glad if I could help, and I've been following the SPARC project with great interest — really appreciate all the great work you’re doing on it.

Regarding the PDOS calculation, there is an equivalent approach where the SPARC wavefunctions are extended to the full real space and then integrated with the atomic orbitals in real space.
In this case, the effect of the k-points is entirely contained in the SPARC wavefunctions, while the atomic orbitals themselves are not affected by k.
Therefore, the normalization factor of the atomic orbital basis should, in my opinion, be independent of k.

In the current implementation, the normalization factor depends on k because it is computed from the sum of multiple atomic orbitals at each grid point.
I think this treatment is not quite appropriate, and this is one of the main reasons why the normalization factor becomes k-dependent in the current script.

[phanish-suryanarayana]

@nichousha233 Thanks for bringing this subtle point. What you are saying is definitely correct. However, in practice, to keep all computations to the simulation cell, the interaction of the wrap around of the atomic orbitals with itself is not neglected (though it should technically be neglected, as you mention). This approximation is commonly adopted, so we have done the same. In the future, maybe we can generalize the code but since our focus is systems on the larger side, this issue should not occur there. Indeed, this issue occurs when the compact support of the atomic orbitals is greater or equal to the cell dimensions.

[1040403998]

Hi @nichousha233, thanks a lot for your interest in the project and for taking the time to share your thoughts!

I think there might be a small misunderstanding regarding what the normalization means in this context. Conceptually, folding the atomic orbitals into the unit cell and extending the KS wavefunctions over the full real-space lattice are simply two mathematically equivalent ways of formulating the same projection. In principle, if both states are normalized consistently over the same domain, they should produce identical overlap values.

However, in practice, when we do not include an additional normalization factor, the computed overlaps can exceed 1 — which is clearly unphysical. This happens because, for small unit cells, the integral domains of the two functions (the extended KS wavefunctions and the localized atomic orbitals) are not normalized on the same measure. In other words, there is a missing normalization factor.

Once we think carefully about how this factor should be defined, performing the normalization after folding the atomic orbitals into the unit cell becomes the more natural and well-defined approach. If we were to extend the wavefunctions instead, the normalization domain — effectively the size of the supercell one chooses — becomes arbitrary, which makes the normalization ambiguous. When we enforce normalization inside the unit cell, the normalization factors naturally become k-dependent, because the Bloch phase changes the distribution of the orbital weight inside the cell.

In principle, one could use the Γ-point normalization factor for all k-points to remove this dependence. But then, at non-Γ k-points, the folded orbitals would be scaled by a factor that no longer matches their own Bloch periodicity, which is inconsistent. Therefore, we keep the k-dependent normalization as the physically meaningful choice.

This is also what QE effectively does — its implementation implicitly includes this unit-cell normalization, even though it is not expressed explicitly. And indeed, this is the only way to ensure that the PDOS projections agree quantitatively with QE’s results.

As Phanish mentioned, all these differences vanish when the unit cell becomes large enough, since the overlap between periodic images becomes negligible. But for small systems, this distinction is necessary to maintain physically valid projection coefficients.

I hope this clarifies what I meant earlier, and I really appreciate the continued discussion — it’s always helpful to have more eyes and perspectives on these details.

[nichousha233]

Hi, I’ve been thinking about this issue over the past few days, and I now agree with your view and the way you’ve handled it. I really appreciate your detailed explanation — I’ve learned a lot from it.

By the way, I think there might be a small issue in this part of the code: the newly generated cell seems to be farther away from the atom than expected.