text:
# Light hadron spectroscopy with two flavors of 𝑂(𝑎) improved dynamical quarks Talk presented by T. Kaneko

## 1 Introduction

One of the major goals of lattice QCD simulation is to confirm the validity of QCD as the theory of strong interaction in the low energy region by comparing its prediction for the hadron spectrum to experiment. Because of huge computational demand to perform simulations with dynamical quarks, however, many works have been forced to neglect effects of dynamical sea quarks. In this quenched approximation, the CP-PACS collaboration found some deviation in the light hadron spectrum from experiment, which motivates us to investigate sea quark effects by performing simulations of QCD with the realistic number of dynamical quark flavors.

The JLQCD collaboration started numerical simulations of dynamical QCD on a supercomputer Hitachi SR8000 model F1, which is newly installed at KEK on March 2000. It has a 1.2 TFlops of peak performance and provides about 400 sustained GFlops for our simulation code. We carry out simulations of QCD with two flavors of dynamical quarks and investigate the chiral extrapolation and finite size effects in the light hadron spectrum, as a step toward studies of QCD with 2+1 flavors (i.e. $`u`$, $`d`$ plus $`s`$ quarks).

## 2 Simulations

We study QCD with two degenerate flavors of dynamical quarks, which are identified with $`u`$ and $`d`$ quarks; the strange quark is treated in the quenched approximation. We employ the $`O(a)`$ improved quark action with the clover coefficient $`c_{\mathrm{SW}}`$ determined non-perturbatively by the ALPHA collaboration. Simulations are performed at $`\beta =5.2`$ and $`c_{\mathrm{SW}}=2.02`$, for which $`a^1`$ evaluated at the physical mass of dynamical $`ud`$ quarks is approximately 2 GeV and hence scaling violation is not expected to be too large. We simulate 5 sea quark masses in the range $`m_{\mathrm{PS}}/m_\mathrm{V}0.8`$–0.6. Since finite size effects could be more important with dynamical quarks, we perform simulations on three lattices with different spatial sizes, $`12^3\times 48`$, $`16^3\times 48`$ and $`20^3\times 48`$. We have accumulated 3000 thermalized HMC trajectories on the $`16^3\times 48`$ lattice, while simulations on other lattices are in progress. Other details of simulation parameters are summarized in Table 1.

## 3 Static quark potential

We calculate the static quark potential on the $`16^3\times 48`$ lattice using the smeared Wilson loops. Since we do not observe any signal of flattening of the potential due to string breaking, potential data are parameterized with the form

$`V(𝐫)=V_0\alpha /|𝐫|g\delta V(𝐫)+\sigma |𝐫|,`$ (1)
where $`\delta V(𝐫)=(G(𝐫)1/|𝐫|)`$ represents a correction to the short-distance Coulomb term calculated with the lattice gluon propagator. Figure 1 shows a plot obtained for the heaviest sea quark mass, where we find that deviation of the potential data from the fit curve is at most 1% in the fitting range of $`𝐫`$. We determine the Sommer scale $`r_0`$  from the parameterization (1). The lattice spacing determined with the condition $`r_0=0.49`$ fm is listed in Table 1.

## 4 Light meson mass measurement

We calculate light hadron correlators for 6 values of the valence quark masses in the range $`m_{\mathrm{PS}}/m_\mathrm{V}`$ 0.80–0.50 for each sea quark mass. Measurements have been completed for the $`16^3\times 48`$ lattice, while only the degenerate case ($`K_{\mathrm{val}}=K_{\mathrm{sea}}`$) is measured on $`12^3\times 48`$ and $`20^3\times 48`$. Most of the analysis results discussed in this talk are based on the $`16^3\times 48`$ lattice data except for the finite size effect in the following.

Figure 2 shows the effective mass of degenerate ($`K_{\mathrm{val}}=K_{\mathrm{sea}}`$) pseudo-scalar meson at the second lightest sea quark mass. Data are shown for three different spatial volumes $`12^3`$, $`16^3`$ and $`20^3`$. We find that the fitted masses on the two larger lattices are in good agreement with each other both for the pseudo-scalar and vector mesons, as shown in Figure 3. This is also observed for other lattices with heavier sea quarks. This suggests that finite size effects are already small for the 163 lattice. We are currently extending similar measurements to the lightest sea quark $`K_{\mathrm{sea}}=0.1355`$, for which finite size effects are expected to be more important.

## 5 Sea quark effect in meson masses

Figure 4 shows the vector meson mass as a function of the pseudo-scalar meson mass squared. Since the effective lattice spacing decreases as the sea quark mass becomes smaller, which complicates the chiral extrapolation if made using the lattice unit, we plot the meson masses normalized by a physical quantity $`r_0`$ measured for each sea quark mass. We also plot the quenched data  in Figure 4, where we clearly see an indication of sea quark effect, i.e. slope in two-flavor QCD is significantly larger than that in quenched QCD. This leads to a larger hyperfine splitting of the strange meson as we shall see later. We also find that our results at different sea quark masses lie almost on one curve. This indicates that the sea quark mass dependence in this quantity is not large in our simulated region of $`K_{\mathrm{sea}}`$ and more precise calculation is needed, particularly for the vector mesons, to see it clearly.

## 6 Chiral extrapolation

We employ two different strategies for the chiral extrapolation. In our main strategy, which we call method (A), we use meson and quark masses normalized by $`r_0`$ throughout the analysis. For instance, the pseudo-scalar meson mass squared is fit to the form

$$(m_{\mathrm{PS}}r_0)^2=A_s(m_{\mathrm{sea}}r_0)+A_v(m_{\mathrm{val}}r_0),$$
(2)
where $`m_{\mathrm{sea}}`$ and $`m_{\mathrm{val}}`$ are respectively sea and average valence quark masses defined through the vector Ward identity (VWI) relation $`m_{q,\mathrm{VWI}}=(1/K1/K_c)/2`$. With this method, data is described very well by a simple linear fit. In the other strategy, method (B), we make the chiral extrapolation using meson and quark masses in the lattice unit. In this case, significant curvature is found in fits for both pseudo-scalar and vector meson masses, and therefore we have to introduce quadratic terms in the fitting function. We find that both fits give consistent results, although the fit error is generally larger for the method (B).

In each method, the physical quark masses $`m_{ud}`$ and $`m_s`$ are determined by tuning meson mass ratios to their experimental values. We use $`m_\pi `$ and $`m_\rho `$ to fix $`m_{ud}`$ and the lattice spacing. The strange quark mass $`m_s`$ is determined from either $`m_K/m_\rho `$ or $`m_\varphi /m_\rho `$.

Our results for the strange vector meson masses ($`K^{}`$ and $`\varphi `$), with pseudo-scalar ($`K`$) used to fix the strange quark mass, are shown in Figure 5 by filled circles (method (A)) or by filled squares (method (B)). For comparison, we also plot the CP-PACS data in the quenched approximation  and the recent two-flavor result . A trend that the sea quark effect pushes up the vector meson masses toward the experimental values is seen in this plot.

## 7 Strange quark mass

The bare strange quark mass is calculated using either of two definitions, i.e. one from VWI and the other from the axial vector Ward identity (AWI), $`2m_{q,\mathrm{AWI}}=_4A_4(t)P(0)^{}/P(t)P(0)^{}`$, where $`A_\mu (t)`$ and $`P(t)`$ are the axial vector current and pseudo-scalar density. We use the improved axial vector current $`A_\mu ^{(imp)}(t)=A_\mu (t)+c_A_\mu P(t)`$, with $`c_A`$ calculated at one-loop . The continuum $`\overline{\mathrm{MS}}`$ quark mass is obtained using one-loop matching  at scale $`\mu =1/a`$ and evolved to $`\mu =2`$ GeV with 3-loop $`\beta `$ function.

Our results for the strange quark mass using the method (A) are

$$m_s(2\mathrm{G}\mathrm{e}\mathrm{V})=\{\begin{array}{cc}\text{94(2) MeV }\hfill & \text{(VWI)},\hfill \\ \text{88(3) MeV }\hfill & \text{(AWI)},\hfill \end{array}$$
(3)
with $`K`$ used as input, or

$$m_s(2\mathrm{G}\mathrm{e}\mathrm{V})=\{\begin{array}{cc}\text{109(4) MeV }\hfill & \text{(VWI)},\hfill \\ \text{102(6) MeV }\hfill & \text{(AWI)},\hfill \end{array}$$
(4)
with $`\varphi `$ used as input. The method (B) gives consistent results.

Only the statistical errors are shown in (3) and (4). The systematic error is more significant, as indicated by the disagreement between VWI and AWI, and the difference between (3) and (4), although the latter is smaller than in the quenched results. One of the most important sources of the systematic error is the use of the one-loop perturbative $`Z`$ factor, for which a naive order counting gives $`O(5\%)`$. In addition, the discretization error of $`O(a^2)`$ or the quenching effect of the strange sea quark could also be important. A recent two-flavor simulation of CP-PACS predicted $`m_s(2\mathrm{G}\mathrm{e}\mathrm{V})=88_6^{+4}`$ MeV or $`90_{11}^{+5}`$ MeV in the continuum limit with the $`K`$ or $`\varphi `$ used as input .

This work is supported by the Supercomputer project No.54 (FY2000) of High Energy Accelerator Research Organization (KEK), and also in part by the Grant-in-Aid of the Ministry of Education (Nos. 10640246, 10640248, 11640250, 11640294, 11740162, 12014202, 12640253, 12640279 and 12740133). K-I.I, T.K and N.Y are supported by the JSPS Research Fellowship.


format: text

text:
# What Future Expects Humanity After the Demographic Transition Time?

## I Introduction

The problem of mankind population growth is one of the global problems of the mankind existence and its prosperity in the future. Will the demographic explosion existing now at the mankind population of the whole world stops and what will be after its stopping? Almost all the phenomenological theories of future of mankind ( see ,, , , ) include the time of demographic transition $`t_0`$. Before this time the growth of mankind population is very large ( hyperbolic or exponential growth), after this time the population growth may practically absent  or may be regulated by driving parameters  that govern the development of mankind (such parameters are present in non-linear open systems (see). Thus the time point $`t_0`$ (the time of demographic transition) is very essential for the future of mankind. What is biological and physical senses of the demographic transition time? The mankind as the whole is a very large system and its behavior can be described by methods of nonlinear dynamic, probability theory or statistical theory of open system using the methods of fractal geometry (, , or power laws ) . Now put the question: how the demographical transition time may be treated from these points of view? There are many different ways to answer on this question. In this paper we consider the demographic transition time as a point of a phase transition, i.e. the point of transition of mankind in a new phase ( in physiological, psychological, behaviour and so on senses). What this phase will be if it may exist? Will be the population of men of the new phase friendly to the old phase of mankind or will they fight with men of the old phase and force out them from their ecological niche? What scenario may exist and be admissible in that case? It seems that if the driving parameters that determine growth of the new phase ( thou it is difficult now to select the main of them, it is the tusk for many research groups ) will the same as now the new phase may diminish the old phase of mankind very soon. We would remind that all phenomenological theories give for the time $`t_0`$ the time in first part of twenty one century. If the treating of time $`t_0`$ as a point of phase transition is correct and has relation to existence of humanity , the embryos of the new phase must appear ( as it well known from theory of phase transitions in physics) before the time of demographic transition , i.e. in very near time. What role they will play in the life of mankind society and what it will be like? Will the new phase consist of international groups of islam terrorist and the fight and wars spread to the whole world? Or it will international groups of corrupted officials that for own prosperity govern the whole world to its death? May be the groups of sexual minorities will became the main part of mankind and diminish it ? Will the new phase characteristics lay in psychological (i.e. behaviour) domain or it lay in the physiological domain ( may be it is incurable genetics diseases or mental disorders reasoned by increasing radiation in our world)? Are the Russian revolution in 1917 year or the Germany fascist state in 1933 year are the examples of appearing of such phases? May be the victory of progressive mankind in the World War $`II`$ and collapse of socialist system in Soviet Union are results only the fact that the demographic time transition is not come? It seems, the search of answers on these questions is the one of main tusks of governments and research demographical gropes of developed countries in near time.
We consider in this paper only one way of mathematical analysis and description of two phases of mankind (the old phase and the new phase) and describe time behaviour of these phases by system of differential equations that are Volterra-Lotka equations for system of predatory-prey.

## II Equations for describing two phases mankind time behaviour

There are two simplest ways (physical and biological) for describing time dependence of two phases (the old and the new) of mankind :
1. In physics the phenomenological theory of phase transitions was developed by Landau (see ). Theory of phase transitions as appearing of order parameters developed by Ginzburg-Landau . These theory may be used for describing time behaviour of two phases. From physical point of view the transition from the old phase of mankind to the new phase may: a) take very short time (as transition of water to ice); b) take long time and has intermediate state when two phases of mankind will co-existent (as in the transition to superconductive state for superconductors of second order);c) be characterized by stretched in time axes phase transition as in some ferroelectrics. The differences of these cases determined by characteristics of the small embryos gropes of new phase and their interaction with the old phase of mankind. If interaction between men of old and new phases will attract the man of old phase to transition to the new phase (we will call such interaction as ”negative” ), the transition time will be small and co-existent of two phases of mankind for the long time period is impossible. Probably the evidences of such transitions are the above examples of Russian revolution of 1917 year and fascist Germany state of 1933 year. If result of interaction of men of the old phase with gropes of men of the new phase will neutral, i.e. men of old phase will not be attracted by embryos of new phase ( ”positive” interaction), amount of new the phase will increase and co-existence of new and old phases is possible. The time of coexistence determined by govern parameters and selected mathematical models.
2. The biological point of view on behavior of mankind allows to use many mathematical models for describing of time behaviour of old and new mankind phases. In this paper we see only simple predatory-prey models where the old phase of mankind will be prey and the new phase will be predatory.

## III Equations of Wolterra-Lotka for old and new phases mankind population

Let us designate the old phase of mankind by $`N_1(t)`$ and the new phase of mankind by $`N_2(t)`$ where $`t`$ is current time. Equations of Volterra-Lotka type is well researched. It is well considered the time behaviour of the predatory and the prey populations (including time behaviour populations near the point of demographic transition $`t_0`$ where fluctuations play important role). So we only write these equations and mark the main consequences of mathematical theories. Let for time $`dt`$ the old phase changes are $`cN_1^\nu dt\gamma N_1N_2dt`$ ( the second member is the result of diminishing character of influences (interactions) between the old and the new phases) and new phase changes are $`\delta N_2dt+bN_1N_2dt`$. The equations for $`N_1(t)`$ and $`N_2(t)`$ may be written in this case as

$$\frac{}{t}N_1(t)=cN_1^\nu (t)\gamma N_1(t)N_2(t)$$
(1)
$$\frac{}{t}N_2(t)=\delta N_2(t)+bN_1(t)N_2(t)$$
(2)
where $`c,b,\delta ,\gamma `$ are constant and defined by selections of considered models. These equations are Volterra equations for $`\nu =1`$ and are Lotka equations for $`\nu =0`$. The humanity phases $`N_1`$ and $`N_2`$ in (1),(2) have links and determine one another.
In the case when populations of the old phase and the new phase are multifractal sets on the time axes (see ) equations (1) and (2) read

$$D_{0,t}^\beta N_1(t)=cN_1^\nu (t)\gamma N_1(t)N_2(t)$$
(3)
$$D_{0,t}^\mu N_2(t)=\delta N_2(t)+bN_1(t)N_2(t)$$
(4)
where $`D_{0,t}^\beta `$ and $`D_{0,t}^\mu `$ are Riemann-Liouville fractional derivatives or generalized Riemann-Liouville derivatives (see)and $`\beta `$ and $`\mu `$ are fractional numbers or functions of govern parameters.
The equations (1) and (2) may be considered as the base equations of rude bioligical models for describing future population of mankind if there are exist two phases of mankind and sets of $`N_1`$ and $`N_2`$ are sets of population of old and new phases . As well known the equations of such type ( type of Volterra-Lotka equations) may be generalized and wrote in more general form ( form of Kolmogorov-Fokker-Plank equations for probability density $`w(N_1,N_2,t)`$ to find in moment $`t`$ the populations $`N_1`$ and $`N_2`$ in the considered system , , or may be considered the case when influences $`N_2`$ on behaviors of $`N_1`$ included only at the bound time intervals and so on). These generalized equations are wide spread and give many known models have used not only in biology (), but in economics (), non linear physics, theoretical chemistry and so on. In the next part of paper we write more detailed description of results of interactions between two mankind phases $`N_11`$ and $`N_2`$ based on well known researches of equations of such type. For case when in above equations there are only one phase $`N_1`$, and $`\nu =0`$ and coefficient $`c`$ is function of time and $`t_0`$ ($`c=c(t,t_0)`$ it is not difficult to chose $`c(t,t_0)`$ in the simple forms and as special cases receive the base equations of work ,, .

## IV How will interact old and new mankind phases?

The mathematical research of solutions characteristics of equations (1)-(2) or more general equations taking into account the fluctuations is well worked out (see for example . As was pointed out these equations may be written for mean values of $`N_1`$ and $`N_2`$ or in the form of govern equations for probability $`w(N_1,N_2,t)`$ to find the system in the moment $`t`$ with values $`N_1`$ and $`N_2`$, or in the form of Kolmogorov- Fokker-Plank equations that are consequences of governor equations. The influences of noises and fluctuations( depending at $`c,b,\delta ,\gamma `$ ) on the behaviour of system may be taken into account too. More complicated equations may be written for case when the phase $`N_2`$ may part of her time have ”a rest time” and in that rest time do not interact with the phase $`N_1`$. So we omit the process of receiving of solutions for different systems of equations describing interaction between old-new phases and only describe main results:
1. Possibility of co-existence between old and new phases with periodical changing of value of $`N_1`$ and $`N_2`$;
2. Possibility for diminishing of one phase by other phase that depends at the initial values;
3. Possibility for calculation of mean values of $`N_1`$ and $`N_2`$ in one of system states;
4. Possibility for calculation of time existence $`t_{ex}`$ for old or new phases if take into account by including in the main equations (or changing the main equation by its Fokker-Plank probability analogies equations) external or internal noises ( additive or multiplicative);
5. Possibility of description of fluctuation characteristics of kinetic variables near bifurcation points;
6. Possibility of researching of noises role in losing stability of stationary state or in different self-excited oscillations states;
7. Possibility of determination of conditions for stochastic regimes appearing in system;
8. Possibility of governing by processes of fighting or co-existence between old and new phases;
9. Possibility of using more complicated models and equations including fractal characteristics of mankind population distribution on the time axis.
10. Possibility of selection of different variants for govern parameters based on experimental sociological dates and so on.
Thus the predatory-prey models and its equations give rich possibilities for researching of the behavior of two phases of mankind. Has this assumption about appearing of two mankind phases relation to problem of future propagation of mankind? The future propagation of mankind will gives answer on this question, but if it assumption is correct it is time now to find the new phase (or phases) and determine (or make attempts to determine) is it one or other phases useful for future of mankind or not useful.

## V Conclusions

As one of advantages to use the equations (1)-(2) for describing demographic problems (with some of them the mankind already has confronted now and may be these problems are consequences of appearing the embryos of new phase of mankind ) we shall stress an opportunity of insert and account in the theory many factors (such as, incurable illnesses, natural cataclysmic etc.) defining future of mankind as a result of influences the control parameters included in the govern parameters of predatory-prey models and fractional dimensions .
There are three main scenario of models based on equations type of (1)-(2):
a) Two phases of mankind co-existent during long time with periodically changes its values and with mean values $`N_1=\delta b^1`$ and $`N_2=c\gamma ^1`$ conserves in bound limits if population fluctuations may be omitted;
b) Mankind population lost its stability as the result of fluctuation influences if the equations of Lotta type used;
c)The old mankind phase has bound time of existence and this time may by calculated if parameters of theory are defined .
The main purpose of this paper was to analyze possibility of introducing of two mankind’s population phases after demographic transition time. Driving parameters included in the considered models (equations of Landau-Ginzburg time theory or equations of predatory-prey type) of the phenomenological theories of the mankind population in the frame of physical Landau or Ginzburg-Landau models or biological predatory-prey models. The predatory-prey type models was chosen asan example. The change of number of mankind (described in the framework of phenomenological theories of the population) can be adjusted by such choice of control parameters with concrete contents of dependencies of them at which the mankind will has stability time and population will grow so slowly, that overpopulation and the problems connected with fighting between two mankind phases will do not arise in the foreseen future. Stress now as was done in  that after demographic transition time also the type of economical development of humankind may be changed .


text:
# Euler, Hipparcos and the five dwarfs : Reconstructing star formation histories in the Local Group

## 1. Introduction

The history of the star formation rate in galaxies is one of the main ingredients required to understand their formation and evolution. Whilst the Roberts’ time scale for gas consumption by star formation gives a rough idea of the evolution of the stellar and gas content, more elaborate models (as summarised, e.g. by Sandage, 1986) can reproduce a wide variety of observables, from the evolution of the bulge-to-disc ratio in the Hubble sequence to the integrated colours as a function of redshift. However, in sharp contrast with these global views where the star formation rate SFR(t) is a monotone function of time, current bottom-up hierarchical scenarios of structure formation predict far more complicated histories, where star formation episodes are related to, if not driven by, halo mergers. The star formation history is, in this context, a more or less faithful reproduction of the merger history, rich in both strong and minor events at redshifts $`z<2`$ or so.

While it is possible to understand some of the observed properties of high redshift galaxies in this framework, a far more direct test can be achieved with the fossil record of the star formation history itself: the colour-magnitude diagrams (CMDs) of the resolved stellar populations in nearby galaxies. Progress in this field has been steady, driven both by ground-based wide-field surveys and by deep, high-resolution, narrow-field HST studies, reaching well below the main sequence turn off point in many satellites of the Local Group (see Aparicio 1998, Mateo 1998 for recent reviews).

## 2. Objective reconstructions of star formation histories

The unprecedented quality of the current CMDs requires powerful tools to invert them in order to derive the star formation history (SFH) which gave rise to the observed distribution of stars in these diagrams. The methods used so far are based on comparisons between the observed CMD and synthetic CMDs computed assuming a given SFH, and then looking for the best matching SFH. For instance Tolstoy & Saha (1996) use a Bayesian likelihood technique, Mighell (1997) a classical $`\chi ^2`$ optimisation, Ng (1998) a Poisson merit function, while Gallart et al. (1999) and Hurley-Keller et al. (1999) minimise a counts in cells statistic for well-defined boxes in the CMD. Alternatively, Dolphin (1997) makes a linear decomposition in terms of fiducial CMDs and solves for the best matching final CMD, in terms of a $`\chi ^2`$ statistic.

One of the problems with most of these techniques is that there is no guarantee that the actual SFH belongs to the set of (usually parametric) SFHs explored within a family of functions defined a priori. For instance, in the case of the Carina dSph, Hurley-Keller et al. (1999) find the best fitting 3-burst solution, leading to results which are at variance with those from Mighell (1997), who uses a non-parametric approach. In addition, it is unclear to what extent these best matching solutions are actually good fits to the observed CMDs. We have therefore developed an entirely new method based on a combination of Bayesian statistics with variational calculus which does not suffer from the limitations listed above. Full details of the method are given in Hernandez et al. (1999, Paper I), and will not be repeated here. Very briefly, the method makes 3 key assumptions: (1) the metallicity of the ensemble of stars is known and has a small dispersion; (2) the initial mass function is given and there are no unresolved binary systems; and (3) both distance and colour excess are known to within some uncertainties. The method then takes as inputs the positions of $`n`$ stars in a colour magnitude diagram, each having a colour $`c_i`$ and luminosity $`l_i`$, with (in this example, uncorrelated) associated errors $`\sigma (c_i)`$ and $`\sigma (l_i)`$ respectively. Using the likelihood technique, we first construct the probability that the $`n`$ observed stars resulted from some function $`SFR(t)`$. This is given by

$$=\underset{i=1}{\overset{n}{}}\left(_{t_0}^{t_1}SFR(t)G_i(t)𝑑t\right),$$
(1)
where

$`G_i(t)={\displaystyle _{m_0}^{m_1}}{\displaystyle \frac{\rho (m;t)}{2\pi \sigma (l_i)\sigma (c_i)}}\mathrm{exp}\left({\displaystyle \frac{D(l_i;t,m)^2}{2\sigma ^2(l_i)}}\right)\mathrm{exp}\left({\displaystyle \frac{D(c_i;t,m)^2}{2\sigma ^2(c_i)}}\right)𝑑m`$ (2)
In this expression $`\rho `$ is the density of stars of mass $`m`$ along the isochrone of age $`t`$, and only depends on the assumed IMF and the set of stellar tracks (and in particular the durations of the different evolutionary phases). The $`D`$ factors are the differences in luminosity and colour of the observed star $`i`$ with respect to the luminosity and colour of a star of mass $`m`$ at time $`t`$. We refer to $`G_i(t)`$ as the likelihood matrix, since each element represents the probability that a given star $`i`$ was actually formed at time $`t`$ with any mass.

Following the discussion of Paper I, we may write the condition that the likelihood has an extremal as the variation $`\delta (SFR)=0`$, allowing a full variational calculus analysis to be used. Developing first the product over $`i`$ using the chain rule, and dividing the resulting sum by $``$, one obtains

$$\underset{i=1}{\overset{n}{}}\left(\frac{\delta _{t_0}^{t_1}SFR(t)G_i(t)𝑑t}{_{t_0}^{t_1}SFR(t)G_i(t)𝑑t}\right)=0$$
(3)
Introducing the new variable $`Y(t)`$ defined as

$$Y(t)=\sqrt{SFR(t)}𝑑tSFR(t)=\left(\frac{dY(t)}{dt}\right)^2$$
(4)
into Equation 3 we can develop the Euler equation to yield

$$\frac{d^2Y(t)}{dt^2}\underset{i=1}{\overset{n}{}}\left(\frac{G_i(t)}{I(i)}\right)=\frac{dY(t)}{dt}\underset{i=1}{\overset{n}{}}\left(\frac{dG_i/dt}{I(i)}\right)$$
(5)
where

$$I(i)=_{t_0}^{t_1}SFR(t)G_i(t)𝑑t$$
(6)
is an integral constraint.

We have now transformed what was an optimisation problem, finding the function that maximises the product of integrals defined by Equation 1, into an integro-differential equation with a boundary condition (at either $`t_o`$ or $`t_1`$) which can be solved by iteration to get, non-parametrically, the function $`SFR(t)`$. Further details on the numerical aspects of the procedure are available in Paper I. It is important to point out that our method has distinctive advantages over other techniques: (1) the variational calculus allows a fully non-parametric reconstruction, free of any astrophysical preconceptions; (2) there is no time-consuming comparisons between CMDs, since the function (not the parameters) $`SFR(t)`$ that maximises the likelihood if solved for directly; (3) the CPU scales linearly with the time resolution required in the reconstruction.

Paper I also presents a detailed study on the influence of the assumptions that were made: the effect of changing the IMF or the presence of unresolved binaries is essentially a normalisation problem, which does not change the overall shape or localisation of a burst of star formation, while a wrong estimate of the metallicity has drastic effects on the position of a burst, a result of the age-metallicity degeneracy. See also Paper I for the effect of photometric uncertainties in the effective time resolution of the reconstructed SFH. Note that since the IMF and the fraction of binaries (and the distribution function of their mass ratios) are unlikely to be measured, an absolute normalisation of the star formation rate cannot be achieved.

## 3. Star formation histories of Local Group dSph satellites

Given the main assumption made so far in order to derive a fully non-parametric variational reconstruction of the SFH, namely a determination of the metallicity of the ensemble of stars, only a few galaxies fulfill this requirement. In addition, in order to minimise any possible systematic effects between different data sets and reduction procedures, as well as crowding/blending effects, we selected from the HST archive an homogeneous sample of five dSph galaxies (Leo I, Leo II, Draco, Carina and Ursa Minor) of the Local Group which were reduced with the same procedures. It is only such an internally consistent data set which allows robust comparisons between different galaxies, with the proviso that the volumes sampled by the WFPC2 field may not be representative of the entire galaxy. Full details of the SFH reconstructions for these galaxies are presented in Hernandez et al. (2000, Paper II).

The dSph UMi can be used as an example of our procedure. Figure 1 gives the WFPC2 CMD we used as input. Note the increase of the 2$`\sigma `$ error bars with decreasing luminosity, and the removal of stars on the horizontal branch (incompatible with the evolutionary phases included in our isochrones). Since incompleteness may be important, and isochrones are degenerate at this level, we also remove stars below $`M_V=5`$, so we are left with $`N=372`$ stars in total. For the fiducial values of distance, metallicity and colour excess indicated (slightly different but consistent with the more recent determinations by Mighell & Burke, 1999), the variational method gives the reconstructed SFH indicated with the black line on the right panel of Fig. 1. This function maximises the likelihood, as defined in Eq. 1, for the given data set. To assess the robustness of this reconstruction, many functions were reconstructed by changing the observational parameters within their error bounds indicated on Fig. 1. The envelope of such solutions is given as the shaded area on the figure, and represents the uncertainties in the reconstruction given by the uncertainties in the observations, so that any function that will be contained within these bounds will maximise $``$, given these errors. The bulk of the stellar population in UMi was therefore formed more than 12 Gyr ago, with a peak at 14 Gyr. Some residual activity may have been present at younger ages, but it is caused by stars bluer than the main turn off point or by blue stragglers. Given the uncertainties in the photometry, we cannot resolve any episodes in the star formation rate within the main burst at 14 Gyr, but we can conclude that the burst lasted less than 2 Gyr (FWHM).

Having found the best solution (in terms of maximising the likelihood) does not guarantee that it is also a good solution. The classical way to assess the goodness of fit, via Fisher’s information matrix, does not apply here (there are no parameters), so we apply again Bayesian techniques to check whether the observed CMD could result from the CMDs produced by the best reconstruction of the SFH. To do this, we do a counts in cells analysis and apply Saha’s W statistic (Saha, 1998), basically

$$W\underset{i=1}{\overset{B}{}}\frac{(m_i+n_i)!}{m_i!n_i!}$$
(7)
with $`B`$ number of cells the CMD is divided into, and $`m_i`$ and $`n_i`$ the number of stars in cell $`i`$ in two CMDs. To compute the distribution function of $`W`$, a series of model-model comparisons are made, that is, synthetic CMDs resulting from realisations of the reconstructed SFH are compared pairwise. To check whether the observed CMD could arise from this distribution, the $`W`$ statistic resulting from comparisons between the observed CMD and a series of CMD realisations are made. If the distribution functions are compatible, at some statistical level, the hypothesis that the observed CMD can be produced by the reconstructed SFH is accepted.

In the case of UMi, the mean and dispersion of the model-model $`W`$ distribution is 47 $`\pm `$ 5, while the data-model distribution gives 44 $`\pm `$ 4, so we can deduce that the synthetic CMDs are compatible at better than 1$`\sigma `$ with the observed CMD. The best solution is also a good solution.

The dSph Leo II (see Fig. 2) presents a similar case, whereby the bulk of stars was formed in a burst slowly rising at 12 Gyr with a maximum around 8 Gyr and an abrupt decrease at 6 Gyr. However in this case, the very sensitive $`W`$ statistic says that the observed CMD is not compatible (at the 2$`\sigma `$ level) with the synthetic diagrams resulting from this best solution. The reason is likely to be the relatively large spread in metallicity, around 0.3 dex, showing the limitation of the method.

The case of the Carina dSph (Fig. 3) is more interesting. There are clearly two sub-giant branches, identified with two epochs of star formation at 8 and 3 Gyr. Given the 2$`\sigma `$ errors indicated in Fig. 3 and simulated CMDs with these errors, we can say that the bursting episodes in Carina lasted more than a 1 Gyr, and that there was no epoch in the 2–8 Gyr interval without star formation activity. The $`W`$ statistic indicates a good solution. It is however puzzling that such an extended star formation history, over 6 Gyr, resulted in such a small metallicity, with a dispersion of 0.2 dex or less. Either these values are only representative of the tip of the giant branch, or the enrichment in metals had a non standard history.

Leo I presents a similar, extended, history of star formation (see Fig. 4). Some activity was present around 10–13 Gyr ago, compatible with the presence of horizontal branch stars recently detected at the NTT (Held et al. 2000) but not present in our HST data. The bulk of the observed stellar populations was formed continuously, with peaks of activity at 8 and 4 Gyr, in agreement with an independent analysis by Gallart et al. (1999). In this case however, the $`W`$ statistic indicates that the maximum likelihood solution is not a good solution at the 2$`\sigma `$ level. Again, a likely explanation is that the dispersion in metallicity (around 0.3 dex) is too large for the method to work. A full bi-variate non-parametric reconstruction of both $`SFR(t)`$ and $`Z(t)`$ is required.

If the Local Group is a representative volume of the Universe, the star formation activity at high redshift was dominated by these dSph that dominate the luminosity function. A properly averaged SFH does not show any significant epoch in the comoving density of SFR (Tolstoy 1999, Paper II). There seems to be no link either between the epoch of the bursts and perigalacticon passages of these satellites.

## 4. The evolution of the SFR in the solar neighbourhood

The evolution of the star formation rate in the Galactic disc is another basic function required to understand the formation of the disc, its chemical evolution and, more generally, the luminosity evolution of spiral galaxies. Previous attempts at determining the SFH in the Galaxy have relied on indirect methods, basically relations between age and some astrophysical property –such as chromospheric activity or metallicity– and then correct with evolutionary models for the stars which have disappeared from the sample. A good example of the complexities of these techniques is provided by Rocha-Pinto et al. (2000b). Here we use the method outlined in § 2. since the solar neighbourhood has a small dispersion in metallicity centred on the solar value. Further details will be found in Hernandez et al. (2000b, Paper III).

The first step is to define a volume-limited sample of well-measured stars, and the ideal catalogue for this task is obviously the Hipparcos catalogue. This provides a direct way to infer the SFH of the solar neighbourhood, as opposed to more indirect methods (e.g. Rocha-Pinto et al. 2000a). Although the completeness of Hipparcos varies both with spectral type and galactic latitude, a cut at $`V=7.9`$ for the sample with parallax errors smaller than 20% provides a reasonable sub-sample, once binaries and variable stars are removed. A typical cut at $`V=7.25`$ can produce an absolute-magnitude limited sample complete in volume with a well understood error distribution. The absolute magnitude limit, say at $`M_V=3.15`$, implies that only stars younger than about 3 Gyr enter the sample, but the kinematical and geometrical corrections are minimised. With such samples, and their error distributions, we simulated synthetic CMDs with two different SFHs: a 3 burst scenario and a constant SFR over 1.5 Gyr. Figure 5 shows that our maximum likelihood variational method reconstructs correctly such star formation histories, even though the number of stars is very small. Note also that the CMDs look very similar, yet they are produced by very different SFHs.

The inversion procedure applied to the actual Hipparcos diagram for the solar neighbourhood is shown on Fig. 6. The superb quality of the data allows us to reconstruct its SFH with the unprecedented resolution of 50 Myr, and a clear pattern emerges: over a roughly low-level constant SFR, there are several distinct peaks, separated by about 0.5 Gyr. As before, the envelope gives the area where any function will maximise the likelihood, when different $`M_V`$ cuts are applied. Note that the envelope increases with age, a reflection of larger uncertainty. This best solution is also a good solution (see Paper III).

A possible interpretation of the quasi-periodic episodes in the star formation rate of the solar neighbourhood may be given in terms of interactions with either spiral arms or the galactic bar. As the pattern speed and the circular velocity are in general different, the solar neighbourhood periodically crosses an arm region, where the increased local gravitational potential might possibly trigger an episode of star formation. In this case, the time interval $`\mathrm{\Delta }t`$ between encounters with an arm at the solar neighbourhood is

$$\mathrm{\Delta }t=\frac{0.22\mathrm{Gyr}}{m}\left(\frac{\mathrm{\Omega }}{29\mathrm{km}\mathrm{s}^1\mathrm{kpc}^1}\right)^1\left|\frac{\mathrm{\Omega }_p}{\mathrm{\Omega }}1\right|^1$$
(8)
where $`m`$ is the number of arms in the spiral pattern. Recent determinations tend to point to large values for the pattern speed (e.g. Mishurov et al. 1979, Avedisova 1989, Amaral and Lépine 1997) close to $`\mathrm{\Omega }_p2324`$ km s-1 kpc-1, which would imply that the regularity present in the reconstructed $`SFR(t)`$ would be consistent with a scenario where the interaction of the solar neighbourhood with a two-armed spiral pattern would have induced the star formation episodes we detect. This is consistent with the current knowledge on the nearby spiral structure (Vallée 1995), and is reminiscent of the explanations put forward to account for the inhomogeneities observed in the Hipparcos velocity distribution function, where well-defined branches associated with moving groups of different ages (Chereul et al. 1999, Skuljan et al. 1999, Asiain et al. 1999) could perhaps be also associated with an interaction with spiral arm(s), although in this case the time scales are much smaller. Of course, other explanations are possible; for example the cloud formation, collision and stellar feedback models of Vazquez & Scalo (1989) predict a phase of oscillatory star formation rate behaviour as a result of a self-regulated star formation régime. Close encounters with the Magellanic Clouds have also been suggested to explain the intermittent nature of the star formation rate, though on longer time scales (Rocha-Pinto et al. 2000b).

## 5. Conclusions

We put forward a maximum likelihood technique, coupled to a variational calculus, which allows the robust, non-parametric reconstruction of the evolution of the star formation rate from the information contained in colour-magnitude diagrams. A full Bayesian analysis is also applied to assess whether the best solutions found are also good fits to the data. Its main limitation, at the moment, is the prior knowledge of the metallicity of the ensemble of stars in a CMD.

For the first time, an objective reconstruction of star formation histories without any a priori or model-dependent information is applied to an homogeneous sample of dwarf spheroidals of the Local Group. We find a wide variety of SFHs, with bursts of activity uncorrelated to any special epoch or event, like perigalacticon passages. Among many other things, this also implies that late accretion was not important in the formation of the Galactic halo (Gilmore et al. 2000).

In the solar neighbourhood observed by Hipparcos, we infer –with an unprecedented resolution of 50 Myr– its star formation history over the past 3 Gyr, finding a surprising regularity of star formation episodes separated by some 0.5 Gyr. A possible explanation is that the solar neighbourhood interacted with two spiral arms or the Galactic bar, triggering star formation at each interaction. These bursts are likely to induce the formation of massive star clusters which slowly dissolve into the galactic disc.

## References

Amaral L.H., Lépine J.R.D., 1997, MNRAS, 286, 885

Aparicio A., 1998, in IAU Symp 192, The stellar content of Local Group galaxies, ed. P. Whitelock & R. Cannon (San Francisco: ASP), 20

Asiain R., Figueras F., Torra J., 1999, A&A, 350, 434

Avedisova V.S., 1989, Astrophys. 30, 83

Chereul E., Crézé M., Bienaymé O., 1999, A&A Suppl. 135, 5

Dolphin, A., 1997, New Ast., 2, 397

Gallart C., Freedman W., Aparicio A., Bertelli G., Chiosi C., 1999, AJ, 118, 2245

Gilmore G., Hernandez X., Valls-Gabaud D., 2000, The Galactic halo: from globular clusters to field stars, 35th Liège Conf., ed. A. Noels et al. in press (astro-ph/9910409)

Held, E.V. et al., 2000, ApJ, 530, L85

Hernandez X., Valls-Gabaud D., Gilmore G., 1999, MNRAS, 304, 705 (Paper I)

Hernandez X., Gilmore G., Valls-Gabaud D., 2000a, MNRAS, 317, 83 (Paper II)

Hernandez X., Valls-Gabaud D., Gilmore G., 2000b, MNRAS, 316, 605 (Paper III)

Hurley-Keller D., Mateo M., Nemec J., 1998, AJ, 115, 1840

Mateo M., 1998, ARA&A, 36, 435

Mighell K.J., 1997, AJ, 114, 1458

Mighell K.J., Burke, C.J., 1999, AJ, 118, 366

Mishurov Y.N., Pavloskaya E.D., Suchkov A.A., 1979, AZh, 56, 268

Ng, Y.K., 1998, A&AS, 132, 133

Rocha-Pinto H.J., Scalo J., Maciel W.J., Flynn C., 2000a, ApJ, 531, L115

Rocha-Pinto H.J., Maciel, W.J., Scalo J., Flynn C., 2000b, A&A, submitted (astro-ph/0001383)

Saha P., 1998, AJ, 115, 1206

Sandage, A. 1986, A&A, 161, 89

Skuljan, J., Hearnshaw, J.B., Cottrell, P.L., 1999, MNRAS, 308, 731

Tolstoy E., 1999, in Dwarf galaxies and cosmology, ed. T.X. Thuan, C. Balkowski, V. Cayatte, J. Tran Thanh Van (Paris: Editions Frontières)

Tolstoy E., Saha A., 1996, ApJ, 462, 672

Vallée, J.P., 1995, ApJ, 454, 119

Vazquez E.C., Scalo J.M., 1989, ApJ, 343, 644


text:
# Perturbative renormalization of the Δ𝐵=2 four-quark operators in lattice NRQCDposter presented by K-I. Ishikawa.

## 1 Introduction

In the lattice calculation of the $`B`$ meson $`B`$-parameters , the use of the effective theory for heavy quark is essential to obtain results with controlled systematic errors. The static approximation in which the heavy quark mass is sent to infinity has been used by many authors , and the next step is to incorporate the correction from finite $`b`$ quark mass. In this direction, we performed an exploratory lattice calculation using the NRQCD action for heavy quark, in which we found that the $`1/m_B`$ correction could be sizable for $`B_B`$  and for $`B_S`$ . A drawback in these calculations was, however, that the perturbative renormalization of the four-quark operators were not available for the NRQCD action and we used the one-loop coefficients for the static action. As a result our results contained large systematic error of $`O(\alpha _s/aM)`$, which is about 10% if evaluated with an order counting at $`\beta `$ = 5.9.

In this paper, we present a perturbative calculations of renormalization constants of the heavy-light four quark operators with the $`O(1/M)`$ lattice NRQCD action. The results connect the missing link in our previous calculations, and we present a reanalysis of $`B_B`$ and $`B_S`$. A full description of our results is published in .

## 2 Operator Definitions

The $`B`$ parameters in the $`B_q\overline{B}_q`$ are defined through

$$\overline{B}_q|O_L^{\overline{\mathrm{MS}}}(\mu )|B_q=\frac{8}{3}M_{B_q}^2f_{B_q}^2B_{L_q}^{\overline{\mathrm{MS}}}(\mu ),$$
(1)
$`\overline{B}_q|O_S^{\overline{\mathrm{MS}}}(\mu )|B_q`$
$`=`$ $`{\displaystyle \frac{5}{3}}M_{B_q}^2f_{B_q}^2B_{S_q}^{\overline{\mathrm{MS}}}(\mu )\left({\displaystyle \frac{M_{B_q}}{\overline{m}_b(\mu )+\overline{m}_q(\mu )}}\right)^2,`$
where a subscript $`q`$ denotes light quark flavor, i.e. $`d`$ or $`s`$. $`B_L`$ is usually called $`B_B`$, but we use a notation $`B_L`$ in this paper to remind that it parameterizes a matrix element of the operator $`O_L`$. In the above equations, the operators $`O_L^{\overline{\mathrm{MS}}}(\mu )`$ and $`O_S^{\overline{\mathrm{MS}}}(\mu )`$ are defined in the continuum theory with the $`\overline{\mathrm{MS}}`$ scheme. Their Dirac structure is

$`O_L`$ $`=`$ $`\overline{b}\gamma _\mu P_Lq\overline{b}\gamma _\mu P_Lq,`$ (3)
$`O_S`$ $`=`$ $`\overline{b}P_Lq\overline{b}P_Lq,`$ (4)
where $`P_L=1\gamma _5`$. The totally anti-commuting convention for $`\gamma _5`$ is assumed in the dimensional regularization.

The matching of these continuum operators is done against the lattice NRQCD. Namely we consider the $`O(1/M)`$ NRQCD action for heavy quark and the $`O(a)`$-improved Wilson fermion for light quark. For gluons the standard plaquette action is used.

For $`O_L`$ and $`O_S`$ the following heavy-light four-quark operators appear in the matching relation:

$`O_L=\left[\overline{b}\gamma _\mu P_Lq\right]\left[\overline{b}\gamma _\mu P_Lq\right],`$
$`O_R=\left[\overline{b}\gamma _\mu P_Rq\right]\left[\overline{b}\gamma _\mu P_Rq\right],`$
$`O_S=\left[\overline{b}P_Lq\right]\left[\overline{b}P_Lq\right],`$
$`O_N=2\left[\overline{b}\gamma _\mu P_Lq\right]\left[\overline{b}\gamma _\mu P_Rq\right]+4\left[\overline{b}P_Lq\right]\left[\overline{b}P_Rq\right],`$
$`O_M=2\left[\overline{b}\gamma _\mu P_Lq\right]\left[\overline{b}\gamma _\mu P_Rq\right]4\left[\overline{b}P_Lq\right]\left[\overline{b}P_Rq\right],`$
$`O_P=2\left[\overline{b}\gamma _\mu P_Lq\right]\left[\overline{b}\gamma _\mu P_Rq\right]+12\left[\overline{b}P_Lq\right]\left[\overline{b}P_Rq\right],`$
$`O_T=(2+N_c)\left[\overline{b}\gamma _\mu P_Lq\right]\left[\overline{b}\gamma _\mu P_Rq\right]`$
$`2(3N_c^22N_c4)\left[\overline{b}P_Lq\right]\left[\overline{b}P_Rq\right],`$ (5)
where $`P_R=1+\gamma _5`$. The heavy quark field $`b`$ has four Dirac components, and it is obtained from the two-component spinor non-relativistic quark (anti-quark), $`Q`$ ($`\chi `$), through the Foldy-Wouthuysen-Tani transformation,

$$b=\left[1\frac{\stackrel{}{\gamma }\stackrel{}{D}}{2M}\right]\left(\begin{array}{c}Q\\ \chi ^{}\end{array}\right).$$
(6)
In the static limit $`M\mathrm{}`$, where algebra is greatly reduced, we consider the $`O(a)`$-improvement of the lattice four-quark operators through $`O(\alpha _sa)`$. In doing that we also define the following dimension seven operators:

$`O_{LD}`$ $`=`$ $`\left[\overline{b}\gamma _\mu P_Lq\right]\left[\overline{b}\gamma _\mu P_L(a\stackrel{}{D}\stackrel{}{\gamma })q\right],`$
$`O_{SD}`$ $`=`$ $`\left[\overline{b}P_Lq\right]\left[\overline{b}P_L(a\stackrel{}{D}\stackrel{}{\gamma })q\right],`$
$`O_{ND}`$ $`=`$ $`2\left[\overline{b}\gamma _\mu P_Lq\right]\left[\overline{b}\gamma _\mu P_R(a\stackrel{}{D}\stackrel{}{\gamma })q\right]`$
$`+4\left[\overline{b}P_Lq\right]\left[\overline{b}P_R(a\stackrel{}{D}\stackrel{}{\gamma })q\right],`$
$`O_{PD}`$ $`=`$ $`2\left[\overline{b}\gamma _\mu P_Lq\right]\left[\overline{b}\gamma _\mu P_R(a\stackrel{}{D}\stackrel{}{\gamma })q\right]`$ (7)
$`+12\left[\overline{b}P_Lq\right]\left[\overline{b}P_R(a\stackrel{}{D}\stackrel{}{\gamma })q\right].`$
## 3 Matching

At one-loop level, the matching of the continuum operator $`O_X^{\overline{\mathrm{MS}}}(\mu )`$ is given by

$`O_X^{\overline{\mathrm{MS}}}(\mu )`$ $`=`$ $`\left[1+{\displaystyle \frac{\alpha _s}{4\pi }}\zeta _{X,X}(\mu ;1/a)\right]O_X^{\mathrm{lat}}(1/a)`$ (8)
$`+{\displaystyle \underset{Y=\mathrm{dim}.6}{}}{\displaystyle \frac{\alpha _s}{4\pi }}\zeta _{X,Y}(\mu ;1/a)O_Y^{\mathrm{lat}}(1/a)`$
$`+{\displaystyle \underset{Z=\mathrm{dim}.7}{}}{\displaystyle \frac{\alpha _s}{4\pi }}\zeta _{X,Z}O_Z^{\mathrm{lat}}(1/a),`$
where the summations with $`Y`$ and $`Z`$ run over dimension six and seven operators respectively. $`\mu `$ is the renormalization scale in the $`\overline{\mathrm{MS}}`$ scheme and $`a`$ is lattice spacing.

The coefficients $`\zeta `$ are determined so that the on-shell scattering amplitudes in the continuum and lattice theories agree with each other. For the external state, zero spatial momentum heavy and light quark suffice to determine the coefficients for dimension six operators. A derivative in terms of spatial momenta have to be considered to obtain the dimension seven coefficients.

### 3.1 $`O_L`$

The matching relation for $`O_L`$ is

$`O_L^{\overline{\mathrm{MS}}}(\mu )`$ $`=`$ $`\left[1+{\displaystyle \frac{\alpha _s}{4\pi }}\zeta _{L,L}(\mu ;1/a)\right]O_L^{\mathrm{lat}}(1/a)`$ (9)
$`+{\displaystyle \frac{\alpha _s}{4\pi }}\zeta _{L,S}(1/a)O_S^{\mathrm{lat}}(1/a)`$
$`+{\displaystyle \frac{\alpha _s}{4\pi }}\zeta _{L,N}(1/a)O_N^{\mathrm{lat}}(1/a)`$
$`+{\displaystyle \frac{\alpha _s}{4\pi }}\zeta _{L,M}(1/a)O_M^{\mathrm{lat}}(1/a)`$
$`+{\displaystyle \frac{\alpha _s}{4\pi }}\zeta _{L,R}(1/a)O_R^{\mathrm{lat}}(1/a)`$
$`+{\displaystyle \frac{\alpha _s}{4\pi }}\zeta _{L,LD}O_{LD}^{\mathrm{lat}}(1/a)`$
$`+{\displaystyle \frac{\alpha _s}{4\pi }}\zeta _{L,ND}O_{ND}^{\mathrm{lat}}(1/a).`$
The coefficient $`\zeta _{L,L}(\mu ;1/a)`$ of the leading operator $`O_L^{\mathrm{lat}}`$ contains logarithmic terms $`6\mathrm{ln}(a^2M_0^2)2\mathrm{ln}(a^2\mu ^2)`$, while others do not have $`\mu `$ dependence. The operator $`O_M`$ does not exist in the matching of the static operator, but it becomes necessary for NRQCD.

Figure 1 shows the heavy quark mass dependence of $`\zeta _{L,X}`$ for the dimension six operators. Filled symbols at $`1/aM_0=0`$ are calculated for the static action. We find a large slope in $`1/aM_0`$ for $`\zeta _{L,L}`$ as shown in the top panel of Figure 1. Most of the mass dependence comes from factorized diagrams, i.e. the four-quark operator can be split into two axial-vector currents and the gluon propagator does not connect these two. As a result, the coefficient $`\zeta _{L,L}`$ is much reduced if we subtract the contributions of two axial vector current, as plotted in the lower panel ($`\zeta _{L,L}2\zeta _A`$). This combination will appear in the analysis of the $`B`$ parameters, when we consider the ratio $`O_L/A_0^2`$.

The coefficients for the dimension seven operators are obtained in the static limit as $`\zeta _{L,LD}`$ = $``$17.20 and $`\zeta _{L,ND}`$ = $``$9.20 .

### 3.2 $`O_S`$

The matching of $`O_S`$ is given by

$`O_S^{\overline{\mathrm{MS}}}(\mu )`$ $`=`$ $`\left[1+{\displaystyle \frac{\alpha _s}{4\pi }}\zeta _{S,S}(\mu ;1/a)\right]O_S^{\mathrm{lat}}(1/a)`$ (10)
$`+{\displaystyle \frac{\alpha _s}{4\pi }}\zeta _{S,L}(\mu ;1/a)O_L^{\mathrm{lat}}(1/a)`$
$`+{\displaystyle \frac{\alpha _s}{4\pi }}\zeta _{S,P}(1/a)O_P^{\mathrm{lat}}(1/a)`$
$`+{\displaystyle \frac{\alpha _s}{4\pi }}\zeta _{S,R}(1/a)O_R^{\mathrm{lat}}(1/a)`$
$`+{\displaystyle \frac{\alpha _s}{4\pi }}\zeta _{S,T}(1/a)O_T^{\mathrm{lat}}(1/a)`$
$`+{\displaystyle \frac{\alpha _s}{4\pi }}\zeta _{S,SD}O_{SD}^{\mathrm{lat}}(1/a)`$
$`+{\displaystyle \frac{\alpha _s}{4\pi }}\zeta _{S,LD}O_{LD}^{\mathrm{lat}}(1/a)`$
$`+{\displaystyle \frac{\alpha _s}{4\pi }}\zeta _{S,PD}O_{PD}^{\mathrm{lat}}(1/a),`$
where $`O_T`$ is new for the NRQCD action. The coefficients $`\zeta _{S,X}`$ are plotted in Figure 2. Like $`\zeta _{L,L}`$ in the matching of $`O_L`$, the leading operator $`O_S`$ has large coefficient, which is reduced by subtracting the factorized contribution, as shown in the upper panel.

The coefficients for the dimension seven operators are $`\zeta _{S,SD}`$ = $``$6.88, $`\zeta _{S,LD}`$ = 2.58 and $`\zeta _{S,PD}`$ = 1.15 in the static limit.

## 4 Reanalysis of $`B`$ Parameters

With these perturbative coefficients, we reanalyze our previous lattice calculation for the $`B`$ parameters, $`B_L`$ and $`B_S`$ .

### 4.1 $`B_L`$

The $`B`$ parameter $`B_L`$ is calculated through

$$B_L^{\overline{\mathrm{MS}}}(\mu )=\underset{X}{}Z_{L,X/A^2}(\mu ;1/a)B_X^{\mathrm{lat}}(1/a),$$
(11)
where $`X`$ runs over $`L`$, $`S`$, $`N`$, $`M`$ and $`R`$. The renormalization constants are

$`Z_{L,L/A^2}(\mu ;1/a)`$ $`=`$
$`1+{\displaystyle \frac{\alpha _s}{4\pi }}\left(\zeta _{L,L}(\mu ;1/a)2\zeta _A(1/a)\right),`$
$`Z_{L,S/A^2}(1/a)`$ $`=`$ $`{\displaystyle \frac{\alpha _s}{4\pi }}\zeta _{L,S}(1/a),`$
$`Z_{L,N/A^2}(1/a)`$ $`=`$ $`{\displaystyle \frac{\alpha _s}{4\pi }}\zeta _{L,N}(1/a),`$
$`Z_{L,M/A^2}(1/a)`$ $`=`$ $`{\displaystyle \frac{\alpha _s}{4\pi }}\zeta _{L,M}(1/a),`$
$`Z_{L,R/A^2}(1/a)`$ $`=`$ $`{\displaystyle \frac{\alpha _s}{4\pi }}\zeta _{L,R}(1/a),`$
and the lattice $`B`$ parameters $`B_{X_q}^{\mathrm{lat}}(1/a)`$ were measured in the numerical simulations through

$`B_{X_q}^{\mathrm{lat}}(1/a)={\displaystyle \frac{\overline{B}_q|O_X^{\mathrm{lat}}(1/a)|B_q}{\frac{8}{3}\overline{B}_q|A_4^{\mathrm{lat}}(1/a)|00|A_4^{\mathrm{lat}}(1/a)|B_q}},`$
in our previous paper . Since the dimension seven operators are not included in the calculation, we could not remove the error of $`O(\alpha _sa\mathrm{\Lambda }_{\mathrm{QCD}})`$.

Figure 3 shows the individual contributions $`Z_{L,X/A^2}(\mu ;1/a)B_X^{\mathrm{lat}}(1/a)`$ for the $`B`$ parameter $`B_L`$. The $`V`$-scheme coupling constant $`\alpha _V(q^{})`$ is used with $`q^{}=2/a`$. Although there is a sizable mass dependence in $`Z_{L,X/A^2}(\mu ;1/a)`$ as we can see from the plot of coefficients (Figure 1), it is canceled by the $`1/m_B`$ corrections in the matrix elements on the lattice $`B_X^{\mathrm{lat}}(1/a)`$, and there is little mass dependence in the products shown in Figure 3.

Our result for $`B_L^{\overline{\mathrm{MS}}}(\mu )`$ is plotted in Figure 4 as a function of $`1/M_P`$ with filled circles. The statistical error shown by solid error bars is very small except for the heaviest heavy quark, while the systematic error given by dotted error bars is more significant. We estimate the systematic errors with an order counting, in which the systematic errors of $`O(\alpha _s^2)`$, $`O(\alpha _sa\mathrm{\Lambda }_{\mathrm{QCD}})`$, $`O(\alpha _s\mathrm{\Lambda }_{\mathrm{QCD}}/M)`$, $`O((a\mathrm{\Lambda }_{\mathrm{QCD}})^2)`$, and $`O((\mathrm{\Lambda }_{\mathrm{QCD}}/M)^2)`$ are considered and added in quadrature. For comparison, we also show a result with one-loop coefficients calculated for the static heavy quark by open circles, for which a large systematic error of $`O(\alpha _s/(aM))`$ is added. The effect of the NRQCD renormalization constants is sizable, but still it is within the estimated systematic error.

The numerical results for $`B_L`$ at physical $`B`$ meson mass are

$`B_{L_d}^{\overline{\mathrm{MS}}}(m_b)`$ $`=`$ $`0.85(3)(11),`$ (13)
$`B_{L_s}^{\overline{\mathrm{MS}}}(m_b)`$ $`=`$ $`0.87(2)(11),`$ (14)
where the first error is statistical and the second is systematic one ($`13\%`$). The SU(3) breaking effect is

$$\frac{B_{L_s}^{\overline{\mathrm{MS}}}(m_b)}{B_{L_d}^{\overline{\mathrm{MS}}}(m_b)}=1.01(1),$$
(15)
where only the statistical error is quoted.

### 4.2 $`B_S`$

Using Eq. (10), $`B_S`$ is obtained as

$$\frac{B_{S_q}^{\overline{\mathrm{MS}}}(\mu )}{_q(\mu )^2}=\underset{X}{}Z_{S,X/A^2}(\mu ;1/a)B_{X_q}^{{}_{}{}^{}\mathrm{lat}}(1/a),$$
(16)
where $`X`$ runs over $`S`$, $`L`$, $`R`$, $`P`$ and $`T`$. In the LHS, we define a ratio of heavy-light matrix elements

$$_q(\mu )=\left|\frac{0|A_4^{\overline{\mathrm{MS}}}|B_q}{0|P^{\overline{\mathrm{MS}}}(\mu )|B_q}\right|.$$
(17)
The renormalization constants are

$`Z_{S,S/A^2}(\mu ;1/a)`$ $`=`$
$`1+{\displaystyle \frac{\alpha _s}{4\pi }}\left(\zeta _{S,S}(\mu ;1/a)2\zeta _A(1/a)\right),`$
$`Z_{S,L/A^2}(\mu ;1/a)`$ $`=`$ $`{\displaystyle \frac{\alpha _s}{4\pi }}\zeta _{S,L}(1/a),`$
$`Z_{S,R/A^2}(1/a)`$ $`=`$ $`{\displaystyle \frac{\alpha _s}{4\pi }}\zeta _{S,R}(1/a),`$
$`Z_{S,P/A^2}(1/a)`$ $`=`$ $`{\displaystyle \frac{\alpha _s}{4\pi }}\zeta _{S,P}(1/a),`$
$`Z_{S,T/A^2}(1/a)`$ $`=`$ $`{\displaystyle \frac{\alpha _s}{4\pi }}\zeta _{S,T}(1/a).`$
For $`B_S`$ we define the lattice $`B`$ parameters as

$`B_{X_q}^{{}_{}{}^{}\mathrm{lat}}(1/a)=`$ (18)
$`{\displaystyle \frac{\overline{B}_q|O_X^{\mathrm{lat}}(1/a)|B_q}{\frac{5}{3}\overline{B}_q|A_4^{\mathrm{lat}}(1/a)|00|A_4^{\mathrm{lat}}(1/a)|B_q}}.`$
The mass dependence of $`B_{S_s}^{\overline{\mathrm{MS}}}/_s`$ is plotted in Figure 4 (filled squares). Open squares are our previous results calculated with the static renormalization constants. The mass dependence is less significant in the results with the correct (NRQCD) renormalization constants.

The numerical results are

$`{\displaystyle \frac{B_{S_d}^{\overline{\mathrm{MS}}}(m_b)}{_d(m_b)^2}}`$ $`=`$ $`1.24(3)(16),`$ (19)
$`{\displaystyle \frac{B_{S_s}^{\overline{\mathrm{MS}}}(m_b)}{_s(m_b)^2}}`$ $`=`$ $`1.24(3)(16),`$ (20)
and

$$\frac{B_{S_s}^{\overline{\mathrm{MS}}}(m_b)_d(m_b)^2}{B_{S_d}^{\overline{\mathrm{MS}}}(m_b)_s(m_b)^2}=1.003(4).$$
(21)
for the SU(3) breaking ratio.

## Acknowledgments

S.H. and T.O. are supported in part by the Grants-in-Aid of the Ministry of Education (Nos. 11740162, 12640279). K-I.I. and N.Y. are supported by the JSPS Research Fellowship.