List of Formulas


  • The list of formulas used in FEFF

The corresponding soubroutines are quoted in the square bracket.


Parameter of angular momentum and quantum number kappa

The quantum number κ is a function of the orbital angular momentum quantum number l and the total angular momentum quantum number j.

κ = lΘ(lj) − (l + 1)Θ(jl)

where Θ(x) is equal to one if x is greater than zero and is equal to zero if x is less than zero. N.B. In the above equation the argument of Θ(x) can never equal zero since l is integral and j is half-integral.

The above equation is equivalent to:

\kappa = -j - 1/2  \qquad  \textrm{ if }\;\; j > l

\kappa = +j + 1/2  \qquad  \textrm{ if }\;\; j < l

Density of States calculations

subroutine ff2rho

if (msapp != 1)

\rho^c_{l}(E) = \rho^c_{l}(E)+\rm{Im}[\chi(E)]\rho^{sc}_l(E)\,

else \rho^c_{l}(E)\,

subroutine fmsdos

\rm{gtr}_l(E) =Tr G_{L,L'}(E)=\sum_{m=-l}^{l} \rm{gg}_{L,L}(E)

FMS calculations

XAS calculations

subroutine setkap(ihole, kinit, linit)

The contribution to the XAS from a given site i\, and orbital angular momentum l\,.

\mu_{li}(E)=\mu{li}^{0\prime}(E)[1+\chi_{li}^\prime(E)]. where \mu_{li}^{0\prime}(E) is the smoothly varying atomic background contributions and \chi_{li}^\prime(E) is the fine-structure or XAFS spectrum. A prime to denote final-state quantities calculated in the presence of a screened core hole.

subroutine ff2xmu

Adds the contributions from each path and absorber, including Debye-Waller factors. Writes down main output: chi.dat and xmu.dat cchi(E) = S02 * TrG(E)


Solving Dirac equations [afovrg]

The spin-orbit eigenfunction (eigenstates of j^2,\,j_z,\,and\, P)

\chi_{km}(\theta,\phi)=\sum_{\sigma}Y_l^{m-\sigma}(\theta,\phi)\phi^\sigma .\langle lm-\sigma\frac{1}{2}\sigma|l\frac{1}{2}jm\rangle

Matrix Elements for Electromagnetic Multipole Transitions (Grant) [xmult]

\langle nkm|\alpha . A|n'k'm'\rangle, A = upexp(iω.r), | ω | = ΔE / c where up is a polarization vector, and ω is the propagation vector.

Simplify the problem by taking ω to define the z-axis \vec{\omega} = \omega e_z, so that we can write the polarization vector u_p = \frac{1}{\sqrt{2}}\left[e_x+i p e_y\right]=-p\xi_p, p=\pm1, where ξ01 − 1 are the usual spherical basis vectors.

The decomposition of eik.r in terms of the irreducible tensors C_0^{(l)} e^{i k.r}=\sum_{l=0}^{\infty}i^l(2l+1)j_l(\omega r) C_0^{(l)}(\Omega), where the angles involved are spherical polar coordinates with ez as polar axis.

The vector operator α is a tensor of order 1. The product of two tensors can be decomposed in terms of the irreducible tensors. \alpha . A=\alpha.u_p \exp(i\omega . r)=\sum_{L=1}^\infty\sum_{l=L-1}^{L+1}a_{lL}(r)\begin{pmatrix}0&X_p^{(l1)L}\\X_p^{(l1)L}&0\end{pmatrix}, where a_{lL}(r)=(-1)^L[L]^{1/2}\,i^l\begin{pmatrix}l&1&L\\0&-p&p\end{pmatrix}j_l(\omega r)

In the form of Wigner-Eckart theorem: \langle nkm|a_{lL}(r)X_p^{(l1)L}\otimes\begin{pmatrix}0&1\\1&0\end{pmatrix}|n'k'm'\rangle=(-1)^{j-m}\begin{pmatrix}j&l&j'\\-m&p&m'\end{pmatrix}R_{k,k'}^{l,L}

In the dipole approximation it is enough to keep only one term (l = 0, L = 1).

Note: the usual notation of atomic physics [k,l,\,...]=(2k+1)(2l+1)\,....

subroutine xmult(k, kp, ls, lb, xm1, xm2)

  • input: k=k, kp=k', ls=l, lb=L\,
  • output: xm1=C_{k,k'}(1), xm2=C_{k,k'}(-1)\,

See Grant eq. 6.30. calculate the factors

\langle k|\alpha\cdot A( l, L)|k' \rangle=(-1)^{j-m}\begin{pmatrix}j & L &j'\\ -m & p & m'\end{pmatrix}R_{k,k'}

R_{k,k'}^{l,L} = \int dr (C_{k,k'}(1)\,P_k Q_{k'}+ C_{k,k'}(-1)\,Q_kP_{k'})j_l(wr)

xm1, xm2 both either real or pure imaginary

\beta=\pm 1 corresponds to the upper(lower) component of Dirac spinor.

  • calculate xm1 (β = 1 in eq.6.30)
  • calculate xm2 (β = − 1 in eq.6.30)

C_{k,k'}^{l,L}(\beta)=\sqrt{6}\left[j,L,j',\lambda,\lambda'\right]^{1/2}(-1)^{\lambda}\begin{Bmatrix}\lambda & l & \lambda' \\ \frac{1}{2} &\frac{1}{2} &1 \\ j & j' & L\end{Bmatrix}\begin{pmatrix}\lambda &l &\lambda' \\ 0 &0 & 0\end{pmatrix}\times\delta(\lambda,j-\frac{1}{2}\alpha\beta)\delta(\lambda',j'+\frac{1}{2}\alpha'\beta)

Performs radial integration for multipole matrix element or central atom absorption [radint]

if ifl = 2

\rm{xirf}=2\times\int_0^\infty dr\,dR_{k,k^\prime}^R(r) \int_0^{r} dr^\prime\,dR_{k,k^\prime}^N(r^\prime)

subroutine xrci ( mult, xm, dgc0, dpc0, p, q, bf, value)

  • output: value

r-dependent multipole matrix element (before r-integration)

dR_{k,k^\prime}^R(r)=P_k Q_{k^\prime}^R \left[C_{k,k^\prime}^{0,1}(-1)j_0(kr)+C_{k,k^\prime}^{2,1}(-1)j_2(kr)\right] +Q_k P_{k^\prime}^R\left[C_{k,k^\prime}^{0,1}(1)j_0(kr)+C_{k,k^\prime}^{2,1}(1)j_2(kr)\right]

dgc0*q* (xm(2)*bf(0) + xm(4)*bf(2)) + dpc0*p * (xm(1)*bf(0) + xm(3)*bf(2))

for double radial integral (use the irregular components of the final state)

dR_{k,k^\prime}^N(r)=P_k Q_{k^\prime}^N \left[C_{k,k^\prime}^{0,1}(-1)j_0(kr)+C_{k,k^\prime}^{2,1}(-1)j_2(kr)\right] +Q_k P_{k^\prime}^N\left[C_{k,k^\prime}^{0,1}(1)j_0(kr)+C_{k,k^\prime}^{2,1}(1)j_2(kr)\right]

Radial Integration Routine

csomm (dr,dp,dq,dpas,da,m,np)

  • output: da

\rm{da}=\int_{r=0}^{dr(np)} (\rm{dp}+\rm{dq})dr^m Modified to use complex p and q. SIZ 4/91

integration by the method of simpson of (dp+dq)*dr**m from 0 to r=dr(np)

dpas=exponential step; for r in the neighborhood of zero (dp+dq)=cte*r**da

r=\exp(-8.8+.05n), dr=.05\times r dn

Wigner 3j symbol [cwig3j]

J=j_1+j_2:\quad \begin{pmatrix} j_1 & j_2 & J\\ m_1 & m_2 & M\end{pmatrix}\equiv\frac{(-1)^{j_1-j_2+M}}{\sqrt{2J+1}}\langle j_1 j_2 m_1 m_2|JM\rangle

function cwig3j (j1,j2,j3,m1,m2,ient)

  • input: j1=j_1, j2=j_2, j3=J, m1=m1, m2=m2, ient\,
  • output: wigner 3j coefficient for integers (ient=1) or semiintegers (ient=2) other arguments should be multiplied by ient

Orthogonality relations

\sum_{m_1=-j_1}^{+j_1}\sum_{m_2=-j_2}^{+j_2}\begin{pmatrix}j_1 & j_2 & j_3\\ m_1 & m_2 & m_3\end{pmatrix}\begin{pmatrix}j_1 & j_2 & j_3'\\ m_1 & m_2 & m_3'\end{pmatrix}=\frac{1}{2j_3+1}\,\delta_{j_3j_3'}\delta_{m_3m_3'}

Composition relation for the spherical harmonics

\int Y_{l_1}^{m_1}Y_{l_2}^{m_2}Y_{l_3}^{m_3} d\Omega = \left[\frac{(2l_1+1)(2l_2+1)(2l_3+1)}{4\pi}\right]^{\frac{1}{2}}\begin{pmatrix}l_1 & l_2 & l_3\\ 0 & 0 & 0\end{pmatrix}\begin{pmatrix}l_1 & l_2 & l_3\\ m_1 & m_2 & m_3\end{pmatrix}

6j symbols [sixj]

J=j+j'+j'':\quad a)\,j'+j=g',\,g'j''=J\quad b)\,j+j''=g'',\,j'+g''=J\quad

\langle j',(jj'')g'';JM|(j'j)g',j'';J'M'\rangle=\delta_{JJ'}\delta_{MM'}\sqrt{(2g'+1)(2g''+1)}(-)^{j+j'+j''+J}\begin{Bmatrix}j' & j & g' \\ j'' & J & g''\end{Bmatrix}

9j symbols [ninej]

J=j_1+j_2+j_3+j_4:\quad a)\,j_1+j_2=J_{12},\,j_3+j_4=J_{34},\,J_{12}+J_{34}=J\quad b)\,j_1+j_3=J_{13},\,j_2+j_4=J_{24},\,J_{13}+J_{24}=J\quad

\langle (j_1 j_2)J_{12},(j_3 j_4)J_{34};JM|(j_1 j_3)J_{34},(j_2 j_4)J_{24};J'M'\rangle=\delta_{JJ'}\delta_{MM'}\sqrt{(2J_{12}+1)(2J_{34}+1)(2J_{13}+1)(2J_{24}+1)}\begin{Bmatrix}j_1 & j_2 & J_{12} \\ j_3 & j_4 & J_{34}\\ J_{13} & J_{24} & J\end{Bmatrix}

developer's resources