Refactor generalized eigensolvers to avoid explicit Cholesky inverse

[QG-phy]

Background

Several code paths solve generalized Hermitian eigenproblems with an overlap matrix using Cholesky orthogonalization. For example, dptb/nn/dftb/dftb_scc.py::SKSCC.eigh_solver currently computes:

chklowt = torch.linalg.cholesky(overlap_mat)
chklowtinv = torch.linalg.inv(chklowt)
H_mat = chklowtinv @ H_mat @ chklowtinv.transpose(1, 2).conj()
eigval, eigvec = torch.linalg.eigh(H_mat)
eigvec = chklowtinv.transpose(1, 2).conj() @ eigvec

This is mathematically valid, but explicitly inverting the triangular Cholesky factor is less efficient and less numerically stable than using torch.linalg.solve_triangular.

Proposed direction

Introduce a shared helper for generalized Hermitian eigenproblems, e.g.

solve_generalized_eigh(H, S=None) -> eigval, eigvec

with a single documented convention for eigenvector shape, likely (nk, norb, nstate).

The helper should replace duplicated logic in:

• dptb/nn/dftb/dftb_scc.py::SKSCC.eigh_solver
• dptb/nn/energy.py::Eigh
• potentially dptb/nn/energy.py::Eigenvalues

Important details

For S = L L†, the standard transformation is:

H' = L⁻¹ H L⁻†
H' c = E c
psi = L⁻† c

When replacing explicit inverse with triangular solves, the right-side transform and eigenvector back-transform need careful treatment. A naïve patch may preserve eigenvalues while returning incorrectly transformed eigenvectors, which can silently affect Mulliken charges and SCC workflows.

Suggested tests

Add tests that compare the new helper against the current inverse-based implementation for:

• batched real Hermitian matrices
• batched complex Hermitian matrices
• eigenvalue agreement
• generalized eigenvector orthonormality: eigvec.conj().T @ S @ eigvec == I
• consistency with existing SCC / Eigh eigenvector shape conventions

Notes

This is not currently blocking PR #330. It is a valid numerical/performance improvement, but should be handled as a focused follow-up with tests rather than as a quick inline replacement.