Programs for Specific Tasks in the World of McPhase

This section contains some further programs which do some important manipulations in magnetism, ...

addj file1.j file2.j:
adds exchange parameters of file2.j to the parameters in file1.j - output is written to stdout, if file2.j is not given, a copy of the input file is saved

Option: -nofcomponents 23
fixes the nofcomponents to 23 by reducing (removing entries) or increasing (by filling with zeroes) the exchange parameter tables
Option: -ni
forces output without indexchange
Option: -sp 0 0 0.1
shift all atomic positions in file 1 by da=0 db=0 dc=0.1
Option: -s 0.2
scales all interactions (including magnetoelastic G) in file1 by 0.2 before adding
Option: -sJ 0.2
scales all interactions (but not magnetoelastic G) in file1 by 0.2 before adding
Option: -rmcomp 9 15
removes components 8 to 15 from result before output
Option: -pd
output a column with distance
Option: -ps
output a column with sublattice index
Option: -pi 5 7 1 2 1 3
print to stdout only interaction tensor (rows 1-2, columns 1-3) of ion 5 with neighbour 7 as matrix
Option: -pG
print to stdout only magnetoelastic G matrices (no interactions)
Option: -v
verbose
If file2.j is not given, a copy of the input file is saved
anisotropy:
program to calculate the magnetic anistropy by doing a mcphas or single ion calculation for different external magnetic field directions and evaluating the expectation value of the magnetic moment.
     usage: anisotropy -h
           anisotropy  T H xn yn zn nofsteps [-r sipffilename Hxc1 Hxc2 ... Hxcnofcomponents]
           anisotropy  T H -p nofthetasteps [-r sipffilename Hxc1 Hxc2 ... Hxcnofcomponents]

     -h           : this (help) message
      T           : temperature in Kelvin
      H           : absolute value of the external magnetic field (T)
      xn,yn,zn    : direction normal to plane, in which the anisotropy
                    should be calculated ... e.g. if you want to
                    calculate the anisotropy in the xy plane, then
                    enter xn yn zn = 0 0 1
      nofsteps    : number of steps to be calculated 
     -p           : calculate polycrystal average
      nofthetastepssteps    : number of theta steps to be calculated for polycrystal average

    option:
    -r sipffilename: filename of single ion parameter file
                      Hxc1,Hxc2,... are the exchange field components (meV)
                     (exchange field is kept constant, external magnetic
                     field is rotated in the anisotropy calculation)

    output files:

    ./results/anisotropy.out  contains anisotropy information
cfsplit:
Calculates via group theory, the multiplicity of crystal field split levels for particular point symmetries. The command is: "cfsplit $\langle$PTGP$\rangle$ $\langle$J$\rangle$" where $\langle$PTGP$\rangle$ is the name of the point group (both Hermman-Maguin and Schoenflies notation are understood), and $\langle$J$\rangle$ is the half-integer value of the total angular momentum J of the lowest multiplet. E.g. "cfsplit Oh 2" gives the familiar cubic Eg/T2g splitting (last line of output). See also: symhmltn which calculates the allowed two-ion exchange terms for a point group.
cif2mcphas [options] $\langle$INPUT.CIF$\rangle$:
Creates a set of McPhase input files (mcphas.j, sipfs) from the crystal structure given in a CIF (Crystallographic Information File), with or without user interaction. CIF input files may be obtained from online repositories such as the Inorganic Crystal Structures Database (ICSD), the Crystallography Online Database (COD) or the American Mineralogist Database (AMCSD), or they may be generated by many crystal structures programs, including FullProf. In fully automated mode, the program is invoked with the name of the CIF input file as: "cif2mcphas $\langle$INPUT.CIF$\rangle$". In this case, the program tries to guess which crystallographic sites contains magnetic or non-magnetic ions, and what the valence (oxidation state) of that ion is. This information may also be given by the user in interactive mode, with the following command: "cif2mcphas $\langle$INPUT.CIF$\rangle$ -i" (e.g. by including the switch -i (or –interactive) before or after the CIF input filename). The online help may be access by invoking the program without arguments or with the -h (–help) switch. Finally, if you do not have a CIF of the structure you want to calculate, you can use the template here:

; Please replace text in square brackets [] with your data
_cell_length_a                     [a_in_Angstrom]
_cell_length_b                     [b_in_Angstrom]
_cell_length_c                     [b_in_Angstrom]
_cell_angle_alpha                  [alpha_in_degrees]
_cell_angle_beta                   [beta_in_degrees]
_cell_angle_gamma                  [gamma_in_degrees]
; The H-M symbol generally needed by cif2mcphas but some spacegroups
; which have only one setting can be identified purely by the number.
; Most orthorhombic or low symmetry spacegroups have multiple centring,
; origin choices or unique axes which will require an H-M symbol.
; If this is the case, and you give just the spacegroup number,
; cif2mcphas may return with an error.
; If you don't use one of these fields, you must comment it out, by
; putting a ";" or "#" at the start of the line.
_symmetry_space_group_name_H-M     [Hermann-Maguin_symbol_with_spaces]
_symmetry_Int_Tables_number        [SpacegroupNumber]
; At the end of the file, please list the non-equivalent
; sites' fractional coordinates of this structure.
loop_
_atom_site_label
_atom_site_fract_x
_atom_site_fract_y
_atom_site_fract_z
[AtomSite] [x_coord] [y_coord] [z_coord]

replacing instances of e.g. [AtomSite] or [x_coord] with the site name and fraction $x$-coordinate, etc. This template can also be generated using the command "cif2mcphas -c" or "cif2mcphas –create". See also: icsdread which downloads and creates a CIF from an entry of the (free) Korean ICSD mirror.

For example, this template may be filled with the following information

_cell_length_a                     3.4891(4)
_cell_length_b                     4.4898(6)
_cell_length_c                     6.0033(8)
_cell_angle_alpha                  90
_cell_angle_beta                   90
_cell_angle_gamma                  90
_space_group_crystal_system        'orthorhombic'
_symmetry_Int_Tables_number        38
_symmetry_space_group_name_H-M     'A m m 2'
loop_
  _atom_site_label
  _atom_site_type_symbol
  _atom_site_fract_x
  _atom_site_fract_y
  _atom_site_fract_z
 Tm1 Tm 0.5000 0.0000   0.3850(2)     
 Ni1 Ni 0.0000 0.0000  -0.0034(3)      
 C1 C   0.0000 0.349(2) 0.1851(18)     
; C-C bond charge as 0 0.5 z as C1
 E1 E 0.0000 0.500(2) C1
; atomic radii  Tm 227 pm C 170pm  Ni 163pm
; C-Ni bonds at Ni+(C1-Ni)*163/(163+170)
 E2 E Nix+(C1x-Ni)*163/(163+170) Ni+(C1-Niy)*163/(163+170)  Ni+(C1-Ni)*163/(163+170)

Note that atomic coordinates may also be expressed as formulas, with positions from previous lines as variables. In the above example the bond charge of bonding electrons is denoted by the symbol E. The position coordinate is defined relative to the position of C1 and Ni. The indices x,y and z refer to the corresponding fractional coordinate and can be omitted if the coordinate to be computed is the same.

Options:

-h
–help
prints help message
-c
- - create
create cif template without any structure information
-i
- - interactive
no automatic guess of pointcharges, magnetic configurations etc.
-pc 3.6
do point charge calculation including charges up to 3.6 Angstroms.
-ch
- - charges
input charges to override defaults, e.g. -pc 5.6 -ch O=-2,Ce=+3,Pd=+2
-sp
- - savepcfile
save the .pc files with list of neighbouring ions in results/ (default is to not save these but to use coordinates internally only). If -pc is not specified, this option has no effect.
-sc
- - savecharges
store pointcharge coordinates in sipf file
-s
- - supercell
creates a super-cell of size #x#x#. (e.g. -s 2x3x4).
-rp
- - readpcfile
don't regenerate point charge coordinates but reread them from the saved .pc files in results/ - this option is intended so that users can modify the charges in the .pc file (e.g. by the program setvariable) and regenerate the CF parameters from the new charges. (You don't have to specify -pc with this option).
-nm
- - nonmagnetic
set nonmagnetic ion overrideing default, e.g. -nm Cu
-m
- - magnetic
set magnetic ion overrideing default, e.g. -m Cu

-so
- - so1ion
force all ions to use so1ion
-ic
- - ic1ion
force all ions to use ic1ion
-ph
- - phonon
force all ions to use phonon
-sf
- - screenfile
apply screening function to charges for calculation of crystal field parameters, e.g. -sf sf.r reads file sf.r, with column 1 r(Å), column 2 screeningfactor for $B_2^m$, column 3 screeningfactor for$B_4^m$, and column 4 screeningfactor for $B_6^m$.

Default is to create sipf files with so1ion module for f-electrons and ic1ion module for d-electrons and phonon module for nonmagnetic ions. The -so and -ic options have no effect unless either -pc or -rp is specified. -so takes precedent over -ic (e.g. when typing cif2mcphas -rp -so -ic file.cif the program will ignore the -ic and force all ions to use so1ion).

Some examples of the syntax are:

cif2mcphas -pc 3.5 ermno3.cif
cif2mcphas -pc 3.5 -sp ermno3.cif
cif2mcphas -rp ermno3.cif
cif2mcphas -pc 3.5 -so ermno3.cif
cif2mcphas -rp -ic ermno3.cif
cpsingleion 10 100 1 file.levels.cef [options]:
By runnning singleion some file.levels.cef is created in folder results. The cpsingleion program may then be used. It calculates the specific heat in the temperature interval 10-100 K with a step width of 1 K. Alternatively a comparison to experimental data can be made by cpsingleion 1 2 cpexp.dat file.levels.cef, where the temperatures are given in column 1 and the experimental specific heat in column 2 of file cpexp.dat. The calculated specific heat is compared to the experimental data and a standard deviation sta is calculated and output is written to stdout. Other quantities can be calculated using the options: -s (calculate entropy (J/molK) instead of cp), -f (calculate free energy (J/mol) instead of cp),-u (calculate magnetic energy (J/mol) instead of cp), -z (calculate partition sum instead of cp)
cpso1ion 10 100 1 [options]:
same as cpsingleion but for output of program so1ion (no file.levels.cef required, takes values from so1ion.out) .
cpic1ion 10 100 1 [options]:
same as cpsingleion but for output of program ic1ion (no file.levels.cef required, takes values from ic1ion.out) .
cpicf1ion 10 100 1 [options]:
same as cpsingleion but for output of program icf1ion (no file.levels.cef required, takes values from icf1ion.out) .
epsdebye Tmax dT Tdebye scale [d1 d2 datafile]:
calculates the phonon induced strain $\epsilon$ using the debye model according to the following formula: $\epsilon=S*T*D(\Theta_{D}/T) $ with $D(z)=\frac{3}{z^3}\int_0^z \frac{x^3 dx}{e^x-1}$ Range is from zero to Tmax in stepwidths dT unless a datafile is given. If a datafile is given, with data column d1 and d2,the strain is calculated for T-values of data column d1 and epsilon is compared to data in column d2 - a standard deviation sta is calculated as a sum of squared deviations. As output the datafile is given, an additional is column added containing the calculated strain epsilon. The datafile has to be sorted according to descending T values !!! output is written to stdout.
extendunitcell n1 n2 n3
program to extend crystallographic unit cell n times in r1 (or r2,r3) direction, meaning take mcphas.j, mcphas.tst and mcdiff.in and generate an extended description of the unit cell n1xr1,n2xr2,n3xr3 put result into extend.j, extend.tst and extend.in
fermicol col T filename
calculates the Fermifunction from energy in column col.
col
column containing energy values (eV) relative to EF
T
temperature (K)
filename
file name
fitfermi T EF fwhm min max filename
fits a (Gaussian convoluted) Fermi function to data in file
T
Temperature (K)
EF
initial value of Fermi Energy (eV)
fwhm
initial vlaue for resolution (eV) (if less than zero the fwhm is not fitted and set to |fwhm|)
min max
energy range of fit (may be less than range of experimental data points)
filename
filename (col 1 is energy in eV and col 2 intensity)
The fermifunction is defined as $f(E)=b+l(E-EF)+[d+k(E-EF)]/(\exp((E-EF)/kT)+1)$, the function $f(E)$ is convoluted with a Gaussian function of given fwhm and the result is compared to experimental data.

output: files can be found in directory results, filename.fit is created with fitted function and parameter values

formfactor *.sipf
program to calculate the neutron magnetic formfactor from the formfactor coefficients in the single ion parameter file *.sipf. If a radial wave function is given in *.sipf then the formfactorand the expectation values of the spherical Bessel functions are evaluated by integration with this radial wave function (see appendix J).
icsdread $\langle$ICSD_ID$\rangle$ $>$ INPUT.CIF
Program to create a CIF (Crystallographic Information File) from the contents of a webpage on the (free) Korean Inorganic Crystal Structures Database (ICSD) mirror (http://icsd.kisti.re.kr/) which is automatically downloaded when the corresponding ICSD ID is given as input. This requires an internet connection. If the webpage corresponding to this ICSD entry was previously downloaded, the filename can be given as input instead and the program will read from this file. Output is sent to the console so should be redirected into a file using the $>$ operator for futher use, for example to set up McPhase input files using cif2mcphas

jjj2j:
transforms file of mcphas.jjj format to mcphas.j format - output is written to stdout

makenn 23.3[options]:
Program to calculate neighbors of atoms given a crystal structure. Note that in order to use makenn you have to set up a working mcphas.j file with the crystal structure. The program makenn takes mcphas.j and creates all neighbours within a sphere of distance 23.3Å, for every neighbour the classical dipole interaction $J_{\alpha\beta}(\vec R =\vec R_j - \vec R_i)=\frac{\mu_0}{4\pi} g_{J_i} g_{J_j} \mu_B^2 \frac{3R_{\alpha}R_{\beta}-\delta_{\alpha\beta}R^2}{R^5}$ (the factor $\frac{\mu_0}{4\pi}$ is necessary in SI units) is calculated and is stored in file makenn.j. If the exchange parameters (and neighbour positions) are not known for your system, you can use this module to generate a list of nearest neighbours and exchange parameters. Currently implemented are not only dipolar interactions, but also exchange interactions via the Bethe-Slater curve or the RKKY model.
option -rkky A(meV) kf(1/A)
calculates the rkky interaction according to $J(R)=A %%@
\cos(2 k_f R)/(2 k_f R)^3$
option -rkky3d A(meV) kx(1/A) ky(1/A) kz(1/A)
calculates the rkky interaction according to $J(R)=A cos(2 \kappa)/(2 \kappa)^3$ with $\kappa^2=k_x^2 R_x^2 + k_y^2 R_y^2 + %%@
k_z^2 R_z^2$
option -rkkz A(meV) kf(1/A)
calculates the rkky interaction according to $J(R)=A %%@
(\sin(2 k_f R)- 2 k_f R \cos(2 k_f R))/(2 k_f R)^4$, scaling A>0, kf should be the Fermi wavevector (see 5.7.35 in [1])
option -rkkz3d A(meV) kx(1/A) ky(1/A) kz(1/A)
calculates the rkky interaction according to $J(R)=A (\sin(2 \kappa)- 2 \kappa \cos(2 \kappa))/(2 \kappa)^4$ with $\kappa^2=k_x^2 %%@
R_x^2 + k_y^2 R_y^2 + k_z^2 R_z^2$
option -kaneyoshi A(meV) D(A) alpha
calculates the Kaneyoshi parametrisation for the Bethe-Slater curve: $J(R)= A [-(R/D)^2+(R/D)^4]exp[-\alpha (R/D)^2]$ with $D$ corresponding to the orbital radius
option -kaneyoshi3d A(meV) Dx(A) Dy(A) Dz(A) alpha
calculates the Kaneyoshi parametrisation for the Bethe-Slater curve: $J(R)= A [-\rho^2+\rho^4]\exp[-\alpha \rho^2]$ with $\rho^2=R_x^2/D_x^2+R_y^2/D_y^2+R_z^2/D_z^2$
option option -bvk filename
for phonon take Born van Karman model with longitudinal and transversal spring constants from file - file format:
#  atom_i_sipf atom_j_sipf bondlength(A) long(N/m) trans(N/m) 
Ce1.sipf  Ce1.sipf  +4.0 200.9 100.0
Ce1.sipf  Ce1.sipf  +4.7 70.9 0.0
mind: into MODPAR2-6 in *.sipf the Einstein-oscillator parameters are written, too. Longitudinal/Transversal springs are described by $c_L$/$c_T$, respectively and the energy is given by
$\displaystyle H$ $\textstyle =$ $\displaystyle \sum_{ij} \frac{c_L(ij)-c_T(ij)}{2\vert\mathbf \mathbf R_{ij}\ver...
...thbf P_j . \mathbf R_{ij})^2
+ \frac{c_T(ij)}{2} (\mathbf P_i - \mathbf P_j )^2$(161)
Fig. 18 shows the Born-van-Karman model definition of longitudinal springs $c_L$ and transversal springs $c_T$. Omit filename to create a sample file with longitudinal springs calculated according to $C_L(ij)=25 \exp(-0.1*r^2_{ij}/\AA^2)$ N/m. Output: file makenn.Cel is created containing just the elastic constants.

option option -cfph [screeningfile.r] [-r]
calculate crystal field phonon interaction: mcphas.j lists magnetic and non magnetic atoms with charges defined in the sipf files by CHARGE= variable. For magnetic atoms the sipf file the variable MAGNETIC=1 has to be set and information about the ion has to be present (IONTYPE etc.). For each magnetic ion a new site is created resembling the magnetic electron charge cloud and this new site is shifted 0.1 Å along $c$ in order to not overlap with the original site. For option [-r] the original site sipf file is replaced automatically to use an sipf file with the MODULE=phonon similar to all the other non magnetic sites (which should have a MODULE=phonon). For the new magnetic site the program pointc is used by makenn with option -d to calculate derivatives $\partial B_{lm}/\partial u$ which are inserted as interaction parameters between MODULE=phonon and MODULE=so1ion sites. In order to use the resulting file resultsmakenn.j a phonon model has to be set up, and the phonon model has to be added to makenn.j, e.g. by program addj. Moreover, magnetic sites sipf files are required, which are created in results/makenn.a*.sipf. Note: A screening file can be used to define distance dependent screening of charges for the pointcharge model calculation format: col1 distance r (Å), col 2 screening factor for $B_2^m$, col 3 for $B_4^m$ and col 4 for $B_6^m$
option option -f [filename]
option option -dm [filename]
read interaction constants from table in file. Use -f for isotropic interactions between momentum $\mathbf J_i$ and $\mathbf J_j$ at positions $i$ and $j$
$\displaystyle \mathcal J \mathbf J_i . \mathbf J_j$ $\textstyle =$ $\displaystyle \mathbf J_i.
\left( \begin{array}{ccc}
\mathcal J & 0 & 0\\
0 & \mathcal J & 0 \\
0 & 0 & \mathcal J \\
\end{array}\right) .\mathbf J_j$(162)
and -dm for Dzyaloshinski Moriya interactions:
$\displaystyle \mathcal D (\mathbf J_i \times \mathbf J_j )$ $\textstyle =$ $\displaystyle \mathbf J_i.
\left( \begin{array}{ccc}
0 & \mathcal D_z & -\mathc...
...al D_x \\
\mathcal D_y & -\mathcal D_x & 0 \\
\end{array}\right) .\mathbf J_j$(163)
To get an sample file use option -f or -dm without a filename.
option -d
puts to the last column the distance of the neighbours (A)

The neigbours of each atom are also stored in seperate files results/makenn.a*.pc, which can be used with the program pointc to evaluate the pointcharge model and calculate crystal field paramaters.

option -npc
do not keep results/makenn.a*.pc files
option -rfunc
create results/makenn.func.r file with distance dependence of rkky and kaneyoshi functions (only for options rkky* and kaneyoshi*)
option -djdx
option -djdy
option -djdz
create files makenn.djdx .djdy .djdz instead of makenn.j, respectively these contain the derivatives of the interaction paramters with respect to displacement of the neighbor in x,y and z direction respectively. These derivatives are useful for the calculation of exchange striction effects.

option -djdeps1
option -djdeps2
option -djdeps3
option -djdeps4
option -djdeps5
option -djdeps6
create files makenn.djdeps1 .djdeps2 ... .djdeps6 instead of makenn.j, respectively these contain the derivatives of the interaction paramters with respect to the strain tensor components epsilon1, epsilon2 ... epsilon6 (in Voigt Notation) as calculated from the kf dependence on strain in the rkky formulas. For a spherical Fermi surface $k_F=(6\pi^2/V)^{1/3}$ and thus $\partial k_f/\partial \epsilon_1=-k_f/3 $. for option
-rkky
$\partial {\mathcal J}/\partial \epsilon_{\alpha} = A[3\cos(2k_fR)+2k_f R \sin(2k_fR)]/24(k_fR)^3$ for $\alpha=1,2,3$.
-rkkz
$\partial {\mathcal J}/\partial \epsilon_{\alpha} = A[(1-(k_fR)^2)\sin(2k_fR)-2k_fR.\cos(2k_fR)]/(12(k_fR)^4)$ for $\alpha=1,2,3$.

For an orthorhombic distorted Fermi sphere the volume is $4\pi k_x k_y k_z/3=8\pi^3/V$ and it follows $dk_x/k_x+dk_y/k_y+dk_z/k_z=V d(1/V)=-dV/V$. Here we have assumed that the strain will not change the effective electron mass and band structure significantly and the change in $k_x,k_y,k_z$ will be only due to the change in electron density which is correlated to a change density of electron states in k-space. In this case $dk_x/k_x=dk_y/k_y=dk_z/k_z=-(d\epsilon_1+d\epsilon_2+d\epsilon_3)/3$ and $\partial k_x/\partial \epsilon_{\alpha}=-k_x/3$ for ${\alpha}=1,2,3$ and similar for $k_y$ and $k_z$. We calculate $\partial J/\partial \epsilon_1=\partial J/\partial \kappa \partial \kappa/\partial \epsilon_1$ and $\partial \kappa/\partial \epsilon_1=\partial \kappa/\partial k_x \partial k_x/\...
...
+\partial \kappa/\partial k_z \partial k_z/\partial \epsilon_1=
= -\kappa/3 $ and similar for $\epsilon_2, \epsilon_3$. Thus for options

-rkky3d
$\partial {\mathcal J}/\partial d\epsilon_{\alpha} = A[3 \cos(2\kappa) + 2\kappa \sin(2\kappa)]/(24\kappa^3)$ for $i=1,2,3$.
-rkkz3d
$\partial {\mathcal J}/\partial d\epsilon_{\alpha} = A[(1-\kappa^2 )\sin(2\kappa)-2\kappa \cos(2\kappa)]/(12\kappa^4) $ for $i=1,2,3$.

These derivatives are useful for the calculation of exchange striction effects.

Note: $xyz$ refers to a right handed Euclidean coordinate system with $y\vert\vert\vec b$, $z\vert\vert(\vec a \times \vec b)$ and $x$ perpendicular to $y$ and $z$.

mcphas2jvx mcphas.j
Creates a JavaView input file (by default results/mcphas.jvx) which shows the atomic positions and exchange interactions in 3D graphics. It requires a mcphas.j type file as input. JavaView can then be used to view the output in 3D using: javaview results/mcphas.jvx. A specific output file may be given using the -o (–output) switch, and for the cluster module the switch -i (–individual) can be used to plot individual atoms within a cluster and the intra-cluster exchange interactions. For example, to create the a non-default file with atomic positions of a cluster calculation for viewing with JavaView, mcphas2jvx mcphas.j -i -o results/cluster.jvx && javaview results/cluster.jvx where the double ampersand means to execute the second (javaview) command if the first is successful.

Figure 31: Example of output of mcphas2jvx+JavaView for a normal setup (left, Bi$_2$Fe$_4$O$_9$ with Fe-Fe exchange interactions shown as lines) and for a cluster setup (right, trimers in LuMnO$_3$ shown as triangles with blue square block indicating trimer centre; intra-trimer exchange shown as purple lines, inter-cluster exchange as black lines). Linewidths are proportional to the strength of the exchange interactions between the linked atoms or clusters.
\includegraphics[angle=0, width=0.4\textwidth]{figsrc/mcphas_jvx.eps} \includegraphics[angle=0, width=0.4\textwidth]{figsrc/trimer_jvx.eps}
pointc[options] ionname$\vert$sipffile charge_and_position$\vert$file.pos:
Program to calculate Crystal field Parameters from Point Charges
example 1: pointc Ce3+ 0.2 4 1 5.3
... calculate Blm (Stevens Parameters) and Llm (Wybourne Parameters) for one point charge of +0.2$\vert e\vert$ in distance x=4 Åy=1 Åz=5.3 Åfrom a Ce$^{3+}$ ion. See equation (3) and appendix E for formulas.
example 2: pointc Ce3+ file.pos
... read several charges+coordinates from file.pos,file format: column 1=charge, column 2-4 = x y z coordinate. (note,program makenn creates useful files for this option from the crystal structure).
example 3: pointc Ce3+ C2.pos 5 6
... same as example 2 but reduced charge model,i.e. B2m calculated, with charges in col 1 of C2.pos,B4m and B6m with charges in col 5 and 6, respectively. Note: if an ion is not implemented, it's parameters can be entered in a single ion property file and pointc is started as
example 4: pointc file.sipf 0.2 4 1 5.3 0.1 0.3
... will use $0.2\vert e\vert$ for $B_2^m$, $0.1\vert e\vert$ for $B_4^m$, $0.3\vert e\vert$ for $B_6^m$. ... the first line of the single ion property file.sipf must be
 #!MODULE=so1ion
The single ion property file must then contain the following information (# denotes comments):
 # file.sipf should contain the following information (# denotes comments):
 # the name of the ion
 IONTYPE=Ce3+
 #Stevens parameters (optional, necessary for output of Blm)
 ALPHA=-0.0571429 BETA=0.00634921 GAMMA=0
 # the radial matrix elements RN=<r^N> in units of a0^N (a0=0.5292 A)
 R2=1.309 R4=3.964 R6=23.31
 #optional radial wave function parameters, for transition metal ions the the values
 #are tabulated in Clementi & Roetti Atomic data and nuclear data tables 14 
 #(1974) 177-478, the radial wave function is expanded as
 # R(r)=sum_p Cp r^(Np-1) exp(-XIp r)(2 XIp)^(Np+0.5)/sqrt(2Np!)
 #rare earth:Freeman&Watson PR127(1962)2058,Sovers J.Phys.Chem.Sol.28(1966)1073
 #e.g. Co2+ is isoelectronic to Fe+, looking at page 422
 #of Clemente & Roetti the parameters are 
 N1=3 XI1=4.95296 C1=0.36301 
 N2=3 XI2=12.2963 C2=0.02707 
 N3=3 XI3=7.03565 C3=0.14777

OUTPUT:

  • stdout ... sipf file with CEF pars,radial matrix elements,Stevens factors, pointcharges and positions (omit with option -o)
  • results\pointc.out ...contains results of convergence when summing up contributions of different neighbours one by one...
  • results\pointc.Blm ... Crystal field parameters $B_l^m$ in Stevens Notation
  • results\pointc.Llm ... Crystal field parameters $L_l^m$ in Wybourne Notation
  • results\pointc.dBlm .. for option -d derivatives $dB_l^m/du$
  • results\pointc.dLlm .. for option -d derivatives $dL_l^m/du$, derivatives are with respect to $u=$neighborposition$(\AA)/a_0$
powdercell2j file:
used to create mcphas.j type file from output of powdercell,output is written to mcphas.j. Example of input file:
  No   name       crystal coordinates          cartesian coordinates
                 x        y        z           x        y        z
  ------------------------------------------------------------------
  1     Sr1    0.3644   0.0000   0.2500     1.0962  -4.1497  -2.7991
  ...

radwavfunc file.sipf:
program to evaluate the radial wave function given by th parametrisation in file.sipf.
reduce_unitcell [options] mcphas.j:
This program checks every atom in the unit cell in file mcphas.j and removes any atom, which is connected to another by a lattice vector. Useful for going from a description in an extended unit cell to the primitive unit cell.
Option: -nofcomponents 23
fixes the nofcomponents to 23 by reducing (removing entries) or increasing (by filling with zeroes) the exchange parameter tables
Option: -ni
forces output without indexchange
Option: -delatoms 1,2,5,7
instead of removing atoms connected by a lattice vector remove atoms number 1,2,5 and 7 from the list and also all interactions with those Note the following approximation is of limited use and NOT recommended, better use -delatoms phon (see below) instead of -delatoms 1,2,...: In case the atoms to be removed have the phonon module, effective interactions are introduced between the remaining atoms, in this way the crystal field phonon interactions with 4f shells can be treated in an effective way: Minimizing the energy by setting zero the derivative of the Hamiltonian (134) with respect to the phonon coordinate $\hat u_i^{\alpha}$ and in the spirit of the Einstein model neglecting bilinear terms (i.e. $K_{\alpha\beta}(ij)$ with $i \neq j$34) in $K_{\alpha\beta}(ii)$, the phonon displacement $u_i^{\alpha}$ can be estimated in a mean field approach by
  $\displaystyle
K_{\alpha\beta}(ii)u_i^{\beta}=\sum_{j\neq i,\gamma} -\Gamma^{...
...mma}(\hat \mathbf J_j)\rangle - G_{mix}^{\gamma\alpha}(i) \epsilon_{\gamma}
$ (164)
From $\overline{K}=\overline{S}^T\overline{\Omega}\overline{S}$ (compare equations (91) ff.) we compute the inverse $\overline{K}^{-1}=\overline{S}\overline{\Omega}^{-1}\overline{S}^T$, or in index notation $K^{-1}_{\alpha\beta}=\sum_{\alpha'=1,2,3}S_{\alpha\alpha'}\Omega^{-1}_{\alpha'\alpha'}S_{\beta\alpha'}$, and compute $u_i^{\alpha}$ from equation (175)

  $\displaystyle
u_i^{\alpha}=\sum_{\alpha'\beta=1,2,3}S_{\alpha\alpha'}\Omega^...
...hat \mathbf J_j)\rangle - G_{mix}^{\gamma\beta}(i) \epsilon_{\gamma} \right]
$ (165)

This expression (176) for $u_i^{\alpha}$ we substitute into the total Hamiltonian (134). Thus we get an effective multipole interaction contribution and an additional magnetoelastic contribution instead of the Einstein term $K_{\alpha\beta}(ii)$ and the crystal field phonon interaction $\Gamma^{\alpha\gamma}(ij) $ and the mixing term $G_{mix}^{\gamma\alpha}(i)$:

$\displaystyle \hat H_{\rm tot}$ $\textstyle =$ $\displaystyle \dots + \frac{1}{2}\sum_{\alpha'\beta\beta',i \neq jj',\gamma\gam...
...O_{\gamma}(\hat \mathbf J_j) \langle O_{\gamma'}(\hat \mathbf J_{j'}) \rangle +$ 
$\displaystyle \dots$ $\textstyle +$ $\displaystyle \sum_{\stackrel{i\neq j,\alpha'\beta\beta'}{\alpha=1-6,\gamma=1,\...
..._{\alpha\gamma=1-6}c_{\rm eff}^{\alpha\gamma}\epsilon_{\alpha}\epsilon_{\gamma}$(166)
Moreover, the terms quadratic in the strain $\epsilon$ lead to a renormalisation of the elastic constants, which become softer because of the relaxation of atomic positions in the unit cell when the crystal is strained (an effect which is not included in the elastic constants in equations (109))

  $\displaystyle
c_{\rm eff}^{\alpha\gamma}=c^{\alpha\gamma}+\sum_{\alpha'\beta\...
...ha'\beta'}S_{\beta\beta'}G_{mix}^{\gamma\beta}(i)}{\Omega_{\beta'\beta'}(ii)}
$ (167)

This effective multipole interaction and magnetoelastic terms are added in the output in case an atom which is removed has the phonon module. Special attention deserves in this context the first term with $j=j'$ which in general for $\gamma=\gamma'$ has to be positive and corresponds to a phonon induced quadrupolar self interaction $\frac{1}{2}\sum_{\alpha'\beta\beta',i \neq j}\frac{S_{\beta\alpha'}S_{\beta'\al...
...alpha'\alpha'}(ii)} O_{\gamma}(\hat \mathbf J_j) O_{\gamma}(\hat \mathbf J_{j})$. Note that this term actually is a sum of $\sum_{\alpha'\beta\beta',i \neq j}\frac{S_{\beta\alpha'}S_{\beta'\alpha'}\Gamma...
...(ii)} O_{\gamma}(\hat \mathbf J_j)\langle O_{\gamma}(\hat \mathbf J_{j})\rangle$. and $-\frac{1}{2}\sum_{\alpha'\beta\beta',i \neq j}\frac{S_{\beta\alpha'}S_{\beta'\a...
...{\gamma}(\hat \mathbf J_j) \rangle\langle O_{\gamma}(\hat \mathbf J_{j})\rangle$. Thus the mean field (37) is computed without the $1/2$ prefactor and the second term is considered in the energy correction term (36) after the end of the mean field loop. The self interactions is always ferroquadrupolar (note that by definition $\Omega_{\alpha'\alpha'}(ii)<0$) and thus enhances the crystal field. The self interaction describes the effect, that when a crystal field is produced by a point charge, then the charge density on the magnetc ion deforms according to the crystal field. This deformation increases the electrostatic force between the point charge and the magnetic ion and the point charge moves accordingly and this movement leads to a further increase of the crystal field. The movement continues until the restoring elastic force is equal to the increase in electrostatic force.

Option: -delatoms phon
instead of removing atoms connected by a lattice vector remove atoms which have a phonon module. This option is recommended to obtain a reduced model containing only magnetic ions. In contrast to the option -delatoms 1,2,5,7 it takes into account the correlated displacements which arise when a stress is applied or a magnetic atom charge density is changed by magnetic field or temperature. The program will perform the following operations
Renormalisation of Elastic Constants
The relaxation of the lattice leads to a reduction of elastic constants with respect to the values obtained by the Born van Karman model equation (109). Effective elastic constants (which correspond to experimental observable elastic constants) can be computed from (134) by applying various external stresses and computing the strain. reduce_unitcell considers in such a calculation only the phonon degrees of freedom (i.e. removes all magnetic ions) and sets nonzero each component of the stress tensor to 0.01 GPa and computes the resulting strains using the mean field loop of the mcphasit program. For the computation the parameters in reduce_unitcell.ini are used, supercells should be specified in order to relax displacements not only with the periodicity of the primitive unit cell as given in the input file. The elastic compliance matrix [64] is obtained by dividing the computed strain by the stress of 0.01 GPa. The elastic constants $c_{\rm eff}^{\alpha\gamma}$ are computed by inverting the matrix of the compliances.
Calculation of effective phonon induced Quadrupolar Self Interaction
A similar self consistent mean field procedure using (134) is adopted to calculate the phonon induced quadrupolar self interaction , initially without any stress/strain to determine the reference energy $U_0$:
  1. From the input parameters all interactions between the magnetic ions are removed
  2. Similar, all the magnetoelastic interactions are removed.
  3. all magnetic interaction operator expectation values are set to zero except one component $\gamma$ for the ion $n$, i.e. $\langle I^i_{\beta}\rangle=\delta_{\beta\gamma}\delta_{in}$.
  4. a mean field loop is performed at zero strain to determine the phonon displacements $\langle \mathbf u_{j} \rangle$ induced by the crystal field phonon interaction $\Gamma^{\alpha\gamma}(jn)$ in equation (133) and to determine the corresponding energy $U'<U_0$.
  5. The phonon induced self interaction is determined by $\mathcal J_{\gamma\gamma}(nn)=-2(U'-U_0)$ and inserted in the output of reduce_unitcell.
  6. this procedure is repeated for all components of all interaction operators of all magnetic ions.
Calculation of effective phonon induced Magnetoelastic Interactions
The next step is to do mean field loops using (134) the same way as for the quadrupolar self interaction (i.e. setting $\langle I^i_{\beta}\rangle=\delta_{\beta\gamma}\delta_{in}$) and setting the strain tensor constant to be $\epsilon_{\alpha}=\delta_{\alpha1} 10^{-6}$. phonon displacements $\langle \mathbf u_{j} \rangle$ and energy $U''$ are calculated self consistently. The phonon induced magnetoelastic interaction $G^{1\gamma}_{\rm phon-induced-cfph}(n)$ is computed using the energy difference $U''-U'$ and the effective elastic energy $E_{el}=\frac{1}{2}c_{\rm eff}^{11}\epsilon_{1}\epsilon_{1}$ by

  $\displaystyle G^{1\gamma}_{\rm phon-induced-cfph}(n)=-\frac{(U''-U')+E_{el}}{\epsilon_{1}}
$ (168)

This phonon induced magnetoelastic interaction $G^{1\gamma}_{\rm phon-induced-cfph}(n)$ is added to the $G^{1\gamma}_{\rm cfph}(n)$ given in the input file and the sum is output by reduce_unitcell. The procedure outlined above for $\epsilon_1$ is repeated for all components of the strain tensor $\epsilon_{\alpha}$ and thus the full set of phonon induced magnetoelastic interactions $G^{\alpha\gamma}_{\rm phon-induced-cfph}(n)$ is obtained.

Calculation of effective phonon induced Quadrupolar Interactions
In a further step the effective phonon induced Quadrupolar Interactions are calculated. In the mean field loop the strain is kept zero. In order to compute $\mathcal J_{\gamma\gamma'}(nn')$ the interaction operator expectation values are set to zero except for the component $\gamma$ for the ion $n$ and $\gamma'$ for the ion $n'$, i.e. $\langle I^i_{\beta}\rangle=\delta_{\beta\gamma}\delta_{in}+\delta_{\beta\gamma'}\delta_{in'}$. Note the self interaction term ( $\mathcal J_{\gamma\gamma}(nn)$ has been computed already). The result of the mean field computation is the energy $U'''$ and the effective phonon induced quadrupolar interaction constants are computed by $\mathcal J_{\gamma\gamma'}(nn')=-(U'''-U_0)-\mathcal J_{\gamma\gamma}(nn)/2-\mathcal J_{\gamma'\gamma'}(n'n')/2$. These interactions are added by reduce_unitcell to the the two ion interactions of the input file.

Option: -delatoms 1:3+-0.75:4+-0.25,2,5,7
removes atoms number 1,2,5,7 for atom 1 the interactions are transferred to 75% to atom 3 and 25% to atom 4 and interactions to atoms 3 and 4 are removed, for 2 5 7 all the interactions with other atoms are removed
Option: -mcdiff
create also mcdiff.in file with reduced unit cell
Option: -v
verbose mode
rotateBlm
    Rotates a set of  crystal field  parameters for  Stevens equivalent
    operators by an azimuthal angle fi about the original z axis and
    a polar angle theta about the new y axis. A right hand axis system is assumed
    and a positive rotation is one which advances a right-hand screw in a
    positive direction along the axis.

    The calculations are  done by means of matrix  multiplication based on
    the method of Buckmaster (phys. stat. sol. a, vol 13,  pp 9, 1972) and
    Rudowicz (J. Phys: Solid State Phys., vol 18, pp 1415, 1985).   

    usage: rotateBlm [-h] [--help] 
              [-i input_file] [--input input_file]
	      [-o output_file] [--output output_file]
              [-v] [--verbose] [-th theta] [-fi fi] [CF parameters]

          
     -h          : this (help) message
     -i in_file  : input CF parameters file in cfield or mcphase formats
     -o out_file : output CF parameters file in mcphase format
     -v          : verbose mode. Will print out input parameters as read.
     -th	 : polar angle theta in degrees
     -fi         : azimuthal angle fi in degrees

    if -i is omitted, the program will  assume the input CF parameters are
          given on the command line in the format: Bkq=x.xx,Bkq=x.xx, etc.
          e.g. $0 B20=0.21,B40=0.0005,B60=0.051,B66=0.626
    negative q parameters such as B_2^{-2},  are specified as:  B22S, with 
          an 'S' at the end, as per the McPhase convention.
    you may also  specify the ion type by a dding another  parameter after
          the CF parameters: e.g. $0 B20=0.21,B40=0.5 Pr3+
    if -o is omitted, the program prints the parameters to standard output.
setup_jqfit[-h] [–help] h k l:
program to setup a fit of exchange parameters in order to reproduce an experimental propagation vector.
     -h          : print help message
      hkl        : Miller indices of propagation vector

    required input files:

    mcphas.j (+ single ion paramter files)
                 :  structural information including all magnetic atoms

    output files:

    mcdisp.par   :  contains propagation vector and list of other hkl to
                    be probed
    mcdisp.mf    :  required input file for mcdisp
    calcsta      :  required input file for simannfit and searchspace
    calcsta.pl.forfit: file with fitparameters for Bethe slater, RKKY fits
    fit.bat      :  batch to start the fit

After running this program you can start immediately a fit of exchange parameters. Edit calcsta.pl.forfit and fit.bat to fine tune the fit according to your needs. During fitting a value of sta $< 1$  indicates, that the maximum of $J(Q)$ is at the propagation vector tau. How much it is below one depends on the magnitude of $J(Q)$ for the competing wavevectors in the list inmcdisp.par.

setup_mcdiff_in [option] T Ha Hb Hc:
setup_mcdiff_in [option] x y:
setup_mcdiff_in [option] out1 out2 out3 ou4 out5 out6 out7:
program to setup mcdiff.in with information on spinconfiguration to be used by program mcdiff. Note, you must have done a mcphas calculation to stabilise a magnetic structure at the desired Temperature/Field. setup_mcdiff_in reads the results of this calculation from results/mcphas.mf and generates an input file mcdiff.in

     -h          : help  message
      T          : Temperature (K)
      Ha,Hb,Hc   : Magnetic Field (T)
      x,y         : x,y values in phasediagram
      out1-out7  : values of external parameters in columns 1-7 of output files mcphas.*

    required input files:

    results/mcphas.sps
                 :  result of a mcphas calculation

    optional input files:

    mcdiff.in
                 :  if present experimental parameters (section 1)
                    and nonmagnetic atoms (section 2) are taken from
                    this file

    output files:

    mcdiff.in    :  required input file for mcdiff

    - after running this program you can start mcdiff to do the calculation
      magnetic diffraction pattern

    options:

    -prefix 001  :  instead of results/mcphas.mf read results/001mcphas.mf
setup_mcdisp_mf [options] T Ha Hb Hc:
setup_mcdisp_mf [options] x y:
setup_mcdisp_mf [options] out1 out2 out3 ou4 out5 out6 out7:

program to setupmcdisp.mf with information on meanfields to be used by program mcdisp. Note, you must have done a mcphas calculation to stabilise a magnetic structure at the desired Temperature/Field. setup_mcdisp_mf reads the results of this calculation from results/mcphas.mf and puts the meanfields into mcdisp.mf.

     -h          : this (help) message
      T          : Temperature (K)
      Ha,Hb,Hc   : Magnetic Field (T)
      x,y        : x,y values in phasediagram
      out1-out7  : values of external parameters in columns 1-7 of output files mcphas.*

    required input files:

    results/mcphas.mf
                 :  result of a mcphas calculation

    output files:

    mcdisp.mf    :  required input file for mcdisp

    - after running this program you can start mcdisp to do the calculation
      of dispersion of excitations or diffuse scattering

    options:

    -prefix 001  :  instead of results/mcphas.mf read results/001mcphas.mf

setup_mcphasjforfit [-h]:
program to setup a fit of exchange parameters by creating mcphas.j.forfit from mcphas.j

     -h          : print help message

    required input files:

    mcphas.j (+ single ion parameter files)
                 :  structural information including all magnetic atoms

    output files:

    mcphas.j.forfit  : all interaction parameters are substituted
                       with parJxxx[0.0,-1e0,1e0,0,1e-6]

    - after running this program you must setup a file calcsta 
      to calculate the standard deviation and then you can start
      a fit by simannfit or searchspace

singleion [options] T[K] Hexta[T] Hextb[T] Hextc[T] Hxc1 Hxc2 Hxc3 ... Hxcnofcomponents
single ion - calculate single ion expectations values $<Ia>, <Ib> $... and transition energies. Options may be used to trigger calculation of magnetic Moment and other Quantities.
 
           T    ..... temperature in Kelvin
           Hext ..... external field in Tesla 
           Hxc... exchange (molecular) field in meV
singleion reads mcphas.j and the singleion parameter files quoted therein and calculates energies, eigenstates, expectation values $<I>$ for the given temperature, external magnetic field Hext and exchange field Hxc (the interaction constants given in mcphas.j are ignored).

for each single ion property file the following files are generated:

   results/file.sipf.levels.cef .. energy levels and eigenstates and <I>
   results/file.sipf.trs ......... transition energies,matrix elements
                                   and (powder) neutron intensities
   results/_file.sipf    ......... input parameter files as read and used by singleion 

options: -nt ......... by default only 5 transition energies are output,
                       if you want more, start e.g. with 
                       option -nt 7 to output 7 transition energies
         -pinit 0.1 .. consider only transitions with population of initial state > 0.1
         -ninit 3  ... consider only transitions from the 3 lowest eigenstates
         -maxE 30  ... consider only transitions with energy lower than 30 meV
         -E        ... output to stdout energy of cf levels instead of transition energy
         -E0       ... output to stdout energy of cf levels relative to ground state ...
         -E0f 0.3 0 2.1 3 4  ... output to stdout energy of cf levels relative to ground state 
                       and calculate sta (sum of squared deviations) of nt (by default 5,
                       see option -nt, which if used must be given prior to -E0f) 
                       transition energies (meV) to the numbers following -E0f 
                       a '0' will exclude this transition from contributing to sta
         -r ion.sipf . do not read mcphas-j but only the single ion
                       parameter file ion.sipf
         -U  ......... calculate energy U, ln of partition sum Z, free energy F
                       instead of <I>
         -t  ......... reads .trs files from previous run (possibly modified by user, 
                       i.e. removing
                       some lines to speed up an approximate calculation of spectra.)
         -Esteps 10 27 for option -d in addition to initial Energy calculate 10 further
                       Energies until 27 meV has been reached
         -Tsteps 10 27 in addition to initial temperature calculate 10 further temperatures
                       until 27K has been reached
         -Hsteps 20 0 0 10 in addition to initial field calculate 20 further external fields
                       until (0 0 10) Tesla has been reached
         -HE ......... in addition to magnetic field also apply electrical field in kV/mm, 
                       i.e. there will be 
                       instead of 3 components  Hexti Hextj Hextk in the command 
                       line 6 components 
                       Hexti Hextj Hextk Eexti Eextj Eextk, similar for Hsteps option 
                       there will be 6 components
                       (therefore: mind that -HE is given before -Hsteps in the command line)
         -opmat 2 .... output operator matrix number n=2 to results/*.opmat
                Operators in results/output op.mat for different values of n:
                n=0                    Hamiltonian
                n=1,...,nofcomponents  operator Matrix In in standard basis
                n=-1,..,-nofomponents  operator Matrix In for Hamiltonian eigenstates 
                                       basis
                n>nofomponents: all operator Matrices (n=0 to n=nofcomponents) in 
                                standard basis
                n<-nofomponents: all operator Matrices (n=0 to n=-nofcomponents) 
                                 in Hamiltonian  eigenstates basis
         -v       .... verbose, output more information on ongoing calculation

  Other Observables:
         -M  ...magnetic moment: calculate expectation values and transition matrix
                       elements for magnetic moment M (muB)instead of I
         -P  ...phonon displacement: calculate phonon displacement in A instead of I
         -pel ..electric dipole moment: calculate electrical dipole moment pel in |e|pm
                       instead of I
         -L  ....orbital momentum: calculate expectation values and transition matrix
                       elements for orbital momentum L
         -S  ....spin: calculate expectation values and transition matrix
                       elements for spin S
         -MQ 0 0 1 ...M(Q),Fourier transform of magnetic moment density:  instead of <I>
                       calculate expectation values, transition matrix elements
                       for M(Q=(0 0 1)/A), the Fourier Transform  of magnetic moment 
                       density M(r) 
         -sx ....spin density: calculate expectation values and transition matrix
                       elements for spindensity coefficients aSx(lm) in expansion 
                       of spindensity-x-component in Ms(r) = sum_lm aS(l,m) R^2(r) Zlm(Omega)
                        E. Balcar J. Phys. C. 8 (1975) 1581
         -sy -sz ......for y and z components use option -sy and  -sz
         -lx ....orbital moment density: calculate expectation values and transition matrix
                       elements for orbital moment density coefficients aLx(lm) in expansion 
                       of orbital moment density-x-component 
                       in Ml(r)=sum_lm  aLx(l,m) F(r) Zlm(Omega)
                       with F(r)==1/r int_r^inf R^2(x) dx,   E. Balcar J. Phys. C. 8 (1975) 1581
         -ly -lz ......for y and z components use option -ly and  -lz

   Susceptibilities:
         -X[observable] 24 0.1  ... calculate susceptibility chi(z) for z=E + i epsilon with
                       E=24 meV and epsilon=0.1 meV. The susceptibility  chi is defined as the
                       derivative of I (or any observable) with respect to the field, e.g. for
                       the observable magnetic moment -XM the susceptibility is chi=dM/dH
The program evaluates equation (267) setting $\hat \mathcal I_{\alpha}$ equal to the magnetic moment operator $\hat M_{\alpha}$ with $\alpha=1,2,3$ for the three cartesian directions x,y and z. Having evaluated the magnetic single ion susceptibility it is possible to use the fluctuation dissipation theorem to calculate the magnetic inelastic single ion neutron scattering cross section in dipole approximation for a polycrystal neglecting form factor, Debye Waller factor by the following relation

  $\displaystyle
S_{\rm mag,poly}^{inel}(Q=0,z)= \frac{2}{3}\frac{(\gamma r_0)^2}{4\pi} \frac{1}{1-exp(- z/k_bT)} \sum_{\alpha=1,2,3}\chi''_{\alpha\alpha}(z)
$ (169)

Note, the prefactor $2/3$ comes from the average of the polarisation factor. Compare (279), (305) and (306).

                       units: standard unit is (unit of observable)^2/meV, e.g. for -XM
                              (muB)^2/meV. you can use the options ... 
         -emu 0.3 0.1 . to obtain the magnetic susceptibility chi in units of emu/mol with
                       constant offset X0=0.3 emu/mol and molecular field constant 
                       lambda=0.1 mol/emu
                       1/(X-X0)=(1/Xcf)-lambda. Xcf obtained the same way as option 
                       -d 0 0 and
                       converting the results to emu/mol by multiplying with factor
                       MU_B MU_B NA/10000=
                       = 0.0578838263*0.55848973464 = 0.0323275227902. 
                       lambda is only applied if Xcf is diagonal.
                       e.g. -XM 0 0 -emu 0 0 will calculate the static magnetic 
                       susceptibility in emu/mol
         -SI ......... to obtain the magnetic susceptibility (-XM) in SI units (Am^2/mol)
         -muBT ......... to obtain the magnetic susceptibility (-XM) in  (muB/Tesla)
         -iX[observable] 24 0.1  ....same as -X, but output inverse susceptibility
                       (works only if off diagonal elements of X are zero)


 Note: for calculating T,H dependencies you can instead of using options -Tsteps 
      or -Hsteps put singleion in a LOOP
      and pipe the result into a file
 ... LOOP linux:   for B in $(seq 0 0.1 14); do singleion 2 $B 0 0 0 0 0; done > results/fielddep.dat
 ... LOOP linux using perl:
perl -e 'for($B=1;$B<14;$B+=0.1){system("singleion 2 ".$B." 0 0  0 0 0");}' > results/sus1Tesla.clc 
 ... LOOP for windows using perl:
perl -e "for($B=1;$B<14;$B+=0.1){system('singleion 2 '.$B.' 0 0  0 0 0');}" > results\sus1Tesla.clc

  .... windows command line: for /L %B in (0,1,14)  do singleion 2 %B 0 0 0 0 0 >> results\fielddep.dat

  .... windows batch file (needed for noninteger numbers):
          @echo off && setlocal ENABLEDELAYEDEXPANSION
          for /L %%I in (0,2,140) do ( set /A W=%%I/10 && set /A "f = %%I %% 10"
          set B=!w!.!f!
          @echo on && singleion 2 0 0 !B! 0 0 0 && @echo off )
          endlocal && @echo on

symhmltn:
Calculates via group theory, the allowed two-ion multipolar exchange terms between two ions whose center has some particular point symmetry. Syntax is: "symhmltn $\langle$PTGP$\rangle$ $\langle$l$\rangle$" where $\langle$PTGP$\rangle$ is the name of the point group (both Hermman-Maguin and Schoenflies notation are understood), and $\langle$l$\rangle$ is the integer multipolar order (l=1 for dipolar exchange, l=2 for quadrupolar, etc.). E.g. "symhmltn Oh 1" gives the familiar Heisenberg exchange terms. See also: cfsplit which calculates the multiplicities of the crystal field split levels of a particular point group.