Simulation of EPR-spectra of randomly oriented samples

Computer Physics Communications 21 (1981) 385 —395 © North-Holland Publishing Company

SIMULATION OF EPR-SPECTRA OF RANDOMLY ORIENTED SAMPLES C. DAUL, C.W. SCHLAPFER Institute of Inorganic Chemistry, University of Fribourg, 1700 Fribourg, Switzerland

B. MOHOS Brown Boveri Company, Baden, Switzerland

J. AMMETER and E. GAMP Institute ofInorganic Chemistry, University of Zurich, 8052 Zurich, Switzerland

Received 27 November 1979; in revised form 29 September 1980

PROGRAM SUMMARY Title of program POWDER

dipole—dipole and/or quadxupole interaction of the electronic and nuclear spins. This program calculates the first derivative of the EPR absorption spectrum of randomly oriented sam-

Catalogue number: ABYG

pies using the following approximations: i) the eigenvalues of the spin Hamiltonian axe given by second order perturbation theory; ii) the intensities of the EPR transitions axe determined by Zeeman interaction only; iii) the paramagnetic species are uniformly or randomly distributed in space. iv) “allowed” transitions (AM 1 = 0) axe calculated only.

Program obtainable from: CPC Program Library, Queen’s

University of Belfast, N. Ireland (see application form in this issue) Computer: CDC 6000; Installation: Rechenzentrum der

ETHZ, 8052 Zurich, Switzerland Operating system: SCOPE 3.4 Programming language used: FORTRAN IV G

Method of calculation

The single crystal spectra for particular orientations are calculateri. They axe summed over all spacial orientations (Simpson rule) and cpnvoluted with a line shape function giving the absorption line. A subsequent numerical derivation yields the 1St derivative spectrum and reduces the random deviations due to the limited number of orientations [1].

High speed storage required: 130 0008 words No. of bits in a word: 60 Overlay structure: none No. of magnetic tapes required: none

Restriction on the complexity of the problem Other peripherals used: card reader, line printer, plotter,

No. of cards in combined program and test deck: 826

In the present version, the stick spectrum is convoluted with a lineshape function after the summation over all orientations has been carried out. This implies that only line widths independent upon orientation and mj can be treated. Further-

Card punching code: BCD

more, since a perturbation calculation is used, it ~c(dipole— is required that ~C(Zeeman) > ~c(hyperfme),~lC(Zeeman)>

disk space

Keywords: EPR, anisotropicg-tensor, anisotropic hyperfme-

dipole) and that ~c(hyperfine)> ~C(quadrupole).

tensor, anisotropic perturbation calculation, numerical integration and differentiation, noise filter.

Ranges between 10 s and 10 mm on the CDC 6000.

Nature of the physical problem

References

The EPR spectra of polycrystalline paramagnetic samples exhibits often complex features due to hyperfine and/or

[1J B. Mohos et al., Anal. Chem. 48 (1976) 231; J. Chem. Phys. 60 (1974) 4633.

Typical running time

385

C. Daul et al. / Simulation of EPR.spectra

386

LONG WRITE-UP 1. Introduction The EPR spectrum of a randomly oriented sample results upon superposition of a large number of randomly oriented single crystals. Thus the simplest way to obtain a powder EPR is to take a single crystal spectrum and to integrate it on the surface of the unit sphere. Such an integration can be done analytically in a restricted number of simple cases [2]. For general purposes this integration must be done numerically and we tried three different types of numerical integrations i.e. Simpson rule, Gaussian quadrature and Monte-Carlo technique. The optimal results were obtained when using the Simpson rule followed by numerical differentiation [1] (EPR spectra are usually recorded as first derivatives). This subsequent differentiation acts as a noise filter and allows to cut down the integration time by a factor of ten.

2. Method of calculation 2.1. Energy levels ofthe spin Hamiltonian We assume that the spin Hamiltonian can be written as ~Cf3(~gH)+~D~+SAI+IQi,

(1)

where fi is the Bohr magneton, ~ the operator of the electronic spin, H the applied magnetic field, g the g-tensor, D a second rank tensor describing the dipole—dipole interaction, A the hyperfine tensor and Q the quadrupole tensor. Al! these tensors are of rank 2 and symmetric. Moreover D and a are traceless. In all calculations we chose the principal axis system of the g-tensor as laboratory systçm. The advantage of this assumption is that the g-tensor is diagonal in this axis system. The application of second order perturbation theory to (1) yields the following energy levels [3] EWS,MI,S,I)=EZe+ED_D+EHF+EQ,

(2)

where S and I are, respectively, the electronic and the nuclear spin andM

5 and Mj the corresponding projection quantum numbers. E~e,EDD, EHF and EQ are respectively the Zeeman, the dipole—dipole, hyperfine and quadrupole contributions, with: E~egI3HMs, E~D= ~[3M~ S(S + 1)] P/f 2g(3H ~ Nfl [4S(S + 1)— 8M~ 1] Ms 2) 2N~+N~’ 2N~det(D)] [2S(S + 1)— 2M~ 1] M ~[Tr(D 5}, 21 M 2) L~][1(1+ 1) MI] Ms EHF = KM5M1 + 2g~1ff[1L2 K 5 + [Tr(A —

—

—

—

—

—

—

—

—

—

—

det(A)[S(S + 1)

—

—

~

M~1M 1)~

EQ[3MI_I(I+1)1P2l_~1j~ {[F~—P~][4I(I+l)—8MI—11M1 2) 2F~+P~ 2F~det(Q)] [21(1+1)— 21WJ 1] M ~[Tr(Q 1, —

—

—

—

C. Daul et al. I Simulation of EPR-spectra

387

where the following abbreviations have been used: 2 eg2e, g g2A1~

eglY’ge,

g2K2

egA2ge,

g2K2L2

egA4ge,

p—1,+land+2,

g2K2P~~ egAQvAge,

v’—l,+land+2,

where e = (sin 0 cos 0, sin 0 sin çb, cos 0) a unit vector describing the direction of the magnetic field H = He. 2.2. Location ofthe absorbtion lines The location of the resonance field of the absorption lines in the presence of an electromagnetic field with frequency w is ~

(3)

Combining (2) and (3) we get H~1fs_~Ms_1,MJ,S,I)=_~JVf[2Ms_1]+ 2~2H{[N22_Nfl 2) X [4S(S ÷1)— 24M5(M5 1)— 9] ~[Tr(D X[2S(S+1)—6MsWs—1)—3]}—.~_KMj —

H

—

2N~+N~ 2N~det(D)] —

—

{[L2 -K2]M~+~[Tr(A2)-L2J[J(J+l)-M7]

+~det(A)(2Mg- 1)M 1}

I

—~

{[P~-pfl[4J(Ii-l)_8MI_1]M1

2g13KM5~lfg 1) 2) 2F~+4 ~-[Tr(Q —

—

—

where H0

=

—

2P~det(Q)] [21(1+ 1)— 2MJ

—

1] M 1},

(4)

hw/gj3 is the centre of symmetry of the spectrum.

2.3. Transition probabilities Let us take an electromagnetic field with2Hf an amplitude Hi I << I(MgLS”ge’1M 2. !HoI in the e’ direction; the transition probability IM~)-~ Ii1f~’) is of then proportional to j3 5’>I Choosing the axis quantization (Sz) parallel to e, the probability of transition due to the H

1 -field is then given by P1(1M5)

-~

1M5

—

1))

{e’g’e’ + 4

-~

(ege)}(M5l~’+IM5 1;). —

g

Integrating out over all possible directions ofH1 we get the probability of the EPR transtition:

P(M5-+M5-_1)=~_.

[Tr(g2)_~(eg4e)][s(s+1)_M.~+Ms1

.

(5)

388

C. Daul et al.

/ Simulation of EPR-spectra

2.4. Spatial integration We choose one integration variable to be u = cos 0 instead of taking simply 0 and the second one as 0. This is so in order to subdivide the surface of the unit sphere into increments (~u)(i~0) of approximately equal size. Thus the integration over the whole (0, 0)-space is given according to Simpson’s rule by 2ir

+1

5 5 —1

2N

0+1 2N +1

f(u, 0) du dO

0

w~wjf(u1,01), ~

r=i

j1

where hu = 1/N is the integration step length in u = cos 0, h = n/NO is the integration step length in 0, (2N0 + 1) and (2N~,+ 1) are the number of integration points in u and 0~respectively, (1 ifiorf

1 or2N+ 1

worw1= ~2ifiorjeven, ~4ifiorf odd. Let us digitalize the magnetic field axis in the following way Hk=HL+Sk,

~

where HL is the lower bound of the magnetic field. If we denote by Hu the upper bound of the magnetic field then the increment ~ is given by ~ = (Hu — HL)/(NptS — 1) i.e. the resolution of the simulated spectrum. Using this notation the EPR-absorption line y (Hk) is obtained by

Y(Hk) =

~f

i,J,MS,MIwiw?(~~i~ 0~,Mg),

H(u1, 0~,Ms,Mj) E [Hk

—

,

H,~

+-~-]

(6)

where the summation in (6) runs by steps of one over all i = 1,2, ..., (2N8 + 1),/ = 1,2, ..., (2Nd, + 1), Mg = —S + 1, —s + 2, ..., Sand M1 = —I, —I + 1, ..., +1 provided that the resonance field H(u1, 0~,Mg, M1) lies in the interval of Hk i.e. if

Hk —~-~H(u1, 0j,Ms,MI)
+-~--,

P(u1, 0~,M~)denotes the intensity of a single crystal absorption line M5 ~ Ms and 01,H(u1, 0j,Ms,Mi) is the corresponding field location.

—

1 with the orientation cos 0

u1

2.5. Convolution with a line shape function The absorption line in eq. (6) corresponds to transitions of rectangular shape with a very narrow width equal to & In order to conform this absorption line to be experimentally observed situation, a convolution of y(Hk) in eq. (6) with an adequate line shape function F(x) is necessary. The program can handle either Gaussian or Lorentzian lineshapes. In practice however, Gaussian lineshapes have proved to be the most adequate. The new absorption line z(Ilk) then becomes

Z(Hk)-~

y(Hk+XJ)~X//w),

(7)

where W is the line width and is incremented in steps of ~. In practice, however, since F(x) is symmetrical and declines very rapidly to zero eq. (7) can be rewritten as

z(Hk)=~Wj[y(Hk+XJ)+y(11k—Xj)1F(X//W),

(8)

C. Daul et al.

where X

1

j~,w

=

~ ~

I Simulation of EPR-spectra

and M— IFIX(Xw/~)+ 1, X

389

2.5 for Gaussian line shape and X = 10 for Lorentzian

lineshapes. 2.6. Numerical calculation of the first derivative The main purpose of this calculation is to eliminate stochastical “noise” due to the finite number of orientations which heve been used. This supplementary numerical work is therefore preferred to the direct method which consists in carrying out the convolution with the first derivative of the lineshape function. The final form of the simulated spectrum T(Hk) is then given by the following transformation:

T(Hk)~wJ[z(Hk+j6)-z(HZ-j~)], where w1 = { ~1’~ ~

(9)

and ~uis chosen equal to IFIX(W/26) + 1.

3. Program description Using eqs. (3) and (4), subroutine SINGLE calculates the resonance field location H(u,PHI,MS,MI) and the transition probab~jtiesTP(u,PHI,MS) for a fixed orientation u = cos 0 and 0 and for all allowed EPR transitions with —S + 1
Subroutines RECHEN and STRICH are called by PLTSP. They calculate the draw particular resonance fields. Other routines like DET calculate the determinant of a 3 X 3 matrix; MSQ calculates the square of a 3 X 3 matrix; QF calculates the quadratic form x- A x; GJRD is a matrix inversion routine. GIRL is used to print a 3 X 3 matrix.

4. Description of common blocks

[SING] G(3) A2(3,3) A4(3,3) DM1(3,3) D(3,3) D2(3,3) QM1(3,3)

: diagonal element of g-tensor : A-tensor squared : A-tensor rised to the 4th power : D-tensor inverted : D-tensor : D-tensor squared : AQ’A

Q(3,3)

: Q-tensor

Q1(3,3) Q2(3,3) TA2 TD2

: AQA 2A :AQ : Tr(A2) : Tr(D2)

390

TQ2 DA DD DQ OMEGA GG TG2

C. Daul et al. /Simulation ofEPR-spectra

: Tr(Q2) : det(A) : det(D) : det(Q) : w the microwave frequency in Hz : g-tensor squared : Tr(G2)

[SPECS] YS(7000) : stick spectrum [SPEC] Y(7000) [HPARI HL • HU NPTS FACT W [RECH] AH(3) RI IFRE NW NNW NSPEC1

: convoluted spectrum

: lower bound of magnetic field axis

: : : :

upper bound of magnetic field axis number of points in the simulated spectrum A, usually equal to 2.5 for Gaussian lineshape, 10 for Lorentzian line width

: diagonal elements of hyperfine tensor : nuclear spin : logical variable for the call of RECHEN : number of different line widths :linewidthindex,NNW1,2,...,NW : number of segments of the calculated spectrum to be plotted.

5. Input description 1st card: Spin definition card NI, NS (215) NI 21 + 1 the nuclear spin multiplicity NS 2S + 1 the electron spin multiplicity 2nd card: Title card (TEXT(I),I=1 ,3) (3A10)

-

3rd card: g-tensor and parameter card G(1), G(2), G(3), OMEGA, HL, HU (6F10.0) G(1), G(2) and G(3) eigenvalues ofg-tensor OMEGA microwave frequency in Hz HL lower field bound in gauss HU

upper field bound in gauss

C. Daul et aL

I Simulation of EFR-spectra

391

4th card: Line width definition card WMIN, WMAX, DW (2F1 0.0) WMIN lower bound of line width WMAX upper bound of line width DW line width increment 5th card: Hyperfine tensor card A(1,1), A(2,2), A(3,3), A(1,2), A(1,3), A(2,3) (6F10.0) (A(I,J) I ~ J are the elements of the real and symmetrical A -tensor expressed in the axis system of the diagonal g-tensor 6th card: Quadrupole interaction tensor card Q(1,1), Q(2,2), Q(3,3), Q(1,2), Q(1,3), Q(2,3) (6F1 0.0) Q(I,J), I ~ J are the elements of the real, symmetrical and traceless Q-tensor expressed in the axis system of the diagonal g-tensor 7th card: Dipole—dipole interaction tensor card D(1,1), D(2,2), D(3,3), D(1,2), D(1,3), D(2,3) (6F10.0) D(I,J), I ~ J are the elements of the real, symmetrical and traceless D-tensor expressed in the axis system of the diagonal g-teflsor 8th card: Numerical parameters card NSPEC, NPRS1, MU, NTHETA, NPHI, CTHETA, FPHI, PPL, FACT (4F10.0) 215,3110 NSPEC NPTS1 MU NTHETA

: : : :

NPHI

CTHETA

: :

FPHI

:

PPL FACT

: :

number of segments of the simulated spectrum to be drawn number of points per segment (default value: 400) number of points used in the first derivative calculation (default value in MU = W/25) number of points used in the spatial integration over cos 0 (default value is NTHETA = SQRT(5OXNPTS)) number of points used in the spatial integration over phi (default value is NPHI = SQRT(5OXNPTS)) upper bound of cos 0 in spatial integration. This option permits one to restrict the summation to the first octant in the case of axial symmetry for example. upper bound of 0 in units of ir (default value is 0.5) number of plotted points per line width (default value is 5) A in note to eq. (8) and depends upon the line shape function used (default value is 2.5).

9th, 10th, cards: Plot parameters cards (HL1(J), HU1 (J), SINT(J), SABS(J), J=1, NSPEC) (4F10.0) HLI(J) Elul (J) ~ field bounds of segment no. i ...

SINT(J) maximal intensity

SABS(J) distance in cm of the spectrum from the bottom of the page

392

C. Daul eta!.

I Simulation ofEPR -spectra

last card: Output option card PRINT, IFRE, IPLOT, LS (Li ,4X,L1 ,4X,Ll ,4X,11) PRiNT logical option to store the simulated spectrum on file IFRE logical option to call RECHEN IPLOT logical option to call PLTSP LS lineshape LS1 Gauss LS2 Lorentz 6. Output description

The print out of POWDER consists of a control copy of the parameters read in followed by a plot of the simulated spectrum. Test run Figs. 1—4 together with the printed test run output show the output from the test run included with the program.

TEST

FREQUENCY LINEI~lIDTH 2~C.D0

293.12

3~UOOOO 6HZ 10~00 GAUSS 546.25

290.37

252.50

I 255.62

268.75

261.87

265.55

FIELD [GAUSS]

286.73

277.25

27q~37

271.55

2811.82

283.15

288.87

~1O’

Fig. 1. Plot of the simulated spectrum; frequency: 9 GHz, linewidth: 10 G, field bound: from 2400 to 2900 G.

288.00

C. Daul et aL

/ Simulation of EPR-spectra

393

TEST

LI FREQUENCY : LINEWIOTH : 235.28

263.12

I

100000 6HZ 10000 GAUSS

2~6.2S

28.37

72.50

245.62

36.75

317.81

310.00

218.72

327.25

32’4.37

327.50

335.82

333.70

336.81

382.00

FIELD [GAUSS) ~iü~ Fig. 2. Same as fig. 1, except field bound: from 2900 to 3400 G.

TEST

L FREQUENCY : LINEkJIIJTH 250.00

0’~3.l2

100000 6HZ 20~00 GAUSS

286.25

282.37

053.50

I 255.62

2SL75

261.87

260.00

268.72

247.55

248.37

FIELD [GAUSS] ~ Fig. 3. Same as fig. 1, except linewidth: 20 G.

277.80

210.62

263.70

286.87

256.85

394

C. Daul et aL

/ Simulation of EPR-spectra

TEST

FREQUENCY LINEWIOTH 250.00

203.12

:

100000 6HZ 2100 GAUSS

286.2S

230.31

302.50

305.62

358.75

3)1.87

315.05

315.12

327.25

328.37

FIELD [GAUSS] 88i0~ Fig. 4. Same as fig. 2., except linewidth: 20 G.

Acknowledgement We thank the Swiss National Science Foundation for financial support.

References [1] B. Mohos et aS., Anal. Chem. 48 (1976) 231; J. Chem. Phys. 60 (1974) 4633. [2]F.K. Kneubühl, 1. Chem. Phys. 33 (1960) 1074.

327.50

330.62

333.75

336.61

350.00

C. Daul et al.

I Simulation of EPR-spectra

395

TEST RUN OUTPUT

2.5 * 1 2.0 * 1 LOSER FIELD 90U’13 ~PPE9 FiflO 80250 NLM9ER OF DATA PCISTS CPEAA FACT3R CFPER 50USD OF COS)THETAI SF5-ES 90U60 IF PHI F1HETA SF1-I LOSE $5155

2

LISE4IIDTH

—

51R1571074

*

2230. SOUlS 4030C GAUSS 160 91350501100.C HZ 2.5 —5.011000 .~ • P1 282 282 OSLiSSIAN

SFECTRUM

5-8~.

1 2

O 8 O

.21522E.Q2 C. 1.

2. .21A3IE4C~. 1.

2 C. C. .2’.203E551

5

7

2

.3Oi2CE—C2 —3. —0.

—0. .111105—73 —2.

—3. —1. .1030115—Il

H0(ISPUT)

C~TE9Tl A 5 9

Z

S

—C.

2.

—2. —C.

—5. —0.

Z C. —7. —2.

C- 059 SOS S

U

5

—5.

—0.

-‘

—0. —3.

—0. —0,

O

Z C. —1.

10.30 20.00

6 5

5)1) 5-52) 6)3) 1-4)4) 801380S

0 5740 0

XX

DY

ZZ

3136.755

3252.073

2679.311

OCS

FIELD 5 74 2

MU

COLCULATEC 4~EI0N8NCE FIELUS 18 lOUIS 15-0L3O1740 SECONO OSCAR EFFACT, EXOLUDOIG

C—IENSSF 11

LINEWI078

.819

.926

.074

3182.655 3149.352 3318.506 3088.638

3015.951 3U53.929 3053.735 35,5. 081

2813.573 2723.676 2635.427 2545.327

OF

PARTIAL SPECTRA

PLOTTED

8S

PLOT N8R.

FRC8

2

2680.008 2900.000

TO 2900.000 3458.000

—

EFFECTS