[CP2K-user] [CP2K:13866] A question about a formula in CP2K source code

[Junting]

Dear Matthias,

Thank you so much for providing me with this insightful work!

However, Eq (3.31) is in Two-electron Integral part, and seems not to 
involve error function. The <a|erf(r)/r|b> term is one-electron integral, 
describing the interaction between electron and core, instead of eletron 
repulsion. I still do not understand how these two are related.

And from intuition, if zetc goes to infinity, ferf should become the same 
as fnuc. But according to source code, they seem quite different, which 
also makes me confused.

Best,

Junting

在2020年9月14日星期一 UTC+8 下午4:15:43<Matthias Krack> 写道：

> Dear Junting
>
>  
>
> Did you compare your results with the work of McMurchie and Davidson 
> <https://www.sciencedirect.com/science/article/pii/002199917890092X> from 
> 1978? The expression for lambda in eq. (3.31) seems to match the 
> implemented ferf for the basic integral.
>
>  
>
> Best 
>
>  
>
> Matthias
>
>  
>
> *From:* c... at googlegroups.com <c... at googlegroups.com> *On Behalf Of *
> Junting
> *Sent:* Montag, 14. September 2020 05:13
> *To:* cp2k <c... at googlegroups.com>
> *Subject:* [CP2K:13866] A question about a formula in CP2K source code
>
>  
>
> Dear CP2K users and developers,
>
>  
>
> I am reading the source code of CP2K, and I have got lots of insights 
> about formulas and algorithms from it. However, I have encountered a 
> question when reading the aobasis/ai_verfc.F file in src directory.
>
>  
>
> The ai_verfc module is used for calculating the matrix <a|erfc(r)/r|b> in 
> all-electron calculation. This term is divided into a nuclear term 
> <a|1/r|b> and an error function term <a|erf(r)/r|b>. In the code prefactors 
> corresponding these two terms are calculated by 
>
>  
>
> fnuc = 2.0_dp*pi*zetp*f0
>
> ferf = 2.0_dp*SQRT(pi**5*zetw)*zetp*zetq*f0
>
>  
>
> and the basic s-orbital integrals are calculated by prefactors and 
> incomplete gamma function
>
>  
>
> t = rcp2/zetp
>
> CALL fgamma(nmax-1, t, f)
>
> vnuc(1, 1, n) = fnuc*f(n-1)
>
>  
>
> t = -f4*rcp2/zetp
>
> CALL fgamma(nmax-1, t, f)
>
> verf(1, 1, n) = ferf*f(n-1)
>
>  
>
> Then the integrals over p or higher orbitals are calculated by a recursion 
> procedure as is discussed in the literature: S. Obara and A. Saika, J. 
> Chem. Phys. 84, 3963 (1986).
>
>  
>
> My question is that according to my derivations, my result for fnuc is 
> consistent with the code, but my ferf is 
>
> ferf = 2.0_dp*pi*zetp*f0*(-f4)**(n-0.5)
>
> which is not consistent with the code, and is dependent on the parameter n 
> in the incomplete gamma function f(n).
>
>  
>
> And it is reasonable that on the limit of zetc going to infinity, the 
> error function term <a|erf(r)/r|b> should be the same as the nuclear term 
> <a|1/r|b>. But according to source code, it seems not to be the case.
>
>  
>
> Therefore I want to ask, is this a potential mistake in the source code, 
> or did I miss anything else in the derivations? I really appreciate it if 
> someone could help me solve this question.
>
>  
>
> Best wishes,
>
> Junting
>
> -- 
> You received this message because you are subscribed to the Google Groups 
> "cp2k" group.
> To unsubscribe from this group and stop receiving emails from it, send an 
> email to cp... at googlegroups.com.
> To view this discussion on the web visit 
> https://groups.google.com/d/msgid/cp2k/89673088-6dba-4f9f-9cf8-464569253964o%40googlegroups.com 
> <https://groups.google.com/d/msgid/cp2k/89673088-6dba-4f9f-9cf8-464569253964o%40googlegroups.com?utm_medium=email&utm_source=footer>
> .
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <https://lists.cp2k.org/archives/cp2k-user/attachments/20200914/e0aface0/attachment.htm>