text:
# Martensitic transition and magnetoresistance in a Cu-Al-Mn shape memory alloy. Influence of ageing.

## I INTRODUCTION

The shape-memory effect Otsuka98  is characteristic of certain alloys (the so-called shape-memory alloys), which exhibit a martensitic transition (MT) from an ordered $`bcc`$ -phase ($`\beta `$-phase) towards a close-packed low-temperature phase ($`M`$-phase). This effect is related to unique thermomechanical properties such as the ability to recover from large permanent deformations produced in the $`M`$-phase by the reverse transition when temperature is increased. During the nineties a great deal of interest has been devoted to the study and development of magnetic shape-memory materials. This interest is mostly due to the possibility of a magnetic control of the shape-memory effect, which has been made evident in the ferromagnetic Ni-Mn-Ga alloy close to the Heusler composition Ni2MnGa Ullakko96 . Furthermore, these materials have been shown to exhibit unexpected pretransitional behaviour Zheludev95 .

The present paper deals with the study of the Cu-Al-Mn alloy. This Hume-Rothery material Ahlers96  shares a number of features with the Ni-Mn-Ga alloy system. It displays the same high-temperature crystallographic structure, a martensitic transformation with associated shape-memory effect and interesting magnetic properties. For Cu-Al-Mn, however, the martensitic transition occurs in a composition region far from the Heusler stoichiometry. For this composition range, the $`\beta `$-phase is only stable at high temperatures, but can be retained at low-temperature by means of suitable cooling. During this cooling the system develops an ordered $`L2_1`$ structure ($`Fm3m`$, Heusler symmetry) in two successive disorder-order transitions: $`A2`$ $`(Im3m)`$ $``$ $`B2`$ $`(Pm3m)`$ at $`T_{c1}`$ and $`B2`$ $``$ $`L2_1`$ at $`T_{c2}`$ Obrado98 . Upon further cooling it undergoes a martensitic transition at a temperature which is strongly composition dependent. This transition has a diffusionless nature which ensures that the atomic distribution of the $`L2_1`$ phase is inherited by the $`M`$-phase. It is worth noting that this feature is common to all Cu-based shape-memory materials Ahlers86 .

Magnetic properties arise from localized magnetic moments at Mn-atoms as occurs in the Ni-Mn-Ga system Webster88 . These magnetic moments are coupled through an oscillating effective interaction (RKKY interaction). ALCHEMI (Atom Location by Channelling Enhanced Microanalysis) experiments Nakanishi93  have shown that Mn atoms are located preferentially in one of the four distinguishable $`fcc`$ sublattices of the $`L2_1`$ structure (the 4b sites in Wyckoff notation). For this configuration, ferromagnetic coupling is dominant and close to the Cu2AlMn composition the system is ferromagnetic. However, for non-stoichiometric alloys the 4b sites are not fully occupied by Mn-atoms and this results in magnetic disorder, which gives rise to different magnetic behaviour depending on the temperature range Prado98 ; Obrado99 . Magnetic clustering has been suggested to be at the origin of the magnetoresistive properties recently reported in Cu-Al-Mn melt-spun ribbons Murthy95 ; Sugimoto98 . An interesting feature is the fact that a phase separation between Cu3Al-rich and Cu2AlMn-rich phases may occur below the $`L2_1`$ ordering line. It is therefore expected that the magnetic and structural properties of the system are sensitive to the temperature history of the material (i.e., ageing). In particular, isothermal annealing at a temperature within the miscibility gap will result in the growth of magnetic clusters. The existence of a miscibility gap was first reported for the Cu3Al $``$ Cu2AlMn pseudobinary composition line Bouchard75  and later confirmed for Cu-rich systems with a composition that slightly deviates from this line Kainuma98 . Interestingly, this composition range includes that for which Cu-Al-Mn displays a martensitic transition. Recently, it has been theoretically shown Marcos01  that the influence of magnetic degrees of freedom on configurational phase stability is at the origin of this phase separation.

The present paper is aimed at experimentally studying the influence of magnetism on the martensitic transformation in Cu-Al-Mn. Since the magnetic properties of this material are expected to be sensitive to ageing, this effect will be studied. In particular, we focus on the study of magnetotransport and magnetic properties through the martensitic transition of bulk samples. The results will contribute to gain a deeper understanding of the magnetoelastic interplay, which will be studied at different levels of coupling.

The paper is organized as follows. In section II, the experimental details are outlined. Section III deals with the experimental results, which are discussed in section IV.

## II EXPERIMENTAL DETAILS

Measurements were performed on a Cu-Al-Mn polycrystal (grain size $`100\mu `$m) prepared by melting pure elements (99.99$`\%`$ purity). The nominal composition of the studied alloy is Cu; 22.8 at$`\%`$ Al; 9.0 at $`\%`$ Mn. From the ingot, rectangular specimens (approximately 12 $`\times `$ 4 mm2 and 0.01 mm thick for resistance measurements or 0.5 mm for calorimetric and magnetic measurements) were first cut with a low-speed diamond saw and mechanically polished down to the desired thickness. All samples were annealed for 10 min at 1080 K and quenched in a mixture of ice and water. This fast cooling avoids precipitation of equilibrium phases, but enables the two ordering transitions to the $`B_2`$ and $`L2_1`$ structures to take place ($`T_{B_2}`$ = 849 K and $`T_{L2_1}`$ = 795 K Obrado98 ). The nominal martensitic transition temperature of as-quenched samples is $`T_M`$ = 157 $`\pm `$ 1 K. For the studied composition the structure of the $`M`$-phase is 18$`R`$ Obrado97 . This structure is monoclinic but it is usually described by a larger (approximately orthorhombic) unit cell containing 18 close-packed atomic planes along the c-axis Lovey87 .

Phase separation between Cu3Al-rich and Cu2AlMn-rich domains occurs under very slow cooling from high temperature. In the as-quenched specimens phase separation can be induced by means of post annealing within the miscibility gap. This process takes place at slow rates so that it is completely negligible at room temperature (even on a month time scale). In all cases, ageing consist in annealing at a temperature $`T_a=`$ 473 K for selected times. We chose this temperature because it is located slightly below the limit of the miscibility gap for the studied composition. Actually, above 525 K no phase separation was detected, while 50 K below, the magnitude of the observed effects is similar, but they occur on a longer time scale.

Four kind of measurements were performed: electrical resistance under applied magnetic field, magnetization, ac-susceptibility and calorimetry. Electrical resistance was measured from 5 K up to 300 K using an ac four-probe method. For magnetization measurements a SQUID magnetometer was used. For calorimetric measurements a highly sensitive and fast response (around 20 s time constant) calorimeter, which was specifically designed for the study of solid-solid phase transitions, was utilized Manosa96 . Calorimetric, magnetization and ac-susceptibility runs were performed through the martensitic transition in the range from 100 K to room temperature. Calorimetric measurements were carried out at a rate of 0.5 K/min, and pairs of data (calorimetric output and temperature) were recorded every 2 seconds. Ac-susceptibility measurements were performed at a frequency $`f`$ = 66 Hz.

## III RESULTS

Electrical resistance measurements were carried out by first cooling the sample down to 5 K and then increasing the temperature in steps up to 300 K. At each plateau, the resistance $`R`$ was measured at different values of the magnetic field $``$ at intervals of 345 Oe from 0 to 10 kOe, and at 2100 Oe intervals within the range from 10 to 50 kOe. The magnetoresistance $`MR`$, defined as the relative change of $`R`$ with $``$, is computed as: $`MR=\left[R(T,)R(T,=0)\right]/R(T,=0)`$. In Fig.1, $`MR`$ is plotted as a function of $``$ at selected temperatures for samples subjected to increasing annealing times at $`T_a`$. In all cases, the alloy exhibits negative magnetoresistance, i.e., the resistance decreases as magnetic field increases and the magnitude of the change is higher at low temperature. Actually, the maximum $`MR`$, of about 7$`\%`$, was obtained at $`T10`$ K and $``$ = 50 kOe for the as-quenched specimen. Overall, the effect of ageing is to reduce the $`MR`$ at high fields. However, at low fields a small increase is observed, which is associated with the material becoming magnetically softer. For instance, the $`MR`$ at $`T=300`$ K and $`=10`$ kOe increases from a practically null value for the as-quenched state to $`0.5\%`$ for the long-term annealed state.

An interesting feature concerns the influence of the martensitic transition on the $`MR`$. This effect is better revealed by plotting the change of $`R`$ with $`T`$ across the MT at selected values of $``$. This is shown in Fig. 2 for an as-quenched sample and for a sample aged for 300 min. The curves obtained are equivalent to those obtained directly by measuring $`R`$ versus $`T`$ under continuous heating at different constant values of $``$ Marcos01b . The large change in resistance in the temperature range between 150 K and 175 K for the as-quenched sample and between 175 K and 200 K for the annealed sample is due to the occurrence of the MT. We will take the temperature difference $`\mathrm{\Delta }T`$ between the maximum and the minimum of the $`R`$ versus $`T`$ curves as an estimation of the spread in temperature of the transition (difference in the starting and finishing transition temperatures). The corresponding change in resistance will be denoted by $`\mathrm{\Delta }R`$. Regardless of the ageing influence, the effect of the magnetic field is mainly to decrease the resistance of the $`M`$-phase (the resistance in the $`\beta `$-phase decreases by a lower amount), which leads to a reduction in $`\mathrm{\Delta }R`$ as the field is increased. In addition, $`\mathrm{\Delta }T`$ is reduced when the field increases (in the as-quenched system this reduction amounts to $``$6 K for an applied field of 50 kOe).

In order to characterize the magnetic response of the system we carried out magnetization ($``$) and ac-susceptibility ($`\chi `$) measurements through the martensitic transition. These measurements were undertaken on samples subjected to increasing annealing times at $`T_a`$. Figure 3(a) shows $``$ versus $`T`$ curves at selected values of the magnetic field. Figure 3(b) shows magnetization versus $`/T`$ curves. It is interesting to correlate the magnetization with the behaviour of the magnetoresistance. Figure 4 shows the $`MR`$ as a function of the square of magnetization $`^2`$ for as-quenched and long-term aged specimens. For each sample, all data for the $`\beta `$-phase scale on a single curve, while data for the $`M`$-phase scale on another curve. In the as-quenched state, the two scaling curves are approximately linear, while remarkable deviations from linearity are evident in the annealed system, particularly for high values of $``$. Moreover, in the annealed state, the scaling (in both $`\beta `$\- and $`M`$-phases) is not as good as in the as-quenched state.

As illustrated in Fig. 3(a), magnetization curves show a significant change $`\mathrm{\Delta }`$ at the MT ($`\mathrm{\Delta }`$ is estimated as the magnetization difference at the transition temperature between extrapolations of the linear behaviour of $``$ versus $`T`$ curves well above and well below the transition). Fig.5 gives $`\mathrm{\Delta }`$ versus $``$ after different annealing times at $`T_a`$. In the as-quenched state, $`\mathrm{\Delta }`$ saturates for a value of $``$ close to 20 kOe; the field needed to reach saturation decreases with increasing annealing time. In all cases, for these values of the field, the magnetization itself is far from having reached saturation. Moreover, the saturation value of $`\mathrm{\Delta }`$ ($`\mathrm{\Delta }_{sat}`$) is strongly dependent on the annealing time. Actually, it decreases with $`t_a`$ as shown in the inset of Fig. 5, reaching, within error, a constant value for $`t_a>`$ 150 min.

As the $`MT`$ is first order, any change of the temperature $`T_M`$ induced by the applied magnetic field $``$ should be accounted for by the Clausius-Clapeyron equation

$$\frac{dT_M}{d}=\frac{\mathrm{\Delta }}{\mathrm{\Delta }S},$$
(1)
where $`\mathrm{\Delta }S`$ is the entropy difference between $`M`$\- and $`\beta `$-phases. We measured $`\mathrm{\Delta }S`$ using calorimetry after different annealing times. In Fig. 6, we present thermal curves obtained during heating runs for different annealing times. These curves show the characteristic noisy structure which is common to many martensitic transitions Clapp ; Levin . Such a noisy structure is enhanced in the as-quenched state. As the annealing time increases, the transformation first becomes smoother and broader, and shifts towards lower temperatures. However, for times longer than 90 min thermal curves partially recover the original noisy character and the transformation shifts towards higher temperatures. The entropy change is determined by a numeric integration of the thermograms Ortin88 . The value obtained for the as-quenched system (1.27 $`\pm `$ 0.02 J/Kmol) is consistent with that reported in Prado95 . The dependence of $`\mathrm{\Delta }S`$ on $`t_a`$ is depicted in Fig. 7. The figure reveals that $`\mathrm{\Delta }S`$ (absolute value) decreases when the annealing time at $`T_a`$ increases, reaching (as occurs with $`\mathrm{\Delta }_{sat}`$) a constant value for $`t_a>`$ 150 min. The relative decrease of $`\mathrm{\Delta }S`$ is about 30$`\%`$, much smaller than the relative decrease $`\mathrm{\Delta }_{sat}`$, which is estimated to be approximately 90$`\%`$.

It is then interesting to compare the values of the derivative $`dT_M/d`$ in the as-quenched and long-term annealed states using equation 1. For the as-quenched system, a value of ($`8\pm 0.4`$) mK/kOe is obtained, while this value reduces to ($`1.10\pm 0.05`$) mK/kOe after 300 min of annealing. The reduction is mainly related to the decrease of $`\mathrm{\Delta }_{sat}`$ with ageing. Actually, the derivative $`dT_M/d`$ provides a good quantification of the dependence of the structural transition temperature on a magnetic field under the assumption that $`\mathrm{\Delta }S`$ is independent of $``$. Such an assumption is based on the fact that the entropy difference between the open and close-packed phases originates from the change of the corresponding vibrational spectrum, as occurs in non-magnetic Cu-based shape-memory alloys Planes01 .

A possible magnetic contribution to $`\mathrm{\Delta }S`$ can be evaluated from the magnetization vs. temperature curves. Actually, at a given temperature $`T`$, the change of entropy of the system with the magnetic field $``$ is given by:

$$\delta S(T,)S(T,)S(T,=0)=_0^{}\left(\frac{}{T}\right)_{}𝑑,$$
(2)
where the thermodynamic (Maxwell) relation:

$$\left(\frac{S}{}\right)_T=\left(\frac{}{T}\right)_{},$$
(3)
has been used. Therefore, application of a magnetic field modifies the entropy change between both the $`M`$\- and $`\beta `$-phases according to:

$$\begin{array}{ccc}\delta \mathrm{\Delta }S& =& \mathrm{\Delta }S()\mathrm{\Delta }S(=0)\hfill \\ & =& _0^{}\left[\left(\frac{^M}{T}\right)_{}\left(\frac{^\beta }{T}\right)_{}\right]𝑑,\hfill \end{array}$$
(4)
where $`\left(\frac{^M}{T}\right)_{}`$ and $`\left(\frac{^\beta }{T}\right)_{}`$ are the derivatives of the magnetization with respect to the temperature in the $`M`$\- and $`\beta `$-phases. The preceding integral was calculated numerically from the magnetization curves above following the procedure given in Pecharsky99 . Results are plotted in Fig.8 as a function of $``$ for different annealing times. This quantity is always negative, which means that the effect of the field is to increase the absolute value of $`\mathrm{\Delta }S`$. However, the increase is very small, and represents, in all cases, less than 1$`\%`$ of the total entropy change. It can therefore be considered negligible for practical purposes. It is, however, interesting to compare this entropy difference with that reported in a Ni51.5Mn22.7Ga25.8 polycrystalline alloy Hu00 . In that case a value of $`\delta \mathrm{\Delta }S4\times 10^3`$ J/K mol is estimated at $`=0.9`$ T. This value is comparable with the value reported here for the Cu-Al-Mn system. Notice that for the Ni-Mn-Ga system, the change in MT temperature with magnetic field is also very weak Planes01 .

An additional, complementary magnetic characterization of the studied system consisted of measuring the ac magnetic susceptibility as a function of the ageing time $`t_a`$. The evolution of the real part of the ac-susceptibility is displayed in Fig. 9. In the as-quenched state, the magnetic susceptibility is small and shows a low-temperature peak associated with the appearence of a glassy magnetic phase. The martensitic transition can be detected on this curve as a very small jump between 160 K and 170 K. The inset of the figure shows a detailed view of the region where the $`MT`$ takes place. In this case, heating and cooling runs are presented so that the transformation hysteresis is revealed.

As ageing time is increased, the susceptibility peak moves towards higher temperatures and becomes simultaneously higher and broader. Besides, the evolution of the martensitic transition temperature with $`t_a`$ follows that of the evolution determined from the calorimetric measurements presented above. For annealing times $`t_a>`$ 60 min the anomalies associated with the freezing temperature and the martensitic (jump) phases overlap. This fact gives rise to a susceptibility that decreases abruptly at an intermediate temperature ($``$ 150 K) as temperature is reduced. With further annealing, the shape of the susceptibility curves remains unchanged, but the value of $`\chi `$ decreases at high temperatures indicating that the system becomes more and more ferromagnetic.

## IV DISCUSSION

In this paper, we have reported results showing the existence of magnetoresistive effects in bulk polycrystalline Cu-Al-Mn samples. Magnetoresistive behaviour had already been observed in $`\beta `$-phase melt-spun ribbons of similar composition Kainuma98 . In the present paper, the investigation has been extended to the bulk system and to the study of magnetoresistance in the martensitic phase. In both $`\beta `$\- and $`M`$-phases, the observed magnetoresistance is typical of heterogeneous alloys Berkowitz92  where magnetoresistance originates from the spin-dependent scattering of conduction electrons on ferromagnetic clusters embedded in a metallic non-magnetic matrix Camley89 . In the studied system, the ferromagnetic clusters are due to the tendency of Mn-atoms to increase their number of third-order neighbours of the same species within the basic $`bcc`$-lattice. At present there is substantial experimental evidence of the existence of these clusters Bouchard75 ; Kainuma98 . X-ray experiments Obrado98  have demonstrated that quenched samples exhibit an $`L2_1`$ structure. In this structure, the Mn-atoms are preferentially located on the 4b sites (Wyckoff notation) of the cubic unit cell of the $`Fm3m`$ structure. However, the distribution of Mn atoms in this sublattice is not completely random. Actually, there are indications of the existence of correlations among the positions of the Mn-atoms within the sublattice Prado98 . While clusters are not strictly well defined in this case, their influence on transport and magnetic properties of the system has to be considered. In particular, they are responsible for the superparamagnetic behaviour of the studied alloy system, as shown in Obrado99 . The clustering tendency is enhanced by ageing within the miscibility gap. Actually, the influence of ageing on the magnetic properties of Cu-Al-Mn shows up by the remarkable effect that ageing has on the $`MR`$ curves; a decrease in the magnitude of $`MR`$ at low temperatures and a slight increase at high temperatures have been observed (see Fig. 1). These modifications are due to an increase in the magnetic softening which takes place after short annealing times. This is clearly a sign of the increase of the ferromagnetic character of the system owing to the growth of magnetic clusters. Such an increase is confirmed by magnetization and susceptibility measurements.

The scaling of $`MR`$ with $`^2`$ was found to be different in both the $`\beta `$\- and $`M`$-phases. This effect mainly arises from the different behaviour of the magnetization in the two phases. Actually, the $``$ vs. $`/T`$ curves also show a tendency to group into two different families corresponding to the $`\beta `$\- and $`M`$-phases, respectively. In the as-quenched state, the $`MR`$ shows a fairly good linear dependence on $`^2`$ for both the $`\beta `$\- and $`M`$-phases. This is the expected behaviour in the absence of magnetic interaction between magnetic clusters Wang94 . However, at high fields the curves for the $`M`$-phase show a clear tendency to deviate from linear behaviour. A similar deviation has previously been observed in heterogeneous alloys and it has been attributed to the progressive field alignment of the disordered spins at the boundaries of the magnetic clusters, which causes a large variation of the $`MR`$ compared with the corresponding change in magnetization Bellouard94 . Such an increase in the slope of the $`MR`$ vs. $`^2`$ curves is usually associated with the existence of high values of the high field susceptibility, as is the case of the as-quenched and short-term annealed samples \[see magnetization curves in Fig. 3(b)\]. After long time annealing, ferromagnetism is responsible for the loss of linearity in the $`MR`$ versus $`^2`$ curves. Actually, the formation of large ferromagnetic clusters whose magnetization saturates at relatively low fields (see Fig. 3(b)) leads to a large change of the total magnetization of the sample and a small variation of $`MR`$ at low fields, since the latter originates from the moment alignment of the small clusters and/or the spins at the boundaries of the ferromagnetic regions, which occurs at higher fields. Consequently, $`MR`$ vs. $`^2`$ shows a quasi-constant regime that corresponds to the range of fields in which the magnetic saturation of large ferromagnetic clusters takes place (see Fig. 4) and a subsequent quasi-linear behaviour due to the magnetoresistance associated with small clusters and boundary spins.

A noticeable effect of the magnetic field on the martensitic transition is the reduction of the temperature range $`\mathrm{\Delta }T`$ in which the transition occurs. A qualitatively comparable effect has been reported in the case of the ferromagnetic Ni2MnGa Ullakko96 , but in this material the reduction is much greater; for that alloy the transition range reduces from about 10 K without an applied magnetic field to less than 2 K for an applied field of 10 kOe. In shape-memory alloys, the transformation is acknowledged to be thermoelastic, which means that the system needs to be continously cooled down (or heated up) in order to increase the transformed fraction of the new phase. The transition is not completed until the temperature is lowered (increased) below (above) a certain value. The free-energy difference between both phases provides the driving force for the transition. Such a driving force is proportional to the temperature difference between equilibrium and the actual temperature to a first approximation equilibrium . The path followed by the system is an optimal path in which accommodation of the transformation elastic strain is almost accomplished. Therefore, $`\mathrm{\Delta }T`$ represents a reasonable measure of the stored elastic energy during the forward transformation Ortin88 . From this point of view, a decrease in $`\mathrm{\Delta }T`$ reflects a lower value of this stored elastic energy. We argue that the physical mechanism behind this effect is that the magnetic field around magnetic clusters breaks the degeneracy of the low-temperature martensitic phase thus favouring the nucleation of these variants with the magnetization easy-axis in the direction of the field. In Cu-Al-Mn, approximately orthorhombic martensitic domains, with the $`c`$-axis oriented along the magnetic field, are expected to have maximum nucleation probability. In the Ni-Mn-Ga alloy such an interpretation is corroborated by the fact that magnetic and mechanical energies required to induce a single variant martensitic structure from a polyvariant crystal are comparable O'Handley98 .

An interesting result concerns the magnetization jump occurring at the martensitic transition. This jump is positive (increase of the magnetization in the forward transition from the $`\beta `$\- to the $`M`$-phase) in the as-quenched state even for small fields. However, for sufficiently long annealing time, when the sample behaves ferromagnetically, this jump is negative for small values of the field and becomes positive above a certain critical field close to the value for which $`\mathrm{\Delta }`$ saturates. This is related to a high magnetic moment of the Mn-clusters and to the strong magnetic anisotropy of the martensitic phase. For fields larger than the critical field, reorientation of martensitic variants enables the Zeeman energy to be minimized and the jump again becomes positive.

The behaviour of magnetic susceptibility provides further understanding of the origin of the interplay between magnetic and structural degrees of freedom. In the as-quenched state, the temperature dependence of the ac-susceptibility is that which is expected for a system with magnetic clusters. Below a freezing temperature (maximum of the $`\chi `$ versus $`T`$ curve) the large magnetic moments of $`L2_1`$-clusters are frozen along the magnetic easy axis. When the ageing time is increased, the peak associated with the occurrence of the frozen magnetic phase shifts towards higher temperatures and becomes broader. This is a consequence of the sensitivity of the susceptibility to the size and shape distribution of the magnetic clusters. For annealing times longer than $``$ 100 min the behaviour of the susceptibility is that expected for a ferromagnetic system: it is quite constant in both $`M`$\- and $`\beta `$-phases. In the temperature range where the $`MT`$ takes place, the increase in the ferromagnetic correlation among clusters explains the evolution of the shape of the hysteresis cycle of the $`MT`$ with an increase of ageing time. For short ageing times, $`\chi `$ increases at the $`MT`$ (forward transition) due to the tendency of the magnetic moments to freeze. In contrast, when the system becomes ferromagnetic, the behaviour of $`\chi `$ at the transition is dominated by the magnetic anisotropy which is larger in the $`M`$\- than in the $`\beta `$-phase.

From a general viewpoint, the results discussed above must be considered as evidence of the coupling between magnetic and structural degrees of freedom occurring at a mesoscopic scale, that is at the length scale of the magnetic/martensitic domains. The coupling also exists at a microscopic level (spin-phonon coupling) which is made clear by the estimated (through Clausius-Clapeyron equation) change of $`T_M`$ with $``$. This coupling is however very weak and is related to the small magnetic contribution to the entropy difference between the $`\beta `$\- and $`M`$-phases. The obtained value is comparable to that estimated for ferromagnetic Ni-Mn-Ga, but is considerably smaller than the value ($`0.1`$ K/kOe) reported for ferrous martensitic alloys Kakeshita85 . Consistently, it has been shown that application of a magnetic field does not significantly modify the phonon branches in Ni-Mn-Ga Manosa01 .

At this point, it is worth noticing that the strong influence of ageing on the entropy difference between $`\beta `$\- and $`M`$-phases may not be related to an increase in the magnetic contribution to the entropy. This effect must be simply ascribed to the fact that Mn-rich regions do not transform martensitically. The important change of $`\mathrm{\Delta }S`$ takes place within the ageing time interval during which the system becomes ferromagnetic. Notice that simultaneously the magnetization change at the transition reaches its saturation value. Therefore, any further evolution of the magnetic properties (as, for instance, revealed by ac-susceptibility measurements) of the system are originated by short-range reordering processes.

To conclude, the results presented prove the existence of both magnetoresistive and magnetoelastic properties in a non-stoichiometric Heusler Cu-Al-Mn alloy with a Mn content of 9 at$`\%`$ and with shape-memory properties. The occurrence of a martensitic transition has enabled the study of the effect of the crystallographic structural change on the above mentioned properties. The magnetoresistive behaviour has been shown to be different for the two structures. With regards to magnetoelastic coupling, it has been found that it mainly occurs at a mesoscopic level between martensitic domains and magnetic clusters.

###### Acknowledgements.

This work has received financial support from the CICyT (Spain), projects MAT2001-3251, MAT2000-0858 and from CIRIT (Catalonia), project 2001SGR00066. J.M. is supported by Direcció General de Recerca (Generalitat de Catalunya).


format: text

text:
# Effect of quantum confinement on exciton–phonon interactions

## Abstract

We investigate the homogeneous linewidth of localized type–I excitons in type–II GaAs/AlAs superlattices. These localizing centers represent the intermediate case between quasi–two–dimensional (Q2D) and quasi–zero–dimensional localizations. The temperature dependence of the homogeneous linewidth is obtained with high precision from microphotoluminescence spectra. We confirm the reduced interaction of the excitons with their environment with decreasing dimensionality except for the coupling to LO phonons. The low–temperature limit for the linewidth of these localized excitons is five times smaller than that of Q2D excitons. The coefficient of exciton–acoustic phonon interaction is $`56`$ times smaller than that of Q2D excitons. An enhancement of the average exciton–LO phonon interaction by localization is found in our sample. But this interaction is very sensitive to the detailed structure of the localizing centers.

The homogeneous linewidth of exciton luminescence is one of the most important features in excitonic dynamics in semiconductors, since it contains directly the information about the interactions between excitons and their environment. During the past two decades, the homogeneous linewidth of excitons in several kinds of quantum well and superlattice systems has been investigated extensively in both time and frequency domains. In the time domain, the excitonic dephasing time was measured from four–wave mixing (FWM), and then the homogeneous linewidth could be deduced.l571797 ; b349027 ; b406442 ; l681006 ; b592215  In the frequency domain, the linewidth was measured directly from photoluminescence,b335512 ; b4610193 ; b485241 ; jap871858  transmission, reflection or absorptionb5116785 ; apl611411  and Raman spectroscopyb501792 . By modeling of experimental data, extensive information about interactions between excitons and acoustic phonons, LO phonons, free carriers and other excitons has been deduced. In these investigations, excitons are quasi–two–dimensional (Q2D). That is, they can move freely in the wells or are localized weakly with a localization energy of several meV. On the other side, the homogeneous linewidth of quasi–zero–dimensional (Q0D) excitons confined in quantum dots, with localization energy of several hundreds of meV, has been studied by spatially resolved measurements.b5315743 ; l744043  The comparison of these two kinds of excitons provides information about the influence of quantum confinement on the interactions between excitons and their environment.

In this paper, we report investigations on homogeneous linewidth of single type–I localized excitons in GaAs/AlAs superlattices which have a global band alignment of type II. The localization energies of these centers are several tens of meV. Thus we can regard these localized excitons as intermediate in dimensionality between Q2D excitons and Q0D excitons. Furthermore, since the investigated centers are found in a small area (1 $`\mu `$m in diameter) of the same sample, we can rule out any artificial effects which come about when comparing different samples. This enables us to discuss the influence of localization on exciton-phonon interactions by comparing these centers, without disturbed by other artificial effects.

The localized excitons are studied by microphotoluminescence ($`\mu `$–PL). The spectral and spatial resolutions are sufficient to detect luminescence from individual localizing centers. Our experimental setup consists of a He flow cryostat with the sample mounted close to a thin window. This allows the use of a microscope objective to image (magnification 20$`\times `$) the excited spot on the sample onto a pinhole. The pinhole defines the spatial resolution. We use a 20–$`\mu `$m pinhole in the present experiments, which corresponds to 1–$`\mu `$m detected area on the sample’s surface. The pinhole is imaged onto the entrance slit of a 0.75–m focal length double grating spectrometer. We use a cooled CCD to record the spectra with a spectral resolution of 30 $`\mu `$eV. The sample is nonresonantly and globally excited by a He–Ne laser. The excitation intensity is about 1 $`\mathrm{W}/\mathrm{cm}^2`$ for all of the temperature–dependent measurements. During the measurement, the temperature of the sample is measured with a diode temperature sensor in good thermal contact. The temperature is stabilized by the He flow and heating to a fluctuation of less than 0.2 K. The measurements are performed in the range of $`780`$ K. We study two samples: (i)140 periods of GaAs(3 nm)/AlAs(2.8 nm) and (ii)140 periods of GaAs(2.3 nm)/AlAs(2.3 nm). Both samples have a type–II band alignment, i.e., the unperturbed conduction-band minimum is in the AlAs layer and the valence–band maximum is in the GaAs layer. Details about the growth and the interface properties of the samples have been reported previously.b571631  The two samples yield quite similar results concerning the exciton–phonon coupling. Thus we will only present data of sample (i).

Figure 1 reviews the luminescence properties of the sample at 20 K. The spatially integrated PL spectrum \[Fig. 1(a)\] is composed of a zero–phonon line at about 1.782 eV and phonon sidebands at the low–energy side. Luminescence intensity maps \[Fig. 1(b)\] show an inhomogeneous distribution of the emission intensity. We can find bright spots of about 1 $`\mu `$m in diameter. This size corresponds to the resolution of the objective in our $`\mu `$–PL system. The actual size of the bright spots was determined to be about $`250300`$ nm \[full width at half maximum (FWHM)\] by scanning near–field optical microscopy with resolution of 100 nm.

Figure 1(c) shows the $`\mu `$–PL spectrum of one of the bright spots. Different from Fig. 1(a), the spectrum is dominated by local emission from the bright spot. On the smooth background, spectrally narrow lines are superimposed.b5915862  Figure 1(d) provides a closer sight of these narrow lines. We have found that the narrow lines observed in Fig. 1(c) could be divided into two groups according to their different temperature behaviors.b5915862  For the lines on the low–energy side of the zero–phonon line, i.e., $`1.751.782`$ eV, the spectral weight shifts red with rising temperature, and their integrated intensity drops. These lines stem from localized type–II states. For the lines in the spectral range of the AlAs LO–phonon replica, i.e., below 1.74 eV, the spectral weight does not change significantly, and their intensities increase exponentially with temperature up to 50 K. We have proved that, although the global band alignment of this sample is type II, the layer thickness fluctuations give rise to local changes in the band alignment toward type I. Recombination of excitons localized in these type–I centers is the origin of the narrow lines in the energy range of $`1.691.74`$ eV. The population mechanism of these localized states has been proved to be electron tunneling from AlAs layers to GaAs layers.b5915862  Since some of these narrow lines are well separated in energy, we can resolve each of them without serious disturbance by adjacent lines. Thus we can analyze the luminescence from single localizing centers in the GaAs layers.

In order to investigate the temperature dependence of the linewidth of excitons localized in the type–I centers, $`\mu `$–PL spectra from several bright spots were measured in the temperature range of $`780`$ K. To check the possible spectral wandering during the integration, we measured the spectra with different integration times. The spectral shape and the linewidth keep unchanged as we vary the integration time in the range of $`50\mathrm{ms}30\mathrm{s}`$. We also measured the sequence of the spectra with an integration time of 50 ms and an interruption time of 1 s. We didn’t find any change of the peak position among these spectra. Thus, the spectral wandering of the sample can be neglected. The narrow lines were fitted by Lorentzian line shapes to obtain the linewidth (FWHM) in each temperature. An example of the narrow lines and the fitting curves is shown in Fig. 2. Recently, Besombes et al.b63155307  found that the line shape of the luminescence from strongly confined CdTe quantum dot deviates from Lorentzian shape with increasing the temperature. The whole spectrum is composed of a zero–phonon line and an additional acoustic–phonon sideband which results from lattice relaxation due to exciton-phonon coupled states. In InAs/GaAs system, both Lorentzianb65033313  and the non-Lorentzianl87157401  line shapes have been observed recently. For GaAs quantum dots, the PL spectrum has been shown to be of Lorentzian line shape.l763005  Also in our experiments on centers of intermediate confinement in GaAs/AlAs superlattices, we find no indications of additional phonon sidebands. The spectral line shape can be well fitted by Lorentzian function up to 80 K (see Fig. 2).

In general, the temperature dependent homogeneous linewidth of the exciton resonance is written as (see, e.g., Ref. b4211218, )

$$\mathrm{\Gamma }_{homo}(T)=\mathrm{\Gamma }_0+\gamma _{AC}T+\gamma _{LO}[\mathrm{exp}(\mathrm{}\omega _{LO}/k_BT)1]^1$$
(1)
where the term linear in temperature is due to exciton scattering with acoustic phonons, and the term nonlinear in temperature is due to interactions with LO phonons. The coefficients $`\gamma _{AC}`$ and $`\gamma _{LO}`$ represent the strength of the exciton–acoustic–phonon interaction and exciton–LO–phonon interaction, respectively. The first term in Eq. (1) is the low–temperature limit of the linewidth. The temperature dependence of the linewidth deduced from the Lorentzian fits was then fitted by Eq. (1), for several well–separated narrow lines. We show one of the fitting results in Fig. 3, as an example. The contributions to the linewidth from acoustic–phonon scattering and LO–phonon scattering are also shown in this figure (short–dashed and dotted lines, respectively). The dashed line represents low–temperature limit of the linewidth. For this narrow line, the fitting parameters are $`\mathrm{\Gamma }_0`$ = 42 $`\mu `$eV, $`\gamma _{AC}`$ = 1.2 $`\mu `$eV/K, and $`\gamma _{LO}`$= 85.7 meV. Since the fitting contains three free parameters, it is necessary to check the sensitivity of these parameters on fitting process. For this purpose, we vary each of these parameters from its optimal value, to check the variation of mean–square deviation between fitting curve and experimental data. We find that a 10 % deviation of $`\mathrm{\Gamma }_0`$, $`\gamma _{AC}`$, or $`\gamma _{LO}`$ from their optimal values corresponds increase of 60 %, 94 % or 580 % in the mean–square deviation, respectively. This result insures the safety of the fitting process. Another possible problem in extracting the parameters from the fitting is whether the temperature range ($`780`$ K) is large enough for an accurate determination of the parameters, especially for the $`\gamma _{LO}`$. In the investigations of the Q2D excitons, the linewidths were measured up to 150 K (Ref. b335512, ) or even room temperaturel681006 ; jap871858 ; b5116785 . In the present study, however, the luminescence of these localized excitons at a temperature above 80 K is too weak for an accurate determination of the linewidth. To check the possible errors introduced by this relatively small temperature range, we redraw the results obtained by the larger temperature range measurements,b335512 ; l681006 ; jap871858 ; b5116785  and fit by Eq. (1) all the data in the whole temperature range or only the data obtained below 80 K, respectively. We find that for all of the results checked, the difference of the parameter values obtained from the two fits using different ranges of data is less than 10 % for $`\gamma _{LO}`$ and less than 3 % for $`\gamma _{AC}`$. This proves the temperature range of $`780`$ K is enough for the determination of the linewidth–temperature curve with a satisfactory accuracy. Furthermore, the data density in the present study (about 50 data points in the temperature range of $`780`$ K) is high enough for an accurate fitting.

In Figure 4(a), we list the fitting results of $`\mathrm{\Gamma }_0`$, $`\gamma _{AC}`$, and $`\gamma _{LO}`$ of the narrow lines analyzed in the present study. In order to distinguish the narrow lines, we define $`\mathrm{\Delta }E`$ as the energy difference between the corresponding narrow line and the peak of zero–phonon line. We show an example of this definition for one narrow line in Fig. 1(c). Such a quantity is temperature independent. We note that $`\mathrm{\Delta }E`$ is close to, but not exactly, the localization energy of the corresponding center, since we cannot regard the zero–phonon–line as the mobility edge exactly. However, such a difference will not influence our discussions. For comparison, we list in Fig. 4(b) the available data of GaAs Q2D excitons deduced from FWM, photoluminescence, or other methods, as a function of the thickness of the GaAs layers. The results of GaAs bulk and superlattices are also listed in this figure, but not included in the calculations of average values, which are shown as dashed lines in Fig. 4(b). In the viewpoint of quantum confinement, we regard the localized excitons investigated here as the intermediate case between Q2D excitons in quantum wells and Q0D excitons in quantum dots. In the following, we will compare the $`\mathrm{\Gamma }_0`$, $`\gamma _{AC}`$, and $`\gamma _{LO}`$ of localized excitons obtained in the present study with that of the other two cases, to discuss the influence of confinement on exciton–phonon interactions in semiconductors.

At first, we discuss the low–temperature limit of the linewidth. We obtain the average value of $`\mathrm{\Gamma }_0`$ to be $`0.057(\pm 0.014)`$ meV for these localized excitons \[dashed line in Fig. 4(a)I\]. This value is five times smaller than the average value of Q2D excitons \[dashed line in Fig. 4(b)I\]. That indicates that the additional in–plane confinement in localizing centers reduces $`\mathrm{\Gamma }_0`$. In order to investigate the contribution of intercarrier scattering to this linewidth, we measured the excitation intensity dependence of the linewidth at 7 K. In the intensity range of $`110\mathrm{W}/\mathrm{cm}^2`$, the linewidth keeps unchanged, while in the range of $`105000\mathrm{W}/\mathrm{cm}^2`$ the linewidth increases slowly with a slope of $`0.010.02`$ $`\mu `$eV/($`\mathrm{W}/\mathrm{cm}^2`$). Due to the complicated population mechanism of these localizing centers (electron tunneling from AlAs layers to GaAs layers), we are not able to relate the excitation intensity to the actual carrier density in the sample. However, we can conclude that the excitation intensity used in the temperature–dependent measurement is quite low, and the intercarrier interaction can be neglected. The intrinsic lifetimes of excitons in quantum wires and quantum dots have been calculated to be several hundreds of picosecounds.l693393  However, recent experiment reveals a 16–ps intrinsic lifetime of excitons in GaAs quantum dots.s2932224  The linewidth obtained in the present study corresponds to a lifetime of 22 ps, which is consistent with this new finding.

Second, we discuss the influence of confinement on exciton–acoustic–phonon interaction. It has been found beforeb349027 ; b501792  and confirmed again recentlyjap871858  that $`\gamma _{AC}`$ of Q2D excitons is smaller than that in bulk GaAs \[Fig. 4(b)II\] . The average value of $`\gamma _{AC}`$ of the localized excitons obtained in the present investigation is $`1.14(\pm 0.24)`$ $`\mu `$eV/K, about $`56`$ times smaller than that of Q2D excitons \[dashed lines in Figs. 4(a)II and 4(b)II\]. The result suggests a reduction of exciton–acoustic–phonon interaction by localization. In quantum wires, such a reduction has been found by direct comparison of free and localized excitons in FWM measurements.b492993  Furthermore, in quantum dots, the homogeneous linewidth has been found to be almost constant up to 50 K.b5315743 ; l744043  These results suggest the extremely small $`\gamma _{AC}`$ for Q0D excitons. The weaker acoustic–phonon interaction with Q0D excitons than with Q2D excitons has also been confirmed in II–VI systems (see, for example, Ref. b63155307, ). The whole evolution discussed above, from bulk via Q2D excitons to localized excitons (this study) and to Q0D excitons, implies strongly that the interaction between exciton and acoustic phonon is steadily reduced by increasing confinement. Such a behavior is consistent with previous theoretical predictions.jlumi44347 ; b428947  In a confined system, the final state of phonon scattering is not always available due to the discrete energy level scheme. Thus, by increasing the quantum confinement, the appearance of the discrete energy levels induces a decrease of the acoustic phonon interaction. We note that when the confinement is so strong that the energy level space is larger then the thermal energy $`k_B`$T, a further increase of the confinement does not further reduce the interaction, since the level space has already been large enough for this bottleneck effect. In this regime, additional effects like lattice relaxationb63155307  can influence the dependence of the acoustic phonon interaction on the quantum dot size. Theoretical calculations revealed an increase, rather than decrease, of the acoustic–phonon coupling when further reducing the size of the quantum dots in this regime.b602638 ; b63155307

The parameters $`\mathrm{\Gamma }_0`$ and $`\gamma _{AC}`$ of the localized excitons obtained here are almost independent of $`\mathrm{\Delta }E`$ in the range of $`0.040.09`$ eV \[Figs. 4(a)I and 4(a)II\]. But, for $`\gamma _{LO}`$, we find a totally different behavior in the same energy range. The values of $`\gamma _{LO}`$ vary in the range of $`30140`$ meV, with no obvious systematic dependence on $`\mathrm{\Delta }E`$. The fluctuations of $`\mathrm{\Gamma }_0`$ and $`\gamma _{AC}`$, which are also obtained in the same fitting process, are all less than 25 %. We attribute the observed scattering in the homogeneous linewidth to an intrinsic feature of the exciton–LO–phonon coupling. In localizing centers, the energy level scheme of excitons is determined by the detailed structure of the center. Due to the monochromatic feature of the LO–phonon dispersion, the exciton–LO–phonon scattering rate depends sensitively on the level scheme. For the center in which the energy level scheme matches the LO–phonon energy well, a strong coupling is observed. In contrast to the LO phonons, the dispersion of the acoustic phonons distributes over a relative wide energy range. Thus, the exciton–acoustic–phonon scattering rate is less sensitive to the detailed structure of the localizing centers. In fact, we do not find the pronounced resonant behavior for the acoustic–phonon coupling \[see Fig. 4(a)II\].

In the strongly confined quantum dots, the explicit size and shape of the localizing potential determines the spatial extension and anisotropy of the electron–hole wave function as well as the electron–hole overlap. This has significant influence on the exciton–phonon interaction.l834654 ; b64241305  In strongly confined CdTe quantum dots, a mixing of the exciton and acoustic–phonon modes, which cannot be described by perturbation treatment, has been proposed.b63155307  That is, the exciton locally distorts the lattice of the dot. This lattice distortion is important for small quantum dots which sizes are comparable with the exciton Bohr radius. For example in II–VI and InAs/GaAs systems, an induced non-Lorentzian broadening have been observed.b63155307 ; l87157401  However, the localizing centers studied here are much larger than the Bohr radius. Thus the distortion is less important, and we do not observe strong deviations from a Lorentzian lineshape even at a temperature of 80 K. For the same reason, the influence of the potential size and shape on the electron–hole wave function is also less pronounced than that in strongly confined quantum dots. So we observe only small variations in the acoustic–phonon coupling strength among these localizing centers with different sizes and shapes.

Despite of the scattering behavior, we can still deduce the enhancement of the exciton–LO–phonon interaction in localized excitons with respect to Q2D excitons. The average value of $`\gamma _{LO}`$ is 71 meV, about five times larger than that of Q2D excitons \[Figs. 4(a)III and 4(b)III\]. The enhancement of exciton–LO–phonon interaction by localization induced by alloying fluctuations in alloy $`\mathrm{GaAs}_{1\mathrm{x}}\mathrm{P}_\mathrm{x}`$ has been found by resonant Raman spectroscopy.b547921  A similar enhancement was also found in GaN quantum wells.apl702882  In those investigations, the LO–phonon replica was used to detect the exciton–phonon interaction. In the present study, we detected luminescence from excitons localized by thickness fluctuations by $`\mu `$-PL. The agreements between different experimental methods as well as the different origins of localization confirm that the additional in–plane confinement on excitons enhances the exciton–LO–phonon interaction.

Up to now, the exciton–LO phonon interaction in quantum dots is still an open problem. In CdTe quantum dots, Besombes et al.b63155307  found that the exciton–LO–phonon scattering is not efficient up to 60 K, while Heitz and co–workersl834654 ; apl773746  observed enhanced exciton–LO–phonon interaction in InAs/GaAs self–organized quantum dots by measuring the phonon–assisted exciton transitions. The Huang–Rhys parameter was found to be five times larger than in bulk InAs. The enhancement was attributed to the quantum confinement and piezoelectric effect. Our result confirms qualitatively the latter finding. According to the extremely sensitive dependence of exciton–LO–phonon interaction on the detailed structures of localizing centers, we suggest that much care should be taken when comparing experimental results of this interaction in different quantum dot samples, since the detailed structure of the dots can be totally different, and this may influence the strength of the interaction to a great extent.

In summary, we have measured the homogeneous linewidth of type–I localized excitons in type–II GaAs/AlAs superlattices using $`\mu `$–PL. These excitons, with a localization energy of several tens of meV, can be regarded as intermediate case between Q2D excitons (free or weakly localized excitons in quantum wells) and Q0D excitons in quantum dots with confinement energy of several hundred meV. The low–temperature limit of the linewidth, $`\mathrm{\Gamma }_0`$, of these localized excitons is found to be five times smaller than that of Q2D excitons. We obtain a $`56`$ times smaller exciton–acoustic–phonon interaction coefficient, $`\gamma _{AC}`$, for the localized excitons with respect to that of Q2D excitons. Together with a comparison of exciton data in bulk and quantum dots, the reduction of exciton–acoustic–phonon interaction by confinement is confirmed. In contrast to the results on $`\mathrm{\Gamma }_0`$ and $`\gamma _{AC}`$, which are independent of localization energy, the coupling to LO phonons, $`\gamma _{LO}`$, shows strong variations. This finding is attributed to the strong influence of the energy level scheme on the exciton–LO–phonon coupling. In average, we confirm an enhancement of exciton–LO–phonon interaction by localization.

We gratefully acknowledge the growth of the excellent samples by W. Braun and K. Ploog, and helpful discussions with K. Kheng. This work was supported by the Deutsche Forschungsgemeinschaft.


text:
# Mixing by Non-linear Gravity Wave Breaking on a White Dwarf Surface

## 1 Introduction

Classical novae result from the ignition (and subsequent explosive thermonuclear burning) of a ($`10^4`$ m) layer of hydrogen-rich material that has accreted from a main sequence companion onto the surface of a white dwarf Truran (1982); Shara (1989); Starrfield (1989); Livio (1994). Observed abundances and explosion energies estimated from observations indicate that there must be significant mixing of the heavier material of the C/O or O/Ne white dwarf into the lighter accreted material (H/He). This mixing is critical because otherwise hydrogen burning would be too slow to reproduce observed nova characteristics in outburst. Further, without this mixing it is difficult to understand the observed abundances of intermediate-mass nuclei in the ejecta. Accordingly, nova models must incorporate a mechanism that will dredge up the heavier white dwarf material (Rosner et al., 2001, and references therein).

A recently proposed mixing mechanism is the breaking of non-linear resonant gravity waves at the C/O surface Rosner et al. (2000); Alexakis, Young, and Rosner (2002); Rosner et al. (2001). The gravity waves, driven by the “wind” of accreted material, can break, forming a layer of well-mixed material. This mixed layer may then be transported upward by convection, thereby enriching the accreted material. Because the length scale of this mixed layer may be very small (much smaller than the length scale of convection), previous precursor simulations have not captured this effect.

In this manuscript, we present a simulation of a wind-driven non-linear gravity wave breaking on the surface of a white dwarf. The simulation was performed with FLASH, a parallel, adaptive-mesh simulation code for the compressible, reactive flows found in many astrophysical environments Fryxell et al. (2000); Calder et al. (2000). This simulation is part of an ongoing study of this mechanism to assess its efficacy for mixing white dwarf material with envelope material.

## 2 Theory and Simulation Details

This work applies the theory of gravity wave generation (originally developed for the air over water interface) (Miles, 1957, and references therein) to a white dwarf with an accreted envelope. The control parameters for this model are $`G=g\delta /U_{\mathrm{max}}^2`$, which is the ratio of potential to kinetic energy of the wind ($`g`$ is the constant gravitational acceleration, $`U_{\mathrm{max}}`$ is the maximum wind speed, and $`\delta `$ is the characteristic length scale of the wind profile) and the ratio of densities at the interface, $`r=\rho _1/\rho _2`$. For this simulation, the domain was $`1.0\times 10^6`$ cm by $`1.0\times 10^6`$ cm, $`g=4.5\times 10^9`$ cm/s, and $`r=10`$, with a density of the white dwarf material at the interface of $`10^4`$ gm/cm3. The white dwarf and accreted envelope materials are modeled as simple $`\gamma =5/3`$ gases with the white dwarf material composed (by mass) of a 50/50 C/O mix and the accreted material composed of a 75/25 H/He mix.

The density and pressure profiles were obtained by integrating the equation of hydrostatic equilibrium

$$\frac{dp}{dy}=\rho g\widehat{𝐤},$$
(1)
which for the case of a compressible, gamma-law gas gives

$$\rho =\rho _i\left[1\left(\gamma 1\right)\frac{g\rho _iy}{P_0\gamma }\right]^{\frac{1}{\gamma 1}}$$
(2)
and

$$P=P_0\left[1\left(\gamma 1\right)\frac{g\rho _iy}{P_0\gamma }\right]^{\frac{\gamma }{\gamma 1}}.$$
(3)
Here $`P_0`$ is the pressure at the interface and $`\rho _i`$ is the density immediately above or below the interface. The wave was created by forcing the interface to be sinusoidal and perturbing the pressure and adding a velocity via a prescription similar to linear theory.

The wind profile is given by

$$U(y)=U_{\mathrm{max}}\left(1e^{y/\delta }\right)$$
(4)
with $`\delta =1.0\times 10^5`$ cm and $`U_{\mathrm{max}}=2\times 10^8\mathrm{cm}/\mathrm{s}`$. The boundary conditions were isothermal hydrostatic on the upper and lower boundaries and periodic on the sides. The simulation was performed with seven levels of adaptive mesh refinement, for an effective resolution of $`512\times 512`$ zones.

## 3 Results

Figure 1 shows the potential energy of the wave and the amount of mixed material during the course of the simulation. The left panel is a plot of potential energy vs. time showing the result from the simulation and the potential energy given by the linear theory as $`4.0\times 10^{26}\mathrm{exp}(140t)`$. The growth rate in this expression (140) comes from the control parameter $`G`$, and $`4.0\times 10^{26}`$ (erg) is the initial potential energy of the wave Alexakis, Young, and Rosner (2002). The two potential energies agree reasonably well until the wave breaks at about 0.04 s. The right panel is a plot of mixed C/O mass per unit area vs. time. The figure demonstrates the dramatic increase in the amount of mixed material that occurs above $`t=0.25`$, the point during the course of the simulation at which the wave begins to break.

Figure 2 shows gray scale images of the log of density with velocity vectors at early times in the simulation. The longest velocity arrows, at the top of the images, correspond to $`2.0\times 10^8`$ cm/s, the maximum wind speed. The left panel shows the initial conditions, and the right panel shows the configuration at $`t`$ = 0.015 s. Visible in the right panel is the development of Kelvin-Helmholtz instability on the upwind (left) side of the wave. Figure 3 shows similar gray scale images of the log of density with velocity vectors at later times in the simulation. By $`t=0.030`$ s (left panel), the wave has begun to break, and the images show vorticity downwind from the crest of the wave. At $`t=0.045`$ s (right panel), there is substantial mixing and obvious vorticity.

## 4 Conclusions

The simulation presented in this manuscript demonstrates the proposed breaking wind-driven non-linear gravity wave mixing mechanism. The results show the development of well-mixed zone just above the surface of the white dwarf. This simulation will be one part of a study of this mixing mechanism investigating effects of wind profiles and speeds. Complete details of the study will appear in Alexakis, et al. (2002).

The expectation is that the breaking of non-linear gravity waves on the surface of the white dwarf will lead to a thin well-mixed layer of material that may then be transported throughout the the envelope by convection. This study investigating the mixing mechanism should provide quantitative information about the mixing rate that will allow for the development of subgrid models that may be applied to multidimensional convection simulations to study the enrichment of the envelope. A preliminary simulation of this kind is also presented in this volume Dursi, et al. (2002).

This work is supported in part by the U.S. Department of Energy (DOE) under Grant No. B341495 to the Center for Astrophysical Thermonuclear Flashes at the University of Chicago. J. W. Truran acknowledges partial support from DOE grant DE-FG02-91ER40606. L. J. Dursi is supported by the Krell Institute CSGF. K. Olson acknowledges partial support from NASA grant NAS5-28524. M. Zingale acknowledges support from the Scientific Discovery through Advanced Computing (SciDAC) program of the DOE, grant number DE-FC02-01ER41176. Additional details about the project and information about requesting a copy of FLASH may be found at http://flash.uchicago.edu.


text:
# Low Temperature Electronic Transport and Electron Transfer through Organic Macromolecules

## Abstract

It is shown that at low temperatures and moderate electron dephasing the electron transmission function reveales a structure containing information about donor/acceptor sites effectively participating in the electron transfer processes and primary pathways of electrons tunneling through molecular bridges in macromolecules. This important information can be obtained as a result of analysis of experimental low temperature current-voltage characteristics for chosen molecules.

It has been estabilished that molecular electron transfer (ET) is essentially a combination of nuclear environment fluctuations and electron tunneling. Due to the large distances between donor and acceptor, ET is mostly provided by intervening molecular bridges, giving rise to a set of intermediate states for the electron tunneling 1 . The expression for the ET rate including both electronic and nuclear factors was first proposed by Marcus 2 ; 3 ; 4  and can be written as follows:

$$K_{ET}=K_{el}K_n\nu _n.$$
(1)
Here, $`K_{el}`$ is the electron transmission coefficient, $`K_n`$ is the nuclear transition or Franck-Condon factor, and $`\nu _n`$ is the effective nuclear vibration frequency.

The enormous size of biological molecules makes calculations of $`K_{el}`$ extremely difficult 5 ; 6 ; 7 ; 8 ; 9 ; 10 ; 11 . In this paper a simple approach is proposed which enables us to show that some intrinsic characteristics of the intramolecular ET such as pathways of tunneling electrons and donor/acceptor coupling to the bridge at different values of tunnel energy $`E`$ can be obtained in experiments on low-temperature electrical transport through the corresponding molecules. These data may be available under moderate electronic dephasing 12 ; 13 ; 14 . A similar approach could also be developed to analyze characteristics of intermolecular ET processes in extended biosystems.

Existing theories of electronic transport through macromolecules are mostly based on a very simple model simulating both electrodes by semi-infinite tight-binding chains attached to the ends of the molecule which is also simulated as a single tight-binding chain of sites 12 ; 13 ; 14 ; 15 ; 16 ; 17 . Here we adopt a resembling model. However, we take into account that donor and acceptor subsystems in realistic macromolecules are complex and include many sites providing effective coupling to the bridge, as well as to the attached electrodes.

To simplify further calculations we assume coherent electron tunneling through the macromolecule to be the predominant transport mechanism. Correspondingly we treat the electron transport through the molecule as a sequence of tunnelings between potential wells. Each well represents one atomic orbital of the donor-bridge-acceptor system. Within this approach, any atom is represented by a set of sites corresponding to its states. For further simplification we assume than intraatomic hopping integrals are smaller that interatomic ones, therefore we consider electron transfer between different sitesas its tunneling between different atoms.

Assume that sites ”i” are associated with those atoms of the bridge which have a strong coupling to the donor, and sites “$`j`$” are associated with atoms interacting with the acceptor. Then we can treat any site “$`i`$” as the origin of a chain along which an electron moves from donor to acceptor, and any site “$`j`$” as the end of the chain. Thus, we have a set of chains (pathways) for an electron moving along the bridge. In further analysis, we concentrate on the case when the chains weakly interact so we can consider them separately. We simulate the donor/acceptor subsystems as sets of semi-infinite tight-binding homogeneous chains. Each chain is attached to a site of donor/acceptor which can be effectively coupled to the bridge. As before we assume that the chains do not interact with each other. The adopted model does not enable us to carry out a proper treatment of coupling of electrodes to the molecule at metal-molecule-metal junctions which is a nontrivial problem 18 ; 19 ; 20 ; 21 . However, it seems reasonable to conjecture that in the molecules with complex donor and acceptor subsystems this process do not strongly depend on details of coupling of the electrodes to the donor/acceptor due to a comparatively large size and complicated structure of these subsystems.

In the following calculations we start from an effective tight-binding Hamiltonian for a single chain included into the bridge:

$$H_{eff}=H_0+H_1+H_𝒟+H_𝒜.$$
(2)
Here first two terms describe the chain itself, and their matrix elements between states $`|k>`$ and $`|l>`$ corresponding to the $`k`$th and $`l`$th sites are given by:

$$(H_0)_{kl}=\alpha _k\delta _{kl};(H_1)_{kl}=𝒱_{kl},$$
(3)
where $`𝒱_{kl}`$= 0 when $`k=l`$ and only states associated with valence electrons are considered. The diagonal matrix elements $`\alpha _k`$ are ionization energies of electrons at sites $`k`$, while $`𝒱_{kl}=𝒱_{lk}`$ includes both direct and exchange energy contributions for an electron to transfer between the $`k`$th and $`l`$th sites. Remainig terms represent self-energy corrections arising due to the coupling of the donor $`(H_𝒟)`$ and acceptor $`(H_𝒜)`$ to the bridge:

$`(H_𝒟)_{ii}=(\mathrm{\Sigma }_𝒟)_i={\displaystyle \underset{m}{}}{\displaystyle \frac{𝒟_{mi}^{2}}{Eϵ_m\sigma _m}},`$ (4)
$`(H_𝒟)_{jj}=(\mathrm{\Sigma }_𝒜)_j={\displaystyle \underset{n}{}}{\displaystyle \frac{𝒜_{jn}^{2}}{Eϵ_n\sigma _n}}.`$ (5)
Here, $`𝒟_{mi}`$ and $`𝒜_{jn}`$ are, respectively, coupling strengths between the $`m`$th donor site or the $`n`$th acceptor site and the $`i`$th or $`j`$th site of the bridge, and $`\sigma _{m,n}=\frac{1}{2}\left\{\theta _{m,n}i\sqrt{4\gamma _{m,n}^2\theta _{m,n}^2}\right\}`$ are the self-energy corrections of the semi-infinite chains attached to the corresponding sites 12 . The parameters $`\theta _{m,n}=Eϵ_{m,n},`$ where $`ϵ_{m,n}`$ and $`\gamma _{m,n}`$ are the ionization energies of electrons at the corresponding donor/acceptor sites, and the nearest-neighbor hopping integrals for the chains. Summation in Eqs. (4) and (5) is carried out over all donor/acceptor sites coupled to the bridge. Due to the presence of the self-energy corrections the eigenvalues of the effective Hamiltonian \[Eq. (2)\] include imaginary parts which give broadening of the bridge energy levels $`E_i.`$ The energy levels are broadened further if include scattering processes in the bridge which are not considered here.

An electric tunneling current $`I`$ flowing from donor to acceptor through the bridge in the presence of a small voltage $`V`$ applied across the bridge has the form 22

$$I=\frac{e}{\pi \mathrm{}}_{\mathrm{}}^{\mathrm{}}𝑑ET(E)[f\left(E\mu _1\right)f\left(E\mu _2\right)],$$
(6)
where $`f(E)`$ is the Fermi function, and the chemical potentials $`\mu _1`$ and $`\mu _2`$ are determined by the equilibrium Fermi energy of the bridge $`E_F`$ and the effective voltage $`V`$ across the bridge 16 , i.e.

$$\mu _1=E_F+\left(1\eta \right)eV;\mu _2=E_F\eta eV.$$

Here, the parameter $`\eta `$ characterizes how the voltage $`V`$ is divided between the two ends of the bridge. The electron transmission function is given by the formula

$$T(E)=2\underset{i,j}{}\mathrm{\Delta }_i|G_{ij}|^2\mathrm{\Delta }_j.$$
(7)
The summation in Eq. (7) is carried over the bridge states $`|i>`$ and $`|j>`$ therefore contributions from all possible pathways are contained here. The quantities $`\mathrm{\Delta }_{i,j}`$ are imaginary parts of the self-energy corrections $`_{𝒟,𝒜},`$ and $`G_{ij}`$ is the matrix element of the Green’s function corresponding to $`H_{eff}:`$

$$𝒢_{ij}=i\left|\left(EH_0H_1H_𝒟H_𝒜\right)^1\right|j.$$
(8)
We see that the dependence of the electron transmission function on energy is determined by the contributions from different donor/acceptor sites, as well as from the Green’s function matrix elements corresponding to different chains included into the bridge subsystem.

We can easily arrive at the first approximation for $`G_{ij}`$ neglecting the broadening of the bridge energy levels. Taking into account only nearest neighbors (NN) and next nearest neighbor (NNN) couplings, we can introduce a new notation for the nonzero hopping integrals, i.e., $`𝒱_{k,k+1}=𝒱_{k1,k}=\beta _k`$ and $`𝒱_{k,k+2}=𝒱_{k2,k}=\gamma _k.`$ Then we can calculate the Green’s function matrix elements for a given chain which makes a start at the site ”i” and ends at the site $`{}_{}{}^{\prime \prime }j_{}^{\prime \prime },`$ applying a diagram technique similar to that described in detail for the (NN) approximation 23 . We present our chain of sites as a graph $`\mathrm{\Gamma }_{ij},`$ using the following rules:

(i) Any site included in the chain corresponds to a graph vertex.

(ii) A term $`E\alpha _k`$ corresponds to a loop attached to the vertex “$`k`$”.

(iii) A hopping integral $`𝒱_{kl}`$ corresponds to an edge originating at the vertex “$`k`$” and ending at the vertex “$`l`$”.

As a result we obtain the graph shown below:

Now we consider all possible cycles for this graph and define the value of a cycle O as the product of values of the edges included into the cycle multiplied by the number of the edges connecting the firstand last sites with a negative sign. When the numbers of the edges associated with the paths from the first site to the last one and back differ, we multiply by the maximum number of edges. When the cycle includes only one vertex, its value is defined as $`E\alpha _k`$ where $`k`$ is the number of the vertex. Then we introduce a cyclic term which is a set of cycles including all sites of the system. The value of the cyclic term is equal to the product of the cyclic values of these cycles. The sum of all possible cyclic terms Or gives us the cyclic value $`\mathrm{\Theta }(\mathrm{\Gamma }_{ij})`$ of the graph $`\mathrm{\Gamma }_{ij}.`$ An expression for the matrix element $`G_{ij}`$ can be obtained in the form:

$$G_{ij}=\underset{\sigma }{}\frac{P_{ij}^\sigma }{\mathrm{\Theta }(\mathrm{\Gamma }_{ij})}.$$
(9)
Here, the summation is carried out over all possible pathways along the chosen chain. The value of $`P_{ij}^\sigma `$ equals the product of all edges along the pathway. When the pathway includes all vertices of the graph $`\mathrm{\Gamma }_{ij}`$ we have

$$P_{ij}=\underset{k=i}{\overset{j1}{}}\beta _k.$$
(10)
When the pathway does not include some vertices, this produces different values for the corresponding $`P_{ij}`$. For example, when the electron isbeing transferred from the site $`l`$ to the site $`l+2`$, passing the intermediate site $`l+1`$, the corresponding edge is equal to $`(E\alpha _{l+1})\gamma _l`$. Consequently, we obtain for the pathway which excludes this vertex:

$$P_{ij}=(E\alpha _{l+1})\gamma _l\underset{k=i,kl,l+1}{\overset{j1}{}}\beta _k.$$
(11)
Following a similar way we can get explicit expressions for all terms in the numerator of Eq. (9). The expression (9) can be significantly simplified if we assume that (NNN) couplings are small compared to (NN) couplings. Then we obtain the following approximation:

$$G_{ij}=\frac{1}{\mathrm{\Theta }_0(\mathrm{\Gamma }_{ij})}\{\underset{k=1}{\overset{j1}{}}\beta _k+\underset{l=i}{\overset{j2}{}}\gamma _l(E\alpha _{l+1})\underset{k=i;kl,l+1}{\overset{j1}{}}\}$$
(12)
where $`\mathrm{\Theta }_0(\mathrm{\Gamma }_{ij})`$ is the cyclic value of the graph calculated within the (NN) approximation.

The Green’s function matrix elements \[Eq. (9)\] have poles which can be detemined by solving the equations

$$\mathrm{\Theta }_{ij}(\mathrm{\Gamma })=0.$$
(13)
These poles correspond to energy states of a given chain. The set of these energy states characterizes the chosen chain (pathway) and different chains have different sets of energy states. When we take into account the broadening of the bridge energy levels, the singularities of the Green’s function matrix elements are converted to the peaks of finite height and width.

To apply the described method of calculation of the matrix elements $`G_{ij}`$ for a real macromolecule we have to determine which atomic orbitals have to be used to arrange a chain between the sites $`{}_{}{}^{\prime \prime }i_{}^{\prime \prime }`$ and $`{}_{}{}^{\prime \prime }j_{}^{\prime \prime }.`$ To choose suitable orbitals we can employ a rigorous quantum mechanical method of tunneling interatomic currents 24  which was alredy successfully used to separate out electron transfer pathways in long-range ET processes 25 . Data concerning the electronic structure of a considered molecule can be obtained within an extended Huckel 25  or some other approximation 26 . As a result we get the set of chains which represent the bridge, and we can use our result \[Eq. (9)\] in further analysis.

It follows from the obtained results \[Eqs. (7) and (9)\] that the electron transmission function contains important information concerning intramolecular electron transfer. Under low temperatures and moderate electronic dephasing $`T(E)`$ can exhibit series of peaks. Their location is determined with contributions from donor/acceptor sites predominating in the ET at a given interval of tunneling energy $`E,`$ and with the energy spectrum of a chain of sites connecting them. When $`E`$ takes on a value close to $`ϵ_{k,l}`$ or to one of the poles of the Green’s function the corresponding term in Eq. (7) can surpass all remaining terms. Thus, for different values of $`E`$ different donor/acceptor sites and different pathways can predominate in the electron transfer through the bridge.

The structure of $`T(E)`$ can be revealed in the low temperature current-voltage characteristics for the electronic transport through the molecular bridges \[Eq. (6)\]. As a result of numerical simulations of electrical transport through porphyrin-nitrobenzene macromolecules it was shown 27  that at low temperatures current-voltage curves can exhibit a steplike behavior which disappears as the temperature raises (see Fig. 1). This originates from a steplike character of Fermi distribution functions at low temperatures along with the ”comblike” structure of $`T(E).`$ At a given voltage $`V`$ the difference of Fermi functions in Eq. (6) takes on nonzero values only in the interior of a certain energy range including $`E_F.`$ Therefore the magnitude of the current $`I`$ at a given voltage is determined by the contributions of peaks of $`T(E)`$ located in this energy range. When the applied voltage increases, this enhances the width of the relevant energy interval. The tunneling current $`I`$ abruptly changes when an extra peak of the electron transmission function appears there. Widths of the plateaus are equal to the distances between adjacent peaks of $`T(E)`$ and magnitudes of sudden changes in the current correspond to the heights and shape of these peaks. At higher temperatures Fermi functions lose their steplike character and the plateous are washed out. Another reason for the structure of $`T(E)`$ (and $`IV`$ curves) to be eroded is the electronic phase-breaking effect which arises due to stochastic fluctuatious of the ion potential. We conjecture, however, that at low temperatures the electronic dephasing effects could be reduced so that the structure of $`T(E)`$ could be revealed.

On these grounds we believe that experiments on electronic transport through macromolecules at low temperatures $`(T1K)`$ could provide an additional information upon characteristics of the electron transfer. Namely, comparing the structure of the electron transmission function reconstructed on the basis of experimental $`IV`$ curves with that obtained as a result of calculations using Eqs. (7) and (9), we can make conclusions concerning actual primary pathways of electrons, as well as sites of the donor/acceptor subsystems involved in the ET process at different values of the tunneling energy. The model adopted in the present work provides us with the results suitable for a quantitative comparison with the results of proposed experiments even avoiding a proper and reliable calculation of the equilibrium Fermi energy of the bridge, and the effective voltage $`V`$. The lack of information about proper values of $`E_F`$ and $`V`$ produces an uncertainty in the location of the origin at the ”voltage” axis in the Fig. 1, so we cannot identify steps of the $`IV`$ curves separately. Nevertheless, changes in the values of $`E_F`$ and $`V`$ do not influence the electron transmission function, therefore the structure of series of the peaks remains fixed. This enables us to identify some series of peaks analyzing sequences of widths of successive steps of the current-voltage curves. Such analysis also can give reasonable estimations for the $`E_F`$ and $`V`$ for the chosen molecule 28 .

Acknowledgments: I thank G.M. Zimbovsky for help with the manuscript.


text:
# Computational complexity arising from degree correlations in networks

## Abstract

We apply a Bethe-Peierls approach to statistical-mechanics models defined on random networks of arbitrary degree distribution and arbitrary correlations between the degrees of neighboring vertices. Using the NP-hard optimization problem of finding minimal vertex covers on these graphs, we show that such correlations may lead to a qualitatively different solution structure as compared to uncorrelated networks. This results in a higher complexity of the network in a computational sense: Simple heuristic algorithms fail to find a minimal vertex cover in the highly correlated case, whereas uncorrelated networks seem to be simple from the point of view of combinatorial optimization.

The last few years have seen a great advance in the study of complex networks albert02 , where the term complex refers to the existence of one or more of the following properties: small world effect watts98 , power-law degree distribution albert02 , and more recently also correlations redner ; pastor01 ; berg ; newman02 . On the other hand, if we focus on the solution of a given task on top of these networks, the term complex is better associated with the time required to solve it, i.e. with its computational complexity garey79 . In this context, a problem is complex if its algorithmic solution time is growing exponentially in the network size. At the core of complex optimization problems one finds the NP-hard class garey79 , where NP stands for non-deterministic polynomial time.

In the case of uncorrelated networks with power law degree distributions we can take profit of the existence of hubs to solve different problems, like destroying the giant component percolation , preventing epidemic outbreaks spreading , and searching searching . The extremely inhomogeneous structure of uncorrelated networks can also be exploited to approximate or even to solve instances of NP-hard problems using heuristic algorithms running in polynomial time. However, the influence of properties like degree correlations or clustering is not clear yet. Recent studies of percolation newman02  and disease spreading boguna02  have shown that degree correlations can quantitatively change, e.g., the transition threshold but qualitatively the results are similar to those obtained for uncorrelated networks.

This changes drastically if we consider hard optimization tasks defined over correlated networks. In this work we study the influence of degree correlations on the computational complexity, and in a more general perspective the relation between the topology of complex networks and the computational complexity of hard problems defined on top of them. For this purpose we generalize the Bethe-Peierls approach to statistical-mechanics models defined on networks with an arbitrary degree distribution and arbitrary degree correlations of adjacent nodes.

The approach is applied to characterize the minimal vertex covers on these graphs. We have chosen this problem for two reasons: It belongs to the basic NP-hard optimization problems over graphs garey79 , and has found applications in monitoring Internet traffic breitbart  and denial of service attack prevention park . Our analytical results are later compared with an approximate solution obtained using a heuristic algorithm. This heuristic fails to find minimal vertex covers in the strongly correlated case, whereas networks with low correlations seem to be simple from the point of view of combinatorial optimization. Within our analytical approach, this change of behavior is associated with replica symmetry breaking (RSB).

Consider the set of undirected graphs with $`N`$ vertices and arbitrary degree distribution $`p_d`$. Following a randomly chosen edge, we will find a vertex of degree $`d+1`$ with probability $`q_d=(d+1)p_{d+1}/c`$, with $`c`$ denoting the average degree. The number $`d`$ of additional edges will be called excess degree. We further assume correlations between adjacent vertices: The probability that a randomly chosen edge connects two vertices of excess degrees $`d,d^{}`$ is given by $`(2\delta _{d,d^{}})e_{dd^{}}`$. The conditional probability that a vertex of excess degree $`d`$ is reached following any edge coming from a vertex of excess degree $`d^{}`$,

$$p(d|d^{})=e_{dd^{}}/q_d^{},$$
(1)
thus explicitly depends on both $`d`$ and $`d^{}`$. Consistency with the degree distribution requires $`_{d^{}=0}^{\mathrm{}}e_{dd^{}}=q_d`$, and $`e_{dd^{}}`$ has to be symmetric. For uncorrelated graphs $`e_{dd^{}}=q_dq_d^{}`$ factorizes. The strongest positive correlations are reached for $`e_{dd^{}}=q_d\delta _{dd^{}}`$ where only vertices of equal degrees are connected.

Let us now consider a general statistical-mechanics model with discrete degrees of freedom defined on vertices, and interactions defined on edges. We use a lattice gas model described by the Hamiltonian

$$\beta H=\underset{i<j}{}J_{ij}w(x_i,x_j)+\mu \underset{i}{}x_i$$
(2)
defined for any microscopic configuration $`x_i=0,1,i=1,\mathrm{},N`$. $`J`$ is the adjacency matrix with entries $`J_{ij}=1`$ if vertices $`i`$ and $`j`$ are adjacent, and $`J_{ij}=0`$ else. The inverse temperature is denoted by $`\beta `$, the chemical potential by $`\mu `$. The interactions $`w(x_i,x_j)`$ are arbitrary, thus including also the case of the ferromagnetic Ising model, $`w(x_i,x_j)=(2x_i1)(2x_j1)`$. The only disorder present in Eq. (2) is given by the edges $`J_{ij}`$. Generalizations to disordered interactions, as present e.g. in spin-glasses, random local fields or non-binary discrete variables are evident. For clarity of the presentation, we restrict ourselves to the simple model given above.

Since the graphs are locally tree-like, the model can be solved by the iterative Bethe-Peierls scheme which becomes exact if only one pure state is present. The free energy can be expressed in terms of simple effective-field distributions acting on vertices of given degree. In the case of multiple pure states this has to be generalized to the cavity approach, see e.g. MePa  for the example of a spin-glass on a Bethe lattice of constant vertex degree. Alternatively, one can apply the replica approach. The simple Bethe-Peierls solution corresponds to the assumption of replica symmetry (RS), whereas the full cavity approach is able to handle also the case of RSB.

Take any edge $`(i,j)`$, i.e. $`J_{ij}=1`$. Let us introduce $`Z_x^{(i|j)}`$ as the partition function of the subtree rooted in $`i`$, with deleted edge $`(i,j)`$, and with $`x_i`$ fixed to the value $`x`$. This partition function can be calculated iteratively,

$`Z_0^{(i|j)}`$ $`=`$ $`{\displaystyle \underset{kj|J_{ik}=1}{}}\left(e^{w(0,0)}Z_0^{(k|i)}+e^{w(0,1)}Z_1^{(k|i)}\right)`$ (3)
$`Z_1^{(i|j)}`$ $`=`$ $`e^\mu {\displaystyle \underset{kj|J_{ik}=1}{}}\left(e^{w(1,0)}Z_0^{(k|i)}+e^{w(1,1)}Z_1^{(k|i)}\right)`$
The effective fields $`h_{(i|j)}=\mathrm{ln}(Z_1^{(i|j)}/Z_0^{(i|j)})`$ are thus determined by the iterative description

$$h_{(i|j)}=\mu +\underset{kj|J_{ik}=1}{}u(h_{(k|i)}).$$
(4)
where $`u(h_{(k|i)})`$ is the effective field induced by $`x_k`$ on site $`i`$, and is given by

$$u(h_{(k|i)})=\mathrm{ln}\left(\frac{e^{w(1,0)}+e^{w(1,1)+h_{(k|i)}}}{e^{w(0,0)}+e^{w(0,1)+h_{(k|i)}}}\right).$$
(5)
The free energy of the system can be written as

$$\beta Nf=\underset{i}{}f_i+\underset{i<j}{}J_{ij}(f_{ij}f_if_j)$$
(6)
where the link contribution equals

$$f_{ij}=\mathrm{ln}\left(\underset{x_i,x_j}{}e^{w(x_i,x_j)+h_{(i|j)}x_i+h_{(j|i)}x_j}\right),$$
(7)
whereas the site contribution

$$f_i=\mathrm{ln}\left(\underset{x_i}{}e^{h_ix_i}\right)$$
(8)
depends on the cavity field

$$h_i=\mu +\underset{k|J_{ik}=1}{}u(h_{(k|i)})$$
(9)
resulting from the influence of all neighbors on vertex $`i`$.

Let us now assume, that the model has only one pure state, which corresponds to the assumption of RS. In this case, the iterative procedure given by Eq. (4) converges to well-defined distributions $`P_d(h)`$ of effective fields $`h_{(i|j)}`$ restricted to vertices of excess degree $`d`$. They are determined by the self-consistency equation

$`P_d(h)`$ $`=`$ $`{\displaystyle _{\mathrm{}}^{\mathrm{}}}{\displaystyle \underset{l=1}{\overset{d}{}}}\left(dh_l{\displaystyle \underset{d_l=0}{\overset{\mathrm{}}{}}}p(d_l|d)P_{d_l}(h_l)\right)`$ (10)
$`\times \delta \left(h\mu {\displaystyle \underset{l=1}{\overset{d}{}}}u(h_l)\right).`$
Please note that, in contrast to the uncorrelated case, we do need the field distributions for all possible excess degrees. In the uncorrelated case, the average of these distributions over $`q_d`$ is sufficient. We may also introduce the analogous distributions $`\stackrel{~}{P}_d(h)`$ of cavity fields $`h_i`$ for vertices of given degree $`d`$ (note: here the full degree is relevant), they can be calculated from $`P_d(h)`$ by

$`\stackrel{~}{P}_d(h)`$ $`=`$ $`{\displaystyle _{\mathrm{}}^{\mathrm{}}}{\displaystyle \underset{l=1}{\overset{d}{}}}\left(dh_l{\displaystyle \underset{d_l=0}{\overset{\mathrm{}}{}}}p(d_l|d1)P_{d_l}(h_l)\right)`$ (11)
$`\times \delta \left(h\mu {\displaystyle \underset{l=1}{\overset{d}{}}}u(h_l)\right).`$
Replacing the sum over vertices and edges in Eq. (6) by the corresponding averages over $`P_d(h)`$ resp. $`\stackrel{~}{P}_d(h)`$, we finally find the free-energy density

$`\beta f`$ $`=`$ $`{\displaystyle \underset{d=0}{\overset{\mathrm{}}{}}}(d1)p_d{\displaystyle _{\mathrm{}}^{\mathrm{}}}𝑑h\stackrel{~}{P}_d(h)\mathrm{ln}\left(1+e^h\right)`$
$`+{\displaystyle \frac{c}{2}}{\displaystyle \underset{d,d^{}=0}{\overset{\mathrm{}}{}}}e_{dd^{}}{\displaystyle _{\mathrm{}}^{\mathrm{}}}𝑑h𝑑h^{}P_d(h)P_d^{}(h^{})`$
$`\mathrm{ln}(e^{w(0,0)}+e^{w(1,0)+h}+e^{w(0,1)+h^{}}+e^{w(1,1)+h+h^{}})`$
The simplest application of this approach is given by the ferromagnetic Ising model. If we look to the ground-states, i.e. to the limit $`\beta \mathrm{}`$, we find that, as to be expected, the global magnetization is determined by the size of the giant component. The existence of a ferromagnetic phase at low temperature is thus related to percolation. The latter was already analyzed in Ref. newman02 .

Another application is given by the vertex cover (VC) problem. It belongs to the basic NP-hard optimization problems garey79  and, therefore, it is expected to require a solution time which is growing exponentially with the graph size. Let us be more precise. Given a graph with vertices $`i\{1,\mathrm{},N\}`$ and edges $`\{(i,j)|1i<jN,J_{ij}=1\}`$, a vertex cover $`V`$ is a subset of vertices, $`V\{1,\mathrm{},N\}`$, such that at least one end-vertex of every edge is contained in $`V`$. So no edge $`(i,j)`$ is allowed to exist with $`iV`$ and $`jV`$. Of course, the set of all vertices forms a trivial VC. The hard optimization problem consists in finding the minimal VC.

Using the hard-sphere lattice-gas representation introduced in WeHa , where $`x_i=1`$ if $`iV`$, and $`x_i=0`$ if $`iV`$, the VC condition can be rewritten as

$$\underset{(i,j)|J_{ij}=1}{}(1x_ix_j)=1$$
(13)
which fits into the above framework by setting $`e^{w(x_i,x_j)}=1x_ix_j`$. The chemical potential can be used to fix the cardinality of the VC, minimal ones are obtained in the limit $`\mu \mathrm{}`$. They correspond to maximal packings in the lattice-gas picture. To perform the limit $`\mu \mathrm{}`$ all fields have to be rescaled as $`h=\mu z`$ WeHa . We obtain

$`P_d(z)`$ $`=`$ $`{\displaystyle _{\mathrm{}}^{\mathrm{}}}{\displaystyle \underset{l=1}{\overset{d}{}}}\left(dz_l{\displaystyle \underset{d_l=0}{\overset{\mathrm{}}{}}}p(d_l|d)P_{d_l}(z_l)\right)`$ (14)
$`\times \delta \left(z1{\displaystyle \underset{l=1}{\overset{d}{}}}\mathrm{max}(0,z_l)\right),`$
which is solved by $`P_d(z)=_{l=1}^+\mathrm{}\rho _l^{(d)}\delta (z+l)`$. This ansatz allows for integer valued fields only, we find a simple relation including only the $`\rho _1^{(d)}`$:

$$\rho _1^{(d)}=\left[\underset{d_1=0}{\overset{\mathrm{}}{}}p(d_1|d)(1\rho _1^{(d_1)})\right]^d$$
(15)
All other $`\rho _l^{(d)}`$ follow easily. The expression inside the parenthesis can be understood as the average probability $`\pi _d`$ that an edge entering a vertex of degree $`d+1`$ carries a constraint, i.e. that it is not yet covered by the neighboring vertex. It thus fulfills the condition

$$\pi _d=\underset{d_1=0}{\overset{\mathrm{}}{}}p(d_1|d)(1\pi _{d_1})^{d_1}.$$
(16)
Having in mind that, due to the limit $`\mu \mathrm{}`$, every vertex with positive $`z`$ is fixed to $`x=1`$, every one with negative $`z`$ has $`x=0`$, we can immediately read off the fraction of vertices belonging to a minimal VC,

$$x_c=1\underset{d=0}{\overset{\mathrm{}}{}}p_d(1\pi _{d1})^{d1}\left(1+\frac{d2}{2}\pi _{d1}\right).$$
(17)
Remember that the last expressions are related to the validity of RS, i.e. to the existence of a single connected cluster of minimal VCs in configuration space. As observed in WeHa , RS is related to the local stability of this solution. In presence of RSB, Eq. (16) has no stable solution. Since it has to be solved by numerical iteration in the general case, an instability prevents the program from convergence and thus provides a precise tool to detect RSB without any RSB calculation.

To see how this works out, we concentrate on networks having equal degree distributions but different correlation properties. We restrict our attention to scale-free graphs with $`p_dd^\gamma `$ for $`d=1,\mathrm{},\mathrm{}`$, with $`\gamma >2`$. For vertex cover, interesting effects are expected to appear for positive correlations, or assortative networks. We therefore consider

$$e_{dd^{}}=q_d\left[r\delta _{d,d^{}}+(1r)q_d^{}\right].$$
(18)
For $`0r1`$, this expression linearly interpolates between uncorrelated ($`r=0`$) and fully assorted ($`r=1`$) networks. Please note that this network is percolated for all $`\gamma `$ as soon as $`r>0`$. In Fig. 1 we show the resulting size of the minimal VCs for different values of $`\gamma `$ as a function of $`r`$. The RS solution breaks at a certain value of $`r`$ that depends on $`\gamma `$. There, the solution-space structure changes drastically, from being unstructured, or RS, in the low correlated case to being clustered, or RSB, for sufficiently high correlations.

To check the consequences of this transition for heuristic optimization algorithms, we have numerically generated scale-free networks with correlations (18), and applied a generalization of the leaf-removal heuristic of Bauer and Golinelli bauer01 . For the special case of correlations given by Eq. (18), the network can be generated using a modification of the Molloy-Reed algorithm MoRe  for random graphs of arbitrary degree distribution. First, each node $`i`$ is assigned a degree $`d_i`$ with probability $`p_{d_i}`$. Then we create a set $``$ containing $`d_i`$ copies of each node $`i`$. Finally, pairs of nodes are connected according to following rule: (i) Select a node $`i`$ in $``$ at random. (ii) With probability $`r`$, select a second node $`j`$ with $`d_j=d_i`$ in $``$ at random; otherwise select an arbitrary node $`j`$ at random. (iii) Connect $`i`$ and $`j`$ and delete both from $``$. This is repeated until $``$ is empty.

Once the network is generated, we construct a VC using a generalized leaf-removal algorithm weigt02  defined as follows: Select a vertex of minimal current degree from the network and cover all its neighbors. The considered vertices and all incident edges are removed from the network. This step is iterated untill the full network is removed. If, for some graph, this algorithm stops without having ever chosen vertices of current degree $`d2`$, the constructed VC is minimal bauer01 . Overestimations may appear if the algorithm is forced to select also vertices of higher degree $`d2`$, where the error can be at most $`d1`$. Thus, summing $`(d1)(1\delta _{d,0})/N`$ over all iteration steps, we get an upper bound $`\mathrm{\Delta }x`$ on the total error made in estimating $`x_c`$ using the above heuristic algorithm. If $`\mathrm{\Delta }x`$ goes to zero in the large-$`N`$ limit, the algorithm has consequently constructed an almost minimal VC.

In Fig. 1 we show the size of VCs found by generalized leaf removal as a function of $`r`$. Up to the RSB point the numerical solutions are close to the analytical values, up to finite-size corrections resulting mainly from a degree cutoff $`d_{max}N^{1/\gamma }`$. Beyond the RSB point we still have a numerical estimate but we cannot be sure that it is optimal. In the inset of Fig. 1 the upper bound on the error is displayed. In the RS region we have $`\mathrm{\Delta }x0`$ for $`N1`$ and, therefore, the heuristic algorithm asymptotically yields the exact value $`x_c`$. However, in the highly correlated region we find a finite $`\mathrm{\Delta }x`$ at any network size, thus the heuristic algorithm fails to find almost minimal vertex covers. Moreover, the point where $`\mathrm{\Delta }x`$ becomes different from zero coincides with the RSB point.

It is interesting to know in which phase realistic networks are. VCs have found applications in monitoring the Internet traffic breitbart  and in denial of service attack prevention park . The analysis of Internet maps has revealed negative (dissasortative) correlations at the autonomous system level pastor01 . Negative correlations are actually common in technological and biological networks newman02 . Hence, the generalized leaf-removal heuristic should output almost optimal VCs in linear time. On the contrary, social networks exhibit positive (assortative) correlations newman02 . In this case VCs can be used to monitor social relations between pairs of individuals but, because of the existence of positive correlations, simple heuristic algorithms may fail to produce near optimal solutions.

To summarize, we have generalized the Bethe-Peierls approach to random networks with degree correlations and analyzed the VC problem as a prototype optimization problem defined over graphs. We have found that uncorrelated power-law networks are simple from the point of view of combinatorial optimization, inhomogeneities of neighboring vertices can be exploited. The introduction of sufficiently large degree correlations leads to RSB and thus to a failure of simple heuristic algorithms. For constructing optimal solutions, complete algorithms including e.g. backtracking have to be used. These, however, result in general in exponential solution times, and thus in a higher algorithmic complexity. Our results point out that optimization problems in many technological and biological networks can be simple due to the strong degree inhomogeneities and negative correlations present on them.

###### Acknowledgements.

We acknowledge fruitful discussions with M. Leone, A. Vespignani and R. Zecchina.