Riemannian gradient descent for Hartree-Fock theory

1. Problema y Contexto

El trabajo aborda el problema de optimización variacional en la teoría de estructura electrónica, específicamente en el marco de la teoría de Hartree-Fock (HF) y la teoría del funcional de la densidad (DFT) de Kohn-Sham.

• El Desafío: El objetivo es minimizar el funcional de energía $\mathcal{E}(\varphi)$ sujeto a restricciones de ortonormalidad de los orbitales moleculares $\varphi = (\varphi_1, \dots, \varphi_N)^T$. Tradicionalmente, los métodos de química computacional (como SCF-DIIS) operan en espacios de dimensión finita o utilizan iteraciones de punto fijo.
• Limitaciones de los Métodos Actuales:
  • Los métodos basados en el producto interno $L^2$ (métricas estándar) presentan dificultades porque el funcional de energía no es diferenciable con respecto a la norma $L^2$ (debido al término de energía cinética).
  • Los algoritmos de punto fijo (SCF) a menudo fallan en converger desde adivinanzas iniciales aleatorias y requieren técnicas de aceleración como DIIS (Direct Inversion in the Iterative Subspace), que pueden ser inestables.
  • La formulación continua en espacios de funciones (en lugar de matrices discretizadas) ha sido menos explorada desde una perspectiva geométrica rigurosa.

2. Metodología Propuesta

El autor propone un marco de optimización Riemanniana formulado directamente en el espacio de Sobolev infinito-dimensional $H^1(\mathbb{R}^3)$, evitando formulaciones distribucionales.

A. Marco Geométrico

• Espacio de Funciones: Se trabaja en el espacio de Sobolev $H^1$, que garantiza que los orbitales tengan energía cinética finita (gradientes cuadráticamente integrables).
• Variedades:
  • Variedad de Stiefel ($St(N)$): Se interpreta el conjunto de orbitales ortonormales como una variedad de Stiefel infinita-dimensional, incrustada en $H^1$ con la métrica $H^1$.
  • Variedad de Grassmann ($Gr(N)$): Para eliminar la invariancia de gauge rotacional (redundancia de grados de libertad), se considera el cociente $Gr(N) = St(N)/O(N)$.
• Métrica: Se utiliza la métrica $H^1$ (definida por $\langle \phi, \psi \rangle_{H^1} = \langle \nabla \phi, \nabla \psi \rangle_{L^2} + \langle \phi, \psi \rangle_{L^2}$) en lugar de la métrica $L^2$. Esto asegura que tanto el funcional de energía como la superficie de restricción sean suaves.

B. Componentes del Algoritmo

El autor deriva expresiones explícitas para los componentes necesarios en la optimización Riemanniana:

1. Gradiente Euclidiano ($\nabla \mathcal{E}$): Se calcula utilizando operadores resolventes $R = (-\Delta + 1)^{-1}$ para manejar la derivada del funcional en el espacio $H^1$.
2. Gradiente Riemanniano ($\text{grad } \mathcal{E}$): Se obtiene proyectando el gradiente euclidiano sobre el espacio tangente de la variedad.
3. Proyección Tangente: Se formula mediante la resolución de ecuaciones de Sylvester para proyectar vectores fuera de la variedad hacia el espacio tangente, manteniendo la ortogonalidad $L^2$.
4. Retracción: Se utiliza la ortogonalización de Löwdin ($\mathcal{R}_\varphi(v) = S^{-1/2}(\varphi + v)$) para mapear pasos de descenso de vuelta a la variedad de Stiefel, garantizando que los orbitales permanezcan ortonormales.
5. Precondicionamiento: Se introduce un precondicionador crítico basado en la inversión del operador de energía cinética (o una aproximación de él). Esto es análogo a la precondición natural en los métodos SCF, pero aplicado al gradiente.
   • El precondicionador es $Prec_\varphi = \text{Proj}_\varphi T^{-1}$, donde $T$ es el operador de energía cinética.
6. Algoritmos de Optimización:
   • Descenso de Gradiente Riemanniano: Con búsqueda de línea tipo backtracking (Armijo).
   • Gradiente Conjugado No Lineal Precondicionado: Incluye transporte de vectores entre espacios tangentes y reinicios tipo Powell para evitar direcciones no deseadas.

C. Discretización

Aunque el marco es continuo, los experimentos numéricos se realizan utilizando multiondas adaptativas (multiwavelets). Esta elección es crucial porque permite evaluar eficientemente las convoluciones de tipo Coulomb y manejar la estructura natural de los gradientes Riemannianos sin depender de bases fijas (como Gaussianas).

3. Contribuciones Clave

1. Formulación en $H^1$: Es la primera vez que se presenta un marco de optimización Riemanniana para HF/DFT directamente en el espacio de Sobolev $H^1$, proporcionando una perspectiva geométricamente consistente y libre de discretización.
2. Derivación Explícita: Se derivan fórmulas cerradas para gradientes, proyecciones y retracciones utilizando operadores resolventes, evitando la necesidad de distribuciones o espacios de Sobolev negativos ($H^{-1}$).
3. Precondicionamiento Físico: Se demuestra cómo incorporar la inversión del operador de energía cinética en el gradiente conjugado, lo que es esencial para la convergencia rápida, imitando la eficiencia de los métodos SCF pero dentro de un esquema de descenso de gradiente.
4. Robustez desde Adivinanzas Aleatorias: Se demuestra que el método converge desde adivinanzas iniciales aleatorias (superposiciones gaussianas aleatorias), algo que los métodos SCF estándar a menudo no logran sin un pretratamiento cuidadoso.
5. Eliminación de Grados de Libertad: Se extiende el método a la variedad de Grassmann para eliminar la redundancia rotacional, mejorando la eficiencia en sistemas grandes.

4. Resultados Numéricos

Los experimentos se realizaron utilizando el software MRChem y la biblioteca VAMPyR con multiondas.

• Sistemas Pequeños (H2, H2He, H2Be, N2):
  • El descenso de gradiente Riemanniano simple (sin precondicionamiento avanzado) superó al método SCF-DIIS en términos de estabilidad, convergiendo desde adivinanzas aleatorias donde DIIS oscilaba o fallaba.
  • El gradiente conjugado precondicionado mostró convergencia rápida (monótona) hacia la solución de referencia.
• Sistemas Grandes y Funcionales DFT:
  • Se probaron moléculas de tamaño medio a grande (desde CO2 hasta Vitamina E y colesterol) utilizando tanto el funcional Hartree-Fock como B3LYP.
  • El método demostró robustez incluso en presencia de ruido numérico inherente a las discretizaciones adaptativas (multiondas).
  • Aunque el método SCF acelerado con KAIN (Krylov-accelerated inexact Newton) es competitivo, el descenso de gradiente Riemanniano mostró curvas de error monótonas y no divergió en casos donde SCF falló debido a criterios de parada demasiado estrictos o inestabilidades.
• Comparación: El método propuesto requiere un número de iteraciones comparable o menor que los esquemas SCF bien precondicionados, pero con una mayor robustez frente a la elección inicial.

5. Significado e Impacto

Este trabajo representa un cambio de paradigma en la optimización de estructura electrónica:

• Independencia de la Discretización: Al formular el problema en el espacio continuo $H^1$, el método es agnóstico a la base de discretización, lo que lo hace ideal para métodos modernos como multiondas, elementos finitos o bases de ondas planas.
• Fundamento Teórico Sólido: Proporciona una base matemática rigurosa para algoritmos de optimización en química cuántica, conectando la teoría de variedades Riemannianas con la física de la materia condensada.
• Potencial Futuro: Abre la puerta a desarrollar métodos de optimización más robustos para problemas de correlación fuerte y sistemas complejos donde los métodos SCF tradicionales tienen dificultades. La capacidad de converger desde adivinanzas aleatorias es particularmente valiosa para la automatización de cálculos de química computacional.

En resumen, Dinvay demuestra que la optimización Riemanniana en $H^1$, combinada con precondicionamiento físico y discretizaciones adaptativas, ofrece una alternativa robusta, eficiente y matemáticamente elegante a los esquemas de campo autoconsistente tradicionales.