MatthewRHermes
commented
1 year ago For testing
pdft_veff, I generally understand what is happening here. Thegen_g_hopgets the gradient, gradient update, hessian_update, and hessian_diag elements, right? Does g_all, in the orbital sector (ieg_all[:ngorb]) contain the Fxy−Fyx where Fxy is the generalized Fock matrix elements given theerisandget_hcore()in the current environment?
Yes. Furthermore, the remaining elements (g_all[ngorb:]) are gradients of the energy with respect to CI vector relaxations towards individual determinants:
$$\langle i|\hat{H}-E_0|0\rangle + \mathrm{h.c.}$$
Regardless, you get the
x0term which is the step in some direction. You then scale that term progressively smaller and checking to see if you get numerical convergence of the difference between the analytical and numerical derivative (summary of these lines). I am a bit confused on whatnp.dot(g_all, x1)represents here mathematically (which is conflated with what exactlyg_allis).
The product of a gradient with a step vector is an approximation to the energy difference after that step, which becomes exact in the limit of small step size with a definite constant of proportionality:
$$ \Delta f = f(1)\delta + \frac{1}{2}f(2)\delta^2 + \cdots $$
$$ \lim_{\delta\to 0} \frac{\Delta f - f(1)\delta}{\Delta f} = \frac{1}{2} \frac{f(2)}{f(1)}\delta $$
In other words, there is an asymptotic plateau where if you halve the step size, you should also halve the relative error of the linear approximation. np.dot (g_all, x1) is $f(1)\delta$ and seminum (x1) is $\Delta f$. The version of this for second derivatives is
$$ \Delta (f') = f(2)\delta + \frac{1}{2}f(3)\delta^2 + \cdots $$
$$ \lim_{\delta\to 0} \frac{\Delta (f') - f(2)\delta}{\Delta (f')} = \frac{1}{2} \frac{f(3)}{f(2)}\delta $$
but now the left and right-hand sides of the expression are both vectors. I compare them by evaluating the norm of the difference vector here.
With that being said, how should I go about testing the
pdft_feffterms numerically? I understand that you wrote theEotOrbitalHessianOperatorand test this, but is there a simpler way of testing without having to copy a bunch of code from the file? You also discuss needing to use a "dressedgen_g_hop" object (here).
I'm sure there is a simpler way to implement because my implementation is incredible spaghetti. But regarding the "dressed gen_g_hop" comment, that is a reminder that what you call the "explicit" contribution to the Hessian is mappable to a CASSCF gradient with a modified Hamiltonian. That's what that entire conditional block of the initializer is about: if incl_d2rho==False, then only the terms in the Hessian-vector product which involve $\mathbf{f}^\mathrm{ot}$ are computed.
Lastly, I am a bit confused on what exactly
paaa_onlydoes for the_ERIS(see these lines). If I understand correctly, you are adding the active space contribution of the HF potential to the HF potential core term? I know this is important since similar terms appear in thegen_g_hopcode and theEotOrbitalHessianOperator
vhf_c by default is
$$ V{pq} = v{pqii}Dii - \frac{1}{2}v{piiq}Dii = \frac{1}{2}v{pqii}D{ii} = v{pqii} = \int v{\Pi} \rho{\mathrm{core}} \frac{\partial\rho}{\partial D_{pq}} $$
If paaa_only, I hijack it to add
$$ V{ai} \mathrel{+}= v{aiuv}D{uv} - \frac{1}{2}v{avui} D{uv} = \frac{1}{2}v{aiuv}D{uv} = \int v{\Pi} \rho{\mathrm{active}} \frac{\partial\rho}{\partial D{ai}} $$
$$ V{aa} \mathrel{+}= v{aauv}D{uv} \cdots = \int v{\Pi} \rho{\mathrm{active}} \frac{\partial\rho}{\partial D{aa}} $$
$$ V{ii} \mathrel{+}= v{iiuv}D{uv+} \cdots = \int v{\Pi} \rho{\mathrm{active}} \frac{\partial\rho}{\partial D{ii}} $$
$$ V{ia} \mathrel{+}= v{iauv}D{uv} \cdots = \int v{\Pi} \rho{\mathrm{active}} \frac{\partial\rho}{\partial D{ia}} $$
since these lead to contributions to the orbital gradient that, by default in PySCF, are computed using explicit two-electron interaction elements of the appropriate index patterns (i.e., of type ppaa and papa), which we have not computed if paaa_only. I think this would actually be a quicker algorithm for doing CASSCF orbital gradients themselves, if CASSCF didn't already need ppaa and papa ERIs anyway for orbital Hessians. But anyway the whole thing is designed so that the gradient that pops out of the CASSCF gen_g_hop code is correct regardless of whether or not paaa_only. But if you want Hessians, you'll need paaa_only=False because you'll need ppaa and papa.
P.S. A bunch of factors of 1/2 are missing in the above five lines; don't take them too seriously. I'm just trying to show the orbital index patterns.
matthew-hennefarth
See how this works to use github to ask questions. So in general, I have written out the kernel functions for generating $\Delta\cdot \mathbf{F}$ which have elements of $$\left[\Delta \cdot \mathbf{F}\right]_{pq} = \int \vec{\rho}^t\Delta\mathbf{f}^\mathrm{ot}\frac{\partial \vec{\rho}}{\partial\gamma{pq}}$$ and for the 4-body term as well here. It uses some generalization of the
veffcode which I extracted/reused. I was able to test that thelazy_kernelandkernelfunctions agree using the same method as theveffcode (I used the ao->mo for the 2-electron terms). What I want to do though is test that the $\Delta\cdot \mathbf{F}$ elements agree when you numerically differentiate the $V{pq}$ and $v{pqrs}$ elements and then contract with a density.For testing
pdft_veff, I generally understand what is happening here. Thegen_g_hopgets the gradient, gradient update, hessian_update, and hessian_diag elements, right? Does g_all, in the orbital sector (ieg_all[:ngorb]) contain the $$F{xy} - F{yx}$$ where $F_{xy}$ is the generalized Fock matrix elements given theerisandget_hcore()in the current environment?Regardless, you get the
x0term which is the step in some direction. You then scale that term progressively smaller and checking to see if you get numerical convergence of the difference between the analytical and numerical derivative (summary of these lines). I am a bit confused on whatnp.dot(g_all, x1)represents here mathematically (which is conflated with what exactlyg_allis).With that being said, how should I go about testing the
pdft_feffterms numerically? I understand that you wrote theEotOrbitalHessianOperatorand test this, but is there a simpler way of testing without having to copy a bunch of code from the file? You also discuss needing to use a "dressedgen_g_hop" object (here).Lastly, I am a bit confused on what exactly
paaa_onlydoes for the_ERIS(see these lines). If I understand correctly, you are adding the active space contribution of the HF potential to the HF potential core term? I know this is important since similar terms appear in thegen_g_hopcode and theEotOrbitalHessianOperatorAny help would be appreciated, and let me know if you would rather me ask these questions via a different method of communication.