text:
# Kinetic model of II-VI(001) semiconductor surfaces: Growth rates in atomic layer epitaxy

## I Introduction

Molecular beam epitaxy (MBE) and its variations continue to attract significant interest as a technique for the production of, for instance, high quality semiconductor films. In addition to its practical relevance, epitaxial growth is highly attractive from a theoretical point of view. It offers many challenging open questions and provides a workshop in which to develop novel tools for the modeling and simulation of non-equilibrium systems in general. Reviews of the experimental techniques and theoretical investigations can be found in, e.g., Refs. pimpinelli, ; corfu, .

II-VI semiconductors, as a promising class of materials, have been the subject of many experimental studies, see Ref. zwosechs,  for a recent overview. Among other aspects, competing surface reconstructions and their interplay with growth or sublimation have been investigated. tatarenko ; neureiter ; faschinger ; hartmann ; cibert ; tatarenko2  As an example we will concentrate on the CdTe(001) surface, in the following, but many of the results should apply to other II-VI compounds as well. The focus will be on atomic layer epitaxy (ALE), a technique which aims at self-regulated layered growth by alternating deposition of the elements in ultra high vacuum. Specific properties of the material system and the attractive clarity of the ALE technique allow for an efficient theoretical modeling of the growth scenario.

In earlier investigations we employed two-dimensional lattice gas models of the CdTe(001) surface epl ; flat  for the investigation of the temperature and flux dependence of surface reconstructions. Chemical bonds and the influence of the underlying crystal were accounted for by effective pairwise anisotropic interactions. Sublimation and ALE growth have been studied in the framework of a simple cubic model,surflett  which accounted only for the most essential features of CdTe. The correct zinc-blende lattice was first considered in Ref. aachen, , but here we extend the model significantly by taking into account the formation of Te-dimers at the surface. More importantly, we will consider weakly bound excess Te, which provides a reservoir of highly mobile atoms. A similar approach with weakly bound arsenic has recently been suggested by Itoh et al.itoh98 ; itoh01  in order to describe the growth on reconstructed GaAs(001) and InAs(001) surfaces.

Within the framework of our model we demonstrate how the surface reconstruction and limited Cd-coverage restrict the ALE growth rate at high temperatures. Our model shows how this limitation is overcome if weakly bound excess Te is present on the surface at low temperatures. We observe a sudden increase of the growth rate with decreasing temperature, in qualitative agreement with experiments.

Due to the increased number of parameters compared with our earlier models a quantitative match with experimental observations is hard to achieve at the present stage of the model. One possibility to overcome this problem is to combine density functional calculations (DFT) with Kinetic Monte Carlo (KMC) simulations as was done e.g. by Grosse et al.grosseprb ; grosseprl  in the context of growth on III-V(001) semiconductor surfaces. However, the available DFT input for II-VI(001) surfaces, especially CdTe(001), does not suffice for a full parametrization of our growth model.

The paper is organized as follows: before introducing the model and method in section III, we summarize essential properties and experimental findings concerning II-VI(001) surfaces in the next section. In section IV we discuss the simulated growth scenario and present the results. The last section summarizes and gives an outlook on perspective projects based on this work.

## II Properties of the material system

A detailed overview of II-VI(001) surfaces, with emphasis on Te-compounds, can be found in Refs. tatarenko, ; cibert, ; tatarenko2, , for instance. Whereas we will mostly refer to CdTe, many of the features mentioned here apply to other compounds, e.g. ZnTe or ZnSe.

CdTe crystallizes in a zinc-blende lattice and consists of alternating layers of Cd and Te which are parallel to the (001) surface of the crystal. In the bulk, each of these layers forms a regular square lattice, cf. Figure 1. Cd atoms are bound to two neighboring Te atoms in the layer below with bonds oriented along the $`[1\overline{1}0]`$ direction, which will be termed $`y`$-direction for short. Accordingly, Te is always bound to two underlying Cd along the $`[110]`$ or $`x`$-direction.

Under vacuum conditions and in absence of particle deposition, the CdTe(001) surface is terminated by a Cd-layer with limited coverage $`\rho ^{\mathrm{Cd}}1/2`$. At low temperatures $`T`$ one finds a dominant $`c(2\times 2)`$ vacancy structure, which corresponds to a checkerboard like occupation of the available square lattice sites. Above $`T570K`$ the surface is dominated by a local $`(2\times 1)`$ ordering where the Cd atoms arrange in rows along the $`[110]`$ or $`x`$-direction which alternate with empty rows. Despite the fact that significant sublimation sets in even below $`570K`$, the temperature driven re-ordering of the Cd-terminated surface has been discussed in terms of a phase transition in effective equilibrium. neureiter ; cibert ; epl ; flat

Electron counting rules,pashley  and density functional theory (DFT) calculations garcia ; park ; gundel  show that the simultaneous occupation of nearest neighbor (NN) sites in the $`y`$-direction should be extremely unlikely within a terminating Cd layer. Such pairs of NN-Cd would lead to an unfavorable local concentration of positive charges. Note that in both the above mentioned arrangements NN pairs of Cd-atoms in $`y`$-direction do not occur, indeed. DFT calculations show furthermore that the $`c(2\times 2)`$ ordering has a slightly lower surface energy than the competing $`(2\times 1)`$ structure. garcia ; park ; gundel

The $`c(2\times 2)`$ ordered Cd-termination can be restored in the presence of a stationary Cd flux even at high $`T`$.tatarenko ; cibert  On the other hand, a steady deposition of Te stabilizes Te-terminated surfaces.tatarenko ; cibert  At relatively small fluxes and high temperatures a coverage $`\rho ^{\mathrm{Te}}1`$ is observed, hence there is no vacancy structure as for Cd-terminated surfaces. The Te atoms arrange in rows of dimers along the $`[110]`$-direction at the surface, corresponding to a $`(2\times 1)`$ reconstruction pattern.

At high flux and low temperatures Te-coverages $`\rho ^{\mathrm{Te}}3/2`$ with a $`(2\times 1)`$ surface symmetry are observed. It has been hypothesized that Te atoms form trimers,tatarenko ; cibert  but the precise mechanism of binding the excess Te to the surface remains unclear.tatarenko2  The evaporation of terminating Te-layers after ending deposition has been studied in the different temperature regimes.tatarenko ; cibert ; tatarenko2  The analysis of the corresponding Arrhenius-plots supports the hypothesis of a very weak binding of excess Te.

As we will demonstrate in the following, the influence of the excess Te is of particular importance for ALE growth. One ALE cycle, i.e. a sequence of one Cd and one Te pulse, could add a complete layer of CdTe to the system, in principle. However, in a wide range of temperatures, one observes a growth rate of at most $`1/2`$ layer per ALE cycle.faschinger ; hartmann  This limitation reflects the restricted coverage of Cd-terminated surfaces. Only for sufficiently low $`T`$ one finds the expected rate of approximately 1 layer per cycle, with a temperature driven, sudden transition from one regime to the other. faschinger ; hartmann  We will explain within the framework of our model how the above mentioned presence of weakly bound Te atoms at the surface allows to overcome the limitation.

## III The model

### III.1 Representation of the crystal and energetics

As a simplification we exclude the incorporation of defects or vacancies into the crystal. Hence, the bulk can be represented in a Solid-On-Solid (SOS) manner, cf. Figure 1. Any atom with more than two chemical binding partners, i.e. with at least one NN in the layer above, is considered immobile. When measuring the surface coverage $`\rho ^{\mathrm{Cd}/\mathrm{Te}}`$ with Cd/Te atoms we count only those atoms which are mobile.

Every Cd–Te bond in the system contributes an energy $`ϵ_c<0`$. The crystal structure is further stabilized by an effective intra-layer NN interaction $`ϵ_b`$ of immobile particles of the same species.

Effective interactions of the mobile particles at the surface are chosen to reproduce known properties of the system, e.g. surface reconstructions. Consequently, mobile Cd and Te are treated on different footings.

As in the planar models of Refs. epl, ; flat,  an infinite repulsion excludes NN pairs of mobile Cd atoms in $`y`$-direction at the surface. The repulsion does not apply to pairs of one mobile and one immobile Cd, if the latter has three or four Te partners. An attractive interaction energy $`ϵ_x<0`$ between mobile Cd-NN in $`x`$-direction and a competing interaction $`ϵ_d<0`$ between diagonal neighbors (NNN) favor the $`(2\times 1)`$ or $`c(2\times 2)`$ arrangement, respectively. As demonstrated in Refs. epl, ; flat, , the temperature dependence of the Cd reconstruction is reproduced qualitatively, if we set $`\left|ϵ_x\right|\stackrel{~}{<}\mathrm{\hspace{0.17em}2}\left|ϵ_d\right|`$.

NN pairs of mobile Te which do not form a dimer contribute an energy $`ϵ_m^{\mathrm{Te}}<0`$, whereas the energy $`ϵ_b^{\mathrm{Td}}<ϵ_m^{\mathrm{Te}}`$ is assigned to a pair which forms a dimer. Note that dimerization introduces additional degrees of freedom to the model without changing the topology of the lattice. Pairwise interactions $`ϵ_x^{\mathrm{Td}}`$ of neighboring dimers in $`x`$-direction favor their arrangement in rows along $`[110]`$, corresponding to a $`(2\times 1)`$ surface symmetry. For NN pairs of one mobile and one immobile Te we assume the same interaction $`ϵ_b`$ as for pairs of incorporated atoms.

In addition to Te particles at regular SOS lattice sites, we include the capture of Te atoms in a weakly bound state, denoted as Te atoms in the following. The energy contribution of each Te atom is $`ϵ^{}`$ where $`\left|ϵ^{}\right|<2\left|ϵ_c\right|`$. This accounts for the fact that the binding of a Te is weaker than that of an isolated mobile atom at a regular lattice site which has two Cd–Te bonds.

Following an earlier hypothesis,tatarenko ; cibert  these weakly bound Te might be associated with atoms which occupy Cd-sites temporarily and hence form a Te-trimer.ahrdiss  But, as already mentioned in section II the precise nature of the Te is still unknown. Here, we abandon the assumption of Te-trimers and instead include Te particles in an effective fashion.

We introduce a reservoir of Te with $`\theta ^{}N`$ particles, where $`N`$ is the number of atoms in one complete bulk layer of Te. The total reservoir occupation $`\theta ^{}`$ is composed of two parts: $`\theta ^{}=\theta _{loc}^{}+\theta _{mob}^{}`$. The quantity $`\theta _{loc}^{}`$ accounts for the fraction of Te which are localized on top of mobile Cd and thus neutralize the infinite repulsion between Cd-NN in y-direction.ahrdiss  This effect will be discussed in more detail in the next section. On the contrary, $`\theta _{mob}^{}`$ corresponds to the remaining Te atoms which are highly mobile and no assumption is made about their precise location and the nature of the binding. In general, during our simulations $`\theta _{loc}^{}\theta _{mob}^{}`$ holds for most of the time.

For the calculation of the Te coverage $`\rho ^{\mathrm{Te}}`$ besides mobile Te atoms at regular lattice sites both localized and mobile Te atoms are taken into account as well, i.e. $`\rho ^{\mathrm{Te}}`$ includes $`\theta ^{}`$ . According to the experimentally observed range of Te-coveragescibert ; tatarenko2  the total reservoir occupation thus is restricted to values $`0\theta ^{}1/2`$. In our model we will exclude larger values by a simple cut-off: if $`\theta ^{}`$ has reached its maximum value all processes leading to an increase of $`\theta ^{}`$ are forbidden, i.e. their rate is set to zero.

We furthermore assume that an attractive pairwise interaction of Te stabilizes the reservoir. In a mean field fashion we represent the latter by an energy contribution $`Nϵ_b^{}\theta _{}^{}{}_{}{}^{2}`$ with parameter $`ϵ_b^{}<0`$.

In summary, the total energy of a given configuration can be written as

$`H`$ $`=`$ $`ϵ_cn_c+ϵ_b\left(n_b^{\mathrm{Cd}}+n_b^{\mathrm{Te}}+n_h^{\mathrm{Te}}\right)+ϵ_xn_x^{\mathrm{Cd}}+ϵ_dn_d^{\mathrm{Cd}}`$ (1)
$`+`$ $`ϵ_m^{\mathrm{Te}}n_m^{\mathrm{Te}}+ϵ_b^{\mathrm{Td}}n_b^{\mathrm{Td}}+ϵ_x^{\mathrm{Td}}n_x^{\mathrm{Td}}`$
$`+`$ $`ϵ^{}\theta ^{}N+ϵ_b^{}\theta _{}^{}{}_{}{}^{2}N`$
which depends on the number of Cd-Te bonds ($`n_c`$), bulk NN pairs of Cd ($`n_b^{\mathrm{Cd}}`$) and Te ($`n_b^{\mathrm{Te}}`$), NN pairs of one mobile and one immobile Te ($`n_h^{\mathrm{Te}}`$), NN and NNN pairs of mobile Cd at the surface ($`n_x^{\mathrm{Cd}}`$ and $`n_d^{\mathrm{Cd}}`$), and NN pairs of mobile Te not forming a dimer ($`n_m^{\mathrm{Te}}`$). The number of Te-dimers is denoted as $`n_b^{\mathrm{Td}}`$, whereas $`n_x^{\mathrm{Td}}`$ counts pairs of neighboring Te-dimers along the $`x`$-direction. The last line in Eq. (1) accounts for the contribution of the Te reservoir which is occupied by $`\theta ^{}N`$ atoms.

Throughout this paper we refer to the following parameter set:

$`ϵ_d=1,ϵ_x=1.95,ϵ_c=6.5,ϵ^{}=5,ϵ_b=0.8,`$
$`ϵ_b^{}=8,ϵ_m^{\mathrm{Te}}=0.6,ϵ_b^{\mathrm{Td}}=1,ϵ_x^{\mathrm{Td}}=0.4`$ (2)
Note that $`\left|ϵ_d\right|`$ only fixes the energy scale in the following. In addition we set $`k_B=\left|ϵ_d\right|=1`$ formally and, hence, measure the temperature $`T`$ on the same scale. Later we will comment on a potential quantitative comparison with experiments.

### III.2 Dynamics and growth

As we are interested in the non-equilibrium dynamics of the system under sublimation or growth conditions, specifying the energy of a given configuration is not sufficient. In fact, only the energy contribution of a few active layers which contain mobile particles will be relevant in the following. The kinetics of the system is governed by the energy barriers which have to be overcome in transitions between the different configurations.

Unfortunately, the available experimental data and first principles results are not sufficient for a systematic fit of our model parameters. In particular, reliable evaluations or estimates of, e.g., diffusion barriers are not available, apart from very few exceptions. Therefore, we do not aim at a precise quantitative description of the CdTe(001) surface at this stage.

Nevertheless, essential features can be deduced from experimental data and physical insight. Important qualitative features of the model turn out to be quite robust against variations of the parameters, as long as they comply with some essential conditions. In our simulations we will use a parameter set which is mainly based on our previous investigations of related models and their comparison with experimental findings.

In our previous studies of simpler models, appropriate choices of the parameter set enabled us to reproduce qualitatively the observed temperature dependent structure of Cd-terminated surfaces in step-flow sublimation. Furthermore, a semi-quantitative match of the re-ordering temperature, macroscopic sublimation rates etc. has been achieved, see Refs. aachen, ; ahrdiss, .

The microscopic processes considered in our Kinetic Monte Carlo (KMC) dynamics are deposition of Cd and Te adatoms, diffusion and desorption of mobile atoms, transitions of Te atoms between regular lattice sites and the reservoir, and desorption of reservoir atoms.

Deposition of Cd and Te atoms occurs with fluxes $`F^{\mathrm{Cd}}`$ and $`F^{\mathrm{Te}}`$, respectively. First, a lattice site is chosen randomly. Then, a search is performed within an incorporation radius $`r_{inc}=2`$ for the site with the lowest height where adsorption is possible.ahrdiss  If this is, for a Cd atom, a site with a neighboring mobile Cd atom in y-direction, adsorption is only possible if there is either a localized Te on top of the Cd neighbor, already, or a mobile Te attaches to the impinging Cd. This occurs with probability $`\theta _{mob}^{}`$. For a Te atom proper adsorption sites include, besides regular Te sites, mobile Cd atoms with no localized Te attached as well as empty Cd sites. In the latter two cases, the Te then becomes a mobile reservoir Te with probability $`12\theta ^{}`$, whereas with probability $`2\theta ^{}`$ deposition is rejected and the particle evaporates. As a consequence, the reservoir occupation is limited to values $`\theta ^{}1/2`$.

The rates for all thermally activated processes are of the Arrhenius form

$$r_i=\nu _i\alpha _i\mathrm{exp}\left(E_i/T\right),$$
(3)
where $`E_i`$ and $`\nu _i`$ denote the activation energy and attempt frequency for process $`i`$, respectively. The prefactor $`\alpha _i`$ accounts for rates which depend explicitly on the reservoir occupation (cf. Table 1). To calculate the activation energy $`E_i`$ we assume the energy of the transition state $`H_{t,i}`$ between the starting ($`s`$) and the final ($`f`$) configuration of the system to be equal to the greater of the energies $`H_{s,i}`$ and $`H_{f,i}`$ of both configurations plus an additional barrier:newman99  $`H_{t,i}=\mathrm{max}\{H_{s,i}+B_{s,i},H_{f,i}+B_{f,i}\}`$. The values of the barrier heights $`B_{s,i}`$ and $`B_{f,i}`$ in general are not the same but depend on the type of process $`i`$. For the activation energies $`E_i=H_{t,i}H_{s,i}`$ one obtains

$$E_i=\mathrm{max}\{B_{s,i},B_{f,i}+\mathrm{\Delta }H_i\}$$
(4)
for diffusion of mobile atoms and transitions of Te atoms between regular lattice sites and reservoir states with the energy change $`\mathrm{\Delta }H_i=H_{f,i}H_{s,i}`$ given by Eq. (1). Diffusion includes moves to available NN or NNN sites in the same layer, as well as hops across terrace edges, cf. Figure 1. For all diffusion steps we set $`B_{s,i}=B_{f,i}=B_o`$ in Eq. (4) with barrier $`B_o=7`$, similar to Refs. aachen, ; ahrdiss, . On the contrary, the assumption that Te atoms in the reservoir are less strongly bound implies $`B_{s,i}=B_o`$ and $`B_{f,i}=B^{}<B_o`$ for the transition of a Te atom into the reservoir and, conversely, $`B_{s,i}=B^{}`$ and $`B_{f,i}=B_o`$ for the opposite process (cf. Table 1).

Desorption of atoms from regular lattice sites or mobile Te requires an activation energy $`E_i=\mathrm{\Delta }H_i`$. We consider only desorption of single atoms, neglecting the fact that Te probably evaporates from the surface in the form of Te2 molecules, in reality.

For all diffusion processes and the desorption of mobile atoms at regular surface sites we assume an attempt frequency $`\nu _i=\nu _o=10^{12}s^1`$, whereas a lower frequency $`\nu _i=\nu ^{}`$ is used for the transition of Te atoms to regular Te states and the evaporation of mobile Te.

This choice is based on the reasonable assumption that Te atoms reside in potential energy minima which are much shallower than those for Te atoms at regular surface sites.zangwill  According to transition state theory the attempt frequency $`\nu ^{}`$ which corresponds to the harmonic oscillation frequency at the local energy minimum will become smaller.tst  A similar argument was used by Itoh itoh01  in order to motivate an attempt frequency of weakly bound arsenic as low as $`10^8s^1`$. In our simulations, we will also set $`\nu ^{}=10^8s^1`$. This choice is further motivated by comparison with experimental observations as discussed in the following section.

The formation or breaking of dimers from NN pairs of mobile Te atoms is taken into account implicitly. A diffusion hop of a Cd atom onto a Te dimer, for instance, requires the breaking of the dimer which has to be considered in the energy balance $`\mathrm{\Delta }H_i`$ of the process. The deposition of a Cd atom will always break all underlying Te dimers involved.

On the other hand, dimer formation is assumed to occur instantaneously, whenever it is possible due to the resulting energy gain. A similar simplification was used in Ref. barnett,  for the description of Si-dimerization on Si(001) surfaces. Hence, if any event results in two un-dimerized Te atoms which are in the correct configuration for forming a dimer pair, a dimer bond will be created immediately. For diffusion and desorption events, dimer formation will be accounted for in $`\mathrm{\Delta }H_i`$.

The existence of a weakly bound state for Te atoms is essential to facilitate growth in our model. Consider a Cd-terminated surface with $`c(2\times 2)`$ reconstruction and maximum coverage $`\rho ^{\mathrm{Cd}}=1/2`$: clearly, deposition of Te at a regular lattice site would violate the SOS condition. On the other hand, the incorporation of additional Cd is impossible due to the infinite repulsion of NN neighbors in $`y`$-direction.

It is plausible to assume, however, that an additional Te bound to a Cd would neutralize the repulsion as it provides three negatively charged dangling bonds. In our model, mobile Cd can hop to empty sites neighboring another mobile Cd in $`y`$-direction, provided there is either already a Te localized on top of the neighboring Cd, or a mobile Te attaches to the diffusing particle. This Te then would become localized itself. The latter scenario is assumed to happen with a probability $`\theta _{mob}^{}`$ which is reflected in the prefactor $`\alpha _i`$ of the corresponding rate (cf. Table 1).

Thus, we assume that mobile Te which are not attached to Cd atoms are available instantaneously everywhere at the surface. Note that after the diffusion step the reservoir occupation $`\theta ^{}=\theta _{loc}^{}+\theta _{mob}^{}`$ remains the same but the ratio $`\theta _{loc}^{}/\theta _{mob}^{}`$ has increased. If at any point of the simulation a localized Te is no longer needed because for example the underlying Cd atom has evaporated, it will become a mobile Te again. Through the mediation of Te atoms, Cd atoms at the surface can coalesce into islands in spite of the NN repulsion, provided the reservoir of Te is not empty.

In our simulations we consider also transitions of Te between regular lattice sites and reservoir states as well as the desorption of mobile Te into the vacuum. In general, the rates of these processes depend on the reservoir occupation, i.e. $`\alpha _i1`$ (cf. Table 1). For instance a mobile Te at the surface may become a mobile Te provided the reservoir is not fully occupied. This happens with a rate according to Eq. (3) with $`\alpha _i=\left(12\theta ^{}\right)`$ which is zero if the reservoir density has reached its maximum value $`1/2`$. Conversely, the transition of a mobile Te to a regular Te state is proportional to the density $`\theta _{mob}^{}`$ of mobile Te atoms. Localized Te atoms which are on top of Cd atoms in order to neutralize the infinite repulsion may become bound at one of the two neighboring regular Te sites, if available, thus reducing the density $`\theta _{loc}^{}`$.

The results presented in the following sections were obtained in simulations with energy parameters according to Eq. (III.1) and the following barriers and frequencies:

$$\nu _o=10^{12}s^1,\nu ^{}=10^8s^1,B_o=7,B^{}=2$$
(5)
## IV Simulations and results

### IV.1 The reservoir of excess Te

Before presenting results of the ALE simulations we discuss the properties of the Te reservoir in greater detail. In our model, the occupation $`\theta ^{}`$ results from the complex interplay of deposition, evaporation and transitions between Te and Te states. However, the basic temperature dependence of $`\theta ^{}`$ can be understood from considering a flat surface in the presence of a steady flux of Te atoms.

After the Te coverage has reached the value $`\rho ^{\mathrm{Te}}=1`$, only deposition directly into the reservoir is possible. If we neglect the desorption of regular Te from the surface and transitions of the type Te $``$ Te, the temporal evolution of $`\theta ^{}`$ is governed by the balance of direct deposition and evaporation of Te. The rate of the latter is a function of $`\theta ^{}`$ itself and one obtains the following simplifying ODE:

$$\frac{d\theta ^{}}{dt}F^{\mathrm{Te}}\mathrm{\Theta }\left[1/2\theta ^{}\right]\theta ^{}\nu ^{}\mathrm{exp}\left[\frac{1}{T}\left(\left|ϵ^{}\right|+\theta ^{}\left|ϵ_b^{}\right|\right)\right]$$
(6)
Whereas the Te flux increases $`\theta ^{}`$, desorption of Te occurs with a frequency proportional to $`\theta ^{}`$ multiplied with the Arrhenius rate for single particle evaporation from the reservoir. Clearly, the filling of the reservoir must saturate at $`\theta ^{}1/2`$. However, the precise nature of the saturation process is not essential for the following arguments and we choose a simple cut-off represented by the Heaviside-function $`\mathrm{\Theta }[x]`$.

Given $`ϵ^{},ϵ_b^{}`$ and $`F^{\mathrm{Te}}`$, the right hand side in Eq. (6) is positive for all $`\theta ^{}[0,1/2]`$ if the temperature is low enough. Hence, the reservoir will be filled in this case. Above a characteristic temperature $`T_{res}`$ the r.h.s. displays two zeros at fairly small values of $`\theta ^{}`$. As a consequence, the differential equation has an attractive fixed point with a stationary value $`\theta ^{}0`$ for $`T>T_{res}`$.

Figure 2 (left) shows $`T_{res}`$ as function of $`F^{\mathrm{Te}}`$ for the choice of parameters $`ϵ^{}=5`$ and $`ϵ_b^{}=8`$. There is a rapid increase of $`T_{res}`$ for very low values of $`F^{\mathrm{Te}}`$ which becomes weaker and weaker for higher fluxes. With $`F=5ML/s`$ we obtain $`T_{res}0.39`$, as an example, and at $`T=0.40`$ the stationary occupation is as low as $`\theta ^{}0.02`$. Hence, starting from $`\theta ^{}=0`$, the reservoir cannot be filled significantly for $`T>T_{res}`$. The behavior predicted from Eq. (6) is in excellent agreement with simulations of the situation in the full model. We would like to mention that in Ref. faschinger,  a simple evaporation model of ALE growth was proposed which produces a flux-dependence of the characteristic temperature which is similar to the one in Figure 2 (left).

The simplified Eq. (6) is also useful in investigating the initial sublimation of Te from a surface with maximum coverage $`\rho ^{\mathrm{Te}}=3/2`$, i.e. $`\theta ^{}=1/2`$. Figure 2 (right) displays the evolution of $`\theta ^{}`$ for $`ϵ^{}=5,ϵ_b^{}=8`$ and three different temperatures as obtained from the numerical integration of (6) with $`F^{\mathrm{Te}}=0`$. The results from simulations of the corresponding situation in the full kinetic model are again in excellent agreement with the simplified description, which neglects evaporation of regular Te and other processes.

One obtains a well-defined characteristic time $`\tau ^{}`$ for the decrease of $`\theta ^{}`$ from $`0.5`$ to, say, $`\theta _\tau ^{}=0.05`$. The precise value of $`\theta _\tau ^{}`$ is of little relevance as $`\theta ^{}(t)`$ decreases to zero very rapidly, cf. Figure 2 (right panel). The temperature dependence of $`\tau ^{}`$ is very well described by an Arrhenius law of the form $`\tau ^{}=\tau _oe^{E_a/T}`$ with macroscopic activation energy $`E_a8.5`$ and prefactor $`\tau _o3.8\times 10^9s`$. In Refs. tatarenko, ; cibert,  experimental data is discussed for the initial sublimation of Te starting from coverage $`3/2`$. There, the macroscopic activation energy is reported to be $`E_a0.96eV`$, providing a potential reference point for attempting a quantitative fit of the model.

Note that the prefactor $`\tau _o10^9s`$ obtained from Fig. 2 (right) is not too different from the time constant $`1/\nu ^{}=10^8s`$ assumed for Te desorption. This should be expected, as the latter is the limiting microscopic process of the scenario. In a naive attempt to extract the prefactor from an analogous Arrhenius plot of experimental data tatarenko  we obtain a value of about the same order, i.e. $`\tau _o10^9s`$. On the contrary, the value from data for Cd desorption tatarenko  corresponds to frequencies which are larger by orders of magnitude. Of course, a reliable determination of prefactors from experimental data for relatively small temperature ranges is questionable if not impossible (and was not attempted by the authors of Ref. tatarenko, ). Nevertheless one might interpret the very rough comparison as another qualitative justification for choosing $`\nu ^{}`$ much smaller than other frequencies in the model.

### IV.2 ALE growth

In order to gain a qualitative understanding of the ALE growth scenario we perform simulations with the parameter set specified in Eqs. (III.1,5). All simulation results presented here were obtained from systems with $`64\times 64`$ Cd and Te atoms in a bulk layer, averaged over $`5`$ independent runs.

We model a situation in which Cd and Te are deposited in alternate pulses of length $`0.9s`$, each with a constant flux $`F^{\mathrm{Cd}}=F^{\mathrm{Te}}=5ML/s`$. The pulses are separated by a dead time of $`0.1s`$, hence, the duration of one complete ALE cycle is $`2s`$ in total. Note that the characteristic time $`\tau ^{}`$ for reservoir sublimation is large compared with the dead time interval for all relevant temperatures. In the results discussed here, the initial surface was perfectly flat and Te-terminated with coverage $`\rho ^{\mathrm{Te}}=1`$, $`\theta ^{}=0`$. Deposition always started with a Cd pulse.

Illustrating surface snapshots of growing surfaces are shown in Figures 3 and 4 for the two different growth regimes discussed in the following. Additional illustrations and mpeg-movies are available from our web pages.web

At fairly high temperatures, growth proceeds according to the following scenario in our simulations:

* The first Cd pulse adds half a layer of Cd to the system. The infinite repulsion of NN in $`y`$-direction prevents coverages $`\rho ^{\mathrm{Cd}}>1/2`$. The predominant arrangement of Cd atoms is in a local $`c(2\times 2)`$ ordering, i.e. a checkerboard pattern. The contribution of the $`(2\times 1)`$ row structure as a thermal excitation increases with higher T.
* After the first Te pulse and the dead time, half a layer of Te has been incorporated. Diffusion of particles leads to the formation of islands upon the Te-terminated surface. Note that the existence of the Te-reservoir state is crucial for this process, as discussed above. Surface Te atoms form dimers, preferentially, and these arrange in rows along the $`x`$-direction, yielding a predominant $`(2\times 1)`$ symmetry of the surface. Although incorporation of Te occurs to a considerable degree through the reservoir state Te, at the end of the dead time one finds $`\theta ^{}0`$, as discussed above.
* The second Cd pulse places mobile atoms on top of existing islands and in between. The pre-dominant local ordering is the $`(2\times 1)`$ pattern, i.e. rows of Cd alternating with empty rows. This is due to the influence of the island edges on the reconstruction, see Refs. aachen, ; surflett, ; ahrdiss,  for a discussion of an analogous effect in layer by layer sublimation. At the end of the dead time, the islands from stage II are still present but now they are Cd-terminated.
* During the second Te pulse, atoms are deposited onto the Cd-terminated islands and in between. The particles rearrange by means of diffusion, leveling off the islands and forming a flat Te-terminated surface. At the end of the dead time, the initial state of the system is restored with one complete layer of CdTe added to the system. In the simulations the filling of gaps is incomplete to a degree which depends on the temperature and the length of the dead time interval.

Hence, at temperatures $`T>T_{res}`$, our model reproduces ALE growth in the same double cycle sequence which was hypothesized and discussed for CdTe and ZnTe in the literature, see Refs. cibert, ; tatarenko2,  and references therein. Figure 5 shows the evolution of the film height with time in the ALE process at $`T=0.44`$, which is well above the characteristic $`T_{res}0.39`$ obtained for $`F^{\mathrm{Te}}=5ML/s`$.

The picture changes qualitatively for temperatures $`T<T_{res}`$. Again, a flat Te-terminated surface with coverage $`\rho ^{\mathrm{Te}}=1`$, $`\theta ^{}=0`$ was prepared, and the deposition started with a Cd pulse. The surface structures obtained after the first Cd and Te pulse are very much alike the ones of stages I and II of the high temperature scenario, since initially the Te reservoir was empty. Note however, that due to a decreased mobility of the atoms the islands on the Te terminated surface are less compact. Since $`T<T_{res}`$, at the end of the Te pulse the surface is not only Te terminated but also the Te reservoir is filled, i.e. $`\theta ^{}1/2`$.

Thus, from this point on the ALE growth scenario proceeds along the following lines in our model:

* The Cd pulse places atoms onto the islands and in between which is analogous to stage III of the high temperature scenario. But now Cd particles at the surface can coalesce, because Te from the reservoir are available to neutralize the NN repulsion of mobile Cd in $`y`$-direction. Therefore the gaps between the islands are leveled out by half a monolayer of Cd atoms and the reservoir which provides also approx. $`1/2`$ monolayer of Te. The surface is then covered by another half layer of Cd atoms arranging preferentially in a $`c(2\times 2)`$ vacancy structure. Thus, the adsorption of one complete monolayer of Cd and a half monolayer of Te which came from the reservoir leads to an increase of the surface height by $`3/4`$ of a complete CdTe layer. At the end of step A) the reservoir is empty.
* The following Te pulse leads to a Te terminated surface with islands, which is identical to the one described above. Again, due to low temperature and short dead time, also the Te reservoir is filled. Hence, another half monolayer of Te is adsorbed and in total one complete layer of CdTe added in only one ALE cycle. At the end of step B) the reservoir has reached its maximum occupation.

As an example, Figure 5 shows the evolution of the mean surface height $`h`$ with time at $`T=0.36`$, i.e. below $`T_{res}`$ for $`F^{\mathrm{Te}}=5ML/s`$. Note that for the first Cd and Te pulse, the curve is essentially identical with that of the high $`T`$ regime, because in both simulations the reservoir was empty, initially. Hence the upper curve of Figure 5 corresponds to the sequence of stages A/B only for times $`t>2s`$. The mean height as displayed in Figure 5 increases by $`3/4`$ of a complete CdTe layer during the Cd pulse A) because $`h`$ takes into account only particles at SOS lattice sites, disregarding the Te reservoir.

From simulations for different temperatures we have obtained the average growth rate per ALE cycle. Figure 6(b) shows that, with increasing $`T`$, the growth rate drops very rapidly in the vicinity of the characteristic $`T_{res}`$ as obtained from the considerations of section IV.1. The comparison with Figure 6(a) shows qualitative agreement with experimental results for an analogous growth scenario.faschinger  Very similar data has been published by Hartmann et al.hartmann  and is also discussed in Ref. cibert, . Note that in both experiments growth rates above 1 layer of CdTe per ALE cycle were observed at very low temperatures. This is presumably due to the formation of Cd crystallites at the surface and cannot be reproduced in our model.tatarenko2

The temperature driven break-down of the growth rate is due to the fact that the reservoir cannot be filled significantly for $`T>T_{res}`$. Simulation results for various values of the flux strongly support this assumption as can be seen from Figure 2 (left): the symbols which indicate the temperature where the growth rate drops in the simulations lie almost exactly on the curve that marks $`T_{res}`$.

In principle, the re-evaporation of Te during the dead time interval could be an alternative cause for the reduced growth rate. Note, however, that even for a flux as high as $`F^{\mathrm{Te}}=10ML/s`$ with $`T_{res}0.41`$, one finds that the characteristic time for reservoir sublimation is $`\tau ^{}4s`$, cf. Figure 2 (right panel). Hence, it appears safe to say that during dead time intervals shorter than, say, $`1s`$, $`\theta ^{}`$ will remain close to the maximum value $`1/2`$, provided all other conditions allow for the filling of the reservoir during a Te pulse.

As a consequence, we expect that the characteristic value of $`T`$ which marks the drop of the growth rate in Figure 6 should be essentially independent of the pulse duration or dead time interval within a wide range of reasonable values.

Comparing the experimental results by Hartmann et al.hartmann  and Faschinger and Sitterfaschinger  for the ALE growth rate one finds that for the latter the temperature where the growth rate drops to approx. $`1/2`$ layers of CdTe per cycle is shifted to higher temperatures by about $`40K`$. Within the framework of our model this can be explained if one takes into account the different experimental conditions. Whereas Hartmann et al. combined particle fluxes of about $`0.5ML/s`$ with $`8s`$ pulse time and $`1s`$ dead times, the flux used by Faschinger and Sitter was $`2ML/s`$ with a pulse duration of $`0.8s`$ and a dead time of $`0.2s`$. As discussed above, the difference of the dead times can be neglected, but for the given fluxes one reads off from Fig. 2 (left) a significant increase of $`T_{res}`$ by about $`9\%`$. As demonstrated, $`T_{res}`$ practically coincides with the temperature where the growth rate drops which consequently is expected to increase by the same amount. The difference of the absolute temperatures observed in Refs. hartmann,  and faschinger,  is about $`8\%`$ (of the lower one) which is indeed comparable with our results.

For very high $`T`$, sublimation of the crystal dominates over the incoming flux and growth becomes impossible. The ALE growth rate vanishes at a temperature which is essentially independent of the properties of Te atoms.

## V Summary and Outlook

In summary, we have presented a detailed lattice gas model of II-VI semiconductor compounds. In contrast to many simpler growth models, it takes into account the correct zinc-blende lattice as well as essential, material specific effects which have been discussed in terms of the CdTe(001) surface, as an example. The model reproduces, for instance, the reconstruction of Cd-terminated surfaces by means of vacancy structures and its dependency on temperature and flux. This includes the temperature driven transition from the $`c(2\times 2)`$ reconstructed surface to a disordered state with local $`(2\times 1)`$ ordering. This re-arrangement of the surface is observed under step flow sublimation, as demonstrated earlier in Refs. aachen, ; surflett, , and a semi-quantitative agreement with experimental data can be achieved.

Here, in addition, the dimerization of Te at the surface is introduced to the model. Thus, the $`(2\times 1)`$ reconstruction of Te-terminated surfaces is reproduced. We demonstrate furthermore that the presence of weakly bound Te atoms is necessary to facilitate growth by, e.g., atomic layer epitaxy. Although the precise nature of such a weakly bound Te state at the surface is unknown, its existence is clearly evident from various experimental observations.cibert ; tatarenko2

By means of Kinetic Monte Carlo simulations of our model we show that the reservoir of Te also explains the experimentally observed temperature dependence of the ALE growth rates. Only below a characteristic temperature a significant amount of Te is present and facilitates growth rates of approximately $`1`$ layer of CdTe per ALE cycle. As in the corresponding experiments, we observe a sudden transition from the self-limited growth rate of about $`1/2`$ layer per ALE cycle at high $`T`$ to the low temperature regime with a rate close to $`1`$ complete layer per cycle.

Our analysis shows that the characteristic transition temperature should depend only weakly on the model (or experimental) parameters within a wide range of reasonable choices. The key dependence is on the particle flux during Cd and Te pulses, respectively, and our model reproduces its effect on the transition qualitatively correct.

As an attempt towards a more quantitative comparison with experiments, one might identify the temperature $`T_{res}0.365`$ for a flux of $`2ML/s`$ in our model (cf. Fig 2) with the experimental value $`T290^oC`$ by Faschinger et al.faschinger  Thus, the energy scale $`ϵ_d`$ may be fixed, and we obtain, for example, a value of $`E_a1.13eV`$ for the macroscopic activation energy for desorption of Te atoms which roughly agrees with the experimentally determined $`0.96eV`$.tatarenko ; cibert  Furthermore, the value of $`T_{res}`$ for our simulations with $`5ML/s`$ translates into $`T_{res}329^oC`$ which agrees well with results of the simple evaporation model from Ref. faschinger, .

We thus believe, that, in principle, our model offers the potential for a quantitative match with experimental data. In simpler versions of the model it was possible to reproduce the temperature dependence of the surface reconstruction as well as macroscopic sublimation rates on a semiquantitative level.aachen  The larger number of parameters in the current extension, however, requires further input from experiment. Moreover, we hope to obtain reliable estimates of microscopic activation barriers from first principle quantum chemical calculations as are already available for III-VI semiconductors.grosseprb ; grosseprl  Such calculations will be extremely useful and should also shed light on the nature of the weakly bound Te states on the surface.

###### Acknowledgements.

We would like to thank S. Tatarenko for useful discussions. This project was supported by the Deutsche Forschungsgemeinschaft through Sonderforschungsbereich 410.


format: text

text:
# References

Jets of modules in noncommutative geometry

G. Sardanashvily

Department of Theoretical Physics, Moscow State University, 117234 Moscow, Russia

E-mail: sard@grav.phys.msu.su

URL: http://webcenter.ru/$``$sardan/

Abstract. Jets of modules over a commutative ring are well known to make up the representative objects of linear differential operators on these modules. In noncommutative geometry, jets of modules provide the representative objects only of a certain class of first order differential operators. As a consequence, a generalization of the standard Lagrangian formalism on smooth manifolds to noncommutative spaces is problematic.

Let $`𝒦`$ be a commutative ring and $`𝒜`$ a commutative $`𝒦`$-ring. Let $`P`$ and $`Q`$ be $`𝒜`$-modules and $`𝒥^r(P)`$ the module of $`r`$-order jets of $`P`$. There is a $`𝒜`$-module isomorphism

$$\mathrm{Hom}_𝒜(𝒥^r(P),Q)=\mathrm{Diff}_r(P,Q),$$
(1)
where $`\mathrm{Diff}_r(P,Q)`$ is the module of $`r`$-order $`Q`$-valued linear differential operators on an $`𝒜`$-module $`P`$ . It follows that the jet module $`𝒥^r(P)`$ is the representative object of the functor $`Q\mathrm{Diff}_r(P,Q)`$.

In particular, let $`𝒜=C^{\mathrm{}}(X)`$ be the $`R`$-ring of smooth real functions on a smooth manifold $`X`$. If $`P`$ is a projective $`C^{\mathrm{}}(X)`$-module of finite rank, it is isomorphic to the module $`Y(X)`$ of global sections of some vector bundle $`YX`$. In this case, $`𝒥^r(P)`$ is the module of global sections of the $`r`$-order jet bundle $`J^rYX`$ of $`YX`$. Furthermore, let $`YX`$ be an arbitrary smooth fibre bundle. An $`r`$-order (not necessarily linear) differential operator on $`YX`$ is defined as a bundle morphism of $`J^rYX`$ to some smooth fibre bundle over $`X`$. As a consequence, Lagrangian formalism on smooth manifolds is conventionally phrased in the jet terms.

If $`𝒜`$ is not commutative, there are different definitions of (linear) differential operators on modules over $`𝒜`$ . Though the notion of a jet can be extended to a module $`P`$ over a noncommutative ring $`𝒜`$, only the first order jet module $`𝒥^1(P)`$ is the representative object of a certain class of first order differential operators on $`P`$.

Let $`𝒜`$ be a commutative $`𝒦`$-ring. Let $`P`$ and $`Q`$ be $`𝒜`$-modules (central bimodules). The $`𝒦`$-module $`\mathrm{Hom}_𝒦(P,Q)`$ of $`𝒦`$-linear homomorphisms $`\mathrm{\Phi }:PQ`$ is endowed with the two different $`𝒜`$-module structures

$$(a\mathrm{\Phi })(p):=a\mathrm{\Phi }(p),(\mathrm{\Phi }a)(p):=\mathrm{\Phi }(ap),a𝒜,pP.$$
(2)
We will refer to the second one as the $`𝒜^{}`$-module structure. Let us put

$$\delta _a\mathrm{\Phi }:=a\mathrm{\Phi }\mathrm{\Phi }a,a𝒜.$$
(3)
The following definitions are equivalent.

Definition 1. An element $`\mathrm{\Delta }\mathrm{Hom}_𝒦(P,Q)`$ is an $`r`$-order $`Q`$-valued differential operator on $`P`$ if $`\delta _{a_0}\mathrm{}\delta _{a_r}\mathrm{\Delta }=0`$ for any tuple of $`r+1`$ elements $`a_0,\mathrm{},a_r`$ of $`𝒜`$.

Definition 2. An element $`\mathrm{\Delta }\mathrm{Hom}_𝒦(P,Q)`$ is a zero order differential operator if $`\delta _a\mathrm{\Delta }=0`$ for all $`a𝒜`$, and $`\mathrm{\Delta }`$ is a differential operator of order $`r>0`$ if $`\delta _a\mathrm{\Delta }`$ for all $`a𝒜`$ is an $`(r1)`$-order differential operator.

Given an $`𝒜`$-module $`P`$, let us consider the tensor product $`𝒜_𝒦P`$ of $`𝒦`$-modules $`𝒜`$ and $`P`$. We put

$$\delta ^b(ap):=(ba)pa(bp),pP,a,b𝒜.$$
(4)
Let us denote by $`\mu ^{k+1}`$ the $`𝒜`$-submodule of $`𝒜_𝒦P`$ generated by elements of the type

$`\delta ^{b_0}\mathrm{}\delta ^{b_k}(\mathrm{𝟏}p).`$
Definition 3. The $`k`$-order jet module $`𝒥^k(P)`$ of a module $`P`$ is defined as the quotient of the $`𝒦`$-module $`𝒜_𝒦P`$ by $`\mu ^{k+1}`$. We denote its elements $`a_kp`$.

The $`𝒦`$-module $`𝒥^k(P)`$ is endowed with the $`𝒜`$\- and $`𝒜^{}`$-module structures

$$b(a_kp):=ba_kp,b(a_kp):=a_k(bp).$$
(5)
There exists the module morphism

$$J^k:Pp\mathrm{𝟏}_kp𝒥^k(P)$$
(6)
of the $`𝒜`$-module $`P`$ to the $`𝒜^{}`$-module $`𝒥^k(P)`$ such that $`𝒥^k(P)`$, seen as an $`𝒜`$-module, is generated by elements $`J^kp`$, $`pP`$.

Theorem 4. Any $`Q`$-valued differential operator $`\mathrm{\Delta }`$ of order $`k`$ on an $`𝒜`$-module $`P`$ factorizes uniquely

$`\mathrm{\Delta }:P\stackrel{J^k}{}𝒥^k(P)Q`$
through the morphism $`J^k`$ (6) and some $`𝒜`$-module homomorphism $`𝔣^\mathrm{\Delta }:𝒥^k(P)Q`$ .

Let us denote $`J:Pp\mathrm{𝟏}p𝒜P`$. The proof is based on the fact that

$$\delta _{b_0}\mathrm{}\delta _{b_k}(𝔣J)(p)=𝔣(\delta ^{b_0}\mathrm{}\delta ^{b_k}(\mathrm{𝟏}p))$$
(7)
for any $`𝔣\mathrm{Hom}_𝒜(𝒜P,Q)`$. The correspondence $`\mathrm{\Delta }𝔣^\mathrm{\Delta }`$ yields the above mentioned isomorphism (1).

Let now a $`𝒦`$-ring $`𝒜`$ need not be commutative. Let $`P`$ and $`Q`$ be two-sided $`𝒜`$-modules. Two-sided $`𝒜`$-modules throughout are assumed to be central bimodules over the center $`𝒵_𝒜`$ of $`𝒜`$. The $`𝒦`$-module $`\mathrm{Hom}_𝒦(P,Q)`$ can be provided with the left $`𝒜`$\- and right $`𝒜^{}`$-module structures (2) and the similar right and left structures

$$(\mathrm{\Phi }a)(p):=\mathrm{\Phi }(p)a,(a\mathrm{\Phi })(p):=\mathrm{\Phi }(pa),a𝒜,pP.$$
(8)
For the sake of convenience, we will refer to $`𝒜𝒜^{}`$ structures (2) and (8) as the left and right $`𝒜𝒜^{}`$ structures, respectively. Let us put

$$\overline{\delta }_a\mathrm{\Phi }:=\mathrm{\Phi }aa\mathrm{\Phi },a𝒜,\mathrm{\Phi }\mathrm{Hom}_𝒦(P,Q).$$
(9)
It is readily observed that $`\delta _a\overline{\delta }_b=\overline{\delta }_b\delta _a`$ for all $`a,b𝒜`$.

We follow the notion of a differential operator in . Let $`P`$ and $`Q`$ be regarded as left $`𝒜`$-modules. Let us consider the $`𝒦`$-module $`\mathrm{Hom}_𝒦(P,Q)`$ provided with the left $`𝒜𝒜^{}`$ module structure (2). We denote $`𝒵_0`$ its center, i.e., $`\delta _a\mathrm{\Phi }=0`$ for all $`\mathrm{\Phi }𝒵_0`$ and $`a𝒜`$. Let $`_0=\overline{𝒵}_0`$ be the $`𝒜𝒜^{}`$ submodule of $`\mathrm{Hom}_𝒦(P,Q)`$ generated by $`𝒵_0`$. Let us consider: (i) the quotient $`\mathrm{Hom}_𝒦(P,Q)/_0`$, (ii) its center $`𝒵_1`$, (iii) the $`𝒜𝒜^{}`$ submodule $`\overline{𝒵}_1`$ of $`\mathrm{Hom}_𝒦(P,Q)/_0`$ generated by $`𝒵_1`$, and (iv) the $`𝒜𝒜^{}`$ submodule $`_1`$ of $`\mathrm{Hom}_𝒦(P,Q)`$ given by the relation $`_1/_0=\overline{𝒵}_1`$. Then we define the $`𝒜𝒜^{}`$ submodules $`_r`$, $`r=2,\mathrm{}`$, of $`\mathrm{Hom}_𝒦(P,Q)`$ by induction as follows:

$$_r/_{r1}=\overline{𝒵}_r,$$
(10)
where $`\overline{𝒵}_r`$ is the $`𝒜𝒜^{}`$ module generated by the center $`𝒵_r`$ of the quotient $`\mathrm{Hom}_𝒦(P,Q)/_{r1}`$.

Definition 5. Elements of the submodule $`_r`$ of $`\mathrm{Hom}_𝒦(P,Q)`$ are said to be the left $`r`$-order $`Q`$-valued differential operators on a two-sided $`𝒜`$-module $`P`$.

An equivalent definition is the following

Definition 6. An element of $`\mathrm{\Delta }\mathrm{Hom}_𝒦(P,Q)`$ is a left $`r`$-order $`Q`$-valued differential operator on an $`𝒜`$-module $`P`$ if it is a finite sum

$$\mathrm{\Delta }(p)=b_i\mathrm{\Phi }^i(p)+\mathrm{\Delta }_{r1}(p),b_i𝒜,$$
(11)
where $`\mathrm{\Delta }_{r1}`$ and $`\delta _a\mathrm{\Phi }^i`$ for all $`a𝒜`$ are left $`(r1)`$-order differential operators if $`r>0`$, and they vanish if $`r=0`$.

The set of these operators is provided with the left and right $`𝒜𝒜^{}`$ module structures.

If $`𝒜`$ is a commutative ring, Definitions 1 and 1 come to Definition 1.

If $`P=Q=𝒜`$, derivations of $`𝒜`$ and their compositions are differential operators on $`𝒜`$ in accordance with Definition 1 .

By analogy with Definitions 1 and 1, one can define right differential operators on two-sided $`𝒜`$-modules as follows.

Definition 7. Let $`P`$ and $`Q`$ be seen as right $`𝒜`$-modules over a noncommutative $`𝒦`$-ring $`𝒜`$. An element $`\mathrm{\Delta }\mathrm{Hom}_𝒦(P,Q)`$ is said to be a right zero order $`Q`$-valued differential operator on $`P`$ if it is a finite sum $`\mathrm{\Delta }=\mathrm{\Phi }^ib_i`$, $`b_i𝒜`$, where $`\overline{\delta }_a\mathrm{\Phi }^i=0`$ for all $`a𝒜`$. An element $`\mathrm{\Delta }\mathrm{Hom}_𝒦(P,Q)`$ is called a right differential operator of order $`r>0`$ on $`P`$ if

$$\mathrm{\Delta }(p)=\mathrm{\Phi }^i(p)b_i+\mathrm{\Delta }_{r1}(p),b_i𝒜,$$
(12)
where $`\mathrm{\Delta }_{r1}`$ and $`\overline{\delta }_a\mathrm{\Phi }^i`$ for all $`a𝒜`$ are right $`(r1)`$-order differential operators.

Definition 1 and Definition 1 of differential operators on two-sided $`𝒜`$-modules are not equivalent, but one can combine them as follows.

Definition 8. Let $`P`$ and $`Q`$ be two-sided modules over a noncommutative $`𝒦`$-ring $`𝒜`$. An element $`\mathrm{\Delta }\mathrm{Hom}_𝒦(P,Q)`$ is a two-sided zero order $`Q`$-valued differential operator on $`P`$ if it is either a left or right zero order differential operator. An element $`\mathrm{\Delta }\mathrm{Hom}_𝒦(P,Q)`$ is said to be a two-sided differential operator of order $`r>0`$ on $`P`$ if it is brought both into the form $`\mathrm{\Delta }=b_i\mathrm{\Phi }^i+\mathrm{\Delta }_{r1}`$, $`b_i𝒜`$, and $`\mathrm{\Delta }=\overline{\mathrm{\Phi }}^i\overline{b}_i+\overline{\mathrm{\Delta }}_{r1}`$, $`\overline{b}_i𝒜`$, where $`\mathrm{\Delta }_{r1}`$, $`\overline{\mathrm{\Delta }}_{r1}`$ and $`\delta _a\mathrm{\Phi }^i`$, $`\overline{\delta }_a\overline{\mathrm{\Phi }}^i`$ for all $`a𝒜`$ are two-sided $`(r1)`$-order differential operators.

One can think of this definition as a generalization of the definition of first order differential operators in noncommutative geometry in .

It is readily observed that two-sided differential operators in Definition 1 need not be left or right differential operators, and vice versa. At the same time, $`𝒜`$-valued derivations of a $`𝒦`$-ring $`𝒜`$ and their compositions obey Definition 1.

Turn now to the notion of jets in noncommutative geometry. One can follow Definition 1 in order to define the left jets of a two-sided module $`P`$ over a noncommutative ring $`𝒜`$. However, the relation (7) does not hold if $`k>0`$. Therefore, no module of left jets is the representative object of left differential operators. The notion of a right jet of $`P`$ meets the similar problem.

Given an $`𝒜`$-module $`P`$, let us consider the tensor product $`𝒜_𝒦P_𝒦𝒜`$ of $`𝒦`$-modules $`𝒜`$ and $`P`$. We put

$`\delta ^b(apc):=(ba)pca(bp)c,`$
$`\overline{\delta }^b(apc):=ap(cb)a(pb)c,pP,a,b,c𝒜.`$
Let us denote by $`\mu ^1`$ the two-sided $`𝒜`$-submodule of $`𝒜_𝒦P_𝒦𝒜`$ generated by elements of the type $`\overline{\delta }^c\delta ^b(\mathrm{𝟏}p\mathrm{𝟏})`$. One can define the first order two-sided jet module $`𝒥^1(P)=𝒜_𝒦P_𝒦𝒜/\mu ^1`$ of $`P`$. Let us denote

$`J:Pp\mathrm{𝟏}p\mathrm{𝟏}𝒜P𝒜.`$
Then the equality

$`\delta _{b_0}\mathrm{}\delta _{b_k}(𝔣J)(p)=𝔣(\delta ^{b_0}\mathrm{}\delta ^{b_k}(\mathrm{𝟏}p\mathrm{𝟏}))`$
holds for any $`𝔣\mathrm{Hom}_{𝒜𝒜}(𝒜P,Q)`$. It follows that $`𝒥^1(P)`$ is the representative object of the functor $`Q\overline{\mathrm{Diff}}_1(P,Q)`$ where the $`𝒦`$-module $`\overline{\mathrm{Diff}}_1(P,Q)`$ consists of two-sided first order differential operators $`\mathrm{\Delta }`$ by Definition 1 which obey the condition $`\overline{\delta }^c\delta ^b\mathrm{\Delta }=0`$ for all $`c,b𝒜`$. They are the first order differential operators defined in .


text:
# Comparing the ensemble mean and the ensemble standard deviation as inputs for probabilistic medium-range temperature forecasts

## 1 Introduction

We consider medium range probabilistic site-specific forecasts of temperatures. There are a number of methods that can be used to produce such forecasts. These methods differ in terms of a) the underlying numerical forecast model, b) the predictors that are taken from the forecast model and c) the statistical methods used to transform these predictors into predictions. For instance, the underlying model could be the ECMWF or the NCEP model (*inter alia*); the underlying predictors could be the ensemble mean and the standard deviation or just the ensemble mean (*inter alia*) and the statistical transform used could be based on the rank histogram or spread regression (*inter alia*).

Much work has gone into trying to understand the sources of medium range forecast skill, and to improve medium range forecasts. Some of this has aimed at determining which of the various numerical models produce the best predictors (point ”a” above). There has also been some investigation into which of the various statistical transforms are the most appropriate (point ”c” above). We will now address a third question but closely related question: what is the relative importance of the different predictors that can be derived from the numerical models (point ”b” above). In particular, we will consider the relative importance of the mean and the standard deviation of ensemble forecast temperatures, when used to make a site-specific probabilistic forecast. We consider *only* the mean and the standard deviation, rather than all the individual ensemble members, since, for the dataset we will use, it has been shown by Jewson (2003b) that there seems to be no useful information in the ensemble beyond the mean and the spread.

The statistical models we use to address this question are all derived from the spread regression model of Jewson et al. (2003). This model takes the mean and the standard deviation from an ensemble forecast and converts them into a probabilistic forecast. This model is particularly useful for addressing the question of relative importance of the different predictors since it allows us to turn each of them on and off rather easily.

We do not think that it is *a priori* obvious which of the ensemble mean and ensemble spread are more useful for making a probabilistic forecast. For instance, one can imagine dynamical systems for which the expectation varies very little, but for which the uncertainty varies a lot. Similarly one can imagine dynamical systems for which the expectation varies a lot, but for which the uncertainty varies little. In the first of these it is likely that forecasts of the expectation are more important, while in the latter it is likely that forecasts of the uncertainty would be more important. We see, therefore, that our question is, in part, a question about the dynamics of the atmosphere. One can also imagine a particular forecast system that, for whatever reason, does very well in predicting the expectation, but very poorly in predicting the uncertainty, or vice versa. From this we see that our question is also a question about specific forecast systems.

In section 2 we describe the data we use for this study. In section 3 we describe the methodology and the statistical models we will use to address this question. In section 4 we present the results of our analysis and in section 5 we discuss the implications.

## 2 Data

We will base our analyses on one year of ensemble forecast data for the weather station at London’s Heathrow airport, WMO number 03772. The forecasts are predictions of the daily average temperature, and the target days of the forecasts run from 1st January 2002 to 31st December 2002. The forecast was produced from the ECMWF model (Molteni et al., 1996) and downscaled to the airport location using a simple interpolation routine prior to our analysis. There are 51 members in the ensemble. We will compare these forecasts to the quality controlled climate values of daily average temperature for the same location as reported by the UKMO.

There is no guarantee that the forecast system was held constant throughout this period, and as a result there is no guarantee that the forecasts are in any sense stationary, quite apart from issues of seasonality. This is clearly far from ideal with respect to our attempts to build statistical interpretation models on past forecast data but is, however, unavoidable: this is the data we have to work with.

Throughout this paper all equations and all values are in terms of double anomalies (have had both the seasonal mean and the seasonal standard deviation removed). Removing the seasonal standard deviation removes most of the seasonality in the forecast error statistics, and partly justifies the use of non-seasonal parameters in the statistical models for temperature that we propose.

## 3 Models

We will address the question of whether the useful information in ensemble forecasts comes to a greater extent from the ensemble mean or from the ensemble spread by considering the site-specific forecasts for temperature at London Heathrow described above. Our forecasts are produced by four *ensemble interpretation models*, which is the name we give to the statistical models that take the non-probabilistic predictors from the ensemble forecasts and convert them into probability density functions for the predicted temperatures. Such models are essential if one wants to produce probability forecasts from ensemble output. The models we choose are ideally suited for answering the question at hand because they allow us to turn on and off the various sources of information in the ensemble very easily and see what impact that has on the skill of the resulting forecast.

The simplest model we consider consists of the use of linear regression between a *single member* of the ensemble and the observed temperature. The regression model serves to correct biases, to adjust the amplitude of the variability, and to predict the uncertainty. We will write this model as:

$$T_iN(\alpha +\beta x_i,\gamma )$$
(1)
where $`T_i`$ is the double anomaly temperature on day $`i`$ and $`x_i`$ is the ensemble member on day $`i`$. The skill of this model gives an indication of the ability of single integrations of the numerical model to predict the distribution of future atmospheric states. It is a baseline that the ensemble forecast should beat if it has any value.

Our second model is the same as the first, but the predictor is the ensemble mean, rather than a single ensemble member. We write this model as:

$$T_iN(\alpha +\beta m_i,\gamma )$$
(2)
where $`m_i`$ is the ensemble mean on day $`i`$. This model should hopefully perform better than the first model.

Our third model is an extension of the first model to include the ensemble standard deviation as a predictor of the uncertainty:

$$T_iN(\alpha +\beta x_i,\gamma +\delta s_i)$$
(3)
where $`s_i`$ is the ensemble standard deviation.

Finally our fourth model combines the second and the third models by using the ensemble mean as a predictor for the temperature, and the ensemble standard deviation as a predictor for the uncertainty around the temperature. This is the spread regression model of Jewson et al. (2003).

$$T_iN(\alpha +\beta m_i,\gamma +\delta s_i)$$
(4)
These models have been designed to allow us to answer the question at hand. For instance, comparison of the results from model one and model two, or between model three and model four, gives the benefit of using the ensemble mean versus a single integration. Comparison between model one and model three, or between model two and model four, gives the benefit of using the ensemble spread versus only using past forecast error statistics.

We fit all of these models using the standard maximum likelihood method from classical statistics. More details of this method as applied to this problem are given in Jewson (2003c).

### 3.1 Scoring the models

In order to compare the results from our four models we need a measure for how good the resulting probabilistic forecast is. To our knowledge only two scores have been suggested in the literature for comparing continuous probabilistic forecasts: the ”ignorance” by Roulston and Smith (2002), and the likelihood, by Jewson et al. (2003). It turns out that these two are, for this case, more or less the same thing. We will therefore use this score to compare our forecasts, in the form of the log-likelihood.

## 4 Results

### 4.1 Parameter values

We first look at the parameter values for our various statistical models. Figure 1 shows the values of alpha for the four models. The two models that are based on the single ensemble member have very similar values of alpha, while the two models that are based on the ensemble mean also have similar values. The values from the models based on the single ensemble member are somewhat smaller.

Figure 2 shows the values of beta for the four models. The values from the four models pair up in the same way as for the values of alpha, with the values for the single model based forecasts being much lower than the values for the ensemble mean based forecasts.

We can understand these values of beta as follows. We can consider a single integration of the numerical model to consist of the ensemble mean plus a noise term. The noise terms are at least somewhat different for the different members of the ensemble. Thus, in forming the ensemble mean the noise terms tend to cancel out and the overall noise is smaller. When we regress a single member onto observed temperatures, the regression coefficient is rather small because of the large noise term. When we regress the ensemble mean onto observed temperatures, the regression coefficient is much larger, because the noise level is lower. The mathematics of this effect have been given in detail by Jewson and Ziehmann (2002).

The value of alpha can be understood in terms of the values of beta. Alpha arises due to a combination of bias, and the beta term. We would expect the bias to be the same for a single ensemble member and the ensemble mean: the differences we see can be explained as being due to the differences in the beta.

Figure 3 shows the values of gamma for the four models. Because gamma plays a slightly different role in the different models, there is no single physical interpretation of what it means. In the two regression models gamma represents the uncertainty: we see that the uncertainty is greater in the model that only uses the single ensemble member, as would be expected. In the spread regression models, gamma represents only a part of the total uncertainty, and can only be interpreted in combination with delta, which is shown in figure 4. Confidence intervals for the value of delta in the spread regression model are given by Jewson et al. (2003), and show that there is significant sampling uncertainty about these estimates. This is not surprising since delta is related to the second moment of the second moment of observed temperatures. This sampling uncertainty is reflected in figure 4 by the jaggedness of the lines. Nevertheless we can see a general trend towards lower values of delta at high leads for both models.

### 4.2 Total Uncertainty

Comparing the gamma and delta parameter values for the four statistical models is somewhat dissatisfactory since, as we have seen, gamma does not have consistent interpretation across the four models. It makes more sense to compare the total uncertainty predicted by the four models. In figures 5 and figure 6 we compare the mean and the standard deviation of this predicted uncertainty respectively. In figure 5 we see that the two models based on the ensemble mean have much lower mean uncertainty than the two models based on the single integration, as would be expected. The ensemble mean gives significantly better forecasts than individual ensemble members, especially at longer leads. Figure 6 shows the variability in the predicted uncertainty for the two spread regression models. The two regression models, both of which ignore the ensemble spread, are not shown because they predict constant levels of uncertainty. We see that the two spread regression models predict roughly the same levels of variability in the uncertainty. There are some differences between the two, but we suspect these differences are mostly due to sampling error.

### 4.3 Likelihood scores

Finally we come to the results that answer the question that we set out to address: what are the relative contributions of the ensemble mean and the ensemble standard deviation to the skill of the final forecast? We present the comparisons in terms of the negative of the log-likehood, for which small values represent a good forecast. We compare the models in pairs to isolate individual effects, as described in section 3.

Figure 7, top left panel, shows a comparison between models one and two. This comparison is designed to illustrate the benefit of using the ensemble mean over using a single ensemble member. We see a very large benefit at all lead times. For instance, the ensemble mean based forecast at lead 6 is roughly as good as the single member based forecast at lead 4. We are left in no doubt that use of the ensemble mean is very beneficial for our probabilistic forecast. The lower left panel shows a comparison between models three and four. This comparison also shows the benefit of using the ensemble mean over using a single ensemble member, but now for the models that also use the ensemble spread. The size of the benefit is roughly the same as before.

The right hand panels show the effects of using the ensemble spread. The top right panel shows the benefit of using the ensemble spread in the models that do not use the ensemble mean, and the lower right panel shows the benefit of using the ensemble spread in the models that do use the ensemble mean. The results are roughly the same: in both cases we see only a very marginal benefit from using the ensemble spread.

Figure 8 now shows the differences that we see in the benefit due to the ensemble mean and the ensemble spread. There are a number of ways one could try and quantify this effect: we choose to look at the numerical differences in the log-likelihood. The top left hand panel shows the differences in the log-likelihood caused by the use of the ensemble mean as a predictor (from comparing models 1 and 2). The values lie between 10 and 70. The top right hand panel shows the differences in the log-likelihood caused by the use of the ensemble spread as a predictor (from comparing models 2 and 4). The values lie between zero and 10. The lower panels show the ratios of these log-likelihood differences, which quantify how much more useful the ensemble mean is as a predictor than the ensemble spread. We see the values for this ratio lie between four and 100, and increase at longer lead times. The ensemble spread is relatively most useful at the shortest lead times, while the ensemble mean is relatively most useful at the longest lead times.

## 5 Discussion

We have investigated the relative contributions of the ensemble mean and the ensemble standard deviation to the skill of probabilistic temperature forecasts. This was done by using four different statistical models to produce the probabilistic forecasts. These models differ in terms of whether they use the ensemble mean or a single ensemble member, and whether they use the ensemble spread or predict uncertainty purely using past forecast error statistics.

The results are very clear: we see that the ensemble mean is *much* more important than the ensemble spread. The difference is least at short leads, where the ensemble mean is only about four times as useful as the ensemble spread (using the particular measure that we have chosen). At medium leads the ensemble mean is between 10 and 20 times as useful as the ensemble spread, and at the longest lead the ensemble mean is roughly 100 times more useful than the ensemble spread.

There are a number of caveats to this study. In particular, we have only used data from one station. Different results would probably be obtained at other stations. We have also only used one year of forecast data, and we have seen that this leads to rather noisy estimates for some of our parameters. Our estimates would be more accurate if we had more data, as long as that data were stationary. Having said that, our main result, which is that the ensemble mean is much more useful than the ensemble standard deviation, seems to be so clear that we doubt whether it would change even if we did repeat this analysis for a number of other stations and with more data.

Finally, we note that have performed all of our analysis in-sample. This is not ideal: if possible, out of sample results are to be preferred. However, Jewson (2003a) has shown that out of sample one cannot detect *any* beneficial effects of ensemble spread whatsoever using this data-set, and so it would not be possible to perform this study if we were to do it out of sample. By using an in-sample analysis, and using statistical models with small numbers of parameters to minimize overfitting and artificial skill, we have been able to shed light on a question that could not be otherwise addressed.


text:
# XMM-Newton confirmation of a warm absorber in the 𝑧=4.4 blazar PMN J0525–3343

## 1 Introduction

PMN J0525–3343 at $`z=4.4`$ is one the the most distant, X-ray bright, radio-loud quasars known (Fabian et al., 2001b). It is one of a handful of such objects that have characteristics typical of blazars with a spectral energy distribution showing peaks in the infra-red and $`\gamma `$-ray regimes, and displaying long-term variability (Fabian et al., 1997, 1998; Moran & Helfand, 1997; Zickgraf et al., 1997; Hook & McMahon, 1998).

PMN J0525–3343 shows spectral flattening at soft X-ray energies, initially seen in observations made with ASCA and BeppoSAX (Fabian et al., 2001b) and now with XMM-Newton. Above $`2\mathrm{keV}`$ the continuum is well fitted by a power-law but this either breaks or is absorbed below $`1\mathrm{keV}`$. A similar effect has also been reported in the $`z=4.72`$ and $`z=4.28`$ blazars, GB 1428+4217 (Boller et al., 2000; Fabian et al., 2001a) and RX J1028.6–0844 (Yuan et al., 2000), respectively. If the soft X-ray deficit is attributed to cold absorption intrinsic to the objects then columns of $`10^{22}10^{23}\mathrm{cm}^2`$ are required. It was not clear, however, if the effect is due to absorption or an underlying break in the continuum.

Spectral flattening has been widely found in lower redshift radio-loud quasars. Fiore et al. (1998) found a systematic decrease with redshift of the spectral slope of soft X-ray emission from radio-loud quasars out to $`z3.9`$. Cappi et al. (1997) reported evidence of absorption for $`1.2<z<3.4`$ radio-loud quasars using ASCA observations and Reeves & Turner (2000) also found systematic spectral flattening for ASCA sources out to a redshift of 4.2.

Absorption has been suggested as the likely explanation but the data quality has been insufficient to distinguish between an intrinsic absorber at the source redshift or an absorption system along the line of sight. The latter is ruled out by statistical arguments (e.g. O’Flaherty & Jakobsen, 1997). The lack of the effect in $`z=23`$ *radio-quiet* quasars would tend to suggest an association with jet activity – either due to an intrinsic absorber associated with the jet or an underlying continuum effect. If the flattening is due to intrinsic absorption there remains much speculation as to its nature. A number of the objects show low levels of UV-optical extinction that are inconsistent with such large columns of cold absorbing material and it has been suggested that the gas may be highly ionized (Fabian et al., 2001b), similar to the ‘warm absorber’ as is common in many Seyfert galaxy nuclei (e.g. Reynolds, 1997; Crenshaw et al., 1999).

Recent XMM-Newton observations (Grupe et al., 2004) of the blazar RX J1028.6–0844 found only marginal evidence for excess absorption, much less than the high column density ($`2\times 10^{23}\mathrm{cm}^2`$; Yuan et al., 2000) implied by the ASCA spectrum. We find that the soft X-ray flattening in PMN J0525–3343 and GB 1428+4217 as detected with ASCA-BeppoSAX and ROSAT-BeppoSAX, respectively, is reduced but robustly detected with XMM-Newton observations (this work; Worsley et al., in preparation).

Alternatively, an underlying break in the continuum emission could explain the flattening; X-ray emission from powerful blazars is thought to be due to inverse Compton scattering by relativistic electrons in the jet plasma. The source photons for the scattering are supplied either by synchrotron emission from the same electron population, or, from an external radiation field of soft photons emitted from the nucleus. A spectral break can be explained in either mechanism by a low-energy cut-off in the electron population, or, in the second process, by a sharp peak in the spectrum of the soft photon field. The association of the break with the jet activity naturally explains the lack of the flattening in radio-quiet quasars.

Here we present the results of a series of XMM-Newton observations of PMN J0525–3343. We confirm the spectral flattening as detected in previous ASCA and BeppoSAX measurements (Fabian et al., 2001b) and conclude that the most likely explanation is the presence of intrinsic warm absorption.

## 2 Identification and the Optical Spectrum

PMN J0525–3343 was identified by Hook & McMahon (1998) in a program to identify radio-loud quasars at high redshift. The quasar was selected as a flat spectrum radio source, with a spectral index $`\alpha =0.1`$ in the Parkes-MIT-NRAO (PMN) Survey at $`4850\mathrm{MHz}`$ and the NRAO VLA Sky Survey (NVSS) at $`1.4\mathrm{GHz}`$ (Wright et al., 1996; Condon et al., 1998). X-ray detection was achieved with the ROSAT High Resolution Imager on 1998 March 4–6, follow-up observations were made with ASCA on 1999 March 1–2 and BeppoSAX on 2000 February 27–28.

An optical spectrum in the range $`40009000\mathrm{\AA }`$ (Péroux et al., 2001) is discussed by Fabian et al. (2001b). Strong emission lines are seen due to the Ly-$`\alpha `$ \+ N v blend along with broad ($`\mathrm{11\hspace{0.17em}500}\mathrm{km}\mathrm{s}^1`$ FWHM) emission due to C iv. A deep C iv absorption doublet ($`1548.195`$, $`1550.770\mathrm{\AA }`$) is seen, redshifted by $`2700\mathrm{km}\mathrm{s}^1`$ with respect to the peak of emission (Fig. 1). Absorption is also seen due to Si iv, N v and Ly-$`\alpha `$. The absence of an optically-thick Lyman-limit system at the quasar redshift indicates that $`N_{\mathrm{H}\text{i}}<3\times 10^{17}\mathrm{cm}^2`$. This constraint, together with the Lyman-$`\alpha `$ absorption line strength indicates the Doppler parameter, $`b`$, has to be at least $`40\mathrm{km}\mathrm{s}^1`$ (Fabian et al., 2001b). We have fitted here the C iv absorption lines and find Doppler parameters of $`b=156\pm 26\mathrm{km}\mathrm{s}^1`$ and $`b=143\pm 29\mathrm{km}\mathrm{s}^1`$ respectively, suggesting an upper limit of $`b182\mathrm{km}\mathrm{s}^1`$. The rest frame equivalent widths of $`3.52\pm 0.60\mathrm{\AA }`$ and $`3.16\pm 0.65\mathrm{\AA }`$ yield the column density in C iv from a curve of growth analysis. Fig. 2 shows the confidence limits on the column density and $`b`$ parameter. Even given the upper and lower limits on the $`b`$ parameter the column is relatively unconstrained at $`10^{1518}\mathrm{cm}^2`$.

## 3 X-ray Observations

XMM-Newton observed PMN J0525–3343 on 2001 February 11 and again on 2001 September 15. The first observation was badly affected by background flaring but the second obtained good data. These results indicated significant differences from the ASCA-BeppoSAX spectra. This prompted us to apply for future observations to search for variability and a series of six more were made between 2003 February 14 and 2003 August 8. All observations were made in full frame imaging mode with the thin 1 filter.

The EPIC data reduction was performed with the Scientific Analysis Software (sas) v5.4.1. Periods of background flaring were excluded and the data were filtered in the standard way to include only the recommended event patterns for spectral analysis and to exclude events from ‘bad’ pixels (Ehle et al., 2003). Table 1 shows the details of the observations. Background-corrected spectra were obtained from the observations, grouped to a minimum of 20 counts in each bin and analysed using xspec v11.2 (Arnaud, 1996).

Spectral fitting to a simple power-law, modified only by galactic absorption ($`N_\mathrm{H}=2.21\times 10^{20}\mathrm{cm}^2`$; Elvis et al., 1994) was made over the full energy range of $`0.312\mathrm{keV}`$ to the data from the PN and MOS cameras simultaneously. The fit was poor (with a reduced chi-squared of $`\chi _\nu ^2=1.109`$; $`\nu =1913`$) with highly curved residuals. Good fits, however, were made over the energy range $`212\mathrm{keV}`$. If the power-law photon index $`\mathrm{\Gamma }`$ is allowed to vary in the $`212\mathrm{keV}`$ range between data sets then we have $`\chi _\nu ^2=1.013`$ ($`\nu =629`$) with values of $`\mathrm{\Gamma }`$ varying from $`1.64`$ to $`1.74`$. Requiring the same value of $`\mathrm{\Gamma }`$ for all observations also results in a good $`212\mathrm{keV}`$ fit, with $`\mathrm{\Gamma }=1.67\pm 0.03`$ ($`\chi _\nu ^2=1.007`$; $`\nu =635`$).

The background-corrected spectrum from revolution 216 was very poor, particularly in the softer energies, and was excluded from further analysis. The revolution 324 observation had a $`210\mathrm{keV}`$ flux of $`9.0\times 10^{13}\mathrm{erg}\mathrm{cm}^2\mathrm{s}^1`$. During the last six observations the flux increases steadily from $`9.210.8\times 10^{13}\mathrm{erg}\mathrm{cm}^2\mathrm{s}^1`$ over a period of $`1\mathrm{month}`$ in the source frame (the error on the flux measurement in a single observation is $`0.4\times 10^{13}\mathrm{erg}\mathrm{cm}^2\mathrm{s}^1`$). For comparison the ASCA $`210\mathrm{keV}`$ flux was $`12\times 10^{13}\mathrm{erg}\mathrm{cm}^2\mathrm{s}^1`$. No statistically significant variability is present in any of the individual lightcurves.

The $`0.12.4\mathrm{keV}`$ flux increases from $`5.4\times 10^{13}\mathrm{erg}\mathrm{cm}^2\mathrm{s}^1`$ (revolution 324) to $`5.76.5\times 10^{13}\mathrm{erg}\mathrm{cm}^2\mathrm{s}^1`$ over the last six observations (flux errors are typically $`0.3\times 10^{13}\mathrm{erg}\mathrm{cm}^2\mathrm{s}^1`$). These can be compared to the ASCA and BeppoSAX values of $`5.2\times 10^{13}\mathrm{erg}\mathrm{cm}^2\mathrm{s}^1`$ and $`3.3\times 10^{13}\mathrm{erg}\mathrm{cm}^2\mathrm{s}^1`$ respectively. The earlier ROSAT count rate corresponds to a value of $`4.2\times 10^{13}\mathrm{erg}\mathrm{cm}^2\mathrm{s}^1`$ (assuming the soft energy photon index measured by ASCA-BeppoSAX).

Fig. 3 shows the ratio of data to the $`\mathrm{\Gamma }=1.67`$ fitted power-law (including Galactic absorption) as observed. Spectral flattening becomes apparent at low energies, with a clear break at $`1\mathrm{keV}`$ and the data/model ratio falling to $`0.5`$ for the softest energies. There appears to be little or no variation in spectral shape between the observations so the PN data from the different revolutions were combined together in the hope of detecting any features in the spectrum. Fig. 4 shows the results (in the rest energy frame of the blazar). A break from the power-law occurs sharply at just below $`5\mathrm{keV}`$. An emission-line-like feature seems to stand out at $`4.5\mathrm{keV}`$ with a possible dip and/or another peak at around $`4\mathrm{keV}`$. None of these correspond to any known spectral features or instrumental artefacts and given the number of data points we do not consider them statistically significant. The MOS camera spectra (Fig. 3) do not corroborate any potential features. Upper limits on the detection of iron K have been estimated using the absorption edge above $`7.1\mathrm{keV}`$. This has an optical depth of $`<0.07`$ (corresponding to a column density of $`<10^{17.6}\mathrm{cm}^2`$) while the K$`\alpha `$ emission line ($`6.40\mathrm{keV}`$) has an equivalent width $`<13\mathrm{eV}`$.

We considered a number of different models to take account of the the soft X-ray spectral flattening. The results are summarised in Table 2. All three models (broken power-law, power-law plus free local absorption and power-law plus free intrinsic absorption) produced statistically acceptable fits with no obvious irregularities in the residuals. The best-fitting models have improvements in chi-squared of between 164 (intrinsic warm/cold absorption) and 152 (broken power-law) when compared to the simple power-law (over the $`0.312\mathrm{keV}`$) range.

The fits show significant variation from the ASCA-BeppoSAX observations which found local/intrinsic total absorption columns of $`20/1800\times 10^{20}\mathrm{cm}^2`$, compared to the $`5.9/130\times 10^{20}\mathrm{cm}^2`$ found here. The broken power-law fits to the ASCA-BeppoSAX spectra also point to a stronger soft energy flattening; they require a break at a higher energy ($`1.72\mathrm{keV}`$) and to a flatter photon index ($`0.9`$). The ratio of the ASCA-BeppoSAX $`(210)/(0.12.4)\mathrm{keV}`$ flux is $`2.3`$ compared to $`1.7`$ for the XMM-Newton data. The ASCA-BeppoSAX observations show a $`210\mathrm{keV}`$ photon index of $`\mathrm{\Gamma }=1.471.68`$ (broken power-law model; 90 per cent confidence) which is consistent with that found here ($`\mathrm{\Gamma }=1.67`$).

To summarise, there is no evidence for variation in the spectral slope of the hard X-ray emission from either the XMM-Newton or earlier observations. The $`210\mathrm{keV}`$ flux would appear to vary on $`6\mathrm{month}`$ time-scales ($`1\mathrm{month}`$ in the blazar frame) with amplitudes $`1020`$ per cent. The $`0.12.4\mathrm{keV}`$ flux also varies but the soft X-ray deficit appears to have lessened in the 1/2 years ($`2/4\mathrm{months}`$ in the blazar frame) between the BeppoSAX-XMM-Newton and ASCA-XMM-Newton observations, characterised by a drop in the intrinsic cold absorber column by an order of magnitude.

## 4 Discussion

The four main potential explanations of the soft X-ray deficit will be examined in turn. Only two of the four scenarios (a break in the continuum or a warm absorber) appear to be self-consistent, the excess intergalactic and cold intrinsic absorber models both appear to be contradicted by other results.

### 4.1 Break in the Underlying Power-Law

An acceptable fit to the data can be obtained by simply allowing a break in the underlying power-law. The best fit was found with a break occurring at $`E_{\mathrm{break}}=0.84\mathrm{keV}`$ where the photon index hardens from $`1.61`$ to $`1.2`$ in the soft X-ray band.

The spectral energy distribution (SED) of blazars is composed of two components; beamed synchrotron emission (peaking in the infra-red to soft X-ray regime) and a high-energy component (peaking in the X-ray up to $`\gamma `$-ray regime), see e.g. Fossati et al. (1998). The likely mechanism for the production of high energy photons is the inverse Comptonization of relatively soft ‘seed’ photons, by relativistic electrons in the jet (Sikora et al., 1994). The two sources for the seed photons are either the locally produced synchrotron emission, or an intense photon field (typically optical/UV) external to the jet, probably radiated from the accretion disc and reprocessed by circumnuclear gas. For the brightest blazars (as in PMN J0525–3343) the latter, the so-called EC (external Compton) mechanism, is thought to dominate over the former, the so-called SSC (synchrotron self-Compton) mechanism.

Spectral breaks have been observed before in X-rays: Concave spectra arise when the observed band spans the steep synchrotron cut-off ($`\mathrm{\Gamma }23`$) before the flat inverse Compton component dominates with $`\mathrm{\Gamma }1.51.8`$. A convex spectrum occurs when the band is dominated by the synchrotron component, the synchrotron rise–peak–cutoff leading to a gradual steepening of the spectral slope (from as flat as $`\mathrm{\Gamma }1.5`$ through the peak at $`\mathrm{\Gamma }2`$ and into the cut-off with $`\mathrm{\Gamma }23`$ or more). PMN J0525–3343 shows a convex spectrum with a *very* flat soft X-ray slope ($`\mathrm{\Gamma }1.2`$). However, the high energy photon index ($`\mathrm{\Gamma }1.61.7`$) is still far too flat to be the synchrotron peak/cut-off. The sharply-flattening convex spectrum is not explicable through any known combination of blazar emission components.

As the hard X-rays observed are compatible with being due to inverse Compton emission, a flattening of the spectrum can be explained by a low-energy cut-off in the electron population of the jet. A cut-off can theoretically produce a maximum flattening to a photon index of up to $`\mathrm{\Gamma }0`$ (Fabian et al., 2001b) and is therefore not inconsistent with the observations here. It would, however, be very difficult to reproduce such a sharp break (the break is very distinct and occurs over an energy span $`<5\mathrm{keV}`$ in the source frame) without a highly energy dependent (and unknown) re-acceleration process. Further difficulties arise because there is no reported evidence of any breaks in nearby radio-loud quasars.

An alternative explanation (applicable only to the inverse Compton emission due to the EC scenario, as is believed to be the case here) requires the seed photon distribution to be sharply peaked at or decrease rapidly below a certain frequency (e.g. Ghisellini, 1996). This mechanism may only be important in the more distant (brighter) radio-loud quasars where the EC emission may dominate over the SSC, hence the lack of an observable break in nearby (weaker) object where the SSC emission could mask the EC spectral decline. As with the alternative mechanism it remains far from clear how to account for the sharpness of the break, even considering the fact that we only see EC photons scattered in a small range of angles. In summary, while we cannot eliminate the broken power-law possibility we consider it very unlikely.

### 4.2 Excess Cold Galactic/Intergalactic Absorption

The flattening can be accounted for by line of sight absorption by cold material with a column of approximately $`N_\mathrm{H}3.7\times 10^{20}\mathrm{cm}^2`$ (in addition to the Galactic level). It is very unlikely that such a large column, localised in the direction of PMN J0525–3343, could be attributable to some kind of Galactic feature, particularly given the large mass of material which would be implied.

Intergalactic absorption systems remain a possibility, in particular damped Ly-$`\alpha `$ systems (DLAs), which can have neutral hydrogen columns of up to several times $`10^{21}\mathrm{cm}^2`$. In order to reproduce the observed absorption the line of sight to the quasar would need to intercept one or more very high column density DLAs of which the chances are low (O’Flaherty & Jakobsen, 1997; Zwaan et al., 1999). Such an explanation, although it cannot be ruled out, is undermined further by the fact that the soft X-ray flattening has been observed in several other high-$`z`$ quasars and has not been seen in $`z=23`$ radio-quiet quasars.

### 4.3 Intrinsic Cold Absorption

For a cold absorber to be viable it must be intrinsic to the source, perhaps arising as the result of some kind of outflow after the initial formation of the central black hole or just because the host galaxy is young and gas-rich.

An intrinsic cold absorber at the source resulted in an acceptable fit to the X-ray data, requiring a column in neutral hydrogen of approximately $`1.3\times 10^{22}\mathrm{cm}^2`$ (Table 2). Such a large amount of obscuring material is completely at odds with the optical observations which imply a column of $`<3\times 10^{17}\mathrm{cm}^2`$ (based upon the lack of a Lyman-limit absorption system; Fabian et al., 2001b). Supersolar metal abundances would appear to unable to account for the factor of $`10^5`$ difference required for $`N_{\mathrm{H}\text{i}}`$. Moreover, the dust in any such column of cold material would be expected to provide strong UV extinction, contrary to observation.

### 4.4 Warm Absorption

A warm (ionized) absorber is able to reproduce both the high X-ray and low optical/UV column densities. Photoionized gas could produce significant photoelectric absorption in the soft X-ray regime. It is likely any dust present in the nuclear regions would have been sublimated by the high source flux or high temperatures, resulting in very low absorption in the optical/UV band.

The optical spectrum shows strong C iv absorption and, whilst the nature of the absorption is unclear, it is still evidence for large quantities of highly ionized material. Indeed, strong C iv features are common in the spectra of Seyfert galaxy nuclei that have X-ray warm absorbers (Crenshaw et al., 1999).

We constructed warm absorber spectral models using the photoionization code cloudy (Ferland et al., 1998). We varied the underlying power-law photon index $`\mathrm{\Gamma }`$, the column density $`N_\mathrm{H}`$ and the ionization parameter $`\xi =L_{210}/nR^2`$ ($`L_{210}`$ is the inferred isotropic luminosity in the $`210\mathrm{keV}`$ rest frame of the blazar, $`n`$ and $`R`$ are the absorber density and the distance from the source respectively).

A grid of spectra were calculated by varying the ionization parameter, column density and photon index. Solar abundances have been assumed. Assuming no spectral break intrinsic to the source the inferred isotropic luminosity, in the $`210\mathrm{keV}`$ blazar frame, is $`L_{210}1\times 10^{47}\mathrm{erg}\mathrm{s}^1`$.

Spectral fits were made to the X-ray spectra. Fig. 5 shows the confidence regions for the warm absorber ionization parameter and column density with a power-law index of $`\mathrm{\Gamma }=1.64`$. The best fit is at $`\mathrm{log}(\xi /\mathrm{erg}\mathrm{cm}\mathrm{s}^1)=0.2`$ and $`\mathrm{log}(N_\mathrm{H}/\mathrm{cm}^2)=22.3`$ but the confidence contours clearly show it was impossible to successfully constrain the parameters. The data is entirely consistent with a warm absorber scenario provided that $`\mathrm{log}(\xi /\mathrm{erg}\mathrm{cm}\mathrm{s}^1)0.8`$ (for a column of $`10^{22.3}\mathrm{cm}^2`$) so that the H i edge has an optical depth $`\tau <1`$. In this case the warm absorber must have $`\mathrm{log}(\xi /\mathrm{erg}\mathrm{cm}\mathrm{s}^1)0.81.2`$ and $`\mathrm{log}(N_\mathrm{H}/\mathrm{cm}^2)=22.322.5`$. These models predict C iv column densities $`10^{14.518.4}\mathrm{cm}^2`$, consistent with the $`10^{1518}\mathrm{cm}^2`$ implied by the equivalent width of C iv line observed in the optical data. The fit of a $`\xi =10\mathrm{erg}\mathrm{cm}\mathrm{s}^1`$ $`N_\mathrm{H}=10^{22.4}\mathrm{cm}^2`$ warm absorber model to the soft X-ray flattening is shown in Fig. 4. The mass of the absorber is $`3\times 10^9R_{\mathrm{kpc}}^2\mathrm{\Omega }M_{}`$, where $`\mathrm{\Omega }`$ is the solid angle subtended by the absorber at the source.

We note the UV absorbers seen in Seyfert galaxies are not the same matter giving rise to the X-ray absorption (e.g. Crenshaw et al., 2003), but are probably dense clouds. Higher resolution spectra of the UV absorption in PMN J0525–3343 is required to determine whether it too is due to separate clouds at different velocities, or to a smooth flow.

It is unfortunate to note that the observed X-ray spectra contain no definite warm absorber signatures e.g. neutral or slightly ionized iron edges, which should occur at observed energies in the $`12\mathrm{keV}`$ range. Given the estimated column and ionization parameter the photoionization models do predict various absorption edges, principally due to oxygen and neon species, but these are in the $`0.10.25\mathrm{keV}`$ range, redshifted out of the observed band. The continuum should recover below the oxygen absorption edge at $`0.135\mathrm{keV}`$ (our frame) which would be a convincing test of the warm absorber model. The MOS data are consistent with such a recovery but the calibration here is so uncertain that further work is required.

## 5 Conclusions

We have investigated the X-ray spectrum of the high redshift blazar PMN J0525–3343 with the XMM-Newton observatory. We confirm the presence of significant soft X-ray spectral flattening, as previously detected in ROSAT, ASCA and BeppoSAX observations. The flattening can be explained by either excess absorption (either intergalactic or intrinsic to the source) or by an underlying break in the expected continuum emission (of a type not seen before). Optical observations, statistical arguments and the lack of any effect in radio-quiet quasars makes any kind of cold absorber seem unlikely – the most plausible explanations remain either a highly ionized warm absorber at the source, or, a spectral break.

A low energy cut-off in the emitting particle distribution responsible for the blazar inverse Compton emission could be a potential mechanism for the spectral break – any such a cut-off would have significant importance in our knowledge of jet physics. The difficulty in this explanation is the lack of the effect in nearby radio-loud quasars, which show no evidence of a similar spectral break (which should be readily detectable at what would be energies $`5\mathrm{keV}`$). An alternative, which would apply to the external Compton (EC) mechanism (as is believed to be important in very bright blazars such as PMN J0525–3343), attributes the flattening to a highly peaked seed photon distribution. Neither process however, is able to easily account for the sharpness of the spectral break observed in PMN J0525–3343 which occurs over a source frame energy range of less than $`5\mathrm{keV}`$.

The presence of strong C iv absorption, the lack of a Lyman-limit absorption system in the observed optical spectrum and low levels of UV reddening makes us favour a warm absorption model. Warm absorbers are a common feature of many observed Seyfert galaxies, where high ($`10^{2123}\mathrm{cm}^2`$) column densities of highly ionized material reside close to the nucleus. We find the X-ray flattening is entirely consistent with a column $`10^{22.322.5}\mathrm{cm}^2`$ of gas with an ionization parameter $`10^{0.81.2}\mathrm{erg}\mathrm{cm}\mathrm{s}^1`$. The equivalent width of the absorption implies C iv column densities in the range $`10^{1518}\mathrm{cm}^2`$, in agreement with the X-ray predicted $`10^{14.518.4}\mathrm{cm}^2`$. A warm absorber also enables Ly-$`\alpha `$ and UV emission from the nucleus to be readily detectable. Further work on the C iv absorber is required to refine the properties of the absorbing column.

## 6 Acknowledgments

Based on observations with XMM-Newton, an ESA science mission with instruments and contributions directly funded by ESA Member States and the USA (NASA). MAW and AKT acknowledge support from PPARC. ACF thanks the Royal Society for support. AC thanks the Italian MIUR and ASI for financial support.