text:
# 1 Introduction

## 1 Introduction

Fragmentation is an irreversible kinetic phenomenon which occurs in many physical and chemical processes. Because of broad range of applications, many recent studies have been carried out to investigate the kinetics of the processes by introducing simple fragmentation models . In  Ziff and MaGrady presented a model of fragmentation in which the rate of break up depends on the size of the fragments. A general discussion of the kinetics of continuous fragmentation processes was given by Cheng and Redner in . In  Krapivsky and Ben-Naim studied the kinetics of random fragmentation of multidimensional objects. In  three models of binary fragmentation were investigated analytically in which at each time step, the largest fragment in the system is broken with some externally tuneable probability. In the second model of that paper, at each time step, there is a fixed rate of fragmentation for the largest fragment and a rate proportional to the inverse of the fragment size for all other fragments smaller than the largest one. The model was solved exactly in the long time limit to reveal a power law distribution with an exponent, which depends on the precise details of the fragmentation process.

Surprisingly, there is a similarity between the result of the stable distribution from the second model of , with some experimental and analytical results of fragment size and mass distributions with power law forms in the shock fragmentation . In  a generalization of the second model of  was carried out. Some models of binary fragmentation were introduced which revealed composite power law distributions for the mass and size of the fragments .

Another behaviour of distribution which is observed in many different processes has a log-normal form . In  a comprehensive study in the log-normal distribution is given. The application of these distributions in many different areas of science is discussed in , such as economics, biology, astronomy, philology, small particle statistics and physical and industrial processes. The log-normal behaviour is known to be able to describe the size distribution well in a wide variety of geological situations , such as that of rocks in a boulder field. In  an analytical and numerical study was carried out on the fragmentation of long thin glass rods. Some stochastic models for one-dimensional brittle fracture were introduced. The models successively extended to describe cascade fracture of long thin glass rods. They are found to give the fragment size distribution close to a log-normal one. In terms of rupture in liquids, experiments in rupture of mercury drops were performed . The result of size distribution for the cumulative number of droplets showed a log-normal distribution at the small falling heights. A transition of the distribution from log-normal to scaling behaviour was observed as the falling height was increased. In  an experimental work on the shock fragmentation of long thin glass rods was carried out by Ishii and Matsushita. The results of fragment size and mass distributions at small falling heights showed a log-normal form for larger fragments and a power law form for smaller fragments. The crossover was seen to be at length scales around the rod diameter. The mass and size distributions for the larger fragments showed a power law form as the falling height was increased.

This paper is organised as follows. In sections 2 and 3 we introduce two new models, which are generalisations of the second model of  and the models of fragmentation in . In the models we introduce a time dependent transition size which produces two regions of fragment sizes with different rates of fragmentation above and below the transition size. A rate of fragmentation proportional to the inverse of the fragment size in the smaller size region produces a power law distribution in that region. A rate of fragmentation combined of two terms, one proportional to the inverse of the fragment size and the other proportional to a logarithmic function of the fragment size yield a log-normal distribution in the larger size region. Finally in section 4 we summarise our conclusions and discuss the application of our results to experimental work on the shock fragmentation of long thin glass rods  and the rupture of mercury drops .

For the rest of this introduction we give a very brief summary of previous work on models of binary fragmentation. In these models, the density of fragments of size $`y`$ at time $`t`$, $`n(y,t)`$, evolves according to

$$n(y,t+\delta t)=n(y,t)\delta tn(y,t)_0^yR(z,yz)𝑑z+\mathrm{\hspace{0.17em}2}\delta t_y^{L(t)}R(y,zy)n(z,t)𝑑z,$$
(1)
where $`R(y,z)`$ is the intrinsic rate that a fragment of size $`(y+z)`$ breaks into fragments of size $`y`$ and $`z`$. The first term on the right hand side is the contribution from those fragments not chosen for the fragmentation in time $`(t,t+\delta t)`$. The second term represents the decrease in the number of fragments of size $`y`$ by the fragmentation into fragments of size $`z`$ and $`(yz)(y>z)`$. The third term represents the increase in the number of particles of size $`y`$ due to the fragmentation of particles of size $`z(>y)`$, such that one of the products is of size $`y`$.

The upper limit of the second integral on the right-hand side is usually set to some fixed value greater than the size of the largest particle in the system. This value is usually $`1`$ or $`\mathrm{}`$ depending on the details of the model. To connect with the analysis of our models in the following sections, we have set the upper limit to $`L(t)`$, the size of the largest particle at time $`t`$. As $`n(y,t)=0`$ for all $`y>L(t)`$, this choice has no effect on our results.

We now turn to a detailed investigation of our models.

## 2 Model A

Firstly, we introduce a time dependent transition size, $`y_m(t)`$, between zero and the size of the largest particle in the system at time $`t`$, $`L(t)`$. The transition size produces two regions of fragment sizes with different rates of fragmentation above and below the transition size. At each time step there is a fixed probability, $`p_1`$, for the fragmentation of the largest particle in the system. A particle of size smaller than the transition size breaks with a probability $`p_2`$ and a rate proportional to the inverse of the fragment size, at each time step. Any other particle of size larger than the transition size can be fragmented, at each time step, either with a probability $`p_3`$ and a rate proportional to the inverse of the fragment size, or with a probability $`p_4`$ and a rate proportional to $`{\displaystyle \frac{\delta t}{y(t)}}\mathrm{ln}\left({\displaystyle \frac{y(t)}{\stackrel{~}{y}(t)}}\right)`$. Here $`y(t)`$ is the size of the fragment at time $`t`$ and $`\stackrel{~}{y}(t)`$ is a size for which the logarithmic rate is zero at time $`t`$. In the model we have $`p_1+p_2+p_3+p_4=1`$. We also assume that $`y_m(t)\stackrel{~}{y}(t)`$, such that for all fragments of sizes larger than the transition size, the logarithmic term in the rate of fragmentation is positive. In order to make the model tractable, we will later choose to have $`y_m(t)`$ proportional to $`L(t)`$.

In the model the distribution function of daughter fragments is uniform. Therefore the rate of breaking of a particle into two smaller pieces is independent of the daughter fragment sizes. With probability $`p_2`$ every particle of size smaller than $`y_m(t)`$ is equally likely to be chosen for the fragmentation. With probabilities $`p_3`$ and $`p_4`$, every particle of size larger than the transition size is equally likely to be selected for the fragmentation with the correspondence rate. We assume that in time $`(t,t+\delta t)`$ the largest size changes from $`L(t)`$ to $`L(t)p_1\delta L`$ and the transition size changes from $`y_m(t)`$ to $`y_m(t)\delta y_m`$.

Then, the fragmentation rate for this model is given by

$`\delta tR(y,z)`$ $`=`$ $`p_1{\displaystyle \frac{\delta L}{L(t)}}\delta (y+zL(t))+p_2{\displaystyle \frac{\delta t}{y+z}}\theta (y_m(t)(y+z))`$
$`+\left[p_3{\displaystyle \frac{\delta t}{y+z}}+p_4{\displaystyle \frac{\delta t}{y+z}}\mathrm{ln}\left({\displaystyle \frac{y+z}{\stackrel{~}{y}(t)}}\right)\right]\theta ((y+z)y_m(t)).`$
The $`\delta `$-function in the first term on the right-hand side ensures that only the largest particle is fragmented and $`{\displaystyle \frac{\delta L}{L(t)}}`$ is the probability of placing the cut in the largest particle in a particular infinitesimal length $`\delta L`$. The second term on the right-hand side of (2) represents the rate of fragmentation for the particles of sizes smaller than $`y_m(t)`$ in time $`(t,t+\delta t)`$. The remaining terms on the right hand side of (2) represent the rate of fragmentation for the fragments of sizes between $`y_m(t)`$ and $`L(t)`$ in time $`(t,t+\delta t)`$. Inserting the rate (2) into the kinetic equation (1) for the fragments in both size regions smaller and larger than $`y_m(t)`$, yields

$`n_1(y,t+\delta t)`$ $`=`$ $`n_1(y,t)p_2\delta tn_1(y,t)+\mathrm{\hspace{0.33em}2}p_1{\displaystyle \frac{\delta L}{L(t)}}n_2(L(t),t)`$ (3)
$`+\mathrm{\hspace{0.17em}2}p_2\delta t{\displaystyle _y^{y_m(t)}}n_1(z,t){\displaystyle \frac{dz}{z}}+2p_3\delta t{\displaystyle _{y_m(t)}^{L(t)}}n_2(z,t){\displaystyle \frac{dz}{z}}`$
$`+\mathrm{\hspace{0.17em}2}p_4\delta t{\displaystyle _{y_m(t)}^{L(t)}}n_2(z,t)\mathrm{ln}\left({\displaystyle \frac{z}{\stackrel{~}{y}(t)}}\right){\displaystyle \frac{dz}{z}}`$
for $`y<y_m(t)`$ and

$`n_2(y,t+\delta t)`$ $`=`$ $`n_2(y,t)p_3\delta tn_2(y,t)p_4\delta t\mathrm{ln}\left({\displaystyle \frac{y}{\stackrel{~}{y}(t)}}\right)n_2(y,t)`$ (4)
$`+\mathrm{\hspace{0.25em}2}p_1{\displaystyle \frac{\delta L}{L(t)}}n_2(L(t),t)+\mathrm{\hspace{0.33em}2}p_3\delta t{\displaystyle _y^{L(t)}}n_2(z,t){\displaystyle \frac{dz}{z}}`$
$`+\mathrm{\hspace{0.25em}2}p_4\delta t{\displaystyle _y^{L(t)}}n_2(z,t)\mathrm{ln}\left({\displaystyle \frac{z}{\stackrel{~}{y}(t)}}\right){\displaystyle \frac{dz}{z}}`$
for $`y>y_m(t)`$. Here $`n_1(y,t)`$ and $`n_2(y,t)`$ are the densities of fragments of size $`y`$ at time $`t`$ which is smaller and larger than $`y_m(t)`$, respectively. The second term on the right hand side of (3) and the second and third terms on the right hand side of (4) represent the decrease in the number of particles of size $`y`$ from fragmentation into smaller particles, with probabilities $`p_2`$, $`p_3`$ and $`p_4`$, respectively. The next term on the right hand sides of (3,4) describes the gain in the number of particles of size $`y`$ from the fragmentation of the largest fragment of size $`L(t)`$. The remaining terms on the right hand sides of (3,4) are the increase in the number of particles of size $`y`$ from the fragmentation of all particles larger than $`y`$ and smaller than $`L(t)`$.

Now we define densities of fragments of length $`y`$ at time $`t`$, which are normalised to a positive constant value A for two regions of sizes smaller and larger than $`y_m(t)`$,

$$g_1(y,t)=A\frac{n_1(y,t)}{{\displaystyle _0^{y_m(t)}}n_1(z,t)𝑑z+{\displaystyle _{y_m(t)}^{L(t)}}n_2(z,t)𝑑z}$$
(5)
for $`y<y_m(t)`$ and

$$g_2(y,t)=A\frac{n_2(y,t)}{{\displaystyle _0^{y_m(t)}}n_1(z,t)𝑑z+{\displaystyle _{y_m(t)}^{L(t)}}n_2(z,t)𝑑z}$$
(6)
for $`y>y_m(t)`$. In time $`tt+\delta t`$ as $`L(t)L(t)p_1\delta L`$ and $`y_m(t)y_m(t)\delta y_m`$ they yield

$$g_1(y,t+\delta t)=A\frac{n_1(y,t+\delta t)}{{\displaystyle _0^{y_m(t)\delta y_m(t)}}n_1(z,t+\delta t)𝑑z+{\displaystyle _{y_m(t)\delta y_m(t)}^{L(t)p_1\delta L}}n_2(z,t+\delta t)𝑑z}$$
(7)
for $`y<y_m(t)`$ and

$$g_2(y,t+\delta t)=A\frac{n_2(y,t+\delta t)}{{\displaystyle _0^{y_m(t)\delta y_m(t)}}n_1(z,t+\delta t)𝑑z+{\displaystyle _{y_m(t)\delta y_m(t)}^{L(t)p_1\delta L}}n_2(z,t+\delta t)𝑑z}$$
(8)
for $`y>y_m(t)`$. Inserting (3), (4) in (7), (8), respectively, gives us

$`g_1(y,t+\delta t)=`$ (9)
$`2p_1{\displaystyle \frac{\delta L}{L(t)}}g_2(L(t),t)+\mathrm{\hspace{0.17em}2}p_2\delta t{\displaystyle _y^{y_m(t)}}g_1(z,t){\displaystyle \frac{dz}{z}}`$
$`+\mathrm{\hspace{0.33em}2}p_3\delta t{\displaystyle _{y_m(t)}^{L(t)}}g_2(z,t){\displaystyle \frac{dz}{z}}+\mathrm{\hspace{0.17em}2}p_4\delta t{\displaystyle _{y_m(t)}^{L(t)}}g_2(z,t)\mathrm{ln}\left({\displaystyle \frac{z}{\stackrel{~}{y}(t)}}\right){\displaystyle \frac{dz}{z}}`$
$`+g_1(y,t)[1p_2\delta t{\displaystyle \frac{p_1}{A}}g_2(L(t),t)\delta L{\displaystyle \frac{p_2}{A}}\delta t{\displaystyle _0^{y_m(t)}}g_1(z,t)dz`$
$`{\displaystyle \frac{p_3}{A}}\delta t{\displaystyle _{y_m(t)}^{L(t)}}g_2(z,t)dz{\displaystyle \frac{p_4}{A}}\delta t{\displaystyle _{y_m(t)}^{L(t)}}g_2(z,t)\mathrm{ln}\left({\displaystyle \frac{z}{\stackrel{~}{y}(t)}}\right)dz]`$
for $`y<y_m(t)`$ and

$`g_2(y,t+\delta t)=`$ (10)
$`2p_1{\displaystyle \frac{\delta L}{L(t)}}g_2(L(t),t)+\mathrm{\hspace{0.17em}2}p_3\delta t{\displaystyle _y^{L(t)}}g_2(z,t){\displaystyle \frac{dz}{z}}`$
$`+\mathrm{\hspace{0.33em}2}p_4\delta t{\displaystyle _y^{L(t)}}g_2(z,t)\mathrm{ln}\left({\displaystyle \frac{z}{\stackrel{~}{y}(t)}}\right){\displaystyle \frac{dz}{z}}+g_2(y,t)[1p_3\delta t`$
$`p_4\delta t\mathrm{ln}\left({\displaystyle \frac{y}{\stackrel{~}{y}(t)}}\right){\displaystyle \frac{p_1}{A}}g_2(L(t),t)\delta L{\displaystyle \frac{p_2}{A}}\delta t{\displaystyle _0^{y_m(t)}}g_1(z,t)𝑑z`$
$`{\displaystyle \frac{p_3}{A}}\delta t{\displaystyle _{y_m(t)}^{L(t)}}g_2(z,t)dz{\displaystyle \frac{p_4}{A}}\delta t{\displaystyle _{y_m(t)}^{L(t)}}g_2(z,t)\mathrm{ln}\left({\displaystyle \frac{z}{\stackrel{~}{y}(t)}}\right)dz]`$
for $`y>y_m(t)`$. In order to solve (9) and (10) we introduce two functions $`F_1(x,t)`$ and $`F_2(x,t)`$ defined by

$$F_1(x,t)=L(t)g_1(xL(t),t)$$
(11)
for $`x<x_m`$ and

$$F_2(x,t)=L(t)g_2(xL(t),t)$$
(12)
for $`x>x_m`$. Here we have defined dimensionless and time independent variables $`x={\displaystyle \frac{y(t)}{L(t)}}`$, $`x_m={\displaystyle \frac{y_m(t)}{L(t)}}`$ and $`\stackrel{~}{x}={\displaystyle \frac{\stackrel{~}{y}(t)}{L(t)}}`$, where $`x`$ changes in the range $`[0,1]`$. $`x_m`$ and $`\stackrel{~}{x}`$ are fixed values and restricted to the same range. Using (9-12) we can obtain two partial differential equations for $`F_1(x,t)`$ and $`F_2(x,t)`$ as

$`{\displaystyle \frac{F_1(x,t)}{t}}=`$
$`\nu x{\displaystyle \frac{F_1(x,t)}{x}}+2\nu F_2(1,t)+\mathrm{\hspace{0.17em}2}p_2{\displaystyle _x^{x_m}}F_1(z,t){\displaystyle \frac{dz}{z}}+2p_3{\displaystyle _{x_m}^1}F_2(x,t){\displaystyle \frac{dx}{x}}`$
$`+\mathrm{\hspace{0.33em}2}p_4{\displaystyle _{x_m}^1}F_2(x,t)\mathrm{ln}\left({\displaystyle \frac{x}{\stackrel{~}{x}}}\right){\displaystyle \frac{dx}{x}}F_1(x,t)[\nu +p_2+{\displaystyle \frac{\nu }{A}}F_2(1,t)`$
$`+{\displaystyle \frac{p_2}{A}}{\displaystyle _0^{x_m}}F_1(x,t)dx+{\displaystyle \frac{p_3}{A}}{\displaystyle _{x_m}^1}F_2(x,t)dx+{\displaystyle \frac{p_4}{A}}{\displaystyle _{x_m}^1}F_2(x,t)\mathrm{ln}\left({\displaystyle \frac{x}{\stackrel{~}{x}}}\right)dx]`$
for $`x<x_m`$ and

$`{\displaystyle \frac{F_2(x,t)}{t}}=`$ (14)
$`\nu x{\displaystyle \frac{F_2(x,t)}{x}}+2\nu F_2(1,t)+\mathrm{\hspace{0.17em}2}p_3{\displaystyle _x^1}F_2(z,t){\displaystyle \frac{dz}{z}}`$
$`+\mathrm{\hspace{0.33em}2}p_4{\displaystyle _x^1}F_2(z,t)\mathrm{ln}\left({\displaystyle \frac{z}{\stackrel{~}{x}}}\right){\displaystyle \frac{dz}{z}}F_2(x,t)[\nu +p_3`$
$`+p_4\mathrm{ln}\left({\displaystyle \frac{x}{\stackrel{~}{x}}}\right)+{\displaystyle \frac{\nu }{A}}F_2(1,t)+{\displaystyle \frac{p_2}{A}}{\displaystyle _0^{x_m}}F_1(x,t)𝑑x`$
$`+{\displaystyle \frac{p_3}{A}}{\displaystyle _{x_m}^1}F_2(x,t)dx+{\displaystyle \frac{p_4}{A}}{\displaystyle _{x_m}^1}F_2(x,t)\mathrm{ln}\left({\displaystyle \frac{x}{\stackrel{~}{x}}}\right)dx]`$
for $`x>x_m`$. Here we have set $`p_1{\displaystyle \frac{\delta L}{L(t)\delta t}}=\nu `$ to be a fixed value, and assumed that $`p_1>0`$. $`\nu `$ is equal to the rate of fragmentation for the largest fragment, and determines the relationship between the real time and the length of the largest fragment at time $`t`$, $`L(t)`$. In the limit $`t\mathrm{}`$, we assume that $`F_1(x,t)`$ and $`F_2(x,t)`$ evolve to time-independent quantities,

$$F_1(x)=lim_t\mathrm{}F_1(x,t)$$
(15)
for $`x<x_m`$ and

$$F_2(x)=lim_t\mathrm{}F_2(x,t)$$
(16)
for $`x>x_m`$. So that, as $`t\mathrm{}`$ then $`{\displaystyle \frac{F_1(x,t)}{t}}0`$ and $`{\displaystyle \frac{F_2(x,t)}{t}}0`$. To solve (2,14) in the long time limit using (15,16), we anticipate solutions as

$$F_1(x)=Bx^\alpha $$
(17)
for $`x<x_m`$ and

$$F_2(x)=C\frac{\mathrm{exp}\left\{\beta \left[\mathrm{ln}\left(\frac{x}{\stackrel{~}{x}}\right)\right]^2\right\}}{x}$$
(18)
for $`x>x_m`$. Substituting (17,18) in (2,14) in the long time limit, using (15,16) give us

$$\alpha =\frac{p_2}{p_3},$$
(19)
$$\beta =\frac{p_4}{2p_3}$$
(20)
and

$$\nu =p_3.$$
(21)
B and C in (17),(18), respectively, are constant values and can be obtained using (2,14) in the long time limit and the continuity condition at $`x_m`$,

$$F_1(x_m)=F_2(x_m).$$
(22)
These give

$$C=\frac{A}{\left({\displaystyle \frac{p_3}{p_3p_2}}\right)\mathrm{exp}\left\{\beta \left[\mathrm{ln}\left({\displaystyle \frac{x_m}{\stackrel{~}{x}}}\right)\right]^2\right\}+\sqrt{{\displaystyle \frac{2\pi p_3}{p_4}}}\left[P\left(\sqrt{{\displaystyle \frac{p_4}{p_3}}}\mathrm{ln}\left({\displaystyle \frac{1}{\stackrel{~}{x}}}\right)\right)P\left(\sqrt{{\displaystyle \frac{p_4}{p_3}}}\mathrm{ln}\left({\displaystyle \frac{x_m}{\stackrel{~}{x}}}\right)\right)\right]}$$
(23)
and

$$B=C\frac{\mathrm{exp}\left\{{\displaystyle \frac{p_4}{2p_3}}\left[\mathrm{ln}\left({\displaystyle \frac{x_m}{\stackrel{~}{x}}}\right)\right]^2\right\}}{x_m^{\left(1\frac{p_2}{p_3}\right)}},$$
(24)
where the second and third terms in the denominator of (23) are a proportion of the normal probability distribution function defined as

$$P(x)=\frac{1}{\sqrt{2\pi }}_{\mathrm{}}^xe^{t^2/2}𝑑t.$$
(25)
We have also assumed that

$$0\alpha =\frac{p_2}{p_3}<1.$$
(26)
Solutions (17,18) using (19,20,23-26) satisfy the normalisation relation

$$_0^{x_m}F_1(x)𝑑x+_{x_m}^1F_2(x)𝑑x=A.$$
(27)
If $`A=1`$, the normalisation is equal to $`1`$. For a choice of $`A`$ as

$$A=\sqrt{\frac{p_4}{2\pi p_3}}\frac{p_3}{p_3p_2}\mathrm{exp}\left\{\beta \left[\mathrm{ln}\left(\frac{x_m}{\stackrel{~}{x}}\right)\right]^2\right\}+P\left(\sqrt{\frac{p_4}{p_3}}\mathrm{ln}\left(\frac{1}{\stackrel{~}{x}}\right)\right)P\left(\sqrt{\frac{p_4}{p_3}}\mathrm{ln}\left(\frac{x_m}{\stackrel{~}{x}}\right)\right),$$
(28)
the solutions (17,18) using (19,20,23,24) yield

$$F_1(x)=\sqrt{\frac{p_4}{2\pi p_3}}\frac{\mathrm{exp}\left\{{\displaystyle \frac{p_4}{2p_3}}\left[\mathrm{ln}\left({\displaystyle \frac{x_m}{\stackrel{~}{x}}}\right)\right]^2\right\}}{x_m^{\left(1\frac{p_2}{p_3}\right)}}x^{\frac{p_2}{p_3}}$$
(29)
for $`x<x_m`$ and

$$F_2(x)=\sqrt{\frac{p_4}{2\pi p_3}}\frac{\mathrm{exp}\left\{{\displaystyle \frac{p_4}{2p_3}}\left[\mathrm{ln}\left({\displaystyle \frac{x}{\stackrel{~}{x}}}\right)\right]^2\right\}}{x}$$
(30)
for $`x>x_m`$. We see that the stable distribution has a power law behaviour in the smaller region, $`x<x_m`$. For a log-normal distribution, the probability of finding a fragment of size between $`x`$ and $`x+dx`$ is given by $`n(x)dx`$, where

$$n(x)=\frac{\mathrm{exp}\left\{\left[\mathrm{ln}\left(x/\stackrel{~}{x}\right)\right]^2/2\sigma ^2\right\}}{\left(2\pi \sigma ^2\right)^{1/2}x},$$
(31)
where $`\stackrel{~}{x}`$ is the mean and $`\sigma `$ is the dispersion of the distribution. Comparing (30) with (31) shows that the stable distribution $`F_2(x)`$ in the larger region exhibits a log-normal behaviour with a mean and a dispersion of the distribution equal to $`\stackrel{~}{x}`$ and $`\sigma =\sqrt{{\displaystyle \frac{p_3}{p_4}}}`$, respectively. To have a positive value for $`A`$ in (28) and therefore physically meaningful distributions, we require $`p_3>p_2`$. This requirement is satisfied in (26).

Now we consider some special cases of the model.

1. There is no fragmentation for the fragments smaller than $`y_m(t)`$, i.e. $`p_2=0`$ and $`p_1+p_3+p_4=1`$. The stable distribution then from (29) yields

$$F_1(x)=\sqrt{\frac{p_4}{2\pi p_3}}\frac{\mathrm{exp}\left\{{\displaystyle \frac{p_4}{2p_3}}\left[\mathrm{ln}\left({\displaystyle \frac{x_m}{\stackrel{~}{x}}}\right)\right]^2\right\}}{x_m}$$
(32)
for $`x<x_m`$ and is the same as (30) for $`x>x_m`$. We see that the distribution is a constant value in the smaller region and has a log-normal form in the larger region.
2. As $`y_m(t)=\stackrel{~}{y}(t)`$ and therefore $`x_m=\stackrel{~}{x}`$, the distributions (29,30) give

$$F_1(x)=\sqrt{\frac{p_4}{2\pi p_3}}x_m^{\left(1\frac{p_2}{p_3}\right)}x^{\frac{p_2}{p_3}}$$
(33)
for $`x<x_m`$ and

$$F_2(x)=\sqrt{\frac{p_4}{2\pi p_3}}\frac{\mathrm{exp}\left\{{\displaystyle \frac{p_4}{2p_3}}\left[\mathrm{ln}\left({\displaystyle \frac{x}{x_m}}\right)\right]^2\right\}}{x}$$
(34)
for $`x>x_m`$. Again the distribution reveals a power law form in the smaller region and has a log-normal behaviour in the larger region.
3. As $`p_40`$ and therefore $`\sigma ^2={\displaystyle \frac{p_3}{p_4}}\mathrm{}`$, the distribution (30) in the larger region yield $`F_2(x)x^{\mathrm{\hspace{0.17em}1}}`$. On the other hand as $`p_40`$ the rate of fragmentation (2) for the fragments of sizes larger than $`y_m(t)`$ is only proportional to the inverse of the fragment size. Consequently, at each time step either the largest fragment breaks with probability $`p_1`$, or a fragment smaller than $`y_m(t)`$ is chosen for fragmentation with probability $`p_2`$, or a fragment larger than $`y_m(t)`$ is selected with probability $`p_3`$, where $`p_1+p_2+p_3=1`$. The stable distributions then exhibit a composite power law behaviour with power laws of exponents equal to $`{\displaystyle \frac{p_2}{\nu }}`$ and $`{\displaystyle \frac{p_3}{\nu }}`$ in the smaller and larger region, respectively. In our case since from (21) we have $`\nu =p_3`$, the exponents reduce to $`{\displaystyle \frac{p_2}{p_3}}`$ and $`1`$ in the two regions. The power law distribution in the smaller region is similar to that from (17,19) and in the larger region has an exponent of $`1`$ as expected. The latter is consistent with one of interesting properties of the log-normal distribution where as the variance $`\sigma ^2`$ is pretty large, then the distribution can be approximated by a power law form with an exponent equal to $`1`$, .

## 3 Model B

In this model again a time dependent transition size, $`y_m(t)`$, produces two regions of fragment sizes with different rates of fragmentation above and below the transition size. At each time step either the largest fragment in the system breaks with a fixed probability, $`p_1`$, or a fragment smaller than the transition size is chosen for the fragmentation with probability $`p_2`$ and a rate proportional to the inverse of fragment size, or another fragment larger than $`y_m(t)`$ is selected with probability $`(1p_1p_2)`$ and a rate proportional to $`{\displaystyle \frac{\delta t}{y(t)}}\left[1+\lambda \mathrm{ln}\left({\displaystyle \frac{y(t)}{\stackrel{~}{y}(t)}}\right)\right]`$. Here $`y(t)`$ is the size of fragment at time $`t`$, $`\lambda `$ is a positive and constant value and $`\stackrel{~}{y}(t)`$ is a size for which the logarithmic term is zero at time $`t`$. In this model we assume that $`y_m(t)<\stackrel{~}{y}(t)`$. For a fragment of size $`y>\stackrel{~}{y}(t)`$ the logarithmic term in the rate of fragmentation is always positive, but a fragment of size $`y_m(t)<y<\stackrel{~}{y}(t)`$ gives a negative value for the logarithmic term. In order to have a positive rate of fragmentation for all particles larger than $`y_m(t)`$, we require $`\left[1+\lambda \mathrm{ln}\left({\displaystyle \frac{y(t)}{\stackrel{~}{y}(t)}}\right)\right]>0`$ which gives $`y(t)>\stackrel{~}{y}(t)e^{\mathrm{\hspace{0.17em}1}/\lambda }`$. If we choose $`y_m(t)`$ to be

$$y_m(t)=\stackrel{~}{y}(t)e^{\mathrm{\hspace{0.17em}1}/\lambda },$$
(35)
then the rate of fragmentation for all fragments of sizes larger than $`y_m(t)`$ is positive.

Consequently, the rate of fragmentation for this model is given by

$`\delta tR(y,z)`$ $`=`$ $`p_1{\displaystyle \frac{\delta L}{L(t)}}\delta (y+zL(t))+p_2{\displaystyle \frac{\delta t}{y+z}}\theta (y_m(t)(y+z))`$
$`+(1p_1p_2){\displaystyle \frac{\delta t}{y+z}}\left[1+\lambda \mathrm{ln}\left({\displaystyle \frac{y+z}{\stackrel{~}{y}(t)}}\right)\right]\theta ((y+z)y_m(t)),`$
where $`y_m(t)`$ is defined as (35). The first two terms on the right hand side of (3) are the same as those in (2) and have the same definitions as before. The last term represents the rate of fragmentation for the fragments of sizes larger than $`y_m(t)`$ in time $`(t,t+\delta t)`$. Using (1,3) for the both size regions give the kinetic equations

$`n_1(y,t+\delta t)=`$ (37)
$`n_1(y,t)p_2\delta tn_1(y,t)+\mathrm{\hspace{0.33em}2}p_1{\displaystyle \frac{\delta L}{L(t)}}n_2(L(t),t)`$
$`+\mathrm{\hspace{0.25em}2}p_2\delta t{\displaystyle _y^{y_m(t)}}n_1(z,t){\displaystyle \frac{dz}{z}}+2(1p_1p_2)\delta t{\displaystyle _{y_m(t)}^{L(t)}}n_2(z,t){\displaystyle \frac{dz}{z}}`$
$`+\mathrm{\hspace{0.25em}2}\lambda (1p_1p_2)\delta t{\displaystyle _{y_m(t)}^{L(t)}}n_2(z,t)\mathrm{ln}\left({\displaystyle \frac{z}{\stackrel{~}{y}(t)}}\right){\displaystyle \frac{dz}{z}}`$
for $`y<y_m(t)`$ and

$`n_2(y,t+\delta t)=`$ (38)
$`n_2(y,t)(1p_1p_2)\delta tn_2(y,t)\lambda (1p_1p_2)\delta t\mathrm{ln}\left({\displaystyle \frac{y}{\stackrel{~}{y}(t)}}\right)n_2(y,t)`$
$`+\mathrm{\hspace{0.25em}2}p_1{\displaystyle \frac{\delta L}{L(t)}}n_2(L(t),t)+\mathrm{\hspace{0.25em}2}(1p_1p_2)\delta t{\displaystyle _y^{L(t)}}n_2(z,t){\displaystyle \frac{dz}{z}}`$
$`+\mathrm{\hspace{0.25em}2}\lambda (1p_1p_2)\delta t{\displaystyle _y^{L(t)}}n_2(z,t)\mathrm{ln}\left({\displaystyle \frac{z}{\stackrel{~}{y}(t)}}\right){\displaystyle \frac{dz}{z}}`$
for $`y>y_m(t)`$. Following the same procedure as that used in the model A, we obtain

$`{\displaystyle \frac{F_1(x,t)}{t}}=`$
$`\nu x{\displaystyle \frac{F_1(x,t)}{x}}+2\nu F_2(1,t)+\mathrm{\hspace{0.17em}2}p_2{\displaystyle _x^{x_m}}F_1(z,t){\displaystyle \frac{dz}{z}}`$
$`+\mathrm{\hspace{0.33em}2}(1p_1p_2){\displaystyle _{x_m}^1}F_2(x,t){\displaystyle \frac{dx}{x}}+\mathrm{\hspace{0.17em}2}\lambda (1p_1p_2){\displaystyle _{x_m}^1}F_2(x,t)\mathrm{ln}\left({\displaystyle \frac{x}{\stackrel{~}{x}}}\right){\displaystyle \frac{dx}{x}}`$
$`F_1(x,t)[\nu +p_2+{\displaystyle \frac{\nu }{A}}F_2(1,t)+{\displaystyle \frac{p_2}{A}}{\displaystyle _0^{x_m}}F_1(x,t)dx`$
$`+{\displaystyle \frac{(1p_1p_2)}{A}}{\displaystyle _{x_m}^1}F_2(x,t)dx+{\displaystyle \frac{\lambda (1p_1p_2)}{A}}{\displaystyle _{x_m}^1}F_2(x,t)\mathrm{ln}\left({\displaystyle \frac{x}{\stackrel{~}{x}}}\right)dx]`$
for $`x<x_m`$ and

$`{\displaystyle \frac{F_2(x,t)}{t}}=`$
$`\nu x{\displaystyle \frac{F_2(x,t)}{x}}+2\nu F_2(1,t)+\mathrm{\hspace{0.17em}2}(1p_1p_2){\displaystyle _x^1}F_2(z,t){\displaystyle \frac{dz}{z}}`$
$`+\mathrm{\hspace{0.25em}2}\lambda (1p_1p_2){\displaystyle _x^1}F_2(z,t)\mathrm{ln}\left({\displaystyle \frac{z}{\stackrel{~}{x}}}\right){\displaystyle \frac{dz}{z}}F_2(x,t)[\nu +(1p_1p_2)`$
$`+\lambda (1p_1p_2)\mathrm{ln}\left({\displaystyle \frac{x}{\stackrel{~}{x}}}\right)+{\displaystyle \frac{\nu }{A}}F_2(1,t)+{\displaystyle \frac{p_2}{A}}{\displaystyle _0^{x_m}}F_1(x,t)𝑑x`$
$`+{\displaystyle \frac{(1p_1p_2)}{A}}{\displaystyle _{x_m}^1}F_2(x,t)dx+{\displaystyle \frac{\lambda (1p_1p_2)}{A}}{\displaystyle _{x_m}^1}F_2(x,t)\mathrm{ln}\left({\displaystyle \frac{x}{\stackrel{~}{x}}}\right)dx]`$
for $`x>x_m`$. Here $`x`$, $`x_m`$, $`\stackrel{~}{x}`$ and $`\nu `$ have the same definitions as before and from (35) we get

$$x_m=\stackrel{~}{x}e^{\mathrm{\hspace{0.17em}1}/\lambda }.$$
(41)
In the long time limit as $`t\mathrm{}`$, we assume that $`F_1(x,t)`$ and $`F_2(x,t)`$ evolve to the same time-independent variables as (15) and (16), respectively. In order to solve (3,3) in this limit we anticipate solutions as

$$F_1(x)=Gx^\gamma $$
(42)
for $`x<x_m`$ and

$$F_2(x)=H\frac{\mathrm{exp}\left\{\eta \left[\mathrm{ln}\left(\frac{x}{\stackrel{~}{x}}\right)\right]^2\right\}}{x}$$
(43)
for $`x>x_m`$. Substituting (42,43) in (3,3) in the long time limit, gives us

$$\gamma =\frac{p_2}{1p_1p_2},$$
(44)
$$\eta =\frac{\lambda }{2}$$
(45)
and

$$\nu =\mathrm{\hspace{0.33em}1}p_1p_2.$$
(46)
G and H in (42,43) are constant values and can be obtained using (3,3) in the long time limit and the continuity condition at $`x_m`$, (22). These give

$$H=\frac{A}{{\displaystyle \frac{1p_1p_2}{1p_12p_2}}\mathrm{exp}\left({\displaystyle \frac{1}{2\lambda }}\right)+\sqrt{{\displaystyle \frac{2\pi }{\lambda }}}\left[P\left(\sqrt{\lambda }\mathrm{ln}\left({\displaystyle \frac{1}{\stackrel{~}{x}}}\right)\right)P\left({\displaystyle \frac{1}{\sqrt{\lambda }}}\right)\right]}$$
(47)
and

$$G=H\frac{\mathrm{exp}\left({\displaystyle \frac{1}{2\lambda }}\right)}{\left(\stackrel{~}{x}e^{\mathrm{\hspace{0.17em}1}/\lambda }\right)^{\frac{1p_12p_2}{1p_1p_2}}},$$
(48)
where we have assumed that

$$0\gamma =\frac{p_2}{1p_1p_2}<1.$$
(49)
The second and third terms in the denominator of (47) are a proportion of the normal probability distribution function with the same definition as (25). Solutions (42,43) using (44-49) satisfy the normalisation relation (27). If $`A=1`$, the normalisation is equal to $`1`$. For a choice of $`A`$ as

$$A=\sqrt{\frac{\lambda }{2\pi }}\frac{1p_1p_2}{1p_12p_2}\mathrm{exp}\left(\frac{1}{2\lambda }\right)+P\left(\sqrt{\lambda }\mathrm{ln}\left(\frac{1}{\stackrel{~}{x}}\right)\right)P\left(\frac{1}{\sqrt{\lambda }}\right),$$
(50)
the solutions (42,43) using (44,45,47,48) yield

$$F_1(x)=\sqrt{\frac{\lambda }{2\pi }}\frac{\mathrm{exp}\left({\displaystyle \frac{1}{2\lambda }}\right)}{\left(\stackrel{~}{x}e^{\mathrm{\hspace{0.17em}1}/\lambda }\right)^{\frac{1p_12p_2}{1p_1p_2}}}x^{\frac{p_2}{1p_1p_2}}$$
(51)
for $`x<\stackrel{~}{x}e^{\mathrm{\hspace{0.17em}1}/\lambda }`$ and

$$F_2(x)=\sqrt{\frac{\lambda }{2\pi }}\frac{\mathrm{exp}\left\{{\displaystyle \frac{\lambda }{2}}\left[\mathrm{ln}\left({\displaystyle \frac{x}{\stackrel{~}{x}}}\right)\right]^2\right\}}{x}$$
(52)
for $`x>\stackrel{~}{x}e^{\mathrm{\hspace{0.17em}1}/\lambda }`$. We see that the stable distribution has a power law behaviour in the smaller region. Comparing (52) with (31) exhibits that the stable distribution $`F_2(x)`$ in the larger region has a log-normal form with a mean of $`\stackrel{~}{x}`$ and a dispersion of the distribution equal to

$$\sigma =\frac{1}{\sqrt{\lambda }}.$$
(53)
It is shown that the log-normal function (52) has a maximum at $`x=\stackrel{~}{x}e^{\mathrm{\hspace{0.17em}1}/\lambda }`$ which we have assumed to be equal to $`x_m`$. This shows that in this model the dimensionless transition size $`x_m`$ is chosen to be equal to the value for which the log-normal distribution in the larger region is maximum. To have a positive value for $`A`$ in (50) and therefore physically meaningful distributions, it is required that $`(1p_1p_2)>p_2`$. This inequality is satisfied from (49).

Now we consider some special cases of the model B.

1. When there is no fragmentation for the fragments smaller than $`y_m(t)`$, i.e. $`p_2=0`$, then the stable distribution from (51) yields

$$F_1(x)=\sqrt{\frac{\lambda }{2\pi }}\frac{1}{\stackrel{~}{x}}e^{1/2\lambda }$$
(54)
for $`x<x_m`$ and is the same as (52) for $`x>x_m`$. The distribution is a constant value in the smaller region and has a log-normal behaviour in the larger region.
2. As $`\lambda 0`$, and therefore from (53) $`\sigma ^2\mathrm{}`$, the situation is then similar to the case $`3`$ of model A and the log-normal distribution (52) in the larger region can be approximated by a power law form with an exponent equal to $`1`$, $`F_2(x)x^{\mathrm{\hspace{0.17em}1}}`$. On the other hand at this limit the rate of fragmentation (3) gives two different probabilities in the two different size regions with a rate proportional to the inverse of the fragment size in each region. The stable distribution then reveal a composite power law with power laws of exponents of $`{\displaystyle \frac{p_2}{\nu }}`$ and $`{\displaystyle \frac{1p_1p_2}{\nu }}`$ in the smaller and larger size regions, respectively. Since in our case from (46) we have $`\nu =(1p_1p_2)`$, the exponents yield $`{\displaystyle \frac{p_2}{1p_1p_2}}`$ and $`1`$, respectively. The result in the smaller size region is equivalent to that from (42,44). In the larger size region, as expected, we get a power law with an exponent of $`1`$ which is consistent with the property of a log-normal distribution with a very large variance.

We see that all the results for the model B can be recovered from the model A by substituting $`p_3`$ by $`(1p_1p_2)`$, $`p_4`$ by $`\lambda (1p_1p_2)`$ and using equation (35).

## 4 Conclusions

In this paper we have studied two models of binary fragmentation. In the models we introduced a time dependent transition size, $`y_m(t)`$, which produced two regions of fragment sizes with different rates of fragmentation above and below the transition size. The rate of fragmentation at each time step is fixed for the largest fragments and is proportional to the inverse of fragment size for the fragments smaller than the transition size. We considered a rate of fragmentation combined of two terms, one proportional to the inverse of the fragment size and the other proportional to a logarithmic function of the fragment size, for the fragments larger than the transition size at each time step. For a size equal to $`\stackrel{~}{y}(t)`$ at time $`t`$, the logarithmic function is zero. In model A we assumed that the transition size, $`y_m(t)`$ is any arbitrary value for which $`y_m(t)\stackrel{~}{y}(t)`$, whereas in the model B we chose a value as (35) for the transition size. The models were then solved exactly in the long time limit to reveal stable time-invariant solutions for the distributions. The results of these distributions for both models A and B exhibited a power law form with an exponent between $`1`$ and zero, in the smaller region. For a special choice of normalisation, the distributions exhibited log-normal behaviours in the larger regions. In this region they revealed one of interesting properties of a log-normal distribution for which as the variance of the log-normal is pretty large, then the distribution is approximated by a power law form with an exponent equal to $`1`$. Special cases of the models with no fragmentation for the smaller fragments were also examined. The stable distributions then exhibited a constant value and a log-normal form in the smaller and larger regions, respectively.

Now we investigate some experimental work on the fragmentation involving log-normal distributions and the transition of the distributions from a log-normal form to a power law one and vice versa. In  the experimental results of fragment size and mass distributions for the long thin glass rods with fixed lengths of $`1500mm`$ and diameters of $`2mm`$ were reported. The rods were dropped horizontally onto a flat hard floor from different heights. At lower falling heights (about $`1m`$ drop), the distributions exhibited a log-normal form for larger fragments and a power law form for smaller ones. The crossover was seen to be at length scales around the rod diameter. This is due to the fact that fragments smaller than the rod diameter, undergo a three-dimensional fracture, whereas the larger fragments are produced by one-dimensional breaking. For a falling height of $`1.20m`$ the cumulative number of fragments on their size were plotted. The data was fitted to a curve calculated by assuming that the size distribution could be described by a log-normal distribution. The fitting gave values of $`\stackrel{~}{l}=30mm`$ and $`\sigma =0.55`$ for the mean and the dispersion of the log-normal distribution. The fitting is excellent for the fragments larger than $`l_c=7mm`$ which is a length scale around the rod diameter. The size distribution of fragments smaller than $`l_c`$ seemed to have a power law form. As the falling height was increased, the cumulative number variation for larger fragments started to show a power law dependence on their size and mass. In  experiments in rupture of mercury droplets were performed. In the experiment at a height of $`h`$, the mercury droplets of about $`2.0mm`$ radius were fallen directly on a glass Petri dish. The results of the cumulative number of the drops versus the drop diameter were plotted. For small falling heights a log-normal behaviour was exhibited. The distribution showed a clear transition from a log-normal form to a scaled one as the falling height $`h`$ was increased.

Our models in the long time limit can be applied to the experimental systems when the fragmentation process is over and the system reaches a stable and time-independent phase. In the models we assumed that the fragments have only one dimension, whereas in real physical systems the fragments have three dimensions. We can therefore only apply our models to the systems for which the fractures occur in one dimension. This assumption seems physically reasonable for a long thin glass rod  for which most of the fractures occur in the length of the rod and the cross-section remains almost fixed. We can also consider this assumption for the mercury droplets  for which at each time step only the radius of the droplet changes and the shape of the fragments remains always spherical.

Now we discuss the similarities between the results from our models and those from the experiments. The experimental results for the size distribution of long thin glass rods at low falling heights exhibited a power law and log-normal behaviour in the smaller and larger size regions, respectively . We got the same behaviour for the distributions in both our models. The data from the experiments for the falling height of $`1.20m`$ revealed that $`l_c<\stackrel{~}{l}`$, where $`l_c`$ is the transition size and $`\stackrel{~}{l}`$ is the mean of the log-normal distribution. In the model B of our work also we have $`x_m<\stackrel{~}{x}`$. We therefore expect model B of our work to describe the experimental results. In model B, the dispersion of the log-normal distribution can be obtained from (53) and the mean and the transition values, $`\stackrel{~}{x}`$ and $`x_m`$, respectively, are any arbitrary values which satisfy the relation (41). Then, by appropriate choices of $`\stackrel{~}{x}`$ and $`\lambda `$ in the model B of our work we can match the mean and the dispersion of the log-normal distribution in the experimental results. Because of (41,53) we should expect a relation as $`l_c=\stackrel{~}{l}e^{\sigma ^2}`$. The values of $`l_c`$ and $`\stackrel{~}{l}`$ from experimental result do not satisfy this relation. As a result, although model B of our work does not match the exact relation between $`l_c`$ and $`\stackrel{~}{l}`$ from experimental result, but it can give an explanation for the transition of distribution from a power law to a log-normal form on the shock fragmentation of long thin glass rods . Also in the smaller region of the model B, the appropriate choices of $`p_1`$ and $`p_2`$ can match the exponent of power law equal to $`{\displaystyle \frac{p_2}{(1p_1p_2)}}`$ from (42,44), with the corresponding exponent of power law distribution from the experimental results. This exponent needs to satisfy the condition (49). For high falling heights a transition of the distribution from log-normal to scaling law was observed in the experiments for the larger rod fragments  and for the mercury droplets . At high heights therefore the distributions from the experiments exhibit a composite and a single power law form for the rod fragments and mercury droplets, respectively. The models in this paper suggest a power law and a log-normal behaviour for the distribution of the smaller and larger fragments, respectively. To explain composite and single power law distributions, we can however use the models in .

To conclude we introduced statistical models with a fixed probability for the breaking the largest particle in the system. All other particles of sizes smaller than the largest size, were divided into two different regions with different rates of fragmentation in each region. In the models therefore the largest fragments are broken with very large probability. This fact together with the form of the fragmentation rate as a function of the fragment size in each region, give a distribution of power law or log-normal form in that region. Our models therefore give a theoretical explanation for the transition of distribution from the scaling to log-normal behaviour, which had been seen before in the experiments.

Acknowledgement

ZT is grateful to the Soudavar Scholarship Fund for the award of a scholarship.


format: text

text:
# Measuring the aspect ratio renormalization of anisotropic-lattice gluons

## 1 Introduction

In principle, an improved action makes it possible to achieve lattice volumes large enough to overcome finite size effects at a computational effort low enough to obtain measurements with good statistical errors. However typical masses for heavy states may be similar to or larger than the inverse lattice spacing on the associated coarse lattice. Propagators for such states then decay too fast in lattice units for accurate measurement.

The problem can be overcome by tuning couplings associated with the “time” direction so that the temporal lattice spacing $`a_t`$ is much smaller than the spatial lattice spacing $`a_s`$ . Refined measurements are then possible while retaining the computational advantages of the improved action and coarse spatial lattice. Such anisotropic actions have already been successfully applied to the glueball spectrum , the spectrum of excitations of the inter-quark potential , heavy hybrids , and the fine structure of the quarkonium spectrum . They are also expected to be a powerful tool in extracting excited-state signals, obtaining high-momentum form factors, pushing thermodynamic calculations to higher temperatures, and in the calculation of transport coefficients.

In this paper, we consider one formulation of the anisotropic action for the pure Yang-Mills sector :

$`S_n=`$
$`\beta {\displaystyle \underset{x,s>s^{}}{}}\chi _0^1\left\{{\displaystyle \frac{5}{3}}{\displaystyle \frac{P_{s,s^{}}}{u_s^4}}{\displaystyle \frac{1}{12}}{\displaystyle \frac{R_{ss,s^{}}}{u_s^6}}{\displaystyle \frac{1}{12}}{\displaystyle \frac{R_{s^{}s^{},s}}{u_s^6}}\right\}\beta {\displaystyle \underset{x,s}{}}\chi _0\left\{{\displaystyle \frac{4}{3}}{\displaystyle \frac{P_{s,t}}{u_s^2u_t^2}}{\displaystyle \frac{1}{12}}{\displaystyle \frac{R_{ss,t}}{u_s^4u_t^2}}\right\}.`$
Here $`s,s^{}`$ run over spatial directions, $`P_{s,s^{}}`$ is a $`1\times 1`$ plaquette and $`R_{s^{}s^{},s}`$ is a $`2\times 1`$ rectangle. The coefficients of these terms are chosen so the action has no $`O(a_s^2)`$ discretization errors in tree-level pertubation theory. $`\chi _0`$ is the anisotropy parameter and is equal to the aspect ratio of the spatial and temporal lattice spacings, $`a_s`$ and $`a_t`$ at tree-level. At higher orders in the perturbative expansion, this aspect ratio recieves quantum corrections, so a renormalized anisotropy determined from a physical process, $`\chi _R`$, differs from $`\chi _0`$ at $`O(\alpha _s)`$. Tadpole improvement (TI) of the perturbative expansion is achieved by tuning the input spatial and temporal link parameters $`u_s`$ and $`u_t`$ for self-consistency at each choice of $`\beta `$ and $`\chi _0`$.

Here we report on the accurate determination of the tadpole parameters in the plaquette and Landau mean-link formulations. (Non-tadpole-improved actions have been studied in ref. .) We compare methods for measuring the renormalized anisotropy and discuss the importance of lattice artifacts. In order that the anisotropic formulation can be used with confidence, physically distinct methods for determining the renormalized anisotropy, $`\chi _R=a_s/a_t`$, should agree. Discrepancies in the results for $`\chi _R`$ give a measure of lattice artifacts and an indication of the importance of further improvements. Our results have been used in applications of anisotropic lattices to heavy quark hybrid states .

In section 2 we report on the determination of the tadpole parameters; in section 3 the renormalized anisotropy measurements are presented, in section 5 we discuss the results and present conclusions.

## 2 Tadpole Parameters

In both the plaquette and Landau mean-link formulations, the tadpole parameters are determined self-consistently and are defined by

$$\begin{array}{ccccccccc}u_s\hfill & =& \left(P_{s,s^{}}\right)^{1/4}\hfill & ,\hfill & u_t\hfill & =& 1\hfill & & \text{plaquette},\hfill \\ & & & & & & & & \\ u_s\hfill & =& U_s_{Landau}\hfill & ,\hfill & u_t\hfill & =& U_t_{Landau}\hfill & & \text{Landau},\hfill \end{array}$$
(2.1)
where the Landau gauge is defined by the field configuration which maximizes

$$F(\{U\})=\underset{x\mu }{}\frac{1}{u_\mu a_\mu ^2}\text{ReTr}\left\{U_\mu (x)\frac{1}{16u_\mu }U_\mu (x)U_\mu (x+\widehat{\mu })\right\}.$$
(2.2)
with respect to gauge transformations. We denote the gauge coupling in the two schemes as $`\beta _L`$ and $`\beta _P`$ respectively.

### 2.1 Plaquette Tadpoles

In the plaquette scheme of Eqn. (2.1), $`u_t=1`$ and $`u_s`$ is determined self-consistently. The expectation value of the spatial plaquette is computed for a range of input parameters $`u_s`$ close to and spanning the self-consistent value, $`u_s^{}`$. A linear interpolation is sufficient to give an accurate value of $`u_s^{}`$. This value is then checked in a further Monte-Carlo simulation. The values of $`u_s^{}`$ for a set of couplings $`\beta _P`$ and tree-level anisotropies $`\chi _0`$ are given in Table 1. For the lattice sizes used in simulation, the plaquette expectation value is found to be independent of the volume at the four-significant-figure level. An extrapolation to infinite volume is unnecessary.

### 2.2  Landau Gauge Fixing

The maximization of $`F(U_\mu ^g)`$ in Eqn. (2.2), where $`U_\mu ^g=g(x)U_\mu (x)g^{}(x+\mu )`$ with respect to a gauge transformation $`\{g(x)\}`$, was carried out using the conjugate-gradient method modified to deal with the group structure of the link elements. At each stage of the iteration the the appropriate conjugate-gradient vector $`\{𝒗(x)\}`$, is computed as a covariant derivative .

$$𝒗(x)=\frac{}{\eta (x)}F(U_\mu ^g)|_{\eta =0},$$
(2.2.1)
where $`g(x)=e^{\eta (x).T}`$ and $`T`$ are the generators of $`SU(3)`$ . The vector $`\{𝒗(x)\}`$ has elements lying in the Lie algebra of $`SU(3)`$ . For each vector a series of group elements is constructed in the associated one-parameter subgroup, each element being a given power of the preceding one:

$$g_1=\mathrm{exp}(ϵ𝒗),g_n=(g_{n1})^p,1<nN.$$
(2.2.2)
The local maximum in the direction $`𝒗`$ is found by evaluating $`F(U_\mu ^{g_p})`$ in descending order from $`p=N`$ to $`p=1`$ . The value of $`ϵ`$, at any stage, is reduced until either the sequence exhibits a maximum or remains constant within a preset tolerance.

For large $`\chi _0`$ the terms including $`U_t`$ dominate the expression for $`F`$ and near to the maximum relatively small changes in $`F`$ correspond to quite large changes in $`u_s`$. Consequently, to be sure that $`u_s`$ as well as $`u_t`$ is accurately calculated, the maximum must be found to a sufficiently high accuracy. The criterion chosen was that $`\delta F/F2\times 10^6`$ where $`\delta F`$ was the accumulated change in $`F`$ for three consecutive iterations of the conjugate-gradient algorithm for which the individual changes in $`F`$ were non-zero. This criterion accounts for the observation that from time to time close to the maximum the change in $`F`$ was zero. We chose $`ϵ=2\times 10^5`$, $`p=4`$, and $`N=6`$.

### 2.3  Landau Gauge Fixed Tadpoles

The self-consistent value of $`𝒖=(u_s,u_t)`$ can be determined by a generalized Newton-Rapheson method or by linear interpolation. The latter method was used for the mean-link tadpoles and was implemented by choosing four input tadpoles $`𝒖_i,i=1,\mathrm{},4`$ and measuring the four corresponding output tadpoles $`𝒖_i^{},i=1,4`$. A linear map is assumed for the incremental vectors $`𝐫_i=(𝒖_i𝒖_1),𝐫_i^{}=(𝒖_i^{}𝒖_1^{}),i=2,3,4`$:

$$𝐫_i=𝑴𝐫_i^{},i=2,3,4,$$
(2.3.1)
where $`𝑴`$ is a $`2\times 2`$ matrix that is determined from the the measured images of $`𝐫_2,𝐫_3`$ and checked for consistency against the measured image of $`𝐫_3`$ . The self-consistent tadpole $`𝒖^{}`$ is then predicted to be

$$𝒖^{}=𝒖_1+(1𝑴)^1(𝐫_1^{}𝐫_1).$$
(2.3.2)
The self-consistency of $`𝒖^{}`$ is then checked computationally. This method was found to be very reliable and only subject to minor adjustments to account for statistical errors.

Moreover, a meaningful statistical error can be assigned to $`𝒖^{}`$ from the image under $`M^1`$ of the statistical error box deduced for $`𝒖^{}`$ from the measurements values of the tadpole parameters. Let the image of $`𝒖^{}`$ under the map $`M`$ be $`𝒖^{}`$ with statistical error $`\delta 𝒖^{}`$. The value $`𝒖^{}`$ is acceptable if in terms of components $`𝒖^{}\delta 𝒖^{}<𝒖^{}<𝒖^{}+\delta 𝒖^{}`$, and the error quoted on $`𝒖^{}`$ is $`\delta 𝒖^{}=M^1\delta 𝒖^{}`$. A typical example is for $`\beta _L=1.8,\chi _0=4`$ on a $`6^3\times 24`$ lattice where

$$𝑴^1=\left(\begin{array}{cc}\hfill 0.341& \hfill 0.227\\ \hfill 0.033& \hfill 0.957\end{array}\right),\delta 𝒖^{}=\left(\begin{array}{c}\hfill 3\times 10^4\\ \hfill 3\times 10^5\end{array}\right)\delta 𝒖^{}=\left(\begin{array}{c}\hfill 1\times 10^4\\ \hfill 2\times 10^5\end{array}\right).$$
(2.3.3)
In practice it is found that the error on $`u_s^{}`$ is typically reduced by a factor of 3 compared with the statistical error from the simulation verifying the self-consistency.

The simulation for the mean-link tadpoles was done on the Hitachi SR2201 computers at the Cambridge High Performance Computing Facility and the Tokyo University Computer Centre. The lattices used were $`L^3\times T`$ with $`T=\chi _0L`$ and $`L=6,8`$ in all cases except one. The high accuracy demanded by the maximization process was very time consuming. Also, many configurations were required to achieve the desired statistical accuracy. Consequently, only one example with $`L=10`$ was done as a check on the finite-size scaling ansatz, $`𝒖^{}(L)=𝒖^{}(\mathrm{})+A/L^2`$ . Typically, for $`10^3\times 40`$ about 400 conjugate-gradient iterations were needed taking about 60 seconds per iteration. In the cases $`L=6,8,10`$ for $`\beta _L=1.8\chi _0=4`$, the fit to the finite-size scaling ansatz was very good as can be seen in Figs. 1 and 2. For other cases the linear extrapolation was assumed to hold. The values measured for $`𝒖^{}(L)`$ and the extrapolation to $`L=\mathrm{}`$ are given in table 2.

## 3 Anisotropy

In this section we report on two methods for determining $`\chi _R`$. The first uses the dispersion relation for the torelon  and the second uses the comparison of the potential measured in the fine and coarse direction using Wilson loops.

### 3.1 The Torelon

The lattice considered is $`S^2\times L\times T`$ with typical values $`S=8`$ in the $`x`$ and $`y`$ directions, $`L=3,4,5`$ in the $`z`$ direction, and $`T=50`$ in the fine or $`t`$ direction. The lattice has periodic boundary conditions and the torelon is created by a Polyakov line that loops around the lattice in the $`z`$ direction. The Polyakov line is associated with a particular point in the $`(x,y)`$-plane and is constructed from links which are covariantly smeared in a manner similar to APE smearing:

$`T(x,y,t)`$ $`=`$ $`Tr{\displaystyle \underset{z=0}{\overset{z=L1}{}}}W_z(x,y,z,t),`$
$`W_z(x,y,z,t)`$ $`=`$ $`(1+{\displaystyle \frac{l^2D^2}{4m}})^mU_z(x,y,z,t),`$ (3.1.1)
where $`D`$ is the appropriate covariant derivative. Typically, $`l=1.5`$, $`m=10`$. The state with momentum $`𝒑=(p_x,p_y)=(n_x,n_y)(2\pi /Sa_s)`$ is

$$T(𝒑,t)=\underset{x,y}{}T(x,y,t)e^{i(p_xx+p_yy)}.$$
(3.1.2)
The torelon propagator

$$G_T(𝒑,t)=\frac{1}{T}\underset{t^{}=0}{\overset{t^{}=T1}{}}T(𝒑,t^{})T^{}(𝒑,t+t^{}),$$
(3.1.3)
is measured for various choices of momentum $`𝒑`$ and a correlated simultaneous fit using SVD decomposition is made to the relativistic dispersion formula:

$`G_T(𝒑,t)`$ $`=`$ $`c(𝒑)e^{E(𝐩)a_t\overline{t}},`$
$`E(𝒑)a_t`$ $`=`$ $`{\displaystyle \frac{a_s\sqrt{𝒑^2+M_T^2}}{\chi _R}},`$ (3.1.4)
where $`M_T`$ is the torelon mass and $`\overline{t}=t/a_t`$. The momenta used were $`𝒏^2=n_x^2+n_y^2=0,1,2,4,5`$ with $`S=8`$. An example of the fit obtained is shown in Fig. 3 for $`\beta _L=1.8,\chi _0=6`$. The fit is over the range $`2t12`$ and has $`\chi ^2/N_{\mathrm{df}}=0.96`$ for 43 d.o.f. The fit is good and gives $`(M_Ta_s)^2=1.85(4),\chi _R=3.61(2)`$. In Fig. 4 the dispersion curves $`E(𝒑^2)`$ versus $`𝒏^2`$ are plotted for $`\beta _L=1.8,\chi _0=4,L=3,4,5`$. From both Figs. 3 and 4 it is clear that rotational invariance is established in the coarse $`xy`$ directions. The good simultaneous fit to the torelon propagators in all cases shows that the continuum dispersion is well satisfied for momenta at any angle to the coordinate axes.

In table 3 we give the full set of results for the torelon anisotropies we have measured. There is no noticeable dependence within errors of $`\chi _R`$ on $`L`$ and the dependence on $`\beta `$ is small for the two values of $`\beta `$ used.

From the values $`M_T(La_s)`$ we can determine the string tension, $`\sigma `$ in units of $`a_s^2`$. We assume that $`M_T(La_s)=\sigma (La_s)La_s`$ where $`\sigma (La_s)`$ is the string tension modified by finite-size corrections : $`\sigma (La_s)=\sigma +D/(La_s)^2`$. A plot of $`\sigma (La_s)a_s^2`$ versus $`1/L^2`$ is shown in Fig. 5 with a linear fit which has $`\chi ^2/N_{\mathrm{df}}=0.05`$ with $`\sigma a_s^2=0.394(3)`$. In another calculation we have determined the absolute value of $`a_t`$ in $`MeV`$ by computing the splitting $`\mathrm{\Delta }M_{PS}=M_\mathrm{{\rm Y}}(1P)M_\mathrm{{\rm Y}}(1S)`$ for the $`\mathrm{{\rm Y}}`$ meson system using $`O(mv^6)`$ NRQCD . In table 4 we give $`a_t^1,a_s^1,\sigma a_s^2`$, the coefficient $`D`$, and $`\sqrt{\sigma }/\mathrm{\Delta }M_{PS}`$. This latter ratio is experimentally close to unity but in ref.  using isotropic $`\beta =6.2`$ UKQCD configurations this ratio was found to be about 1.25 which, from table 4, agrees with our finding except possibly for the $`(\beta =1.7,\chi _0)`$ lattice which has the most coarse lattice spacings. It is generally accepted that the discrepancy is due to quenching and the significant outcome is that our results from the anisotropic lattice agree well with those of a spatially finer isotropic lattice of $`a=0.06`$fm, indicating that we have correctly reproduced the physics expected for this comparison on lattices with $`a_s=0.250.30`$fm. A naive estimate for the coefficient $`D`$ is obtained by calculating the contribution to $`M(L)`$ from the zero-point fluctuations of the torelon regarded as a periodic flux tube. The outcome is $`D=\pi /31.05`$ which is compatible with our fit of $`D1.35(5)`$. The same calculation predicts that there is no $`L`$-independent constant contribution to $`M(L)`$ and this is verified by our fits.

The torelon dispersion relation was also studied on the plaquette-tuned parameters. For these simulations, lattices of extent $`(8\times 10)\times 4\times N_t`$ were used, to enable us to investigate a wider range of momenta combinations. At momenta close to the cut-off, we expect the discretization errors to be $`O(a_s^4p^4,\alpha _sa_s^2p^2)`$, thus in order not to contaminate our determination of $`\chi _R`$ with large discretization errors, we first determined the range of momenta over which a good correlated fit to the continuum dispersion relation could be made. The data for different momentum ranges from simulations at $`\beta _P=2.1,\chi _0=6`$ were tested and these results are shown in Table 5. Both the measured anisotropy renormalizations and the torelon rest energies determined are all consistent within statistical precision up to the highest momentum measured. For the subsequent computations of the anisotropy, presented in Table 6, the largest momentum used was the (1,1) data (shown in bold in Table 5).

### 3.2 The Sideways Potential

In the sideways potential method , a coarse direction, $`z`$, on the anisotropic lattice is chosen to be the time direction. There are then two types of spacelike direction, coarse and fine. The measurement of the potential between static quarks with a separation lying in the plane of coarse links is compared with the measurement of the potential when the separation lies along the line of fine links. The demand that the two measurements yield the same function of physical distance determines the renormalized anisotropy. Points in the coarse-coarse plane are denoted by $`\stackrel{}{x}=(x,y)`$ and points in the fine direction by $`t`$ where $`x,y,z,t`$ are integers.

We measure appropriate spatial Wilson loops $`W_{ss}(\stackrel{}{x},z)`$ and also loops using the fine direction $`W_{ts}(t,z)`$. We define

$$V_s(\stackrel{}{x},z)=\mathrm{log}\left(\frac{W_{ss}(\stackrel{}{x},z)}{W_{ss}(\stackrel{}{x},z+1)}\right),V_t(t,z)=\mathrm{log}\left(\frac{W_{ts}(t,z)}{W_{ts}(t,z+1)}\right).$$
(3.2.1)
As $`z\mathrm{}`$, $`V_s(\stackrel{}{x},z)V_s(|\stackrel{}{x}|)`$, and $`V_t(t,z)V_t(t)`$ where $`V_s(|\stackrel{}{x}|)`$ and $`V_t(t)`$ are the two versions of the interquark potential. For a physical distance $`r`$ we have $`|\stackrel{}{x}|a_s=ta_t=r`$ . We therefore estimate the renormalized anisotropy $`\chi _R`$ by tuning it so that $`V_s(|\stackrel{}{x}|)=V_t(t/\chi _R)`$, where the right side is evaluated by means of linear interpolation between the values measured at integral $`t`$. It is implicit in the method that there is effective rotational invariance in the $`\stackrel{}{x}`$-plane. We find that if we exclude potentials at the smallest distances, $`|\stackrel{}{x}|=1,\sqrt{2}`$, then the values of $`\chi _R`$ associated with different directions in the $`\stackrel{}{x}`$-plane generally agree within errors. The agreement is particularly good if the links are smeared in an appropriate manner. The results for the renormalized anisotropy, $`\chi _R`$, for both mean-link and plaquette schemes are shown in table 7.

An alternative approach to making the comparison is to fit the measured potentials with the forms

$`a_sV_s(\stackrel{}{x})=a_sV_0+\sigma a_s^2x+{\displaystyle \frac{e}{x}}`$
$`a_sV_t(t)=a_sV_0+\sigma a_sa_tt+{\displaystyle \frac{a_se}{a_tt}}.`$ (3.2.2)
The renormalized anisotropy, $`\chi _R`$, is then determined from the ratio of the coefficients of the linear terms in the two cases. $`\chi _R`$ can in principle be determined from the ratio of the coefficients of the coulombic terms however such an estimate depends on short distance effects and is inherently more sensitive to discretization errors.

This approach was tested on an $`8^3\times 48`$ lattice at $`\beta _P=2.3,\chi _0=6`$. Prior to measurement, the lattice is blocked to an $`8^4`$ volume by thinning the time-slices of spatial links and replacing the temporal links with the product of six fine links connecting the appropriate time-slices. The $`t`$ and $`z`$ axes are then interchanged, since as before we want to use a coarse direction for Euclidean decay. Then a standard APE-smearing algorithm is applied to links on the new time slices. These degrees of freedom are then used to measure the potential in the $`(x,y)`$ and $`z`$ axes in the standard fashion. Fig. 6 shows the results of this simulation. The best fits of this data to Eqn. (3.2.2) for a range of static-source separations were computed and the string tension results are shown in Table 8. In Fit C, the coulomb term was poorly resolved and was thus fixed to zero. The use of these different ranges allowed us to investigate the systematic uncertainty and discretization errors in determining $`\chi _R`$ from this method. The three results are consistent within $`2\%`$. The anisotropy measurement from the full range of separations (Fit A, shown in bold in Table 8) is in good agreement with the torelon dispersion result of Table 6.

Both the above approaches yielded results reasonably consistent with each other and with the toleron results at the 3% level. Discrepancies can easily be explained by the presence of discretization errors. Because the flux tube generating $`V_s(|\stackrel{}{x}|)`$ has a coarse-fine cross-section and that generating $`V_t(t)`$ a coarse-coarse cross-section we anticipate, as in the case of the toleron (eq.(4.5), that some discretization error in $`\chi _R`$ that is $`O((\alpha _sa_s/R)^2)`$ will remain. Since $`a_s0.3\mathrm{fm}`$, $`R0.7`$ to $`1\mathrm{fm}`$, and $`\alpha _s0.3`$, we expect discretization errors to be around 5%.

Because we are concerned to make the long distance behavior consistent in both the fine and coarse directions it is advantageous to use Wilson loops of the largest possible spatial extent. However, in practice, the statistical errors on large Wilson loops grow exponentially with separation. In order to achieve acceptable errors we were restricted to using loops of size 2 or 3 in coarse lattice units. This is to be compared with the torelon, where the practical size is 3 or 4 lattice units.

## 4  Discretization Errors

One way of viewing discretization errors is to regard them as arising from the absence of the correction term in the lattice action that yields continuum results. The corresponding term in the action density can be presumed to be a local (redundant) operator. The locality of the operator implies that while there may be a correction to the mass per unit length of the torelon, there is no correction that directly alters the asymptotic proportionality of the torelon mass and its length. The appropriate length scale against which to measure lattice artifacts therefore is the radius $`R`$ of the toleron cross-section or equivalently, $`\mathrm{\Lambda }_{QCD}^1`$ or $`\sqrt{\sigma }`$ where $`\sigma `$ is the string tension. We expect therefore that the $`O(a_s^2)`$ errors in $`\chi _R`$ are proportional to $`(a_s/R)^2`$ .

A persuasive plausibility argument for this behavior is as follows. The torelon state of momentum $`𝒑`$ and mass $`M_T`$ is an eigenstate of the transfer matrix $`T`$ with an eigenfunctional labelled by $`𝒑`$ and $`M_T`$ and eigenvalue $`\mathrm{exp}(E(𝒑^2)a_t/\chi _R)`$. Suppose an operator is added to the action which corrects for the $`O(\alpha _sa_s^2)`$ errors. This operator will be local on the scale of the lattice and the effect on $`E(𝒑^2)`$ can be estimated using first-order perturbation theory:

$$\delta \frac{a_sE(𝒑^2)}{\chi _R}=\alpha _sf(\overline{𝒑}^2,\overline{M}_T),$$
(4.1)
where $`\overline{𝒑}^2=𝒑^2a_s^2`$, $`\overline{M}_T=M_Ta_s`$, and $`f(\overline{𝒑}^2,\overline{M}_T)`$ is the dimensionless matrix element of the added operator which is proportional to $`a_s^2`$ by construction. The question is what scale balances the dimension of $`a_s^2`$? We find

$`{\displaystyle \frac{\delta \chi _R}{\chi _R^2}}`$ $`=`$ $`\alpha _s\varphi ^{}(0,\overline{M}_T),`$ (4.2)
$`{\displaystyle \frac{\delta \overline{M}_T^2}{2\chi _R}}`$ $`=`$ $`\alpha _s(\varphi (0,\overline{M}_T)\overline{M}_T^2\varphi ^{}(0,\overline{M}_T)),`$ (4.3)
where $`\varphi (\overline{𝒑}^2,\overline{M})=a_sE(𝒑^2)f(\overline{𝒑}^2,\overline{M}_T)`$ and $`\varphi ^{}=\varphi /(\overline{𝒑}^2)`$. It must be that the change in the action by a local operator corresponds to a change in the string tension, $`\sigma `$. Hence, since $`M_T=\sigma L`$, it follows that $`\delta \overline{M}_T^2\overline{M}_T^2`$ as this is the only parameter depending on $`L`$, and from eqn. (4.3) this implies that

$$\varphi (0,\overline{M}_T)\overline{M}_T^2,\varphi ^{}(0,\overline{M}_T)\overline{M}_T\text{-independent constant}.$$
(4.4)
It is conceivable that an accidental cancellation between the two terms in eqn. (4.3) would allow different behavior to be inferred for $`\varphi `$ and $`\varphi ^{}`$ but these results must hold true for all possible local perturbations, not just the particular one that eliminates $`O(\alpha _sa_s^2)`$ errors. We consider this kind of cancellation to be unlikely. Substituting this behavior for $`\varphi ^{}`$ into Eqn. (4.2) we find that

$$\frac{\delta \chi _R}{\chi _R^2}C\frac{a_s^2}{R^2},$$
(4.5)
where $`R`$ is a typical length scale associated with the torelon which cannot be $`M_T^1`$, for example, the flux tube radius. Alternatively, $`R^1`$ can be taken to be $`\sqrt{\sigma }`$, $`\mathrm{\Lambda }_{QCD}`$.

## 5  Conclusions

In this paper we have investigated various methods for measuring the renormalized aspect ratio for pure QCD on an anisotropic lattice. In the main we used a tadpole-improved action with the Landau gauge mean-field definition for the tadpole parameters, but we also included results for the action with tadpoles defined by the mean plaquette. The object of the investigation was to assess the consistency of different methods in order to judge the effectiveness of the improvement scheme. In principle, measurements of the anisotropy from different physical probes should agree close to the Euclidean-symmetric continuum limit. The bare anisotropy $`\chi _0`$ is renormalized by the effect of operators which are irrelevant in the neighborhood of the fixed point controlling the continuum limit. From a renormalization group (RG) view point the location of this fixed point is ambiguous up to redefinitions of the RG transformation used to locate it. This ambiguity is due to redundant operators  which have no effect on physical observables of the continuum theory. Consequently, continuum actions with different $`\chi _R`$ can differ only by redundant operators since they must correspond to the same continuum physics. Thus, while $`\chi _R`$ can be changed by tuning $`\chi _0`$, an action with no lattice artifacts must give rotationally invariant physical results and, consequently, any physical method for measuring $`\chi _R`$, such as the ratio of two physical observables, must give the same answer. In as much as this is not the case the differences will give an estimate for the effect of lattice artifacts and the necessity for improvement.

The spatial lattice spacing $`a_s`$ was measured in a separate NRQCD simulation by fitting the $`1P1S`$ mass-splitting $`\mathrm{\Delta }M_{PS}`$ for bottomonium. This estimate can be compared with the string tension, deduced from a $`D/L^2`$ extrapolation to $`\mathrm{}`$ of the torelon mass. The coefficient $`D`$ was found to be $`1.38`$ which is in tolerable agreement with $`\pi /3`$ predicted from analysis of flux-tube fluctuations . The ratio $`\sqrt{\sigma }/\mathrm{\Delta }M_{PS}`$ is $`1.2`$ in agreement with earlier NRQCD analyses . The departure from the experimentally observed value of $`1`$ is attributed to quenching. The values for both the mean-link and plaquette tadpoles were found self-consistently. This was very resource intensive in the mean-field case since it required a very accurate gauge fixing to Landau gauge and the self-consistent iteration is in the 2D space of $`(u_s,u_t)`$ which requires additional effort. Also, for given $`\chi _0`$ the tadpoles showed an $`L`$ dependence and hence required an extrapolation to $`L=\mathrm{}`$. This in turn requires accurate data for a good fit. In contrast, the plaquette tadpole was easily found and there was no discernible $`L`$-dependence.

The torelon method deduces $`\chi _R`$ from a fully correlated fit to the dispersion relation. Statistical errors are produced by the fit. The sideways potential approach requires some method for matching the coarse and fine potentials and the error analysis is more complex because of both systematic and statistical errors and because the signal rapidly decreases as the loop size increases. Two methods were used to extract $`\chi _R`$. The first compared the loop predictions for the coarse and fine potentials (3.2.1) and deduced $`\chi _R`$ from the rescaling of $`t`$ needed for them to agree. This method was applied to actions with mean-field improvement. Typical loops were of edge length 2-3 in units of the coarse spacing $`a_s`$. There is a clear difference in results compared with the torelon computation of about $`34\%`$. In the second approach the fine direction is first blocked by a scaling of $`\chi _0`$ to give a lattice that is approximately isotropic and then the potential is fitted to a standard form (3.2.2) and $`\chi _R`$ deduced by requiring that the linear slopes agree. The latter method has the advantage that it depends only on the long range structure of the potential and excludes the short range coulomb part. However, it does require large computing resources to extract a reliable signal at large separations, in this case up to 7 lattice units. This technique was applied mainly to actions with plaquette improvement. Although less extensively investigated, the results are consistent with the torelon computation.

The observed discrepancies can be easily attributed to residual lattice artifact effects. For the torelon we expect $`\delta \chi _R/\chi _RC\alpha _sa_s^2/R^2`$ which can be sizeable although $`C`$ is unknown. However, it should also be remarked that neither the torelon nor the second sideways potential method include the coulomb part of the static potential. Because of its short-range nature this part will be the most sensitive to the effect of lattice artifacts. In this case agreement of prediction for $`\chi _R`$ is good. This may merely indicate that $`\chi _R`$ is being deduced from the properties of the flux tube in both cases and so they should agree in that artifact effects are the same. However, it is encouraging that different methods based on the long-range properties of the action are consistent and it is worth noting that the torelon flux tube was 3–5 lattice units in length whereas it was necessary to extend the potential method to 7 lattice units. Simulation times are correspondingly reduced for the torelon. Unfortunately, it was not possible to apply the second method to data gathered using the mean-field action because of limitations in resources.

We do not find any evidence that the mean-link tuned action is superior to the plaquette tuned action in suppressing lattice artifacts, and it should be noted that the computation required to determine the mean-link tadpoles for the former case is very time consuming. However, there is other evidence that the mean-link action is superior. It is known to give smaller scaling errors in the NRQCD charmonium hyperfine splitting , better hadron mass scaling , and to give a clover coefficient for SW fermions (with the Wilson gauge action) that agrees more closely  with the non-perturbatively determined value . Using our configurations we present evidence that it also gives better rotational invariance for the static quark potential. In Fig. 7 we show the deviation of the potential from a fit to the standard form $`V(r)=b+c/r+\sigma r`$. Comparing the mean-link TI data with the plaquette TI data, it is clear that mean-link TI gives smaller deviations from the rotationally invariant fit. In table 9 we present some of the features of Fig. 7 in a quantitative form. As expected, mean-link TI gives better rotational invariance at $`r=\sqrt{3}`$, and it causes the two $`r=3`$ potentials ($`(2,2,1)`$ vs $`(3,0,0)`$) to agree better.

As mentioned in the introduction, there are many applications for anisotropic lattices. It would be valuable to measure anisotropies for a wider variety of lattice spacings and bare anisotropies, using the methods we have investigated. It would also be very useful to calculate perturbative formulae for the anisotropy and lattice spacing as a function of the bare parameters of the action, $`\beta `$ and $`\chi _0`$ . Finally, it will be necessary to repeat this tuning process for improved light quark actions. Some of the requisite perturbative calculations have already been performed .

## 6 Acknowledgements

This work is supported in part by funds provided by the U.S. Department of Energy (D.O.E.) under cooperative research agreement #DF-FC02-94ER40818. Part of the calculation was performed on the SP-2 at the Cornell Theory Center, which receives funding from Cornell University, New York State, federal agencies, and corporate partners. Calculations were carried out on the Hitachi SR2201 computers at the University of Cambridge High Performance Computing Facility and the Tokyo University Computer Centre.