bfk - Inelastic neutron-scattering from RE ions in a crystal field including damping effects due to the exchange interaction with conduction electrons


This is based on an extension of the theory published by Klaus W. Becker, Peter Fulde and Joachim Keller in Z. Physik B 28,9-18, 1977 "Line width of crystal-field excitations in metallic rare-earth systems" (for details see appendix N). Here we give an introduction to the computer program bfk for the calculation of the neutron scattering cross section. The computer program bfk is written by J. Keller, University of Regensburg.

Description of the program:

The program calculates the dynamical susceptibility and the neutron scattering cross-section of single RE ions in the presence of crystal fields and Landau damping due to the exchange interaction with conduction electrons.

It needs the following input-files (not all are needed for all tasks)

1. A file containing the information about the RE ion: Type of ion, number of CF-levels, energy eigenvalues and eigenstates. The date are extracted from the input-file by reading the information contained in lines starting with $\char93 !$ or blanks, see the attached example.

2. File with the formfactor data for RE ion

3. File with a list of ($\vec Q,\omega$)-values, for which the calculation shall be performed

4. A parameter-file containing the names of the files with the formfactor, the table with the ($\vec Q,\omega$)-values, the energy range, scattering direction etc., see the attached example.

5. The value of the coupling constant $g=j_{ex}N(0)$, the temperature, the mode of calculation, the form of the out-put, the name of the file with the CEF-data, the name of the parameter file are provided by the commandline, which is used to start the program.

The program consists of a number of modules and subroutines which are briefly described in the following:

1. Modules CommonData, MatrixElements, FormfactorPreparation

These modules contain definitions of global variables and arrays used in the program and in different subroutines. FormfactorPreparation also contains the subroutine FormfactorTransformation which transforms an input-file with formfactor data into a file with formfactor values for equidistant Q-values. and the function Formfac to calculate the formfactor at arbitrary Q-values.

2. Subroutine ReadData

Subroutine to read-in data needed to calculate the dynamical susceptibility and the neutron scattering cross-section.

It reads the commandline, containing the coupling $g=j_{ex}N(0)$, the temperature $T$ (in Kelvin), mode of calculation (see below), form of out-put, name of the file with RE data, name of the parameter-file (containing also the name of the file with the formfactor data). The information about the RE ion is transferred into a workfile cefworkfile.dat for inspection and use in the following runs. The data contained in the parameter-file are stored in the file bfkdata.dat. The latter two have to be given only in the first run. If they are left-out in the following runs, the are assumed to be unchanged.

3. Subroutine Matrixelements

a) Calculates angular momentum matrices jjx, jjy, jjz for the crystal-field eigenstates (2-dim arrays, dimension Ns x Ns). The three directional components are also stored in the 3-dimensional array jjj(3,Ns,Ns).

b) Calculates Boltzmann-factors $p(n)$. A cut-off in the exponent $\beta
E(n)$ is introduced such that Boltzmann factors with large negative exponents are set equal to zero.

c) Defines a set of transitions $\nu$ between states n1 and n2, stored in two 1-dim arrays v1($\nu$), v2($\nu$). If both Boltzmann factors of the two states involved are zero, this transition is eliminated from the set of allowed transitions.

d) Calculates static suscepibilities $P(\mu)$ for the standard basis operators $K_{n,m}$ for the allowed transitions.

e) All these reults are stored in a file bfkmatrix.dat for examination, if something goes wrong.

4. MatrixInversionSubroutine

adapted from Numerical Recipes, to be used for the inversion of the complex matrix $\Omega_{\nu\mu}$. Called by 5.

5. Subroutine Relmatrix

Calculates the matrix relaxation function $\Phi_{\mu,\nu}(\omega)$ for the set of dynamical variables obtained from the standard basis operators for a given energy (freqency) $\omega$.

6. Subroutine Suscepcomponents

Calculates the different components of the dynamical susceptibility

\begin{displaymath}
\chi^{\alpha,\beta}(\omega) = \beta\sum_{\mu \nu}(J^\alpha)...
...\delta_{\mu\nu}-\omega \Phi_{\mu\nu}(\omega)\bigr]J^\beta_\nu
\end{displaymath}

and

\begin{displaymath}
Im \chi^{\alpha\beta}(\omega)/(1-\exp(-\beta\omega))
\end{displaymath}

7. Function Scatfunction(Q,$\omega$)

Calculates

\begin{displaymath}
S(\vec Q,\omega)= \sum_{\alpha\beta}\bigl(\delta_{\alpha\be...
...de Q^\beta\bigr) Im \chi^{\alpha\beta}/(1-\exp(-\beta\omega))
\end{displaymath}

8. Subroutine OutputResults

Here the results for the dynamical susceptibility, the scattering function and the differential neutron scattering cross section for different scattering geometries are calculated, and the results written into files bfkm.res for different scattering-modes m=0-6, which are written into the subdirectory /results. Depending on the value of ms=1,2 the new results over-write the previews results or append.

Depending on the number m=0-6 (3. entry of the commandline) the following results are calculated.

mode=0: all nine components $\chi^{\alpha\beta}(\omega)$ of the complex dynamic susceptibility are calculated for $Npoint$ equidistant energies $\omega$ between $emin$ and $emax$.

mode=1: the diagonal components of $Im \chi^{\alpha\alpha}(\omega)/\tanh{\beta\omega/2}$ are calculated and the frequency integral is compared with the sum-rule

\begin{displaymath}
\sum_\alpha {1\over\pi}\int d\omega {Im
\chi^{\alpha\alpha}\over \tanh(\beta\omega/2)}
= J(J+1)
\end{displaymath}

The sum-rule sometimes is not very well fulfilled, since within this approximation the Landau damping does not fall-off fast enough at large energies

mode=2: The scattering function

\begin{displaymath}
S(\vec Q,\omega)=\sum_{\alpha,\beta} (\delta_{\alpha, \beta...
...ta){ \chi''_{\alpha \beta}(\omega)\over 1-\exp(-\beta\omega)}
\end{displaymath}

is calculated for a given set of values for energy loss $\omega$ and scattering vector $\vec Q$ contained in a file specified in the parameterfile.

mode=3: The 9 different components of the scattering-function

\begin{displaymath}
(r_0g_JF(\vec Q)/2)^2{1\over \pi} S^{\alpha\beta}(\vec Q,\o...
...)= {Im \chi^{\alpha,\beta}(\omega)\over
1-exp(-\beta\omega)}
\end{displaymath}

are calculated.

mode= 4-6: the neutron scattering cross section

\begin{displaymath}
{d^2\sigma \over d\Omega d E'} = {k' \over
k}S(\vec Q,\omega)
\end{displaymath}

with

\begin{displaymath}
S(\vec Q, \omega)=({r_0\over 2}g F(Q))^2{1\over \pi }
\sum...
...
{Im \chi^{\alpha,\beta}(\omega)\over 1-exp(-\beta \omega)}
\end{displaymath}

is calculated for different scattering geometries: In mode 4 the direction of the wave vector $\vec k$ and the energy $E=k^2/2m$ of the incident beam is fixed. The direction of the scattering wave vector $\vec Q$ is fixed, but the length of $\vec Q$ is variable. The wave vector of the scattered particles is $\vec k'=\vec k-\vec Q$, their energy is $E'=k'^2/2m$ and the energy loss is $\omega =E-E'$. In mode 5 the direction of the wave vectors $\vec k$ and $\vec k'$ of the incoming and scattered beam are fixed, while the energy $E'$ of the scattered beam is variable. In mode 6 the energy $E$ of the incident particles is variable and the energy $E'=k^2/2m$ of the scattered particles fixed.

How to run the program:


The translated program is started with a command-line like

bfk 0.1 10 0 1 prlevels.cef paramfile.par

with the following structure:

name of the program: bfk; coupling constant g; temperature T (in K); type of calculation: mode =1...6; type of output: mst=1 overwrite, mst=2 append new results; name of file with RE ion data; name of parameter file.


The last two entries can be skipped in later runs, if they are not changed.

The mode number mode = 1 ...6 refers to the subject of calculation. The output- number mst=1,2 refers to the type of output-storage.


The file with RE data should have the form produced by sol1on (see the attached example). The lines starting with numbers or blanks contain information, the lines starting with # are commentaries, the lines starting with #! also carry information.


The parameterfile contains additional parameters needed to run the program: energy range and number of energy values. Energies of incident or scattered particles, direction of incident or scattered particles.

mode=4: E energy of incident particles, k11,k12,k13 direction if incident particles (vector with arbitrary length), k21,k22,k23 direction of scattered particles.

mode=5: E energy of incident particles, k11,k12,k13 direction if incident particles (vector with arbitrary length), k21,k22,k23 direction of scattering vector $\vec Q$.

mode =6: E energy of scattered particles, k11,k12,k13 direction if incident particles (vector with arbitrary length), k21,k22,k23 direction of scattered particles.


The parameterfile also contains the namme of a file with a list of scattering vectors $\vec Q$ and energy loss ($\omega$) needed for mode 2,3

Finally it contains the name of a file with the formfactor of the ion.


J. Keller, May 2013

Martin Rotter 2017-01-10