text:
# Electronic Transport through YBa2Cu3O7-δ Grain Boundary Interfaces between 4.2 K and 300 K

## Abstract

The current-induced dissipation in YBa2Cu3O7-δ grain boundary tunnel junctions has been measured between 4.2 K and 300 K. It is found that the resistance of 45 (100)/(110) junctions decreases linearly by a factor of four when their temperature is increased from 100 K to 300 K. At the superconducting transition temperature $`T_\mathrm{c}`$ the grain boundary resistance of the normal state and of the superconducting state extrapolate to the same value.

Soon after the discovery of superconductivity in La2-xBaxCuO4  it was realized that the normal-state properties of the high-$`T_\mathrm{c}`$ cuprates differ significantly from those of the low-$`T_\mathrm{c}`$ superconductors. The most prominent and still unexplained normal state properties of the high-$`T_\mathrm{c}`$ compounds are the linear temperature dependence of the in-plane resistivity in the optimally doped compounds , spin and charge inhomogeneities , and the pseudogap present in the underdoped regime .

The complex electronic behavior of the high-$`T_\mathrm{c}`$ materials is also reflected in the characteristics of their interfaces. For $`T<T_\mathrm{c}`$, for example, large-angle grain boundaries are excellent tunnelling barriers , a feature unknown from conventional superconductors. Several mechanisms responsible for the electronic transport properties of grain boundaries and other interfaces have been identified for the temperature range below $`T_\mathrm{c}`$, such as the d$`_{x^2y^2}`$ dominated order parameter symmetry, structural effects and space charge layers . In fact, grain boundaries operated below $`T_\mathrm{c}`$ have been used to reveal fundamental properties of the high-$`T_\mathrm{c}`$ cuprates, such as the existence of Cooper pairs  or the unconventional order parameter symmetry .

As is the case for $`T<T_\mathrm{c}`$, we expect valuable information on the interfaces of the cuprates, or even on the bulk materials themselves, to be contained in the tunneling characteristics of grain boundaries for $`T>T_\mathrm{c}`$. However, only very few data are available for this temperature regime . This lack of data arises from the difficulty to measure the grain boundary current-voltage (I($`V_{\mathrm{gb}}`$)) characteristics free of voltages produced by the bridges that are needed to contact the interfaces. By utilizing a difference-technique to subtract these unwanted voltages, we have now succeeded in measuring the current-voltage characteristics of YBa2Cu3O7-δ grain boundaries for $`4.2\text{K}T300\text{K}`$. These studies provide evidence that the normal state resistance of 45 (100)/(110) tilt grain boundaries is not influenced by the onset of the superconducting transition. The experiments furthermore reveal that for $`T>T_\mathrm{c}`$ the resistance of these grain boundaries decreases linearly with increasing temperature.

The experiments were performed with bicrystalline YBa2Cu3O7-δ films grown by pulsed laser deposition at 760C in $`0.25`$ mbar of O2 to a typical thickness of 40–50 nm. The SrTiO3-substrates contained symmetric and (100)/(110)-asymmetric 45 -tilt grain boundaries, specified to within 1. After deposition, the samples were cooled during two hours in an O2-pressure of $`400`$ mbar. Because the deposition geometry was tuned to optimize the homogeneity of the samples, $`T_\mathrm{c}`$ varied by less than 0.5 K across the wafers, while the thickness variations were smaller than 10 %. For the four-point measurements, gold contacts were structured photolithographically before the YBa2Cu3O7-δ films were patterned by etching in a H3PO4 solution, or by dry etching with an Ar ion-beam. The YBa2Cu3O7-δ bridges had widths and lengths between the voltage probes of $``$ 2.5 $`\mu `$m and 30 $`\mu `$m, respectively (see Fig. 1), and their $`T_\mathrm{c}`$ varied less than 200 mK. In addition, several samples were patterned into Wheatstone-bridge configurations, of which each arm consists of a meander line containing 23 straight sections with widths of $``$ 6 $`\mu `$m. Two of the meander lines were patterned such that each of the 23 sections crosses the boundary once.

The technique to measure precisely the grain boundary resistance $`R_{\mathrm{gb}}`$ below $`T_\mathrm{c}`$ is based on measuring the I($`V_{\mathrm{gb}}`$) characteristics of the Josephson junctions while applying a microwave field large enough to suppress their critical current. For this purpose, we used low frequency microwaves ($`<`$ 10 GHz) to avoid inaccuracies resulting from the Shapiro steps.

All measurements were conducted in shielded rooms using motor-driven dip-sticks. The samples were cooled at 4.2 K in liquid Helium, at 77 K in liquid nitrogen, and at other temperatures via He oder N2 vapor and thermal conduction across the sample holder. Immersion into the cooling liquids did not affect the characteristics, providing evidence that self-heating is insignificant under the experimental conditions.

For $`T>T_\mathrm{c}`$ the transport properties of the interfaces cannot be measured with a straightforward four-point technique because the voltage $`V_{\mathrm{meas}}`$ generated by a bridge straddling the boundary is composed of two contributions: a) the grain boundary voltage $`V_{\mathrm{gb}}`$ and b) the voltages $`V_{\mathrm{g1}}`$ and $`V_{\mathrm{g2}}`$ caused by the resistive grains inside the bridge. This problem can be solved by subtracting from $`V_{\mathrm{meas}}`$ the grain contributions $`V_{\mathrm{g1}}`$ and $`V_{\mathrm{g2}}`$. For this, $`V_{\mathrm{g1}}`$ and $`V_{\mathrm{g2}}`$ are obtained by approximating them with the voltages $`V_{\mathrm{g1}}^{}`$ and $`V_{\mathrm{g2}}^{}`$ generated by independent intragrain bridges that are measured simultaneously. Not surprisingly, test measurements revealed variations of $`V_{\mathrm{g1}}^{}`$ and $`V_{\mathrm{g2}}^{}`$ of several percent as a function of the bridge location. Therefore, the intragrain bridges were placed as closely as possible to the boundary bridge (see Fig. 1). The desired measurement accuracy also determined the optimum size of the bridges. On the one hand the bridges were patterned to be as small as possible so that they could be located close to the grain boundary. On the other hand, their size was chosen to be large enough for the patterning-induced scatter of the bridge aspect ratios to be insignificant. In samples optimized this way, the voltages across the two bridges $`V_{\mathrm{g1}}^{}`$ and $`V_{\mathrm{g2}}^{}`$ were identical to within 1 % at 100 K and to 0.3 % for $`T>200`$ K. Therefore the grain boundary voltage equals with the same accuracy $`V_{\mathrm{gb}}V_{\mathrm{meas}}(V_{\mathrm{g1}}^{}+V_{\mathrm{g2}}^{})/2`$.

As introduced by Mathur et al.  for studies of La1-xCaxMnO3 bicrystals, and similar to the work described in Ref. , we also patterned several samples into Wheatstone-bridges. This approach has the advantage that the generated voltages are huge and that $`V_{\mathrm{g1}}^{}`$ and $`V_{\mathrm{g2}}^{}`$ are subtracted from $`V_{\mathrm{meas}}`$ already during the measurement. Having performed measurements with such bridges, we checked their balance by photolithographically cutting the Wheatstone-configurations, to individually assess the resistances of the meander lines. Because these studies revealed $`T`$ dependent balancing errors of 1–10 %, we preferred the three-bridge approach for measurements of the $`R_{\mathrm{gb}}(T)`$) characteristics, while Wheatstone-bridges were chosen when a large signal was required, as was the case for some studies of $`I(V)`$ characteristics.

To gain insight into the electronic transport mechanisms, we measured the $`I(V)`$ characteristics between 4.2 K and 300 K by using Wheatstone-bridges. For $`T>T_\mathrm{c}`$ the characteristics are non-linear on this large voltage scale, in particular below 150 K (see Fig. 2). Analyzing the non-linearity, we first note that the microstructure of the grain boundaries is inhomogeneous down to the unit cell level . The $`I(V)`$ curves are therefore generated by large numbers of microstructurally different channels connected in parallel. The behavior of the averaged junctions is consistent with the one of a back-to-back Schottky contact as predicted by the band bending model . If tentatively described by a Simmons-fit , the averaged $`I(V)`$ characteristics correspond to heights and widths of hypothetical effective junction barriers of $``$ 100 meV and 1–2 nm. Thus, these data suggest unusually small barrier heights. They are, in particular, much smaller than the energy scale of grain boundary built-in potentials measured by electron-holography of $``$ 2 eV .

At small voltages, $`V_{\mathrm{gb}}`$ 10 mV, the non-linearity of the $`I(V)`$ curves is insignificant, and they are characterized by their ohmic resistance, which can be obtained by selecting an appropriate voltage ($`V_{\mathrm{cr}}`$) or current ($`I_{\mathrm{cr}}`$) criterion $`V_{\mathrm{cr}}=I_{\mathrm{cr}}\times R_{\mathrm{gb}}<`$ 10 mV. For $`T<T_\mathrm{c}`$, due to the microwave irradiation, the $`I(V)`$ characteristics are linear, too, even over a larger voltage range. The only exception is given by 45 boundaries at $`T<40`$ K, which will be discussed below. Applying the three-bridge technique, we were therefore able to deduce the $`R_{\mathrm{gb}}(T)`$ dependence of a given grain boundary by measuring in one temperature sweep from 4.2 K to 300 K just one boundary bridge plus the two intragrain reference bridges. Hereby one voltage or current criterion, typically $`I_{\mathrm{cr}}=100\mu `$A, is used for all temperatures.

The $`R(T)`$-dependence of a 45 symmetric grain boundary measured this way is presented in Fig. 3. The resistance-area product $`R_{\mathrm{gb}}A`$ (77 K) $``$ $`1\times 10^8`$ $`\mathrm{\Omega }`$cm2 compares well with literature values . Three remarkable features are displayed by the temperature dependence of the resistance. First, the resistance reaches a maximum at $``$ 30 K (see also Fig. 4). This maximum is not displayed by symmetric 24 and 36 boundaries measured as reference samples. Second, around $`T_\mathrm{c}`$ the resistance develops a distinct peak structure. This peak, which is not associated with the microwave irradiation, was found to vary from sample to sample, and to change with time over weeks. Except for this peak, the grain boundary resistance of the normal state and of the superconducting state extrapolate within the measurement accuracy of $`2\%`$ to the same value $`R_{\mathrm{gb}}`$ ($`T_\mathrm{c}`$) = 10.7 $`\mathrm{\Omega }`$. This temperature dependence is characteristic for all samples we have studied, except for few that displayed at $`T_\mathrm{c}`$ a resistance step of $`5\%`$. As this step was found to be not reproducible and to increase for a given sample over weeks, we consider it to be an artifact resulting from chemical reactions or diffusion. Third, between 100 K and 300 K the resistance decreases with increasing temperature by a factor of four. The temperature dependence is hereby remarkably linear, with a barely noticeable positive curvature (Fig. 3, upper inset). This is in striking contrast to the linear resistance increase of the adjacent YBa2Cu3O7-δ grains (Fig. 3, lower inset). In the following, we will discuss the three features and their implications for the understanding of the interfaces, beginning with the superconducting regime.

At low temperatures, the resistance of 45 grain boundaries changes non-monotonically with $`T`$, and even depends on the value of the tunneling current (see Fig. 4). This maximum and its current bias dependence are consistent with the formation of a zero-bias anomaly by the faceted 45 junctions. At small bias currents $`V_{\mathrm{gb}}`$ is of the order 1 mV, and therefore sufficiently small to probe the temperature dependent shape of the zero-bias conductance peak in the superconducting energy gap. This zero-bias conductance peak is usually presumed to arise from Andreev bound states generated by the d$`_{x^2y^2}`$ order parameter symmetry of YBa2Cu3O7-δ .

The peak structure at 89 K provides evidence that at $`T_\mathrm{c}`$ the resistance of the grain boundary bridge is much larger than the resistance of the intragrain-bridges. The observed height of the peak shown in Fig. 3 is consistent with a $`T_\mathrm{c}`$ difference of the bridges of 150 mK. The shape of the peak is influenced by the differences in the fluctuations of the intragrain bridges and the one containing the tunnel junction. Unfortunately we cannot extract data on this effect from the peak, as the spatial distribution of $`T_\mathrm{c}`$ is unknown.

The grain boundary resistance values above and below $`T_\mathrm{c}`$ extrapolate to the same value at the transition temperature. The transition into the superconducting state is therefore found not to affect the dissipation in the junctions at $`T_\mathrm{c}`$. This observation cannot be reconciled with the predictions for surface charging effects in the hole superconductivity scenario introduced by Hirsch and coworkers . In this scenario, charges are separated in the superconducting state in the vicinity of a grain boundary, but not above $`T_\mathrm{c}`$. This charge separation is expected to change the electronic properties of the interfaces when passing through $`T_\mathrm{c}`$ . The grain boundary data do not show this effect.

For $`T>T_c`$, the grain boundary resistance decreases linearly. This behavior cannot be explained by elastic tunneling from the vicinity of a well-defined Fermi surface across standard Schottky contacts with $`T`$-independent barrier heights and widths . There, a nonlinear and weak $`R_{\mathrm{gb}}(T)`$ characteristic is expected for the parameter range derived above. A more detailed elastic tunneling model has to include $`T`$ dependent Schottky barriers and to consider the spatial inhomogeneity of the interface. To significantly alter the weak and nonlinear $`T`$ dependence both effects demand the thermal energy $`k_BT`$ to exceed the energy scales of the barrier for large parts of the interface areas. Hereby a very special energy profile of the barriers is required to yield the linear $`R_{\mathrm{gb}}(T)`$ behavior over a large temperature range. This situation seems unphysical to us. We note that inelastic tunneling via resonant impurity channels  cannot account for the linear temperature characteristic either. Considering the problems which standard semiconductor models have in describing the observed interface characteristics, we now turn our attention to the effect of the electronic correlations present in the cuprates.

As pointed out by Miller and Freericks  electronic correlations strongly affect the barrier properties of high-$`T_c`$ interfaces in the superconducting state. This is also true in the normal conducting state. As suggested by the band bending model , the barrier is an underdoped region, with the possible formation of antiferromagnetic fluctuations or even of local magnetic moments. The resultant magnetic scattering is expected to increase the interface resistance with decreasing temperature. Since the density of the magnetic moments is not known, the $`R_{\mathrm{gb}}(T)`$ dependence cannot be determined yet. Besides correlations in the barrier, electronic correlations in the bulk may also influence the transport through the interface. What would be expected if, due to the grain boundary charging, the YBa2Cu3O7-δ layers adjacent to the grain boundary were underdoped, characterized by a pseudo gap and, possibly, by preformed pairs? The spectral densities of the electronic excitations on both sides of the interface control the tunneling process. It is only for the 45 grain boundaries under consideration, that the nodal and antinodal directions of the pseudogap are directly coupled, and the tunneling current is most strongly reduced. As the pseudo gap diminishes for increasing temperature, the boundary conductivity is expected to increase. Here we note that the $`R_{\mathrm{gb}}(T)`$ curve is characterized by one temperature scale, $`T_0`$, the extrapolated zero-resistance intercept. In the samples investigated, $`T_0350`$ K, which corresponds to the pseudo gap energy scale.

In summary, it is found that the resistance of 45 -tilt YBa2Cu3O7-δ grain boundary junctions decreases linearly with increasing temperature for $`T>T_\mathrm{c}`$. In the vicinity of $`T_\mathrm{c}`$ the grain boundary resistance and thus the dissipation are indistinguishable for $`T<T_\mathrm{c}`$ and $`T>T_\mathrm{c}`$, in contradiction to predictions of the theory of hole superconductivity. The $`I(V)`$ characteristics between 4.2 K and 300 K are nonlinear, suggesting a tunneling barrier with an effective height and width of $``$ 0.1 eV and $``$ 1–2 nm, respectively. The linear $`R_{\mathrm{gb}}(T)`$ dependence is proposed to be caused by correlation controlled tunneling.

The authors gratefully acknowledge helpful discussions with Y. Barash, H. Bielefeldt, M. Blamire, U. Eckern, P. J. Hirschfeld, C. Laschinger, J. Ransley, A. Rosch, D. G. Schlom, and M. Siegel. This work was supported by the DFG through the SFB 484 and by the BMBF via project 13N6918A.

Schneider et al., Figure 1

Schneider et al., Figure 2

Schneider et al., Figure 3

Schneider et al., Figure 4


format: text

text:
# A model of transport nonuniversality in thick-film resistors

\[

## Abstract

We propose a model of transport in thick-film resistors which naturally explains the observed nonuniversal values of the conductance exponent $`t`$ extracted in the vicinity of the percolation transition. Essential ingredients of the model are the segregated microstructure typical of thick-film resistors and tunneling between the conducting grains. Nonuniversality sets in as consequence of wide distribution of interparticle tunneling distances.

PACS numbers: 72.60.+g, 64.60.Fr, 72.80.Tm

\]

Thick-film resistors (TFRs) are glass-conductor composites based on RuO2 (but also Bi2Ru2O7, Pb2Ru2O6, and IrO2) grains mixed and fired with glass powders. Besides the widespread use of TFRs in pressure and force sensor applications, their transport properties are of great interest also for basic research. The percolating nature of transport in TFRs has been reported since long time and now it is well documented. As shown in Fig.1 where we reports a selection of previously published data on different TFRs, the conducting phase concentration $`x`$ dependence of the conductance $`G`$ of TFRs follows a percolating-like power-law equation of the form:

$$G=G_0(xx_c)^t,$$
(1)
where $`G_0`$ is a prefactor, $`x_c`$ is the critical concentration below which $`G`$ vanishes and $`t`$ is the transport critical exponent. The values of $`G_0`$, $`x_c`$, and $`t`$ which best fit the experimental data are reported in the inset of Fig.1.

According to the standard theory of transport in isotropic percolating systems, $`G_0`$ and $`x_c`$ depend on microscopic details such as the microstructure and the mean value of the junction resistances connecting two neighbouring conducting sites, while, unless the microscopic resistances have a diverging distribution function (see below), the critical exponent $`t`$ is universal, i. e., it depends only upon the lattice dimensionality $`D`$. For $`D=3`$, random resistor network calculations predict $`t=t_02.0`$, in agreement with various granular metal systems, or other disordered compounds.

As it is clear from Fig.1, TFRs have values of $`t`$ ranging from its universal limit $`t2.0`$ (filled squares, Ref.) up to very high values like $`t5.0`$ (filled diamonds, Ref.) or even higher. Despite of their clear percolating behavior, TFRs do not fulfil therefore the hypothesis of universality common to many systems and instead belong to a different, quite vast, class of materials which display nonuniversal transport behavior, that is a regime where the transport exponent $`t`$ depends on microscopic details (microstructure etc.). Typical examples of nonuniversal systems are carbon-black–polymer composites, and materials constituted by insulating regions embedded in a continuous conducting phase.

Despite that TFRs have been historically among the first materials for which transport nonuniversality has been reported, the microscopic origin of their universality breakdown has not been specifically addressed so far. In this letter we show that the cross-over between universality and nonuniversality reported in Fig.1 can be explained within a single model whose basic features are the peculiar microstructure of TFRs and the tunneling processes between conducting grains.

Before describing our model for TFRs, let us first recall the mathematical requisites for universality breakdown in random resistor networks. Consider a regular lattice of sites and assign to each neighbouring couple of sites a bond which has finite conductance $`g`$ with probability $`p`$ and zero conductance with probability $`1p`$. The resulting conductance distribution function is then:

$$\rho (g)=ph(g)+(1p)\delta (g),$$
(2)
where $`\delta (g)`$ is the Dirac delta-function and $`h(g)`$ is the distribution function of the finite bond conductances. For well behaved distribution functions $`h(g)`$, conductivity is universal and follows Eq.(1) with $`t=t_0`$. Instead, as first shown by Kogut and Straley, if $`h(g)`$ has a power law divergence for small $`g`$ of the form:

$$\underset{g0}{lim}h(g)g^\alpha ,$$
(3)
and $`\alpha `$ is larger than a critical value $`\alpha _c`$, then transport is no longer universal. Renormalization group analysis predicts in fact that

$$t=\{\begin{array}{cc}t_0\hfill & \text{if}\nu +1/(1\alpha )<t_0\hfill \\ \nu +1/(1\alpha )\hfill & \text{if}\nu +1/(1\alpha )>t_0\hfill \end{array},$$
(4)
where $`\nu 0.88`$ is the correlation-length exponent for a three dimensional lattice. By using $`t_02.0`$ we obtain therefore $`\alpha _c0.107`$. Equations (3) and (4) have been shown to arise from a system of insulating spheres embedded in a continuous conducting material (swiss-cheese model), and from a tunneling-percolation model with highly fluctuating tunneling distances. Here we show that Eq.(3) \[and consequently Eq.(4)\] arises naturally from a simple representation of TFRs in terms of their microstructure and elemental transport processes.

Let us start by considering the highly non-homogeneous microstructure typical of TFRs. These systems are constituted by a mixture of large glassy particles (typically with size $`L`$ of order $`1`$-$`3`$ $`\mu `$m) and small conducting grains of size $`\mathrm{\Phi }`$ varying between $`10`$ nm up to $`200`$ nm. In this situation, the small metallic grains tend to occupy the narrow regions between the much larger insulating zones leading to a filamentary distribution of the conducting phase. A classical model to describe such a segregation effect was proposed already in the 1970’s by Pike. This model treats the glassy particles as cubes of size $`L\mathrm{\Phi }`$ whose edges are occupied by chains of adjacent metallic spheres of diameter $`\mathrm{\Phi }`$. Such chains define channels (bonds) which form a cubic lattice spanning the whole sample. Let us assume for the moment that a bond has probability $`p`$ of being occupied by a fixed number $`n+1`$ of spheres and probability $`1p`$ of being empty. To each couple of adjacent spheres we assign an inter-sphere conductance $`\sigma _i`$ ($`i=1,\mathrm{},n`$). A random resistor network can be therefore defined as in Eq.(2) where $`h_n(g)`$ is the distribution function of the total channel conductance $`g`$ of $`n`$ conductances $`\sigma _i`$ in series:

$$g^1=\underset{i=1}{\overset{n}{}}\frac{1}{\sigma _i}.$$
(5)
The high values of piezoresistance (i.e., the strain sensitivity of transport) typical of TFRs, and the low temperature dependence of transport strongly indicate that the main contribution to the overall resistance stems from tunneling processes between neighbouring metallic grains. Hence, if the centers of two neighbouring metallic spheres are separated by a distance $`r`$, then the intergrain tunneling conductance $`\sigma `$ is approximatively of the form:

$$\sigma =\sigma (r)\sigma _0e^{2(r\mathrm{\Phi })/\xi },$$
(6)
where $`\sigma _0`$ is a constant which we set equal to the unity, $`\xi 1/\sqrt{V}`$ is the tunneling factor and $`V`$ is the intergrain barrier potential. Let us make the quite general assumption that the centers of the spheres are set randomly along the channel, so that the distances $`r`$ change according to the distribution function $`P(r)`$ for a set of impenetrable spheres arranged randomly in a quasi one dimensional channel. By following Ref., $`P(r)`$ can be calculated exactly and its explicit expression is:

$$P(r)=\frac{1}{a_n\mathrm{\Phi }}e^{(r\mathrm{\Phi })/(a_n\mathrm{\Phi })}\mathrm{\Theta }(r\mathrm{\Phi })$$
(7)
where $`a_n=(1+L/n\mathrm{\Phi })\mathrm{\Phi }/2`$ is the mean inter-sphere distance and $`\mathrm{\Theta }`$ is the step function. By combining Eq.(6) with Eq.(7) the distribution $`f(\sigma )`$ of the inter-sphere conductances is then:

$$f(\sigma )=_0^{\mathrm{}}𝑑rP(r)\delta [\sigma \sigma (r)]=(1\alpha _n)\sigma ^{\alpha _n},$$
(8)
where

$$\alpha _n=1\frac{\xi /2}{a_n\mathrm{\Phi }}.$$
(9)
To obtain the distribution function $`h_n(g)`$ of the occupied channels, we first note that Eq.(5) implyies that $`g`$ is dominated by the minimum inter-sphere conductance $`\sigma _{\mathrm{min}}`$ among the set of $`n`$ conductances in series. Hence the small-$`g`$ limit of $`h_n(g)`$ is just the distribution function $`\stackrel{~}{f}`$ of $`\sigma _{\mathrm{min}}`$:

$$\stackrel{~}{f}(\sigma _{\mathrm{min}})=nf(\sigma _{\mathrm{min}})\left[1_{\sigma _{\mathrm{min}}}^1𝑑\sigma _{\mathrm{min}}f(\sigma _{\mathrm{min}})\right]^{n1},$$
(10)
which, from Eq.(8) and by setting $`g\sigma _{\mathrm{min}}`$, leads to:

$$\underset{g0}{lim}h_n(g)n(1\alpha _n)g^{\alpha _n}.$$
(11)
The conducting bond distribution function behaves therefore as Eq.(3) so that for $`\alpha _n>\alpha _c0.107`$ transport universality breaks down and $`t>t_0`$. For $`L=1`$ $`\mu `$m, $`\mathrm{\Phi }=10`$ nm, and $`\xi =1`$ nm this is achieved already for $`n<90`$, i. e., slightly less then the maximum number $`L/\mathrm{\Phi }=100`$ of spheres which can be accommodated inside a channel.

Our model of universality breakdown in TFRs can be readily generalized to describe more realistic situations. For example, the number $`n`$ of spheres inside the occupied channels can vary according to a given distribution. It is also straightforward to rewrite Eq.(7) in order to describe cases in which the diameter $`\mathrm{\Phi }`$ of the spheres is not fixed, or to let the size of the insulating grains to change by assigning a distribution function for $`L`$. It is then possible to have different scenarios and, more importantly, to obtain a crossover from transport universality ($`t=t_0`$) to nonuniversality ($`t>t_0`$) within the same framework. This reminds the experimental situation reported for TFRs and summarized in Fig.1.

Let us now comment on the capability of other existing theories to describe transport universality breakdown in TFRs. At a first glance, the swiss-cheese model of Ref. is a natural candidate since the large values of $`L/\mathrm{\Phi }`$ typical of many TFRs may lead to an effective continuous conducting phase filling the voids between the large glassy grains. However, there are examples in which nonuniversality has been reported for TFRs with $`L/\mathrm{\Phi }`$ only of order $`510`$, a value probably too small to be compatible with the swiss-cheese picture. Even more problematic are the cases for which $`t=t_02.0`$ has been measured for TFRs with $`L/\mathrm{\Phi }100`$ (see for example Fig.1), while the swiss-cheese model would have predicted $`t>t_0`$. Regarding instead the model proposed by Balberg, for which transport is dominated by random tunneling processes in a percolating network, it is important to point out that it was defined by using a phenomenological distribution function for the nearest-neighbour particle distances very similar to our Eq.(7). Balberg argured that such form of $`P(r)`$ is a reasonable compromise between the distribution function of spheres randomly placed in three dimensions, and the effect of interactions between the conducting and insulating phases. Instead we have shown that Eq.(7), and consequently the power-law divergence of $`h(g)`$, is a straightforward outcome of the quasi one-dimensional geometry of the conducting channels in the segregation model of TFRs. Despite that our model has been formulated specifically for TRFs, nevertheless it could be applied also to other segregated disordered compounds for which tunneling is the main mechanism of transport and nonuniversality has been reported.

In conclusion, we have proposed a simple tunneling-percolation model capable of describing the observed transport universality breakdown in TFRs. Essential ingredients of the theory are the segregated structure, modelled by quasi-one dimensional channels occupied randomly by the conducting particles, and intergrain tunneling taking place within the channels.

This work is part of TOPNANO 21 project n.5557.2.


text:
# SUSY GUT Models of 𝜈 Mass and 𝜇→𝑒𝛾

## 1 SUSY GUT See-saw Models of $`\nu `$ Mass

There are three great advantages that SUSY GUT see-saw models have over other models of neutrino mass and mixing. The first advantage is that neutrino masses are explained in essentially the same way that the other known fermion masses are explained, namely the existence of both left- and right-handed varieties of fermion that have a Dirac coupling to each other through a Yukawa interaction with a Higgs doublet. GUTs, with the notable exception of $`SU(5)`$, generally require that right-handed neutrinos exist. Thus they generally predict that neutrino masses should exist. Moreover, GUTs are very well motivated on other grounds. By contrast, in many other approaches to explaining neutrino masses they arise from some postulated new physics that is not strongly motivated on other grounds, such as R-parity-violating interactions in SUSY, singly charged scalar fields (as in the Zee mechanism), or triplet Higgs .

The second advantage of the SUSY GUT see-saw mechanism is that it explains in an order of magnitude sense the size of neutrino masses. In the see-saw formula, the Dirac mass of a neutrino is expected to be related to the mass of the corresponding up-type of quark. So, if one estimates the mass $`m_3`$ of the third family neutrino, by using the top quark mass for the Dirac neutrino mass and the GUT scale for the Majorana mass of the right-handed neutrino, one obtains

$$m_3m_{Dirac}^2/M_R(174GeV)^2/(2\times 10^{16}GeV)1.5\times 10^{12}GeV.$$
(1)
On the other hand, if neutrino masses are hierarchical, then the mass $`m_3`$ can be determined from the mass-squared splitting observed in atmospheric neutrino oscillations:

$$m_3\sqrt{m_3^2m_2^2}=\sqrt{\delta m_{atm}^2}6\times 10^{11}GeV.$$
(2)
The agreement is gratifyingly good. By contrast, most other mechanisms for generating neutrino mass rely upon new physics that involves parameters that are very poorly if at all constrained by experiment. For example, we have little idea what the magnitudes of the coefficients of R-parity violating terms would be, or the couplings of the Zee scalar, or the couplings and masses of Higgs triplets. The neutrino masses in such models can come out to be anywhere within a huge range.

The third advantage of SUSY GUT see-saw models of neutrino mass is that they can be very predictive, since unified symmetries can precisely relate the Dirac neutrino mass matrix, the mass matrices of the quarks, and the mass matrix of the charged leptons.

Although SUSY GUT see-saw models are extremely well motivated, they face what at first sight appears to be a serious problem in explaining why some neutrino mixing angles are large. Since quarks and leptons are unified in GUTs, there was the expectation for a long time that leptonic mixing angles would turn out to be small like the CKM angles. Indeed, that was true in GUT models constructed before experimental evidence began to emerge in the mid-1990s that the atmospheric angle might actually be large.

More specifically, many models assume hierarchical mass matrices, i.e. the elements get smaller going up and to the left of diagonal elements. For instance, consider the classic idea proposed independently by Wilczek-Zee, Weinberg, and Fritzsch in 1977 , whereby the down quark mass matrix has the form

$$M_D=\left(\begin{array}{cc}0& ϵ\\ ϵ& 1\end{array}\right)m,ϵ1.$$
(3)
This gives $`\mathrm{tan}\theta ϵ`$ and $`m_d/m_sϵ^2`$, where $`\theta `$ is the angle required to diagonalize $`M_D`$. In fact, it gives exactly $`\mathrm{tan}\theta =\sqrt{m_d/m_s}`$. If we ignore the contribution to the Cabibbo angle coming from diagonalizing the up quark mass matrix $`M_U`$, then one reproduces the famous empirical relation $`\mathrm{tan}\theta _C\sqrt{m_d/m_s}`$.

There are countless models in the literature where the small mass ratios ($`m_em_\mu m_\tau `$, $`m_um_cm_t`$, and $`m_dm_sm_b`$) and small CKM angles are explained by assuming that the mass matrices of the charged leptons ($`M_L`$), up quarks ($`M_U`$), and down quarks ($`M_D`$) are hierarchical.

The trouble is that in unification based on $`SO(10)`$ or larger groups the Dirac mass matrix of the neutrinos ($`M_N`$) is related to these other mass matrices, and thus might be expected to be hierarchical also. However, if $`M_N`$ is hierarchical, then by the see-saw formula $`M_\nu =M_NM_R^1M_N^T`$ it would typically be the case that the light neutrino mass matrix $`M_\nu `$ would be hierarchical too, leading to small leptonic mixing angles.

For instance, if the Dirac mass matrix of the neutrinos has the hierarchical form

$$M_N\left(\begin{array}{ccc}ϵ^2& ϵϵ^{}& ϵ^{}\\ ϵϵ^{}& ϵ^2& ϵ\\ ϵ^{}& ϵ& 1\end{array}\right)m,$$
(4)
where $`ϵ`$ and $`ϵ^{}`$ are small parameters, and the Majorana mass matrix $`M_R`$ has a “random” form, i.e. all its elements being of the same order, then the matrix $`M_\nu `$ is easily seen to have the same hierarchical form as $`M_N`$, and the leptonic mixing angles will be small.

There is now strong evidence that two of the three leptonic mixing angles are large: $`|U_{e2}|\mathrm{sin}\theta _{sol}0.5`$, and $`|U_{\mu 3}|\mathrm{sin}\theta _{atm}0.7`$. There are two main ways to explain these large angles within the framework of the see-saw mechanism in SUSY GUTs. These have been called “Single Right-Handed Neutrino Dominance” or SRND , and “Lopsided Models” .

(A) Single Right-Handed Neutrino Dominance

(i) $`\nu _{R2}`$ or $`\nu _{R1}`$ Dominates in See-saw Formula.

This idea is best explained through a simple example. Suppose that the matrices $`M_N`$ and $`M_R`$ have the following forms:

$$M_N\left(\begin{array}{ccc}\eta & \eta & \eta \\ \eta & ϵ& ϵ\\ \eta & ϵ& 1\end{array}\right)m,M_R^1\left(\begin{array}{ccc}0& 0& 0\\ 0& 1& 0\\ 0& 0& 0\end{array}\right)m_R,$$
(5)
where $`\eta ϵ1`$. Then it is clear that when these are inserted into the see-saw formula the matrix $`M_R^1`$ approximately projects out the second column of $`M_N`$, giving

$$M_\nu \left(\begin{array}{c}\eta \\ ϵ\\ ϵ\end{array}\right)\left(\begin{array}{ccc}\eta & ϵ& ϵ\end{array}\right)m^2/m_R\left(\begin{array}{ccc}\eta ^2& \eta ϵ& \eta ϵ\\ \eta ϵ& ϵ^2& ϵ^2\\ \eta ϵ& ϵ^2& ϵ^2\end{array}\right)m^2/m_R.$$
(6)
One sees that as a result of the similar magnitudes of the 22 and 32 elements of $`M_N`$, the elements of the 23 block of $`M_\nu `$ are all of the same order, leading to a large atmospheric angle. Moreover, the fact that $`M_R^1`$ has nearly the form of a projection matrix of rank one, leads the 23 block of $`M_\nu `$ to have a nearly “factorized form”, i.e. to be approximately a rank one matrix. This gives the hierarchy between the neutrino masses $`m_2m_3`$.

(i) $`\nu _{R3}`$ Dominates in See-saw Formula.

Here the idea is that the matrices $`M_N`$ and $`M_R`$ have the forms:

$$M_N\left(\begin{array}{ccc}& & ϵ\\ & & 1\\ & & 1\end{array}\right)m,M_R^1\left(\begin{array}{ccc}0& 0& 0\\ 0& 0& 0\\ 0& 0& 1\end{array}\right)m_R,$$
(7)
where the dashes represent much smaller entries. This yields

$$M_\nu \left(\begin{array}{ccc}ϵ^2& ϵ& ϵ\\ ϵ& 1& 1\\ ϵ& 1& 1\end{array}\right)m^2/m_R.$$
(8)
Here $`M_N`$ has a large off-diagonal element. This kind of model is similar in some respects to the lopsided models about to be discussed. Indeed, in lopsided models, if one goes to a basis where $`M_L`$ is diagonal, the matrix $`M_N`$ acquires the form shown in Eq. (7).

(B) Lopsided Models

The basic insight behind lopsided models is that $`SU(5)`$ unification unifies left-handed down quarks with right-handed charged leptons, and right-handed down quarks with left-handed charged leptons. Consequently, what are related to the CKM angles by $`SU(5)`$ unification are the mixing angles of the right-handed leptons, not the observed neutrino angles, which of course are mixings of the left-handed leptons. Thus there is really no discrepancy in the fact that all of the CKM angles are small and some of the neutrino angles are large.

Suppose that the mass matrix of the charged leptons $`M_L`$ is highly left-right asymmetric. In particular, suppose that its 23 element is of the same order as its 33 element, while its 32 element is smaller by a factor $`ϵ1`$. Now, $`SU(5)`$ relates $`M_L`$ to the transpose of $`M_D`$. Indeed, in minimal $`SU(5)`$ one has $`M_L=M_D^T`$. This is just due to the fact, already mentioned, that $`SU(5)`$ relates down quarks and charged leptons of opposite chirality. Therefore, in the situation we are considering, if there is an $`SU(5)`$ symmetry governing the structure of the mass matrices, the matrix $`M_L`$ will have its 32 element and 33 element being of the same order, while its 23 element will be a factor of order $`ϵ`$ smaller. Notice, then, that the 23 elements, which are the ones that control the mixing of the second and third families of left-handed fermions (the 32 elements control the corresponding right-handed mixing) will be large for the charged leptons but small for the down quarks. Thus, $`U_{\mu 3}`$ could be big and $`V_{cb}`$ small, as observed.

A highly predictive $`SO(10)`$ example of a lopsided model is given in Ref. 5. The Dirac mass matrices have the following forms:

$$\begin{array}{cc}M_U=\left(\begin{array}{ccc}\eta & 0& 0\\ 0& 0& ϵ/3\\ 0& ϵ/3& 1\end{array}\right)m_U,\hfill & M_D=\left(\begin{array}{ccc}0& \delta & \delta ^{}\\ \delta & 0& ϵ/3\\ \delta ^{}& \sigma +ϵ/3& 1\end{array}\right)m_D,\hfill \\ & \\ M_N=\left(\begin{array}{ccc}\eta & 0& 0\\ 0& 0& ϵ\\ 0& ϵ& 1\end{array}\right)m_U,\hfill & M_L=\left(\begin{array}{ccc}0& \delta & \delta ^{}\\ \delta & 0& \sigma +ϵ\\ \delta ^{}& ϵ& 1\end{array}\right)m_D,\hfill \end{array}$$
(9)
where $`\eta \delta ,|\delta ^{}|ϵ\sigma 1`$. It is the large lopsided entry $`\sigma `$ that leads to large atmospheric neutrino mixing. The large solar mixing in this model arises from the matrix $`M_\nu `$. In this model there are six effective Yukawa operators, and seven real parameters ($`m_U/m_D`$, $`\sigma `$, $`ϵ`$, $`\delta `$, $`|\delta ^{}|`$, arg$`\delta ^{}`$, and $`\eta `$) to explain 14 measured quantities (8 mass ratios of quarks and charged leptons, 4 CKM parameters, and 2 of the neutrino mixing angles).

In the model just discussed, the large atmospheric angle comes from the charged lepton mass matrix, while the large solar angle comes from the neutrino mass matrix. Another possibility is that both large angles come from the charged lepton mass matrix, as happens in “doubly lopsided models” . In doubly lopsided models the entire third column of the matrix $`M_L`$ is large. The large 13 element gives the large solar angle, while the large 23 element gives the large atmospheric angle.

## 2 The Problem of $`\mu e\gamma `$

It is well known that SUSY GUT models with high scale SUSY breaking can lead to flavor-changing processes in the lepton sector . It is easy to see why. Consider the superpotential of the MSSM augmented to include the right-handed neutrinos:

$$\begin{array}{ccc}W& =& (Y_U)_{ij}\widehat{Q}_i\widehat{U}_j\widehat{H}_2+(Y_D)_{ij}\widehat{Q}_i\widehat{D}_j\widehat{H}_1\hfill \\ & +& (Y_N)_{ij}\widehat{L}_i\widehat{N}_j\widehat{H}_2+(Y_L)_{ij}\widehat{L}_i\widehat{E}_j\widehat{H}_1\hfill \\ & +& \frac{1}{2}(M_R)_{ij}\widehat{N}_i\widehat{N}_j+\mu \widehat{H}_1\widehat{H}_2.\hfill \end{array}$$
(10)
In the basis where $`Y_L`$ and $`M_L`$ are diagonal, the matrix $`Y_N`$ will in general have off-diagonal elements. Let the slepton masses be degenerate at a scale $`M_{}`$ that is at or above the unification scale $`M_G`$: $`m_{ij}^2(M_{})=m_0^2\delta _{ij}`$. The off-diagonal elements in $`Y_N`$ will induce off-diagonal elements in $`m^2`$ through the renormalization group running between $`M_{}`$ and the scale of the right-handed neutrino masses, which we will call $`m_R`$. The induced off-diagonal elements are given in the leading logarithm approximation by

$$(\delta m^2)_{ij}\frac{3m_0^2+A^2}{8\pi ^2}\underset{\mathrm{}}{}(Y_N)_i\mathrm{}^{}(Y_N)_j\mathrm{}\mathrm{ln}(M_{}/m_R).$$
(11)
The off-diagonal elements $`(\delta m^2)_{12}`$ and $`(\delta m^2)_{23}`$ will lead, when inserted into one-loop diagrams, to the decays $`\mu e\gamma `$ and $`\tau \mu \gamma `$ respectively. These effects have been calculated by several groups . One sees from the previous equation that if $`(Y_N)_{13}`$ and $`(Y_N)_{23}`$ are large at scales above $`m_R`$, there will be a significant contribution to $`(\delta m^2)_{12}`$, i.e. the mixing of $`\stackrel{~}{\mu }_L\stackrel{~}{e}_L`$, and thus to $`\mu e\gamma `$.

Let us consider the lopsided model shown before. One starts in the original basis with

$$Y_N\left(\begin{array}{ccc}\eta & 0& 0\\ 0& 0& ϵ\\ 0& ϵ& 1\end{array}\right)h_t,M_L=\left(\begin{array}{ccc}0& \delta & \delta ^{}\\ \delta & 0& \sigma +ϵ\\ \delta ^{}& ϵ& 1\end{array}\right)m_D,$$
(12)
where $`h_t`$ is the top quark Yukawa coupling, which is of order one. The matrix $`M_L`$ can be diagonalized in stages, the first being to eliminate the large 23 element $`\sigma +ϵ1`$ that is the source of the large atmospheric neutrino mixing angle. Rotating in the 23 plane of the left-handed lepton doublets by an angle $`\theta _{23}\theta _{atm}`$ one obtains

$$Y_N\left(\begin{array}{ccc}\eta & 0& 0\\ 0& ϵs_{atm}& s_{atm}\\ 0& ϵc_{atm}& c_{atm}\end{array}\right)h_t,M_L=\left(\begin{array}{ccc}0& \delta & \delta ^{}\\ \overline{\delta }& ϵs_{23}& 0\\ \overline{\delta ^{}}& ϵc_{23}& 1/c_{23}\end{array}\right)m_D.$$
(13)
We see that there is a large element $`(Y_N)_{23}`$ directly stemming from the large atmospheric mixing. The remaining steps in the diagonalization of $`M_L`$ will involve a rotation in the 12 plane of the left-handed lepton doublets to eliminate the 12 element $`\delta `$. This rotation is in fact just analogous to the Cabibbo rotation of the quarks, so we will call the rotation angle the “leptonic Cabibbo angle”, and denote it $`\theta _C^{\mathrm{}}`$. Just as $`\mathrm{tan}\theta _C\sqrt{m_d/m_s}`$, one expects in general that $`\mathrm{tan}\theta _C^{\mathrm{}}\sqrt{m_e/m_\mu }`$. Indeed, in this particular model $`\mathrm{tan}\theta _C^{\mathrm{}}\sqrt{m_e/m_\mu }`$. The result of this rotation on the matrix $`Y_N`$ is to give:

$$(Y_N)_{23}\mathrm{cos}\theta _C^{\mathrm{}}\mathrm{sin}\theta _{atm}h_t,(Y_N)_{13}\mathrm{sin}\theta _C^{\mathrm{}}\mathrm{sin}\theta _{atm}h_t.$$
(14)
As a result, in lopsided models the amplitude for $`\mu e\gamma `$ is proportional to $`\frac{1}{2}\mathrm{sin}2\theta _C^{\mathrm{}}\mathrm{sin}^2\theta _{atm}h_t^2\frac{1}{2}\sqrt{m_e/m_\mu }`$. That is, it is not suppressed by small Yukawa couplings or small mixing angles. In Ref. 10 the rate for $`\mu e\gamma `$ was calculated in the illustrative model that we have been discussing. There it was found that for $`\mathrm{tan}\beta =10`$, $`m_{1/2}=250`$ GeV, $`A_0=m_0`$, and $`\mu `$ positive the branching ratio for this process varies between $`2\times 10^9`$ (for $`m_0=100`$ GeV) and $`3\times 10^{11}`$ (for $`m_0=1`$ TeV) . The experimental limit is $`1.2\times 10^{11}`$ . Calculations that go beyond the leading logarithm approximation have been done , and these suggest that the branching ratio for $`\mu e\gamma `$ computed in Ref. 10 may be underestimated by about an order of magnitude. The lopsided models thus have a serious issue to confront.

However, t is not only lopsided models that face this issue. Consider SRND models in which the right-handed neutrino of the third family dominates the see-saw formula. Such models, as we see from Eq.(7), have $`(Y_N)_{23}\mathrm{sin}\theta _{atm}`$ in the original basis. And, since one typically expects that there will be a leptonic analog of the Cabibbo angle, the diagonalization of $`M_L`$ should lead to the relations in Eq. (14). In other words, the problem arises whether the large atmospheric neutrino mixing is the result of a large 23 element in $`M_L`$ or in $`M_N`$.

We see that there are four ingredients of the problem of excessive $`\mu e\gamma `$. (1) High scale SUSY breaking. The foregoing discussion assumed that at some high scale above $`m_R`$ there are nonzero slepton masses. (2) SO(10) or a larger group. This assumption relates $`(Y_N)_{33}`$ to $`(Y_U)_{33}h_t1`$. (3) Large $`\theta _{atm}`$ from large 23 element of either $`M_L`$ (lopsided) or $`M_N`$ ($`\nu _{R3}`$-dominated see-saw). This assumption leads to the result that $`(Y_N)_{23}(Y_N)_{33}1`$ in the basis where $`M_L`$ is diagonal. And (4) Significant leptonic Cabibbo mixing. This yields $`(Y_N)_{13}\mathrm{tan}\theta _C^{\mathrm{}}(Y_N)_{23}`$ in the basis where $`M_L`$ is diagonal.

## 3 Ways Out of the $`\mu e\gamma `$ Problem

Corresponding to the four ingredients that go into making the $`\mu e\gamma `$ problem serious, there are at least four approaches to avoiding or alleviating this problem. We will discuss them in the reverse order.

(4) Small leptonic Cabibbo angle. In several recent papers , models have been proposed that give acceptable $`\mu e\gamma `$ by making $`\theta _C^{\mathrm{}}`$ small. One example is one of the types of models analyzed by T. Blazek and S.F. King in Ref. 14. In this SRND type of model $`Y_N`$ has the following form (in the basis where $`M_L`$ is diagonal):

$$Y_N\left(\begin{array}{ccc}& \lambda ^n& \lambda ^m\\ & \lambda ^n& 1\\ & \lambda ^n& 1\end{array}\right)f,$$
(15)
where $`m2`$, $`\lambda 0.15`$, and $`f1`$. The matrix $`M_R`$ is assumed to be approximately diagonal. This form for $`Y_N`$ gives acceptably small $`\mu e\gamma `$. The reason is simple. As comparison of Eqs. (14) and (15) shows, these models implicitly assume that $`\mathrm{tan}\theta _C^{\mathrm{}}\lambda ^m\lambda ^20.02`$. This is to be compared to the typical value of order $`\sqrt{m_e/m_\mu }0.07`$ that arises in predictive SUSY GUT models. This gives, then, at least an order of magnitude suppression of $`\mu e\gamma `$.

Another example is provided by a model presented by Bi in Ref. 10. In this model, which is of the “doubly lopsided” type, the matrix $`M_N`$ is approximately diagonal, while $`M_L`$ has the form

$$M_L=\left(\begin{array}{ccc}0& \delta & \sigma \\ \delta & 0& 1\sigma \\ 0& ϵ& 1\end{array}\right)m,$$
(16)
where $`\sigma 0.58`$, $`ϵ0.12`$, and $`\delta 0.00077`$.

This model is quite similar in spirit to the lopsided model described above. Indeed, if one “rotates away” the large 13 element one gets a matrix of very similar form to Eq. (12). Going through the diagonalization procedure, one easily sees that the leptonic Cabibbo angle comes out to be of order $`\delta /ϵ0.01`$.

The trouble is that in a GUT model the matrix $`M_L`$ is intimately related to the matrix $`M_D`$. Thus, if this model were to arise from a GUT one would have to explain why $`\mathrm{sin}\theta _C0.21`$, while $`\mathrm{sin}\theta _C^{\mathrm{}}0.01`$. Unfortunately, the model as presented in Ref. 10 gives only the leptonic matrices, and so is not a complete unified model of quark and lepton masses and mixings.

To see that it is difficult to make the leptonic Cabibbo angle small in a predictive SUSY GUT model of quark and lepton masses, let us see what happens if we attempt to do this by modifying the model of Ref. 5 described earlier. To allow the leptonic Cabibbo angle to be small let us add some more parameters. In particular, take the matrices $`M_L`$ and $`M_D`$ to have the forms:

$$M_L=\left(\begin{array}{ccc}0& \delta \overline{\delta }& \delta ^{}\overline{\delta }^{}\\ \delta +\overline{\delta }& 0& \sigma +ϵ\\ \delta ^{}+\overline{\delta }^{}& ϵ& 1\end{array}\right)m_D,M_D=\left(\begin{array}{ccc}0& \delta +\frac{1}{3}\overline{\delta }& \delta ^{}+\frac{1}{3}\overline{\delta }^{}\\ \delta \frac{1}{3}\overline{\delta }& 0& \frac{1}{3}ϵ\\ \delta ^{}\frac{1}{3}\overline{\delta }^{}& \sigma +\frac{1}{3}ϵ& 1\end{array}\right)m_D.$$
(17)
Here the new parameters $`\overline{\delta }`$ and $`\overline{\delta }^{}`$ (both complex) have been added, so that the model has 11 instead of 7 real parameters. These forms can be obtained from fairly simple effective operators.

It is clear that the leptonic Cabibbo angle will be suppressed if the 12 element of $`M_L`$ is sufficiently small. This can be achieved here through a cancellation between $`\delta `$ and $`\overline{\delta }`$. Let us define a suppression parameter $`\alpha `$ by $`\theta _C^{\mathrm{}}=\alpha \sqrt{m_e/m_\mu }`$. Thus, in the model of Ref. 5 one has that $`\overline{\delta }=\overline{\delta }^{}=0`$ and $`\alpha =1`$.

Interestingly, one finds that although we have added several parameters, the model is almost as predictive as before. In fact, it gives all the same predictions, essentially, except that the prediction for $`V_{ub}`$ is altered. One finds

$$\frac{V_{ub}}{V_{cb}V_{us}/V_{cs}}\frac{1}{\sigma ^2}\left[1+\sqrt{\sigma ^2+1}(2s_1s_2\alpha )\right],$$
(18)
where $`s_1`$ and $`s_2`$ are phases. The model of Ref. 5 corresponds to the case $`\alpha =1`$, $`s_1=s_2e^{i\theta }`$, and $`\sigma \sqrt{3}`$., giving $`V_{ub}/(V_{cb}V_{us}/V_{cs})\frac{1}{3}+\frac{2}{3}e^{i\theta }`$, which gives an excellent fit to the data. (It gives a circle in the complex plane of radius $`\frac{2}{3}`$ centered at the point $`(\frac{1}{3},0)`$. This passes through the experimentally allowed region.) However, if one tries to suppress the leptonic Cabibbo angle, that means that $`\alpha 1`$, giving $`V_{ub}/(V_{cb}V_{us}/V_{cs})\frac{1}{3}+\frac{4}{3}e^{i\theta }`$, which is inconsistent with data by a large margin, since experimentally $`|V_{ub}/(V_{cb}V_{us}/V_{cs})|0.4`$, while the prediction obviously implies that $`|V_{ub}/(V_{cb}V_{us}/V_{cs})|1.0`$.

Another way to modify the model of Ref. 5 to allow small leptonic Cabibbo angle is by changing the entries for the first family as follows:

$$M_L=\left(\begin{array}{ccc}0& \delta \eta & 0\\ \delta ^{}+\eta & 0& \sigma +ϵ\\ 0& ϵ& 1\end{array}\right)m_D,M_D=\left(\begin{array}{ccc}0& \delta ^{}+\frac{1}{3}\eta & 0\\ \delta \frac{1}{3}\eta & 0& \frac{1}{3}ϵ\\ 0& \sigma +\frac{1}{3}ϵ& 1\end{array}\right)m_D.$$
(19)
Here again one sees that the possibility exists of a suppression of the leptonic Cabibbo angle through a cancellation between $`\delta `$ and $`\eta `$. The forms in Eq. (19) give a reasonable fit to the masses of the first family fermions, and to the Cabibbo angle, but yields the prediction

$$\frac{V_{ub}}{V_{cb}V_{us}/V_{cs}}\frac{1}{\sqrt{1+\sigma ^2}}\left[2\alpha s_1\right]s_2,$$
(20)
where $`\alpha `$ is defined as before, and $`s_1`$ and $`s_2`$ are phases. For $`\alpha 1`$ (i.e. small leptonic Cabibbo angle) one has then that $`V_{ub}/(V_{cb}V_{us}/V_{cs})`$ is a complex number of unit magnitude, a result completely inconsistent with experiment.

There are other ways to modify the model of Ref. 5 to make it more flexible and able to accomodate a small leptonic Cabibbo angle, besides the two examples just shown. However, if the models remain predictive, they tend to give a bad fit to some quantity, usually $`V_{ub}`$. Moreover, the suppression of the leptonic Cabibbo angle is achieved generally by an artificial cancellation, as in these two examples. This does not prove that predictive SUSY GUT models of quark and lepton masses with small leptonic Cabibbo angle cannot be constructed, but it suggests that constructing such models is not a trivial matter.

(3) $`\theta _{atm}`$ does not arise from large 23 elements of $`M_L`$ or $`M_N`$. One can imagine a number of ways that a large atmospheric angle could arise without there being a large 23 element in either $`M_L`$ or $`M_N`$. One way is through SRND with the first or second family right-handed neutrino dominating as in Eq. (5). The main question is whether this idea can be implemented in predictive SUSY GUT models. There may be some difficulty in doing this, since the form of $`M_N`$ shown in Eq. (5) may clash with the form of $`M_U`$ needed to obtain $`m_c/m_tV_{cb}`$. (Having a form for $`M_U`$ that is similar to the form for $`M_N`$ given in Eq. (5) would give $`m_c/m_tϵ`$ and $`V_{cb}ϵ`$.)

A second possibility is that a large atmospheric angle might arise from a non-see-saw mechanism. That is it could be that $`M_\nu =M_NM_R^1M_N^T+M_\nu ^{nonseesaw}`$. Then $`M_N`$ and $`M_L`$ could both have small 23 elements. However, such a scenario sacrifices most of the advantages of SUSY GUT see-saw models enumerated in the first section of this paper. One would have to invoke new ad hoc physics to generate $`M_\nu ^{nonseesaw}`$. One would have to explain why the non-see-saw mechanism happened to give neutrino masses that were of the same magnitude (roughly speaking) as the see-saw mechanism. And one would give up any predictivity for the neutrinos’ masses and mixings.

A third possibility is that in fact both $`M_L`$ and $`M_N`$ have large 23 elements, but the flavor-changing effects are suppressed by a cancellation between them. An example of this is the model of Babu, Pati, and Wilczek . In that model, $`(Y_N)_{23}/(Y_N)_{33}=\sigma ϵ`$ and $`(Y_L)_{23}/(Y_L)_{33}=\eta ϵ`$, where it happens that $`\sigma ϵ=0.395`$ and $`\eta ϵ=0.436`$. Thus, in the basis where $`M_L`$ ($`Y_L`$) is diagonal, the matrix $`Y_N`$ has a small 23 element by cancellation. This cancellation does not occur for the atmospheric angle, however, since the form of $`M_R`$ is such that $`M_\nu `$ has a small 23 element. It should be noted that the same cancellation between $`\sigma `$ and $`\eta `$ happens for the quark sector and makes $`V_{cb}`$ small. One drawback to this model is that the cancellation is a numerical accident, though perhaps another version of the model could be found where this accident is explained.

(2) An SU(5) or SU(5)-like model. The group $`SO(10)`$ relates $`Y_N`$ to $`Y_U`$ and tends to predict that $`(Y_N)_{33}h_t`$, but $`SU(5)`$ leaves $`(Y_N)_{33}`$ free to be much smaller than $`h_t`$. Unfortunately, $`SU(5)`$ is much less predictive for quark and lepton masses than $`SO(10)`$ or larger groups.

An interesting possibility is the following. Suppose that $`SO(10)`$ breaks to $`SU(5)`$ at some scale $`M_{10}M_{}`$ (the scale at which slepton masses are degenerate). In the $`SO(10)`$ theory above $`M_{10}`$ let the Dirac Yukawa coupling matrix of the neutrinos be $`Y_N^{(0)}`$, which couples the $`\nu `$ and $`N^c`$ in the three families ($`16_i`$) to each other. Let there also be $`SO(10)`$-singlet fermions $`N_i`$ in each family which get masses from the terms $`M_{ij}N_iN_j+Y_{ij}N_i^cN_j\mathrm{\Sigma }(\overline{16})`$, where both $`M`$ and $`\mathrm{\Sigma }`$ are of order $`M_{10}`$. Then below $`M_{10}`$ the remaining $`SU(5)`$-singlet fermions (i.e. “right-handed neutrinos”) are linear combinations of $`N_i^c`$ and $`N_i`$. Thus, the effective Dirac Yukawa coupling matrix of the neutrinos below $`M_{10}`$ is $`Y_NY_N^{(0)}[(Y\mathrm{\Sigma })^1M]Y_N^{(0)}\lambda `$. The parameters can easily be such that the matrix $`\lambda `$ has small elements. Then the elements of $`Y_N`$, and in particular the 33 element can be small. This would suppress flavor changing effects due to running between $`M_{10}`$ and $`m_R`$, as one sees from Eq. (11).

However, above the scale $`M_{10}`$ there are large Dirac Yukawa couplings for the neutrinos, since one expects that $`(Y_N^{(0)})_{33}h_t`$ by $`SO(10)`$ symmetry. Thus, if $`M_{}M_{10}`$, one will still get contributions to flavor-changing slepton masses, but with $`\mathrm{ln}(M_{}/m_R)`$ in Eq. (11) replaced by $`\mathrm{ln}(M_{}/M_{10})`$, which could easily be of order one. In this way, one suppresses the rate for $`\mu e\gamma `$ by one or two orders of magnitude.

(1) SUSY breaking occurs at a low scale. The whole discussion of the flavor changing problem above is predicated upon the assumption that there is high-scale SUSY breaking, so that there are slepton masses at a scale $`M_{}`$ that is large compared to the right-handed neutrino mass scale $`m_R`$. In schemes such as gauge-mediated SUSY breaking, there would not be a problem with $`\mu e\gamma `$ such as we have been analyzing.

In conclusion, it would seem that most ways to suppress $`\mu e\gamma `$ in the context of high-scale SUSY breaking only do so to a rather limited extent. If $`SO(10)`$ is broken to $`SU(5)`$ at a scale near that at which slepton masses are degenerate, it gives a suppression of only an order of magnitude or two in the rate. The leptonic Cabibbo angle is difficult to suppress. The SRND models in which the right-handed neutrino of the first or second family dominates could give a strong suppression, but it is doubtful whether realistic and predictive schemes based on this idea can be constructed; so far they have not been.

By far the cleanest way to get rid of the problem is by some low-scale SUSY breaking. Perhaps this is also the most plausible possibility. There is, even if we ignore large neutrino mixings and GUTs, a well-known problem with flavor changing processes in models with high-scale SUSY breaking. None of the proposed solutions to this “SUSY Flavor Problem” have so far been very convincing. One can regard the problem of $`\mu e\gamma `$ in SUSY GUT see-saw models as simply a particularly acute form of a more general problem with high-scale SUSY breaking. It is reasonable to suppose that whatever mechanism resolves the other flavor-changing problems of SUSY may also resolve this one.


text:
# The masses of 𝐷_{𝑠𝑗}^∗(2317) and 𝐷_{𝑠𝑗}^∗(2463) in the MIT bag model

## I Introduction

A new scalar meson with the $`c\overline{s}`$ flavor content and a mass around 2.32 GeV has been recently discovered by the BaBar collaboration auber  as a narrow resonance in the $`D_s\pi ^0`$ decay channel and confirmed by the CLEO collaboration cleo . This state aroused a particular interest because of a plain conflict with the theoretical calculations which predict the mass about $`100300`$ MeV above its experimental value godfrey . This discrepancy opened a possibility for more refine explanations of the new state barnes . However, it is not clear, weather the difference between the predictions and experiments rises a real issue or it is rather a problem of the models. Besides all that, there are successful explanations of the meson spectrum based on the chiral models nowak ; bardeen1 ; feldmann ; bardeen2 ; nowak2  under the assumption of the $`SU(3)_L\times SU(3)_R`$ chiral symmetry. The QCD sum rules results dai  also agree with the experiments within the theoretical uncertainty.

This short note is, as a matter of fact, a supplement to the papers devoted to the spectroscopy of hadrons which contain one heavy quark shuryak ; izatt ; sad . In these papers the MIT bag model chodos  was adapted to heavy-light systems ($`Q\overline{q},Qqq`$ were $`Q=c,b`$ and $`q=u,d,s`$) and masses of mesons and barions were calculated for the lowest energy states $`J^P=0^{},1^{},(\frac{1}{2})^+,(\frac{3}{2})^+`$. Using the same model and the same parameters (with justified changes in the quark masses) one can find the values of 2.29 GeV for the $`0^+`$ and $`2.46`$ for the $`1^+`$ $`c\overline{s}`$ states.

## II Mass calculation

The mass of the bag can be written in the form chodos2 :

$$M=E_Q+E_q+E_{bag}+E_0+E_M+E_E$$
(1)
where $`E_{Q(q)}`$ is the energy of the heavy(light) quark, $`E_{bag}=\frac{4}{3}\pi BR^3`$ is the volume energy of the spherical bag of the radius $`R`$, $`E_0=Z/R`$ is the Casimir energy and the last two terms describe the energies of the color-magnetic and color-electric interactions between quarks. Let us first notice that for the heavy-light systems the center of mass of the bag is well defined. It coincides with the heavy quark position and the corrections due to the c.m. motion can be neglected. The energy of the heavy quark is approximated by the formula $`E_Q=m_Q+p_{light}^2/2m_Q`$, where $`m_Q`$ is the mass of the heavy quark and the second term is the non-relativistic kinetic energy. The total momentum of the bag in the c.m. frame should be zero thus the values of the heavy and light quark momenta are equal. Let us notice that the non-relativistic kinetic energy term in $`E_Q`$ was neglected in the previous calculations izatt ; sad . This is justified if one wishes to study the lowest energy states. However for exited states, the contribution of the kinetic energy of the heavy quark is more important and should be included. The energy of the color-magnetic interactions between the heavy and light quarks is equal to izatt ; sad

$$E_M=\frac{\lambda A\alpha _s}{m_QR^2}\stackrel{}{\sigma }_Q\stackrel{}{\sigma }_q$$
(2)
where $`\lambda =\{1,2\}`$ for barions and mesons respectively, $`A`$ is calculated from the light quark wave-function and depends on its mass. The constant

$$\alpha _s=\frac{0.47}{\mathrm{log}(1+1/\mathrm{\Lambda }R)},\mathrm{\Lambda }=0.42\text{GeV}$$
(3)
is inspired by the QCD running coupling formula sad ; donoghue . The spin-spin interaction term $`\stackrel{}{\sigma }_Q\stackrel{}{\sigma }_q=\{3,1\}`$ for scalar and vector mesons respectively. There are two sources which contribute to the color-electric energy: the interactions between quarks and their self-interactions. For mesons with one heavy quark it can be written as izatt :

$$E_E=\frac{4\alpha _s}{3R}\left(C+\frac{1}{2}\right)$$
(4)
where $`C`$ can be calculated and depends on the mass of the light quark. The second term in the parenthesis comes from the self-interaction of the heavy quark. The heavy quark self-interaction is divergent and the infinite term is omitted in (4)izatt . The radius of the bag is determined by minimizing the formula $`E^{}=E_{bag}+E_0+E_q`$ with respect to $`R`$111One can include the kinetic part of the heavy quark energy in $`E^{}`$. The change in the radius can be absorbed by the changes of the parameters of the model with no important impact on the spectrum. The discussion of the interaction treatment for $`R`$ determination can be found in izatt ..

The value of the bag constant $`B^{1/4}=0.1383`$ and the Casimir parameter $`Z=0.73`$ are taken from the paper sad . The masses of the quarks $`m_u=m_d=0`$, $`m_s=0.198`$ GeV, $`m_c=1.824`$ GeV, $`m_b=5.21`$ GeV were determined to reproduce the masses of the mesons $`D,D_s,B`$. Using model parameters in the $`0^{},1^{}`$ sector one finds the light quark energy $`E_q=2.04/R`$ for the $`u,d`$ quarks and $`E_s=\sqrt{m_s^2+(2.35)^2/R^2}`$ for the strange quark. The radius of the non-strange meson is $`R=4.1`$ GeV-1 and for the strange meson $`R_s=4.2`$ GeV-1. The color-magnetic interaction parameter is $`A(m_q=0)=0.52`$ and $`A(m_s)=0.5`$ whereas for the color-electric part one finds $`C(m_q=0)=0.52`$ and $`C(m_s)=0.62`$. These numbers together with the equations (1 \- 4) give the spectrum written in the Table 1.

For the scalar and pseudovector mesons, where the light quark is in the $`1P_{1/2}`$ state, the light quark energy is $`E_q=3.81/R`$ for the $`u,d`$ quarks and $`E_s=\sqrt{m_s^2+(3.91)^2/R^2}`$ for the strange quark. The radii are $`R=5.09`$ GeV-1 for the non-strange and $`R_s=5.08`$ GeV-1 for the strange mesons. The color-magnetic interaction parameter is $`A(m_q=0)=0.87`$ and $`A(m_s)=0.84`$ whereas the color-electric parameter takes values $`C(m_q=0)=0.88`$ and $`C(m_s)=0.79`$. The masses of the $`0^+,1^+`$ states are given in the Table 2.

The mass splitting relations $`\mathrm{\Delta }M_{1^+1^{}}=\mathrm{\Delta }M_{0^+0^{}}`$ nowak ; bardeen1  (where $`\mathrm{\Delta }M_{1^+1^{}}=M_{1^+}M_1^{}`$ and so on) and $`\mathrm{\Delta }M_{1^+0^+}=\mathrm{\Delta }M_{1^{}0^{}}`$ bardeen2 ; nowak2  should hold in the chiral limit. Indeed the experiment finds cleo  $`\mathrm{\Delta }M_{1^+1^{}}\mathrm{\Delta }M_{0^+0^{}}=1.2\pm 2.1`$ MeV and $`\mathrm{\Delta }M_{1^+0^+}\mathrm{\Delta }M_{1^{}0^{}}=1.6\pm 2.1`$ MeV for the $`c\overline{s}`$ mesons. The MIT bag model fails to reproduce the mass splitting equalities at this level of accuracy. This can be expected because the MIT bag model does not handle the chiral symmetry in the appropriate way. In spite of that the model predicts that the following relations hold (at the level of $`1020`$ MeV) independently of the flavor of the light quark:

$`\mathrm{\Delta }M_{1^+1^{}}(c\overline{q})`$ $``$ $`\mathrm{\Delta }M_{1^+1^{}}(c\overline{s}),\mathrm{\Delta }M_{0^+0^{}}(c\overline{q})\mathrm{\Delta }M_{0^+0^{}}(c\overline{s})`$ (5)
$`\mathrm{\Delta }M_{1^+0^+}(c\overline{q})`$ $``$ $`\mathrm{\Delta }M_{1^+0^+}(c\overline{s}),\mathrm{\Delta }M_{1^{}0^{}}(c\overline{q})\mathrm{\Delta }M_{1^{}0^{}}(c\overline{s})`$
where $`q=u,d`$. These relations follow from the $`SU(3)_L\times SU(3)_R`$. However this symmetry is broken by the non-zero, strange quark mass and it is not clear why it should work so well. Let us stress that from the relations (5) only the last one is confirmed by the experiments at the accuracy of $`12`$ MeV. Other relations are the assumptions based on the $`SU(3)_L\times SU(3)_R`$ chiral symmetry or follow from the model calculations.

Finally, to be clear, let us notice that the considered MIT bag model joins several features of works izatt ; sad . The parameters $`B,Z`$ as well as the description of the coupling $`\alpha _s`$ comes from the paper sad . The treatment of the infinity in the electric heavy quark interaction follows paper izatt . The heavy quark kinetic term was also added. The quark masses are adjusted anew to accommodate the changes in the model.

## III Conclusions

The result of the previous paragraph, at the face value, can be regarded as a success of the MIT bag model. The description of the meson spectrum(similar for barions shuryak ; izatt ; sad ) at the level of 4-5 per cent is already at the edge of the reliability of the model.

However the more careful scrutiny leads to additional questions. Let us first notice that the quark masses are not treated on an equal footing. The $`m_u,m_d,m_s`$ are current masses whereas $`m_c,m_b`$ are constituent masses. This follows from the treatment of the heavy quark self-interaction in the color-electric energy. The infinite term is absorbed to the mass of the heavy quark. Secondly, the coupling constant at the energy scale $`1/R=0.250.1`$ GeV takes values $`\alpha _s=1.01.2`$. This is certainly not a small parameter. Indeed one can try to fit a new set of the bag model parameters with current quark masses and a small coupling constant. First, the infinite term could be cut-off by the heavy quark Compton wave-length $`1/m_Q`$. Secondly, the scale of the interaction between the light and heavy quarks may be rather set by the $`m_Q`$ or geometrical mean $`\sqrt{m_QE_q}`$. If done, then the coupling constant can be reduced to the value of $`\alpha _s0.30.4`$ and the heavy quarks get their current masses. But in this situation, the mass of $`D_{sj}^{}(2317)`$ immediately falls in the range $`2.442.55`$ GeV. The simultaneous requirement of the current heavy quark masses and a small coupling constant can not be reconciled with a good fit to the whole meson spectrum222For example $`B^{1/4}=0.1383`$ GeV, $`Z=1.85`$, $`m_u=m_d=0`$, $`m_s=0.143`$ GeV, $`m_c=1.482`$ GeV, $`m_b=4.6`$ GeV gives a relatively good fit to the lowest energy hadron spectrum. The energy scale for the heavy-light quarks interactions is set to $`\sqrt{m_QE_q}`$. The infinite term is regularized by the Compton wave-lenght of the heavy quark. The author tried different sets of parameters with reasonable assumptions about the energy scale of the interactions and regularization methods with similar success..

Despite of these drawbacks the MIT bag model is capable to predict the $`0^\pm ,1^\pm `$, $`D`$ mesons spectrum with the accuracy of few per cent. One can also find that the mass splitting relations (5) are independent of the flavor of the light quark. These relations follow directly from the chiral symmetry $`SU(3)_L\times SU(3)_R`$. However this symmetry is broken by the strange quark mass and it is not obvious why the relations (5) should hold.

In the further work one can try to include the light quark in the $`1P_{3/2}`$ state. However the calculation of the $`1P_{3/2}`$ state is more involved and require reconsideration of the model assumptions (e.g. sphericity of the bag) thus demanding much more extensive elaboration then this short note assumed.

Acknowledgements. I would like to thank Leszek Motyka and Kacper Zalewski for useful discussions. This work was supported by the Polish State Committee for Scientific Research, (KBN) grant no. 2 P03B 09322.


text:
# H2 in the interstitial channels of nanotube bundles

## I INTRODUCTION

Theoretical and experimental research in carbon nanotubes (CN’s) is nowadays a very active field in which new avenues are opening continuously, all related to its particular structure at the nanometer-scale level. Besides their singular electrical and mechanical properties, CN’s are able to adsorb some atoms and molecules offering a real possibility for quasi-one-dimensional environments. These are mainly the inner part of the nanotubes and the interstitial channels (IC’s) among them due to their natural bundle arrangement. When CN’s are formed, they appear as long capped cylinders which only can adsorb a given substance if their caps are removed by chemical means. Instead, IC’s are always present offering very narrow channels, distributed according to an hexagonal lattice, and with a sizeable adsorption capability in comparison with adsorption in planar graphite. Recently, Talapatra et al. talapatra  studied experimentally the binding energies of Xe, CH4, and Ne on close-ended nanotube bundles and concluded that none of them is adsorbed in the intersites. Probably, the size of these atoms and molecules is too large compared with the radius of the interchannel ($`r_{\text{IC}}3`$ Å) to fit into them. However, the situation for lighter species like He or H2 seems different. There are several experiments which claim that both of them can be adsorbed in the IC’s due to their small radius.teizer ; williams ; ye  H2 is certainly the most interesting system since a CN bundle could be one of the best options to reach the target energy densities for a lightweight hydrogen-storage system usable for transportation. Dylon ; dylon2

In a recent theoretical work, Calbi et al. Cole1  have studied the adsorption of H2, and other molecules and atoms, in the IC’s of a CN bundle. Introducing in the formalism a possible dilation of the bundle, by means of a harmonic approximation, they concluded that in equilibrium the bundle can swell. The dilation is observed in all cases, with different intensity depending on the particular system. In the geometry there considered, corresponding to (10,10) tubes in the standard terminology, and when the gas adsorbed is H2 that swelling amounts approximately to 1 %. This is the result of a delicate balance between two competing effects when the lattice is expanded: a decrease in the CN bundle cohesion energy and an increase in the H2-nanotube binding energy. Obviously, this prediction for H2 acquires special relevance since it would imply a significant increase of the CN bundle storage power. However, at present there is not any experimental confirmation of this possible adsorption-induced dilation. Probably, even if such effect is indeed present, its manifestation could be hardly observed due to the present experimental uncertainties. Raman spectroscopy seems a promising method for a better insight. In a recent work, this method was applied successfully by Williams et al. williams  to discriminate the different adsorption sites of H2 in CN bundles.

In the present work, we present a diffusion Monte Carlo (DMC) study of H2 adsorbed in IC’s of a (10,10) CN bundle. In the past, and using the same methodology, we characterized the ground-state properties of 4He and H2 inside a single CN nano4 ; nano2  and of 4He adsorbed in an IC.nanoic  In this paper, our aim is to determine the equation of state of H2 using microscopic theory, with special attention to the possibility of an adsorption-induced bundle dilation. To this end, the influence on swelling of inputs like the nanotube radius, the molecule-molecule and CN-molecule interactions, the geometry of the intersite, and the transversal degree of freedom are thoroughly discussed. Our results agree somehow with the ones of Ref. Cole1,  but the magnitude of a possible swelling is seen to be comparatively smaller or even inexistent.

The rest of the paper is organized as follows. In Section II, the influence of the several parameters of the model on a possible dilation of the bundle is studied using a one-dimensional approximation. The accuracy of this model is then tested in Section III by means of a three-dimensional DMC calculation. Finally, Section IV comprises a brief summary and the main conclusions.

## II ONE-DIMENSIONAL APPROXIMATION: TESTING THE MODEL

As in Ref. Cole1, , we study H2 adsorbed in the IC of a (10,10) CN bundle. The bundle is disposed in such a way that, in the plane perpendicular to the CN long axis, a triangular lattice is formed. The lattice constant, i.e., the distance between adjacent centers, is $`D_0=17`$ Å. $`D_0`$ corresponds to the equilibrium position and there is an overall agreement on its value. Therefore, we consider it as a fixed parameter in all the calculations.

Having established the value of $`D_0`$, the radius of a single nanotube is directly related to the C-C minimum separation between two neighbors. In a theoretical calculation, Tersoff and Ruoff Tersoff  concluded that this distance is 3.4 Å, in agreement with the most accepted value for the Lennard-Jones $`\sigma _{\text{C-C}}`$. In fact, this value coincides with the experimentally accepted distance between graphite sheets. Considering this C-C distance, the radius of a (10,10) tube is 6.8 Å. It is worth mentioning that this radius is the same one obtains constructing the same nanotube by rolling up a graphitic surface with a C-C distance of 1.42 Å. However, there is not a general agreement on those values. Recently, Charlier et al. euro  carried out a density functional calculation of a CN bundle and obtained an equilibrium geometry corresponding to a C-C intertube distance of $``$ 3.2 Å. The CN radius is then 6.9 Å, a value which was used in Ref. Cole1, . At present and to our knowledge, there is not experimental data on the C-C distance which could help to choose between different theoretical models. Therefore, we consider in the present calculation two possible CN radius, 6.8 and 6.9 Å.

The interstitial channel between three adjacent nanotubes in a (10,10) bundle can lodge a hard sphere with radius 2.9 or 3 Å, corresponding to CN’s with radius 6.8 and 6.9 Å, respectively. Since the $`\sigma `$ parameter for the C-H2 Lennard-Jones potential is around 3 Å, one can reasonably assume a one-dimensional (1D) approximation for hydrogen adsorption in IC’s. In this approximation, the H2 energy per particle ($`e=E/N`$) can be written as Cole1

$$e(\lambda ,D)=ϵ_0(D)+e_{\text{1D}}(\lambda )+h(\lambda ,D)+\frac{1}{\lambda }\frac{3}{4}k(DD_0)^2,$$
(1)
with $`\lambda `$ the linear density of the H2 molecules and $`D`$ the lattice parameter. In Eq. (1), $`ϵ_0(D)`$ is the binding energy of a single molecule in the IC, $`e_{\text{1D}}(\lambda )`$ is the 1D H2 energy per particle, $`h(\lambda ,D)`$ corresponds to the interaction energy with H2 in neighboring IC’s, and the last term takes into account the cohesion energy of the bundle when it is dilated from the equilibrium lattice distance $`D_0`$. In the following, we analyze the several terms entering Eq. (1) to disentangle which are the relevant inputs influencing a possible dilation of the bundle.

Both $`ϵ_0(D)`$ and $`e_{\text{1D}}(\lambda )`$ are calculated using the diffusion Monte Carlo (DMC) method.les ; boro4he  The DMC method solves in a stochastic way the $`N`$-body Schrödinger equation in imaginary time for the wave function $`f(𝑹,t)=\psi (𝑹)\mathrm{\Psi }(𝑹,t)`$,

$`{\displaystyle \frac{f(𝑹,t)}{t}}`$ $`=`$ $`D\mathbf{}^2f(𝑹,t)+D\mathbf{}(𝑭f(𝑹,t))`$
$`+(E_\text{L}(𝑹)E)f(𝑹,t),`$
with $`\psi (𝑹)`$ a trial wave function introduced for importance sampling. In Eq. (II), $`D=\mathrm{}^2/2m`$, $`E_\text{L}(𝑹)=\psi (𝑹)^1H\psi (𝑹)`$ is the local energy, and $`𝑭(𝑹)=2\psi (𝑹)^1\mathbf{}\psi (𝑹)`$ is the so-called drift force which guides the diffusion movement to regions where $`\psi (𝑹)`$ is large.

The trial wave function $`\psi `$ used in the calculation of 1D H2 is a Jastrow one,

$$\psi ^{\text{1D}}(𝑹)=\underset{i<j}{}f(r_{ij}),$$
(3)
with $`f(r_{ij})=\mathrm{exp}\left(0.5\left(b/r_{ij}\right)^5\right)`$ a McMillan two-body correlation factor. The parameter $`b`$ in $`f(r)`$ has been determined by means of a variational Monte Carlo (VMC) optimization. Near the equilibrium density, the optimal value is $`b=3.996`$ Å and it increases gradually with the density $`\lambda `$ (at $`\lambda =0.358`$ Å-1, $`b=4.026`$ Å).

The trial wave function used in the calculation of $`ϵ_0(D)`$ includes two-body correlations with the three nanotubes surrounding the IC. They are also of McMillan type,

$$\psi ^{\text{IC}}(𝑹)=\underset{n=1}{\overset{3}{}}\mathrm{exp}\left[\frac{1}{2}\left(\frac{a}{r_n}\right)^5\right],$$
(4)
$`r_n`$ being the distance of the hydrogen molecule to the center of any of the three tubes. The parameter $`a`$, optimized using VMC, varies from $`a=22.5`$ Å for zero dilation to $`a=21.5`$ Å for $`DD_0`$ = 0.2 Å. On the other hand, the optimal values for the parameter $`a`$ present a negligible dependence with the particular C-H2 potential chosen in the calculation. The IC is so narrow that one can approximate the sum of individual CN-molecule potentials by a new one which is only a function of the radial distance $`r`$ to the center of the IC. This simplified model is obtained by an azimuthal average of the three individual potentials.Cole1  In this case, the trial wave function is simpler than the previous model (4). We have chosen a gaussian

$$\psi _\text{c}(𝑹)=\mathrm{exp}(cr_n^2),$$
(5)
the parameter $`c`$ varying from 7.6 Å-2 to 11.5 Å-2, depending on the radius of the tube and on the dilation. The greater the radius, the smaller the value of $`c`$.

Para-hydrogen is spherical to a large extent. As usual in most of molecular hydrogen calculations, the intermolecular interaction is considered purely radial and described by the Silvera and Goldman (hereafter SG) model.Silvera  In addition, and to make comparisons with previous work (Ref. Cole1, ), we have also made some calculations with the potential proposed by Kostov, Cole, Lewis, Diep, and Johnson.kostov  This H2-H2 potential (hereafter KCLDJ) incorporates, in a rather crude way, three-body corrections to the pair potential coming from the triplets H2-C-H2.

Much more critical for estimating a possible bundle dilation is the model for the CN-H2 interaction. As in our previous work on the equation of state of H2 adsorbed inside a nanotube,nano2  we use the CN-H2 potential proposed by Stan and Cole.Stan  This potential results from an average over all the C-H2 interactions between the C of an infinite CN and a single molecule located at a radial distance $`r`$. The interaction is thereby independent of $`z`$ (corrugation effects are neglected) and can be used for adsorption both inside and outside the nanotube. The resulting CN-H2 potential depends explicitly on the $`\sigma `$ and $`ϵ`$ parameters of the Lennard-Jones potential between C and H2. At present, there is not an overall agreement about which are the best set of parameters ($`\sigma `$,$`ϵ`$) describing this interaction. In order to study the influence of this choice on the calculations, two different sets have been studied. The first one, used for example in Refs. Cole1,  and uptake, , is derived with the Lorentz-Berthelot combining rules

$`\sigma _{\text{C-H}_2}`$ $`=`$ $`{\displaystyle \frac{\sigma _{\text{H}\text{2}\text{-H}\text{2}}+\sigma _{\text{C-C}}}{2}}`$ (6)
$`ϵ_{\text{C-H}_2}`$ $`=`$ $`\sqrt{ϵ_{\text{H}\text{2}\text{-H}\text{2}}ϵ_{\text{C-C}}}.`$
Using $`\sigma _{\text{C-C}}=3.4`$ Å, $`ϵ_{\text{C-C}}=28`$ K, $`\sigma _{\text{H}\text{2}\text{-H}\text{2}}=3.05`$ Å, and $`ϵ_{\text{H}\text{2}\text{-H}\text{2}}=37.0`$ K the first set is $`ϵ_{\text{C-H}_2}=32.2`$ K and $`\sigma _{\text{C-H}_2}=3.23`$ Å(hereafter LB). The second option is to consider the optimal parameters which describes the interaction of H2 with planar graphite. This set from Wang, Senbetu, and Woo wang  is probably more realistic than the first one and it has also been used in the past.williams ; nano2 ; sc  The values are $`ϵ_{\text{C-H}_2}=42.8`$ K and $`\sigma _{\text{C-H}_2}=2.97`$ Å(hereafter WSW). Notice the sizeable difference between both sets which will generate significant differences in the final results.

Dilation and compression of the lattice constant in the CN bundle have an energetic cost. Following Refs. Mizel,  and Cole1, , that contribution is assumed to be harmonic in the displacement around the equilibrium position $`D_0`$ (see Eq. 1). This approximation is expected to hold for small $`DD_0`$ values but the uncertainty in the real value of $`k`$ is rather large. From compression modes measured in graphite by Nicklow et al.,nicklow  Mizel et al. Mizel  estimated a value $`k=1740`$ K Å-3. The constant $`k`$ can also be obtained from the bulk modulus $`B`$. Considering that any deformation of the tube takes place only in a direction perpendicular to its long axis, $`k=2/\sqrt{3}B`$.Tersoff  From the experimental measure of $`B`$ by Tang et al. Tang  ($`B=41.66`$ GPa) one obtains $`k=3015`$ K Å-3, a figure much larger than the one from Ref. Mizel, . Part of the difference can somehow emerge from the fact that the experimental data for $`B`$ was obtained for CN’s with wider radius (7.04 Å). In the results presented below, we have used $`k=1740`$ K Å-3 everywhere but the influence on the final results of the particular $`k`$ value is analyzed in selected cases.

Finally, the term $`h(\lambda ,D)`$ in Eq. (1) sums the interaction energy between the H2 molecules of a given IC and the ones of its neighboring channels. This contribution may be readily obtained assuming a mean field approximation,

$$h(\lambda ,D)=\frac{\lambda }{2}_{\mathrm{}}^{\mathrm{}}𝑑xV_{\text{H}\text{2}\text{-H}\text{2}}\left(\sqrt{x^2+D^2}\right),$$
(7)
which considers the neighboring channels as uniform arrays without correlation effects that could modify the H2 kinetic energy. The calculation in Eq. (7) is extended to the nearest neighbors, the next-nearest neighbors, and then on up to the desired accuracy. Obviously, $`h(\lambda ,D)`$ changes when the bundle swells. In Ref. yo, , we have verified the high accuracy of Eq. (7) by comparing its estimation with an exact DMC calculation.

The influence of the CN radius and the C-H2 potential in the equation of state of adsorbed H2 is clearly observed in the value of $`ϵ_0`$. Table I contains results for the H2 binding energy considering two CN radius, 6.8 and 6.9 Å, the two C-H2 potentials discussed above, LB and WSW, and two geometries for the IC, the real one with the triangular composition of the tubes surrounding the IC (T) and the radial (R) model from an azhimutal average. The narrowness of the IC makes the differences between the T and R results not to be larger than 15 %. The dominant effect is unquestionably the C-H2 potential. The significant differences between the pairs ($`\sigma `$,$`ϵ`$) for the two models generate a big discrepancy in the value of $`ϵ_0`$. The largest binding energies correspond to the most accurate WSW potential due to its larger parameter $`ϵ`$. It is also noticeable, and will be discussed below in connection with dilation, the different behavior of $`ϵ_0`$ with the CN radius: it strongly depends on the radius for LB whereas remains nearly unchanged for WSW. This last feature is a consequence of the larger $`\sigma `$ value of the LB potential which effectively reduces the space available inside the IC to accommodate the H2 molecule.

In order to determine if a dilation of the bundle is energetically preferred, a series of calculations using Eq. (1) has been carried out. As commented before, our aim is to quantitatively determine the influence of the potentials, the CN radius, the elastic constant $`k`$, and the geometry model for the IC on a possible dilation. The direct output of Eq. (1) is the energy per particle as a function of the linear density $`\lambda `$ and $`D`$. For a given dilation $`DD_00`$, with fixed $`D_0=17`$ Å, $`e(\lambda )`$ can be obtained and from it the equilibrium point ($`\lambda _0`$, $`e_0`$) corresponding to zero pressure. Illustrative outputs of this procedure are shown in Figs. 1 and 2. Fig. 1 contains the H2 equation of state in the IC for values of the lattice dilation $`DD_0`$ = 0.18, 0.17, 0.16, and 0.15 Å. All the curves have been calculated assuming a CN radius of 6.9 Å, the SG H2-H2 potential, $`k=1740`$ K Å-3, the LB C-H2 potential, and a T channel. The lowest energy per particle at the minimum is obtained for a dilation 0.16 Å, nearly coincident with the result reported in Ref. Cole1,  (0.166 Å) obtained with the H2-H2 KCLDJ potential and an R channel. If the C-H2 potential is the WSW model, the H2 equation of state changes dramatically and the dilation becomes much smaller. This is shown in Fig. 2 which differs from Fig. 1 in the C-H2 potential, the other inputs being unchanged. The curves now correspond to dilations $`DD_0`$ = 0.06, 0.05, and 0.04 Å. The minimum corresponds to a swelling of 0.05 Å. The equilibrium density is $`\lambda _0=0.232`$ Å-1, slightly inferior to the one obtained with the LB potential, $`\lambda _0=0.268`$ Å-1. The magnitude of the energy per particle is very different in the two cases: going from LB to WSW the energy increases in a factor two. The major part of this increment comes from the differences in the binding energies of a single molecule $`ϵ_0`$ (Table I).

The full set of results combining CN radius, C-H2 and H2-H2 interactions, and IC geometry is reported in Table II. $`\lambda _0`$ corresponds to the H2 equilibrium density at the energetically preferred swelling, reported also in the Table. Due to the large uncertainty in the value of $`k`$, we only report the influence of its change ($`3015`$ K Å-3 instead of $`1740`$ K Å-3) in the most favorable cases for the bundle swelling. Inspection of the data shows that neither the H2-H2 potential nor the detailed form of the IC, R or T, are relevant in the swelling process. Maintaining fixed the radius to a value 6.9 Å and the LB C-H2 interaction, the dilation comes out 0.16-0.17 Å, in overall agreement with Calbi et al.. Cole1  However, if the LB interaction is substituted by the WSW one, the swelling observed is reduced by a factor of three (from 0.17 Å to 0.05 Å), and the equilibrium density $`\lambda _0`$ decreases by 15 %.

The second part of Table II contains a similar analysis for a CN radius of 6.8 Å. A general trend arising from the comparison between the first half of the Table (6.9 Å) and the second one (6.8 Å) is the increase of the dilation with the CN radius when the rest of the inputs are kept fixed. This seems quite obvious since, with a smaller nanotube and the same CN-CN distance, the room for adsorption increases. Thus, the swelling decreases from 0.17 to 0.11 Å with the LB C-H2 interaction whereas it turns absolutely negligible in case of using the probably most accurate WSW potential. On the other hand, $`\lambda _0`$ is observed to systematically decrease when the radius is moved from 6.9 to 6.8 Å. Finally, the influence of the value of $`k`$ appears in the Table marked with an asterisk for a particular case in which the dilation is one of the largest values. As it could be expected, $`(DD_0)`$ decreases if $`k`$ increases: using $`k=3015`$ K Å-3 instead of 1740 K Å-3 the dilation is reduced by 30 %. Similar reductions would be obtained in the other cases.

## III FULL DMC CALCULATION

In the preceding Section it has been assumed that radial and longitudinal degrees of freedom of H2 inside IC’s are not coupled. The correlations between H2 molecules were purely longitudinal and the interaction with the surrounding walls was solved only for the one-body problem. In this Section, we check the validity of that approximation by making an exact three-dimensional DMC calculation of H2 in the IC. The CN bundle cohesion term and the mean field contribution $`h(\lambda ,D)`$ are summed up to the DMC energy as in Eq. (1).

The trial wave function for importance sampling is written as

$$\psi (𝑹)=\psi ^{\text{1D}}(𝑹)\psi ^{\text{IC}}(𝑹),$$
(8)
where $`\psi ^{\text{1D}}(𝑹)`$ and $`\psi ^{\text{IC}}(𝑹)`$ are the same than in the 1D calculation, Eqs. (3) and (4), respectively. The optimal values for the parameters entering Eq. (8) are the same than the ones in the preceding Section. The DMC calculation has been carried out for a selected case which we consider contains the most reasonable model. The CN radius is 6.8 Å and the IC is of T type. The H2-H2 interaction is SG and the parameter set for C-H2 is the WSW one.

In Fig. 3, the total energy per particle of H2 in the IC is shown as a function of $`\lambda `$ and for different bundle lattice values, $`DD_0`$ = 0, 0.05, 0.10, and 0.15 Å. The minimum of the energy is achieved with $`D=D_0`$, what implies an absence of dilation in the bundle of nanotubes. This is in agreement with the corresponding case in Table II, marked by a dagger. In the full 3D calculation the equilibrium density is $`\lambda _0=0.2184(12)`$ Å-1, to be compared with $`0.2194(5)`$ Å-1 from Table II. The second value is within the error bar of the first one. The energy per hydrogen molecule at the equilibrium density is also very similar in both cases, $`e_0=1101.43(1)`$ K and $`e_0^{\text{1D}}=1102.31(2)`$ K. This clear overlap means that the 1D model (1) is a very good approximation to the physics of H2 inside the interchannels of a bundle of nanotubes around the equilibrium density.

Finally, in Fig.4 the comparison between both calculations is extended to higher densities. In the energy scale, the binding energy $`ϵ_0=1096.02`$ K (Table I) is now subtracted. Equation (1) provides un upper bound to the exact result which is very close to the exact energies near $`\lambda _0`$. However, its quality worsens when the density increases due to the emergence of the radial degree of freedom which is frozen in the 1D model.

## IV SUMMARY AND CONCLUSIONS

The equation of state of H2 adsorbed in the IC’s of a CN bundle has been calculated using microscopic theory. By means of the DMC method, which provides exact results for given model potentials, the possible dilation of the bundle induced by adsorption has been carefully analyzed. The most complete analysis has been done in the 1D model.Cole1  In spite of its simplicity, this approach (1) has proven to be very accurate when compared with a full 3D calculation, specially near the equilibrium density. Playing with the different alternatives for the parameters entering into the calculation the influence of each one has been established. Summarizing, the results show that the critical one is the particular C-H2 interaction. The deeper and probably more realistic WSW potential reduces the swelling predicted by the LB model in a significant way. The same can be concluded about the C-C interaction, represented here by the $`k`$ parameter, but to a lesser extent. On the other hand, sizeable changes in the H2-H2 interaction do not modify the physics of the dilation process. Our results do not exclude completely the possibility of some lattice dilation but if it exists its value seems clearly smaller than the figures quoted in Ref. Cole1, .

A key point in the discussion is the issue of what C-H2 potential is more realistic. To this end, one should compare theoretical data with experimental results. However, the latter are really scarce. Recently, Vilches and collaborators vilches1 ; vilches2  reported data on H2 isotherms on bundles of closed-capped carbon nanotubes for several temperatures. At very low coverages, they deduced the existence of H2 adsorption on the ridges created between two nanotubes in the outer part of the bundle and/or in the IC’s. For these sites, they estimated a binding energy $``$ 700 K, 1.5 times larger than the one for molecular hydrogen on graphite. Comparison with Table I rules out the case with R= 6.9 Å and the LB potential: its $`ϵ_0`$ is too small. On the other hand, the value for R = 6.8 Å and the LB potential seems to fit perfectly well the experimental result. However, there is a problem: the binding energy 700 K is supposed to be an average of the energies of hydrogen adsorbed on the ridges and in the IC’s. Then, the real binding energy in IC’s should be greater than the experimental finding, fitting well with the results obtained with the WSW potential.

###### Acknowledgements.

We are grateful to Ari Mizel, Oscar E. Vilches, and Milton W. Cole for helpful discussions. This work was supported by DGI (Spain) Grant N0 BFM2002-00466 and BQU2001-3615-C02-01, and Generalitat de Catalunya Grant N0 2001SGR-00222.