text:
# 1 Introduction

## 1 Introduction

Top-quark and in general heavy-flavour production physics (see, for a review, ) is currently one of the main fields of investigation in both experimental and theoretical particle physics. At the next Run II of the Tevatron accelerator and, in the future, at the LHC  and at $`e^+e^{}`$ Linear Colliders , the production of a large amount of $`t\overline{t}`$ pairs will allow accurate studies of the top quark properties and an improved measurement of its mass.

The measurement of these properties will rely by-and-large on the accuracy of the theoretical modeling of the exclusive properties of the final states, for example jet distributions. Use of inclusive parton-level calculations, although resulting in the exact inclusion of the next-to-leading-order (NLO) contributions , is not sufficient to fully describe the effects associated with the large logarithms corresponding to soft or collinear parton emission. Calculations based on the fragmentation function formalism , where collinear logarithms can be resummed up to the NLO, are on the other hand too inclusive to allow a complete study of the final states. Monte Carlo event generators  are the best possible tool to perform the resummation of these enhanced logarithms, to simulate multiple radiation in high energy processes, and to provide a description of the hadronization transition leading to the final observable particles.

Inclusion of these higher-order corrections by the QCD event generators is done however in the soft/collinear approximation. Furthermore, the Monte Carlo evolution suppresses entirely emission of radiation inside some regions of the physical phase space (‘dead zones’) corresponding to hard and large-angle parton radiation. These regions are unfortunately sometimes crucial from the experimental point of view. Emission inside these dead zones can be performed using the exact amplitudes by following the method discussed in . This method has been applied in the past to jet production in $`e^+e^{}`$ annihilation , in deep inelastic scattering  and, more recently, to the description of top decays  and of Drell–Yan processes . Alternative approaches have also been proposed in the literature, see e.g. ref. .

In  a marked impact of matrix-element corrections to top decays was found for $`e^+e^{}`$ interactions slightly above the threshold for $`t\overline{t}`$ production. In this paper we wish to perform a similar analysis for top production and decay at hadron colliders and investigate the effect of the implemented hard and large-angle gluon radiation on jet observables and on the top mass reconstruction. As far as the top mass is concerned, we shall consider final states with leptons and $`J/\psi `$ since the LHC experimentalists claim it is a favourite channel, with a systematic error no larger than 1 GeV , and we shall give more details on the analysis and the preliminary results presented in .

In Section 2 we briefly review the method applied in  to implement matrix-element corrections to the HERWIG description of top decays. In Section 3 we shall study phenomenologically-relevant jet observables at the Tevatron and at the LHC and investigate the impact of the improved treatment of top decays. In Section 4 we shall discuss the method of reconstructing the top mass by using final states with leptons and $`J/\psi `$ and the effect of matrix-element corrections on the top mass measurement. In Section 6 we shall make some concluding remarks and comments on possible further improvements of the study here presented.

## 2 Matrix-element corrections to simulations of top decays

In the HERWIG Monte Carlo event generator, the top quark decay $`tbW`$ is performed in the top rest frame, as discussed in . The top quark cannot emit soft gluons in its decay stage as it is at rest, while the $`b`$ quark is allowed to radiate in the cone $`0<\theta _g<\pi /2`$, $`\theta _g`$ being the soft-gluon emission angle relative to the direction of the $`b`$ quark. The subsequent parton shower is performed following the prescription of the angular ordering . The $`W`$ hemisphere, corresponding to $`\pi /2<\theta _g<\pi `$, is completely empty and the soft phase space is not therefore entirely filled by the HERWIG algorithm. In  the authors showed that neglecting the ‘backward’ gluon radiation correctly predicts the total energy loss, however problems are to be expected when dealing with angular distributions.

In  matrix-element corrections to the HERWIG simulation of top decays have been implemented: the missing phase space is populated according to a distribution obtained from the calculation of the exact first-order matrix element (hard correction), the shower in the already-populated region is corrected by using the exact amplitude any time a hardest-so-far emission is encountered in the evolution (soft correction).

Since the HERWIG dead zone includes part of the soft singularity, matrix-element corrections to top decays are not a straightforward extension of the method applied in . The soft singularity has been avoided by setting a cutoff on the energy of the gluons which are radiated in the backward hemisphere by the $`b`$ quark. As shown in , and discussed later on in this paper, the sensitivity to this cutoff is however very small.

In  $`e^+e^{}`$ interactions at a centre-of-mass energy $`\sqrt{s}=360`$ GeV, slightly above the threshold for $`t\overline{t}`$ production, were considered. This is an ideal phenomenological environment to test the impact of the implemented corrections, with most of the radiation being associated to the top-decay stage. Three-jet events were analysed and a remarkable impact of matrix-element corrections was found when comparing different versions of HERWIG111The latest public version HERWIG 5.9 had some errors in the treatment of top decays, which have been corrected in the intermediate version 6.0. In order to estimate the impact of matrix-element corrections, we therefore compare the versions 6.1 and 6.0. The results of HERWIG 6.1 , the new version which includes also the improved treatment of top decays, were also compared to the ones obtained by the exact $`𝒪\left(\alpha _S\right)`$ matrix-element calculation of the process $`e^+e^{}t\overline{t}\left(bW^+\right)\left(bW^{}\right)\left(g\right)`$ . While the authors of  had found serious discrepancies when comparing the exact $`𝒪\left(\alpha _s\right)`$ results with the ones of HERWIG before matrix-element corrections, good agreement was found in  after matrix-element corrections in the region of large energies and angles, where fixed-order calculations are reliable.

In  it was also shown that, although the fraction of events generated in the dead zone varies from about 2% to 4% when the cutoff changes from 5 GeV to 1 GeV, the dependence of phenomenological distributions on its chosen value is pretty negligible after one applies typical experimental cuts on the jet transverse energy $`E_T>10`$ GeV and on the invariant opening angle between jets $`\mathrm{\Delta }R>0.7`$. The value $`E_{\mathrm{min}}=2`$ GeV was then chosen as the default value and this value will be kept throughout this paper as well. In the following sections, we shall consider hadronic production of $`t\overline{t}`$ pairs and analyse the effect of matrix-element corrections to top decays on jet observables and on the top mass reconstruction at the Tevatron and at the LHC.

## 3 Jet activity in dilepton $`𝒕\overline{𝒕}`$ events

We start our analysis by considering inclusive jet observables. To emphasize the effect of the matrix element corrections to top decays, we confine ourselves to the case of leptonic decays of both $`W`$’s in each event. The most likely hard and well-separated jets are then those from the $`b`$ and the $`\overline{b}`$. Extra radiation from either the initial state (ISR) or the top decay may give rise to extra jets (final state radiation from the produced $`t\overline{t}`$ pair tends to be small because of the large top mass). Interesting observables which may show the effects of the new HERWIG treatment of top decays are related to these extra jets, and in particular to the one with the largest value of $`E_T`$, the ‘third jet’.

We cluster jets according to the inclusive version of the $`k_T`$ algorithm , setting a radius parameter $`R=0.5`$ at the Tevatron and $`R=1`$ at the LHC222Such a choice is due to the relation between the $`R`$ parameter of the $`k_T`$ algorithm and the radius $`R_{\mathrm{cone}}`$ of the cone algorithm $`R_{\mathrm{cone}}0.75\times R`$  and to the fact that the Tevatron experimentalist run a pure cone algorithm with $`R_{\mathrm{cone}}=0.4`$. For the LHC we shall nevertheless stick to the recommended value $`R=1`$.. We set a cutoff on the transverse energy of the resolved jets $`E_T>10`$ GeV, and we study the following inclusive distributions: transverse energy $`E_T`$ and rapidity $`\eta _3`$ of the third jet, the minimum invariant opening angle $`\mathrm{\Delta }R=\sqrt{\mathrm{\Delta }\varphi ^2+\mathrm{\Delta }\eta ^2}`$ among the three hardest jets, the threshold variable $`d_3`$ of the $`k_T`$ algorithm, according to which all events are forced to be three-jet-like. Finally, we consider the number of jets $`n_{\mathrm{jets}}`$ that pass the 10 GeV cut in transverse energy.

We start by comparing the Monte Carlo jet distributions before and after matrix-element corrections. To make the comparison more realistic, we normalize the plots to the expected integrated luminosities (2 fb-1 at $`\sqrt{s}=2`$ TeV for the Tevatron Run II, and 10 fb-1 for 1 year of LHC low-luminosity running 333Using total cross-sections of 7 pb and 830 pb for Tevatron and LHC , and assuming a conservative 1% overall efficiency and BR for the dilepton final states, this correponds to approximately 150  and 8$`\times 10^4`$  events, respectively.), and smear the contents of a given bin with $`N`$ events according to a Gaussian distribution with average $`N`$ and standard deviation $`\sqrt{N}`$. Furthermore, to partially account for detector effects we smear the value of the reconstructed observables with a 10% resolution. The final jet distributions are plotted in Figs. 1-5. 444With no smearing, the $`\mathrm{\Delta }R`$ distributions at the Tevatron and at the LHC would have shown a sharp cutoff for $`\mathrm{\Delta }R=R`$, as predicted by the $`k_T`$ algorithm we have been using, but nevertheless the 10% smearing allows a small fraction of the jet events to have even $`\mathrm{\Delta }R<R`$.

All observables show very small changes due to the matrix-element corrections. To quantify the residual differences, we performed a Kolmogorov–Smirnov test on the distributions, and evaluated the number of events necessary to establish the difference between the 6.0 and 6.1 distributions at 95% CL. The results of this test are shown in table 1. We see that at both machines a number of events of the order of $`10^310^4`$ is sufficient to see a difference in the shape of the two histograms. We therefore conclude that the effects of matrix-element corrections to top decays are in principle detectable at the LHC, where about $`8\times 10^4`$ events per year are expected in the dilepton channel, while the foreseen statistics are too low at the Tevatron.

We conclude this Section by pointing out that the small differences detected in the case of inclusive third-jet observables are largely a consequence of the small fraction of events for which the matrix-element corrections are applied. We evaluated in fact that for the Tevatron energy only approximately 7% of events that contain a third jet required the evaluation of the matrix-element corrections to the top decays (the number is 6% for the LHC). We expect that larger effects and differences will appear once the matrix-element corrections to the $`t\overline{t}`$ production process is included. This work is under way, along similar lines to vector boson production .

## 4 Top mass reconstruction

One observable for which the higher-order corrections described in this paper may induce relevant changes is the invariant mass distribution of $`J/\psi +\mathrm{}`$ pairs from $`tbW`$ decays, with the $`J/\psi `$ coming from the decay of a $`b`$ hadron, and the isolated lepton from the decay of the $`W`$. This distribution was recently suggested  as a potential way of measuring the top mass during the high-luminosity phase of the LHC with experimental accuracies better than 1 GeV. The dominant source of experimental uncertainty is the statistical one, due to the strong suppression from the small branching ratios 555Ref.  estimates a sample of approximately $`10^3`$ events in one year of high luminosity running at the LHC, $`=10^5`$ $`\mathrm{pb}^1`$, for a production cross section $`\sigma _{\mathrm{LHC}}(t\overline{t})=833`$ pb and a total branching fraction $`B=5.3\times 10^5`$, using the current expectations for tracking and reconstruction efficiencies. . The shape of the $`m_{J/\psi +\mathrm{}}`$ distribution can be compared with a template of shapes parametrised by the top mass, and the value of $`m_t`$ can therefore be fitted. The Lorentz invariance of the observable makes it completely independent of the details of the $`t\overline{t}`$ production mechanism, and of the structure of the ISR. Given that the spectrum of the leptons from the $`W`$ decay is known very well, the dominant theoretical uncertainty comes therefore from the predictions for the spectrum of the $`J/\psi `$. This is obtained from the convolution of the energy spectrum of $`b`$ hadrons in the top decay, with the spectrum of $`J/\psi `$’s in the $`b`$-hadron decay. In principle this second element can be measured with high accuracy over the next few years at the $`B`$-factories 666Although small corrections are expected to be needed for this application, due to the different composition of $`b`$-hadrons in top decays relative to that in the decays of the $`\mathrm{{\rm Y}}(4S)`$..

We shall therefore concentrate here on the problem of the $`b`$-hadron spectrum in top decays, and investigate how the top mass measurement is affected by the matrix-element corrections to top decays in HERWIG. For simplicity, we shall analyse the $`m_B\mathrm{}`$ spectra (instead of the $`m_{J/\psi +\mathrm{}}`$ ones) obtained by running HERWIG with and without matrix-element corrections to top decays.

The Tevatron statistics will be too low to use this channel as a probe of the top mass. Nevertheless we shall present results for the Tevatron as well, to show that indeed the details of the production mechanism (which is mainly $`q\overline{q}t\overline{t}`$ at the Tevatron and $`ggt\overline{t}`$ at the LHC) have no impact on the top mass determination. Since in $`q\overline{q}`$ annihilation the $`t\overline{t}`$ pair is always produced in a colour octet state, while in $`gg`$ fusion it may come either in a colour singlet or colour octet state, this comparison will indicate that non-perturbative corrections to the $`m_B\mathrm{}`$ mass spectrum are very weakly dependent on the details of the colour-neutralization model in the MC.

We generate $`t\overline{t}`$ samples using HERWIG 6.0 and 6.1 and plot the $`m_B\mathrm{}`$ spectra for different values of $`m_t`$. We then evaluate the average value and the standard deviation of our distributions, the differences in the average values $`m_B\mathrm{}`$ and the corresponding statistical errors. In Fig. 6 the $`m_B\mathrm{}`$ spectra are plotted for $`m_t=175`$ GeV, at the Tevatron and at the LHC, before and after matrix-element corrections, while in Fig. 7 one can find the distributions at the LHC, according to HERWIG 6.1, for $`m_t=171`$ GeV and $`m_t=179`$ GeV. In tables 2 and 3 we summarize the results of our statistical analysis.

We observe a systematic shift of about $`800900`$ MeV towards lower values of $`m_B\mathrm{}`$ after matrix-element corrections to top decays. Furthermore, the results at the Tevatron are the same as the ones at the LHC, to within 150 MeV. Using the Kolmogorov–Smirnov test, we find that about $`N_{\mathrm{eff}}6000`$ reconstructed final states are sufficient to verify that the shapes of the 6.1 and the 6.0 distributions are different at the confidence level of 95%; one year of high-luminosity run would nevertheless allow one to distinguish the two distributions at 70% confidence level.

We studied the dependence of our results on the chosen infrared cutoff for the energy of the gluons emitted in the dead zone, which was discussed in Section 2. We find a negligible variation: at the LHC, and for $`m_t=175`$ GeV, we find $`m_B\mathrm{}=93.78,\mathrm{\hspace{0.17em}93.54}`$ and 93.47 GeV for $`E^{\mathrm{cutoff}}=1,\mathrm{\hspace{0.17em}2}`$ and 5 GeV, respectively. The distributions in Fig. 8, obtained for $`m_t=175`$ GeV and different values of $`E^{\mathrm{cutoff}}`$, are essentially identical.

If we set a cut $`m_B\mathrm{}>50`$ GeV on the invariant mass, to reduce the sensitivy to low-mass tails possibly affected by backgrounds, we find the results summarized in tables 4 and 5 for the Tevatron and the LHC respectively. Once we cut off part of the spectrum, it is to be expected that we find higher values for $`m_B\mathrm{}`$ and lower values for the differences $`m_B\mathrm{}^{6.0}m_B\mathrm{}^{6.1}`$, which are indeed now of about $`400500`$ MeV. Once again the shifts at the Tevatron and at the LHC are of similar size, within 150 MeV.

In order to evaluate the impact of the found discrepancies on the top mass, we perform a linear fit of the $`m_B\mathrm{}`$ distribution as a function of $`m_t`$, by means of the least square method. We find, after considering all the $`m_B\mathrm{}`$ values:

$`6.1:m_B\mathrm{}`$ $`=`$ $`0.563m_t5.087\mathrm{GeV},ϵ\left(\mathrm{GeV}\right)=0.046\left(\mathrm{Tevatron}\right);`$ (1)
$`6.0:m_B\mathrm{}`$ $`=`$ $`0.568m_t5.139\mathrm{GeV},ϵ\left(\mathrm{GeV}\right)=0.023\left(\mathrm{Tevatron}\right);`$ (2)
$`6.1:m_B\mathrm{}`$ $`=`$ $`0.568m_t6.004\mathrm{GeV},ϵ\left(\mathrm{GeV}\right)=0.057\left(\mathrm{LHC}\right);`$ (3)
$`6.0:m_B\mathrm{}`$ $`=`$ $`0.559m_t3.499\mathrm{GeV},ϵ\left(\mathrm{GeV}\right)=0.052\left(\mathrm{LHC}\right);`$ (4)
where $`ϵ`$ is the mean square deviation in the fit. We see that the linear fit is very good, and well within the required accuracy. The obtained fits are plotted in Figs. 9 and 10 at the Tevatron and at the LHC respectively. The error bars correspond to the statistical errors found on $`m_B\mathrm{}^{6.0}m_B\mathrm{}^{6.1}`$, which, as can be seen from the figures, are significantly lower than the difference due to the implementation of matrix-element corrections to top decays. This means that the impact on the top mass is a physical effect and not just the result of statistical fluctuations.

Using these fits, and inverting the relation between $`m_B\mathrm{}`$ and $`m_t`$, we find that for a given value of $`m_B\mathrm{}`$ consistent with the range $`171\stackrel{<}{}m_t\stackrel{<}{}179`$, the values of $`m_t`$ extracted using the two versions of HERWIG differ by 1.5 GeV, both at the Tevatron and at the LHC. This is a rather large value, competitive with the expected systematic error at the LHC, indicating that such corrections are relevant, and must be applied.

After setting a cut of 50 GeV on the $`m_{Bl}`$ spectra, we obtain the following fits:

$`6.1:m_B\mathrm{}`$ $`=`$ $`0.525m_t+6.125\mathrm{GeV},ϵ\left(\mathrm{GeV}\right)=0.049\left(\mathrm{Tevatron}\right);`$ (5)
$`6.0:m_B\mathrm{}`$ $`=`$ $`0.524m_t+6.765\mathrm{GeV},ϵ\left(\mathrm{GeV}\right)=0.022\left(\mathrm{Tevatron}\right);`$ (6)
$`6.1:m_B\mathrm{}`$ $`=`$ $`0.526m_t+5.974\mathrm{GeV},ϵ\left(\mathrm{GeV}\right)=0.026\left(\mathrm{LHC}\right);`$ (7)
$`6.0:m_B\mathrm{}`$ $`=`$ $`0.520m_t+7.578\mathrm{GeV},ϵ\left(\mathrm{GeV}\right)=0.040\left(\mathrm{LHC}\right).`$ (8)
The differences in the top mass extraction between the 6.0 and 6.1 versions from a given value of $`m_B\mathrm{}`$ are now reduced to 1 GeV, smaller than for the fully inclusive distribution, but still significant relative to the overall accuracy goal of 1 GeV.

To conclude our study, we wish to investigate the dependence of our results on the hadronisation model used by HERWIG. We start by studying the invariant mass distribution of the lepton with the $`b`$ quark, as opposed to the $`b`$-hadron. We do this by considering the HERWIG final state at the end of the perturbative evolution, just before the non-perturbative gluon-splitting phase which precedes the formation of the colour-singlet clusters, and the eventual hadronisation.

In Fig. 11 we plot the 6.0 and 6.1 distributions of the invariant mass of the $`b`$-lepton system for $`m_t=175`$ GeV and the 6.1 ones at hadron- and parton-level. The results at parton level for different values of the top mass can be found in table 6. We find that, at fixed $`m_t`$, the average values of the parton-level invariant masses $`m_b\mathrm{}`$ are larger than the ones after the hadronisation of the $`b`$ quark, with differences between the 6.0 and the 6.1 version of the order about $`600700`$ MeV.

After a linear fit, we find the following relations:

$`6.1:m_b\mathrm{}`$ $`=`$ $`0.640m_t10.256\mathrm{GeV},ϵ\left(\mathrm{GeV}\right)=0.050.`$ (9)
$`6.0:m_b\mathrm{}`$ $`=`$ $`0.646m_t10.528\mathrm{GeV},ϵ\left(\mathrm{GeV}\right)=0.033`$ (10)
Given the values of the slopes, the shift in the extracted top mass between the two versions is now of about 1.0 GeV, compared to the 1.5 GeV after hadronisation.

As a whole, the non-perturbative contribution to the shape of the $`b`$-hadron spectra is therefore pretty important, as is already known in the case of $`Z^0`$ decays. Confidence in the accurate description of the non-perturbative phase should be gained from the study of the $`B`$-fragmentation function in $`Z^0`$ decays. Our study of the impact of matrix element corrections at the Tevatron and at the LHC suggests that non-perturbative corrections do not depend significantly on the production mechanism. This we also checked by performing a similar study in the case of $`e^+e^{}`$ production. We therefore expect that once a tuning of the non-perturbative $`b`$ fragmentation function in HERWIG is achieved, using for example the latest high-precision SLD results , the results can be extended to the study of the $`B\mathrm{}`$ mass spectrum in top decays. Given the size of the mass shift induced by hadronisation (of the order of 8 GeV, as found by comparing Table 6 and 3), a control over the fragmentation function at a level better than 5% will have to be achieved in order to maintain this contribution to the theoretical systematics well below $`𝒪`$(1 GeV). This should be possible, since the current size of the experimental uncertainty on $`x_B`$ is at the level of $`1`$% .

## 5 Conclusions

We have studied $`t\overline{t}`$ events in the dilepton channel at the Tevatron and at the LHC using the new version of the HERWIG Monte Carlo event generator, which includes matrix-element corrections to the description of top decays. We considered observables involving the third hardest jet in transverse energy, to enhance possible effects of the implemented corrections.

We have found that the distributions obtained before and after matrix-element corrections are rather similar. A Kolmogorov test, which compares the shapes of two distributions, allows however to detect differences at the 95%CL with the large statistics available at the LHC. The statistics of the Run II at the Tevatron are however not sufficient.

We have also investigated the reconstruction of the top mass by looking at final states with leptons and a $`J/\psi `$. These final states are an excellent candidate for an experimental determination of the top mass with systematic errors in the range of 1 GeV. While the main production mechanisms of top quarks are different at the Tevatron and at the LHC, we found equivalent results in the two cases, indicating that this method of reconstruction of the top mass is not sensitive to the details of the colour-neutralisation model.

We considered the spectra of the invariant mass $`m_B\mathrm{}`$, where the $`B`$ meson comes from the hadronisation of the $`b`$ quark produced in $`tbW`$ and $`\mathrm{}`$ is the charged lepton from the decay $`W\mathrm{}\nu `$, for different values of $`m_t`$ and obtained that the implementation of matrix-element corrections to top decays results in a shift of about 1.5 GeV on the top mass if one is able to reconstruct the whole $`m_B\mathrm{}`$ spectrum and of about 1 GeV after setting the cut $`m_B\mathrm{}>`$ 50 GeV. The shifts we found are physical effects related to the inclusion of hard and large-angle gluon radiation in the Monte Carlo shower, since they are much larger than the statistical errors on them. Analyses at the parton level have shown an impact of a similar magnitude.

It will be now very interesting to compare the new HERWIG results for the invariant-mass distributions with the ones obtained after performing a next-to-leading order calculation of the process $`tbWg`$, and convoluting the result with the fragmentation function for the hadronisation of the $`b`$ quark into a $`B`$ meson taken from LEP and SLD data. In order to perform such a comparison in detail, the HERWIG cluster model used to simulate the hadronisation process will have to be tuned to fit that data.

Furthermore, although in this paper we have concentrated our analysis on the top mass reconstruction in final states with leptons and $`J/\psi `$, it will be worthwhile redoing the Tevatron analysis to determine $`m_t`$ using the HERWIG parton shower model, provided with matrix-element corrections.

We finally recall that, though we feel safely confident that the new version of HERWIG will be a trustworthy event generator for the purpose top decays, top production is still performed in the soft/collinear approximation, with dead zones in the phase space which are needed to be filled. Matrix-element corrections to top production are in progress and may have an impact on jet observables and on the top mass measurement at the Tevatron as well.

## Acknowledgements

We wish to acknowledge Bryan Webber for useful discussions on these and related topics. We thank A. Kharchilava for several discussions on the measurement of the top mass using the $`\psi \mathrm{}`$ final states. G.C. is also grateful to the R.A.L. Theory Group and to the CERN Theory Division for hospitality during some of this work.

## References

1. S. Frixione, M.L. Mangano, P. Nason and G. Ridolfi, Heavy Quark Production, hep-ph/9702287, in Heavy Flavours II, ed. A.J. Buras and M. Lindner, World Scientific, p.609;
S. Frixione, M. L. Mangano, P. Nason and G. Ridolfi, Nucl. Phys. B431 (1994) 453.
2. M. Beneke et al., hep-ph/0003033, to appear in Proceedings of 1999 CERN Workshop on Standard Model Physics (and more) at the LHC, G. Altarelli and M.L. Mangano eds.
3. E. Accomando et al. \[ECFA/DESY LC Physics Working Group Collaboration\], Phys. Rept. 299, 1 (1998) \[hep-ph/9705442\].
4. M. L. Mangano, P. Nason and G. Ridolfi, Nucl. Phys. B373, 295 (1992);
S. Frixione, M. L. Mangano, P. Nason and G. Ridolfi, Phys. Lett. B351, 555 (1995) \[hep-ph/9503213\].
5. M. Cacciari, M. Greco and P. Nason, JHEP 9805, 007 (1998) \[hep-ph/9803400\];
F. I. Olness, R. J. Scalise and W. Tung, Phys. Rev. D59, 014506 (1999) \[hep-ph/9712494\].
6. G. Marchesini and B.R. Webber, Nucl. Phys. B310 (1988) 461. G. Marchesini et al., Comput. Phys. Commun. 67 (1992) 465.
7. T. Sjöstrand, Comp. Phys. Comm. 46 (1987) 367.
8. F. E. Paige and S. D. Protopopescu, BNL-29777; H. Baer, F. E. Paige, S. D. Protopopescu and X. Tata, ISAJET 7.48, hep-ph/0001086.
9. M.H. Seymour, Comput. Phys. Commun. 90 (1995) 95.
10. M.H. Seymour, Z. Phys. C56 (1992) 161.
11. M.H. Seymour, Matrix Element Corrections to Parton Shower Simulation of Deep Inelastic Scattering, contributed to 27th International Conference on High Energy Physics (ICHEP), Glasgow, 1994, Lund preprint LU-TP-94-12, unpublished.
12. G. Corcella and M.H. Seymour, Phys. Lett. B442 (1998) 417.
13. G. Corcella and M.H. Seymour, Nucl. Phys. B565 (2000) 227.
14. J. André and T. Sjöstrand, Phys. Lett. B442 (1998) 417
G. Miu and T. Sjöstrand, Phys. Lett. B449 (1999) 313.
15. A. Kharchilava, Phys. Lett. B476 (2000) 73.
16. G. Corcella, hep-ph/9911447, On the Top Mass Reconstruction Using Leptons, contributed to the U.K. Phenomenology Workshop on Collider Physics, Durham, U. K., 19-24 September 1999.
17. G.Marchesini and B.R. Webber, Nucl. Phys. B330 (1990) 261
18. G. Marchesini and B.R. Webber, Nucl. Phys. B238 (1984) 1.
19. G. Corcella, I.G. Knowles, G. Marchesini, S. Moretti, K. Odagiri, P. Richardson, M.H. Seymour, B.R. Webber, hep-ph/9912396.
20. L.H. Orr, T. Stelzer and W.J. Stirling, Phys. Rev. D56 (1997) 446
21. S. Catani, Yu.L. Dokshitzer, M.H. Seymour and B.R. Webber, Nucl. Phys. B406 (1993) 187
M.H. Seymour, Nucl. Phys. B421 (1994) 545.
22. S.D. Ellis and D.E. Soper, Phys. Rev. D48 (1993) 3160.
23. R. Bonciani, S. Catani, M.L. Mangano and P. Nason, Nucl. Phys. B529 (1998) 424.
24. The CDF II Detector, Technical Design Report, FERMILAB-Pub-96/390-E.
25. ATLAS detector and physics performance, Technical Design Report, Volume I, CERN LHCC 99-14.
26. K. Abe et al. \[SLD Collaboration\], hep-ex/9912058.


format: text

text:
# DIRECT OBSERVATION OF THE FOURTH STAR IN THE ZETA CANCRI SYSTEM

## 1 Introduction

$`\zeta `$ Cancri is a visual triple stellar system which has been studied as such for over 200 years. The system consists of three bright stars, all dwarfs of close to solar type with visible magnitudes near 6. The orbital planes of the system are not far from the plane of the sky. There is a close pair AB separated by about 1 arcsec with an orbital period 58.9 years, and a third star C about 6 arcsec away with an orbital period estimated at 1100 years. Star C itself has proper-motion irregularities that were first remarked upon by Sir John Herschel (1831) and were interpreted by Otto Struve (1874) in terms of orbital motion in a period on the order of 20 years around an unseen companion. Subsequent observations and discussions have corroborated Struve’s interpretation, and have recently been reinforced by a spectroscopic orbit of C by Griffin (2000), who has also described in some detail the interesting observational history of the system. The situation concerning the three visual orbits in which the $`\zeta `$ Cnc stars are involved (those of AB, CD, and AB–CD) have been summarized most recently by Heintz (1996). All three have appreciable eccentricity, with that of CD being lowest at 0.12.

The principal question that has remained outstanding concerning the $`\zeta `$ Cnc system is the nature of the companion to star C, since its orbit indicates a companion of comparable mass and a separation of about 0.4 arcsec. Griffin (2000) recounts some of the earlier claims regarding the nature of star D, and it is clear that it must be considerably fainter than C at visible wavelengths. Suggestions have consequently been made that it is a white dwarf, or else a close pair of M stars which would be quite faint in the visible. Griffin concluded that the question remained open, but that a short exposure by a large telescope with adaptive optics might provide the definitive evidence. Accordingly, the system was observed during a run with the Canada France Hawaii telescope by one of us (JBH) in January 2000, and followed up by further data from the same system a month later by another of us (FM). As we describe below, we discovered star D immediately and easily.

## 2 Observations and data

The observations of $`\zeta `$ Cnc were carried out during runs of the adaptive optics (AO) camera of the Canada France Hawaii telescope (CFHT), with the near-infrared camera KIR (see Rigaut et al 1998). They were obtained during gaps in other programs, using filters that happened to be in the camera. Thus, as described below, the observing program was not particularly well planned, but it occupied a total exposure time of under one minute! The field of view of 35 easily includes all of the stars in the system. The brightness of the stars made the AO correction extremely good, with FWHM close to the diffraction limit at all wavelengths. The filters used were generally narrow-band, to avoid detector saturation and the need for very short exposures. Table 1 lists the details of the observations.

The first observations by one of us (JBH) in January 2000 used the AB pair as guide star. That proved not to be a happy choice: the two bright stars are separated by less than 1′′ and the guiding system ‘hunted’ back and forth between the two, producing double images of all the stars present, all with attendant diffraction rings! The ‘PSF’ pattern also changed depending on where in the field of view the guide star(s) were, during the dither pattern used. In spite of the complex images produced, it was clear that a companion star to C (star D) was present, and it was possible to measure its position and brightness by subtracting the complex PSF. The process was fortuitously made simpler by the fact that star D is quite well separated from C along the same direction as the separation of AB (and hence the double PSF).

Since it was clear that better results would be obtained by guiding on star C, which is effectively single at the visible wavelengths used for guiding, and in view of the possibility of measuring the relative motions of the stars over an interval of several weeks, further observations were obtained in February by one of us (FM). The AO system produced very well-corrected images from the bright guide star, with FWHM of 0′′.08, 0′′.10, and 0′′.12 at $`J`$, $`H`$, and $`K`$; they are shown in Figure 1.

Since the stars are so bright and the exposures so short, there was no need for flat fields or sky subtraction in order to measure the star fluxes and positions. Small changes seen in the diffraction patterns are due to the shortness of the exposure times. The residual seeing halo due to any incompleteness of the AO correction is resolved into individual speckles, as is typical for short exposures on AO systems. However, we established that measurements from images dithered to different places on the detector produced values that differed by amounts too small to matter in this investigation.

## 3 Measurements and discussion

For each of the observing runs, the positions and fluxes of the stars were determined, using all the images taken. Fluxes were measured by the IRAF task ‘imexam’, and also by summing the signal from sections of image and subtraction of the off-star signal levels. Star positions were determined by imexam, and also by visual inspection of contour plots over the central pixels of each star image.

The camera was removed and re-installed between the two observing runs. Thus, it was necessary to check the orientation of the images for small rotations. Since systems AB and CD moved in their orbits over the interval, we used the position angle of BC as the fiducial, since its predicted change in the interval is only 0.03 (Heintz 1996). That check showed that a marginal rotation between runs of 0.05$`\pm `$0.08 had occurred. This correction has been applied to the position-angle differences we report.

Our main interest is in relative fluxes between the system stars, so no absolute calibration was applied, except to check that the overall colors between visible and IR are as expected for the bright stars, which are all of almost the same spectral type. Table 2 shows the measurements made of fluxes and Table 3 shows the position and angle measurements.

Table 2 compares the color differences expected for main-sequence stars of various likely types, from Bessell & Brett (1988). They have been corrected from the standard $`JHK`$ values to account for the different bandpass centres of the Fe II and Br$`\gamma `$ filters. The spectral type of star C is quoted as G5 in the Bright Star Catalogue (Hoffleit 1982) but Griffin (2000) has argued that it is based upon a misconception and that G0 is more likely. As A and B are bright and overlapping, while C and D have different fluxes, the color differences C$``$D are more reliable than A$``$D or B$``$D.

As a figure of merit, we sum the differences between model and observed colors in Table 2. For B$``$D, that indicates that M2 is better than M0 for D, and for C$``$D G0+M2 fits better than G0+M0 or G5+M2. Thus, we conclude from the color differences that C and D have spectral types G0 and M2, respectively, with an uncertainty only on the order of one spectral subtype.

The mass ratio C/D is close to unity — Heintz (1996) gives masses of 0.99 and 0.93 for C and D — which does not accord with the mass ratio of single main-sequence stars of types G and M. However, if D is itself a binary system of two stars close to spectral type M2, the mass would be 0.78, and for an M0 pair it would be 0.94 (see Allen 1973). From the mass function of 0.05 found by both Griffin and Heintz, such a combined mass corresponds to an orbital inclination of about 35.

The angle changes noted in Table 3 are fully consistent with the published orbits of Heintz (1996). Note that the orientation of our detector is based on adopting the predicted angle of BC from Heintz’s orbits. The ‘predictions’ in the table are derived from Heintz’s tables and periods, assuming a mass ratio near unity for CD. The observed separation and direction of D from C, and their changes, indicate that to be the case. Monitoring of the CD visual orbit with further AO images will enable much better estimates to be made of the mass ratio.

Given that the stars are all at the same distance, absolute magnitudes can be used to predict the magnitude difference between C and D in the $`V`$ band as 4.5 for two M2 stars and 3.5 for two M0 stars, compared with G0 (Allen 1973). Applying the standard ($`VK`$) colors to the difference we measure between stars C and D we derive $`V`$ magnitude differences of 4.1 and 3.5 for pairs of M2 and M0 stars, respectively. Thus, the colors and relative fluxes are fully consistent with stars D being an unresolved binary of two M stars. A challenge for future observations will be to resolve D at visual wavelengths and measure its $`V`$ magnitude, and even to demonstrate (directly and/or spectroscopically) that it is double.

We are grateful to David Crampton for suggesting the observations to the CFHT observer after a conversation following a colloquium on $`\zeta `$Cnc given by one of us (RFG).

References

Allen C.W., 1973, Astrophysical Quantities, Athlone Press, London (third edition).

Bessell M.S., & Brett J.M., 1988, PASP, 100, 1134

Griffin R.F., 2000, Observatory, 120, 1

Heintz W.D., 1996, AJ, 111, 408

Herschel J.F.W., 1831, Mem.RAS, 5, 171

Hoffleit D., 1982, The Bright Star Catalogue, Yale University Observatory

Rigaut F., et al, 1998, PASP, 110, 152

Struve O., 1874, Comptes Rendus, 79, 1463

Caption to Figure 1

Images of $`\zeta `$ Cnc guided on star C. N is 2.8 left of vertical, and E to the left. Filters are $`J`$ (top left, 1.2$`\mu `$), Fe II cont (top right, 1.7$`\mu `$), Br$`\gamma `$ (lower left, 2.2$`\mu `$). Lower right is detail of CD pair with Fe II cont filter. CD separation is 0′′.32.


text:
# Critical currents in Josephson junctions, with unconventional pairing symmetry: 𝑑_{𝑥²-𝑦²}+𝑖𝑠 versus 𝑑_{𝑥²-𝑦²}+𝑖𝑑_{𝑥𝑦}

## I Introduction

One of the main questions in the research activity on high-$`T_c`$ superconductors nowadays is the identification of the order parameter symmetry and its underlying mechanism. The most possible scenario is that the bulk pairing state has a $`d_{x^2y^2}`$-wave character. Theoretical calculations, suggest that an imaginary $`s`$-wave component which breaks the time reversal symmetry is induced in some cases, wherever the $`d_{x^2y^2}`$-wave order parameter varies spatially such as near a vortex, or near the surface . Also the observation of fractional vortices on a triangular grain boundary in YBa2Cu3O7 by Kirtley $`et`$ $`al.`$ , may indicate a possible violation of the time-reversal symmetry near grain boundary. Theoretical explanation of this experiment is given by Bailey et. al. in Ref.  where they study a triangular grain boundary in $`d`$-wave superconductors. They conclude that under the assumption of $`d`$-wave symmetry, the flux at the edges of this triangle can take the values $`\pm \mathrm{\Phi }_0/2`$, which does not agree with the experiment. However under the assumption of $`d_{x^2y^2}+is`$-wave symmetry an intrinsic phase shift $`\varphi _c(x)`$ exists in each triangle edge. In turn the phase $`\varphi (x)`$ must change in order to connect the different values of $`\varphi _c`$ in each segment. This arrangement leads to fractional vortices or antivortices at each three corners, in agreement with the experiment.

Another pairing state which breaks the time reversal symmetry is the $`d_{x^2y^2}+id_{xy}`$-wave. Patches of complex $`d_{xy}`$ components are induced around magnetic impurities at low temperatures in a $`d_{x^2y^2}`$-wave superconductor forming a phase coherent state as a result of tunneling between different patches . Violation of parity and time reversal symmetry occurs in this state. Also on the high field region, $`HH_{c2}`$ the $`d_{x^2y^2}`$-wave state can be perturbed by the external filed, producing a $`d_{x^2y^2}+id_{xy}`$ state in the bulk .

The observation of the splitting of the zero energy peak in the conductance spectra at low temperatures indicates that a secondary component is induced which violates locally the time reversal symmetry . Theoretical explanation based on surface-induced Andreev states, has been proposed . Recently the field dependence of this splitting has been observed in the tunneling spectra of YBCO . This observation is consistent with a $`d_{x^2y^2}+is`$ surface order parameter as well as a $`d_{x^2y^2}+id_{xy}`$ bulk order parameter. Another question which can be asked is to what extend, the observation of a symmetric magnetic interference pattern in the corner junction experiments  is an identification of $`d_{x^2y^2}`$-wave symmetry, or could also imply a $`d_{x^2y^2}+id_{xy}`$ pairing state also? In this work we propose a phase sensitive experiment based on the Josephson effect, which may be used to distinguish which symmetry ($`d_{x^2y^2}+is`$ or $`d_{x^2y^2}+id_{xy}`$) the order parameter has near the surface. We study the static properties of a frustrated junction which is made of two one-dimensional junctions, of $`d_{x^2y^2}+ix`$-wave, (where $`x=d_{xy}`$ or $`x=s`$) and $`s`$-wave superconductors. By introducing an extra relative phase in one part of this junction, the above junction can be mapped into the corner junctions experiments. We examine the spontaneous flux and the critical current modulation of the vortex states with the junction orientation angle $`\theta `$, the magnitude of the secondary component $`n_s`$, and the magnetic field $`H`$. In each case we derive simple arguments which are useful to discriminate between the time reversal symmetry broken states. For example, when the orientation is exactly such that the lobes of the dominant $`d_{x^2y^2}`$-wave order parameter points towards the junction interface the magnetic interference pattern is symmetric (asymmetric) when the secondary order parameter is $`x=d_{xy}(s)`$. This is verified for small junctions as well as in the long junction limit, and can be used to distinguish between broken time reversal symmetry states.

The rest of the paper is organized as follows. In Sec. II we discuss the Josephson effect between a superconductor with broken time reversal symmetry and an $`s`$-wave superconductor. In Sec. III the geometry of the corner junction is discussed. In Sec. IV we present the results for the magnetic flux of the spontaneous vortex states in corner junctions with some intrinsic magnetic flux. In Sec. V the parameters which can modulate the spontaneous flux and the critical currents are considered. In Sec. VI a connection with the experiment is made. Finally, a summary and discussion are presented in the last section.

## II Josephson effect between two superconductors with mixed wave symmetry

We consider the junction shown in Fig. 1(a), where two superconductors ($`A`$ in the region $`z>t`$ and $`B`$ in the region $`z<0`$), are separated by an intermediate layer. We assume that each superconductor has a two component order parameter. The order parameter for each component $`k(k=1,2)`$ in the superconductors, can be written as

$$n_k=\{\begin{array}{cc}\stackrel{~}{n}_k^Ae^{i\varphi _k^A},\hfill & z>t\hfill \\ \stackrel{~}{n}_k^Be^{i\varphi _k^B},\hfill & z<0\hfill \end{array}.$$
(1)
Here $`\varphi _k^{A(B)}`$ is the phase of the order parameter $`n_k`$ in superconductor $`A(B)`$. Then phenomenological Ginzburg-Landau theory is used to calculate the supercurrent density given by

$$J=\underset{k,l=1}{\overset{2}{}}J_{ckl}\mathrm{sin}(\varphi _k^B\varphi _l^A),$$
(2)
where

$$\begin{array}{cc}J_{c11}=\hfill & (2e\mathrm{}/m_a^{}d)\stackrel{~}{n}_1^A\stackrel{~}{n}_1^B\hfill \\ J_{c21}=\hfill & (2e\mathrm{}/m_\nu ^{}d)\stackrel{~}{n}_1^A\stackrel{~}{n}_2^B\hfill \\ J_{c12}=\hfill & (2e\mathrm{}/m_\nu ^{}d)\stackrel{~}{n}_2^A\stackrel{~}{n}_1^B\hfill \\ J_{c22}=\hfill & (2e\mathrm{}/m_b^{}d)\stackrel{~}{n}_2^A\stackrel{~}{n}_2^B\hfill \end{array},$$
(3)
$`m_a^{},m_\nu ^{},m_b^{}`$ are the effective masses that enter into the Ginzburg-Landau equations. In the following these masses are taken equal to an effective mass $`m^{}`$.

We restrict to the case where $`B`$ is $`s`$-wave. In this case $`\stackrel{~}{n}_1^B=0`$, and $`\stackrel{~}{n}_2^B`$=constant. We define $`\varphi =\varphi _2^B\varphi _1^A`$, as the relative phase difference between the two superconductors. We consider the case where the intrinsic phase difference within superconductor $`A`$ is $`\varphi _2^A\varphi _1^A=\pi /2`$. Then the order parameter in $`A`$ is complex and breaks the time reversal symmetry. The supercurrent density can be written as:

$$J(\varphi )=\stackrel{~}{J}_c\mathrm{sin}(\varphi +\varphi _c),$$
(4)
with

$$\stackrel{~}{J}_c=\sqrt{J_1^2+J_2^2},$$
(5)
$$\varphi _c=\{\begin{array}{cc}\mathrm{tan}^1\frac{J_2}{J_1},\hfill & J_1>0\hfill \\ \pi +\mathrm{tan}^1\frac{J_2}{J_1},\hfill & J_1<0\hfill \end{array},$$
(6)
where $`J_1=J_{c21}`$, $`J_2=J_{c22}`$. The Josephson critical current density $`\stackrel{~}{J}_c`$ is scaled in units of $`J_{c0}=\frac{e\mathrm{}}{m^{}d}`$. Two special cases are the following:

i) For $`d_{x^2y^2}+is`$ wave case the magnitude of the $`d_{x^2y^2}`$-wave component in (1) is $`\stackrel{~}{n}_1^A=n_{10}\mathrm{cos}(2\theta )`$ , where $`\theta `$ is the angle of the crystalline $`a`$-axis of superconductor $`A`$ with the junction interface. The magnitude of the secondary order parameter in superconductor $`A`$ is $`\stackrel{~}{n}_2^A=n_{20}=0.1n_{10}`$.

ii) For $`d_{x^2y^2}+id_{xy}`$ wave case, the magnitude of the $`d_{x^2y^2}`$-wave component in (1) is given by $`\stackrel{~}{n}_1^A=n_{10}\mathrm{cos}(2\theta )`$, while the $`d_{xy}`$ wave component is $`\stackrel{~}{n}_2^A=n_{20}\mathrm{sin}(2\theta )`$, where $`n_{20}=0.1n_{10}`$. This order parameter can occur in the following way: The order parameter magnitude for the $`d`$-wave state $`\mathrm{\Delta }_0(\theta )=\mathrm{\Delta }_0\mathrm{cos}(2\theta )`$ is an equal admixture of pairs with orbital moment $`L_z=\pm 2`$, and can be written as $`\mathrm{\Delta }_0(\theta )=(\mathrm{\Delta }_0/2)[\mathrm{exp}(2i\theta )+\mathrm{exp}(2i\theta )]`$. In the presence of perturbation such as (ferromagnetically) ordered impurity spins $`S_z`$ the coefficients of $`L_z=\pm 2`$ components will shift linearly in $`S_z`$ with opposite signs. The final state will be $`\mathrm{\Delta }_0(\theta )\mathrm{\Delta }_0(\theta )+iS_z\mathrm{\Delta }_1(\theta )`$, where $`\mathrm{\Delta }_1(\theta )=\mathrm{sin}(2\theta )`$. The strength of the secondary component is proportional to the perturbation $`S_z`$.

## III The corner junction geometry

We consider the corner junction shown in Fig. 1(b), between a superconductor with broken time reversal symmetry at the surface and an s-wave superconductor. If the angle of $`a`$-axis with the interface in the $`x`$-direction is $`\theta `$, then the corresponding angle in the $`z`$-direction will be $`\pi /2\theta `$. We map the two segments each of length $`L/2`$ where $`L=10\lambda _J`$ of this junction into a one-dimensional axis. In this case the two dimensional junction can be considered as being made of two one dimensional junctions described in Sec. II connected in parallel. Their characteristic phases $`\varphi _{c1}`$ and $`\varphi _{c2}`$ depend upon the angle $`\theta `$. We call this junction frustrated since the two segments have different characteristic phases $`\varphi _{c1},\varphi _{c2}`$. The fabrication details of corner junctions or superconducting quantum interference device (SQUID), between sample faces at different angles can be found in Ref. .

The superconducting phase difference $`\varphi `$ across the junction is then the solution of the sine-Gordon (s-G) equation

$$\frac{d^2\varphi (x)}{dx^2}=\stackrel{~}{J}_c\mathrm{sin}[\varphi (x)+\varphi _c(x)]I^{ov},$$
(7)
with the boundary conditions

$$\frac{d\varphi }{dx}|_{x=0,L}=H.$$
(8)
The length $`x`$ is scaled in units of the the Josephson penetration depth given by

$`\lambda _J=\sqrt{{\displaystyle \frac{\mathrm{}c^2}{8\pi edJ_{c0}}}},`$

where $`d`$ is the sum of the $`s`$-wave, and mixed wave $`\lambda _{ab}`$ penetration depths plus the thickness of the insulator layer. The relative phase $`\varphi _c(x)`$ is $`\varphi _{c1}(\varphi _{c2})`$ in the left (right) part of the junction. The external magnetic field $`H`$, scaled in units of $`H_c=\frac{\mathrm{}c}{2ed\lambda _J}`$ is applied in the $`y`$ direction, which is considered small compared to $`\lambda _J`$. The bias current per unit length $`I^{ov}`$ in the overlap geometry is scaled in units of $`\frac{c}{4\pi }H_c`$, and is uniformly distributed along the entire $`x`$ axis of the junction.

We can classify the different solutions obtained from Eq. (7) with their magnetic flux content

$$\mathrm{\Phi }=\frac{1}{2\pi }(\varphi _R\varphi _L),$$
(9)
where $`\varphi _{R(L)}`$ is the value of $`\varphi `$ at the right(left) edge of the junction, in units of the flux quantum $`\mathrm{\Phi }_0=\frac{hc}{2e}`$.

## IV Spontaneous vortex states

Firstly let us consider the case where the two one-dimensional junctions of $`d_{x^2y^2}+ix`$-wave where $`x=s`$ or $`x=d_{xy}`$, and $`s`$-wave supeconductors, each of length $`L/2`$, described in Sec. II are uncoupled. Then for $`0<x<L/2`$ the stable solutions for the s-G equation are $`\varphi (x)=\varphi _{c1}+2n_1\pi `$, where $`n_1=0,\pm 1,\pm 2,\mathrm{}`$, while for $`L/2<x<L`$ the stable solutions for the s-G equation are $`\varphi (x)=\varphi _{c2}+2n_2\pi `$, where $`n_2=0,\pm 1,\pm 2,\mathrm{}`$, where $`\varphi _{c1}`$, $`\varphi _{c2}`$, are the relative phases in each part of the junction due to different orientations. These solutions are plotted in Fig. 2, for $`n_i=0,1`$, $`i=1,2`$ as a function of the orientation angle $`\theta `$. When the frustrated junction is formed, and we consider the above junctions in parallel, the phase $`\varphi `$ is forced to change around $`x=L/2`$, to connect these stable solutions. This variation of the phase $`\varphi `$, along the junction describes the Josephson vortices. The flux content of these states (in units of $`\mathrm{\Phi }_0`$) is

$$\mathrm{\Phi }=[\varphi (L)\varphi (0)]/2\pi =(\varphi _{c2}+\varphi _{c1}+2n\pi )/2\pi ,$$
(10)
where the $`n`$-value ($`n=n_1n_2=0,\pm 1,\pm 2,\mathrm{}`$) distinguishes between solutions with different flux content. We will concentrate to solutions called modes with the minimum flux content i.e., $`n=0,1,1`$. Their magnetic flux in terms of $`\varphi _{c1},\varphi _{c2}`$ is shown in table I. Generally the flux content is fractional i.e. is neither integer nor half-integer, as a consequence of the broken time reversal symmetry of the problem.

In the actual numerical simulations, the stable solutions of the sine-Gordon equation in the left(right) part of the junction are taken as the initial conditions for the iteration procedure. For example for the $`n=0`$ solution the phase $`\varphi (x)`$ is taken $`\varphi (x)=\varphi _{c1}`$ $`(\varphi _{c2})`$ in the left (right) part of the junction, as an initial condition and then is iterated until convergence. Besides if we take as initial condition , $`\varphi (x)=\varphi _{c1}`$, in the left side, and $`\varphi (x)=2\pi \varphi _{c2}`$ in the right side, the final state of the system, after the iteration procedure, is the solution which we call $`n=1`$, with negative magnetic flux, and not exactly opposite to $`n=0`$. We comment here that the solutions after the iteration procedure have smooth variation as a function of the position, as opposed to the step function variation of the initial conditions.

For the $`00`$ junction $`\varphi _{c1}=\varphi _{c2}=0`$, and the flux becomes $`\mathrm{\Phi }=n`$, so we say that the flux is quantized in integer units of $`\mathrm{\Phi }_0`$. In this case, there exist solutions with flux $`\mathrm{\Phi }=\mathrm{},1,0,1,\mathrm{}`$ . These solutions, when $`n0`$ are stabilized by the application of an external magnetic field. In the case of a junction with some spontaneous flux, at least for the modes with lower flux content, the external field is not necessary since the spontaneous magnetization state is stable.

In the case of $`0\pi `$ junction, where the intrinsic phase in the right (left) part of the junction is $`\varphi _{c2}=\pi `$ ($`\varphi _{c1}=0`$), the stable solutions of the s-G equation are $`\varphi (x)=2n\pi `$ for the left part, while $`\varphi (x)=\pi (2n+1)`$ for the right part of the junction. In this case a $`0\pi `$ junction is formed. The corresponding flux becomes $`\mathrm{\Phi }=(n+1/2)\pi `$, and the particular values of $`n=0`$, $`n=1`$ give the half vortex and antivortex solutions, with opposite fluxon content, $`\mathrm{\Phi }=0.5`$ and $`\mathrm{\Phi }=0.5`$ respectively.

## V Magnetic flux and critical current modulation

In the following we will describe three parameters which can alter the spontaneous flux and the critical currents of the vortex states described in the previous section, in a corner junction between a superconductor with time reversal broken symmetry and an $`s`$-wave superconductor. These include the orientation angle $`\theta `$, the magnitude of the secondary order parameter $`n_s`$, and the magnetic field $`H`$. In each parameter separately we will point out the differences between the $`d_{x^2y^2}+is`$-wave, and $`d_{x^2y^2}+id_{xy}`$-wave.

### A Junction orientation

For the $`d_{x^2y^2}+is`$-wave case, we consider first the situation where $`\theta `$ is varied from $`0`$ to $`\pi /2`$. In Fig. 3 we plot the spontaneous magnetic flux versus the angle ($`\theta `$) for the different modes $`n=0,1,1`$ in the corner junction geometry. As we can see the magnetic flux changes with orientation. For angle $`\theta `$ close to $`0`$ or $`\pi /2`$ the spontaneous modes existing at $`H=0`$ are separated by an integer value of the magnetic flux. This is also the case in the pure $`s`$-wave superconductor junction problem. The difference is that the modes are found displaced to fractional values of magnetic flux, contrary to the $`s`$-wave case where the magnetic flux takes on integer values at $`H=0`$. In particular the vortex solution in the $`n=0`$ mode (solid line) contains less that half a fluxon for $`\theta =0`$, and as we increase the angle $`\theta `$ towards $`\pi /4`$ it continuously reduces its flux, i.e. it becomes flat exactly at $`\theta =\pi /4`$ and then it reverses its sign and becomes an antivortex with exactly opposite flux content at $`\theta =\pi /2`$ from that at $`\theta =0`$. In addition we have plotted in Fig. 4a the phase distributions for the mode $`n=0`$ in different orientations $`\theta =0,\pi /4,\pi /2`$. The transition form the vortex to the antivortex mode as the orientation changes is clearly seen in this figure. Note that the solutions in this mode remain stable for all junction orientations. This is seen in Fig. 5 where we plot the lowest eigenvalue ($`\lambda _1`$) of the linearized eigenvalue problem as a function of the angle $`\theta `$ . We see that $`\lambda _1>0`$, denoting stability for all values of the angle $`\theta `$ in this mode.

Let as now examine the solution in the $`n=1`$ mode, (dotted line in Fig. 3). We see that at $`\theta =0`$ it has negative flux, which in absolute value is more than $`\mathrm{\Phi }_0/2`$ and as we increase the angle $`\theta `$ it decreases its flux to a full antifluxon when the orientation is slightly greater than $`\pi /4`$ and than to flux greater than $`\mathrm{\Phi }_0`$ when $`\theta `$ reaches $`\pi /2`$. As seen in Fig. 5 this solution becomes unstable at a point to the left of $`\theta =\pi /4`$ (point $`\iota `$) due to the abrupt change of flux at this angle. More strictly the instability sets in due to the competition between the slope of the phase at the edges of the junction and at the junction center as the angle $`\theta `$ approaches the value $`\pi /4`$. At this point the slope competition makes the antivortex unstable. This is seen in Fig. 4b) (dotted line) where the phase distribution for the $`n=1`$ mode solution is plotted at the point where the instability starts i.e. $`\theta =0.242\pi `$.

Finally the solution in the $`n=1`$ mode contains more than one fluxon at $`\theta =0`$ and is clearly unstable. It becomes stable at an angle slightly on the right of $`\theta =\pi /4`$, (point $`\nu `$ in Fig. 5) where the flux varies more smoothly, see $`\theta =0.258\pi `$ in Fig. 4c. At $`\theta =\pi /2`$ it contains more than $`\mathrm{\Phi }_0/2`$ in flux. We expect a time reversal broken symmetry state like $`d_{x^2y^2}+is`$ to be characterized by either the solution in the fractional vortex or antivortex mode, because due to the different character of these solutions a change from one variant to the other would demand the application of an external current or magnetic field and in this sense it would cost additional energy. So since these states are stable in external perturbations, once the system is prepared in one of these it will remain to that state.

In general we see that for each value of $`\theta `$ there exist in the junction a pair of stable solutions which when applying an external bias current will lead to observable critical currents. In Fig. 6 we plot the overlap critical current per unit length $`I_c^{ov}`$ as a function of $`\theta `$, at $`H=0`$, for the $`n=0,1,1`$-mode solutions, in the $`d_{x^2y^2}+is`$-wave case. In the overlap geometry the current is distributed in the entire $`x`$-axis. In the calculations we have taken into account that the Josephson critical current density $`\stackrel{~}{J}_c`$ has a characteristic variation with the orientation. We find that for a given orientation it is possible for the junction current density to vary in the way that several modes with different critical currents can exist. In Fig. 7 we plot the current density when the total current is maximum, for different modes, and orientations, which will give us information about the actual shapes of the vortices. Let us consider the situation where the junction contains a solution in the mode $`n=0`$, at $`\theta =0`$, when the net current is maximum. The spatial variation of $`\varphi `$ is described by a fractional vortex which is displaced around the value $`\varphi =\pi `$, from the corresponding distribution at zero current which is around $`\pi /2`$ (see Fig. 4a). The current density distribution as seen in Fig. 7a (solid line) at the maximum current is flat above unit with a small variation around the junction center giving rise to the large value on the net current, seen in Fig. 6. Also at $`\theta =\pi /4`$ the flat phase distribution corresponding to the $`n=0`$ solution at zero current is displaced towards the value $`\varphi =\pi `$ when applying an external current. The corresponding current distribution seen in Fig. 7a, (dotted line) is straight line and the net current is small for this orientation. For the $`n=1`$ solution at the point where the instability sets in i.e. $`\theta =0.242\pi `$, the current density distribution is symmetric around zero as seen in Fig. 7b (dotted line) and carries zero net current at this point. Thus the instability occurs just before the angle where a full antifluxon enters the junction. A slightly different situation occurs in the magnetic interference pattern of a pure $`s`$-wave superconductor junction  where, the net current is zero at the magnetic field where a full fluxon or antiluxon enters the junction, in the no flux $`0`$-mode. At the point $`\theta =\pi /2`$, of the $`n=1`$-mode the junction contains more than one fluxon causing the characteristic oscillations in the current density around the junction center as seen in Fig. 7b (dashed line). This reduces the critical current for this orientation.

For the $`d_{x^2y^2}+id_{xy}`$ pairing symmetry state, we plot in Fig. 8a) the flux content for the $`n=0,1,1`$, versus the angle $`\theta `$. Note the half integer or multiplies value of $`\mathrm{\Phi }`$ at $`\theta `$ close to $`0`$ or $`\pi /2`$. For this grain orientation the magnetic flux is only sensitive to the real part of the order parameter, which has a sign change but does not break time-reversal symmetry. In the $`d_{x^2y^2}+is`$-wave state the order parameter is complex for all junction orientations and breaks the time-reversal symmetry. Close to $`0`$ or $`\pi /2`$ the flux is fractional. The flux quantization at $`\theta =0`$ can be used to discriminate between these states.

In Fig. 8b) we plot the critical current per unit length evolution with the grain angle $`\theta `$ in the $`d_{x^2y^2}+id_{xy}`$-wave state. Close to $`\theta =0`$ we see that the $`I_c^{ov}`$ for the $`n=0,1`$ solutions, coincide. This happens also at $`\theta =\pi /2`$ for the $`n=0,1`$ solutions. In these orientations the order parameter becomes pure real and does not break the time-reversal symmetry. As a result the critical current at these angles is the same as in a junction with pure $`d`$-wave symmetry. At $`\theta =\pi /4`$ the order parameter is pure imaginary and has the same magnitude for both pairing states. As a consequence for $`\theta =\pi /4`$, the critical currents for both junctions coincide. Also the unstable part of the $`n=1`$ branch, in the $`I_c`$ vs $`\theta `$ is almost the same for the two symmetry states, due to the small difference in the flux, compared with the large flux content of the solutions in this region.

### B Magnitude of the secondary order parameter

In the above calculations the magnitude of the secondary order parameter is small compared to the dominant (i.e. $`n_{20}=0.1n_{10}`$). However the maximum fraction of the secondary component that has been observed in phase coherent experiments employing different materials, geometries, and techniques is up to $`25\%`$ of the dominant . This triggered our interest to study the magnetic flux and also the critical currents as a function of the strength ($`n_s`$) of the secondary order parameter, where the magnitude of the dominant order parameter $`n_d`$ is also varied in a way that $`n_s+n_d=1`$. When $`n_s=0`$ only the $`d_{x^2y^2}`$-wave order parameter is present, while when $`n_s=1`$ only the $`s`$-wave order parameter appears. This situation can be realized for example near the surface where the $`d_{x^2y^2}`$-wave order parameter is suppressed and the $`s`$-wave order parameter is enhanced. The result is presented in Fig. 9a) and 9b) for the $`d_{x^2y^2}+is`$-wave case at $`\theta =0`$. We see that when the secondary component is absent (i.e. $`n_s=0`$) the picture of the $`d_{x^2y^2}`$-wave state is reproduced. The same picture also holds for the $`d_{x^2y^2}+id_{xy}`$-wave state at $`\theta =0`$, since the order parameter for the $`d_{x^2y^2}+id_{xy}`$-wave state at $`\theta =0`$ is real not breaking the time-reversal symmetry. So for $`\theta =0`$ the magnetic flux and the critical current for the $`d_{x^2y^2}+id_{xy}`$-wave state would not show any change with the variation of the secondary order parameter $`d_{xy}`$. As $`n_s`$ is increasing the modes $`n=0`$ and $`n=1`$ are no more degenerate, in the sense that their flux deviates from the value $`\mathrm{\Phi }_0/2`$ and $`\mathrm{\Phi }_0/2`$ respectively and also their critical currents are no longer equal. The mode $`n=0`$ has larger critical current because it has smaller flux content in absolute value. For values of $`n_s`$ close to unity, the different modes contain integer magnetic flux, as in the junction between $`s`$-wave superconductors, and also their critical currents have the same values as in the perfect junction problem. The conclusion is that the larger the secondary component is in a sample the easier is to be detected in a flux measurement experiment.

### C Magnetic field

We now examine the influence of the magnetic field on the spontaneous vortices for broken time reversal symmetry pairing states. In Fig. 10 we plot the magnetic flux at zero current versus the magnetic field $`H`$ for the $`d_{x^2y^2}+id_{xy}`$-wave state at $`\theta =0`$. In the pure $`s`$-wave superconductor junction there is no overlap between different modes in the magnetic flux, and each mode has magnetic flux which is more than $`n\mathrm{\Phi }_0`$ and less than $`(n+1)\mathrm{\Phi }_0`$. In this problem due to spontaneous magnetization the range of the modes is different and in some cases overlapping, and the labeling is with a single index $`n`$, corresponding to the pure $`s`$-wave superconductor junction (n,n+1) mode . Moreover the range in magnetic flux of each mode is displaced compared to the pure $`s`$-wave superconductor junction problem by an amount which corresponds to the intrinsic flux. Also we have the existence of stable vortex states i.e. $`n=0,1`$, together with the unstable ones i.e. $`n=1`$, $`2`$ in a large interval of the magnetic field, which is almost the same. The $`n=2`$ mode extends to zero magnetic field, and the reason we didn’t examined this mode in Sec. IV is because the stability analysis shows negative eigenvalues for all the range of junction orientations, at $`H=0`$. In the long $`s`$-wave junction the extremum of the mode $`(0,1)`$ in $`H`$ is the critical field for one fluxon (antifluxon) penetration from the edges, \[denoted by $`H_{cr}`$ $`(H_{cl})`$, for the right (left) edge\], and is equal to $`2(2)`$. The solution for the phase at these extremum values of the field becomes unstable because the value of the phase at the junction edges reaches a critical value. In the problem of a junction with some spontaneous flux, we consider here, the range of the corresponding mode $`0`$ in $`H`$ is significantly broadened and also the instability at the boundaries sets in due to different reasons. In particular the instability occurs due to the interaction of the flux entering from the junction edges, when the magnetic field reaches the critical value $`H_{cr}(H_{cl})`$, with the spontaneous flux at the center. Similar features are encountered in the problem of flux pinning from a macroscopic defect in a conventional $`s`$-wave junction.

We now examine the magnetic-interference pattern for the two symmetries where the bias current enters in the overlap geometry. In the $`d_{x^2y^2}+id_{xy}`$-wave case, where $`\theta =0`$, this pattern has a symmetric form as we can see from Fig. 11(a). This is because this result is only sensitive to the real part of the order parameter, which has a sign change but does not break time-reversal symmetry. For the angle $`\theta =0.5`$ where the order parameter has a finite imaginary part and breaks the time-reversal symmetry this pattern becomes asymmetric and the ”dip” appears to a value of flux slightly different than zero. Note that the asymmetry refers mainly to the modes $`n=0`$, and $`n=1`$. The other modes are not influenced much due to their higher flux content. Also the critical current is suppressed compared to the case where $`\theta =0`$ as can be seen in Fig. 11(b), due to a drop in $`J_c`$.

In the $`d_{x^2y^2}+is`$-wave symmetry, in the limit where $`\theta 0`$, the order parameter is complex and the pattern is asymmetric as can be seen in Fig. 11c, for the angle $`\theta =0`$. This is in agreement with our previous work for the inline current input for a junction with $`d_{x^2y^2}+is`$ symmetry . There it was found that the pattern is asymmetric for lengths as long as $`L=10\lambda _J`$. For angles close to $`\pi /4`$, the magnetic interference pattern is similar with the $`d_{x^2y^2}+id_{xy}`$-state. This is because the $`\mathrm{sin}(2\theta )`$ dependence of the $`d_{xy}`$ component is almost unity. This is seen in Fig. 11d where we present the variation of the critical current per unit length versus the enclosed flux for $`\theta =0.5`$, and the symmetry state is $`d_{x^2y^2}+is`$.

In the short junction limit $`L<\lambda _J`$ the same argument can be applied without any explicit reference to fractional vortex and antivortex solutions. However as we found in our previous work , both $`n=0`$ and $`n=1`$ (there $`f_{va},f_a`$) exist, with reduced flux content, in this limit as a continuation of the corresponding solutions in the large junction limit. In this case the external applied magnetic field becomes equal to the self field, and the maximum current can be calculated analytically ,

$$\frac{I_m(\mathrm{\Phi })}{I_{m0}}=\left|\frac{\mathrm{sin}(\pi \mathrm{\Phi }/2\mathrm{\Phi }_0)\mathrm{cos}[\pi \mathrm{\Phi }/2\mathrm{\Phi }_0+(\varphi _{c2}\varphi _{c1})/2]}{\pi \mathrm{\Phi }/2\mathrm{\Phi }_0}\right|.$$
(11)
As we see at $`\theta =0`$ for the $`d_{x^2y^2}+id_{xy}`$-wave case, the relation $`\varphi _{c2}\varphi _{c1}=n\pi `$ holds, and the magnetic interference pattern becomes symmetric, while for the $`d_{x^2y^2}+is`$, this difference is a fraction of $`\pi `$ and the pattern is asymmetric. However as we increase the junction length, we expect this symmetric pattern for the $`d`$-wave order parameter to be continued. This symmetry in the large junction limit, is described more effectively by the assumption of the $`n=0,1`$ solutions which give a symmetric magnetic interference pattern as we presented. Also the $`n=1`$ solution extends to values for the magnetic flux, where the $`n=0`$ solution is absent. Eliminating one of them will break the symmetry of the diagram.

## VI Experimental relevance

The symmetric pattern with a minimum at zero applied field observed in corner junction experiments between YBCO and Pb at $`\theta =0`$ has been interpreted as an indication of $`d_{x^2y^2}`$-wave symmetry.  This result refers to short junctions where the junction size is much smaller than the Josephson penetration depth. However as we found here these experimental data are also consistent with an order parameter with $`d_{x^2y^2}+id_{xy}`$ pairing symmetry at $`\theta =0`$.

Also the critical current $`I_c`$ versus the magnetic flux $`\mathrm{\Phi }`$ of a SQUID, consisting of two planar Josephson junctions on the faces of YBCO superconducting crystal, connected by a loop of a second superconductor, for $`\theta =0`$ or $`\theta =\pi /2`$ is found shifted by $`\mathrm{\Phi }=0.5\mathrm{\Phi }_0`$ and has a minimum at $`\mathrm{\Phi }=0`$ (instead of a maximum as in a SQUID involving conventional $`s`$-wave superconductors or the edge SQUID in which both junctions are on the same crystal face) but is still symmetric. This result has been attributed to an order parameter with $`d_{x^2y^2}`$-wave symmetry. However the theoretical analysis done by Beasley et al.  shows that it is also consistent with an order parameter with $`d_{x^2y^2}+id_{xy}`$-pairing symmetry at $`\theta =0`$.

In both cases of SQUID and corner junction the symmetric pattern observed at $`\theta =0`$ rules out the $`d_{x^2y^2}+is`$-wave pairing state where the order parameter is complex everywhere resulting in an asymmetric $`I_c`$ versus $`\mathrm{\Phi }`$ pattern for all angles $`\theta `$. However the small asymmetry (less than $`2\%`$) observed at $`\theta =0`$ in some experiments can be attributed to various complicating factors e.g. fluxon trapping as will be discussed latter in this section.

The experiment proposed here to resolve ambiguity between $`d_{x^2y^2}+id_{xy}`$ and $`d_{x^2y^2}`$ at $`\theta =0`$, is to execute the same experiments using SQUID or corner junction at an angle between sample faces $`\theta `$ between $`0`$ and $`\pi /2`$. Our theory predicts symmetric(asymmetric) pattern for the $`d_{x^2y^2}`$-wave ($`d_{x^2y^2}+id_{xy}`$)-wave pairing state for the corner junction case. This kind of experiment has already been done in the case of SQUID geometry.  The tunneling directions are defined lithographically and patterned by ion milling of a $`c`$-axis oriented film. A YBCO thin film is patterned into a circle with a series of Nb-Au-YBCO edge junctions at orientations spaces every $`7.5^{}`$. The measurement of the $`I_c`$ vs $`\theta `$, which probes mainly the magnitude of the order parameter has an angular anisotropy, indicating an anisotropic order parameter. Also the execution of this experiments is not easy due to the difficulty in cleaning, polishing a crystal at angle $`\theta `$, between $`0`$ and $`\pi /2`$.

Also in an experiment analogous to the corner junction Miller, Ying et. al.  used frustrated thin-film tricrystal samples to probe the pairing symmetry of YBCO. They found a minimum in the $`I_c`$ vs the externally applied flux $`\mathrm{\Phi }_e`$ diagram at $`\mathrm{\Phi }_e=0`$ in the short junction limit and a maximum at $`\mathrm{\Phi }_e=0`$ for a wide junction where the junction length is much larger than the $`\lambda _J`$. However for a wide junction the correct quantity to be compared should be the total flux $`\mathrm{\Phi }`$ which involves contribution both from the externally applied flux and the intrinsic flux. Also in the tricrystal magnetometry experiments on half-flux quantum Josephson vorticies one can only observe spontaneous magnetization of $`\mathrm{\Phi }_0/2`$ in a frustrated geometry only in the large junction length limit .

There is a number of complicating factors in the interpretation of the experiments involving corner junctions that could lead to an asymmetric ($`I_c`$ vs $`\mathrm{\Phi }`$) pattern even for $`\theta =0`$. These are the asymmetry of the junction (meaning that the critical current of the two junction faces are not equal). This will only cause the dip to be shallower and will maintain the symmetry of the $`I_c`$ vs $`\mathrm{\Phi }`$ diagram. Also these experiments are influenced by the sample geometry and the effect of flux trapping i.e. there can be vortices trapped between the planes of the cuprate superconductors that could affect the $`I_c`$ vs $`\mathrm{\Phi }`$ diagram. In the corner junction case, it creates an asymmetry in the flux modulation curves. However these flux trapping effects are not sufficiently large to change the qualitative interpretation of these experiments.

## VII Conclusions

We studied numerically the possible spontaneous vortex states that may exist in a corner junction between a superconductor with time reversal symmetry broken, (i.e. $`d_{x^2y^2}+id_{xy}`$ or $`d_{x^2y^2}+is`$), and an $`s`$-wave superconductor, in the long junction limit. We studied separately three parameters which can be used to modulate the spontaneous flux. These are the magnetic field $`H`$, the interface orientation $`\theta `$, and the magnitude of the subdominant order parameter $`n_s`$. We pointed out the differences between time reversal broken states under these modulation parameters.

We found that in flux modulation experiments involving superconductors with some spontaneous flux the range in magnetic flux of each mode is displaced compared to the case of a pure $`s`$-wave superconductor junction by an amount which corresponds to the intrinsic flux. In particular when the magnetic field $`H`$ is considered as the modulation parameter, the range in $`H`$ of the lower fluxon modes is significantly broadened compared to the $`s`$-wave case, and the instability at the boundary values of the field sets in due to the interaction of the flux entering from the junction edges with the intrinsic flux. In any case, for each value of the parameter which changes the flux, the modes are separated by a single flux quantum.

We also derived some simple arguments to discriminate between the different pairing states that break the time reversal symmetry. For the $`d_{x^2y^2}+id_{xy}`$-wave pairing state, the junction orientation where $`\theta =0`$ i.e. the lobes of the dominant $`d_{x^2y^2}`$-wave order parameter are at right angles for the corner junction, give flux quantization condition $`\mathrm{\Phi }=n\mathrm{\Phi }_0/2`$ as in the $`d_{x^2y^2}`$-wave state, which is different from the corresponding flux quantization for the $`d_{x^2y^2}+is`$-wave pairing state, at $`\theta =0`$, which is $`\mathrm{\Phi }=(n/2+f)\mathrm{\Phi }_0`$, where $`f`$ is a small quantity. These different conditions provide a way experimentally to distinguish between time reversal broken symmetry states. Note that since the magnitude of the secondary order parameter is small compared to the dominant, the detection of time reversal broken states requires a very precise measurement of the spontaneous magnetic flux.

Also we showed that the magnetic interference pattern at $`\theta =0`$ is symmetric (asymmetric) for the $`d_{x^2y^2}+id_{xy}`$ ($`d_{x^2y^2}+is`$), and this also can be used to probe which symmetry the order parameter has, at least where the junctions are formed. We expect our findings, for the magnetic field dependence of the critical current, to hold even in the short junction limit, where the most experiments on corner junctions have been performed .

## VIII Acknowledgments

One of us N.S. is grateful to A.V. Balatsky, J. Betouras for useful discussions that led to this article. Also N.S. would like to acknowledge the ESF/FERLIN programme for partial support to participate to conferences.