text:
# The 𝐵→𝑋_𝑠𝛾𝛾 decay in the standard model and the general two-Higgs-doublet model

## I Introduction

The flavor changing neutral current (FCNC) inclusive transition $`BX_s\gamma `$ has been a subject of great interest during the past few years. The basic theoretical framework is the standard model (SM) at scales $`m_w`$ or $`m_t`$. QCD short distance corrections  are then incorporated via the renormalization group technique to yield the low energy effective Hamiltonian valid at scale $`O(m_b)`$ which is relevant for B decay processes.

As is well known, the rare decay $`BX_s\gamma `$ is theoretically very clear, measured at CLEO  with increasing accuracy and in remarkable agreement with theoretical estimates. Great progress in both the theoretical calculation  and the experimental measurement  for $`BX_s\gamma `$ decay enable one to constrain the new physics models, such as the two-Higgs-doublet model (2HDM) , the minimal supersymmetric standard model (MSSM)  and the Technicolor models . In light of these developments, it is natural to consider other inclusive channels which as a whole may separate out contributions from various operators in the effective Hamiltonian.

The rare flavor changing neutral current inclusive transition $`BX_s\gamma \gamma `$, which is expected to be $`10^3`$ smaller in its branching ratio relative to $`BX_s\gamma `$, has attracted new attention in view of the planned experiments at KEK and SLAC B-factories and existing accelerators, which may test branching ratios as low as $`10^8`$. Just like the decay $`BX_s\gamma `$, the decay $`BX_s\gamma \gamma `$ is relatively clean after some proper precaution to take out effects due to strong resonance, for instance, the effect of $`\eta _c`$ at its peak to the two photon spectrum . The studies about these rare decay modes will then provide further opportunities in testing the whole technology of weak decays, or better yet in providing some clues of new physics.

There have been some theoretical investigations for the process $`BX_s\gamma \gamma `$, which corresponds at the quark level to the transition $`bs\gamma \gamma `$. Theoretical calculations in the SM were firstly done on the basis of pure electroweak theory  and subsequently improved to include the leading order (LO) QCD effects . Nearly all the studies of $`BX_s\gamma \gamma `$ and $`BX_s\gamma `$ decays are based on the free decay of $`b`$ quark and can be justified from the heavy quark effective theory (HQET). According to the argument in Ref., the HQET corrections to these rare decay modes are strongly suppressed by powers of $`(\mathrm{\Lambda }/m_b)^2`$. Furthermore, Choudhury and Yao  have shown that even when$`(\mathrm{\Lambda }/m_c)^2`$ and $`(\mathrm{\Lambda }/m_c)^4`$ corrections are included, the overall HQET corrections to the free quark results of $`BX_s\gamma \gamma `$ are only a few percent and can be safely neglected. In this paper, we will present our results at free quark level.

Recently, we estimated the new physics contributions to the inclusive decay $`bqg^{}q\overline{q^{}}q^{}`$ with $`q\{d,s\}`$ and $`q^{}\{u,d,s\}`$ in the model III: the third-type two-Higgs-doublet model , and found that the new physics contributions can be significant . In this paper, we calculate the new physics corrections to the branching ratio and differential distributions of the rare decay $`BX_s\gamma \gamma `$ induced by the new gluonic and electroweak charged-Higgs penguin diagrams in the model III.

The plan of this paper is as follows: in the next section, we give a brief review of the process $`BX_s\gamma \gamma `$ in the SM at LO approximation, present the relevant formulae and discuss some general characteristics of this process. In Sec. III, we will investigate in detail this rare decay in the model III and present some numerical results. In the last section, concluding remarks are added.

## II $`BX_s\gamma \gamma `$ in the SM

By the use of an extension of Low’s low energy theorem , or alternatively, by applying the equation of motions , the most general effective Hamiltonian that describes radiative $`bs`$ decays with up to three emitted gluons or photons is given by ***Strictly speaking, operators $`(\overline{s}_\alpha \gamma ^\mu Lu_\alpha )(\overline{u}_\beta \gamma _\mu Lb_\beta )`$ and $`(\overline{s}_\alpha \gamma ^\mu Lu_\beta )(\overline{u}_\beta \gamma _\mu Lb_\alpha )`$ should be added. Here we neglect them for their small coefficients.

$`H_{eff}`$ $`=`$ $`{\displaystyle \frac{4G_F}{\sqrt{2}}}V_{tb}V_{ts}^{}{\displaystyle \underset{i=1}{\overset{8}{}}}C_i(\mu )O_i,`$ (1)
where the current-current, penguin, electroweak- and chromo-magnetic dipole operators are of the form

$`O_1`$ $`=`$ $`(\overline{s}_\alpha \gamma ^\mu Lc_\beta )(\overline{c}_\beta \gamma _\mu Lb_\alpha ),`$ (2)
$`O_2`$ $`=`$ $`(\overline{s}_\alpha \gamma ^\mu Lc_\alpha )(\overline{c}_\beta \gamma _\mu Lb_\beta ),`$ (3)
$`O_3`$ $`=`$ $`(\overline{s}_\alpha \gamma ^\mu Lb_\alpha ){\displaystyle \underset{q=u\mathrm{}b}{}}(\overline{q}_\beta \gamma _\mu Lq_\beta ),`$ (4)
$`O_4`$ $`=`$ $`(\overline{s}_\alpha \gamma ^\mu Lb_\beta ){\displaystyle \underset{q=u\mathrm{}b}{}}(\overline{q}_\beta \gamma _\mu Lq_\alpha ),`$ (5)
$`O_5`$ $`=`$ $`(\overline{s}_\alpha \gamma ^\mu Lb_\alpha ){\displaystyle \underset{q=u\mathrm{}b}{}}(\overline{q}_\beta \gamma _\mu Rq_\beta ),`$ (6)
$`O_6`$ $`=`$ $`(\overline{s}_\alpha \gamma ^\mu Lb_\beta ){\displaystyle \underset{q=u\mathrm{}b}{}}(\overline{q}_\beta \gamma _\mu Rq_\alpha ),`$ (7)
$`O_7`$ $`=`$ $`{\displaystyle \frac{e}{16\pi ^2}}\overline{s}_\alpha \sigma ^{\mu \nu }(m_bR+m_sL)b_\alpha F_{\mu \nu },`$ (8)
$`O_8`$ $`=`$ $`{\displaystyle \frac{g_s}{16\pi ^2}}\overline{s}_\alpha \sigma ^{\mu \nu }(m_bR+m_sL)\lambda _{\alpha \beta }^ab_\beta G_{\mu \nu }^a.`$ (9)
where $`\alpha `$ and $`\beta `$ are color indices, $`a=`$1, $`\mathrm{}`$,8 labels $`SU(3)_C`$ generators, and $`L,R`$=$`(1\gamma _5)/2`$, while $`F_{\mu \nu }`$ and $`G_{\mu \nu }^a`$ denote the QED and QCD field strength tensors, respectively.

The Wilson coefficients $`C_i(\mu )`$ in Eq.(1) are process independent and their renormalization is determined only by the basis of operators $`O_i`$. In our calculation of the leading order matrix elements of $`bs\gamma \gamma `$ shown in Fig.17, we find no regulation scheme dependence enters through the new class of penguin diagrams with two external photon. This means that it is suffice to use LO regularization-scheme-independent Wilson coefficients $`C_{1,\mathrm{},6}`$, $`C_{7\gamma }`$ and $`C_{8g}`$ and we need not consider the matrix elements due to the insertion of $`O_5`$ and $`O_6`$ into the one-photon penguin diagrams. In Table I, we present the numerical values of the Wilson coefficients at the scale $`\mu =O(m_b)`$  These Wilson coefficients are numerically well consistent with those given in Ref.. .

At the leading order, the amplitude for the decay $`b(p)s(p^{})+\gamma (k_1)+\gamma (k_2)`$ can be expressed as

$`A`$ $`=`$ $`{\displaystyle \frac{ie^2G_F}{\sqrt{2}\pi ^2}}V_{tb}V_{ts}^{}\overline{u}_s(p^{})[F_2W_2^{\mu \nu }+F_5(m_bW_{5,b}^{\mu \nu }R+m_sW_{5,s}^{\mu \nu }L)`$ (11)
$`+F_7W_7^{\mu \nu }]u_b(p)ϵ_\mu (k_1)ϵ_\nu (k_2),`$
with

$`F_2`$ $`=`$ $`\left[N_cC_1(\mu )+C_2(\mu )\right]Q_u^2\kappa _c`$ (15)
$`+C_3(\mu )\left\{N_c\left[Q_d^2(\kappa _d+\kappa _s+\kappa _b)+Q_u^2(\kappa _u+\kappa _c)\right]+Q_d^2(\kappa _s+\kappa _b)\right\}`$
$`+C_4(\mu )\left\{Q_d^2\left[(\kappa _d+\kappa _s+\kappa _b)+Q_u^2(\kappa _u+\kappa _c)\right]+N_cQ_d^2(\kappa _s+\kappa _b)\right\}`$
$`\left[N_cC_5(\mu )+C_6(\mu )\right]\left[Q_d^2(\kappa _d+\kappa _s+\kappa _b)+Q_u^2(\kappa _u+\kappa _c)\right],`$
$`F_5`$ $`=`$ $`\left[C_5(\mu )+N_cC_6(\mu )\right]Q_d^2,`$ (16)
$`F_7`$ $`=`$ $`C_7(\mu )Q_d,`$ (17)
$`W_2^{\mu \nu }`$ $`=`$ $`i\left\{{\displaystyle \frac{1}{k_1k_2}}\left[k_1^\nu ϵ_{\mu \rho \sigma \lambda }\gamma ^\rho k_1^\sigma k_2^\lambda k_2^\mu ϵ_{\nu \rho \sigma \lambda }\gamma ^\rho k_1^\sigma k_2^\lambda \right]+ϵ_{\mu \nu \rho \lambda }\gamma ^\rho (k_2k_1)^\lambda \right\}L,`$ (18)
$`W_{5,q}^{\mu \nu }`$ $`=`$ $`{\displaystyle \frac{1}{m_q^2}}\left(iϵ_{\mu \nu \lambda \sigma }k_1^\lambda k_2^\sigma \gamma _5+k_1k_2g^{\mu \nu }k_1^\nu k_2^\mu \right)\left(12\kappa _q\right)+4\left(g^{\mu \nu }{\displaystyle \frac{k_1^\nu k_2^\mu }{k_1k_2}}\right)\kappa _q,`$ (19)
$`W_7^{\mu \nu }`$ $`=`$ $`{\displaystyle \frac{1}{2}}[{\displaystyle \frac{1}{2pk_2}}\overline{)}k_1\gamma ^\mu (m_bR+m_sL)(\overline{)}p\overline{)}k_2+m_b)\gamma ^\nu `$ (21)
$`+{\displaystyle \frac{1}{2p^{}k_2}}\gamma ^\nu (\overline{)}p\overline{)}k_1+m_s)\overline{)}k_1\gamma ^\mu (m_bR+m_sL)],`$
where $`k_1`$ and $`k_2`$ are four-momenta of the emitted photons, $`N_c=3`$ denotes the number of colors, $`Q_u=2/3`$ and $`Q_d=1/3`$ are the up-type and down-type quark electric charges, and the factor $`\kappa _q`$ is defined as

$`\kappa _q`$ $`=`$ $`{\displaystyle \frac{1}{2}}+{\displaystyle \frac{1}{z_q}}{\displaystyle _0^1}{\displaystyle \frac{dx}{x}}\mathrm{ln}\left(1z_qx+z_q^2x^2\right)`$ (22)
$`=`$ $`\{\begin{array}{cc}\frac{1}{2}\frac{2}{z_q}\left(\mathrm{arctan}\sqrt{\frac{z_q}{4z_q}}\right)^2\hfill & \mathrm{if}z_q<4,\hfill \\ \frac{1}{2}+\frac{1}{z_q}\left[\frac{\pi ^2}{2}+2\left(\mathrm{ln}\frac{\sqrt{z_q}+\sqrt{z_q4}}{2}\right)^22i\pi \mathrm{ln}\left(\frac{\sqrt{z_q}+\sqrt{z_q4}}{2}\right)\right]\hfill & \mathrm{otherwise}.\hfill \end{array}`$ (25)
Here $`z_q=2k_1k_2/m_q^2`$. The masses appearing in $`W_7^{\mu \nu }`$ arise from operator $`O_7`$ while those in $`W_{5,q}^{\mu \nu }`$ are internal quark masses in the loops (see Fig.1). In our numerical analysis, we take all these quark masses as current mass. The factor $`\kappa _q`$ in Eq.(25) is the loop integral function  and its absorption part reflects intermediate quark threshold effect. In getting the expressions of $`W_i^{\mu \nu }`$ that are essential to obtain the matrix elements of $`bs\gamma \gamma `$ decay presented in , we used the following identities:

$`\gamma _\mu \gamma _\lambda \gamma _\sigma `$ $`=`$ $`g_{\sigma \lambda }\gamma _\mu +g_{\mu \lambda }\gamma _\sigma g_{\mu \sigma }\gamma _\lambda +iϵ_{\mu \rho \sigma \lambda }\gamma ^\rho \gamma _5,`$ (26)
$`\gamma _{\mu _1}\gamma _{\mu _2}\gamma _{\mu _3}\gamma _{\mu _4}`$ $`=`$ $`iϵ_{\mu _1\mu _2\mu _3\mu _4}\gamma _5g_{\mu _3\mu _4}g_{\mu _1\mu _2}+g_{\mu _2\mu _4}g_{\mu _1\mu _3}g_{\mu _2\mu _3}g_{\mu _1\mu _4}+g_{\mu _3\mu _4}\gamma _{\mu _1}\gamma _{\mu _2}`$ (28)
$`g_{\mu _2\mu _4}\gamma _{\mu _1}\gamma _{\mu _3}+g_{\mu _2\mu _3}\gamma _{\mu _1}\gamma _{\mu _4}+g_{\mu _1\mu _4}\gamma _{\mu _2}\gamma _{\mu _3}g_{\mu _1\mu _3}\gamma _{\mu _2}\gamma _{\mu _4}+g_{\mu _1\mu _2}\gamma _{\mu _3}\gamma _{\mu _4}.`$
The square amplitude summed over spins and polarizations is then given by

$`|A|^2`$ $`=`$ $`{\displaystyle \frac{1}{4}}\left({\displaystyle \frac{e^2G_F}{\sqrt{2}\pi ^2}}V_{tb}V_{ts}^{}\right)^2m_b^4\{|F_2|^2A_{22}+|F_5|^2A_{55}+|F_7|^2A_{77}+2Re[F_7F_2^{}]A_{27}`$ (30)
$`+2Re[F_5F_2^{}(12\kappa _s)]A_{25}^s+2Re[F_7F_5^{}]A_{57}+2Re[F_5F_2^{}(12\kappa _b)]A_{25}^b\},`$
where the quantities $`A_{ij}`$ denote the contractions between the tensor $`W_i^{\mu \nu }`$ and $`W_j^{\mu \nu }`$, and the first two terms and the third term of $`|A|^2`$ arise from one-particle irreducible (1PI) diagrams and one-particle reducible (1PR) diagrams, respectively. In order to give the explicit expressions of $`A_{ij}`$, we introduced the following notations,

$`s`$ $`=`$ $`{\displaystyle \frac{2k_1k_2}{m_b^2}},t={\displaystyle \frac{2pk_2}{m_b^2}},u={\displaystyle \frac{2pk_1}{m_b^2}},rho_1={\displaystyle \frac{m_s^2}{m_b^2}},`$ (31)
$`\rho `$ $`=`$ $`1ut+s={\displaystyle \frac{p^2}{m_b^2}},\overline{s}={\displaystyle \frac{s}{1\rho }},\overline{t}={\displaystyle \frac{t}{1\rho }},\overline{u}={\displaystyle \frac{u}{1\rho }},`$ (32)
where $`m_s`$ is $`s`$ quark mass appearing in $`W_{5,7}^{\mu \nu }`$ and $`p^{}`$ is the momentum of the outgoing $`s`$ quark. In this framework, coefficients $`A_{ij}`$ are then given by

$`A_{22}`$ $`=`$ $`2[(1\rho )^2(1+\rho )s],`$ (33)
$`A_{55}`$ $`=`$ $`\left\{16|\kappa _b|^2+|(12\kappa _b)s+4\kappa _b|^2+\rho _1[16|\kappa _s|^2+|(12\kappa _s){\displaystyle \frac{s}{\rho _1}}+4\kappa _s|^2]\right\}(1s+\rho )`$ (35)
$`+16Re\left\{8\rho _1\kappa _b\kappa _s^{}+s[\kappa _b2(1+\rho _1)\kappa _b\kappa _s^{}+\rho _1\kappa _s^{}]\right\},`$
$`A_{25}^{b,s}`$ $`=`$ $`\pm s(1\rho s),`$ (36)
$`A_{27}`$ $`=`$ $`\left[2(1+\rho _1)s+{\displaystyle \frac{(\rho +\rho _1)s^2+2(\rho _1\rho )ts}{t(st)}}+{\displaystyle \frac{(\rho +\rho _1)s^2+2(\rho _1\rho )us}{u(su)}}\right],`$ (37)
$`A_{57}`$ $`=`$ $`Re\{8(\kappa _b+\rho _1\kappa _s)s[4\rho _1(\kappa _b+\kappa _s)+s[(12\kappa _s)+\rho _1(12\kappa _b)]]`$ (39)
$`\times [{\displaystyle \frac{s^2}{t(st)}}+{\displaystyle \frac{s^2}{u(su)}}]4\left[(\kappa _b+\rho _1\kappa _s)(\rho \rho _1)({\displaystyle \frac{s^2}{t(st)}}+{\displaystyle \frac{s^2}{u(su)}})\right]\},`$
$`A_{77}`$ $`=`$ $`(1+\rho _1)\left[(1\rho )A_{77}^1+A_{77}^2\right]+A_{77}^3,`$ (40)
$`A_{77}^1`$ $`=`$ $`{\displaystyle \frac{1}{\overline{t}}}[1+\overline{u}+{\displaystyle \frac{2\overline{u}(\overline{u}2)}{(1\overline{u})}}\overline{t}+{\displaystyle \frac{2\overline{u}1}{1\overline{u}}}\overline{t}^2]+(\overline{t}\overline{u}),`$ (41)
$`A_{77}^2`$ $`=`$ $`{\displaystyle \frac{2}{\overline{t}^2}}\left[1{\displaystyle \frac{1+\rho _1}{1\overline{u}}}\overline{t}+{\displaystyle \frac{\rho _1}{(1\overline{u})^2}}\overline{t}^2\right]`$ (43)
$`+(\rho _1\rho )({\displaystyle \frac{2\overline{s}}{\overline{t}(1\overline{u})}}{\displaystyle \frac{\overline{s}\overline{u}}{(1\overline{u})^2}}{\displaystyle \frac{\overline{u}^2}{(1\overline{u})^2}})+(\overline{t}\overline{u}),`$
$`A_{77}^3`$ $`=`$ $`2{\displaystyle \frac{s}{\overline{t}\overline{u}}}t\left\{(1+\rho _1)(2+\overline{u}\overline{t})+{\displaystyle \frac{\rho _1}{1\rho }}\left[1{\displaystyle \frac{2(1+\rho _1)\overline{t}\overline{u}}{(1\overline{t})(1\overline{u})}}\right]\overline{s}\right\}`$ (45)
$`+{\displaystyle \frac{2s}{tu}}(\rho _1\rho ){\displaystyle \frac{(1+\rho _1)(2s+\overline{s})\rho _1(\overline{s}^22\overline{t}\overline{u})+2(2\rho )(1+\rho _1)\overline{t}\overline{u}}{(1\overline{u})(1\overline{t})}},`$
Unlike Ref., we here distinguish different $`s`$ quark contributions which enter in the amplitude via $`\rho `$ and $`\rho _1`$ respectively. This is crucial to our results, as described below. Note that in the $`A_{77}`$ term of the square amplitude, there exists infrared (IR) divergence when one integrates the amplitude over the physical phase space. This divergence can be canceled out when one considers $`O(\alpha _e)`$ virtual corrections to the $`bs\gamma `$ amplitude . In order to calculate the physical rate of interest that is free of divergence, we have to impose a cut on the energy of each photon, which will naturally correspond to the experimental cut imposed on the minimum energy of detectable photons.

In numerical calculations, we use the input parameters listed in Table II and take the following cuts:

$`E_\gamma >100MeV,E_s>600MeV,\theta >{\displaystyle \frac{\pi }{9}},`$ (46)
where $`E_\gamma `$ and $`E_s`$ denote the energy of photon and that of of the outgoing mesons respectively, and $`\theta `$ is the angle between any two outgoing particles. The first constraint is required to avoid $`IR`$ divergence while the last constraint is to exclude photons that are emitted too close to each other or to the outgoing s quark. For the mass of $`s`$ quark, cares must be taken. In principle, constituent mass should be used in the phase space integration, while the masses appeared in $`W_7^{\mu \nu }`$ and $`W_5^{\mu \nu }`$ should be the current mass. This means that $`\rho `$ and $`\rho _1`$ appeared in Eqs.(33-45) should be

$`\rho ={\displaystyle \frac{m_s^2(constituent)}{m_b^2}}{\displaystyle \frac{m_K^2}{m_b^2}},\rho _1={\displaystyle \frac{m_s^2(current)}{m_b^2}}.`$ (47)
Following Ref., the branching ratio of $`BX_s\gamma \gamma `$ decay can be written as

$`(BX_s\gamma \gamma )\left[{\displaystyle \frac{\mathrm{\Gamma }(bs\gamma \gamma )}{\mathrm{\Gamma }(bcl\nu _l)}}\right]^{th}\times (BX_cl\nu _l)^{expt}.`$ (48)
Using the input parameters presented in Table II, setting $`\mu =m_b`$ and taking the cuts given in Eq.(46), we find numerically that

$`(BX_s\gamma \gamma )4.6\times 10^7,A_{FB}(BX_s\gamma \gamma )=0.79,`$ (49)
where $`A_{FB}`$ is the forward-backward asymmetry of the $`BX_s\gamma \gamma `$ decay.

If we use a common $`s`$ quark mass $`m_s=0.5`$ GeV in the numerical calculations as Ref. did, we find that

$`(BX_s\gamma \gamma )=3.8\times 10^7,A_{FB}(BX_s\gamma \gamma )=0.76,`$ (50)
which agree well with the results presented in Ref., where $`(BX_s\gamma \gamma )=3.7\times 10^7`$ and $`A_{FB}=0.78`$. For $`m_s=0.15`$ GeV, however, we find that

$`(BX_s\gamma \gamma )=5.1\times 10^7,A_{FB}(BX_s\gamma \gamma )=0.81.`$ (51)
From above numerical results, we find that the $`m_s`$-dependence is weak for $`A_{FB}`$, but relatively strong for the branching ratio.

Among all the contributions to the decay rate, the contribution from 1PR diagrams is predominant, larger than $`97\%`$ of the total, which is due to the cancellation between $`C_1`$ and $`C_2`$ and the QCD enhancement of $`|C_7|`$. Numerical results also show that the branching ratio is sensitive to the cuts we have imposed. For example, if we demand $`\theta `$ in Eq.(46) larger than $`\pi /6`$, the branching ratio will reduce to $`3.6\times 10^7`$. Moreover, because of the scale dependence of Wilson coefficients, there exists a $`25\%`$ theoretical uncertainty at leading order approximation, as is the case for $`BX_s\gamma `$.

In Fig.17, we present the normalized differential distribution $`(1/\mathrm{\Gamma })d\mathrm{\Gamma }/ds`$ versus $`s`$ ( solid curve ). Comparing with the corresponding result given in , we find that the peak due to $`c\overline{c}`$ threshold effect is smeared and the average invariant mass of the two photons is lowered when QCD corrections are added. In Figs.(17-17) , the normalized differential distribution $`(1/\mathrm{\Gamma })d\mathrm{\Gamma }/d\mathrm{cos}\theta _{\gamma \gamma }`$ versus $`\mathrm{cos}\theta _{\gamma \gamma }`$ (solid curve in Fig.17) and the spectrum of the two photons, defined as the photon with lower energy (solid curve in Fig.17) and the photon with higher energy (solid curve in Fig.17), are plotted.

From these four figures, one can understand the kinematics of the process as follows: an energetic $`s`$ quark with mean energy around $`2.0GeV`$ tends to be emitted, compensated by the harder of the two photons; while the less energetic photon tends to go in the direction of the $`s`$ quark. This topology is typical of a bremsstrahlung event of the $`s`$ quark.  Such kinematics are very useful for us to separate the short-distance (SD) signal from the long-distance (LD) background, which may come from the channel $`BX_s\eta (\eta ^{})X_s\gamma \gamma `$ . For example, by demanding $`s0.3`$ and $`\mathrm{cos}\theta _{\gamma \gamma }0`$, about $`75\%`$ of the signal is remained while almost all the LD background are removed.

In Figs.(17-17), the contribution from 1PI diagrams (dashed curve) and that from the interference between 1PI and 1PR (dot-dashed curve) are also plotted. From these figures, we can see that except in some marginal area, the contributions from 1PI diagrams and the interference are much smaller than those coming from 1PR diagrams. Taking this in mind, one can infer that though the four normalized differential distributions are calculated at leading order approximation, they are insensitive to the variation of scale. The reason is that, to a good approximation, a common factor $`|C_7(\mu )|^2`$ can be extracted from both $`\mathrm{\Gamma }`$ and $`d\mathrm{\Gamma }`$ and will disappear in the ratio of them. This feature becomes more evident when the cuts of $`s0.3`$ and $`\mathrm{cos}\theta _{\gamma \gamma }0`$ are imposed, which removes nearly all the contributions of $`O_1O_6`$ (see Fig.17). In fact, a more general statement is that: as long as the new physics effects beyond the SM appear only in the matching of the Wilson coefficients of the standard effective operator basis, which is so in the MSSM , and $`|C_7|`$ is not suppressed as required by $`bs\gamma `$, these four observable distributions are insensitive to the variation of scale $`\mu `$ and the new physics.

## III $`BX_s\gamma \gamma `$ In Model III

### A Model description

The two-Higgs-doublet model (2HDM)  is the simplest extension of the SM . During past years, the models I and II have been studied extensively in literature and tested experimentally and the model II has been very popular since it is the building block of the MSSM. In Ref., the authors studied the $`BX_s\gamma \gamma `$ decay in the SM and models I and II in the basis of pure electroweak theory , and found that the branching ratios in models I and II can be appreciably different from that in the SM. In this paper, we focus on estimating the new physics effects on the $`BX_s\gamma \gamma `$ decay in the framework of the third type of 2HDM, usually known as the model III . In the model III, no discrete symmetry is imposed and both up- and down-type quarks may couple with either of the two Higgs doublets. The Yukawa coupling in quark sector in this case is

$`_Y=\eta _{ij}^U\overline{Q}_{iL}\stackrel{~}{\varphi _1}U_{jR}+\eta _{ij}^D\overline{Q}_{iL}\stackrel{~}{\varphi _1}D_{jR}+\xi _{ij}^U\overline{Q}_{iL}\stackrel{~}{\varphi _2}U_{jR}+\xi _{ij}^D\overline{Q}_{iL}\stackrel{~}{\varphi _2}D_{jR}+h.c,`$ (52)
where $`\varphi _1`$ and $`\varphi _2`$ are the two Higgs doublets and $`\eta _{ij}^{U,D}`$ and $`\xi _{ij}^{U,D}`$ are the Yukawa couplings. As described in , in order to let $`\varphi _1`$ correspond to generate fermion masses while $`\varphi _2`$ to introduce new interactions, one can choose the following parameterization of the Higgs doublets

$`\varphi _1={\displaystyle \frac{1}{\sqrt{2}}}\left[\left(\begin{array}{c}0\\ v+H^0\end{array}\right)+\left(\begin{array}{c}\sqrt{2}\chi ^+\\ i\chi ^0\end{array}\right)\right],\varphi _2={\displaystyle \frac{1}{\sqrt{2}}}\left(\begin{array}{c}\sqrt{2}H^+\\ H^1+iH^2\end{array}\right),`$ (59)
where $`v=(\sqrt{2}G_F)^{1/2}=246GeV`$. After the rotation that diagonalizes the mass matrices of quark fields and that of the Higgs doublets, we can get the Lagrangian that is relevant for our following discussions. The main features of model III are as follows :

* FCNC may exist at the tree level. The neutral and the charged flavor changing couplings are related by:

$`\xi _{charge}^U=\xi _{neutral}^UV_{CKM},\xi _{charge}^D=V_{CKM}\xi _{neutral}^D,`$ (60)
where $`V_{CKM}`$ is the ordinary mixing matrix between down-type quarks, and

$`\xi _{neutral}^{U,D}=(V_L^{U,D})^1\xi ^{U,D}V_R^{U,D}.`$ (61)
* Like the models I and II, there are also five Higgs bosons in model III: the light and heavy CP-even neutral Higgs boson $`h^0`$ and $`\overline{H}^0`$, one CP-odd neutral Higgs boson $`A^0`$ and a pair of charged Higgs bosons $`H^\pm `$. The transformation relation between $`(H^0,H^1,H^2)`$ in Eq.(59) and the mass eigenstates $`(\overline{H}^0,h^0,A^0)`$ can be found in .

In the following of this paper, we will parameterize $`\xi _{neutral}^{U,D}`$ as Ref.

$`(\xi _{neutral}^{U,D})_{ij}=\lambda _{ij}^{U,D}{\displaystyle \frac{\sqrt{m_im_j}}{v}},`$ (62)
and treat $`\lambda _{ij}^{U,D}`$ as basic free parameters.

### B Experimental constraints

There is a considerable interest in the constraint of the parameter space of the 2HDM, especially in model III, since the FCNC may appear at the tree level. Compared with the SM, the additional free parameters of model III are the masses of the additional Higgs bosons and the coupling constants $`\lambda _{ij}`$. In this subsection, we summarize the main constraints on these parameters from direct searches at LEP experiments, and from the measurements of $`F^0\overline{F}^0`$ mixing with $`F^0=K^0,B_d^0`$.

Let’s firstly turn to mass constraints. The LEP Working Group for Higgs boson searches  has recently reported excess of events that might be due to the production of a neutral Higgs boson weighing about $`115`$ GeV, and placed a lower mass limit for the SM Higgs boson $`M_{H^0}>113.5`$ GeV at $`95\%C.L.`$. For the $`h^0`$ and $`A^0`$ Higgs bosons of the MSSM, the new $`95\%C.L.`$ limit is $`M_{h^0}>88.3`$ GeV and $`M_{A^0}>88.4`$ GeV respectively. But one should note that these constraints on the masses of $`h^0`$ and $`A^0`$ are not applicable in the 2HDM’s because the coupling of $`h^0`$ to $`Z^0Z^0`$ and to $`A^0Z^0`$ go like $`\mathrm{sin}\alpha `$ and $`\mathrm{cos}\alpha `$( $`\alpha `$ is a free parameter of the model), respectively. For a very small mixing angle $`\alpha `$ and with $`A^0`$ ($`h^0`$) sufficiently heavy no lower limit on $`M_{h^0}`$ ($`M_{A^0}`$) can be set from the LEP data. For the charged Higgs bosons in the 2HDM’s, the $`95\%C.L.`$ limit is $`M_{H^\pm }>78.6`$ GeV.

Indirect constraints on Higgs masses come from $`R_b`$ and $`\rho `$, which have been measured at LEP. According to studies in Ref., if one requires $`R_b^{theory}`$ to be within $`1\sigma `$ deviation from $`R_b^{exp}`$ and $`\rho `$ within $`2\sigma `$ deviation from $`\rho ^{exp}`$, the preferred range of Higgs boson masses in model III are $`80GeV<M_{h^0}M_{A^0}<120GeV`$ and $`M_{H^\pm }(180220)GeV`$. On the other hand, if one allows $`R_b^{theory}`$ to vary within $`2\sigma `$ errors of $`R_b^{exp}`$, a charged Higgs mass of a few hundreds GeV is feasible. In addition, unitary condition requires these masses to be less than $`1TeV`$.

As for the couplings $`\lambda _{ij}^{U,D}`$, a lot of processes have to be analyzed to give reasonable constraints, as have been done in Refs.. In Table III, we list current constraints on these couplings along with the processes from which constraints have been or will be placed. It should be noted that despite current constraints on flavor changing couplings $`\lambda _{bs}^D`$ and $`\lambda _{tc}^U`$ are rather loose, some much more stringent constraints will be imposed in the future experiments.

For simplicity, we set (a) all the FC couplings $`\lambda _{ij}^{U,D}(ij)=0`$; (b) $`\lambda _{ii}^U=\lambda _1`$ for $`i=u,c,t,`$ and $`\lambda _{ii}^D=\lambda _2`$ for $`i=d,s,b`$; and (c) $`\lambda _{1,2}`$ are real numbers. The advantage of such a setting is the new contributions to Wilson coefficients discussed below come only from the charged Higgs penguin diagrams and the neutral Higgs bosons are irrelevant to our discussion. The model III therefore differs from the model II only in the couplings of the charged Higgs boson to fermions. The Feynman rules for the $`\overline{U_i}H^+D_j`$ and $`\gamma H^+H^{}`$ vertex are

$`\overline{U_i}H^+D_j:`$ $`{\displaystyle \frac{i}{2}}\left[V_{CKM}\xi _{ij}^D(1+\gamma _5)\xi _{ij}^UV_{CKM}(1\gamma _5)\right],`$ (63)
$`\gamma H^+H^{}:`$ $`ie(p_2p_1)^\mu ,`$ (64)
For more details about the couplings between the Higgs bosons and quarks or gauge bosons, one can see Ref..

Another very important source of constraint is the process $`BX_s\gamma `$. In the SM, the theoretical prediction of $`(BX_s\gamma )`$ is $`(2.8\pm 0.8)\times 10^4`$ at LO approximation , and $`(3.28\pm 0.33)\times 10^4`$ at NLO approximation . The $`95\%C.L.`$ limit from CLEO measurement  is $`2\times 10^4<(BX_s\gamma )<4.5\times 10^4`$. One can see that the theory and the data are in remarkable agreement, which implies that only a narrow parameter space of new physics can survive. Since the aim of this paper is to estimate the new physics effect on the process $`BX_s\gamma \gamma `$, we will vary the input parameters of 2HDM in a rather large range allowed by

$`{\displaystyle \frac{2.0\times 10^4}{2.8\times 10^4}}<{\displaystyle \frac{_{III}^{LO}(BX_s\gamma )}{_{SM}^{LO}(BX_s\gamma )}}|{}_{\mu =m_b}{}^{}<{\displaystyle \frac{4.5\times 10^4}{2.8\times 10^4}}.`$ (65)
### C Calculation of $`BX_s\gamma \gamma `$

In model III, the operator bases presented in Eq.(9) are insufficient for the process $`BX_s\gamma \gamma `$ and we need following additional operators when all the FC couplings are zero : Here we neglect the operators $`O_{9,10}`$ and $`O_{9,10}^{}`$ introduced in  because in our case, the corresponding Wilson coefficients keep zero at any scale.

$`O_1^{}`$ $`=`$ $`(\overline{s}_\alpha \gamma ^\mu Rc_\beta )(\overline{c}_\beta \gamma _\mu Rb_\alpha ),`$ (66)
$`O_2^{}`$ $`=`$ $`(\overline{s}_\alpha \gamma ^\mu Rc_\alpha )(\overline{c}_\beta \gamma _\mu Rb_\beta ),`$ (67)
$`O_{3,5}^{}`$ $`=`$ $`(\overline{s}_\alpha \gamma ^\mu Rb_\alpha ){\displaystyle \underset{q=u\mathrm{}b}{}}(\overline{q}_\beta \gamma _\mu (R,L)q_\beta ),`$ (68)
$`O_{4,6}^{}`$ $`=`$ $`(\overline{s}_\alpha \gamma ^\mu Rb_\beta ){\displaystyle \underset{q=u\mathrm{}b}{}}(\overline{q}_\beta \gamma _\mu (R,L)q_\alpha ),`$ (69)
$`O_7^{}`$ $`=`$ $`{\displaystyle \frac{e}{16\pi ^2}}\overline{s}_\alpha \sigma ^{\mu \nu }(m_bL+m_sR)b_\alpha F_{\mu \nu },`$ (70)
$`O_8^{}`$ $`=`$ $`{\displaystyle \frac{g_s}{16\pi ^2}}\overline{s}_\alpha \sigma ^{\mu \nu }(m_bL+m_sR)\lambda _{\alpha \beta }^ab_\beta G_{\mu \nu }^a.`$ (71)
For the evaluation of Wilson coefficients, we need their initial values with standard matching computations. Denoting the Wilson coefficients in the SM with $`C_i^{SM}(m_W)`$ and those from the additional charged Higgs contribution (see Fig.17) with $`C_i^H(m_W)`$, we have the initial values of the Wilson coefficients for the first set of operators \[ Eq. (9) \]

$`C_i^{III}(m_W)`$ $`=`$ $`C_i^{SM}(m_W)+C_i^H(m_W),`$ (72)
where

$`C_{1\mathrm{}6}^H(m_W)`$ $`=`$ $`0,`$ (73)
$`C_7^H(m_W)`$ $`=`$ $`{\displaystyle \frac{1}{2}}\left[\lambda _1^2f_1(y)+\lambda _1\lambda _2f_2(y)+{\displaystyle \frac{m_s^2}{m_t^2}}\lambda _2^2f_1(y)\right],`$ (74)
$`C_8^H(m_W)`$ $`=`$ $`{\displaystyle \frac{1}{2}}\left[\lambda _1^2g_1(y)+\lambda _1\lambda _2g_2(y)+{\displaystyle \frac{m_s^2}{m_t^2}}\lambda _2^2g_1(y)\right],`$ (75)
and the explicit expressions of $`C_i^{SM}(m_W)`$ can be found in the literature . For the primed Wilson coefficients we have

$`C_i^{III}(m_W)=C_i^H(m_W),`$ (76)
and $`C_i^H`$ are given as

$`C_{1\mathrm{}6}^H(m_W)`$ $`=`$ $`0,`$ (77)
$`C_7^H(m_W)`$ $`=`$ $`{\displaystyle \frac{1}{2}}\left[{\displaystyle \frac{m_sm_b}{m_t^2}}\lambda _2^2f_1(y)+{\displaystyle \frac{m_s}{m_b}}\lambda _1\lambda _2f_2(y)+{\displaystyle \frac{m_s}{m_b}}\lambda _1^2f_1(y)\right],`$ (78)
$`C_8^H(m_W)`$ $`=`$ $`{\displaystyle \frac{1}{2}}\left[{\displaystyle \frac{m_sm_b}{m_t^2}}\lambda _2^2g_1(y)+{\displaystyle \frac{m_s}{m_b}}\lambda _1\lambda _2f_2(y)+{\displaystyle \frac{m_s}{m_b}}\lambda _1^2g_1(y)\right].`$ (79)
In our calculations we neglect the contributions of the internal $`u`$ and $`c`$ quarks and define functions $`f_1(y)`$, $`f_2(y)`$, $`g_1(y)`$, and $`g_2(y)`$ as

$`f_1(y)`$ $`=`$ $`{\displaystyle \frac{y(75y8y^2)}{72(y1)^3}}+{\displaystyle \frac{y^2(3y2)}{12(y1)^4}}\mathrm{ln}[y],`$ (80)
$`f_2(y)`$ $`=`$ $`{\displaystyle \frac{y(5y3)}{12(y1)^2}}+{\displaystyle \frac{y(3y+2)}{6(y1)^3}}\mathrm{ln}[y],`$ (81)
$`g_1(y)`$ $`=`$ $`{\displaystyle \frac{y(y^2+5y+2)}{24(y1)^3}}+{\displaystyle \frac{y^2}{4(y1)^4}}\mathrm{ln}[y],`$ (82)
$`g_2(y)`$ $`=`$ $`{\displaystyle \frac{y(y3)}{4(y1)^2}}+{\displaystyle \frac{y)}{2(y1)^3}}\mathrm{ln}[y],`$ (83)
where $`y=m_t^2/m_{H^\pm }^2`$. Numerical analysis show that $`f_1(y)`$ and $`f_2(y)`$ have the following properties:

$`f_1(y)<0,f_2(y)>0,f_2(y)>3f_1(y),`$ (84)
and $`g_i(y)`$ have similar properties.

Using the renormalization group, we can get the LO Wilson coefficients at any scale. The Wilson coefficients $`C_i^{III}(\mu )`$ can be obtained from  with $`C_i^{SM}(m_W)`$ replaced by $`C_i^{III}(m_W)`$, while $`C_i^{III}(\mu )`$ are given by§§§Because the strong interaction preserve chirality, the first set of operators (Eq.(9)) cannot mix with its chirality flipped counterparts, the primed operators, and the anomalous dimension matrices of the two separate set of operators are the same and do not overlap . This means the evaluation of $`C_i^{}(\mu )`$ are the same as in the SM.

$`C_{1\mathrm{}6}^{III}(\mu )`$ $`=`$ $`0,`$ (85)
$`C_7^{III}(\mu )`$ $`=`$ $`\eta ^{\frac{16}{23}}C_7^{III}(m_W)+{\displaystyle \frac{8}{3}}(\eta ^{\frac{14}{23}}\eta ^{\frac{16}{23}})C_8^{III}(m_W),`$ (86)
where $`\eta =\alpha _s(m_W)/\alpha _s(\mu )`$.

The square amplitude summed over spin and polarizations can be written as

$`|A|_{III}^2=|A|^2+|A^{}|^2+2|A^{}A^{}|,`$ (87)
where the first and second term comes from $`O_1O_8`$ and $`O_1^{}O_8^{}`$ respectively and the third term is the interference between them. The explicit expression of $`|A|^2`$ is the same as that in Eq.(30) with $`C_i`$ being replaced by $`C_i^{III}`$, and $`|A^{}|^2`$ and $`|A^{}A^{}|`$ are given by

$`|A^{}|^2`$ $`=`$ $`{\displaystyle \frac{1}{4}}\left({\displaystyle \frac{e^2G_F}{\sqrt{2}\pi ^2}}V_{tb}V_{ts}^{}\right)^2m_b^4|F_7^{}|^2A_{77},`$ (88)
$`|A^{}A^{}|`$ $`=`$ $`{\displaystyle \frac{1}{4}}\left({\displaystyle \frac{e^2G_F}{\sqrt{2}\pi ^2}}V_{tb}V_{ts}^{}\right)^2m_b^4{\displaystyle \frac{m_s}{m_b}}\left\{|F_2^{}F_7^{}|A_{27}^{}+|F_5^{}F_7^{}|A_{57}^{}+|F_7^{}F_7^{}|A_{77}^{}\right\},`$ (89)
where

$`F_7^{}`$ $`=`$ $`C_7^{III}Q_d,`$ (90)
$`A_{27}^{}`$ $`=`$ $`2s{\displaystyle \frac{(1+\rho _1)s^2}{t(st)}}+s(\rho _1\rho )({\displaystyle \frac{1}{t}}{\displaystyle \frac{1}{st}})+(tu),`$ (91)
$`A_{57}^{}`$ $`=`$ $`Re[8(\kappa _b+\kappa _s)s[4(\kappa _b+\rho _1\kappa _s)+2s(1\kappa _b\kappa _s)]{\displaystyle \frac{s^2}{t(st)}}`$ (93)
$`+4(\kappa _b+\kappa _s)(\rho _1\rho ){\displaystyle \frac{u}{st}}]+(ut),`$
$`A_{77}^{}`$ $`=`$ $`2[(1\rho )A_{77}^1+A_{77}^2]2s(2+{\displaystyle \frac{2}{t}}+{\displaystyle \frac{2}{u}}+{\displaystyle \frac{2(1\rho _1s)}{tu}}`$ (95)
$`{\displaystyle \frac{(\rho _1\rho )(3+5\rho _1)(1+\rho _1)s}{(st)(su)}}{\displaystyle \frac{(1\rho _1)(s^2+4\rho _1tu)}{(st)(su)tu}}).`$
In Eq.(89-95), $`m_s`$ is the current mass of $`s`$ quark and the explicit expressions of $`A_{77},A_{77}^1`$ and $`A_{77}^2`$ have been given in Eqs.(40-43). For $`|A|_{III}^2`$, two characters should be noted. Firstly, the interference term between $`O_i`$ and $`O_i^{}`$ is suppressed by $`m_s/m_b`$ and therefore, to a good approximation, can be neglected. Secondly, the total 1PR contributions are proportional to $`(|F_7|^2+|F_7^{}|^2)A_{77}`$. Since $`|F_7|^2+|F_7^{}|^2`$ is not suppressed, as required by $`BX_s\gamma `$ in 2HDM , the 1PR contribution is then predominant, which is very similar with the case in the SM. As we discussed in section II, these characters mean that the four differential distributions defined in last section cannot have a large deviation from SM predictions and we do not expect to see any new physics signal by the measurement of these four differential distributions. Indeed, we have checked our numerical results in a vast parameter space constrained by Eq.(65), and find it correct.

We are now ready to calculate the new physics corrections to the decay rate of $`BX_s\gamma \gamma `$. In order to reduce the effects of theoretical and experimental uncertainties, we define two ratios

$`R_{\gamma \gamma }(\mu )`$ $`=`$ $`{\displaystyle \frac{_{BX_s\gamma \gamma }^{III}(\mu )}{_{BX_s\gamma \gamma }^{SM}(\mu )}},`$ (96)
$`R_\gamma (\mu )`$ $`=`$ $`{\displaystyle \frac{_{BX_s\gamma }^{III}(\mu )}{_{BX_s\gamma }^{SM}(\mu )}},`$ (97)
where the four branching ratios denote the LO theoretical predictions in the SM and model III, respectively.

In the numerical calculation, we consider the following three typical scenarios.

* Scenario I: $`|\lambda _1|=|\lambda _2|=\lambda `$ and $`m_{H^\pm }`$ is about a few hundred GeVs. This scenario can be further divided into two limiting cases:

+ Case-1: $`\lambda _1=\lambda _2=\lambda `$.
    + Case-2: $`\lambda _1=\lambda _2=\lambda `$ or $`\lambda _2=\lambda _1=\lambda `$, where $`\lambda `$ is positive.
* Scenario II: $`\lambda _2=2/\lambda _1=\lambda `$ and $`m_{H^\pm }`$ is about a few hundred GeVs. By setting The $`\lambda `$ defined in this paper is different from that in Ref. by a factor of $`\sqrt{2}`$. $`\lambda =\sqrt{2}\mathrm{tan}\beta `$, one reproduces the numerical results of the Model II, a very popular version of the two-Higgs-doublet models.
* Scenario III: $`|\lambda _2||\lambda _1|`$ and $`m_{H^\pm }`$ is about $`200`$ GeV, which is favored by current experiments.

Let us firstly concentrate on Scenario I. From the explicit calculations, we find that the new physics contribution in case-1 tends to increase the value of $`C_7`$. In Figs. 17 and 17, we show the $`\lambda `$ dependence of the ratio $`R_{\gamma \gamma }`$ in case-1 for fixed $`m_{H^\pm }=300GeV`$ and assuming $`\mu =2m_b`$ (dashed curve), $`m_b`$ (solid curve) , $`m_b/2`$ (dotted curve), respectively. In determining the range of $`\lambda `$, we have put the constraint of Eq.(65) and required $`C_7`$ to be negative and positive respectively.

For the case-1, one can see from Figs. 17 and 17 that

* The $`\mu `$ dependence of the ratio $`R_{\gamma \gamma }`$ is weak for a negative $`C_7`$: we have $`0.7R_{\gamma \gamma }1.0`$ within the range of $`1/2\mu /m_b2`$.
* The $`\mu `$ dependence of the ratio $`R_{\gamma \gamma }`$ is strong for a positive $`C_7`$: we have $`0.3R_{\gamma \gamma }2.7`$ within the range of $`1/2\mu /m_b2`$; which implies that the new physics contribution to the branching ratio $`BX_s\gamma \gamma `$ can be significant.
* In the same region of $`\lambda `$ as specified in Figs. 17 and 17, the numerical results show that the ratio $`R_{\gamma \gamma }/R_\gamma `$ only varies within a very small range of $`1\pm 0.02`$ for $`1/2\mu /m_b2`$. This fact reflects the dominance of 1PR contribution to the decay rate of $`BX_s\gamma \gamma `$, and a strong correlation between the branching ratios of $`BX_s\gamma \gamma `$ and $`BX_s\gamma `$.

The main feature of case-2 is that the new physics contribution tends to decrease the value of $`C_7`$, as a result, coefficient $`C_7`$ keeps negative for all the values of $`\lambda `$. In Fig. 17, we show the dependence of $`R_{\gamma \gamma }`$ on $`\lambda `$ in case-2 for fixed $`m_{H^\pm }=300GeV`$ and $`\mu =2m_b`$ (dashed curve) , $`m_b`$ (solid curve), $`m_b/2`$ (dotted curve). It is easy to see that the $`\mu `$ dependence of $`R_{\gamma \gamma }`$ is weak. Study of the ratio $`R_{\gamma \gamma }/R_\gamma `$ gives the same conclusion as in case-1.

In Fig. 17, we plot the $`m_{H^\pm }`$ dependence of the ratio $`R_{\gamma \gamma }`$ in Scenario I. As one can see from this figure, $`R_{\gamma \gamma }`$ tends to approach 1 when the charged-Higgs boson becomes heavier, which is just the decoupling behaviour of the new heavy particle.

Secondly, we consider the Scenario II. By setting $`\lambda =\sqrt{2}\mathrm{tan}\beta `$, one reproduces the result of the popular model II: the second type of the two-Higgs-doublet models. Wilson coefficient $`C_7`$ keeps negative for all the values of $`\lambda `$. In model II, two additional free parameters will enter into our calculation: the charged Higgs mass $`m_{H^\pm }`$ and the ratio $`\mathrm{tan}\beta =v_2/v_1`$ with $`v_2`$ and $`v_1`$ are the vacuum expectation value of the two Higgs doublets $`\varphi _1`$ and $`\varphi _2`$.

In Fig.17, we plot the $`\mathrm{tan}\beta `$ and $`\mu `$ dependences of the ratio $`R_{\gamma \gamma }`$ in model II for fixed $`m_{H^\pm }=500GeV`$. In this figure, the three curves correspond to $`\mu =2m_b`$ (dashed curve), $`m_b`$ (solid curve), and $`m_b/2`$ (dotted curve), respectively. One can see that (a) the ratio $`R_{\gamma \gamma }`$ is larger than $`1.32`$ in the whole considered range of $`\mathrm{tan}\beta `$ and $`\mu `$, and (b) the ratio $`R_{\gamma \gamma }`$ has a weak $`\mathrm{tan}\beta `$ dependence for $`\mathrm{tan}\beta 3`$, but sensitive to the variation of scale $`\mu `$.

In Fig.17, we plot the $`\mathrm{tan}\beta `$ and $`m_{H^\pm }`$ dependences of the ratio $`R_{\gamma \gamma }`$ in model II for fixed $`\mu =m_b`$. In this figure, the three curves correspond to $`\mathrm{tan}\beta =1`$ (solid curve), $`10`$ (dashed curve), and $`100`$ (dotted curve), respectively. One can see that the ratio $`R_{\gamma \gamma }`$ is insensitive to the variation of $`\mathrm{tan}\beta `$.

Finally, we consider the Scenario III. The main characteristic of this Scenario is that $`C_7`$ increases with increasing $`\lambda _1`$ for fixed $`\lambda _2`$. In Figs. 17 and 17, we illustrate the $`\lambda _1`$ and $`\mu `$ dependences of the ratio $`R_{\gamma \gamma }`$ for fixed $`m_{H^\pm }=200`$ GeV and $`\lambda _2=50`$. In these two figures, the three curves correspond to $`\mu =2m_b`$(dashed curve), $`m_b`$ (solid curve), and $`m_b/2`$ (dotted curve), respectively. In determining the range of $`\lambda _1`$, we have considered the constraint of Eq.(65) and required $`C_7`$ to be negative and positive for Figs. 17 and 17, respectively. Again, the ratio $`R_{\gamma \gamma }/R_\gamma `$ only varies within a very small range of $`1\pm 0.02`$ for $`1/2\mu /m_b2`$.

In Fig. 17, we show the $`\lambda _2`$ and $`\mu `$ dependences of the ratio $`R_{\gamma \gamma }`$ for fixed $`m_{H^\pm }=200`$ GeV and $`\lambda _1=0.2`$. In this figure, the three curves correspond to $`\mu =2m_b`$(dashed curve), $`m_b`$ (solid curve), and $`m_b/2`$ (dotted curve), respectively. In determining the range of $`\lambda _2`$, we have required Eq.(65) to be satisfied. It is easy to see that the ratio $`R_{\gamma \gamma }`$ has a strong (moderate) dependence on $`\mu `$ ($`\lambda _2`$). The future NLO calculation of $`R_{\gamma \gamma }`$ in the SM and 2HDM’s will decrease the $`\mu `$ dependence. Study of the ratio $`R_{\gamma \gamma }/R_\gamma `$ also shows strong correlations between the branching ratios of $`BX_s\gamma `$ and $`BX_s\gamma \gamma `$ decays.

In Figs. 17 and 17, we plot $`m_{H^\pm }`$, $`\lambda _1`$ and $`\lambda _2`$ dependences of the ratio $`R_{\gamma \gamma }`$ in model III and assuming $`\mu =m_b`$. In Fig. 17, the soild and dashed curve refers to $`\lambda _1=0.018`$ and $`0.010`$, respectively. When choosing the value of $`\lambda _1`$, we keep $`C_7`$ to be negative and require $`|C_7|>|C_7^{SM}|`$ ( solid curve) and $`|C_7|<|C_7^{SM}|`$ (dashed curve), respectively. One can find from Fig. 17 that the ratio $`R_{\gamma \gamma }`$ tends to approach 1 when charged Higgs boson is becoming heavier. In Fig. 17, the three curves correspond to $`\lambda _1=0.2`$, $`\lambda _2=40`$(dashed curve), $`50`$ (soild curve), and $`60`$ (dotted curve), respectively. To determine the range of $`m_{H^\pm }`$, we have required Eq.(65) to be satisfied and $`C_7>0`$. Numerical results also show that for a much larger $`m_{H^\pm }`$, $`C_7`$ is driven to be negative.

## IV conclusion

In this paper, based on the low-energy effective Hamiltonian, we calculated the rare decay $`BX_s\gamma \gamma `$ in the SM and the general two-Higgs-doublet model with the restriction $`\lambda _{ij}^{U,D}=0`$ for $`ij`$. We focused on the estimation of the new physics contributions to the branching ratio $`(BX_s\gamma \gamma )`$ and the differential distributions $`(1/\mathrm{\Gamma })d\mathrm{\Gamma }/ds`$, $`(1/\mathrm{\Gamma })d\mathrm{\Gamma }/d\mathrm{cos}\theta _{\gamma \gamma }`$ and $`(1/\mathrm{\Gamma })d\mathrm{\Gamma }/dX_\gamma `$.

Within the considered parameter space allowed by currently available data, we found the following:

* In model III, the prediction of the branching ratio $`(BX_s\gamma \gamma )`$ ranges from one third to three times of the SM prediction, but it is highly correlated with the corresponding theoretical prediction of $`(BX_s\gamma )`$.
* In model II, the new physics enhancement to the branching ratio $`(BX_s\gamma \gamma )`$ can be as large as $`(3050)\%`$ with respect to the SM prediction.
* In the SM and model III, the contribution from 1PR diagrams is dominant and hence those four observable differential distributions are insensitive to the variation of scale $`\mu `$ as well as possible new physics corrections considered in this paper.
* Although the process $`BX_s\gamma \gamma `$ provides many new physical observables, it is not a better process to detect new physics than $`BX_s\gamma `$ because of the smallness of its decay rate and the long-distance background.

At the next-to-leading order, the coefficients of $`O_3O_6`$ and $`O_1^{}O_6^{}`$ may get enhanced. However, due to their subleading feature, we do not expect drastic deviation from the LO predictions.

## ACKNOWLEDGMENTS

The authors acknowledge the support by the National Natural Science Foundation of China under Grant Nos. 19575015, 19775012 and 10075013, and the Excellent Young Teachers Program of Ministry of Education, P.R.China.


format: text

text:
# Introduction.

## Introduction.

G. Ellingsrud et C. Peskine ont démontré () que, excepté pour un nombre fini de composantes du schéma de Hilbert, une surface lisse de $`𝐏_4`$ est de type général. D’autre part, les seuls exemples de surfaces lisses irrégulières de $`𝐏_4`$ connus sont d’irrégularité égale à 1 ou 2. Ces surfaces sont construites comme images réciproques par des endomorphismes d’un nombre fini de surfaces irrégulières connues. Et il est conjecturé de longue date que l’irrégularité des surfaces de $`𝐏_4`$ est bornée. Dans , Peskine conjecture que, excepté pour un nombre fini de composantes du schéma de Hilbert, la section hyperplane générale d’une surface lisse de $`𝐏_4`$ est linéairement complète (cf. §1). Il suggère aussi que la section hyperplane générale d’une surface obtenue comme image réciproque par un endomorphisme de $`𝐏_4`$ d’une surface irrégulière, est linéairement complète, excepté peut-être pour un nombre fini de familles de surfaces irrégulières. Nous montrons ici que toute section hyperplane d’une surface image réciproque d’une surface par un endomorphisme (qui n’est pas un automorphisme) est linéairement complète. Ce qui, vue la nature des surfaces irrégulières connues, va dans le sens de la conjecture de Peskine. La démonstration repose sur une décomposition de certaines images directes de fibrés (1.1) et conduit en fait à un énoncé général sur les sections hyperplanes de variétés images réciproques de sous-variétés lisses par un endomorphisme de l’espace projectif (2.2). La dernière section regroupe quelques remarques sur l’adjonction de ces variétés.

## 1 Décomposition de $`\pi _{}(𝒪_X^{}(lH^{}))`$.

Soit $`n`$ un entier $`1`$. Soient $`𝐏`$ et $`𝐏^{}`$ des espaces projectifs sur C de dimension $`n`$ et $`H`$ (resp. $`H^{}`$) une section hyperplane de $`𝐏`$ (resp. $`𝐏^{}`$). Nous dirons qu’un sous-schéma $`X`$ d’un espace projectif $`𝐏`$ est linéairement complet si le morphisme de restriction : $`\mathrm{H}^0(𝐏,𝒪_𝐏(H))\stackrel{}{}\mathrm{H}^0(X,𝒪_X(H))`$ est surjectif.
Soient $`𝐏^{}\stackrel{\pi }{}𝐏`$ un morphisme d’image non réduite à un point, donc fini, surjectif et plat, $`X𝐏`$ un sous-schéma de degré $`d`$ et $`X^{}=\pi ^1(X)𝐏^{}`$ son image réciproque. Si $`\pi ^{}(𝒪_𝐏(H))=𝒪_𝐏^{}(kH^{})`$, $`X^{}`$ est de degré $`d.k^{\mathrm{codim}\left(X\right)}`$. En particulier $`\pi `$ est de degré $`k^n`$.

Lorsque $`X`$ est localement Cohen-Macaulay, on désignera par $`\omega _X`$ son faisceau dualisant, et, dans ce cas, $`\pi `$ étant fini et plat, $`X^{}`$ est aussi localement Cohen-Macaulay et son faisceau dualisant est $`\omega _X^{}=\pi ^!(\omega _X)`$ (où $`\pi ^!()`$ est défini par $`\pi _{}\pi ^!()=\underset{¯}{\mathrm{Hom}}(\pi _{}(𝒪_X^{}),)`$, cf.  ou \[4, III.7\]). On souhaite étudier les faisceaux $`𝒪_X^{}(lH^{})`$ et, lorsque $`X`$ est localement Cohen-Macaulay, $`\omega _X^{}(lH^{})`$, et en particulier leur cohomologie. La proposition suivante et l’annulation des images directes supérieures de $`\pi `$ nous permettent de ¡¡travailler sur $`X`$¿¿.

###### Proposition 1.1

Le fibré $`\pi _{}(𝒪_X^{}(lH^{})),l\mathrm{}`$, se décompose naturellement :

$$\pi _{}(𝒪_X^{}(lH^{}))=\underset{d\mathrm{}}{}E_{l,d}𝒪_X(dH)$$

$`E_{l,d}`$ est le conoyau de la multiplication :

$$\mathrm{H}^0\left(𝒪_𝐏(H)\right)\mathrm{H}^0\left(\pi _{}(𝒪_𝐏^{}(lH^{}))((d1)H)\right)\stackrel{}{}\mathrm{H}^0\left(\pi _{}(𝒪_𝐏^{}(lH^{}))(dH)\right),$$

non nul exactement lorsque $`\left[l/k\right]d\delta (n,k,l)`$, avec $`\delta (n,k,l)=n+1+\left[\frac{n+1+l}{k}\right]`$.

* Démonstration de 1.1 Comme $`\pi `$ est affine, $`\pi _{}(𝒪_X^{}(lH^{}))=\pi _{}(𝒪_𝐏^{}(lH^{}))_{|X}`$, on peut donc supposer que $`X=𝐏`$. Or $`\pi _{}(𝒪_𝐏^{}(lH^{}))`$ est un fibré sur $`𝐏`$ et

$$\mathrm{H}^i\left(\pi _{}(𝒪_𝐏^{}(lH^{}))(dH)\right)=\mathrm{H}^i\left(𝒪_𝐏^{}((l+k.d)H^{})\right).$$

Le terme de droite s’annule pour $`0<i<n`$ et tout entier $`d`$. Donc $`\pi _{}(𝒪_𝐏^{}(lH^{}))`$ est une somme de fibrés en droites. Un tel fibré est engendré en degré $`\delta `$ si et seulement si

$$\mathrm{H}^n\left(\pi _{}(𝒪_𝐏^{}(lH^{}))((\delta n)H)\right)=0,$$

ce qui est équivalent à $`l+k(\delta n)>n1`$. Notons $`V_{l,d}=\mathrm{H}^0\left(\pi _{}(𝒪_𝐏^{}(lH^{}))(dH)\right)`$ (la composante homogène de degré $`l+kd`$ de l’anneau gradué de $`𝐏^{}`$). Le diagramme commutatif suivant (où $`l+kd0`$ et les flèches horizontales (resp. verticales) sont induites par la multiplication des sections sur $`𝐏`$ (resp. $`𝐏^{}`$)) :

montre que $`E_{l,d}=0`$ entraîne $`E_{l,d+1}=0`$ . D’où la décomposition annoncée.

* 1.2 Lorsque $`\pi `$ n’est pas un automorphisme (i.e. $`k>1`$), on vérifie que $`\mathrm{H}^0\left(_X^{}(H^{})\right)=0`$, donc $`X^{}`$ est non dégénéré.

* 1.3 D’autre part si $`\mathrm{h}^1\left(_X(dH)\right)=0\mathrm{pour}\mathrm{tout}d<s`$ alors

$$\mathrm{h}^1\left(_X^{}(lH^{})\right)=0\mathrm{pour}\mathrm{tout}l<sk.$$

En effet, si $``$ est un faisceau cohérent sur $`X`$, on a

$$\mathrm{H}^i\left(\pi ^{}()(lH^{})\right)\mathrm{H}^i\left(\pi _{}(𝒪_X^{}(lH^{}))\right).$$

En particulier si $`\mathrm{h}^0\left(𝒪_X\right)=1`$ et si $`\pi `$ n’est pas un automorphisme, $`X^{}`$ est linéairement complet. Remarquons cependant que $`\mathrm{h}^1\left(_X(dH)\right)=0\mathrm{pour}\mathrm{tout}d`$ si et seulement si $`\mathrm{h}^1\left(_X^{}(lH)\right)=0\mathrm{pour}\mathrm{tout}l`$.

* 1.4 Si $`X`$ est lisse on sait d’après le théorème de transversalité de Kleiman (voir  ou \[4, III.10.8, p. 273\]) que, pour un automorphisme linéaire général $`\sigma `$, $`\pi ^1(\sigma X)`$ est lisse. Donc si $`X`$ est lisse et en position générale, $`X^{}`$ est lisse et $`\mathrm{h}^i\left(𝒪_X^{}\right)=\mathrm{h}^i\left(𝒪_X\right),i<\mathrm{dim}\left(X\right)`$.

Dans toute la suite on supposera que $`\pi `$ n’est pas un automorphisme (i.e. $`k>1`$).

## 2 La section hyperplane.

Soient $`X𝐏`$ une variété lisse connexe, $`s`$ une section non-nulle de $`𝒪_𝐏^{}(H^{})`$ et $`\mu _s`$ la multiplication par $`s`$ : $`𝒪_X^{}\stackrel{\mu _s}{}𝒪_X^{}(H^{})`$.

###### Proposition 2.1

$`\mathrm{H}^i\left(\mu _s\right)`$ est injective pour $`i<\mathrm{dim}\left(X\right)`$.

Comme d’après 1.2 et 1.3 l’image réciproque de $`X`$ est linéairement complète et non dégénérée, on en déduit :

###### Corollaire 2.2

Soient $`X𝐏`$ une varieté lisse de dimension au moins égale à 2, $`\pi :𝐏^{}\stackrel{}{}𝐏`$ un endomorphisme qui n’est pas un automorphisme et $`X^{}=\pi ^1(X)`$. Alors pour tout hyperplan $`H^{}`$ coupant $`X^{}`$ proprement, la section $`X^{}H^{}`$ est linéairement complète. En particulier, si $`X^{}`$ est irréductible, toute section hyperplane de $`X`$ est linéairement complète.

* Démonstration de 2.1 D’après 1.1, $`\pi _{}(\mu _s)`$ se décompose sous la forme: $`\pi _{}(\mu _s)=\left(\mu _s^{p,q}\right)_{p,q0}`$ avec :

$$\mu _s^{p,q}:E_{0,q}𝒪_X(qH)\stackrel{}{}E_{1,p}𝒪_X(pH)$$

Remarquons que $`\mu _s^{\mathrm{0,0}}`$ s’identifie à :

$$𝒪_X\stackrel{s1}{}\mathrm{H}^0\left(𝒪_𝐏^{}(H^{})\right)𝒪_X$$

En appliquant le théorème d’annulation de Kodaira on voit que, pour $`i<\mathrm{dim}\left(X\right)`$, $`\mathrm{H}^i\left(\pi _{}(\mu _s)\right)`$ se réduit à :

$$\mathrm{H}^i\left(𝒪_X\right)\stackrel{s1}{}\mathrm{H}^0\left(𝒪_𝐏^{}(H^{})\right)\mathrm{H}^i\left(𝒪_X\right)$$

qui est évidemment injective.

Sans supposer $`X`$ lisse, on démontre de la même manière :

###### Proposition 2.3

Soit $`j`$ un entier. Si $`\mathrm{H}^i\left(𝒪_X(dH)\right)=0`$ pour tout entier $`i<j`$ et tout entier $`d`$ tel que $`0<d\delta (n,k\mathrm{,0})`$ alors $`\mathrm{H}^i\left(\mu _s\right)`$ est injective pour tout entier $`i<j`$.

## 3 L’application canonique.

On suppose toujours que $`X`$ est lisse, connexe et en position générale. Si $`0l<k`$ et $`\delta _l=\delta (n,k,l)`$, on a $`0<\delta _l<n+1`$. D’après 1.1 :

$$\pi _{}(\omega _X^{}(lH^{}))=\underset{¯}{\mathrm{Hom}}(\pi _{}(𝒪_X^{}(lH^{})),\omega _X)=(E_{l\mathrm{,0}}^{}\omega _X)\mathrm{}(E_{l,\delta _l}^{}\omega _X(\delta _lH))$$

Un élément non nul de $`E_{l,\delta _l}^{}`$ définit naturellement un morphisme non nul :

$$\pi ^{}(\omega _X(\delta _lH))\stackrel{}{}\omega _X^{}(lH^{})$$

* 3.1 Si $`l=0`$, on obtient une factorisation :

La flèche verticale de gauche est une projection.

* 3.2 Si $`l=1`$ et $`\delta _1\mathrm{dim}\left(X\right)`$, on sait (cf. \[1, 8.8.5, p. 239\]) que $`\omega _X(\mathrm{dim}\left(X\right)H)`$ est engendré par ses sections sauf si $`(X,𝒪_X(H))=(𝐏_m,𝒪_{𝐏_m}(1))`$. Donc, hormis ce cas, $`\mathrm{H}^0\left(\omega _X^{}(H^{})\right)`$ est $`0`$. Si $`s`$ est une section non nulle de $`\omega _X^{}(H^{})`$, le système linéaire $`s\mathrm{H}^0\left(𝒪_X^{}(H^{})\right)\mathrm{H}^0\left(\omega _X^{}\right)`$, défini une application rationnelle qui est un plongement hors du lieu des zéros de $`s`$ et donc l’application canonique de $`X^{}`$ est birationnelle.

En particulier :

###### Proposition 3.3

Soient $`S𝐏_4`$ une surface lisse et $`\pi :𝐏_4\stackrel{}{}𝐏_4`$ un endomorphisme qui n’est pas un automorphisme. Alors, si $`S`$ est en position générale, $`S^{}=\pi ^1(S)`$ est lisse et son fibré canonique est très ample, sauf si $`S`$ est un plan et $`k=2`$.

* 3.4 La surface $`S^{}`$ est donc de type général, sauf dans le cas d’exception de la proposition, $`S^{}`$ étant alors une surface de Del Pezzo.

* Démonstration de 3.3 Comme $`N_S^{}𝐏^{}=\pi ^{}(N_S𝐏)`$, alors $`\omega _S^{}\pi ^{}(\omega _S^{})^{}=(\omega _𝐏^{}\pi ^{}(\omega _𝐏)^{})_{|S^{}}=𝒪_S^{}(5(k1)H^{})`$. Donc $`\omega _S^{}=\pi ^{}(\omega _S(3H))𝒪_S^{}((2k5)H^{})`$. Cela démontre le résultat puisque, sauf si $`S`$ est un plan, $`\omega _S(2H)`$ est engendré par ses sections.

Guillaume Jamet

Université Pierre et Marie Curie

Institut de Mathématiques

Analyse Algébrique, Case 82

Tour 46–00, $`3^e`$ étage

4, place Jussieu

F-75252 PARIS Cedex 05

Tél : 01 44 27 26 90

Fax : 01 44 27 26 87

e-mail : jamet@math.jussieu.fr


text:
# The evolution of core and surface magnetic field in isolated neutron stars

## 1 Introduction

Since the process of the generation of neutron star (NS) magnetic fields (MFs) and, hence, their initial structure, strength and localization is still under discussion, we here intend to investigate the evolution of NS MFs which penetrate the entire star. Though there are arguments that the NS MF may be confined to the crustal layer, where it has been generated by thermoelectric effects during the early hot period of the NS’s life when the temperature gradients in the crust are immense (Blandford, Applegate & Hernquist 1983; Urpin, Levshakov & Yakovlev 1986; Wiebicke & Geppert 1996), one has also to consider the possibility that the NS MF permeates both the crust and the core of the star. This could be the field structure from the very beginning of the NS’s life, either according to the simple model of flux conservation during the collapse or due to a very efficient dynamo action in the core of the proto–NS as described by Thompson & Duncan (1993).
In a recent paper (Konenkov & Geppert 2000, hereafter KG00) we considered the flux expulsion from the superfluid core of a NS. Since the transition from the normal to the superfluid state of the matter in the cores of NSs may occur rather early after the NS’s birth (Page 1998), we assumed that the protons in the core form a superconductor of type II (Baym, Pethick & Pines 1969) and the magnetic flux is concentrated in an array of proton flux tubes (fluxoids). There are several forces acting on to the fluxoids (for details see Sec. 2), driving them outward into the normal conductive crust of a NS, where the magnetic flux may decay ohmically. Considering an initial field structure where the MF at the surface and in the core has the same strength, we showed that the field evolution in the crust leads to a deceleration of the flux expulsion which is the more pronounced the stronger the initial surface MF is. It turned out, that the main force, which is responsible for the flux expulsion from the core, is the buoyancy force.
We have also confirmed the result of Ding, Cheng & Chau 1993 (hereafter DCC), that not the entire magnetic flux is expelled from the core, but a certain part of it remains there for eternity. This results in a nonvanishing residual field $`B_{res}`$, which can be in the range of $`10^710^{10}`$ G, depending on the model parameters.
It is well known that the equation of state (EOS), which describes the state of matter in the core region of the NS, has a crucial influence on the crustal MF decay in isolated NSs (Urpin & Konenkov 1997; Page, Geppert & Zannias 2000). The EOS determines the compactness of the NS and, hence, the scale length of the MF in the crust, as well as the cooling history of the NS. While a softer EOS causes an smaller scale which, in turn, leads to a more rapid field decay, a softening of the EOS decelerates also the cooling. In a warmer crust the MF will be dissipated faster too. Since the EOS influences also - via the moment of inertia - the spin–down which is one important process for the core flux expulsion, we will consider in this paper the field evolution for NSs modelled with three different EOS, covering the whole range from very stiff to soft ones.
However, the assumption that the field strengths at the NS’s surface and in its core are initially of the same order of magnitude is not the only conceivable one. There are good reasons to assume that the field strengths at these positions differ considerably. This could be, e.g., a consequence of a field amplification during the collapse by flux conservation which results in a flux permeating the whole star with, say, $`10^9`$ G and a very efficient field generation after the formation of the NS by the thermoelectric instability in the crustal layers, producing there, say, $`10^{12}`$ G. Thus, we want to investigate here the effect of different initial ratios $`B_{p0}/B_{c0}1`$, that is, initial magnetic configurations, for which the crust contains the bulk of the magnetic flux and only a tiny fraction of it is anchored in the core.
Another process which determines the initial MF structure in the new–born NS is the post core–collaps accretion of fall–back matter after the supernova explosion. If this accretion is hypercritical the ram pressure overwhelms the pressure of a possibly existent MF and submerges it. The depth of submergence depends on the total amount of accreted matter and on the EOS of the NS matter (for details see Geppert, Page & Zannias 1999). Given the parameters as estimated for SN 1987A (Chevalier 1989), a MF generated by dynamo action in the proto–NS or amplified simply by flux conservation during the collaps would be submerged down to the crust–core boundary or even into the core. The rediffusion of that submerged MF would last at least $`10^8`$ years, thus, the NS born in SN 1987A will not appear as a pulsar in the near future, unless field generation processes in the crust will transform thermal into magnetic energy relatively fast. However, also fall–back much weaker than observed in SN 1987A will cause a certain submergence and, hence, a delayed switching on of pulsars (Muslimov & Page 1995). The discrepancy between the real age ($`1.710^5`$ years) and active age ($`1.610^4`$ years) of the pulsar B1757-24 (Gaensler & Frail 2000) can be explained at least partly by the assumption that the MF of this pulsar was submerged and rediffused during $`10^5`$ years (see Sec. 3). Generally, if there was initially some field generation in the whole NS, the post–supernova accretion onto the new born NS will lead to a field structure with $`B_{p0}/B_{c0}1`$, i.e. to a field structure having a much weaker MF strength at the surface than in the core. It is our aim to consider here the effect of different submergence depths onto the flux expulsion from the core and its consequences for the early NS evolution.
The paper is organized as follows. In Section 2 we give a short description of the model. It will be based on the model of the interplay of forces that determine the flux expulsion from the NS’s core presented in detail in KG00. The numerical results for the different initial field structures are presented in Section 3 and Section 4 is devoted to the discussion and conclusions.

## 2 Description of the model

We calculate the velocity of the flux carrying proton tubes (fluxoids) using the model described in detail in KG00. Namely, fluxoids are supposed not to be rigidly tied with the neutron vortices (Srinivasan et al. 1990, Jahan-Miri & Bhattacharya 1994), but to move under the common action of the buoyancy force, the drag force and the force exerted by the neutron vortices.
The superfluid core participates in the rotation of the NS by forming an array of quantized neutron vortices. The radial velocity of vortices at the crust-core boundary is determined by

$$v_n=\frac{R_c\dot{\mathrm{\Omega }}_s}{2\mathrm{\Omega }_s},$$
(1)
where $`\mathrm{\Omega }_s`$ is the averaged angular velocity of the core superfluid, $`R_c`$ is the radius of the NS core. In this paper we consider the evolution of isolated NSs. We assume, that the rotational evolution is determined by energy losses due to magneto–dipole radiation:

$$P\dot{P}=\frac{2B_p(t)^2R^6}{3Ic^3},$$
(2)
where $`B_p`$ is the surface field strength at the magnetic pole, $`R`$ is the radius of the NS, $`I`$ its moment of inertia, $`c`$ is the speed of light and the magnetic axis is perpendicular to the rotational one.

The buoyancy force, acting on the unit length of a fluxoid, is given by (Muslimov, Tsygan 1985):

$$f_b=\left(\frac{\mathrm{\Phi }_0}{4\pi \lambda }\right)^2\frac{1}{R_c}\mathrm{ln}\left(\frac{\lambda }{\xi }\right),$$
(3)
where $`\mathrm{\Phi }_0=2\times 10^7`$ G cm2 is the quantum of the magnetic flux, $`\lambda `$ is the London penetration depth, and $`\xi `$ is the proton coherence length, which is less than $`\lambda /\sqrt{2}`$ for a superconductor of type II. The buoyancy force is always positive, i.e. directed outward.

The drag force per unit length of a fluxoid, resulting from the interaction of the normal electrons in the core with the magnetic field of a fluxoid (Harvey, Ruderman, Shaham 1986, Jones 1987), is given by

$$f_v=\frac{3\pi }{64}\frac{n_ee^2\mathrm{\Phi }_0^2}{E_F\lambda }\frac{v_p}{c},$$
(4)
where $`n_e`$ is the number density of electrons in the core, which is assumed to be about 5% of the number density of neutrons, $`e`$ is the elementary charge, $`E_F`$ the Fermi energy of the electrons, and $`v_p`$ is the velocity of fluxoid. For the Fermi energy and electron number density we take the values determined by the density at the crust–core interface, $`\rho _c=210^{14}`$g/cm3. When $`v_p>0`$ (fluxoids are moving outward), one has $`f_v<0`$, i.e. $`f_v`$ hampers the expulsion of the flux from the core of NS.

The force, acting on the unit length of a fluxoid exerted by the neutron vortices, is given by (DCC)

$$f_n=\frac{2\mathrm{\Phi }_0\rho r\mathrm{\Omega }_s(t)\omega (t)}{B_c(t)},$$
(5)
where $`\omega =\mathrm{\Omega }_s\mathrm{\Omega }_c`$ is the lag between the rotational velocity of the core superfluid, $`\mathrm{\Omega }_s`$, and the rotational velocity of both the crust and the charged component of the core, $`\mathrm{\Omega }_c`$, while $`B_c`$ denotes the averaged core magnetic field. The observed spin period $`P`$ of a NS is related to $`\mathrm{\Omega }_c`$ by $`P=2\pi /\mathrm{\Omega }_c`$. The maximum lag $`\omega _{cr}`$, which can be sustained by the pinning force, defines also the maximum force, which can be exerted by the vortices on to the fluxoids (see DCC and Jahan–Miri 2000 for a more detailed discussion). According to DCC,

$`\omega _{cr}=8.7\times 10^2x_p\alpha _gr_6^1\left({\displaystyle \frac{m_pm_p^{}}{m_p}}\right)\left({\displaystyle \frac{m_p^{}}{m_p}}\right)^{1/2}\times `$
$`\times B_{c12}^{1/2}\mathrm{ln}\left({\displaystyle \frac{\lambda }{\xi }}\right)\mathrm{sin}(2\chi )\mathrm{rad}\mathrm{s}^1,`$ (6)
where $`x_p`$ is the fractional concentration of protons, $`\alpha _g`$ is a numerical factor of order of the unity, $`r_6`$ is the distance from the NS rotational axis in $`10^6`$ cm, $`m_p`$ is the mass of the proton, $`m_p^{}`$ its effective mass, $`B_{c12}`$ is the core magnetic field strength in units of $`10^{12}`$ G, and $`\chi `$ is the angle between rotational and magnetic axis. We assume hereafter $`x_p=0.025`$, $`\alpha _g=1`$ and $`m_p^{}=0.8m_p`$

Vortices can either move outward faster than fluxoids (forward creeping, $`\omega =\omega _{cr}`$), or the velocities of both kinds of flux tubes can coincide (comoving,$`\omega _{cr}<\omega <\omega _{cr}`$ ), or neutron vortices can move slower than fluxoids (reverse creeping, $`\omega =\omega _{cr}`$). The vortex acting force can be both positive or negative, depending on the sign of $`\omega `$, i.e. $`f_n`$ can either promote the expulsion or impede it.

The sum of the powers of the forces, acting onto the fluxoids in the core, is equal to the Poynting flux through the surface of the NS core:

$$\underset{fluxoids}{}(f_b+f_n+f_v)v_pdl=\frac{c}{4\pi }\underset{S_{core}}{}\left[\stackrel{}{E}\times \stackrel{}{B}\right]𝑑\stackrel{}{S}_{core},$$
(7)
where the summation runs over all fluxoids. The integral on the l.h.s. of (7) is taken over the length of each fluxoid, while the integral on the r.h.s. of (7) is performed over the surface of the core; the normal vector of this surface is directed into the core. We consider the evolution of a poloidal field, which is supposed to be dipolar outside the NS. Introducing a vector potential $`A=(0,0,A_\varphi )`$ and the Stokes stream function $`A_\varphi =S(r,t)\mathrm{sin}(\theta )/r=B_{p0}R^2s(r,t)\mathrm{sin}(\theta )/r`$, where $`r`$, $`\theta `$, $`\varphi `$ are spherical coordinates, $`s(r,t)`$ the normalized stream function, and $`B_{p0}`$ the initial magnetic field strength at the magnetic pole, one can express the field components as

$$B_r=B_{p0}R^2\frac{s}{r^2}\mathrm{cos}\theta ,B_\theta =\frac{B_{p0}R^2}{2}\frac{\mathrm{sin}\theta }{r}\frac{s}{r}.$$
(8)
By means of the variable $`s`$ the integral on the r.h.s. of the equation (7) can be rewritten as:

$$\frac{c}{4\pi }\underset{S_{core}}{}\left[\stackrel{}{E}\times \stackrel{}{B}\right]𝑑\stackrel{}{S}_{core}=\frac{B_{p0}^2R_c^2}{6}\frac{s(R_c,t)}{t}\frac{s(R_c,t)}{r}.$$
(9)
Assuming the MF in the core to be uniform, one can introduce the total core forces

$$F_{n,b,v}=f_{n,b,v}4R_c/3N_p,$$
(10)
where $`4R_c/3`$ is the mean length of the fluxoid and $`N_p=\pi R_c^2B_c/\mathrm{\Phi }_0`$ is the number of fluxoids. It is worth mentioning that $`F_nB_c^{1/2}`$ both in the forward and in the reverse creeping regimes, but $`F_b,F_vB_c`$.

One can now introduce the crustal force

$$F_{crust}=\frac{B_{p0}^2R^4}{6}\frac{1}{v_p}\frac{s(R_c,t)}{t}\frac{s(R_c,t)}{r}.$$
(11)
When $`s(R_c,t)/r`$ is positive and $`s(R_c,t)/t`$ is negative the crustal force is negative, i.e. it prevents the roots of the fluxoids to move towards the magnetic equator. The crustal force is proportional to $`B_c^2`$.

Now one can rewrite equation (7) as:

$$F_n+F_b+F_v(v_p)+F_{crust}(v_p)=0,$$
(12)
which determines the velocity of fluxoid $`v_p`$.

The equation of the balance of electromagnetic energy in the region outside NS core reads:

$$\underset{V_{\mathrm{crust}}}{}\frac{j^2}{\sigma }𝑑V+\frac{d}{dt}\underset{V}{}\frac{B^2}{8\pi }𝑑V+\frac{c}{4\pi }\underset{S_{core}}{}\left[\stackrel{}{E}\times \stackrel{}{B}\right]𝑑\stackrel{}{S}_{core}=0.$$
(13)
The first integral is restricted to the crust, the second has to be taken over all space excluding NS’s core, and the third is a surface integral over the core. Therefore, the crustal force determined by the equation (11) coincides with the crustal force determined by the equation (8) from KG00, however, a smoother numerical procedure can is achieved with the recent version.

The evolution of the magnetic field in the solid crust is governed by the induction equation without a convective term:

$$\frac{\stackrel{}{B}}{t}=\frac{c^2}{4\pi }\times \left(\frac{1}{\sigma }\times \stackrel{}{B}\right).$$
(14)
It can be rewritten in terms of the Stokes stream function as:

$$\frac{s}{t}=\frac{c^2}{4\pi \sigma }\left(\frac{^2s}{r^2}\frac{2s}{r^2}\right),$$
(15)
and should be supplemented by the boundary conditions at the surface ($`r=R`$) and at the crust-core boundary ($`r=R_c`$). While the outer boundary condition is given by $`R\frac{s}{r}=s`$ at $`r=R`$, the inner one is time-dependent and determined by equation (12).

The conductivity in the regions of interest in the crust is determined by the scattering of electrons at impurities and phonons. The phonon conductivity is dependent on the temperature, and, hence, on the cooling history of the NS. The impurity conductivity is temperature independent and is determined by the concentration of the impurities $`Q`$. The relaxation times for both electron scattering processes depend strongly on the density and on the chemical composition of the crust. For the phonon conductivity we use the numerical data given by Itoh et al. 1993, for the impurity conductivity we apply the analytical expression derived by Yakovlev & Urpin 1980.
The evolution of the magnetic field in the superconducting core is driven only by the motions of the fluxoids:

$$\frac{\stackrel{}{B_c}}{t}=\times \left(\stackrel{}{v}_p\times \stackrel{}{B_c}\right).$$
(16)
This equation can be rewritten in terms of the normalized stream function $`s`$ by:

$$\frac{s}{t}=v_p\frac{s}{r}.$$
(17)
The assumption of the homogeneity of the magnetic field in the core leads to the following ansatz for the fluxoid velocity

$$v_p(r,t)=\alpha (t)r,$$
(18)
and for $`s(r,t)`$ in the core

$$s(r,t)=\frac{B_c(t)r^2}{B_{p0}R^2}.$$
(19)
With this, the solution of equation (17) describes the evolution of the core MF at the crust–core boundary as

$$B_c(t)=B_{c0}\mathrm{exp}\left(2_0^t\alpha (t^{})𝑑t^{}\right).$$
(20)
Note, that in equations (15) and (17) of KG00 there are two misprints. However, the numerical computation was based on the correct equations. We use the procedure described in KG00 to find the velocity of the fluxoids and, consequently, the magnetic and rotational evolution of NS.

## 3 Numerical results

Although the induction equation is linear, the problem, as shown in Section 2, is a nonlinear one in terms of the magnetic field, and has a time and field dependent boundary condition at the crust–core interface. Additionally, equation (12) which determines the flux expulsion velocity is nonlinear in $`v_p`$ too. For the calculations described below we use the standard cooling scenario considered by Van Riper (1991) for the neutron star model based on a stiff (Pandharipande & Smith 1975), medium (Friedman–Pandharipande 1981) and a soft (Baym, Pethick & Sutherland 1971) EOS. The cooling histories for those EOS are taken from Van Riper (1991).
The main influence of the choice of the EOS on to the MF evolution comes via the different density profiles which determine the scale lengths of the MF in the corresponding NS models.
In all our model calculations we apply the initial rotational period $`P_0=0.01`$ s, a value which is generally accepted for new–born NSs.

### 3.1 The effect of different EOS

In order to extract the effects of different EOS on to the flux expulsion from the NS’s core we apply for all models the initial surface and core magnetic field strengths $`B_{p0}=B_{c0}=210^{12}`$ G, an impurity parameter $`Q=0.1`$ and the standard cooling scenario.
We consider three qualitatively different EOS and use the corresponding density profiles calculated for a $`1.4M_{}`$ NS (for details see Urpin & Konenkov (1997) and references therein). The soft EOS is represented by the BPS model with a radius of $`7.35`$km and a crust thickness of $`310`$m. The corresponding values for the medium (FP model) and stiff (PS model) EOS are $`10.61`$km , $`940`$m, and $`15.98`$km, $`4200`$m, respectively. The effect of increasing stiffness of the EOS describing the state of the core matter is at least threefold: i) since $`P\dot{P}B_p^2R^4/M`$ the increasing radius of the NS $`R`$ leads to a more efficient spin–down for a given surface MF $`B_p`$ and NS mass $`M`$; ii) the larger scale of the MF and iii) the faster cooling of the NS cause a deceleration of the crustal field decay.
In Figure 1 we show the evolution (from top to bottom) of the MF, velocities of vortices and fluxoids, forces and rotational periods of NS models based on a very soft (BPS) and a very stiff (PS) EOS. As in the case of the medium FP EOS (see KG00), at the beginning the dominating expulsive force is $`F_n`$. This is because in the core of the magnetized and rapidly spinning NS there exists a large number of fastly outward moving neutron vortices, which act on the fluxoids. As shown in the third panel of Figure 1, during that relatively short early period of $`t<10^410^5`$ years, $`F_n`$ is balanced by $`F_v`$, while the other forces in eq.(12) are negligible. However, because of the preceding NS spindown, the number of vortices in the core (and, hence, $`F_n`$) decreases, and $`F_b`$ becomes the main force, which drives the fluxoids outward. In KG00 we have shown, that in case of an intermediate EOS, the timescale of the expulsion of the magnetic field from the core $`\tau _e`$ is determined by the balance of $`F_b`$ and $`F_{crust}`$. One can see that, although $`F_b`$ and $`F_{crust}`$ are still the dominant forces governing the expulsion timescale both for NSs with stiff and soft EOS, in the case of the BPS model one can not neglect $`F_n`$ in the balance of forces. NSs with softer EOS spin down slower than NSs with stiffer EOS. In $`10^6`$ years a NS with the soft EOS spins down to only 0.19 s, while a NS with the stiff EOS spins down up to 1.26 s. Thus, $`F_n`$, which is proportional to the number of vortices, i.e. inversely proportional to the spin period of the NS, plays a more important role in the dynamics of fluxoids for NSs with a soft EOS. The evolution of their surface MF follows closely the core field decrease, because the diffusion timescale of the crustal currents is smaller than the expulsion timescale of the core field. The expulsion of the flux in case of the BPS model stops at about $`3.510^7`$ years, leaving behind a residual field $`B_{res}510^{10}`$ G, since then the vortex acting force balances $`F_b`$ and prevents further expulsion. Therefore, the slight spindown results in a relatively short residual spin period $`P_f0.5`$ s and a high strength of the residual MF.
In the case of a NS model with a stiff EOS, the expulsion timescale appears to be much shorter than the dissipation timescale of the field in the thick crust. When the MF in the core has been reduced by $`4`$ orders of magnitude (after $`210^8`$ years), the surface field is still $`10^{12}`$ G. Thus, the evolution of the surface MF in NSs with stiff EOS is insensitive to the details of the physics of the superfluid core (but, of course, is sensitive indirectly via the cooling history of NS which is extremely sensitive to the occurrence of superfluidity). In order to check the sensitivity of the results with respect to another a priori unknown parameter of our models, we varied $`Q`$ in the range $`0.011`$ and found, that the surface MF begins to decay for any $`Q`$ when the core flux is almost completely expelled.
Since the NS with stiffening of the EOS spins down more effectively (because of the larger radius and the longer decay–timescale of the surface MF), at the end of the evolution a smaller number of neutron vortices remains in the core which can hold there only a smaller number of fluxoids than in case of a softer EOS NS model. Thus, the residual field becomes weaker with stiffening of the EOS.

### 3.2 The case of $`B_{c0}<<B_{p0}`$

There exist some observational and theoretical evidences that the initial core field might be weaker than the surface one. Thus, Chau, Cheng & Ding (1992) concluded from an analysis of glitch observations in the Vela pulsar that the MF at the crust–core boundary is in the order of $`10^810^9`$ G while the surface MF is about three to four orders of magnitude larger. Such a difference in the field strengths at the surface and at the crust–core boundary could be produced by the thermoelectric instability (Blandford, Applegate & Hernquist 1983, Urpin, Levshakov & Yakovlev 1986, Wiebicke & Geppert 1996), which may transform thermal into magnetic energy effectively during the very early period of NS’s life in its outermost crustal layers.
In Figure 2 we show the initial field configurations in terms of the Stoke’s stream function $`s(r,t=0)`$ for $`B_{p0}=210^{12}`$ G and $`B_{c0}=10^9,10^{10},`$ and $`10^{11}`$ G. The calculations have been performed for an initial period $`P_0=0.01`$ s, an impurity parameter $`Q=0.1`$ and an initial density where the crustal field matches the inner one, $`\rho _0=10^{13}`$g cm-3. In order to show clearly the effect of different $`B_{p0}/B_{c0}`$, we will not vary $`Q`$ and $`\rho _0`$ as well as present the results for the medium FP EOS and standard cooling only.
We show in Figure 3 the evolution of the field components, of the forces acting on to the fluxoids in the core, of the velocities of both the fluxoids $`v_p`$ and of the neutron vortices $`v_n`$, and the rotational evolution in terms of the rotational period $`P(t)`$.
Since $`F_{crust}B_c^2`$, $`F_b,F_vB_c`$, but $`F_nB_c^{1/2}`$, a relative decrease of the field strength in the core reduces $`F_{crust}`$, $`F_b`$ and $`F_v`$ much faster than $`F_n`$, thereby enforcing the role the vortices play for the dynamics of the fluxoids for the weaker $`B_c`$. As seen from the third panels of Figure 3 apart from the early stages of evolution ($`t=10^3\mathrm{}10^4`$ years) the fluxoids are completely in the comoving regime, and their motion is totally determined by those of the neutron vortices, and the velocity of fluxoids at $`t=0`$ is the greater, the smaller $`B_{c0}`$ is. Consequently, a MF configuration with $`B_{c0}<<B_{p0}`$ results in a “spin–down induced” flux expulsion as described by Jahan-Miri & Bhattacharya (1994) for accreting NSs, and the influence of the crustal force on the dynamics of fluxoids is negligible.
The decay of the surface MF almost fully coincides with that in case of the purely crustal magnetic configurations (Urpin & Konenkov 1997, Page et al. 2000). Only in the case of $`B_{c0}=10^{11}`$ G, the value of $`B_p`$ during the plateau phase at $`t=10^510^7`$ years is slightly higher than in the case of $`B_{c0}=10^{10}`$ or even $`10^9`$ G. However, the most remarkable difference is that the surface MF decays not down to zero, but to a nonzero residual value. Since $`v_p=v_n`$ after $`t>10^4`$ years the ceasing of the spin–down of the NS, when $`B_p`$ has been decayed down to the weak $`B_c`$, results in a ceasing of the further field expulsion from the core. Thus, in the case of $`B_{c0}=10^{11},10^{10}`$ and $`10^9`$ G, $`B_{res}10^{10},210^9`$ and $`1.510^8`$ G, respectively.
It is worth mentioning that the choice of a stiffer EOS leads to a more extended phase of forward creeping. For such stiff EOS and relatively strong core fields ($`B_{c0}/B_{p0}0.1`$) after about $`10^8`$ years even the reverse creeping stage can be entered.

### 3.3 The effect of different submergence depths

A possible explanation for the discrepancy between the estimated rates of supernovae and of pulsar births in the galaxy (Frail 1998) is the delayed switch–on of a pulsar. This feature has been discussed by Muslimov & Page (1995), who considered a relatively shallow submergence of the initially present MF. However, the immediate consequence of a supernova (type Ib or II) explosion will be a fall back of matter on to the new born NS (Colpi, Shapiro & Wasserman, 1996). In case this fall back accretion is hypercritical the ram pressure of the back–falling matter is larger than the pressure of a, say, $`10^{12}`$ G initially present MF and the accretion flow is purely hydrodynamic; conditions as could be estimated for the SN 1987A (see Chevalier 1989). Applying the time dependence of the fall back accretion rates given by Colpi et al. (1996), Geppert, Page & Zannias (1999) calculated the submergence of a NS MF for different EOS as a function of the total amount of accreted matter $`\mathrm{\Delta }M`$. It became clear that for values as inferable for SN 1987A (Chevalier 1989) the MF has been submerged so deeply that the created NS will not shine up as a radiopulsar for more than $`10^8`$ years.
When the fall back accretion ceases the submerged field will start the rediffusion towards the surface of the NS. Once the NS is an isolated one during all its life, for given conductive properties of the layers above the field maintaining electrical currents and a given EOS determining the cooling and the scalelength of the field, the submergence density $`\rho _{sub}`$ down to which the MF has been buried by fall back accretion is the only parameter which defines the duration of that rediffusion process $`\tau _{red}`$. Evidently, the larger $`\mathrm{\Delta }M`$ the larger $`\rho _{sub}`$ and, hence, the larger $`\tau _{red}`$.
In order to show clearly the effect of the different $`\rho _{sub}`$ we will use the same representation as in the preceding section for the medium EOS with $`Q=0.1`$. In that way Figure 4 shows the evolution of the NS for $`\rho _{sub}=10^{12}`$, and $`10^{14}`$g cm-3, which corresponds roughly to an total amount of accreted matter of about $`3\times 10^5`$, and $`10^2M_{}`$, respectively.
The process of flux expulsion from a NS with initially submerged MF is qualitatively different from that of a standard NS. Since at the beginning of the NS’s life, after the fall back has finished, the surface MF $`B_{p0}`$ is practically zero, there is no spin–down and, hence, the velocity of the vortices $`v_n`$ is zero too. Thus, the neutron vortices, because being fixed in number and location, will not force the fluxoids to move outward, but will impede them to be expelled by the buoyancy force. The immobility of the fluxoids and the initial smallness of the crustal forces results in a situation where the balance of forces is given only by $`F_b+F_n=0`$, where the vortex acting force $`F_n`$ is just compensating the buoyancy, while $`F_v`$ and $`F_{crust}`$ are many orders of magnitude smaller than $`F_n`$, $`F_b`$.
When the field would remain buried there were no chance to expell the core flux because the number of vortices remains constant in a stationary rotating NS. However, depending on $`\rho _{sub}`$ the submerged crustal field rediffuses towards the surface thereby causing an increasing braking of NS’s rotation. As shown in Figure 4 that the spin–down onsets the later the deeper the crustal field has been submerged. The initially immobile neutron vortices start to move outward and the proton fluxoids, tightly bounded to them by $`F_b`$ follow their motion towards the crust, i.e. the fluxoids are in the comoving state. With increasing $`v_p`$ the drag force $`F_v`$ rises and, since the crustal force is still negligible, the balance of forces is determined by $`F_b+F_n+F_v=0`$. During that early epoch the core field does not change which results in a constant and positive $`F_b`$. Since, on the other hand both $`F_v`$ and $`F_n`$ are negative, the increasing absolute value of drag force reduces the absolute value of the vortex acting force. By the vertical line with mark $`t_1`$ we indicated the moment, when $`F_n=0`$; $`t_1`$ is $`310^2`$ years for $`\rho _{sub}=10^{12}`$g cm-3 and $`610^4`$ years for $`\rho _{sub}=10^{14}`$gcm-3. Exactly at that moment there is no interaction between the fluxoids and the vortices and the balance of forces is given by $`F_b+F_v=0`$.
However, the rediffusion of the crustal field continues beyond $`t=t_1`$ and the velocity $`v_n`$ (and hence $`v_p`$) increases too. The change in the sign of $`F_n`$ reflects the change in the state of flux expulsion: while early on the vortices hampered the movement of the fluxoids, now the positive $`F_n`$ supports the flux expulsion, though it is still $`v_n=v_p`$. The process of increasing $`F_v`$ and $`F_n`$ proceeds until $`F_n`$ reaches its maximal value, determined by $`F_n=F_n(\omega _{cr})`$. This moment is indicated in Figure 4 by $`t_2`$, which is $`610^2`$ years for $`\rho _{sub}=10^{12}`$g cm-3 and $`710^4`$ years for $`\rho _{sub}=10^{14}`$gcm-3. At that point the fluxoids can no longer follow the motion of the vortices, which cut now through the fluxoids, hence, the forward creeping regime has been reached. Close to this moment, the crustal force becomes comparable to $`F_b`$, however, both are small in comparison with $`F_n`$ and $`F_v`$. Later on, when the surface MF $`B_p`$ becomes comparable to the core field $`B_c`$, the transition from a “submerged” to a “standard” field configuration has been performed. The velocity of fluxoids is no longer determined by $`F_n`$ and $`F_v`$, but instead by $`F_{crust}`$ and $`F_b`$ , a situation which is the same as described in KG00 and the evolution goes through the phases of the forward creeping, the comoving, and the reverse creeping regime. It turns out, that the long–term magnetic and rotational evolution does not depend crucially on $`\rho _{sub}`$; it affects mainly the NS behaviour in its early life and may lead to the phenomenon that a young NS is recognized to be an old one because its active age, $`\tau _a=P/2/\dot{P}`$, is apparently high. However, especially for very deep submergence of the field, the pulsar may also turn to look younger than its real age when the rediffusion progresses (see Geppert et al 1999). The rediffusion of the initially submerged MF may lead to some interesting metamorphoses of ages of NSs: the young NSs look older or younger than they really are. For example, from the magnetic evolution of the NS at the right panel of Figure 4, one can estimate at $`t=130000`$ years an active age of about $`50000`$ years, and a braking index $`n2\ddot{P}P/\dot{P}^21.4`$ . The pulsar B1757-24 in SNR G5.4-1.2 has an active age of $`16000`$ years and a braking index $`n<1.33`$, while its real age was estimated from the morphology of the SNR to be $`170000`$ years (Gaensler & Frail 2000). Obviously, the assumption of MF submergence in NSs can at least partly explain the discrepancy between true and active age, and the occurence of a small braking index for this pulsar (Muslimov & Page 1995).
The other process which generates a strong surface field after the NS’s birth, the thermoelectric instability, may not account for those apparent discrepancies in the ages, because this instability acts effectively at most during the first few thousand years of the NS’s life.
Note however, that with increasing $`\rho _{sub}`$, $`B_{res}`$ increases while $`P_{res}`$ decreases. An increase of $`\rho _{sub}`$ from $`10^{12}`$g cm-3 to $`10^{14}`$gcm-3 causes an increase of $`B_{res}`$ by a factor of $`4`$ and decreases $`P_{res}`$ by a factor of $`2`$. This can be understood from the fact that for larger $`\rho _{sub}`$ the spin–down starts later. During that longer period a small amount of core flux has been expelled and dissipated in the inner crust as well as the NS cooled down, thereby enhancing the conductivity in the surface layers which decelerates the rediffusion. Thus, in the process of rediffusion the surface MF can not reach the same high value as for the corresponding case with shallow submergence. Therefore, the spin–down in case of the deep submerged MF is less efficient, $`P_{res}`$ is smaller, which reduces the efficiency of the flux expulsion from the core, resulting eventually in a larger $`B_{res}`$.

### 3.4 Shear stresses induced by flux expulsion

In this paper the crust of the NS is assumed to be solid and not deformed by Maxwell (and others) stresses, i.e. the crustal MF evolves only through Ohmic dissipation and does not affect the elastic properties of the crust. The opposite case was considered by Ruderman (1991a,b) and Ruderman, Zhu and Chen (1998), who investigated the consequences of the changing core MF configuration for the spin history, especially glitches. Moving the roots of the fluxoids at the crust–core interface, crustal strains are built up and relaxed by large–scale crust–cracking events. From our calculations we can estimate the shear stress at the crust–core boundary. The total force, acting upon the crust by the expelled flux is balanced by the force, acting from the crust on the fluxoids, i.e. $`F_{crust}`$. The maximum of $`\left|F_{crust}\right|`$ is reached at that stage of evolution, when the balance of $`F_b`$ and $`F_{crust}`$ determines the velocity of fluxoids. It seems to be a rather common situation, that the maximum force, which acts on the crust, does not exceed greatly the buoyancy force, which acts on the fluxoids from the beginning of the NS evolution. This force is about $`10^{37}`$ dyn, as one can see from Figure 1, for $`B_{c0}=210^{12}`$ G. The shear stress exerted at the base of the crust is given by

$$S_{crust}(t)\frac{\left|F_{crust}(t)\right|}{R_c\mathrm{\Delta }R},$$
(21)
where $`\mathrm{\Delta }R`$ is the thickness of the crust and $`\left|F_{crust}(t)\right|`$ is evaluated by use of equation (11). In Figure 5 we compare $`S_{crust}`$, calculated for $`B_{c0}=210^{12}`$G and $`Q=0.1`$, with the maximum shear stress the crust can sustain before yielding, as estimated by (Ruderman 1991a),

$$S_{max}\frac{\mathrm{\Delta }R}{R}\mu \theta _{max}<\mathrm{\hspace{0.33em}3}10^{26}\mathrm{dyn}/\mathrm{cm}^2,$$
(22)
where $`\theta _{max}`$ is the elastic strain limit, $`\theta _{max}<\mathrm{\hspace{0.33em}10}^2`$ for a $`10^610^8`$ years old NS, and $`\mu <\mathrm{\hspace{0.33em}10}^{30}`$ dyn cm-2 is the shear modulus of the crust. It is seen that the shear stresses reach their maximum after $`10^410^6`$ years and become negligible when the magnetic flux is expelled from the core. Clearly, the probability for crust–cracking increases with a softening of the EOS. Therefore, while for NSs with a soft EOS early crust–cracking events, as glitches, may occur, their appearance in young ($`t<10^5`$ yrs) NSs with stiffer EOS hints to an initial field strength larger than $`210^{12}`$ G and/or to a smaller yield strain ($`10^4`$) of the young NS crust.

## 4 Discussion and conclusion

In order to investigate the magnetic and rotational evolution of isolated NSs whose MF penetrates both its core and crust, we solved the system of equations which describe the expulsion of the MF from the superfluid core and its subsequent Ohmic dissipation in the solid crust. Thereby we assumed that the spin–down is driven by magneto–dipole radiation.
In a recent paper (KG00) we considered this process for a NS whose state of matter is described by a medium EOS and discussed the effects of different initial field strengths and impurity concentrations in the crust. For initial field strengths larger than $`10^{12}`$ G, our study yielded a decay time for the surface field $`>10^7`$ years. In the present paper we considered the effect of different EOS on the MF evolution and its consequences. We studied also qualitatively different initial MF structures.
In the case of a configuration, where the initial MF strengths at the surface and in the core are of the same order of magnitude (modelled in this paper by $`B_{p0}=B_{c0}`$) we have found, that the main force, which drives the fluxoids outward, is the buoyancy force $`F_b`$. The expulsion timescale is determined by the balance of $`F_b`$ and $`F_{crust}`$, which depend mainly on the core field strength and on the crustal conductivity, and is not affected by the spin history of the NS. The surprising result is that the expulsion timescale depends only weakly on the EOS. Since stiffer EOS results in a thicker crust, one would expect, that $`F_{crust}`$, which depends on the rate of dissipation of crustal currents, counteracts the expulsion of the magnetic flux from the core of a NS with the stiffer EOS for a longer time. However, the total core forces are proportional to $`R_c^3`$ (see eq. (10)), i.e., the total buoyancy force is about $`6`$ times stronger for the PS-model than for the BPS model. The interplay of these effects results in rather similar flux expulsion timescales for different EOS.
However, the evolution of the observable surface MF differs substantially for soft and stiff EOS. For the soft BPS and medium FP (see KG00) EOS the diffusion timescale of the crustal field appears to be much shorter than the core flux expulsion timescale. Thus, the evolution of the surface and core MFs are tightly connected. Moreover, for the soft and the medium EOS the decay timescale of the surface MF depends on its initial strength. In the case of the thick crust as modelled by the NS with the PS EOS we have found, that the diffusion of the crustal field lasts longer than the expulsion of the core field by approximately two orders of magnitude. That means, in case of a stiff EOS the surface field evolution is not sensitive to the details of physics of the forces acting upon the fluxoids in the core, and does not depend strongly on the initial MF strength. Hence, for a stiff EOS the evolution of the surface MF coincides with that, considered by Konar & Bhattacharya (1999), for an initially expelled flux. This result is valid for initial surface MFs in the range of $`10^{11}10^{13}`$ G, and for impurity parameters in the range of $`0.01<Q<1`$.
It is generally accepted that the MF of normal radiopulsars decays weakly during their lifetime (Bhattacharya et al. 1992, Hartman et al. 1996). Unfortunately, the results obtained from our scenario of field evolution do not allow to select or to reject some EOS as being not in accordance with observational facts. All models lead to decay timescales of the surface field comparable or exceeding the radiopulsars lifetime.
As in the case of the medium EOS (see KG00), not the entire magnetic flux is expelled from the core, but some part remains there for eternity, resulting in a residual MF of NSs. This can be important for the explanation of the long-living MFs of millisecond pulsars. However, for those NSs a more complex analysis is required, which takes into account all the effects of accretion occuring in binary systems.
For a “crustal” MF configuration, when initially almost the whole magnetic flux is confined to the NS’s crust, we have found that the neutron vortices play a very important role for the dynamics of the fluxoids. The velocity of the fluxoids appears to be equal to that of the vortices during almost the whole evolution, so that the flux expulsion can be really called “spin–down induced”(Srinivasan et al. 1990, Jahan Miri & Bhattacharya 1994). However, the evolution of the surface MF coincides almost completely with that of a purely crustal field as considered by Urpin & Konenkov (1997), and Page et al. (2000). Correspondingly, all conclusions of those papers are applicable to this magnetic configuration, namely the conclusion, that models of NSs, based on the medium and stiff EOS with standard cooling yield a satisfactory agreement with observations.
For a NS with an initially submerged crustal field we found that the submergence affects the evolution only during the early periods of the NS’s life, as long as the crustal field is diffusing back to the surface. Then, the further magnetic and rotational evolution of the NS becomes similar to the “standard” evolution as described in KG00. The only reminder of the submergence episode is the relatively weak correlation of the residual field strength with the submergence depth and its anti–correlation with the final rotational period.
We also calculated the shear stresses at the bottom of the crust, which arise when the footpoints of the fluxoids are moved along the crust–core interface. We found that total force, which acts on the crust from the moved fluxoids (being equal to the force, which acts from the crust on to the fluxoids, i.e. $`F_{crust}`$), does not exceed strongly the buoyancy force, which acts onto the fluxoids from the very beginning of the evolution. It can be strong enough to break the crust, especially when a soft EOS decribes the state of the core matter and the crust is relatively thin.

## Acknowledgment

The work of D.K. was supported by a scholarship of the Alexander von Humboldt–Stiftung and by RFBR grant 00-02-04011. We are indebted to the M. Reinhardt and O. Gnedin for carefully reading the manuscript and for suggestions to improve it.

## References

Baym G., Pethick C., Pines D. 1969, Nature, 224, 673

Baym G., Pethick C., Sutherland P.G., 1971, ApJ, 170, 299

Bhattacharya D., Wijers R.A.M.J., Hartman J. W., Verbunt F. 1992 A&A 254, 1990

Blandford R., Applegate J., Hernquist L., 1983, MNRAS, 204, 1025

Chau H. F., Cheng K. S., Ding K. Y. 1992 ApJ 399, 213

Chevalier R. 1989 ApJ 346, 847

Colpi M., Shapiro S.L. & Wasserman I. 1996 ApJ 470, 1075

Ding K. Y., Cheng K. S., Chau H. F. 1993 ApJ 408, 167

Frail D. A. 1998, in The Many Faces of Neutron Stars, eds a. Alpar, R. Buccheri, & J. van Paradijs, (Kluwer Academic Press: Dordrecht), p. 179

Friedman B., Pandharipande V. 1981, Nucl. Phys. A, 361, 502

Gaensler B., Frail D. 2000, Nature 406, 158

Geppert U., Page D., Zannias T. 1999, A&A 345, 847

Hartman J., Bhattacharya D., Wijers R., Verbunt F. 1996, A&A, 322, 477

Harvey J., Ruderman M., Shaham J. 1986, Phys. Rev. D. 33, 2084

Itoh N., Hayashi H., Koyama Y. 1993, ApJ 418, 405

Jahan-Miri M., Bhattacharya, D. 1994, M.N.R.A.S. 269, 455

Jahan-Miri M. 2000, ApJ 532, 514

Kluzniak W. 1998, ApJ 509, L37

Konar S., Bhattacharya D. 1999, M.N.R.A.S. 308, 795

Konenkov D., Geppert U. 2000, M.N.R.A.S. 313, 66

Muslimov A., Page D., 1995, ApJ, 440,L77

Muslimov A., Tsygan A. 1985, Ap&SS 115, 41

Page D. 1998, in Neutron Stars and Pulsars, eds N. Shibazaki, N. Kawai, S. Shibata, & T. Kifune, (Universal Academy Press: Tokyo), p. 183.

Page D., Geppert U., Zannias T. 2000, A&A 360, 1052

Pandharipande V. Smith R. 1975, Nucl. Phys. A. 237, 507

Ruderman M. 1991a, ApJ, 366, 261

Ruderman M. 1991b, ApJ, 382, 576

Ruderman M., Zhu T., Chen K. 1998, ApJ 492, 267

Srinivasan G., Bhattacharya D., Muslimov A., Tsygan A., 1990 Curr. Sci., 59, 31

Thompson Ch., Duncan R. 1993, ApJ, 408, 194

Urpin V., Konenkov D. 1997, M.N.R.A.S. 292, 167

Urpin V., Levshakov S., Yakovlev D. 1986, M.N.R.A.S. 219, 703

Van Riper, K. 1991, ApJS 75, 449

Wiebicke H.-J., Geppert U. 1996, A&A 309, 203

Yakovlev D., Urpin V. 1980, Sov. Astron. 24, 303