The influence of phosphorus content on magnetic anisotropy in ferromagnetic (Ga, Mn) (As, P)/GaAs thin films

For the past two decades, research in diluted magnetic semiconductors (SCs) has attracted tremendous interest because of their great potential for different spintronic applied devices [1–3]. Ferromagnetic SCs might combine storage, detection, logic and communication capabilities on a single chip to produce multifunctional devices such as magnetic sensors, spin injectors, and magnetic memories. Among them, the Ga_1−xMn_xAs compounds have received a lot of attention (see [4] and references therein). These alloys exhibit many interesting magnetic properties, such as a relatively high Curie temperature [5] and complex magnetic anisotropy (MA) that determines the easy magnetization axes in the absence of an external field [6–10]. Since information is stored in the direction of magnetization in magnetic storage devices, switching between several easy axes and precession of the magnetization have been considered for information storage and/or processing [11, 12].

As a consequence of the SC band structure (BS) of Ga_1−xMn_xAs, the MA is strongly dependent on Mn and carrier concentrations, temperature, and strain [10, 13, 14]. This dependence allows the use of various static or quasi-static methods to intentionally manipulate the magnetic easy axes. The control of the desired direction of the easy magnetization (in-plane or out-of-plane) has been achieved by (i) changing the Mn concentration [6, 15], (ii) changing the carrier concentration by hydrogenation [16], (iii) varying the temperature [6, 17], (iv) choosing different substrates and therefore different growth-induced strain in the magnetic layer [17, 18], (v) bonding the sample to a piezo-transducer and applying a voltage to change the strain [19–21], and (vi) inducing strain relaxation by lithography patterning [12, 22]. Furthermore dynamical modification of the MA by strain was achieved by using short strain pulses to induce a transient change of the MA in order to launch magnetization precession [23, 24].

The prediction of MA under various conditions of strain, Mn and carrier concentrations and its comparison with experimental results has received a lot of attention in Ga_1−xMn_xAs. Strain control of MA in Ga_1−xMn_xAs_1−yP_y, the focus of this paper, has so far been explored much less extensively. Introducing a few per cent of phosphorus at the arsenic sites in Ga_1−xMn_xAs layers modifies the epitaxial strain continuously from compressive to tensile independently of the Mn doping level [25–27]. One can thereby tune the magnetic easy axis from in-plane to out-of-plane simply by increasing the phosphorus content. Moreover, these layers contain much fewer defects than Ga_1−xMn_xAs layers designed to produce normal-to-plane magnetic easy axis, such as Ga_1−xMn_xAs grown on a relaxed (Ga, In)As buffer layer. Indeed, the self-organization of magnetic domains in periodic patterns, impossible to reach in (Ga, Mn)As/(Ga, In)As, was observed in (Ga, Mn)(As, P)/GaAs [28]. These quaternary layers also have other important advantages over Ga_1−xMn_xAs. For example, it has been predicted that the decrease in the lattice constant of Ga_1−xMn_xAs_1−yP_y should lead to an increase in the Curie temperature and the exchange integral J_pd [29], which was partly confirmed by experiments [28].

As far as we know, although several works were devoted to the study of Ga_1−xMn_xAs_1−yP_y [25–31], the phosphorus composition dependence of the Ga_1−xMn_xAs_1−yP_y band parameters and the link between MA and the BS has not been theoretically investigated. Let us note that one cannot simply map onto Ga_1−xMn_xAs_1−yP_y the calculations performed for Ga_1−xMn_xAs. These generally use the band parameters of GaAs. Here one needs to input the BS of the GaAs_1−yP_y alloy.

In a recent article [32], we have investigated the effect of epitaxial strain on the MA in a ferromagnetic (Ga, Mn)(As, P) thin film with a typical composition of Mn 7% and P 7%. One of the remaining open questions is the influence of phosphorus fraction on the electronic and magnetic properties of (Ga, Mn)(As, P).

This paper is organized as follows. Section 2 presents the calculation of the energy BS and band parameters of strained Ga_1−xMn_xAs_1−yP_y/GaAs within a 40-band k⋅p model using the mean-field Zener model [13]. The influence of subsequent approximations to the electronic structure (6-band k⋅p model) is discussed. The carrier masses and Luttinger parameters are determined as a function of the phosphorus fraction. In section 3, we present a quantitative comparison of the anisotropy parameters obtained from ferromagnetic resonance (FMR) experiments on a set of samples with similar Mn concentration and varying P concentration with our model calculation of the MA using the BS determined in section 2. Section 4 of this work is devoted to the conclusion and perspectives.

2.1. Methodology

In order to obtain the dependence of the MA as a function of phosphorus fraction, the initial step is the BS calculation of Ga_1−xMn_xAs_1−yP_y bulk SCs using k⋅p theory. For this purpose, we use a 40-band k⋅p Hamiltonian first introduced by Saïdi et al [33]. It takes into account the (s, p, d) levels both in the valence bands (VBs) and conduction bands (CBs). We denote by H₄₀ the Hamiltonian matrix. Our 40 × 40 k⋅p Hamiltonian involves matrix elements between the ten $({\Gamma }_{3\mathrm{d}}^{\mathrm{C}}+{\Gamma }_{5\mathrm{d}}^{\mathrm{C}})$, ten $({\Gamma }_{3\mathrm{d}}^{\mathrm{V}}+{\Gamma }_{5\mathrm{d}}^{\mathrm{V}})$, six Γ_5C, six Γ_5V, two Γ_1C, two Γ_1V, two Γ_1u, and two Γ_1q wavefunctions. A full expression of the 40 × 40 Hamiltonian matrix is explicitly given in [34]. This model depends on 10 Γ-centered eigenvalues, 18 spin–orbit coupling coefficients, and 30 matrix elements. In [33], the energy gaps, the momentum matrix elements, and the spin–orbit interaction coefficients used within this model are given in detailed form.

When the crystal is strained, additional terms must be introduced in the k ⋅p Hamiltonian [35–37]. In the case of a biaxially strained Ga_1−xMn_xAs_1−yP_y layer grown on a [001]-oriented cubic GaAs buffer, the biaxial strain tensor can be written as a diagonal tensor with elements _xx,_yy,_zz [10, 13].

In this work, we have adopted the notations of Glunk et al [9] for the relaxed lattice parameter and strain. Hence, ${\epsilon }_{z z}=\frac{{a}_{\perp }-{a}_{\mathrm{rel}}}{{a}_{\mathrm{rel}}},{\epsilon }_{x x}={\epsilon }_{y y}$ and ${\epsilon }_{z z}=-2\frac{{C}_{1 2}}{{C}_{1 1}}{\epsilon }_{x x}$. The relaxed lattice parameter a_rel of a biaxially strained layer on a [001]-oriented substrate is given by

$$\begin{eqnarray} &&{a}_{\mathrm{rel}}=\frac{2{C}_{1 2}}{{C}_{1 1}+2{C}_{1 2}}{a}_{\parallel }+\frac{{C}_{1 1}}{{C}_{1 1}+2{C}_{1 2}}{a}_{\perp }, \end{eqnarray} \tag{ 1 }$$

where a_∥ and a_⊥ are the lateral and vertical lattice parameters of the layer. C₁₁ and C₁₂ are elastic stiffness constants of bulk GaAs_1−yP_y. They are related to the phosphorus content y by (in 10¹¹ dyn cm⁻²) [38]

$$\begin{eqnarray} &&\begin{array}{@{}c@{}} \displaystyle {C}_{1 1}=1 2.3 7+2.1 8 9 5 y+0.2 8 4 5{y}^{2}\\ \displaystyle {C}_{1 2}=5.2 6+0.9 3 1 0 y+0.1 2 1 0{y}^{2}. \end{array} \end{eqnarray} \tag{ 2 }$$

Using equation (1), the vertical strain _zz can be written as

$$\begin{eqnarray} &&{\epsilon }_{z z}=\mathrm{lm}\frac{1-\frac{{C}_{1 1}}{({C}_{1 1}+2{C}_{1 2})}}{1+\frac{{C}_{1 1}}{({C}_{1 1}+2{C}_{1 2})}\mathrm{lm}}, \end{eqnarray} \tag{ 3 }$$

where $\mathrm{lm}=\frac{({a}_{\perp }-{a}_{\parallel })}{{a}_{\parallel }}$ is the lattice mismatch determined by high-resolution x-ray diffraction.

The characteristics of the samples studied in this paper are summarized in table 1. The seven Ga_0.93Mn_0.07As_1−yP_y samples are selected in order to compare layers with the same nominal Mn concentration (7%) and different phosphorus concentrations. The samples were grown by low-temperature molecular-beam epitaxy on GaAs(001) substrates. The Mn and P concentrations were calibrated using Vegard's law, based on the lattice constants of a series of GaAs_1−yP_y references samples grown under the same conditions. After growth, the samples were systematically annealed under nitrogen atmosphere at 250 ° C for 1 h [25]. The data listed in table 1 clearly show that the strain changes from compressive (_zz > 0) to tensile (_zz < 0) when the phosphorus composition increases in the 50 nm thick Ga_1−xMn_xAs_1−yP_y epilayer.

Table 1.  Ga_1−xMn_xAs_1−yP_y sample parameters: the phosphorus concentration [P], the lattice mismatch lm, the corresponding strain component _zz, the saturation magnetization M_S at T = 4 K, the effective Mn concentration x_eff, and the splitting parameter of the valence subbands B_G.

Sample	A	B	C	D	E	F	G
[P] (%)	0	2.6	3.4	4.3	5.6	7	11.3
lm (ppm)	4990	280	0	−1200	−2320	−3580	−6620
zz (%)	0.23	0.013	0	−0.055	−0.11	−0.16	−0.31
MS (emu cm−3)	37.6	37.7	37.3	40.75	42.6	39.2	35
xeff (%)	3.69	3.67	3.63	3.96	4.14	3.80	3.49
BG (meV)	−26.0	−26.3	−26.0	−28.5	−30.0	−27.5	−24.7

The Bir–Pikus H_S Hamiltonian matrix used in this work is given explicitly in [34]. Because the influence of strain on p-type CB is unknown, we have taken the p-type CB hydrostatic deformation potential (a_Γ5C) and the p-type CB shear deformation potential (b_Γ5C) equal to zero. Moreover, we do not use the Bir–Pikus Hamiltonian for the s-type VB, which is at about −12.5 eV from the top of the VB, because (i) no parameters are available for this band and (ii) its influence on the CB and on the VB of interest is very small. Thus, we consider a_Γ1V = 0. In the same way there is no need to apply the strain to the d-type VBs and CBs.

The presence of Mn ions in Ga_1−xMn_xAs_1−yP_y is taken into account by including the s, p–d exchange interactions. Hence, H_exc = H_pd + H_sd describes the exchange coupling between the localized d-electrons and the delocalized holes (p) and electrons (s). H_pd describes the p–d hybridization, which yields an effective exchange interaction between the hole spin s and the Mn spin S (S = 5/2) carrying a magnetic moment −Sgμ_B, where g is the Landé factor and μ_B the Bohr magneton. The exchange interaction can be further simplified by the mean-field approximation where the virtual crystal approximation is made [2, 13]. The conduction and valence bands exchange splittings are obtained from the exchange energies N₀α and N₀β, where N₀ is the concentration of cation sites and α and β are the s–d, and p–d exchange integrals, respectively. In order to obtain the values of N₀β we assume ${N}_{0}\sim 4/{a}_{\mathrm{rel}}^{3}$ and β =− 54.5 meV nm³. a_rel varies from 5.6661 to 5.6359 Å and N₀β from −1.20 to −1.22 eV for y from 0 to 11.3%. Since there is no report concerning experimental or theoretical determination of the s–d exchange energy N₀α the value 200 meV, as observed for II–VI diluted magnetic SCs [39], has been adopted in this work. Owing to the presence of Mn atoms in interstitial sites the effective concentration of Mn atoms participating in the ferromagnetism is lower than the total one [40]. x_eff is deduced from the low-temperature saturation magnetization M_S as x_eff = M_S/(N₀Sgμ_B). Table 1 clearly shows that only about 50% of Mn atoms contribute to the magnetization for all samples, revealing that the incorporation of Mn in interstitial sites is not modified by alloying with phosphorus. One should also note that, in the Zener model, the exchange splitting of the valence subbands B_G = A_FβM_S/6gμ_B is an important parameter for the BS and the MA. The values used in our calculations are given in table 1. The A_F parameter describes the enhancement of spin density of states by the carrier–carrier interactions. The best fit with experimental results is obtained with A_F = 1.4 (close to the value 1.2 used in [13]).

According to the k ⋅p effective Hamiltonian theory, the total Hamiltonian H of the system used in this work has the general form: H = H_kp + H_S + H_exc. Let us note that we have not taken into account the bandgap renormalization arising from the high density of carriers in the (III–Mn)V compounds since we are interested in MA and not in optical transitions.

2.2. Results

In this section we present and discuss the results of the application of our k ⋅p model to the electronic structure calculations for bulk Ga_1−xMn_xAs_1−yP_y SC. All these calculations were carried out at the relaxed lattice constant. An example of the BS calculated with the full Hamiltonian H for x_eff = 0.0369 and y = 0 is shown in figure 1. The electronic BS of this compounds shows features similar to that of bulk GaAs. The band diagram is well reproduced on a width of about 18 eV. The extended k⋅p method is known to be valid up to 6 eV above and 12 eV below the top of the VB in four directions, namely ΓX,ΓL,ΓK and XU.

Zoom In Zoom Out Reset image size

One of the important effects of the strain and exchange is observable in the VB near-Γ region. In figures 2(a)–(c), we show a magnified BS near k = 0, for three different phosphorus fractions (a) y = 0%, (b) y = 3.4%, and (c) y = 11.3%, corresponding to the regime of compressive, zero, and tensile strain, respectively. It is clear that the degeneracy between heavy-holes (HHs) and light-holes (LHs) at the center of the Brillouin zone is lifted. Moreover, the band energies shift evidently near the Γ point, while there is a smaller shift on other VB regions. The influence of the P concentration is mostly evidenced along the [111] direction (figure 2(c)).

Zoom In Zoom Out Reset image size

Generally, due to the presence of the strain, the BS of SCs is altered, which changes the lattice constant and reduces the symmetry of the crystal. It also modifies the energy gaps, and lifts degeneracy. The influence of subsequent approximations on the electronic structure is considered. In figure 2(d) we show the comparison of the VB dispersion obtained from the 40-band model (solid line) and that from the 6-band model (dashed line) for the case of Ga_0.9631Mn_0.0369As over the [001] and [111] directions. This figure shows a satisfactorily agreement between the two models only up to k ∼ 0.1 Å⁻¹ measured from the Γ edge. For a Fermi energy in (Ga, Mn)As up to ≈100 meV the 6-band model is probably sufficient. For higher Fermi energies (doping level above ≈5 10²⁰ cm⁻³) the deviation of the 6-band k⋅p model BS from the real BS might lead to an inaccurate estimation of the MA. Let us emphasize that, within the 6-band model, one needs to specify the Luttinger parameters in order to take into account the interaction of the VB with the remote bands. This can be achieved for (Ga, Mn)As by taking the parameters of GaAs, but these parameters are not available for (Ga, Mn)(As, P). In the framework of the 40-band model the coupling between bands is described in terms of the matrix elements and band edges, which are deduced from the experiments. Once the BS is well reproduced with the 40-band model, we then calculate the Luttinger parameters [33] which can be used within the 6-band model. This is an important advantage of the 40-band model compared to the 6-band model.

To investigate the effect of phosphorus fraction on the electronic properties of alloys, we also performed calculations of energy gaps at high-symmetry points in the Brillouin zone. The bandgap energy at Γ point, ${\mathbf{E}}_{\mathrm{g}}^{\Gamma }$, is defined as the difference between the lowest conduction band and the highest VB in the presence of strain and exchange. The variation of ${\mathbf{E}}_{\mathrm{g}}^{\Gamma }$ as a function of the phosphorus content in the quaternary alloy layers Ga_1−xMn_xAs_1−yP_y/GaAs is displayed in figure 3. This variation results from four contributions: (i) the increase of the gap of the bulk Ga(As, P) with P concentration (alloying effect), (ii) the down-shift of the gap resulting from the exchange splitting of the valence and conduction bands (almost constant here since N₀β and x_eff have small variations with y), (iii) the effect of the hydrostatic part of the strain, which increases the gap for compressive strain and decreases it for the tensile one, thus opposing the alloying effect, and finally (iv) the effect of the biaxial strain on the VB, which gives a variation similar to the hydrostatic component but with a smaller amplitude. The increase in the gap seen in figure 3 results from the alloying effect (i). The plateau near y ≈ 4% results from the cancellation of strain contributions (iii) and (iv) at zero strain together with a slightly increased exchange down-shift (13% larger for y = 6% than for 0% and 11.3%).

Zoom In Zoom Out Reset image size

We proceed to understand the role of phosphorus fraction on the effective masses and Luttinger parameters. For that purpose we calculate the effective masses for the s-type CB (m_e(Γ_6C)) and for p-type VBs (m_hh,m_lh) using the total Hamiltonian H. m_e(Γ_6C) was obtained from the second derivative of s-type CB energy with respect to the wavevector around the Γ point. As can be seen in figure 4(a), the electron effective mass of quaternary Ga_1−xMn_xAs_1−yP_y exhibits a non-monotonous dependence on the phosphorus composition. Indeed, m_e(Γ_6C) increases mildly with decreasing compressive biaxial strain, while it decreases rapidly with tensile strain. This non-monotonous behavior reflects a competition between alloying and strain effects. Let us note that the variation of the Mn composition has no effect on the masses since the exchange splitting does not change the band curvature.

Zoom In Zoom Out Reset image size

Generally speaking, these observations are consistent with previous experiments and a theoretical study of other bulk SCs [41]. In figures 4(b)–(c), we show the effective masses for both the HH and LH bands along the [001], [110] and [111] directions as a function of phosphorus fraction. m_hh and m_lh have the same behavior along all three directions. Similarly to the electron mass the HH mass exhibits a non-monotonous behavior with P concentration. It decreases slightly under decreasing compressive strain and increases mildly with tensile strain. The LH mass decreases monotonously with increasing y but the slope of m_lh(y) changes when going from compressive to tensile strain. Note finally that the m_hh effective mass is the largest along [111] direction, while m_lh is the largest along the [001] direction.

An accurate description of the VB effective masses makes our results reliable to calculate the Luttinger parameters γ_ℓ (ℓ = 1–3) from the following expression [42]

$$\begin{eqnarray*} &&\begin{array}{@{}l@{}} \displaystyle {\gamma }_{1}=1/2(1/{m}_{\mathrm{lh}}^{0 0 1}+1/{m}_{\mathrm{hh}}^{0 0 1}) \\ \displaystyle {\gamma }_{2}=1/4(1/{m}_{\mathrm{lh}}^{0 0 1}-1/{m}_{\mathrm{hh}}^{0 0 1}) \\ \displaystyle {\gamma }_{3}=1/4(1/{m}_{\mathrm{lh}}^{0 0 1}+1/{m}_{\mathrm{hh}}^{0 0 1})-1/2{m}_{\mathrm{hh}}^{1 1 1} \end{array} \end{eqnarray*}$$

and shows that the behavior of the Luttinger parameters is a consequence of the competition between the alloying effect and the strain effect. In figure 5, we show the dependence of the Luttinger parameters on y and compare them to data from the literature for unstrained GaAs_1−yP_y alloys. The full lines show the y-variation obtained by a linear interpolation between the GaAs and GaP values [43]. The decrease in γ_ℓ (ℓ = 1–3) with increasing y (full lines) can be explained by an important alloying effect. For our calculated values (symbols), the effect of strain differs depending on the regime of strain considered. Indeed, comparatively to the unstrained case, the compressive strain leads to a decrease of the Luttinger parameters while the tensile strain leads to an increase of γ_ℓ.

Zoom In Zoom Out Reset image size

Magnetic anisotropy has attracted much attention because it is essential not only for the fundamental understanding of the microscopic origins of ferromagnetism in diluted magnetic SCs, but also it is of central importance for applications of these materials, since it determines the easy axis of magnetization, the anisotropy fields and the magnetic domain wall width and specific energy. Let us recall that the MA is the dependence of the free-energy density F of a magnetic system on the orientation of the magnetization direction m = (M/M) with respect to the symmetry axes. Recent studies of the ferromagnetic phase in Ga_1−xMn_xAs_1−yP_y epilayers demonstrated the existence of rich MA [25, 26, 30]. In order to provide a concise description of the MA, the normalized quantity F_M = F/M is considered instead of F. As given in [44–46], there are three contributions to F_M: (i) the contribution F_M,1(m) that arises from spin–orbit coupling in the VB in a cubic crystal with exchange and biaxial strain, (ii) the contributions F_M,2(m) caused by the demagnetization field and a uniaxial in-plane contribution along $[1\overline{1}0]$, whose microscopic origin is still under debate [6, 47–50], and (iii) the contribution F_M,3(m) from the single-ion anisotropy [46, 51], which is small and usually neglected. The free-energy functional F_M is then given by F_M(m) = F_M,1(m) + F_M,2(m) + F_M,3(m). We consider a phenomenological description of the free energy in terms of anisotropy constants and the direction cosines of the magnetization vector m_x,m_y, and m_z with respect to the main crystallographic directions 〈100〉. We restrict our calculations to zero temperature and zero magnetic field. In this approximation and taking only up to the fourth-order term for the m_i components, the free-energy density for [001]-oriented layers can be written as

$$\begin{eqnarray*} {F}_{\mathrm{M}}(\mathbf{m})&=&{B}_{2\perp }{m}_{z}^{2}+{B}_{4\parallel }({m}_{x}^{4}+{m}_{y}^{4})+{B}_{4\perp }{m}_{z}^{4}+{B}_{\mathrm{ d}}{m}_{z}^{2}\nonumber\\ &&+~\frac{{B}_{\overline{1}1 0}}{2}({m}_{x}-{m}_{y})^{2}+{B}_{2\perp }^{\mathrm{si}}{m}_{z}^{2}\nonumber\\ &&+~{B}_{4\parallel }^{\mathrm{si}}({m}_{x}^{4}+{m}_{y}^{4})+{B}_{4\perp }^{\mathrm{si}}{m}_{z}^{4}, \end{eqnarray*}$$

where ${B}_{\mathrm{d}}=\frac{{\mu }_{0}M}{2}$. The first three anisotropy components represent contributions to F_M,1. The next two contribute to F_M,2, and the last three to F_M,3. The parameters B_ℓ denote the magnetic anisotropy constants. $({B}_{2\perp }+{B}_{2\perp }^{\mathrm{si}})$ is the perpendicular uniaxial anisotropy field resulting from the biaxial strain introduced during the molecular-beam epitaxy growth of GaMnAs layers on substrates with a different lattice constant. $({B}_{4\parallel }+{B}_{4\parallel }^{\mathrm{si}})$ and $({B}_{4\perp }+{B}_{4\perp }^{\mathrm{si}})$ are constants representing the in-plane and the perpendicular cubic anisotropy which reflects the non-equivalence of the 〈100〉 and 〈110〉 directions. In the case of a perfect cubic crystal, the symmetry requires B_4∥ = B_4⊥ and ${B}_{4\parallel }^{\mathrm{si}}={B}_{4\perp }^{\mathrm{si}}$. Since growth on a lattice-mismatched substrate introduces tetragonal distortion and breaks the purely cubic symmetry of the crystal lattice, the in-plane and perpendicular cubic anisotropy constants need to be distinguished. B_d is the demagnetizing anisotropy field, arising from shape anisotropy. ${B}_{\overline{1}1 0}$ or B_2∥ is the in-plane uniaxial anisotropy field that reflects the non-equivalence of $[1\overline{1}0]$ and [110] directions.

In our theoretical approach, we have calculated the F_M,1 contribution. The numerical procedure to determine the magnetic anisotropy is described in [9, 32]. Using the BS calculated in section 2, F₁ is obtained first by summing over all energy eigenvalues within the four spin-split, HH and LH Fermi surfaces. The resulting energy density is then normalized by the saturation magnetization M. We use the sample parameters specified in table 1. Due to the presence of charge compensating defects depending on growth conditions in Ga_1−xMn_xAs_1−yP_y samples, the exact determination of the Fermi level in ferromagnetic Ga_1−xMn_xAs_1−yP_y is difficult. Hence, for reasons of simplicity and lack of data on the hole concentration, we have assumed that the Fermi level has a constant value equal to −110 meV from the top of the VB for all samples except sample G. In order to justify this approximation, we calculate the hole concentration p_theo for the chosen Fermi energy, and compare it to the estimated hole concentration p_est. The values of p_est are obtained by assuming that the reduction of the hole concentration is only due to the presence of interstitial Mn that act as double donors. We have x_eff = x − 2x_i, where x_i is the interstitial Mn content. Consequently, the hole concentration is p_est = N₀(x − 3x_i) = N₀(3x_eff − x)/2 [52]. The values of p_theo and p_est are listed in table 2 for the studied samples. We can see that taking the Fermi level as constant is a fairly reasonable approximation. Let us note that sample G (y = 11.3%) is close to the transition from the metallic to the impurity band conduction regime, as deduced from the temperature dependence of the resistivity. As argued in [2], our model might nevertheless be applicable. Indeed, using a slightly smaller Fermi energy (−81 meV), in agreement with the estimated carrier concentration p_est (table 2), we will see that a good agreement between theoretical and experimental MA fields is obtained.

In order to compare the microscopic theory with the phenomenological description of the MA, we consider the anisotropic part ΔF_M,1 of the F_M,1 contribution. When the reference magnetization direction is taken along [100],ΔF_M,1 is given by the following expression

$$\begin{eqnarray} \Delta {F}_{\mathrm{M},1}&=&{F}_{\mathrm{M},1}(\mathbf{m})-{F}_{\mathrm{M},1}(\mathbf{m}=[1 0 0])\nonumber\\ &=&{B}_{2\perp }{m}_{z}^{2}+{B}_{4\parallel }({m}_{x}^{4}+{m}_{y}^{4}-1)+{B}_{4\perp }{m}_{z}^{4}. \end{eqnarray} \tag{ 4 }$$

For m rotated in the (001) and in the (010) plane, we can write equation (4) as [9]:

$$\begin{eqnarray} &&\Delta {F}_{\mathrm{M},1}(\varphi )={B}_{4\parallel }({\cos }^{4}(\varphi )+{\sin }^{4}(\varphi )-1) \end{eqnarray} \tag{ 5 }$$

and

$$\begin{eqnarray} \Delta {F}_{\mathrm{M},1}(\theta )&=&{B}_{2\perp }{\cos }^{2}(\theta )+{B}_{4\parallel }({\sin }^{4}(\theta )-1)\nonumber\\ &&+~{B}_{4\perp }{\cos }^{4}(\theta ). \end{eqnarray} \tag{ 6 }$$

Once ΔF_M,1 is calculated numerically in the microscopic model as a function of the crystalline orientations (θ,φ) for the different phosphorus fractions, the phenomenological description of MA, namely equations (5) and (6), can be fitted to the resulting angular dependences by using B_2⊥,B_4∥, and B_4⊥ as fit parameters.

The results of the microscopic calculations are shown by full symbols in figures 6(a) and (b). The full lines show least-squares fits to the calculated data using equations (5) and (6). We observe that the fit curves are in very good agreement with the microscopic results. This good fit clearly demonstrates that the microscopic theory can be well described by B_2⊥,B_4∥, and B_4⊥. Therefore a realistic description of the F_M,1 part of the free energy only needs terms up to the fourth order in m. We proceed to understand the role of the phosphorus fraction on the magnetic properties of Ga_1−xMn_xAs_1−yP_y quaternary alloys. For that purpose we calculate ΔF_M,1(φ) and ΔF_M,1(θ) for three representative P concentrations: y = 11.3% (tensile strain), y = 3.4% (≈ zero strain), and y = 0.0% (compressive strain), respectively. Figure 6(a) shows that in Ga_1−xMn_xAs_1−yP_y there is a stronger fourfold cubic in-plane anisotropy constant B_4∥ under tensile strain than under compressive strain. Such a behavior has also been shown in the case of Ga_1−xMn_xAs ferromagnetic SC [9, 53].

Zoom In Zoom Out Reset image size

From figure 6(b), we observe that in the regime of tensile strain ΔF_M,1 exhibits a minimum for m oriented along [001] (by symmetry also along $[0 0\overline{1}]$). For compressive strain, the minima occur for m oriented along $[1 0 0],[\overline{1}0 0]$, [010], and $[0\overline{1}0]$. Thus, the theoretical model correctly describes the well-known experimental fact that the magnetically hard axis along [001] in compressively strained layers turns into an easy axis in tensile strained layers, as observed in [26]. Note, however, that for the regime of zero strain, the easy axis is not aligned along any high-symmetry direction. This result is due to the values approximately equal taken by the cubic anisotropy parameters (B_4∥,B_4⊥). These observations clearly point out that the easy axis of magnetization is sensitive to the variation of the phosphorus fraction or, in other words, the preferential axes of the magnetization are closely related to the strain conditions.

Applying the phenomenological procedure described above to the whole set of Ga_1−xMn_xAs_1−yP_y/GaAs layers under study, the anisotropy parameters B_2⊥,B_4∥, and B_4⊥ were determined as a function of phosphorus fraction.

Let us now consider the single-ion anisotropy contribution F_M,3. The uniaxial term ${B}_{2\perp }^{\mathrm{si}}$ can be estimated from the D parameter of the single-ion Hamiltonian used in the analysis of electron paramagnetic resonance experiments (EPR) [46, 51] as ${B}_{2\perp }^{\mathrm{si}}=\frac{D S}{g{\mu }_{\mathrm{B}}}$. Extrapolating the variation of D obtained from EPR for very low Mn concentrations [46] would give ${B}_{2\perp }^{\mathrm{si}}=9~\mathrm{mT}$ for Sample A. This is 5% of the value of B_2⊥. The ${B}_{2\perp }^{\mathrm{si}}$ contribution might even be much smaller since the extrapolation was found to overestimate the D value by a factor of ≈4 for a concentration as low as 0.0027 [51]. The value of D is not known for the same Mn spin state in GaP. However D is proportional to the spin–lattice coupling coefficient G₁₁, which is at least three times smaller in GaP than in GaAs [54]. Therefore we expect that the ${B}_{2\perp }^{\mathrm{si}}$ contribution in (Ga, Mn)(As, P) (for the same absolute value of strain) would be even smaller. There might also be a cubic contribution from the single-ion anisotropy. From the a parameter obtained from EPR experiments [46, 51], a ${B}_{4}^{\mathrm{si}}$ contribution can be estimated as $\frac{a{S}^{3}}{6 g{\mu }_{\mathrm{B}}}$, which amounts to 4 mT. Given the smallness of the single-ion contribution to anisotropy, in the following the experimental results will be compared to the carrier contributions B_2⊥,B_4∥, and B_4⊥.

The experimental magnetic anisotropy constants were obtained by a standard X-band FMR spectroscopy with a 9 GHz microwave source. The angular dependence of the FMR spectra were measured by rotating the static magnetic field in two crystallographic planes: (110) and (001) named 'out-of-plane' and 'in-plane configuration', respectively. These two sets of angular variations of the resonance field enable us to determine the anisotropy constants and associated anisotropy fields [44].

We have taken the experimental data at T = 4.2 K in order to compare with our calculations performed at zero temperature. Firstly, we focus on the behavior of the B_2⊥ parameter. In figure 7(a), B_2⊥ is plotted versus phosphorus fraction y. The microscopic calculations are in very good agreement with experimental data. However, the results of B_2⊥ deduced from the microscopic model are found to be slightly different from the experimental data for y = 3.4% and y = 4.3%. The microscopic model predicts that B_2⊥ is linearly proportional to the uniaxial deformation _zz. Indeed, a zero value of the uniaxial out-plane parameter is found for zero strain (y = 3.4%). However a non-zero value is obtained experimentally. This discrepancy might originate from the experimental difficulty in determining the strain in this range giving a large error bar at zero strain: ±600 ppm. Additionally, it might be necessary to include higher order terms in F_M,1 for the analysis of experimental results when B_2⊥,B_4∥, and B_4⊥ are close to zero. For the range y = 0–7% an almost linear decrease of B_2⊥ is evidenced. For larger concentrations 7.0% < y < 11.3%,B_2⊥ depends less strongly on y. This result can be explained by the behavior of Luttinger parameters in the range 7.0% < y < 11.3% and the lower Fermi energy for y = 11.3%.

Zoom In Zoom Out Reset image size

The fourth-order parameters B_4∥ and B_4⊥ are presented in figures 7(b) and (c), respectively. We can see that the experimental data of B_4∥ are very well reproduced by the microscopic calculations. The B_4∥ values obtained using the 40-band k⋅p model in (Ga, Mn)(As, P) are in much better agreement with experimental values than in ternary (Ga, Mn)As using the 6-band k⋅p model [9]. We believe that the more accurate description of the MA originates from a better description of the bands with our 40-band k⋅p model in the range used for the calculation of the MA, in particular from a good determination of the Luttinger parameters. Besides, the assumption that consists in taking the Fermi level as a constant can be considered as relevant because of the almost perfect agreement between experiment and theory.

Unlike the results of B_4∥, those of B_4⊥ (figure 7(c)) show large variations between the microscopic values and the experimental ones. This discrepancy becomes less important for samples which present close experimental values of B_4∥ and B_4⊥ (samples B, D, E, and G). However, the predictions of microscopic calculations could not explain the large variation observed in samples A, C and F. Moreover, we can notice that the change of sign for B_4⊥ is not reproduced by microscopic calculations. The same problem has also been reported in [9] for Ga_1−xMn_xAs. At the moment, the reason for this discrepancy is not yet understood. We cannot exclude the possibility that the change of sign of B_4⊥ is an artifact of the experimental method for determining B_4⊥ [9, 44]. In samples with high B_2⊥ coefficient the experimental determination of B_4⊥ might be unreliable [10].

To summarize, we have presented a theoretical survey of the magnetic properties of ferromagnetic Ga_1−xMn_xAs_1−yP_y/GaAs epilayers using the mean-field approximation and the 40-band k⋅p model, in which the BS and band parameters are calculated. The present k⋅p description accurately reproduces the overall band diagram, as well as the band shifts and the carrier effective masses versus phosphorus content. We have considered the effect of the phosphorus composition on the MA. The present microscopic model has been validated through an accurate set of comparisons with experimental results. In particular, the values of B_2⊥ and B_4∥ anisotropy parameters are very well reproduced.

In the next step, it would be instructive to examine the effect of the introduction of shear deformation in our microscopic calculations. Firstly, the anisotropy coefficient B_2∥, which may be associated with magnetostriction or may result from anisotropic Mn distribution [6, 10, 50], plays an important role in the orientation of the in-plane magnetic easy axes. This in-plane uniaxial MA has the same symmetry as the one that would arise from an _xy shear strain although never evidenced experimentally. Secondly, the recent demonstration of magnetization precession triggered by quasi-transverse picosecond strain pulses emphasizes the importance of shear deformation for magnetization manipulation [55].

The 40-band k ⋅p model also offers the possibility to calculate the optical properties, in particular the dielectric permittivity and the spectral dependence of the magneto-optical Kerr effect, which is of great importance for the optimization of optical pump–probe experiments and magneto-optical imaging of magnetic domains in these ferromagnetic SCs.

The authors thank Laura Thevenard from INSP for useful discussions. This work was performed in the framework of the MANGAS project (ANR 2010-BLANC-0424-02).