ELF doesn't rely on third-order derivatives of electron density.

It is worth to note that Multiwfn is able to analytically calculate 3-order derivative tensor of orbital wavefunction, see "subroutine orbderv" in function.f90. Because electron density can be easily evaluated based on wavefunctions of occupied orbitals, realizing 3-order derivative of electron density in Multiwfn is straightforward.

In fact, in the "subroutine calchessmat_lapl" in function.f90, the 3-order derivative tensor of electron density (the array "rhotens3") is just the intermediate quantity for evaluating gradient of Laplacian of electron density, you can directly extract out this part of code.

Dear Tian Lu, best regards.
We are interested in the computing of an electron localization function (LF) within the Multiwfn code. Such LF depends on  Third-order derivatives of the density and also orbital basis. My question is how can I incorporate the calculation of these derivatives in a suitable way within the Multiwfn code?