text:
# 1 Introduction

## 1 Introduction

Despite the success of the Standard Model (SM), it is far from complete. For instance, the origin of the fermion mass hierarchy and issues such as baryon asymmetry are yet to be resolved. We would expect new physics to contribute to these areas. At the more fundamental level, nobody knows of any deep explanation as to why there should only be three generations of quarks & leptons. For these reasons, and others, extensions of the three-generation SM are being investigated. Here, we investigate the possibility of incorporating a fourth-generation into a GUT model.
As is well known, all fourth-generation models must adhere to certain experimental constraints, the first of which stems from precise measurements of the decay characteristics of the $`Z`$-boson performed at LEP. This has set a lower bound of $`M_F\frac{M_Z}{2}`$ on any non-SM particles that couple to the $`Z`$-boson. Ignoring for the moment the unnatural hierarchy emerging within the neutrino sector, we assume a Dirac mass $`M_N\left(\frac{M_Z}{2}\right)`$ for the heavy neutral lepton. We label the fourth-generation explicitly with the notation

$$Q_4=\left(\begin{array}{c}T\\ B\end{array}\right)_4T_4^cB_4^c;L_4=\left(\begin{array}{c}N\\ E\end{array}\right)_4N_4^cE_4^c$$
(1)
Physics beyond the Standard Model is severely constrained by precision electroweak data111We assume $`\left|V_{tb}\right|^2`$, $`\left|V_{TB}\right|^21`$ and $`\left|V_{Tb}\right|^2`$, $`\left|V_{tB}\right|^21`$ so that contributions to the $`Zb\overline{b}`$ decay can be ignored.. Assuming SM contributions, fits to LEP data give the radiative correction parameter222S is the well known radiative correction parameter (weak isospin symmetric), normalized to zero for the SM with $`M_{Higgs}=`$ 100 GeV. $`S=0.04\pm 0.11`$. For a heavy ($`M_Z`$) degenerate fourth-generation we obtain $`\mathrm{\Delta }S=\frac{2}{3\pi }`$ and so is ruled out at $`99.8\%`$ $`C.L`$. However, Maltoni et al.  have shown that particles with mass $`M\frac{M_Z}{2}`$ give drastically different contributions to $`S`$. It is possible for a heavy neutrino $`N`$ with mass $`M_N\frac{M_Z}{2}`$ to cancel the contributions from the heavy $`T,B,E`$ with the SM solution333Throughout this paper, we shall employ the notation that upper case letters ($`M_x`$) denote pole masses and lower case letters ($`m_x`$) denote running masses.

$$M_E>M_N;M_T>M_B$$
(2)
$$M_N\frac{M_Z}{2};M_B=M_B^{min}$$
(3)
$$(M_EM_N)3(M_TM_B)60GeV$$
(4)
where the $`N`$ must be relatively stable to avoid detection (i.e. mixing matrix elements $`V_{Ne,\mu ,\tau }<10^6`$). In Eq. (3), $`M_B^{min}`$ is the minimum mass of the B quark as allowed by experimental searches. An extra generation can be accomodated below the $`1\sigma `$ level (or even two generations at $`1.5\sigma `$)444Updated fits to recent LEP data show that an extra generation with $`M_N=M_U`$ 180 GeV; $`M_D=M_E`$ 130 GeV gives a $`\chi ^2`$ fit for the fourth-generation as good as that of the three-generation case . Perturbative unification is not consistent with such large masses.. However, it was shown by Gunion et al.  that the fourth-generation charged and neutral leptons must be relatively light ($`M_{E,N}\frac{M_Z}{2}`$) in order to stay in the perturbative regime below the GUT scale. Although they worked within a supersymmetric framework, similar results are expected to hold in the SM . Setting $`M_E\frac{M_Z}{2}`$ to achieve perturbative unification is not consistent with the above SM fits (Eq. (4)) to the precision data. Therefore, we require the cancellations to arise from another sector and so we consider a supersymmetric theory. Indeed, it is shown in  that light supersymmetric particles can also effect the fit to precision data, allowing for a fourth-generation below the $`2\sigma `$ level. In particular, neutralinos ($`\stackrel{~}{\chi }_1^0`$) and charginos ($`\stackrel{~}{\chi }_1^\pm `$) with masses less than $`M_Z`$ can provide the correct sign contributions, whilst being consistent with current experimental limits if nearly degenerate 555 $`\stackrel{~}{\chi }_1^0`$ and $`\stackrel{~}{\chi }_1^\pm `$ are defined as the lightest mass eigenstates in the neutralino and chargino sectors respectively..
In this paper we investigate the possibility of consistently incorporating a fourth-generation into a $`N=1`$ R-parity conserving supergravity model. We assume a structure akin to that of the minimal supersymmetric standard model (MSSM3), adding a complete chiral fourth-generation and its associated SUSY partners (the so called MSSM4). Specifically we require $`(i)`$ perturbative values for all Yukawa couplings at energies up to the GUT scale and $`(ii)`$ gauge coupling unification. These two constraints will be termed collectively as perturbative unification.
Having satisfied precision data fits, it remains for us to provide solutions that $`(i)`$ unify perturbatively at the GUT scale and $`(ii)`$ evade the direct experimental searches performed at $`CDF`$. In section 2 we investigate the specific decay channels of the fourth-generation quarks and ensure we can provide consistency with the experimental direct searches. The solution we present requires the introduction of light SUSY particles ($`\stackrel{~}{\chi }_1^0,\stackrel{~}{\chi }_1^\pm ,\stackrel{~}{B}`$), so that the two-body decays $`T\stackrel{~}{B}\stackrel{~}{\chi }_1^+`$ and $`B\stackrel{~}{B}\stackrel{~}{\chi }_1^0`$ are kinematically allowed. These will always dominate over the one-loop FCNC decays $`BbZ`$ and two-generation decays $`BcW^{}`$ that traditional searches have looked for, whilst at the same time suppress the decay rate $`TbW^+`$. We note a light, degenerate chargino/neutralino pair is just what is needed to provide the necessary cancellations in precision data.
In section 3 we discuss the influence of specific fourth-generation and SUSY particle masses on precision data fits in more detail and point out new features that appear when trying to satisfy all constraints. In section 4 we proceed with a renormalisation group study of the four generation MSSM. Under the assumption of a common mass scale for the coloured/weak sparticle spectrum ($`M_{col}`$/$`M_{wk}`$ respectively) we will derive upper bounds on $`M_{col}`$. This bound stems from the fact that if we are to achieve perturbative values for the top quark Yukawa coupling constant $`h_t`$ up to the GUT scale, then the effects of the coloured sparticle spectrum must be included in the running of the strong coupling $`\alpha _3`$ at an early stage. Finally, in section 5 we investigate the lightest Higgs mass in the MSSM4.

## 2 Experiment and a Fourth Generation

To reduce the number of free parameters we take a common mass scale for the fourth-generation quarks $`M_Q`$ and for the leptons $`M_L`$. It is understood, however, that in reality a small non-degeneracy $`\mathrm{\Delta }M_Q=M_TM_B`$; $`\mathrm{\Delta }M_L=M_EM_N>0`$ exists as will be emphasized in the following text. We begin our discussion by considering the leptonic sector. Under the assumption that the mixing between ($`E,N`$) and the first three generations is negligible, the decay $`ENW^{}`$ will be dominant. Current experimental limits searching for $`ENW^{}`$ from $`e^+e^{}E^+E^{}`$ production have been performed by the OPAL/L3 collaborations up to the kinematic limit $`M_E`$ 100 GeV . However, to be consistent with perturbative unification we require $`M_E\frac{M_Z}{2}`$. Therefore, in order to evade experimental bounds, the mass difference $`\mathrm{\Delta }M_L`$ must be less than $``$ 5 GeV so as to result in a virtual $`W^{}`$ whose decay products are too soft to be triggered. Regarding the heavy neutrino, OPAL/L3 have set the bound $`M_N>`$ 70 - 90 GeV based on the search for $`NlW^{}`$ ($`l=e,\mu `$ or $`\tau `$) provided the mixing matrix elements satisfy $`V_{Ne,\mu ,\tau }>10^6`$ . However, we have already assumed that the mixing of the fourth generation leptonic sector is in fact negligible ($`V_{Ne,\mu ,\tau }<10^6`$); so the neutrino is stable enough to leave the detector and, in this case, only the DELPHI bound of $`M_N>`$ 45 GeV from measurements of the $`Z`$-width applies. For the purpose of the numerical study we take $`M_L=`$ 46 GeV.
Turning to the quark sector, the requirement of perturbative unification places strict upper limits on the masses of the $`T`$ and $`B`$ and they certainly must be below the top quark whose mass is $`M_t=174\pm `$ 5.1 GeV. Experimental searches for the fourth-generation quarks have mainly concentrated on the $`B`$-quark where CDF search for

$$e^+e^{}B\overline{B}b\overline{b}ZZ(b\overline{b})(l^+l^{})(q\overline{q})$$
(5)
The current bound excludes $`M_B<`$ 199 GeV assuming the branching ratio $`BR(BbZ)1`$. The search is also sensitive to other decay modes as long as $`BR(BbZ)`$ is large enough to provide the leptonic trigger. For instance, the decays $`Bbh`$ and $`BcW^{}`$ will be triggered, since the hadronic decay of the $`h`$ and $`W`$ is kinematically similar to those of the $`Z`$. One possible escape might come about if we note that $`hN\overline{N}`$ would compete with $`hb\overline{b}`$ for $`M_hM_Z`$ (dominating for larger $`M_h`$) and would provide an invisible signature. For $`M_hM_Z`$ we would also expect $`BR(Bbh)BR(BbZ)`$ . However, CDF still exclude a $`B`$-mass in the range 104 GeV $``$ 152 GeV, assuming $`BR(BbZ)50\%`$ and no sensitivity to the other decay modes . For the region $`M_B<`$ 104 GeV, we require $`BR(h\text{invisible})>`$ 30% in order to evade the B search. This comes from the fact that $`BR(Zq\overline{q})=`$ 70% and that we also assume the acceptance for the two processes

$$e^+e^{}B\overline{B}b\overline{b}Z\{Z\text{or}h\}(b\overline{b})(l^+l^{})\{q\overline{q}\text{or}b\overline{b}\}$$
(6)
are equal. However, this scenario is inconsistent with the latest bounds from LEPII  which exclude $`M_h<`$ 106 GeV if $`BR(h\text{invisible})>`$ 30%. Moreover, we have yet to account for the $`T`$ decays which turn out to be highly constrained if we assume $`SM`$-like processes. The channel $`TbW^+`$ is prohibited for the obvious reason that the $`T`$-quark would have been picked up in the $`CDF`$ search for the top quark. Although we might assume that $`TBW^{}`$ is dominant by suppressing $`V_{Tb}`$ we must notice that, since the $`T\overline{T}`$ production cross-section is similar to that of $`B\overline{B}`$, we would effectively expect double the event rate on the $`B`$-quark search. Taking this into account would further strengthen existing bounds on the $`B`$-quark mass.
From the ideas presented so far, we might conclude that the fourth-generation with perturbative unification is not consistent with experimental bounds on the ($`T,B`$) quarks. We have not, however, considered the possibility of light SUSY particles providing decay channels for ($`T,B`$). In this situation one can constrain the light (i.e. $`<M_Z`$) neutralino/chargino pair, which is already required by fits to precision data, in order to provide the following two-body decays

$$T\stackrel{~}{B}_1\stackrel{~}{\chi }_1^\pm ;B\stackrel{~}{B}_1\stackrel{~}{\chi }_1^0$$
(7)
where $`\stackrel{~}{B}_1`$ is the lightest mass eigenstate, defined by the amount of mixing between the L/R weak eigenstates; $`\stackrel{~}{B}_1=\stackrel{~}{B}_L\mathrm{cos}\theta _{\stackrel{~}{B}}+\stackrel{~}{B}_R\mathrm{sin}\theta _{\stackrel{~}{B}}`$. Ensuring both decay channels are kinematically accessible, combined with the fact that perturbative unification requires that $`M_QM_Q^{max}`$ 108 GeV666$`M_Q^{max}`$ is defined as the maximum allowed value of $`M_Q`$. This is obtained by ensuring perturbative consistency in renormalization group calculations (see section 4). (see Table 1), places severe restrictions on the allowed spectrum. In the next section we combine the experimental constraints with those obtained from precision fits to obtain the allowed fourth-generation spectrum.

## 3 Precision Measurements and a Fourth Generation

It is difficult to provide bounds from precision data without a fully consistent study taking into account exact particle masses, any light SUSY spectra present and mixings between different flavours. However, Maltoni et al.  pointed out that a highly degenerate neutralino/chargino pair can provide the necessary contributions needed to cancel that of the fourth-generation, whilst at the same time being consistent with LEP bounds. Specifically they require

$`M_{\pi ^+}\stackrel{<}{}\mathrm{\Delta }M_{\stackrel{~}{\chi }}\stackrel{<}{}3GeVM_{\stackrel{~}{\chi }}`$ $``$ $`60GeV`$ (8)
where we define the notation $`\mathrm{\Delta }M_{\stackrel{~}{\chi }}=M_{\stackrel{~}{\chi }_1^\pm }M_{\stackrel{~}{\chi }_1^0}`$; $`M_{\stackrel{~}{\chi }}=M_{\stackrel{~}{\chi }_1^0}M_{\stackrel{~}{\chi }_1^\pm }`$. Looking at their results, we see that the magnitude of the contribution to the fitted parameters from this sector are highly dependent on $`M_{\stackrel{~}{\chi }}`$. Deviations from $`M_{\stackrel{~}{\chi }}`$ 60 GeV larger than + 30 GeV / - 5 GeV are ruled out at the $`2\sigma `$ level. This scenario is inconsistent with perturbative unification, since $`M_{\stackrel{~}{\chi }}`$ 60 GeV with $`M_{\stackrel{~}{B}_1}>M_{\stackrel{~}{\chi }}`$ would require $`M_Q>`$ 120 GeV in order to retain the decay channels in Eq. (7)777We are assuming here that the lightest supersymmetric particle is the neutralino $`\stackrel{~}{\chi }_1^0`$.. One possible solution presents itself, however, if we assume that the lightest supersymmetric particle (LSP) is in fact the fourth-generation sneutrino $`\stackrel{~}{N}`$ which is stable due to R parity. Assuming the leptonic mixing angle between the third- and fourth-generation is non-zero, the $`\stackrel{~}{\chi }_1^0`$ would decay invisibly via $`\stackrel{~}{\chi }_1^0\upsilon _\tau \stackrel{~}{N}`$ . The following masses are possible:
$`M_{\stackrel{~}{\chi }}`$ $``$ $`60GeV`$    $`M_{\stackrel{~}{B}_1}`$ $``$ $`45GeV`$    $`M_{\stackrel{~}{N}_1}`$ $`<`$ $`M_{\stackrel{~}{B}_1}`$        $`M_Q`$ $`=`$ $`106GeV`$    $`M_L`$ $`=`$ $`46GeV`$  (9)   $`M_t`$ $`=`$ $`170GeV`$
where we assume that the mixing angles $`\theta _{\stackrel{~}{B}/\stackrel{~}{N}}`$ are such that $`\stackrel{~}{B}_1/\stackrel{~}{N}_1`$ decouple from the Z-boson at tree-level, and therefore do not contribute to its total width. The fourth-generation bottom squark will decay via $`\stackrel{~}{B}_1c\tau \stackrel{~}{N}_1`$. Such a decay involves the factor $`|V_{Bc}V_{N\tau }|^2`$, leading to a long lifetime. Current bounds looking for hadronizing sbottom quarks exclude the range 5 GeV $`M_{\stackrel{~}{B}_1}`$ 38 GeV, if the mixing angle $`\theta _{\stackrel{~}{B}}=`$ 1.17 rad where the production cross-section is minimized . We note this corresponds to the mixing angle at which $`\stackrel{~}{B}_1`$ decouples from the Z-boson, as required from the Z-width constraint.
Although the sparticle masses seem to be contrived, we note that they can be obtained from reasonable assumptions about the supersymmetric sector. In order to obtain $`M_{\pi ^+}\stackrel{<}{}\mathrm{\Delta }M_{\stackrel{~}{\chi }}\stackrel{<}{}`$ 3 GeV, we require the hierarchy

$$|\mu |M_1M_2$$
(10)
where $`|\mu |`$ is the Higgs mixing parameter. $`M_1`$ and $`M_2`$ are the electroweak gaugino masses (for a review of supersymmetry see ). This structure can occur naturally within specific supergravity models . The final mass eigenstates $`M_{\stackrel{~}{X}_{1,2}}`$ ($`X=B,N`$) depend on the appropriate squark/slepton soft mass terms and scalar tri-linear couplings that are present in the soft supersymmetric breaking sector, and the $`\mu `$ term. Here, we assume it is possible to constrain these free parameters to give the required spectrum without contradicting the model in any way. However, more detailed studies taking into account the high-energy behaviour and renormalization group improvement needs to be performed, in order to establish consistency with radiative electroweak breaking and other constraints such as the absence of charge/colour breaking etc.
We note that the fits to precision data performed in  constrained the fourth-generation masses to satisfy $`M_N=M_T`$ and $`M_E=M_B`$, whereas we consider the case $`M_TM_B>M_NM_E`$. From Eqs. (5.9–5.11, 5.13) in , we can see that as long as we satisfy the conditions $`|M_{T,N}M_{B,E}|M_Z`$, the two scenarios differ only in the universal contributions to the radiative parameters (see Maltoni, Ph.D. thesis  for terminology). In particular, we would require an extra (positive) universal contribution in order to exactly reproduce the fits. This could arise if we were to find large SU(2) breaking within the $`(\stackrel{~}{t},\stackrel{~}{b})`$ sector. This has been shown to provide the correct type of contribution that is needed to compensate between the two scenarios .

## 4 Renormalization Group Study of the MSSM4

Here we investigate the effect of the fourth-generation on the evolution of couplings to the GUT scale, where we require gauge coupling unification. This places upper limits on the masses of the extra particles to ensure their Yukawa couplings run perturbatively to the unification scale $`M_U`$ ($`h^2(\mu )4\pi `$, $`M_Z\mu M_U`$). Starting at the low-energy scale $`M_Z`$, the electroweak gauge couplings $`\alpha _1(M_Z)`$, $`\alpha _2(M_Z)`$ are fixed through the relations $`\frac{1}{\alpha _i(M_Z)}=\frac{3}{5}\frac{\left(1\mathrm{sin}^2\theta _W\right)}{\alpha _{em}(M_Z)};\frac{\mathrm{sin}^2\theta _W}{\alpha _{em}(M_Z)}`$ for $`i=1,2`$ respectively. The strong coupling $`\alpha _3(M_Z)`$ is taken from the Particle Data Group (PDG)  to be $`0.1181\pm 0.002`$. We take the best fit values for $`\alpha _{em}(M_Z)`$/$`\mathrm{sin}\theta _W`$ from the PDG. In principle, one should extract the $`Z`$-pole couplings assuming the full MSSM4, thereby accounting for the fourth-generation fermions and light SUSY spectra in a fully consistent way888We account for the one-loop leading logarithmic corrections from the SUSY sector when running the RGE by employing the step-function approach . This procedure is accurate in the limit of heavy sparticles but fails for masses $`\stackrel{~}{M}<M_Z`$ where both logarithmic and finite corrections will influence the extraction of the couplings.. We have also performed our study in the $`\overline{MS}`$ scheme although it is the $`\overline{DR}`$ scheme that is consistent with supersymmetry. However, differences between the $`\overline{DR}`$ and $`\overline{MS}`$ schemes are not significant at the low-energy scale as compared to other uncertainties .
In our analysis we neglect all Yukawa couplings from the first three generations, except that of the $`t`$-quark whose mass we take to be $`M_t=`$ 170 GeV. As is typical with four-generation models, we require small values of $`\mathrm{tan}\beta `$ (the ratio of the Higgs vev’s) so as to avoid $`h_B(M_Z)`$$`𝒪`$($`\sqrt{4\pi }`$. Once all couplings at $`M_Z`$ have been fixed, we integrate up in energy scale using the two-loop renormalization group equations (RGE). The one-loop leading logarithmic threshold corrections from the SUSY sector are accounted for in the numerical procedure. We select the point where $`\alpha _1(\mu )=\alpha _2(\mu )`$ as the unification scale $`M_U`$ with coupling $`\alpha _U(M_U)`$. Any deviation in $`\alpha _3(M_U)=\alpha _U(M_U)`$, which we parameterize as $`\delta =\frac{\alpha _3(M_U)\alpha _U(M_U)}{\alpha _U(M_U)}`$, can arise from either of two sectors. On the one hand, we have the $`\overline{MS}`$ vs. $`\overline{DR}`$ mismatch, experimental errors in $`\alpha _{em}(M_Z)`$/$`\mathrm{sin}\theta _W`$ and the variations in the best fit values of $`\alpha _{em}(M_Z)`$ and $`\mathrm{sin}\theta _W`$ as the fourth-generation and light SUSY particles ($`\stackrel{~}{M}<M_Z`$) are included. However, more importantly, assuming no intermediate scales, high-energy threshold corrections from specific GUT/string models can provide contributions to $`\delta `$. Following , we note that these corrections (for particular models) can be large. We conservatively assume that unification is consistent if $`|\delta |`$ 5 %.
The SUSY threshold corrections in the MSSM4 are important as they can influence whether or not a particular set of parameters ($`M_t,M_T,M_B,M_N,M_E,\mathrm{tan}\beta `$) will retain perturbative consistency to the GUT scale. This can be observed analytically if we write the one-loop leading logarithmic correction to the strong coupling from the SUSY sector

$$\frac{1}{\alpha _3^+(M_Z)}\frac{1}{\alpha _3^{}(M_Z)}=\frac{b_3^{MSSM4}b_3^{SM4}}{2\pi }\mathrm{ln}\left(\frac{M_{col}}{M_Z}\right)$$
(11)
Here $`M_{col}`$ represents an effective threshold scale and $`b_3^{MSSM4/SM4}`$ represents the one-loop beta function contribution to the strong coupling in the MSSM4/SM4 respectively. This correction is implemented at the scale $`M_Z`$ and accounts for the coloured sparticles with masses $`M_{col}>M_Z`$ in the running of the strong coupling. $`\alpha ^+`$/$`\alpha ^{}`$ represents the renormalized gauge coupling just above/below the scale $`M_Z`$ and $`\alpha _3^{}(M_Z)`$ is fixed to be 0.118 from experimental input. From Eq. (11) we see that the higher the scale $`M_{col}`$, the lower $`\alpha _3^+(M_Z)`$. Writing the one-loop RGE for the third-generation top Yukawa coupling $`h_t`$, accounting for the $`t,T,B`$ quarks and the strong coupling, we have

$`{\displaystyle \frac{1}{y_t}}{\displaystyle \frac{dy_t}{d\mathrm{\Lambda }}}`$ $`=`$ $`{\displaystyle \underset{\delta =t,T,B}{}}a_{t\delta }y_\delta {\displaystyle \frac{16}{3}}\alpha _3`$ (12)
$`a_{t\delta }`$ $`=`$ $`(6,3,1)`$
where $`y_\delta =\frac{h_\delta ^2}{4\pi }`$ and $`\mathrm{\Lambda }=\frac{\mathrm{ln}\mu }{2\pi }`$. We can see that the corrections from the sparticle sector (see Eq. (11)) are constrained, since we require a large $`\alpha _3^+(M_Z)`$ in order to compensate for the contributions from $`h_{T,B}`$ in Eq. (12) that are driving $`h_t`$ to non-perturbative values. In particular, the Yukawa couplings $`h_{t,T,B}`$ are fixed at the low-energy scale through the relations

$$h_{t,T}(m_{t,T})=\frac{m_{t,T}(m_{t,T})}{\upsilon \mathrm{sin}\beta }h_B(m_B)=\frac{m_B(m_B)}{\upsilon \mathrm{cos}\beta }$$
(13)
where $`\upsilon =`$ 174 GeV is the Higgs vacuum expectation value and $`m_{t,T,B}(m_{t,T,B})`$ are the quark running masses at the scale $`m_{t,T,B}`$. For a given $`\mathrm{tan}\beta `$, if $`M_{t,T,B}`$ is sufficiently large so that the RG running of $`h_t`$ is tending towards non-perturbative values, then we obtain an upper bound $`M_{col}^{max}`$, above which the initial $`\alpha _3^+(M_Z)`$ is too small to counteract the effect of the fourth-generation couplings when solving the RGE in Eq. (12).
In practice, we perform a full numerical study, accounting for the threshold corrections from $`t,T,B`$ and all sparticles by changing the $`\beta `$-functions and using the step-function approach in the running of the gauge couplings . We assume separate degeneracies amongst the coloured and weak SUSY spectrum

$`M_{wk}`$ $`=`$ $`M_{\stackrel{~}{L}}=M_{\stackrel{~}{l}}=M_{\stackrel{~}{H}}=M_{\stackrel{~}{W}}=M_H`$ (14)
$`M_{col}`$ $`=`$ $`M_{\stackrel{~}{Q}}=M_{\stackrel{~}{q}}=M_{\stackrel{~}{g}}`$ (15)
where $`l`$ and $`q`$ denote the leptons and quarks of the first three-generations respectively. We assume 100 GeV $`M_{col}`$ 2 TeV and we restrict the values of $`M_{wk}`$ to be 500 GeV, 1 TeV or 2 TeV. Of course, considering separate degeneracies amongst the coloured and weak sparticle spectra is an approximation. In fact, our model demands the introduction of light sparticles to ensure consistency with experimental searches for the fourth-generation and precision data bounds, thus providing significant deviations from degeneracy. Nevertheless, looking at Figure 1, we can clearly see that the effective mass ($`M_{col}`$) of the coloured spectrum is bounded from above, in the process of retaining perturbative consistency. The effect of $`M_{wk}`$ is clearly weak since $`\alpha _1(M_Z)`$ and $`\alpha _2(M_Z)`$ are an order of magnitude lower than the strong coupling $`\alpha _3(M_Z)`$. In particular, although we expect the sparticles to deviate in mass around the effective scales $`M_{col}`$ and $`M_{wk}`$, we can still conclude that as long as all the sparticles have masses less than 2 TeV, then perturbative unification is consistent with $`(M_t,M_Q,M_L)=`$ (170 GeV, 106 GeV, 46 GeV) for $`1.63\mathrm{tan}\beta 1.67`$.

## 5  The Higgs Sector of the MSSM4

In this section, the effect of the fourth-generation on the lightest Higgs mass ($`M_h`$) is investigated. We shall employ the one-loop effective potential with contributions from the top / stop and fourth-generation fermions / sfermions. Since we are at low $`\mathrm{tan}\beta `$, we can ignore the contributions from the third generation bottom / sbottom masses. At tree-level, the Higgs sector of minimal supersymmetry models is fixed by two parameters $`(i)\mathrm{tan}\beta `$ and $`(ii)M_A`$ (the pseudoscalar Higgs-boson mass). We set $`\mathrm{tan}\beta =`$ 1.64 from the requirement of perturbative unification (see Figure 1). Assuming the pseudoscalar Higgs-boson mass ($`M_A`$) is degenerate with the heavy Higgs sector ($`M_H`$) we set $`M_A=`$ 500 GeV. Accounting for one-loop terms, the Higgs mass $`M_h`$ becomes dependent on ($`M_{t,\stackrel{~}{t}}`$; $`M_{T,\stackrel{~}{T}}`$; $`M_{B,\stackrel{~}{B}}`$; $`M_{N,\stackrel{~}{N}}`$; $`M_{E,\stackrel{~}{E}}`$). In general, the sparticle mass matrix squared in the L/R basis can be written as

$$M_{\stackrel{~}{i,j}}^2=\left(\begin{array}{cc}M_{S_{i,j}}^2+M_{i,j}^2+\mathrm{\Delta }_{D_{i,j}}& M_{i,j}\left(A_{i,j}+\mu \mathrm{cot}\beta \right)\\ M_{i,j}\left(A_{i,j}+\mu \mathrm{cot}\beta \right)& M_{\overline{S}_{i,j}}^2+M_{i,j}^2+\mathrm{\Delta }_{\overline{D}_{i,j}}\end{array}\right)$$
(16)
where $`\mathrm{\Delta }_{D_{i,j}}`$, $`\mathrm{\Delta }_{\overline{D}_{i,j}}`$ represents the $`D`$-term contributions; $`i=t,T,B`$ and $`j=E,N`$. Since large values of $`|\mu |`$ are required by radiative electroweak breaking and also to preserve the hierarchy in Eq. (10), we take $`|\mu |=`$ 500 GeV. The sparticle mass eigenvalues squared, denoted by $`(M_{\stackrel{~}{i,j}})_{1,2}^2`$, are obtained by diagonalizing the matrix in Eq. (16). Since the Higgs mass is dependent on the mixing in the sparticle sector, it is insufficient to consider a fully degenerate SUSY spectrum as in the approximation used in section (4). To account for this, we set the soft supersymmetric breaking terms for the coloured squarks to be $`M_{S_i}=M_{\overline{S}_i}=M_{col}^{soft}`$ which we vary over the range 100 GeV $`M_{col}^{soft}`$ 1 TeV. For the sleptons we take $`M_{S_j}=M_{\overline{S}_j}=M_{wk}^{soft}=`$ 500 GeV. We then randomly vary the two mixing parameters $`A_i(=A_{t,T,B})`$; $`A_j(=A_{E,N})`$ and retain the maximum value returned for the lightest Higgs mass $`M_h^{max}`$. As shown in previous studies of the Higgs sector, this will occur for $`X_{i,j}\pm \sqrt{6}M_{S_{i,j}}`$ where $`X_{i,j}=A_{i,j}+\mu \mathrm{cot}\beta `$. Large mixings will induce light sparticles ($`<M_Z`$) which are required in the MSSM4 to provide the decay channels to the fourth-generation quarks. The minimum Higgs mass is obtained when $`X_{i,j}`$ 0.
In the MSSM4, for $`M_h`$ 100 GeV, the channel ($`h\text{invisible}`$), where invisible represents ($`N\overline{N},E\overline{E},\stackrel{~}{N}\overline{\stackrel{~}{N}},\stackrel{~}{\chi }\overline{\stackrel{~}{\chi }}`$), will open and would dominate over conventional $`hb\overline{b}`$ rates. Exclusion limits will now come from the missing energy search ($`e^+e^{}ZhZ+\text{missing energy}`$) that currently sets the lower bound at 114.4 GeV  assuming a SM-like Higgs and $`BR(h\text{invisible})=`$ 1. We therefore need to check, for our constrained set of parameters ($`M_{t,\stackrel{~}{t}}`$; $`M_{T,\stackrel{~}{T}}`$; $`M_{B,\stackrel{~}{B}}`$; $`M_{N,\stackrel{~}{N}}`$; $`M_{E,\stackrel{~}{E}}`$; $`\mathrm{tan}\beta `$), that $`M_h^{max}`$ lies above this experimental lower bound.
We find a large region 370 GeV $`M_{col}^{soft}`$ 1000 GeV where $`M_h^{max}>`$ 114.4 GeV, although some mixing in the third-generation stop and/or fourth-generation squark sector is required. The region where $`\mathrm{tan}\beta 1.64`$, although excluded in the MSSM3 scenario from experimental lower limits , is allowed in the MSSM4 due to the extra one-loop contributions that arise from $`T,\stackrel{~}{T}`$; $`B,\stackrel{~}{B}`$; $`E,\stackrel{~}{E}`$ and $`N,\stackrel{~}{N}`$. The lightest Higgs-boson would then decay invisibly with an upper bound on its mass of $`M_h`$ 140 GeV 999We have not accounted for two-loop terms that have been shown to correct the Higgs-mass $`M_h`$ by up to -10 GeV in the MSSM3 ..

## 6 Conclusions

We have seen that it is possible to incorporate a fourth-generation into a GUT model, requiring the existence of supersymmetric particles in order to provide the necessary cancellations in precision fits. The fourth-generation masses are then tightly constrained with typical values of $`M_Q=`$ 106 GeV; $`M_L=`$ 46 GeV. To retain perturbative consistency to the unification scale, assuming all sparticles have masses less than 2 TeV, we constrain $`\mathrm{tan}\beta `$ to lie in the interval 1.63 $`\mathrm{tan}\beta `$ 1.67. We find $`M_{col}^{soft}`$ 370 GeV 101010$`M_{col}^{soft}`$ represents the scale of the soft breaking parameters in the $`\stackrel{~}{t},\stackrel{~}{T},\stackrel{~}{B}`$ sector. from the requirement that the lightest Higgs mass must be above its present experimental lower bound, although some mixing between the L/R weak eigenstates in the third and/or fourth-generation squark sector is required. In order to provide decay channels to the fourth-generation quarks, it might be that the LSP is the sneutrino $`\stackrel{~}{N}_1`$ with a mass $`M_{\stackrel{~}{N}_1}`$ 40 GeV. This has implications for dark matter constraints. We would also expect a light $`\stackrel{~}{B}_1`$ with $`M_{\stackrel{~}{B}_1}`$ 40 - 50 GeV along with a degenerate neutralino/chargino pair $`M_{\stackrel{~}{\chi }}=`$ 55 - 65 GeV. The SUSY spectrum needed to satisfy all these constraints cannot be obtained from MSUGRA scenarios with universal parameters at the unification scale.

## Acknowledgements

We would like to thank David Sutherland for discussions and useful remarks.


format: text

text:
# Helium-Star Mass Loss and its Implications for Black-Hole Formation and Supernova Progenitors

## 1 Introduction

A helium star is the naked core of a star that has lost its H-rich envelope, as a result of either a strong stellar wind or binary interaction. In a very massive star (initial mass $`M_\mathrm{i}\text{ }>40`$$`\mathrm{M}_{}`$) the stellar wind is strong enough to remove the envelope before or during the central He-burning phase of evolution. Such stars can thus leave single naked He-burning stars with masses larger than about 15 $`\mathrm{M}_{}`$. Less massive helium stars can be produced by mass transfer in a close binary system, if the primary (more massive) component of a binary has $`M_\mathrm{i}>23`$$`\mathrm{M}_{}`$ and the orbital dimensions are such that Roche lobe overflow (RLOF) starts during the main sequence (case A mass transfer) or during the Hertzsprung gap or the first giant branch (case B). The remnant of mass transfer will then be an almost naked helium star with $`M>0.32`$$`\mathrm{M}_{}`$ (the minimum mass for helium burning) orbiting a more massive main sequence star (van den Heuvel 1994). This allows quite a large range in orbital periods as well as initial masses, and helium stars in binaries are therefore expected to be quite common.111The fact that not many such systems are known can be understood by considering that the main sequence star completely dominates the spectrum at optical wavelengths (Pols et al. 1991). Here we will consider helium stars that are massive enough to undergo core collapse and end their lives as neutron stars or black holes, i.e. $`M_{\mathrm{He}}\text{ }>22.5`$$`\mathrm{M}_{}`$. This requires initial masses of at least $`810`$$`\mathrm{M}_{}`$.

### 1.1 Mass Loss and the Formation of Black Holes

Helium stars with $`M\text{ }>510`$$`\mathrm{M}_{}`$ have strong mass loss themselves; they are identified with hydrogen-free Wolf-Rayet (WR) stars of type WN or WC. Mass loss strongly influences the evolution of these stars, as well as their final masses and fate (Langer 1989). Unfortunately, the observationally derived mass-loss rates for WR stars are very uncertain. Evolution models for He stars computed with the mass-loss parametrisations suggested by Langer (1989) or Hamann, Koesterke, & Wessolowski (1995) result in strong mass convergence: even for the largest initial masses the final mass before core collapse is no more than 3 to 4 $`\mathrm{M}_{}`$ (e.g. see Woosley, Langer, & Weaver 1995; Wellstein & Langer 1999). However, recent WR wind models that take into account the inhomogeneous structure of the wind have led to a downward revision of the mass-loss rate by a factor of 3 to 5 (Hamann & Koesterke 1998; Nugis & Lamers 2000). In fact, the Nugis & Lamers rate is nearly an order of magnitude smaller than the Hamann et al. (1995) rate for WN stars.

Whether the helium star leaves a neutron star or a black hole remnant will be determined to a large extent by its final mass before core collapse, or rather, by its final core mass. However, the outcome depends on the details of the explosion mechanism and the limiting mass is very uncertain (e.g. see Fryer et al. 2002). Important observational constraints come from X-ray binaries with low-mass companions (LMXB) in which the dynamically inferred mass of the compact star exceeds 3 $`\mathrm{M}_{}`$, the maximum possible neutron star mass. A strict lower limit to the black hole (BH) mass is set by the mass function, which in several LMXB is at least 6 $`\mathrm{M}_{}`$ (e.g. Casares, Charles, & Naylor 1992; McClintock et al. 2001). In a few cases, the inferred BH mass is very likely to be $`\text{ }>10`$$`\mathrm{M}_{}`$ (McClintock 1998; Orosz et al. 2001). In the evolutionary scenarios for the formation of BH–LMXB the immediate progenitor of the black hole is a naked He star. Because the black holes in these systems can hardly have accreted any mass (King & Kolb 1999), the pre-explosion mass must have exceeded the BH mass, probably by a substantial amount if the collapse was accompanied by a supernova explosion. Clearly, if all naked He stars reach final masses of only $`34`$$`\mathrm{M}_{}`$ these facts cannot be accounted for. One possible solution is that the progenitors of the observed systems result from case C mass transfer (i.e. RLOF started after central He exhaustion) rather than case A/B (Brown, Lee, & Bethe 1999). In that case the naked He star has already gone through core He burning and is close to core collapse when it forms, so there is insufficient time to lose a significant amount of mass in a stellar wind. Whether this scenario can explain all (or any) of the observed BH-LMXBs depends critically on the range of initial masses and orbital periods for which case C mass transfer is possible, which in turn depends quite sensitively on uncertain details of stellar evolution models. However, in the light of the revised WR mass-loss rates, it is worthwhile to reconsider whether case A/B mass transfer cannot after all produce massive black holes. Nelemans & van den Heuvel (2001), using a simple analytic estimate, suggest that this is indeed possible. In this paper we present full-scale evolutionary calculations incorporating the Nugis & Lamers (2000) mass-loss rate in order to investigate this question.

### 1.2 Supernovae of Types Ib and Ic

The amount of mass loss from a helium star, and hence its final mass, also has important consequences for the type of supernova (SN) explosion it produces. Type Ib supernovae show helium lines in their spectra, and their connection to the core collapse of helium stars is quite straightforward. Type Ic supernovae show little or no evidence for helium, and their progenitors are less obvious. However, the similarity of their late-time spectra indicates that the progenitors of SN Ib and SN Ic are related, and in fact small amounts of helium may be present in SNe Ic (Filippenko, Barth, & Matheson 1995). Hence helium cores that have been stripped of all or most of their helium layers are the best candidates for SNe Ic. As with the formation of helium stars themselves, either strong stellar-wind mass loss or mass transfer in a binary may be responsible for this additional stripping. If the initial mass is large enough, the strong WR wind can expose the C-rich core and remove most of the helium. In binary systems that have gone through case A/B mass transfer, a second phase of mass transfer is possible from the helium star to its companion (case BB mass transfer), if the orbit is close enough and the helium star is not too massive ($`M\text{ }<6`$$`\mathrm{M}_{}`$).

In this paper we attempt to identify likely progenitors of Type Ic supernovae and of black holes in X-ray binaries in the light of revised mass-loss rates for Wolf-Rayet stars. To this end we present evolutionary calculations of massive helium stars, losing mass both through a stellar wind and by mass transfer to a nearby compact binary companion. We describe the calculations and their initial conditions in Section 2, and the relevant results in Section 3. The implications are discussed in Section 4.

## 2 Evolutionary Calculations and Initial Conditions

The stellar evolution calculations on which this paper is based are described in detail by O.R. Pols (in preparation) and Dewi et al. (2002). The initial configuration at the start of each calculation is a homogeneous star of almost pure helium with a solar fraction of heavier elements ($`Y=0.98`$, $`Z=0.02`$) which we evolve through central helium burning and through carbon burning. We consider initial helium star masses $`M_{\mathrm{He},\mathrm{i}}`$ between 2 and 32 $`\mathrm{M}_{}`$. The underlying assumption is that such a He star has formed by case A or case B mass transfer in a binary, and we consider two extreme cases.

In the first case (I) we assume that the companion is a main sequence star and that mass transfer was conservative, which leads to a widening of the orbit (van den Heuvel 1994). For this case we computed the evolution of the He star in isolation, i.e. not considering the possibility of a second (case BB) mass transfer phase as the He star expands. This will be correct except perhaps for He stars less massive than 2.5 $`\mathrm{M}_{}`$ which can attain radii of 100 $`\mathrm{R}_{}`$ or more. In the calculations we apply the stellar-wind mass-loss rate of Nugis & Lamers (2000), which depends on the luminosity $`L`$ and on the surface abundances, thus:

$$\dot{M}=1.0\times 10^{11}(L/\mathrm{L}_{})^{1.29}Y^{1.7}Z^{0.5}\mathrm{M}_{}/\mathrm{yr}.$$
(1)
The composition dependence causes mass loss to accelerate when C and O are exposed at the stellar surface, which has interesting consequences. This set of single-star models applies both to sufficiently wide case A/B binaries, and to actually single He stars with $`M_{\mathrm{He},\mathrm{i}}\text{ }>15`$$`\mathrm{M}_{}`$, which can be formed by stellar-wind mass loss. These calculations are described by O.R. Pols (in preparation).

In the second case (II) we assume that the companion is a neutron star (NS) in a close orbit, and we compute the non-conservative mass transfer that ensues when the He star fills its Roche lobe. We assume that the NS accretes up to its Eddington limit and that the excess mass is lost from the system with the specific orbital angular momentum of the NS. Such a system can form as a result of the spiral-in of the neutron star in the envelope of a massive star evolving off the main sequence in an initially wide orbit, i.e. as a remnant of a Be/X-ray binary (van den Heuvel 1994). We consider a range of orbital periods for the He-star + NS binary, consistent with the (very uncertain) periods expected after such a spiral-in (Dewi et al. 2002). This case can be considered as an advanced evolutionary stage of case I, i.e. after the first-formed He star has undergone core collapse and become a neutron star. In these calculations we apply a somewhat different mass-loss rate, namely one fourth of the $`L`$–dependent rate proposed by Hamann et al. (1995). This rate is larger (by roughly a factor of two) than the Nugis & Lamers rate, but there is no real discrepancy with case I because the final mass in case II is determined by RLOF and not by the stellar wind for the masses considered (up to 6.6 $`\mathrm{M}_{}`$; more massive He stars expand very little and usually avoid RLOF). We take full account of the orbital evolution resulting from non-conservative RLOF, the stellar wind, and gravitational wave radiation (for details see Dewi et al. 2002).

We note that other cases are possible, e.g. non-conservative case B mass transfer leading to a He-star + main sequence binary in a fairly close orbit. In such a system we also expect the He star to undergo case BB mass transfer, but this time to a main sequence star and — in all likelihood — conservatively. The result would probably be intermediate between cases I and II, but this should be borne out by actual calculations.

## 3 Results: Final Masses and Helium Amounts

In this section we concentrate on the final configurations resulting from the calculations: the stellar mass $`M_{\mathrm{He},\mathrm{f}}`$ and the amount of helium left in the envelope $`\mathrm{\Delta }M_{\mathrm{He},\mathrm{f}}`$ just before core collapse,222The values given refer to the end of the calculations, which in most cases extended well into carbon shell burning. During the very short remaining time until core collapse, the masses cannot change significantly and can indeed be considered as final. and we discuss their implications for supernovae and black hole formation. For other aspects of the models we refer to the papers by O.R. Pols (in preparation) and Dewi et al. (2002).

### 3.1 Case I

The final masses of the single He-star models (i.e. remnants of conservative case A/B mass transfer) are plotted in Figure 1. At the top of the figure, the initial zero-age main sequence (ZAMS) mass has been indicated, under the assumption that the initial mass of the He star equals the core mass at helium ignition. This is approximately correct for case B, but underestimates the ZAMS mass if the He star formed by case A. The relation $`M_{\mathrm{He},\mathrm{i}}=0.098M_{\mathrm{ZAMS}}^{}{}_{}{}^{1.35}`$ has been used to estimate the ZAMS mass (Hurley, Pols, & Tout 2000, equation 44). Final masses up to 13 $`\mathrm{M}_{}`$ are reached for the most massive stars considered. Although the final–initial mass relation levels off, no mass convergence is apparent. We see that for case B binaries starting with $`M_{\mathrm{ZAMS}}>27`$$`\mathrm{M}_{}`$, final masses in excess of 6 $`\mathrm{M}_{}`$ are reached. This is consistent with the largest observed black hole masses if no additional mass is lost during the collapse of the black hole. On the other hand if we assume, for example, that 25 per cent of the mass is ejected when the black hole forms, the ZAMS mass needs to be greater than 40 $`\mathrm{M}_{}`$. Under this assumption it is still possible to obtain even a 10 $`\mathrm{M}_{}`$ black hole after case B mass transfer, starting from a 32 $`\mathrm{M}_{}`$ He star (or $`M_{\mathrm{ZAMS}}70`$$`\mathrm{M}_{}`$).

A slight break in the initial–final mass relation is apparent between $`M_{\mathrm{He}}=10`$ and 12 $`\mathrm{M}_{}`$. This is caused by the composition dependence in the Nugis & Lamers mass-loss rate which increases with $`Z`$, the fraction of heavy elements. For initial mass up to 10 $`\mathrm{M}_{}`$, the products of He burning are never exposed to the surface, while for larger masses the surface becomes carbon enriched (these stars make the transition from WN to WC stars). As a result the mass loss accelerates when this transition is made.

In Figure 2 we show the final mass of 4He left in the envelope prior to core collapse, as a function of initial He star mass. For case I two regimes can be distinguished. For $`M_{\mathrm{He},\mathrm{i}}\text{ }<10`$$`\mathrm{M}_{}`$ more than 1 $`\mathrm{M}_{}`$ of He remains in the envelope, while for $`M_{\mathrm{He},\mathrm{i}}\text{ }>12`$$`\mathrm{M}_{}`$ only $`0.4`$$`\mathrm{M}_{}`$ of He is left. The transition is quite sharp, for the same reason as indicated above: mass loss accelerates when the star becomes a WC star and removes a large fraction of the remaining helium. The transition is smoother if mass loss does not depend on composition but only on luminosity, as in the calculations of Wellstein & Langer (1999, their Figure 6). Following the arguments in that paper, it is tempting to identify progenitors of Type Ib supernovae with $`M_{\mathrm{He},\mathrm{i}}\text{ }<10`$$`\mathrm{M}_{}`$ and those of Type Ic SNe with $`M_{\mathrm{He},\mathrm{i}}\text{ }>12`$$`\mathrm{M}_{}`$. Although the latter still have a substantial amount of helium, this may not show up in the SN spectrum if it is not mixed with the radioactive 56Ni (Woosley & Eastman 1997). Such Type Ic SN progenitors would be massive ($`>7`$$`\mathrm{M}_{}`$), and probably leave black hole remnants, although in view of the large uncertainties in the core-collapse mechanism, neutron star remnants cannot be excluded. Although rather massive black holes may form by direct collapse (Fryer 1999), in at least one LMXB there is strong evidence that the formation of the black hole was indeed accompanied by a supernova explosion (Israelian et al. 1999). These explosions can possibly be identified with hypernovae such as SN 1998bw, which was of Type Ic and associated with a $`\gamma `$-ray burst. A massive, almost bare CO core in combination with a large explosion energy can explain the bright and slowly declining lightcurve and broad spectral features of such hypernovae (Iwamoto et al. 2000).

### 3.2 Case II

Helium stars in close orbits around a NS companion initially also lose mass in a stellar wind. When non-conservative RLOF starts, much higher mass-loss rates are reached, up to $`10^4`$$`\mathrm{M}_{}`$/yr, because mass transfer occurs on the thermal timescale of the He star (see Dewi et al. 2002 for more details). At the end of the calculations most of the envelope has been transferred and lost from the binary system, unless RLOF starts when the He star is already close to carbon burning. As a consequence, the final masses, shown in Figure 1, depend somewhat on the initial orbital period of the He-star + NS system, but are usually close to the CO core mass except for the widest orbits and most massive He stars. In all cases they are substantially smaller (1.5 – 3 $`\mathrm{M}_{}`$) than in Case I, where only the stellar wind operates. For $`M_{\mathrm{He},\mathrm{i}}>6.6`$$`\mathrm{M}_{}`$ RLOF becomes dynamically unstable, but such systems are rare because He stars of such masses hardly expand. The small final mass implies core collapse will result in neutron star formation. Unless the explosion disrupts the binary, the remnants of these systems are therefore double NS binaries which, if close enough to merge in a Hubble time, are candidate progenitors of $`\gamma `$-ray bursts.

The final mass of He in the envelope, as shown in Figure 2, is also much smaller than in case I but, like the final stellar mass, depends on the initial He star mass and on the orbital period. In most cases it is $`\text{ }<0.2`$$`\mathrm{M}_{}`$, and for the shortest orbital periods, $`P0.1`$ day, it can be as small as 0.04 $`\mathrm{M}_{}`$. This is much less than is achieved for massive He stars by stellar wind only (case I, see above). It is also less than in the conservative case BB mass-transfer models by Wellstein & Langer (1999). We conclude that non-conservative case BB mass transfer to a compact companion is the most efficient way to produce almost bare CO cores prior to explosion. These stars will almost certainly produce a Type Ic supernova. The small progenitor mass would result in a relatively faint, fast declining lightcurve. Such a model was first suggested by Nomoto et al. (1994) as the progenitor of the Type Ic SN 1994I, and shown to match the observed lightcurve (Iwamoto et al. 1994).

## 4 Summary and Conclusions

If the reduced and composition-dependent mass-loss rate for WR stars of Nugis & Lamers (2000) is adopted, the following picture emerges. Binary components with $`M\text{ }>35`$$`\mathrm{M}_{}`$ that form WR stars through case A or B mass transfer, as well as single stars massive enough to form WR stars directly, reach final He-star masses in excess of 7 $`\mathrm{M}_{}`$ and have only a small amount of He ($`0.4`$$`\mathrm{M}_{}`$) left in their envelopes. Such stars probably leave black hole remnants, and the final masses are large enough to be consistent with the observed BH masses in X-ray binaries. Unless such black holes form by direct collapse, these stars very likely undergo a Type Ic supernova explosion which can possibly be identified with bright, slowly declining hypernovae such as SN 1998bw.

Less massive stars in binaries undergoing case A/B mass transfer, but wide enough to avoid case BB mass transfer, have smaller final masses and at least 1 $`\mathrm{M}_{}`$ of He left in their envelope. These stars undergo a Type Ib supernova and leave either neutron stars or, possibly, black holes for stars at the upper end of the mass range.

Binary components with initial mass $`\text{ }<20`$$`\mathrm{M}_{}`$ and in close enough orbits undergo case BB mass transfer. Independently of the adopted WR mass-loss rate, this results in even smaller final masses and a smaller amount of He left in the envelope. The most efficient way to end up with an almost He-free star prior to core collapse is by case BB mass transfer to a neutron star in a very close orbit (i.e. RLOF from the initial secondary component of the binary, after the primary has already collapsed). For small initial masses and orbits with $`P0.1`$ day it is possible to have as little as 0.04 $`\mathrm{M}_{}`$ of He left in the envelope. Such stars almost certainly explode as a SN Ic, and have very small pre-SN masses (1.5 to 3 $`\mathrm{M}_{}`$). They can be identified with faint, fast declining Type Ic supernovae like SN 1994I.

In summary, our most important conclusions are:

* Adoption of the Nugis & Lamers mass-loss rate leads to final He-star masses after case A/B mass transfer that are consistent with the observed black hole masses. There is no need to resort to case C mass transfer to explain the observed BH binaries, although such a scenario remains a possibility.
* Massive binary evolution leads to two distinct classes of SN Ic progenitors: one with large pre-explosion mass, $`M>7`$$`\mathrm{M}_{}`$, formed from initial masses $`M>35`$$`\mathrm{M}_{}`$, which probably leave BH remnants and can possibly be identified with hypernovae, and the other with small pre-explosion mass, 1.5 – 3 $`\mathrm{M}_{}`$, formed by case BB mass transfer in binaries with $`M1220`$$`\mathrm{M}_{}`$, which leave NS remnants.

The latter conclusion is similar to that of Wellstein & Langer (1999), but in their models massive SN Ic progenitors have much smaller pre-explosion masses, $`34`$$`\mathrm{M}_{}`$.

We note that other consequences of the mass-loss rate of WR stars also need to be explored, in particular for the properties of WR stars themselves. Their luminosity and abundance distributions and the number ratio of WN to WC stars all depend on the adopted mass-loss rate. This is beyond the scope of this paper, but we emphasise that in order to draw definite conclusions about WR mass loss and BH formation, all these aspects have to be considered in conjunction.

## Acknowledgments

It is a pleasure to thank Norbert Langer, Ed van den Heuvel, Gijs Nelemans, and Gert-Jan Savonije for stimulating discussions and valuable suggestions. Jasinta Dewi acknowledges sponsorship by NWO Spinoza Grant 08-0 to E.P.J. van den Heuvel.

## References

Brown, G.E., Lee, C.H., & Bethe, H.A. 1999, NewA, 4, 313

Casares, J., Charles, P.A., & Naylor, T. 1992, Nature, 355, 614

Dewi, J.D.M., Pols, O.R., Savonije, G.J., & van den Heuvel, E.P.J. 2002, MNRAS, in press (astro-ph/0201239)

Filippenko, A.V., Barth, A.J., & Matheson, T. 1995, ApJ, 450, L11

Fryer, C.L. 1999, ApJ, 522, 413

Fryer, C.L., Heger, A., Langer, N., & Wellstein, S. 2002, ApJ, submitted (astro-ph/0112539)

Hamann, W.R., & Koesterke, L. 1998, A&A, 335, 1003

Hamann, W.R., Koesterke, L., & Wessolowski, U. 1995, A&A, 299, 151

Hurley, J.R., Pols, O.R., & Tout, C.A. 2000, MNRAS, 315, 543

Israelian, G., Rebolo, R., Basri, G., Casares, J., Martin, E.L. 1999, Nature, 401, 142

Iwamoto, K., Nomoto, K., Hoflich, P., Yamaoka, H., Kumagai, S., & Shigeyama, T. 1994, ApJ, 437, L115

Iwamoto, K., et al. 2000, ApJ, 534, 660

King, A.R., & Kolb, U. 1999, MNRAS, 305, 654

Langer, N. 1989, A&A, 220, 135

McClintock, J.E. 1998, in Accretion Processes in Astrophysical Systems, AIP Conf. Proc. 431, ed. S.S. Holt & T.R. Kallman (New York: AIP), 290

McClintock, J.E., et al. 2001, ApJ, 555, 477

Nelemans, G., & van den Heuvel, E.P.J. 2001, A&A, 376, 950

Nomoto, K., Yamaoka, H., Pols, O.R., van den Heuvel, E.P.J., Iwamoto, K., Kumagai, S., & Shigeyama, T. 1994, Nature, 371, 227

Nugis, T., & Lamers, H.J.G.L.M. 2000, A&A, 360, 227

Orosz, J.A., et al. 2001, ApJ, 555, 489

Pols, O.R., Coté, J., Waters, L.B.F.M., & Heise, J. 1991, A&A, 241, 419

van den Heuvel, E.P.J. 1994, in Interacting Binaries, Saas-Fee Advanced Course 22, ed. H. Nussbaumer & A. Orr (Berlin: Springer-Verlag), 263

Wellstein, S., & Langer, N. 1999, A&A, 350, 148

Woosley, S.E., & Eastman, R.G. 1997, in Thermonuclear Supernovae, NATO ASI Series C, Vol. 486, ed. P. Ruiz-Lapuente, R. Canal, & J. Isern (Dordrecht: Kluwer), 821

Woosley, S.E., Langer, N., & Weaver, T.A. 1995, ApJ, 448, 315


text:
# Discrete Time Leads to Quantum-Like Interference of Deterministic Particles

## 1 Introduction

It is well known that historically the results of experiments with elementary particles were interpreted as the evidence of the impossibility to provide deterministic, classical-like, description of motion of these objects. The main attitude in the development of quantum theory was deeper and deeper understanding of the fact that quantum randomness has fundamental, irreducible, character – in the opposite to classical randomness. Randomness in classical statistical mechanics can be reduced to uncertainty of initial conditions. The evolution of probabilistic density described by Liouville’s equation can be reduced to deterministic evolution of an ensemble of particles described by Hamiltonian equations on the phase space. As we have just mentioned, such a picture was considered as totally inadequate to quantum situation. The formulation by N. Bohr the complementarity principle was the culmination of anti-deterministic development of views to experiments with elementary particles. The collection of these views (originated by Bohr, Heisenberg and many others, see e.g. , ) is now days known as the orthodox Copenhagen interpretation of quantum mechanics. By this interpretation it is (even in principle) impossible to provide deterministic description of motion of elementary particles. Thus we could not reproduce statistical results of quantum experiments by using classical statistical mechanical approach: deterministic equations for trajectories of e.g. electron that reproduce statistical behaviour described by a wave function. This idea appeared already in letter’s exchange between Heisenberg and Bohr directly after the publication of the famous Heisenberg paper . By discovering dynamical equations that describe physical observables, e.g. position and momentum observables, W. Heisenberg claimed that deeper, ontic, description of quantum systems is even in principle impossible. This viewpoint was strongly supported by the discovery of Heisenberg’s uncertainty principle. This principle was used by N. Bohr as the basis for starting great changes in philosophy of physics, resulted by the complementarity principle.

We would like to notice, see  for the details, that, in fact, Heisenberg-Bohr conclusions were not logically justified. By creating a mathematical formalism for physical observables, Heisenberg matrice-mechanics, we do not prove the impossibility of ontic deterministic description. Neither by referring to uncertainty principle. If this principle would be interpreted statistically, see e.g. , , then this would be simply a relation for dispersions of two random variables. It seems very doubtful that such a relation should imply such strong restrictions on mathematical description of reality as e.g. impossibility to describe trajectories of e.g. electrons in configuration or phase space.

Now days we can definitely say that Heisenberg-Bohr conclusions were not totally justified, since we have e.g. Bohmian formalism ,  or stochastic electrodynamics, see e.g. ,  on results related to this paper. For example, the Bohmian formalism provides the deterministic description of trajectories of e.g. electron. Statistical results given by quantum formalism can be reproduced on the basis of this deterministic picture. However, by some reasons Bohmian theory is commonly considered as unacceptable. It seems that the appearance of this ontic model for quantum mechanics does not change essentially the orthodox Copenhagen orientation of quantum community. By ourself we do not have definite point of view to the validity of Bohmian model. In any case it could be used as an argument against Heisenberg-Bohr conclusions.

We think that there might be created other, non-Bohmian, ontic models reproducing probabilistic results given by quantum formalism. In particular, it might be that some of these models could be local. Of course, the reader may argue that there are Bell’s arguments. These arguments imply that local deterministic description is impossible, see e.g. , . However, recent investigations, see e.g.  (papers of Accardy and Regoli, Ballentine, De Muynck, Gudder, Volovich),  –, demonstrated that Bell’s conclusions were not totally justified. It seems that experimental violations of Bell’s inequality need not be interpreted as arguments against local realism.

‘Local realism’ is the standard terminology used in Bell-discussions. Of course, the use of such a terminology was a consequence of EPR-discussions. We would not like to tell about realism. Of course, we do not deny the existence of independent physical reality. But we understood well that all our models are simply approximations of this reality. We would never create the model that would be totally adequate to physical reality. We prefer to speak about deterministic models that provide some mathematical picture of reality. The crucial point of this consideration is that, in principle, a mathematical model of space-time need not be identified with the ‘continuous’ model given by real numbers. When we say deterministic description, we do not mean that this is some model based on differential equations in the continuous real space. For instance, it can be some discrete deterministic model or a $`p`$-adic model. The latter models were intensively studied in the connection to string theory, see, for instance, -. At the moment we do not say the $`p`$-adic program of reconstruction of theoretical physics was totally successful. In any case it is an interesting (deeply developed) alternative approach to mathematical modelling of physical phenomena. However, it is not the subject of our present consideration.

Regarding to space-time models, we, finally, remark that locality in the most of ‘local realism’ discussions also is considered in the old fashioned, real continuous, form. It would be natural to extend this 18th-20th century approach and consider not only real locality, but e.g. locality on a discrete space or $`p`$-adic locality. As having large experience of the work with $`p`$-adic numbers, we can tell that there is crucial difference between e.g. real and $`p`$-adic locality.

As the reader understood, we are looking for local deterministic models in above mentioned extended meaning that could reproduce probabilistic behaviour described by quantum theory.

We recall that one of the most distinguishing features of ‘quantum’ probabilistic behaviour is the appearance of interference structures created by ensembles of elementary particles. In fact, the impossibility to provide deterministic description of such a phenomenon was the main reason to create the orthodox Copenhagen interpretation. Thus deterministic models reproducing interference pictures are of the large interest as at least simulating quantum behaviour.

We remark that according to the orthodox Copenhagen interpretation it is (even in principle) impossible to create such models. In particular, we can not describe ‘self-interference’ of e.g. electron in the two slit experiment without to use wavelike arguments. In particular, there is a rather common opinion that there is crucial difference between classical and quantum probabilistic rules for addition of probabilities of alternatives:

$$P=P_1+P_2$$
(1)
$$P=P_1+P_2+2\sqrt{P_1P_2}\mathrm{cos}\theta .$$
(2)
However, recently it was demonstrated in series - of papers that the difference between ‘classical’ and ‘quantum’ probabilities can be explained in classical probabilistic terms by taking into account contextual dependence of probabilities involved into quantum interference. Moreover, in – there was presented the idea that quantum formalism was merely the pure mathematical discovery of such calculus of probabilities depending on contexts, complexes of experimental physical conditions. Of course, everybody knows that contextualism was originally incorporated into quantum theory by N. Bohr. Thus, in fact, contextual probabilistic approach is nothing else than mathematical (classical probabilistic) formalization of Bohr’s contextualism. We also remark that ‘classical probabilistic’ should not identified with some concrete probabilistic model; in particular, conventional measure-theoretical model, Kolmogorov, 1933 . For us, ‘classical probabilistic’ is a frequency description. There can be various mathematical models for such a description. Of course, quantum probabilistic calculus, Hilbert space calculus, also gives a frequency description. However, here we start with Hilbert space that appears without direct relation to frequencies; such a Hilbert space description (with corresponding interpretation of superposition principle) looks merely as the mathematical model for wave phenomena. In – we start with frequencies and reproduce complex wave amplitudes and the Hilbert space structure as a consequence of contextual dependence of probabilities.

Understanding of the contextual structure of quantum probabilities implies that, on one hand, we need not follow to orthodox Copenhagen 111In particular, there can be created local deterministic models for quantum statistics.; on the other hand, quantum-like probabilistic behaviour need not be related only to experiments with elementary particles. In particular, we can try to obtain interference-like effects in experiments with macroscopic systems by taking into account context dependence. In our paper  we presented numerical experiment for macroscopic charged balls that can be considered as the direct analogue of the two slit experiment. We found the interference effect and introduced corresponding complex waves (of course, of probabilities).

In this paper we study the deterministic model for scattering of charged particles on the charged screen with a single slit. The resulting interference picture has a nontrivial minimum-maximum distribution. This interference picture has no relation to ‘self-interference’ of particles, no wave-structure is involved into considerations. The basic source of interference is the discrete time scale used in our mathematical model: instead of Newton’s differential equations, continuous time evolution, we consider difference equations, discrete time evolution. Interference effect disappears as the time discreetness parameter goes to zero that could be interpreted as an analog of quantum decoherence.

The common viewpoint might be that we study just a discrete approximation to continuous Newton model. It is supposed that the latter model gives the right picture of ‘classical physical reality.’ However, we think that continuous Newton’s model is just an approximation of physical reality. The right picture is given by discrete difference equations. Hence the contradiction between statistical description provided by quantum formalism and Newtonian mechanics could not be considered as a contradiction between quantum and classical (deterministic) physics. Such a contradiction, that typically discussed in quantum literature, should be interpreted as simply a consequence of the use of an approximation, namely continuous Newtonian mechanics, instead of the use of the adequate model, namely discrete model with some level of discretization depending on an experimental context.

## 2 Discrete Time in Newton’s Equations

Classical particles are believed to obey the well known Newton’s equation

$$𝑭=m\ddot{𝒓}$$
(3)
Here we modify this equation to produce an interference picture similar to quantum interference. We introduce a parameter of time discreetness $`\tau `$ described below.

Let us rewrite the second order differential equation (3) as a system of first order differential equations, we have

$$\begin{array}{cc}\hfill 𝑭& =m\dot{v}\hfill \\ \hfill v& =\dot{r}\hfill \end{array}$$
(4)
In the system (4) the derivatives assume the continuousness of time. Let us now introduce a discreetness parameter $`\tau `$.

$$\begin{array}{cc}\hfill 𝑭& =m\frac{v(t+\tau )v(t)}{\tau }\hfill \\ \hfill v(t+\tau )& =\frac{r(t+\tau )r(t)}{\tau }\hfill \end{array}$$
(5)
In the limit of $`\tau 0`$ (5) is equivalent to (3) and (4).

In the model described below we consider particles which move obeying the system (5) where the force is produced by a charged screen.

Please note that in our model the coordinate space is left continuous, although it would be interesting to consider it on the discrete coordinate space, i.e. on the lattice.

## 3 The Model

We consider a scattering on the single slit (Fig.1). Uniformly charged round particles are emitted at point $`e`$ (emitter) with fixed velocity with angles evenly distributed in the range $`(\pi /2,\pi /2)`$. Each particle interacts with the uniformly charged flat screen $`S_1`$. The charge distribution on the particle and the screen stays unchanged even if the particle comes close to the screen. Physically this is a good approximation when the particle and the screen are both made of dielectric. There is a rectangular slit in the screen (on the Fig.1 the slit is perpendicular to the plane of the picture). Particles pass through the slit in screen $`S_1`$ and gather on screen $`S_2`$. We are interested in the particle distribution on the second screen.

Now let us write the laws of motion for the particles. The force affecting the particle is given by the Coulomb’s law

$$𝑭_i=\underset{𝒟_i}{}\frac{q\sigma }{|𝒓^{}|^2}\frac{𝒓^{}}{|𝒓^{}|}𝑑s$$
(6)
where $`𝒓^{}`$ is a vector from an element on the screen to the particle, $`q`$ is charge of the particle, $`\sigma `$ is charge density on the screen, i.e. charge of a unit square. We integrate over the surface of the screen, the integration region $`𝒟_i`$ is plane of the screen except the split.

Projecting equation (6) to $`xy`$-plane, where $`x`$ and $`y`$ denotes horizontal and vertical coordinates of the particle respectively we get

$$\begin{array}{cc}\hfill F_x& =q\sigma \underset{\mathrm{\Gamma }}{}𝑑y^{}\underset{}{}𝑑z^{}\frac{x}{(x^2+(yy^{})^2+z_{}^{}{}_{}{}^{2})^{3/2}}\hfill \\ \hfill F_y& =q\sigma \underset{\mathrm{\Gamma }}{}𝑑y^{}\underset{}{}𝑑z^{}\frac{yy^{}}{(x^2+(yy^{})^2+z_{}^{}{}_{}{}^{2})^{3/2}}\hfill \end{array}$$
(7)
where $`2R`$ is the height of the slit, $`F_x`$ and $`F_y`$ denote the projections of the force $`𝑭`$ to $`x`$ and $`y`$ axes, and the integration region

$$\mathrm{\Gamma }=(\mathrm{},R)(R,+\mathrm{})$$
(8)
Integrating the rhs of (7) we get

$$\begin{array}{cc}\hfill F_x& =2q\sigma \left(\pi +\mathrm{arctan}\frac{yR}{x}\mathrm{arctan}\frac{y+R}{x}\right)\hfill \\ \hfill F_y& =q\sigma \mathrm{ln}\frac{x^2+(Ry)^2}{x^2+(R+y)^2}\hfill \end{array}$$
(9)
We take the following initial values

$$\begin{array}{cc}\hfill x(0)& =D\hfill \\ \hfill y(0)& =0\hfill \end{array}\begin{array}{c}\hfill \dot{x}(0)=v_0\mathrm{cos}\alpha \\ \hfill \dot{y}(0)=v_0\mathrm{sin}\alpha \end{array}$$
(10)
where angle $`\alpha `$ is a random variable uniformly distributed in $`[0,2\pi )`$. The constant parameters $`v_0`$ and $`D`$ are initial velocity and distance between emitter and the screen.

Particles are emitted at point $`e`$ (see Fig.1), move obeying affected by force (9) pass through slit in the screen $`S_1`$ and gather on the screen $`S_2`$. Having points where particles hit the screen $`S_2`$ we compute frequencies with which particles appear on screen $`S_2`$ as a function of coordinates on the screen, we call this function a particle distribution. We are interested in computing the particle distribution over a vertical line on screen $`S_2`$ with $`z=0`$. That is why we consider a motion only in the $`xy`$-plane and initial values (10) do not contain $`z`$-coordinate.

The second screen was separated with cells of equal size, the diameter of a particle. The number of particles which hit into each cell was calculated and interpreted as a particle distribution. The details of numeric computations are given in the appendix.

## 4 Conclusion

In this note we have shown that a quantum-like interference picture could appear as a statistical effect of deterministic particles, i.e. having trajectories and obeying deterministic equations, if one introduces a discrete time. The nature of the resulting interference picture (particle distribution, see Fig. 2-5 in appendix) does not follow from the geometry of force field, but is strongly attached to the discreetness parameter $`\tau `$.

The described behavior stays without contradiction with a contextual approach to quantum probabilities. It would be interesting to investigate the scattering on the two slit screen.

## 5 Acknowledgments

One of the authors, Andrei Khrennikov, would like to thank I would like to thank S. Albeverio, L. Accardi, L. Ballentine, V. Belavkin, E. Beltrametti, W. De Muynck, R. Gill, D. Greenberger, S. Gudder, T. Hida, A. Holevo, P. Lahti, D. Mermin, A. Peres, J. Summhammer, for (sometimes extremely critical) discussions on probabilistic foundations of quantum mechanics.

References

1. N. Bohr, Phys. Rev., 48, 696-702 (1935).

2. W. Heisenberg, Physical principles of quantum theory. Chicago Univ. Press, 1930.

3. W. Heisenberg, Zeits. für Physik, 33, 879 (1925).

4. A. Khrennikov, Växjö interpretation of quantum mechanics. Preprint quant-ph/0202197.

5. H. Margenau, Phil. Sci., 25, 23 (1958).

6. L. E. Ballentine, The statistical interpretation of quantum mechanics, Rev. Mod. Phys., 42, 358–381 (1970).

7. D. Bohm (1951), Quantum theory, Prentice-Hall. Englewood Cliffs, New-Jersey.

8. D. Bohm and B. Hiley, The undivided universe: an ontological interpretation of quantum mechanics. Routledge and Kegan Paul, London, 1993.

9. G. Cavalleri, Nuovo Cimento B, 112, 1193 (1997) .

10. A. Zecca and G. Cavalleri , Nuovo Cimento B, 112, 1, (1997); Cavalleri G. and Tonni E., ”Discriminating between QM and SED with spin”, in C. Carola and A. Rossi, The Foundations of Quantum Mechanics (Hystorical Analysis and Open Questions) (World Sceintific Publ., Singapore), p.111, 2000.

11. J.S. Bell, Rev. Mod. Phys., 38, 447–452 (1966). J. S. Bell, Speakable and unspeakable in quantum mechanics. Cambridge Univ. Press, 1987.

12. J.F. Clauser , M.A. Horne, A. Shimony, R. A. Holt, Phys. Rev. Letters, 49, 1804-1806 (1969); J.F. Clauser , A. Shimony, Rep. Progr.Phys., 41 1881-1901 (1978). A. Aspect, J. Dalibard, G. Roger, Phys. Rev. Lett., 49, 1804-1807 (1982); D. Home, F. Selleri, Nuovo Cim. Rivista, 14, 2–176 (1991). H. P. Stapp, Phys. Rev., D, 3, 1303-1320 (1971); P.H. Eberhard, Il Nuovo Cimento, B, 38, N.1, 75-80(1977); Phys. Rev. Letters, 49, 1474-1477 (1982); A. Peres, Am. J. of Physics, 46, 745-750 (1978). P. H. Eberhard, Il Nuovo Cimento, B, 46, N.2, 392-419 (1978); J. Jarrett, No s, 18, 569 (1984).

13. Proc. Conf. ”Foundations of Probability and Physics”, Växjö-2000, Sweden; editor A. Khrennikov, WSP, Singapore, 2001.

14. I. Pitowsky, Phys. Rev. Lett, 48, N.10, 1299-1302 (1982); Phys. Rev. D, 27, N.10, 2316-2326 (1983); S.P. Gudder, J. Math Phys., 25, 2397- 2401 (1984); A. Fine, Phys. Rev. Letters, 48, 291–295 (1982); P. Rastal, Found. Phys., 13, 555 (1983). W. Muckenheim, Phys. Reports, 133, 338–401 (1986); W. De Baere, Lett. Nuovo Cimento, 39, 234-238 (1984); 25, 2397- 2401 (1984).

15. A. Yu. Khrennikov, Interpretations of Probability. VSP Int. Sc. Publishers, Utrecht, 1999.

16. A. Yu. Khrennikov, A perturbation of CHSH inequality induced by fluctuations of ensemble distributions. J. of Math. Physics, 41, 5934-5944(2000).

17. V. S. Vladimirov, I. V. Volovich, and E. I. Zelenov, $`p`$-adic Analysis and Mathematical Physics, WSP, Singapore, 1994.

18. P. G. O. Freund and E. Witten, Adelic string amplitudes. Phys. Lett. B, 199, 191-195 (1987).

19. P. H. Frampton and Y. Okada, $`p`$-adic string $`N`$-point function. Phys. Rev. Lett. B , 60, 484–486 (1988).

20. A. Yu. Khrennikov, $`p`$-adic valued distributions in mathematical physics. Kluwer Acad. Publ., Dordrecht, 1994.

21. A. Yu. Khrennikov, Non-Archimedean analysis: quantum paradoxes, dynamical systems and biological models. Kluwer Acad. Publ., Dordreht, 1997.

22. A. Yu. Khrennikov, Linear representations of probabilistic transformations induced by context transitions. J. Phys.A: Math. Gen., 34, 9965-9981 (2001).

23. A. Yu. Khrennikov, Origin of quantum probabilities. Proc. Conf. ”Foundations of Probability and Physics”, Växjö-2000, Sweden; editor A. Khrennikov, p. 180-200, WSP, Singapore (2001).

24. A. Yu. Khrennikov, Ensemble fluctuations and the origin of quantum probabilistic rule. Rep. MSI, Växjö Univ., 90, October (2000).

25. A. N. Kolmogoroff, Grundbegriffe der Wahrscheinlichkeitsrechnung. Springer Verlag, Berlin (1933); reprinted: Foundations of the Probability Theory. Chelsea Publ. Comp., New York (1956).

26. A. Yu. Khrennikov, Ya. I. Volovich, Numerical experiment on interference for macroscopic particles. Preprint quant-ph/0111159, 2000.

## Appendix

Below we present a sample trajectories plots for two different discretness parameters $`\tau `$ (Fig.4,5). And the corresponding particle distributions (Fig.2,3).

Parameters of the model where the following $`D=5,d=25,R=5,v_0=12`$ the particle and the screen $`S_1`$ (Fig.1) had the opposite charges, i.e. $`q\sigma =1`$ and the radius of a particle $`r=0.2`$.

To produce a particle distribution pictures with a trustable precision about $`10^7`$-$`10^8`$ trajectories where computed. To produce such a large amout of computations even for modern stations we used a parallel Sun-UltraSPARC 4-processor station located at Växjö University. The program was implemented using GNU C++ (g++).

Since the computation time is proportional to $`1/\tau `$ and the computations where performed simultaneously the total number of computed trajectories for $`\tau =0.05`$ (Fig.4) and $`\tau =0.01`$ (Fig.3) differs approximately five times.


text:
# Subharmonic gap structure in d-wave superconductors

## Abstract

We present a self-consistent theory of current-voltage characteristics of d-wave/d-wave contacts at arbitrary transparency. In particular, we address the open problem of the observation of subharmonic gap structure (SGS) in cuprate junctions. Our analysis shows that: (i) the SGS is possible in d-wave superconductors, (ii) the existence of bound states within the gap results in an even-odd effect in the SGS, (iii) elastic scattering mechanisms, like impurities or surface roughness, may suppress the SGS, and (iv) in the presence of a magnetic field the Doppler shift of the Andreev bound states leads to the splitting of the SGS, which is an unambiguous fingerprint of d-wave superconductivity.

Although the d-wave scenario is emerging as the new paradigm in superconducting cuprates Scalapino1995 ; Harlingen1995 ; Tsuei2000 , there is still a lacking consensus with regards to current transport in SIS junctions Kashiwaya2000 . By now it is clear that tunneling in unconventional superconductors is a phase-sensitive technique. An important consequence of an order parameter showing a sign change in different momentum directions is the formation of Andreev bound states at zero energy confined to the surface Buchholtz1981 ; Hu1994 . One of the consequences of the presence of these bound states is the appearance of a zero-bias conductance peak in the current-voltage characteristics. This peak has been observed in different experiments using SIN-type junctions Kashiwaya1995  and SIS grain boundaries Alff1998 . However, there is still a controversy related to the observation of the subharmonic gap structure (SGS) in cuprate SIS junctions. Devereaux and Fulde suggested that the Andreev scattering in SNS contacts could be used to identify the symmetry of the order parameter Devereaux1993 . Based on their analysis the SGS revealed in different experiments in YBCO contacts were interpreted as an evidence of the existence of a well-defined nonvanishing gap Polturak1993 ; Ponomarev1995 . After more detailed but non-self consistent analysis, different authors have concluded that the SGS are weak in d-wave superconductors, since averaging over the anisotropic gap washes out any prominent features Hurd1997 ; Lofwander2001 . However, there are also more recent experimental reports of SGS in YBCO edge Josephson junctions Engelhardt1999 ; Nesher1999 , and even in the c-axis tunneling of different cuprate contacts Ponomarev2001 . Auerbach and Altman proposed an alternative interpretation of the appearance of pronounced SGS in Ref. \[Nesher1999, \], claiming that this structure is an indication of magnon tunneling that can be explained in the context of the SO(5) theory Auerbach2000 .

In this Communication we present a self-consistent theoretical analysis of the ac Josephson effect in d-wave superconductors. Using quasiclassical methods we determine the local electronic properties of the superconductors in the interface region. This includes effects on the order parameter (OP) profile and on the local density of states (DOS) by pair breaking caused both by quasiparticle scattering off the interface and off homogeneously distributed impurities in the crystals Matsumoto1995 ; Buchholtz1995 ; Fogelstrom1997 ; Barash1997 ; Poenicke1999 ; Fogelstrom-preprints . The I-V characteristics of the d-wave/d-wave SIS-junction are computed using the local surface Green’s functions containing all relevant information of superconducting electrodes and solving the appropriate boundary conditions for a point contact Cuevas2001 . We shall show that, not only does the SGS survive after averaging for all crystal misorientations, but there are qualitatively new features in the I-V characteristics that can only be captured by a truly self-consistent calculation. For instance, there is an even-odd effect in the SGS which originates from processes connecting the zero-energy states with Andreev states formed close to the gap edges. Such an even-odd effect in the differential resistance has been reported in YBCO edge-junctions by Nesher and Koren Nesher1999 . These states are missed in a non-self consistent calculation where the suppression of the OP close to the interface is neglected. We also show that the SGS in d-wave contacts can be tuned with a magnetic field. The Doppler shift of the Andreev bound states Fogelstrom1997 ; Aprili1999  leads to a very peculiar splitting of the SGS, which provides a clear signature of d-wave superconductivity.

In the limit of not too low interface transparency, $`D0.1`$, the main feature of the I-V characteristics of conventional SIS contacts is the appearance of SGS, which consists of a series of maxima in the conductance at voltages $`eV=2\mathrm{\Delta }/n`$, where $`\mathrm{\Delta }`$ is the superconducting gap and $`n`$ is an integer number. This peculiar structure is due to the occurrence of multiple Andreev reflections (MAR) and by now well established both theoretically Averin1995  and experimentally Scheer1997  in the case of s-wave superconductors. In this work, our goal is to extend the analysis of the SGS to the case of superconducting cuprates. For this purpose, we consider a voltage biased contact, consisting of two $`d_{x^2y^2}`$ superconductors separated by a single interface of arbitrary transparency. The order parameter on side $`i`$, $`i=L,R`$, is rotated by $`\alpha _i`$ with respect to the surface normal, and we denote junction type by the relative crystal orientations as $`d_{\alpha _L}`$-$`d_{\alpha _R}`$. There are several experimental realizations of this system, among which the bicrystal grain-boundary junctions are ideal examples Alff1998 . We carry out the calculation of the current following the formulation introduced by two of the authors in Ref. \[Cuevas2001, \], and refer the reader to this work for all the technical details. In the case of a constant bias voltage, $`V`$, one can show that the current oscillates in time with all the harmonics of the Josephson frequency, i.e. $`I(t)=_mI_me^{im\varphi (t)}`$, where $`\varphi (t)=\varphi _0+(2eV/\mathrm{})t`$ is the time-dependent superconducting phase difference. We concentrate ourselves in the analysis of the dc current, denoted from now on as $`I`$. Furthermore, we assume that the interface conserves the momentum of the quasiclassical trajectories, which allows us to write the current as a sum over independent trajectory contributions: $`I=\frac{1}{2}_{\pi /2}^{\pi /2}𝑑\widehat{𝒑}_FI(\widehat{𝒑}_F)\mathrm{cos}(\widehat{𝒑}_F)`$, where $`\widehat{𝒑}_F`$ defines the Fermi surface position. For the angular dependence of the transmission coefficient we use the expression $`D(\widehat{𝒑}_F)=D\mathrm{cos}^2(\widehat{𝒑}_F)/[1D\mathrm{sin}^2(\widehat{𝒑}_F)]`$, resulting from a $`\delta `$-like potential. Here $`D`$ is the transmission for the trajectory perpendicular to the interface.

Let us start by analyzing the case of a symmetric $`d_0`$-$`d_0`$ junction in the clean limit. In this case, the order parameter is constant up to the surface, and there are no bound states for any trajectory. The I-V characteristics of a single trajectory, $`I(\widehat{𝒑}_F)`$, coincide with those of isotropic s-wave superconductors, and exhibit a pronounced SGS at $`eV=2\mathrm{\Delta }(\widehat{𝒑}_F)/n`$ \[Averin1995, \]. Since the different trajectories see different gaps, the relevant question is whether the SGS survives after averaging. The answer can be seen in Fig. 2, where both the current and differential conductance are shown for different transmissions. One can clearly see a SGS at $`eV=2\mathrm{\Delta }/n`$, although it is more rounded than in the s-wave case. Notice also that the SGS are better defined as minima in the conductance at the subharmonic voltages, i.e. as maxima in the differential resistance. At this point a natural question arises: why is the maximum gap the energy scale revealed in the SGS? The idea is that the jump of the current at the opening of a new MAR scales with the gap of the corresponding trajectory. Therefore, the trajectories with the largest gaps dominate the contribution to the SGS.

Let us now consider the case of a $`d_{\pi /4}`$-$`d_{\pi /4}`$ junction. In this case, assuming specular quasiparticle scattering at the interface, an Andreev bound state forms at zero energy for every trajectory Hu1994 . This implies that the surface acts as a pair-breaker Buchholtz1995  and the gap is depressed in the vicinity of the interface, vanishing exactly at the barrier. This order-parameter profile induces the appearance of new bound states at the gap edges for some trajectories Barash1997 , which are also visible in angle-averaged DOS, see Fig. 1a. As we show in Fig. 1a these gap singularities are not present in the non-self consistent calculation. This fact has important consequences in the I-V characteristics. As is shown in Ref. \[Cuevas2001, \], assuming a constant order parameter the most prominent feature of the I-Vs is the appearance of SGS at $`eV=\mathrm{\Delta }/n`$, instead of in $`eV=2\mathrm{\Delta }/n`$ as in conventional superconductors. Its origin can be understood as follows. Inside the gap the current is dominated by multiple Andreev reflections (MAR). In this case there are two types of MAR processes: (a) those which connect the bound states with the gap edges and (b) the usual ones connecting the gap edges. The first ones give rise to the SGS at $`eV=\mathrm{\Delta }/n`$, while the second could give rise to the series $`eV=2\mathrm{\Delta }/n`$. However, for the non-self consistent calculation the DOS at the gap edges is rather small (see Fig. 1a), and thus, at the opening of this second type of processes their probability is small. In the self-consistent case the main consequence of the presence of the gap singularities is the recovery of the odd terms in the SGS series $`eV=2\mathrm{\Delta }/n`$. The odd terms appear now due to the enhancement of the probability of the MARs connecting the gap edges. However, the bound states at the gap edges do not appear for every trajectory, which weakens the SGS due to these processes. In this sense, in Fig. 3 one can clearly see a difference between the even and odd maxima of the conductance. This even-odd effect was reported in Ref. \[Nesher1999, \]. This result show clearly the relevance of the self-consistency, even in the ideal clean case.

In d-wave superconductors the order parameter is very sensitive to scattering from nonmagnetic impurities and surface roughness. In particular, it is known that these elastic scattering mechanisms provide an intrinsic broadening for the zero-energy bound states (ZEBS) Poenicke1999 ; Fogelstrom-preprints . For the case of Born scatterers this broadening is $`\sqrt{\mathrm{\Gamma }\mathrm{\Delta }}`$, where $`\mathrm{\Gamma }=1/2\tau `$ is the effective pair-breaking parameter locally at the surface. This is illustrated in Fig. 1b for the case of bulk impurities. The natural question now is: what is the influence of the broadening of the ZEBS in the SGS? In Fig. 4 we show the differential conductance of a $`d_{\pi /4}`$-$`d_{\pi /4}`$ junction for different values of the bulk-impurity scattering rate. Notice that as the elastic scattering rate increases the SGS disappears, which can be understood as follows: the increase of density of states in the gap region enhances the probability of single-quasiparticle processes, producing the subsequent reduction of the probability of the Andreev processes, which in turn leads to the suppression of the SGS.

Probably the clearest signature of d-wave superconductivity in the SGS is its evolution with magnetic field. It is known that a magnetic field perpendicular to the $`ab`$-plane, $`𝐇=H\widehat{𝐳}`$, leads to a Doppler shift in the continuum excitations given by $`𝐯_f𝐩_s`$, where the condensate momentum is $`𝐩_s=(e/c)A(x)\widehat{𝐲}`$, with $`A`$ the self-consistently determined vector potential Fogelstrom1997 . As shown in Fig. 1c, this means that the Andreev bound states are shifted to an energy which, in the limit of a large ratio $`\lambda /\xi _0`$, can be estimated to be $`ϵ_b(\widehat{𝒑}_F)=(e/c)v_fH\lambda \mathrm{sin}\widehat{𝒑}_F`$, $`\lambda `$ being the $`ab`$-plane penetration depth. We shall use a natural field scale set by a screening current of order the bulk critical current, $`H_0=c\mathrm{\Delta }/ev_f\lambda `$, which is of the order of a Tesla Fogelstrom1997 . As can be seen in Fig. 1d, the screening currents flow parallel to the interface and in opposite directions in both electrodes, which means that the trajectory resolved DOS of the left and right superconductors are shifted by $`2ϵ_b(\widehat{𝒑}_F)`$ relative to each other. This shift modifies the threshold voltages of MARs starting and ending in different electrodes, leading to the splitting of the peaks with an odd order $`n`$ in the SGS. On the contrary, since the magnetic field produces a rigid shift of the spectrum, the threshold voltages of those MARs starting and ending in the same electrode are not modified. This means that the positions of the peaks with an even order $`n`$ in the SGS remain unchanged. This fact gives rise to a rich SGS in the trajectory resolved current, which consists in conductance peaks at the following voltage positions: (a) $`[2\mathrm{\Delta }(\widehat{𝒑}_F)\pm 2ϵ_b(\widehat{𝒑}_F)]/n`$, with $`n`$ odd, due to MARs connecting the gap edges of the left and right electrodes (the signs $`\pm `$ correspond to electron and hole processes), (b) $`[\mathrm{\Delta }(\widehat{𝒑}_F)\pm 2ϵ_b(\widehat{𝒑}_F)]/n`$, with $`n`$ odd, due to MARs connecting the bound states and the gap edges of different electrodes, (c) $`2\mathrm{\Delta }(\widehat{𝒑}_F)/n`$, with $`n`$ even, due to MARs connecting the gap edges of the same electrode, (d) $`\mathrm{\Delta }(\widehat{𝒑}_F)/n`$, with $`n`$ even, due to MARs connecting the bound states and the gap edges of the same electrode. Finally, there is also structure at $`eV=2ϵ_b(\widehat{𝒑}_F)/n`$, with $`n`$ odd, due to MARs connecting the bound states.

After doing the angle average most of these features are still clearly visible. This is illustrated in Fig. 5 where we show the differential conductance of a $`d_{\pi /4}`$-$`d_{\pi /4}`$ junction with transmission $`D=0.4`$ for different values of the magnetic field. Starting at large voltages, the weak structure seen at $`2\mathrm{\Delta }`$ splits with applied field. Around $`eV=\mathrm{\Delta }`$ there are both type (b) and type (c) processes leading to a maximum at $`eV=\mathrm{\Delta }`$, unaffected by the applied field, as well as a field-shift of the dip just above $`\mathrm{\Delta }`$. The field dependence of the differential conductance is most clearly resolved at larger biases, $`eV\mathrm{\Delta }/2`$, as the marks of the various processes begin to overlap at small bias. What is important to note is that it is still the bulk maximum gap that gives the characteristic energy scale for the SGS. The effect of the Doppler shift on the SGS or the differential conductance is only prominent in junctions with a sizable misorientation. For junctions close to the $`d_0`$-$`d_0`$ case, the main contribution to the SGS comes from trajectories close to perpendicular incidence, i.e. with $`\mathrm{sin}\widehat{𝒑}_F0`$ and thus having a vanishing Doppler shift.

In summary, we have presented a self-consistent analysis of the I-V characteristics of d-wave/d-wave contacts at high transparencies. We have shown that it is possible to observe SGS for all crystal misorientations. We predict that the presence of bound states inside the gap gives rise to two new qualitative effects in the SGS: (a) an even-odd effect, and (b) a very peculiar splitting in an external magnetic field. These features are unique hallmarks of the d-wave scenario, and we hope that our analysis will trigger off a more detailed experimental study of the SGS in cuprate Josephson junctions.

This work has been supported by the EU TMR Network on Dynamics of Nanostructures, the CFN supported by the DFG, and by the Swedish research council.