text:
# S-matrix for a spin 1/2 particle in a Coulomb +Scalar potential

## Abstract

The S-matrix for a spin 1/2 particle in the presence of a potential which is the sum of the Coulomb potential $`V_c=A_1/r`$ and a Lorentz scalar potential $`V_s=A_2/r`$ is calculated.

PACS numbers:03.65.Pm,11.80.-m

e-mail:vaidya@if.ufrj.br

1.Introduction

We consider a spin $`\frac{1}{2}`$ particle in a potential which is the sum of the Coulomb potential $`V_c=A_1/r`$ and a Lorentz scalar potential $`V_s=A_2/r`$.The scalar potential is added to the mass term in the Dirac equation and may be interpreted as an effective position dependent mass.If the scalar potential is assumed to be created by the exchange of massless scalar mesons, it has the form $`V_S=\frac{A_2}{r}`$.

Exact solutions for the bound states in this mixed potential can be obtained by separation of variables . In this paper we consider the scattering problem for such a potential which does not seem to have been treated in the literature. We calculate the phase shifts by the conventional technique and show how the scattering problem can also be solved algebraically.

This paper is organized as follows.In section 2 we separate variables in the Dirac equation obtaining the radial equations,in section 3 we solve the radial equations for the scattering problem and calculate the phase shifts,in section 4 we apply an algebraic technique to obtain the phase shifts.Section 5 contains the conclusions.

2. Separation of variables in the Dirac equation

The time independent Dirac equation in the presence of the mixed potential may be written as

$`(𝜶𝐩+\beta (M{\displaystyle \frac{A_2}{r}})(E+{\displaystyle \frac{A_1}{r}}))\mathrm{\Psi }(𝐱)=\mathrm{𝟎}`$ (1)
where $`𝐩=i\frac{}{𝐱}`$ and $`r=|𝐱|`$. To separate variables we write $`\mathrm{\Psi }(𝐱)`$ in terms of two component spinors

$$\mathrm{\Psi }=\left(\begin{array}{c}\varphi \\ \chi \end{array}\right).$$
(2)
The two component angular solutions are eigenfunctions of $`J^2,J_z,L^2,S^2`$ and are of two types

$$\varphi _{}^{(+)}{}_{j,m}{}^{}=\left(\begin{array}{c}(\frac{l+\frac{1}{2}+m}{2l+1})^{\frac{1}{2}}Y_{l,m\frac{1}{2}}\\ (\frac{l+\frac{1}{2}m}{2l+1})^{\frac{1}{2}}Y_{l,m+\frac{1}{2}}\end{array}\right),$$
(3)
for $`j=l+\frac{1}{2}`$ and

$$\varphi _{}^{()}{}_{j,m}{}^{}=\left(\begin{array}{c}(\frac{l+\frac{1}{2}m}{2l+1})^{\frac{1}{2}}Y_{l,m\frac{1}{2}}\\ (\frac{l+\frac{1}{2}+m}{2l+1})^{\frac{1}{2}}Y_{l,m+\frac{1}{2}}\end{array}\right),$$
(4)
for $`j=l\frac{1}{2}`$.

In the above basis one can verify that

$$J^2\varphi _{}^{(\pm )}{}_{j,m}{}^{}=j(j+1)\varphi _{}^{(\pm )}{}_{j,m}{}^{}.$$
(5)
$$(1+𝝈𝐋)\varphi _{}^{(\pm )}{}_{j,m}{}^{}=\kappa \varphi _{}^{(\pm )}{}_{j,m}{}^{},$$
(6)
where $`\kappa =\pm (j+\frac{1}{2})`$ for $`j=l\frac{1}{2}`$. Next,we put

$$\mathrm{\Psi }=\left(\begin{array}{c}\frac{iG_{lj}}{r}\varphi _{}^{l}{}_{jm}{}^{}\\ \frac{F_{lj}}{r}𝝈𝐧\varphi _{}^{l}{}_{jm}{}^{}\end{array}\right),$$
(7)
where $`𝐧=\frac{𝐫}{r}`$ and

$$\varphi _{}^{l}{}_{jm}{}^{}=\varphi _{}^{(\pm )}{}_{j,m}{}^{},$$
(8)
for $`j=l\pm \frac{1}{2}`$. Defining the Dirac operator

$$K=\gamma ^0(1+𝚺𝐋)$$
(9)
we have

$$K\mathrm{\Psi }=\kappa \mathrm{\Psi }$$
(10)
where we have used the relation

$$[1+𝝈𝐋,𝝈𝐧]_+=0.$$
(11)
Next,using the relations

$$𝝈𝐩\frac{f(r)}{r}\varphi _{}^{l}{}_{jm}{}^{}=\frac{i}{r}(\frac{df}{dr}+\frac{\kappa f}{r})𝝈𝐧\varphi _{}^{l}{}_{jm}{}^{},$$
(12)
and

$$𝝈𝐩𝝈𝐧\frac{f(r)}{r}\varphi _{}^{l}{}_{jm}{}^{}=\frac{i}{r}(\frac{df}{dr}\frac{\kappa f}{r})\varphi _{}^{l}{}_{jm}{}^{},$$
(13)
we get the radial equations

$`{\displaystyle \frac{dG_{lj}}{dr}}+{\displaystyle \frac{\kappa }{r}}G_{lj}(E+M{\displaystyle \frac{A_2}{r}}+{\displaystyle \frac{A_1}{r}})F_{lj}`$ $`=`$ $`0,`$
$`{\displaystyle \frac{dF_{lj}}{dr}}{\displaystyle \frac{\kappa }{r}}F_{lj}+(EM+{\displaystyle \frac{A_2}{r}}+{\displaystyle \frac{A_1}{r}})G_{lj}`$ $`=`$ $`0.`$ (14)
Next, let $`\rho =kr`$ where $`k^2=E^2M^2`$ Then the radial equations take the form (we omit the indices $`l,j`$)

$`{\displaystyle \frac{dG}{d\rho }}`$ $`=`$ $`({\displaystyle \frac{E+M}{k}}+{\displaystyle \frac{A_1A_2}{\rho }})F{\displaystyle \frac{\kappa }{\rho }}G,`$
$`{\displaystyle \frac{dF}{d\rho }}`$ $`=`$ $`({\displaystyle \frac{EM}{k}}+{\displaystyle \frac{A_1+A_2}{\rho }})G+{\displaystyle \frac{\kappa }{\rho }}F.`$ (15)
3. Direct calculation of the phase shifts

In this section we solve the radial equation by using a technique similar to that used by Lin . Defining the functions $`u(\rho )`$ and $`v(\rho )`$ by

$`G`$ $`=`$ $`{\displaystyle \frac{1}{2}}e^{i\rho }(u+v),`$
$`F`$ $`=`$ $`{\displaystyle \frac{i}{2}}({\displaystyle \frac{EM}{E+M}})^{\frac{1}{2}}e^{i\rho }(uv).`$ (16)
we get

$`{\displaystyle \frac{du}{d\rho }}`$ $`=`$ $`{\displaystyle \frac{iE\gamma }{k\rho }}u+{\displaystyle \frac{1}{\rho }}(\kappa +{\displaystyle \frac{iM\gamma ^{}}{k}})v,`$
$`{\displaystyle \frac{dv}{d\rho }}`$ $`=`$ $`({\displaystyle \frac{\kappa }{\rho }}+i{\displaystyle \frac{M\gamma ^{}}{k\rho }})ui{\displaystyle \frac{\gamma E}{k\rho }}v2iv,`$ (17)
where $`\gamma =A_1+\frac{M}{E}A_2`$ and $`\gamma ^{}=A_1+\frac{E}{M}A_2`$. Eliminating $`v`$ we get

$$\rho \frac{d^2u}{d\rho ^2}+(1+2i\rho )\frac{du}{d\rho }+(\frac{2E\gamma }{k}\frac{\lambda ^2}{\rho })u=0,$$
(18)
where $`\lambda =\frac{\kappa }{|\kappa |}(\kappa ^2A_{1}^{}{}_{}{}^{2}+A_{2}^{}{}_{}{}^{2})^{\frac{1}{2}}`$. Let

$$u(\rho )=\rho ^{|\lambda |}w(\rho ).$$
(19)
Then $`w(\rho )`$ satisfies the equation

$$\rho \frac{d^2w}{d\rho ^2}+(2|\lambda |+1+2i\rho )\frac{dw}{d\rho }+2i(|\lambda |\frac{iE\gamma }{k})w=0.$$
(20)
Restoring the indices $`l,j`$ we get

$$u_{lj}=a_{lj}\rho ^|\lambda |\mathrm{\Phi }(|\lambda |\frac{iE\gamma }{k},2|\lambda |+1,2i\rho ),$$
(21)
where $`a_{lj}`$ is a constant and $`\mathrm{\Phi }(a,b,z)`$ is the confluent hypergeometric function. Using equation (17)it is easy to show that

$$v_{lj}=a_{lj}\frac{|\lambda |\frac{iE\gamma }{k}}{\kappa +\frac{iM\gamma ^{}}{k}}\rho ^{|\lambda |}\mathrm{\Phi }(|\lambda |\frac{iE\gamma }{k}+1,2|\lambda |+1,2i\rho ).$$
(22)
For $`r\mathrm{}`$ we get

$`G^{out}`$ $`=`$ $`{\displaystyle \frac{a_{lj}}{2^{|\lambda |+1}}}e^{\frac{E\gamma \pi }{2k}}e^{ikr+i\frac{E\gamma }{k}\mathrm{ln}2kr}{\displaystyle \frac{\mathrm{\Gamma }(2|\lambda |+1)}{\mathrm{\Gamma }(1+|\lambda |+i\frac{E\gamma }{k})}},`$
$`G^{in}`$ $`=`$ $`{\displaystyle \frac{a_{lj}}{2^{|\lambda |+1}}}e^{\frac{E\gamma \pi }{2k}}e^{ikri\frac{E\gamma }{k}\mathrm{ln}2kr}{\displaystyle \frac{\mathrm{\Gamma }(2|\lambda |+1)}{\mathrm{\Gamma }(1+|\lambda |i\frac{E\gamma }{k})}}\times `$ (23)
$`e^{i\pi |\lambda |}`$ $`{\displaystyle \frac{|\lambda |i\frac{E\gamma }{k}}{\kappa +i\frac{M\gamma ^{}}{k}}}.`$
Thus the phase shifts are given by

$$e^{2i\delta _{l,j}(k)}=\frac{\kappa +i\frac{M\gamma ^{}}{k}}{|\lambda |i\frac{E\gamma }{k}}\frac{\mathrm{\Gamma }(1+|\lambda |i\frac{E\gamma }{k})}{\mathrm{\Gamma }(1+|\lambda |+i\frac{E\gamma }{k})}e^{i\pi \lambda }.$$
(24)
4. Algebraic calculation of the phase shifts

In this section we apply an algebraic technique to calculate the phase shifts.This technique was used earlier by Alhassid, Gürsey and Iachello  for the relativistic Coulomb problem.

With $`\mathrm{\Phi }=\left(\begin{array}{c}G_{lj}\\ F_{lj}\end{array}\right)`$ we get

$`[{\displaystyle \frac{d}{dr}}+{\displaystyle \frac{1}{r}}(\kappa \rho _3+A_2\rho _1iA_1\rho _2)M\rho _1iE\rho _2]\mathrm{\Phi }=0,`$ (25)
where $`\rho _i`$ are the Pauli matrices.

The potential matrix $`\mathrm{\Lambda }=A_2\rho _1iA_1\rho _2+\kappa \rho _3`$ may be diagonalized by using the result

$$e^{i\beta \rho _2}e^{\alpha \rho _1}\mathrm{\Lambda }e^{\alpha \rho _1}e^{i\beta \rho _2}=\lambda \rho _3,$$
(26)
where

$`\mathrm{tanh}2\alpha `$ $`=`$ $`{\displaystyle \frac{A_1}{\kappa }},`$
$`\mathrm{tan}2\beta `$ $`=`$ $`{\displaystyle \frac{A_2}{\kappa ^{}}},`$
$`\kappa ^{}`$ $`=`$ $`ϵ(\kappa )(\kappa ^2A_{1}^{}{}_{}{}^{2})^{\frac{1}{2}},`$
$`\lambda `$ $`=`$ $`ϵ(\kappa )(\kappa _{}^{}{}_{}{}^{2}+A_{2}^{}{}_{}{}^{2})^{\frac{1}{2}},`$ (27)
where $`ϵ(\kappa )=\frac{\kappa }{|\kappa |}`$. Defining

$$\mathrm{\Phi }=e^{\alpha \rho _1}e^{i\beta \rho _2}\widehat{\mathrm{\Phi }},$$
(28)
we get

$`[{\displaystyle \frac{d}{dr}}+{\displaystyle \frac{\lambda }{r}}\rho _3+i𝐤𝝆]\widehat{\mathrm{\Phi }}=0,`$ (29)
where

$`ik_1`$ $`=`$ $`{\displaystyle \frac{EA_1A_2M\kappa _{}^{}{}_{}{}^{2}}{\kappa ^{}\lambda }},`$
$`ik_2`$ $`=`$ $`iE{\displaystyle \frac{\kappa }{\kappa ^{}}},`$
$`ik_3`$ $`=`$ $`{\displaystyle \frac{(MA_2+EA_1)}{\lambda }}.`$ (30)
One can verify that

$$E^2=𝐤^2+M^2.$$
(31)
Multiplying on the left by the operator $`\frac{d}{dr}\frac{\lambda }{r}\rho _3i𝐤𝝆`$ we get the second order equation

$$[\frac{d^2}{dr^2}\frac{\lambda }{r^2}\rho _3\frac{\lambda ^2}{r^2}+\frac{2E\gamma }{r}+k^2]\widehat{\mathrm{\Phi }}=0,$$
(32)
where

$$\gamma =A_1+\frac{M}{E}A_2.$$
(33)
For $`r\mathrm{}`$ omitting the $`\frac{1}{r}`$ terms in equation (25) we have two types of free particle solutions which satisfy the conditions

$`(M\rho _1+iE\rho _2)\mathrm{\Phi }^{in}`$ $`=`$ $`ik\mathrm{\Phi }^{in},`$
$`(M\rho _1+iE\rho _2)\mathrm{\Phi }^{out}`$ $`=`$ $`ik\mathrm{\Phi }^{out},`$ (34)
where $`k^2=E^2M^2`$.The explicit form of these solutions is

$`\mathrm{\Phi }^{in}`$ $`=`$ $`g^{in}\left(\begin{array}{c}1\\ \frac{ik}{E+M}\end{array}\right)e^{ikr},`$
$`\mathrm{\Phi }^{out}`$ $`=`$ $`g^{out}\left(\begin{array}{c}1\\ \frac{ik}{E+M}\end{array}\right)e^{ikr}.`$ (35)
To calculate the scattering matrix we need to find the ratio $`\frac{g^{out}}{g^{in}}`$.

Next,equations (29)and (30) give,in the limit $`r\mathrm{}`$

$$\frac{g^{out}}{g^{in}}=\frac{\widehat{g}^{out}}{\widehat{g}^{in}}\frac{ik_{}S_{11}ik_3S_{12}ikS_{12}}{ik_{}S_{11}ik_3S_{12}+ikS_{12}},$$
(36)
where $`S=e^{\alpha \rho _1}e^{i\beta \rho _2}`$. After considerable effort we get

$`(ik_{}S_{11}ik_3S_{12})`$ $`=`$ $`{\displaystyle \frac{E+M}{A_1A_2}}(\kappa +\lambda ){\displaystyle \frac{(A_1(\kappa ^{}+\lambda )A_2(\kappa ^{}+\kappa ))}{4\kappa ^{}\lambda \mathrm{cosh}\alpha ,\mathrm{cos}\beta }}`$
$`ikS_{12}`$ $`=`$ $`ik{\displaystyle \frac{(A_1(\kappa ^{}+\lambda )A_2(\kappa ^{}+\kappa ))}{4\kappa ^{}\lambda \mathrm{cosh}\alpha \mathrm{cos}\beta }}.`$ (37)
Thus we get

$$\frac{ik_{}S_{11}ik_3S_{12}ikS_{12}}{ik_{}S_{11}ik_3S_{12}+ikS_{12}}=\frac{\kappa +\lambda +\frac{i}{k}(E\gamma M\gamma ^{})}{\kappa +\lambda \frac{i}{k}(E\gamma M\gamma ^{})},$$
(38)
where we have used

$$E\gamma M\gamma ^{}=(EM)(A_1A_2).$$
(39)
Finally, using the relation

$$\frac{\kappa i\frac{M\gamma ^{}}{k}}{\lambda i\frac{E\gamma }{k}}=\frac{\lambda +i\frac{E\gamma }{k}}{\kappa +i\frac{M\gamma ^{}}{k}},$$
(40)
we get

$$\frac{g^{out}}{g^{in}}=\frac{\widehat{g}^{out}}{\widehat{g}^{in}}\frac{\kappa i\frac{M\gamma ^{}}{k}}{\lambda i\frac{E\gamma }{k}}.$$
(41)
One may calculate the phase shifts by an algebraic technique by observing that the differential equation for $`\widehat{g}`$ obtained from equation (32) is of the same form as the radial equation for the two-dimensional hydrogen atom with the Hamiltonian

$$H=\frac{𝐩^2}{2}\frac{\alpha }{r}.$$
(42)
The latter problem has been treated by Alhassid,Engel and Wu .They use the fact that the problem has the dynamical symmetry under the $`SO(2,1)`$ group so that the asymptotic states obey recursion relations under the action of the generators of the $`SO(2,1)`$ algebra.On the other hand, the asymptotic states are linear superposition of eigenstates of the Casimir operator of the algebra of the Euclidean group $`E(2)`$.As a consequence one gets a recursion relation for the ratio of the out and in amplitudes which can be solved. Making appropriate changes in their results:$`\alpha 2E\gamma `$ and $`m\lambda +\frac{1}{2}`$, we get

$$\frac{\widehat{g}^{out}}{\widehat{g}^{in}}=e^{i\pi \lambda }\frac{\mathrm{\Gamma }(1+\lambda i\frac{E\gamma }{k})}{\mathrm{\Gamma }(1+\lambda +i\frac{E\gamma }{k})}.$$
(43)
Using equation (43) in equation (41) we have the same phase shifts as those of equation (24).

5. Conclusion

In conclusion we have shown that the phase shifts for generalized Dirac equation can be found algebraically. The results are in agreement with the conventional calculation.The advantage of the algebraic method is that we work in the $`r\mathrm{}`$ limit instead of solving for the radial functions and then taking the limit.

References

1. W. Greiner, Relativistic Quantum Mechanics,Springer Verlag,New York,(1990).

2 S. I. Ikhadir, O.Mustafa,R. Sever,Hadronic J.16, 57,(1993).

3. Q.G.Lin Phys. Lett. A 260(1999) 17

4. Y. Alhassid,F. Gürsey and F. Iachello J. Phys.A:Math. Gen.

22(1989) L947

5.Y. Alhassid,J. Engel and J. Wu Phys. Rev. Lett. 53(1984) 17


format: text

text:
# Optimal Cooperation and Submodularity for Computing Potts’ Partition Functions with a Large Number of States

## Introduction

Spin models with discrete symmetry are often used to describe order-disorder transitions, such as in magnetic systems or in absorbed monolayers, etc, for a review see in. In two-dimensional regular lattices and with homogeneous couplings there are some models (Ising model, critical Potts model, eight-vertex model, etc.) with an asymptotically exact solution for the partition function in the thermodynamic limit. In the presence of quenched, that is time independent, disorder this type of exact solutions are scarce which has greatly hampered our understanding about collective phenomena in disordered systems, for example about the properties of spin glasses.

For disordered systems, such as for spin models with random couplings and/or fields one usually performs an exact computation on finite lattices and the results are then extrapolated. An important limitation of this finite lattice method is the available sizes, which can be treated by different computational methods. For most of the discrete spin models a computation of the ground state represents a hard optimization problem when there is frustration and it is even much more complex to compute the partition function.

In this respect the Ising model, which is equivalent to the Potts model with $`q=2`$ states, is an exceptional problem since the ground (minimum energy) states can be found in polynomial time using the solution of max cut problems in planar graphs that involves matching theory   . The full partition function of the Ising problem can also be computed for lattices of bounded genus . Another example of application of optimization theory to lattice statistical mechanics is the Random Field Ising Model where the ground state can be found in polynomial time in any dimension  . A discussion about optimal cuts in statistical mechanics can be found in  and a general review about the application of combinatorial optimization methods in statistical physics can be found in.

For the $`q=3`$-state Potts model, when the number of states is three, the computation of the minimum energy of the system involves (NP-)hard combinatorial problems such as minimizing the 3-cut in a graph, and no efficient procedure is known for computing the partition function.

If the number of states can be arbitrary large (and given as part of the input) the partition function has an exponential number of terms, and the problem may look hopeless. In this paper we show that this is not the case: when the number of states tends to infinity, the partition function can be determined up to an integer multiplicative factor, and some average properties of the system can be inferred. This evaluation goes trough the minimization of a submodular function. Note that in  a simulated annealing method and an approximate combinatorial method have been used. The solution we give here is both exact and much more efficient.

Even though there are ready-made algorithms for general submodular function minimization    , the resources they require are far beyond the possibilities of the present application. While the existing first implementations of general submodular function minimization may require too much time already to solve the optimization problem on a lattice containing a dozen of sites (vertices of the input graph corresponding to the variables of the general problem), we want and will solve problems involving more than two hundred thousand ‘sites’. For this we use methods inspired by     that fit the particularity of the application and develop an algorithm which is tailored to graphic submodular functions. Besides the computational advantages of this approach, we are also aiming at conceptual simplicity and elegance. Readers who seek for more mathematical background and connections can consulte .

In Section 1 we describe the problem of physics and show how it reduces to a function maximization. In Section 2 we formulate the mathematical problem, the preliminaries and the broader context of our algorithm; in Section 3 and 4 we show how to construct an optimal solution and in Section 5 we state the algorithm. Finally, in Section 6 we analyze the practical implementation of the method and display a few computational results. The paper is self-contained, the algorithms and their correctness will be fully proved.

## 1 The problem in terms of physics.

For a review of the huge amount of work devoted to the Potts model see . In this section we recall the Potts model and show why its partition function (up to an integer multiplicative factor) is determined, if the number of ‘states’ tends to infinity by the minimum value of a submodular function.

A lattice of spins is given at each site of which is a variable. We denote the number of sites by $`n`$. Each variable $`\sigma _1,\mathrm{},\sigma _n`$ can take values in $`\{0,1\mathrm{},q1\}`$ for a given integer $`q`$. Pairs of neighboring sites on the lattice are called bonds, the number of bonds by $`m`$.

Each configuration $`\sigma =(\sigma _1,\mathrm{},\sigma _n)`$ has an energy

$$E(\sigma )=\underset{ij}{}K_{ij}\delta _{\sigma _i\sigma _j},$$
(1)
where $`\sigma _i\{0,\mathrm{},q1\}`$ the sum runs over all bonds $`ij`$, $`K_{ij}\text{IR}_+`$ is a given non-negative weight of bond $`ij`$, and $`\delta _{ab}`$ is the Kronecker symbol ($`1`$ if $`a=b`$ and $`0`$ otherwise). For mathematicians: a lattice is a graph, sites are vertices, bonds are edges. The aim is to compute or approximate the partition function

$$Z(K_{ij}:ij\text{ is a bond})=\underset{\sigma }{}\mathrm{exp}(E(\sigma ))=\underset{\sigma }{}\underset{ij}{}e^{K_{ij}\delta _{\sigma _i\sigma _j}},$$
(2)
where the summation runs over all assignments of values $`\sigma =(\sigma _1,\mathrm{},\sigma _n)`$, and the product over all bonds $`ij`$ (the inverse temperature has been included in the couplings $`K_{ij}`$).

We follow . Note that $`e^{K\delta }=1+(e^K1)\delta `$ for $`\delta \{0,1\}`$ ; introduce $`\nu _{ij}=e^{K_{ij}}1`$ and expand the product of sums:

$$\underset{ij}{}e^{K_{ij}\delta _{\sigma _i\sigma _j}}=\underset{ij}{}(1+(e^{K_{ij}}1)\delta _{\sigma _i\sigma _j})=\underset{F}{}\underset{ijF}{}\nu _{ij}\delta _{\sigma _i\sigma _j},$$

where the summation runs over all subsets of bonds $`F`$. Substituting this to (2) we get:

$$Z(K_{ij})=\underset{\sigma }{}\underset{F}{}\underset{ijF}{}\nu _{ij}\delta _{\sigma _i\sigma _j}$$

$$=\underset{F}{}\underset{\sigma }{}\underset{ijF}{}\nu _{ij}\delta _{\sigma _i\sigma _j}=\underset{F}{}q^{c(F)}\underset{ijF}{}\nu _{ij},$$
(3)
where the sum runs over all subsets $`F`$ of bonds, and where $`c(F)`$ is the number of connected components of $`F`$ on the set of all sites, counting also the isolated sites among the components. (By convention if $`F`$ is empty then $`_{ijF}\nu _{ij}=1`$.) We got the last equality by counting the number of different $`\sigma =(\sigma _1,\mathrm{},\sigma _n)`$ for which the product is nonzero; since it is nonzero if and only if $`\sigma (i)=\sigma (j)`$ for every $`ijF`$, that is, if and only if $`\sigma `$ is constant on every connected component of $`F`$, all possible $`\sigma `$ can be enumerated by choosing an element of $`Z_q`$ for every connected component of $`F`$ independently. Therefore we have exactly $`q^{c(F)}`$ such $`\sigma `$.

Note that this sum has $`2^m`$ terms, and that in (3) $`q`$ does not need to have an integer value. (This provides a way of defining Potts’ model for non-integer values .) Clearly, $`K_{ij}0`$ is equivalent to $`\nu _{ij}0`$, and we can introduce a new set of variables $`\alpha _{ij}`$ with

$$\nu _{ij}=q^{\alpha _{ij}},$$

and the partition function becomes

$$Z=\underset{F}{}q^{c(F)+_{ijF}\alpha _{ij}}.$$
(4)
Finally one introduces the function

$$f(F)=c(F)+\underset{ijF}{}\alpha _{ij}$$
(5)
so that

$$Z=\underset{F}{}q^{f(F)}.$$
(6)
As pointed out in ref , while $`q`$ tends to infinity the sum in (6) is asymptotically equal to $`Z=Nq^f^{}`$ where $`f^{}=max_{F\{1,\mathrm{},n\}}f(F)`$ and $`N`$ is the number of optimal sets. Note that if the weights $`w_{ij}`$ are arbitrary reals, the degeneracy $`N`$ is likely to be $`1`$.

We will analyze the combinatorial optimization problem of finding this maximum. The problem of finding the degeneracy of the optimal set will not be addressed. In the analysis of the Potts model magnetic properties are related to geometrical properties of the optimal sets $`A`$, and specifically to their fractal dimension (if any). We assume that, on average and for large enough lattice, this fractal dimension is the same for all sets, and therefore finding a single optimal set is sufficient for our purpose. We will present methods which apply to any set of weights $`\{\alpha _{ij}\}`$ in an arbitrary lattice (graph) and provide a simply stated and efficient algorithm.

## 2 The problem in terms of graphs.

For basic graph theoretic notions, terminology and notation we refer to . Given a graph $`G=(V,E)`$ on $`n`$ vertices with weight function $`w:E(G)\text{IR}`$, we want to solve the following problem:

$$\text{Maximize}\{f_{G,w}(A)=c_G(A)+\underset{eA}{}w(e):AE(G)\}$$
$`(POTTS),`$
where $`c_G(A)`$ is the number of connected components of the graph $`G(A)=(V,A)`$. Note a slight change of notation with Eq. 5 where an edge is denoted by $`ij`$ the two corresponding adjacent vertices. We can suppose $`0<w(e)<1`$ for all $`eE(G)`$, because if $`w(e)=0`$, then deleting the edge, if $`w(e)1`$, then contracting it (identifying its two endpoints) we get an equivalent problem.

When no confusion is possible we will simply use $`f`$ or $`f_G`$ for $`f_{G,w}`$, $`c`$ for $`c_G`$. This function that has to be maximized in order to solve (POTTS) has a crucial and well-known property (see for instance , Exercise 6.2). For any two subsets $`A`$ and $`B`$ of edges of $`G`$ :

$$f(A)+f(B)f(AB)+f(AB).$$
$`(SUPER)`$
Indeed, it is clear that $`_{eA}w(e)+_{eB}w(e)=_{eAB}w(e)+_{eAB}w(e)`$, so we only have to show that $`c(A)+c(B)c(AB)+c(AB)`$; that is $`c`$ satisfies the (SUPER) property. For the sake of completeness we include a proof. Proceed by induction on $`|AB|+|BA|`$. If this number is 0 then the statment is obvious. If not, then there exists say $`eAB`$. By induction $`c(Ae)+c(B)c(ABe)+c(AB)`$. Deleting an edge increases the number of connected components by at most 1 and so we are done unless the equality holds and $`c(Ae)=c(A)`$ and $`c(ABe)=c(AB)+1`$. But this is impossible since $`c(Ae)=c(A)`$ means that the two endpoints of $`e`$ are in the same connected component of $`G(Ae)`$, and then also in the same component of $`G(ABe)`$.

A function $`g`$ defined on the set of subsets of a set $`S`$, with values in IR is said to be submodular if $`g(A)+g(B)g(AB)+g(AB)`$ for all subsets $`A`$ and $`B`$ of $`S`$. The (SUPER) property shows that the function $`f`$, whose minimization solves (POTTS), is submodular.

By now several polynomial algorithms     are known to minimize a submodular function, therefore we can already claim that there exists a polynomial algorithm solving our problem (POTTS). However, these general algorithms are not efficient enough for our instances. Fortunately, the specific properties of our problem make possible the use of a considerably simpler and quicker algorithm. Several more specialized algorithms were already known (see  for an extended presentation) and in particular Cunningham  worked out a specialized algorithm to the “optimal attack problem” :“The weight of an edge represents the effort required by an attacker to destroy the edge, and the attacker derives a benefit for each new component created by destroying edges $`\mathrm{}`$”. (POTTS) is equivalent to the “optimal attack problem” but with a ‘complementary’ viewpoint : the weight of an edge represents the benefits of cooperation between two vertices (say, researchers), furthermore there is a unit support for each component (say, for each research project). So there is a loss of support when two components unite. Cooperate optimally!

The algorithm described below has its roots in  and is similar to , . It has the best asymptotic worst case complexity, the same as , and is probably the simplest, both conceptually and in the use of computational resources of the implementation. Sophisticated ingredients were necessary for finding this solution, but we don’t need to make explicit use of them. The interested reader may find more details on that subject in .

We return now to the (POTTS) problem. An important observation is the following: if $`A`$ is a set of edges, and $`X`$ is the vertex set of a connected component of $`G(A)`$, then adding to $`A`$ edges induced by $`X`$, that is with both endpoints in $`X`$, increases the value of $`f`$. Thus in an optimal solution $`A^{}`$ of (POTTS) each connected component of $`A^{}`$ contains all the edges of $`G`$ it induces. Let us call Potts partition any partition of the vertices such that each class induces a connected subgraph of $`G`$.

All we have to determine is : a Potts partition maximizing the total weight of edges induced by the classes plus the number of classes.

As a consequence, following , during all the procedure it will be sufficient to consider subsets of vertices instead of subsets of edges.

A set $`F`$ of edges optimizing (POTTS) will be called optimal, its associated Potts partition $`𝒫_G(F)=\{X_1,\mathrm{}X_k\}`$ will be called optimal as well. For any subset $`F`$ of edges let $`w(F)=_{eF}w(e)`$ and for any subset $`X`$ of vertices let $`w(X)=_{eE(X)}w(e)`$ where $`E(X)`$ is the set of edges with both extremities in $`X`$.

## 3 Extension of a solution

In this section we will see how any set of edges maximizing (POTTS) for a graph $`G`$ minus a vertex, may be completed so that it maximizes (POTTS) for the graph $`G`$ itself.

Let $`G`$ be a graph with a weight function $`w`$ on the edges, let $`x`$ be a vertex of $`G`$ and let $`G^{}`$ be the subgraph of $`G`$ obtained by deleting $`x`$. The edges of $`G^{}`$ keep their weight. Let $`F^{}E(G^{})`$ be an optimal solution of (POTTS) for $`G^{}`$.

The following Lemma is an easy consequence of the (SUPER) property.

###### Lemma 1

There exists an optimal solution of (POTTS) for $`G`$ which contains all edges of $`F^{}`$.

Proof: Let $`F`$ be any optimal solution of (POTTS) for $`G`$. By the (SUPER) property one has in $`G`$:

$$f_G(FF^{})f_G(F)+f_G(F^{})f_G(FF^{}).$$

Both $`F^{}`$ and $`FF^{}`$ are subsets of edges of $`G^{}`$ and since $`F^{}`$ is optimal in $`G^{}`$, $`f_G(F^{})=1+f_G^{}(F^{})1+f_G^{}(FF^{})=f_G(FF^{})`$, whence :

$$f_G(FF^{})f_G(F).$$

Since $`F`$ is optimal in $`G`$ by hypothesis, we get that $`FF^{}`$ is optimal too, and the Lemma is proved. $`\mathrm{}`$

As a consequence of Lemma 1 we can obtain an optimal solution $`F^{}`$ for $`G`$ by adding edges to $`F^{}`$. Let $`X_1,\mathrm{},X_k`$ be the vertices of the connected components of $`G^{}(F^{})`$, then $`𝒫_G(F^{})=\{X_1,\mathrm{},X_k,\{x\}\}`$. Each connected component in $`G(F^{})`$ will be either an element of $`𝒫_G(F^{})`$ or the union of at least two elements of $`𝒫_G(F^{})`$. The next Lemma shows how the value of a Potts partition is affected when a subset of connected components are put together.

###### Lemma 2

Let $`𝒫`$ and $`𝒫^{}`$ be two Potts partitions of $`G`$ such that $`𝒫=(𝒫^{}𝒲)\{X_i;X_i𝒲\}`$ for some $`𝒲𝒫`$. Then

$$f_G(𝒫)=f_G(𝒫^{})(|𝒲|1w(E(𝒲))),$$

where $`E(𝒲)`$ denotes the set of edges of $`G`$ joining vertices belonging to two different sets in $`𝒲`$.

In particular if $`𝒫^{}`$ is optimal then $`|𝒲|1w(E(𝒲))0`$ and if $`𝒫`$ is optimal then $`|𝒲|1w(E(𝒲))0`$.

Proof: Replacing in $`𝒫`$ the sets of $`𝒲`$ by their union decreases the cardinality of $`𝒫`$ by $`|𝒲|1`$, on the other hand let $`A`$ and $`A^{}`$ be the sets of edges induced by the classes respectively of $`𝒫`$ and $`𝒫^{}`$, $`A=A^{}E(𝒲)`$ and so the weight of the edges augment by $`w(E(𝒲))`$. The equality is then proved and the rest follows directly. $`\mathrm{}`$

We are now ready to prove the following theorem:

###### Theorem 1

For any $`𝒲𝒫_G(F^{})`$ containing $`\{x\}`$ and minimizing $`|𝒲|1w(E(𝒲))`$, the set $`F^{}=F^{}E(𝒲)`$ is optimal for (POTTS) in $`G`$.

Proof: By Lemma 1 we know that there exists an optimal solution $`F^{}`$ of (POTTS) for $`G`$ which contains $`F^{}`$. Let $`𝒲𝒫_G(F^{})`$ such that $`W=_{X_i𝒲}X_i`$ is an element of $`𝒫_G(F^{})`$. If $`𝒲`$ doesn’t contain $`\{x\}`$ then, since $`F^{}`$ is optimal in $`G^{}`$, we get by Lemma 2 that $`|𝒲|1w(E(𝒲))0`$ (this value is the same in $`G`$ and $`G^{}`$). On the other hand $`F^{}`$ is optimal in $`G`$ which implies that $`|𝒲|1w(E(𝒲))=0`$, but then $`F^{}E(𝒲)`$ is also optimal for $`G`$. Hence there exists an optimal solution $`F^{}`$ of (POTTS) for $`G`$ containing $`F^{}`$ and such that any element of $`𝒫_G(F^{})`$ not containing $`x`$ is already in $`𝒫_G(F^{})`$.

So any $`𝒲𝒫_G(F^{})`$ containing $`\{x\}`$, such that $`W=_{X_i𝒲}X_i`$ induces a connected subgraph of $`G`$, and minimizing under these conditions $`|𝒲|1w(E(𝒲))`$, will provide an optimal solution $`F^{}E(𝒲)`$. (Notice that this minimum is $`0`$ since for $`𝒲=\{\{x\}\}`$ we get $`0`$.)

In order to prove the theorem it remains only to show that we may skip the connectivity constraint : we show now that any subset $`𝒲𝒫_G(F^{})`$, $`\{x\}𝒲`$, minimizing $`|𝒲|1w(E(𝒲))`$ corresponds to a connected subgraph in $`G`$. Assume not, and let $`𝒲_1,\mathrm{},𝒲_l`$ ($`l2`$) be the subsets of $`𝒫_G(F^{})`$ corresponding to the vertices of the connected components induced by $`W=_{X_i𝒲}X_i`$ in $`G`$, $`𝒲_1`$ being the one containing $`\{x\}`$. We have

$$|𝒲|1w(E(𝒲))=|𝒲_1|+\mathrm{}+|𝒲_l|1w(E(𝒲_1))\mathrm{}w(E(𝒲_l))$$

For $`i2`$, $`\{x\}𝒲_i`$ and so as already noticed $`|𝒲_i|w(E(𝒲_i))1`$, which implies $`|𝒲|1w(E(𝒲))>|𝒲_1|1w(E(𝒲_1))`$ and so we get a contradiction to our assumption on $`𝒲`$. $`\mathrm{}`$

At this point we see that any way to find $`𝒲`$ will provide a constructive algorithm for getting an optimal solution of (POTTS) : we start with a solution for a small subgraph (for example a one vertex subgraph) and then add the vertices one by one, computing at each step a new extended solution. Fortunately an optimal $`𝒲`$ can be computed and in the next section we are going to show how.

We remark that the value $`|𝒲|1w(E(𝒲))`$ doesn’t depend on the subgraphs induced by the subsets $`X_i`$ in $`𝒲`$, so we may ignore these and work in a possibly smaller graph. To be precise : the result of shrinking the pairwise disjoint sets $`X_1,\mathrm{},X_k`$ of $`𝒫_G(F^{})`$ in $`G`$ is the graph $`\mathrm{shr}(G)=(\mathrm{shr}(V),\mathrm{shr}(E))`$, where $`\mathrm{shr}(V)=\{x,x_1,\mathrm{},x_k\}`$ where $`x_1,\mathrm{},x_k`$ are distinct new vertices, and the function $`\mathrm{shr}:V\{x,x_1,\mathrm{},x_k\}`$ is defined with $`\mathrm{shr}(v):=x_i`$ if $`vX_i`$ and $`\mathrm{shr}(x):=x`$; the image of an edge $`e`$ with extremities $`x`$ and $`y`$ is $`\mathrm{shr}(e)`$ with extremities $`\mathrm{shr}(x)\mathrm{shr}(y)`$, moreover, if this is a loop ($`\mathrm{shr}(x)=\mathrm{shr}(y))`$, it is deleted; edges keep their weight, that is $`w_{\mathrm{shr}}(\mathrm{shr}(e))=w(e)`$; sets of vertices or of edges are replaced by the image sets. There is a one to one correspondence between the subsets $`𝒲`$ of $`𝒫_G(F^{})`$ containing $`\{x\}`$ and subsets $`W`$ of vertices of $`\mathrm{shr}(G)`$ containing $`x`$ and for any such $`W=\mathrm{shr}(𝒲)`$ one has $`|𝒲|1w(E(𝒲))=|W|1w_{\mathrm{shr}}(W)`$ (as defined in the preceding chapter $`w_{\mathrm{shr}}(W)`$ is the sum of weights of the edges in the subgraph of $`\mathrm{shr}(G)`$ induced by $`W`$). So $`𝒲`$ will be optimal if and only $`W=\mathrm{shr}(𝒲)`$ minimizes $`|W|1w_{\mathrm{shr}}(W)`$.

## 4 A network flow model

Given a graph $`H`$ with weight function $`w`$ on the edges and $`WV(H)`$, let $`b(W)=|W|1w(W)`$. From the previous chapter it is clear that any way of solving the following problem will provide an algorithm solving (POTTS) : given a vertex $`x`$ in $`H`$, find $`W^{}`$ a subset of vertices of $`H`$ containing $`x`$ such that $`b(W^{})=min(b(W);WV(H),xW)`$. This problem is solved both in Picard-Queyranne  and Padberg-Wolsey , using a network flow model that will be described below. (The work ,  solve much more complex problems involving several network flow computations – the subroutine below is the adaptation of an auxiliary procedure.)

We first give some definitions and wellknown facts. Given a directed graph and a subset $`X`$ of its vertices, $`\delta (X)`$ will denote the set of arcs leaving $`X`$, and is called a cut; if $`sX`$, $`tX`$ it is an $`(s,t)`$-cut; $`\delta (X,Y)`$ denotes the set of arcs oriented from $`X`$ to $`Y`$. A function $`c`$ with nonnegative value on the arcs is called a capacity function. If $`F`$ is a set of arcs then $`c(F):=_{eF}c(e)`$. A network is a directed graph with capacity function. By a well known theorem of Ford and Fulkerson  a minimum $`(s,t)`$-cut in a network may be found by computing a maximum flow in this network. Very efficient maximum flow algorithms are known .

We now describe in several steps a network $`(D,c)=N(H,w)`$ that will be associated to $`H`$ and $`w`$:

1. $`V(D):=V(H)\{s,t\}`$ where $`s,t`$ are distinct new vertices;

2. Define for all $`uV(H)`$:

$$p(u):=1/2\underset{vV(H),uvE(H)}{}w(uv),\text{ and then }$$

* if $`p(u)>1`$ add an arc $`(s,u)`$ of capacity $`c(s,u)=p(u)1`$,
* if $`p(u)<1`$ then add an arc $`(u,t)`$ of capacity $`c(u,t)=1p(u)`$.

3. to each edge $`uv`$ of $`E(H)`$ we associate the arcs $`(u,v)`$ and $`(v,u)`$ of capacities $`c(u,v)=c(v,u)=\frac{1}{2}w(uv)`$.

Let $`WV(H)`$, $`\delta (\{s\}W)`$ is an $`(s,t)`$-cut of $`N(H,w)`$ and reciprocally any $`(s,t)`$-cut of $`N(H,w)`$ corresponds to a subset of vertices of $`H`$. There is a close relationship between the capacities of the $`(s,t)`$-cuts of $`N(H)`$ and the values of $`b`$ on subsets of vertices of $`H`$ :

$$c(\delta (\{s\}W))=|W|w(W)+K=b(W)+K+1,$$
$`(CUT)`$
where $`K:=c(\delta (s))=_{(s,x)A}c(s,x)`$.

Indeed, let us see how the capacity of the cuts changes if we start with the set $`\{s\}`$ inducing an $`(s,t)`$-cut of capacity $`c(\delta (s))=K`$ and then ‘add’ to it the vertices of $`W`$ one by one:

The contribution of adding $`v`$ to the side of $`s`$ is $`1p(v)`$, for either it decreases from $`p(v)1`$ to $`0`$ (this happens if $`p(v)>1`$ see the first case of ‘2.’ in the construction of $`N(H,w)`$), or it increases from $`0`$ to $`1p(v)`$. Thus the contribution of these arcs is:

$$\underset{vW}{}1p(v)=|W|w(W)1/2\underset{xyE(H),xW,yW}{}w(xy).$$

On the other hand the contribution of the arcs between $`W`$ and $`V(H)W`$ is clear: at the beginning it is zero, and at the end it is:

$$c(W,V(H)W)=1/2\underset{xyE(H),xW,yW}{}w(xy).$$

The change comparing to $`K`$ is provided by the sum of the two contributions, which is $`|W|w(W)`$, and so equation (CUT) has been proved.

From (CUT) we get that any $`(s,t)`$-cut $`C=\{s\}W`$ containing the special vertex $`x`$ and of minimum capacity will correspond to our goal : a subset $`W`$ of vertices of $`H`$ containing $`x`$ and minimizing $`b`$. Such a cut is easy to obtain : shrink $`\{s,x\}`$ or equivalently, add to $`N(H,w)`$ an arc of infinite capacity from $`s`$ to $`x`$, and compute in the new network a minimum $`(s,t)`$-cut : this cut will contain $`x`$ and will be minimum among $`(s,t)`$-cuts containing $`x`$.

Let us note that this method can be easily generalized to minimize any modular shift of $`w(W)`$. Adding $`|W|1`$ is just a particular choice of a modular function. That is, one can minimize $`(_{wW}m(w))w(W)`$, where $`m:V\text{IR}`$. Indeed, for generalizing from $`|W|`$ to $`m(W)`$ one only has to write $`m(x)`$ instead of $`1`$ in both cases of Step 2 in the construction of N(H,w).

## 5 The Algorithm

In this section we state the ‘Optimal Cooperation Algorithm’ whose validity is a consequence of the results of the two preceding chapters.

At each iteration a subset $`UV`$ and a partition $`𝒫`$ of $`U`$ will be at hand. In Step 0 we give trivial initial values; in Step 1, we choose an arbitrary vertex $`u`$ to be added to $`U`$ and through the following steps we compute a subset $`𝒲`$ of $`𝒫\{u\}`$ providing a partition of $`U\{u\}`$ consisting of the sets contained in $`𝒫`$ but not in $`𝒲`$ plus a new class containing all the other vertices (that is : the new class is the union of all sets in $`𝒲`$). The set of edges corresponding to this new partition is the union of the set of edges corresponding to the previous one and the edges between the sets in $`𝒲`$.

OPTIMAL COOPERATION ALGORITHM

INPUT: A graph $`G=(V,E)`$, and a weight function $`w`$: $`E(0,1)`$, where $`(0,1):=\{t\text{IR}:0<t<1\}`$.

OUTPUT: A partition $`𝒫^{}`$ optimizing (POTTS) : the set of edges $`A^{}`$ induced by the sets in $`𝒫^{}`$ maximizes $`\{f_{G,w}(A):=c_G(A)+_{eA}w(e):AE\}`$.

0. $`U:=\mathrm{}`$, $`𝒫:=\mathrm{}`$.

Do $`n`$ times consecutively Steps 1 to 5, and then define $`𝒫^{}=𝒫`$:

1. Choose a vertex $`uVU`$.

2. Define $`H`$ and $`w_H`$ as the result of shrinking the classes of $`𝒫`$ in $`G(U\{u\})`$ with weight function $`w`$ restricted to the edges of $`G(U\{u\})`$.

3. Construct $`N(H,w_H)`$, add an arc of ”infinite” capacity from $`s`$ to $`u`$ and using any maximum flow algorithm compute, in the so obtained network, a minimum $`(s,t)`$-cut defined by the set of vertices $`C=\{s,u,x_1,\mathrm{},x_k\}`$ where each $`x_i`$ is a vertex of $`H`$ corresponding to a set $`X_i𝒫`$.

4. Define $`R:=\{u\}X_1\mathrm{}X_k`$

5. Redefine $`U`$ and $`𝒫`$: $`U:=U\{u\}`$, $`𝒫:=(𝒫\{X_1,\mathrm{},X_k\})\{R\}`$;

END;

Remarks:

\- ”Infinite” means big enough to avoid a minimum cut containing this arc: for example the sum of the capacities of all arcs.

\- Note that $`C`$ can be equal to $`\{s,u\}`$.

\- The graphs $`H`$ can also be constructed iteratively with only one shrinking in each iteration, by adding $`u`$ and then shrinking the set $`C\{s\}`$.

\- The choice for $`u`$ is completely free; this freedom could be used for making the computations easy.

###### Theorem 2

The output $`𝒫^{}`$ is an optimal partition for $`(G,w)`$.

Proof: At step 0, $`𝒫:=\mathrm{}`$ is an optimal solution for the subgraph of $`G`$ induced by $`U:=\mathrm{}`$. Assume now that $`𝒫`$ is an optimal solution for the subgraph of $`G`$ induced by $`U`$, and show that after applying steps 1 to 5 the new partition is optimal for $`U\{u\}`$. From the preceding chapter we know that $`W:=C\{t\}=\{u,x_1,\mathrm{},x_k\}`$ is a subset of vertices of $`H`$ containing $`u`$ minimizing $`|W|w_H(W)1`$. But $`H`$ and $`w_H`$ are obtained by shrinking the classes of $`𝒫`$ in $`(G(U\{u\}),w)`$, hence $`𝒲=\{u,X_1,\mathrm{},X_k\}`$ is a subset of $`\{u\}𝒫`$ minimizing $`|𝒲|1w(E(𝒲))`$ and by Theorem 1 the partition $`(𝒫\{X_1,\mathrm{},X_k\})\{R\}`$ is then optimal $`\mathrm{}`$

Notice that the ”Optimal Cooperation algorithm” consists merely of $`n1`$ network flow computations and $`n1`$ shrinking.

## 6 Implementation and evaluation

In this section we discuss some aspects of the practical implementation of the optimal cooperation algorithm. Here we consider a two dimensional square lattice where the weights, $`w`$, follow a bimodal distribution:

$$P(w)=p\delta (ww_1)+(1p)\delta (ww_2)$$
(7)
where $`\delta (x)=1`$ if $`x=0`$ and $`\delta (x)=0`$ otherwise. We restrict ourself to the symmetric case $`p=1/2`$. As explained in  in the infinite lattice limit there is a second-order phase transition in the model at $`w_1+w_2=1`$. In the strongly disordered phase, $`w_1+w_21`$ the optimal set of edge is $`A=\mathrm{}`$, while in the strongly ordered phase, $`w_1+w_21`$ it is $`A=E`$. In between, the optimal set $`A`$ has a non-trivial structure and at the transition point it is a fractal reflecting the critical properties of the system.

As explained in the previous section, for a graph with $`n`$ vertices the algorithm consists in $`n1`$ steps. At each step a new site is incorporated to the current graph for which the optimization has already been performed. There is a freedom of the order of the added sites, which can be used to speed up the computation, as well as to obtain the optimal sets of intermediary lattices. For example from the computation of the optimal set for an $`L\times L`$ square lattice with free boundary conditions one can also obtain the optimal sets of $`L1`$ minors of size ranging from 1 to $`L1`$ with free boundary conditions, which is useful for finite size analysis.

The “engine” of the algorithm is the max flow solver, which is called at each step. We have used the Goldberg and Tarjan algorithm . Another important point concerns the building of the network (point 3 in the Optimal Cooperation Algorithm). It is possible to build it explicitly as another data structure than the lattice. Alternatively one can increase up to a sufficiently large value all the capacities of the edges of a spanning tree of each connected component as mentioned in the previous paragraph. In general the cpu time required to explicitly shrink the connected components is larger than the cpu time needed to simply increase some capacities. On the other hand, keeping all vertices will increase the size of the network in which the max flow has to be found. If all the capacities are large ($`w_1+w_2>1`$, that is we are in the ordered phase) the network will have only two sites at each step, and therefore the effort to build the network will be largely compensated by the speed at which the maximum flow is found. In the opposite extreme, when the optimal set tends to be empty, it is questionable which of the two strategies is the best. We have found that explicitly shrinking the connected components is always the most efficient. As a matter of fact the computation time is about ten times larger at the critical point $`w_1+w_2=1`$ than at the ordered phase ($`w_1+w_2>1`$) and about three times larger than at the disordered phase ($`w_1+w_2<1`$).

As an illustration we display in the table below some typical cpu times needed to find an optimal set on a periodic square lattice of various size, $`L`$, at the critical point $`w_1+w_2=1`$. From these results one can extract the empirical rule that in our algorithm the computation time grows approximately as the fourth power of the linear size $`L`$, thus as the second power of the number of sites. These cpu times refer to a PC with a pentium III, 800 megahertz processor.

In conclusion in this paper we have presented and implemented a combinatorial optimization algorithm which solves the optimal cooperation problem in polynomial time. The procedure, which can be used to obtain the partition function of the $`q`$-state random bond Potts model in the large $`q`$-limit is, up to our knowledge, the first application of submodular function optimization in (statistical) physics. Detailed results about the $`q`$-state random bond Potts model, such as properties of phase transitions in two- and three-dimensional lattices will be presented in a separate publication.

Acknowledgment: The authors are highly indebted to Maurice Queyranne for answering a query by e-mail with pointers to , . Realizing that these provide only a separation algorithm we made a second query about sharper algorithms that optimize on (GRAPHIC). We are thankful to Tom McCormick for leading us to Baïou, Barahona and Mahjoub’s paper , which also opened the doors to other relevant work: Cunningham’s ‘Attack and Reinforcement’ and Barahona’s algorithm. Finally we thank further electronic discussions with Francisco Barahona, further useful comments of Maurice Queyranne concerning the manuscript, and András Frank for his courses on submodular functions.

F.I.’s work has been supported by the Hungarian National Research Fund under grant No OTKA TO34183, TO37323, MO28418 and M36803, by the Ministry of Education under grant No. FKFP 87/2001 and by the EC through Centre of Excellence (No. ICA1-CT-2000-70029).


text:
# Report on GR16, Session A3: Mathematical Studies of the Field Equations

## 1 The Horowitz-Myers Soliton

According to a conjecture by Maldacena, supported by a growing body of evidence, there is a correspondence between string theory in Anti–deSitter (AdS) spacetime and a conformal field theory (CFT) on the boundary of AdS. Stability considerations motivated by the AdS/CFT conjecture led to the conjecture that the boundary of asymptotically AdS space–times should be connected. The Riemannian version of this was proved by Witten and Yau and generalized by Cai and Galloway. The non–supersymmetric version of the AdS/CFT correspondence, leads to consideration of asymptotically locally AdS spacetimes with conformal boundary such that one spatial direction is compactified on a circle. If we allow nontrivial topology at infinity, the positive mass theorem for general relativity is not valid in general. In case of negative cosmological constant, $`\mathrm{\Lambda }<0`$, a partial Cauchy surface in an asymptotically locally AdS $`n+1`$-dimensional space–time may have any orientable $`n1`$–manifold as boundary at infinity. The Kottler and Nariai metrics give explicit examples of static AdS metrics on manifolds with nontrivial topology at infinity. These examples all have horizons.

Horowitz and Myers have found a static solution of the Einstein equations with negative cosmological constant, $`\mathrm{\Lambda }<0`$, which is asymptotically locally AdS and contains no horizon. I will refer to this spacetime as the HM soliton. The HM soliton in $`n+1`$ dimensions is given by the metric

$$ds^2=r^2dt^2+\frac{1}{V(r)}dr^2+V(r)d\varphi ^2+r^2\underset{i=1}{\overset{n2}{}}(dy^i)^2$$
(1)
where $`V(r)=\frac{r^2}{\mathrm{}^2}\left(1\frac{r_0^n}{r^n}\right)`$, $`\mathrm{}^2=\frac{n(n1)}{2\mathrm{\Lambda }}`$, and $`r_0`$ is a constant. In order to avoid an orbifold singularity, $`\varphi `$ must be identified with period $`\beta _0=\frac{4\pi \mathrm{}^2}{nr_0}`$. After making the $`y^i`$ coordinates periodic as well (with arbitrary period), the resulting spacetime has conformally compact Cauchy surfaces with flat torus conformal boundary.

Computing the mass of the soliton using a standard method gives the negative answer $`E=\frac{\beta r_0^2}{4G_5\mathrm{}}`$. The mass at infinity of metrics asymptotic to the soliton is related via the AdS/CFT correspondence to the negative Casimir energy of a non–supersymmetric gauge theory on the boundary. If a non–supersymmetric version of the AdS–CFT conjecture holds, then stability considerations imply that the soliton must be the lowest energy state with the given asymptotic conditions. This leads to a new positive mass conjecture,

###### Conjecture 1 (Horowitz and Myers)

Spacetimes satisfying a dominant energy condition with $`\mathrm{\Lambda }<0`$, and with the asymptotic behavior of the soliton (1), have mass bounded from below by $`E`$. The soliton is the unique lowest mass solution for all spacetimes in its class.

This conjecture is supported by perturbation calculations up to second order. A further piece of supporting evidence was provided by the recent work of Galloway, Surya and Woolgar, described by Greg Galloway in his talk at session A.3. Galloway et al. proved that the HM soliton is the unique static metric in a class of metrics with the given asymptotic conditions. In order to formulate the main result of Galloway et. al., let $`(\mathrm{\Sigma },h,N)`$ be a solution of the static vacuum field equations $`R_{ab}=N^1_a_bN+\frac{2\mathrm{\Lambda }}{n1}h_{ab},\mathrm{\Delta }N=\frac{2\mathrm{\Lambda }}{n1}N,`$ which is conformal to a compact manifold with boundary $`(\stackrel{~}{\mathrm{\Sigma }},\stackrel{~}{h})`$ with conformal factor $`\stackrel{~}{N}=1/N`$ (with $`\stackrel{~}{N}=0`$ on $`\stackrel{~}{\mathrm{\Sigma }}`$. Then $`(\mathrm{\Sigma },h,N)`$ is said to satisfy condition (S) if the level sets $`\stackrel{~}{N}=c`$ near infinity are weakly convex in the sense that the principal curvatures are either all positive or all negative. The mass is computed as an integral over the conformal boundary of the mass aspect of $`(\mathrm{\Sigma },h,N)`$, defined as (specializing the Ashtekar–Magnon conformal mass to the static case) $`\mu =^{n2}\stackrel{~}{R}/x^{n2}|_{x=0}`$ (where $`x=\stackrel{~}{N}`$ and $`\stackrel{~}{R}`$ is the scalar curvature function of $`(\stackrel{~}{\mathrm{\Sigma }},\stackrel{~}{h}))`$.

###### Theorem 1.1

Let $`(\mathrm{\Sigma },h,N)`$ be a static spacetime as above, satisfying condition (S). Suppose in addition, that (a) The conformal boundary of $`(\mathrm{\Sigma },h,N)`$ is the same as that of the $`H\&M`$ soliton, (b) The mass aspect $`\mu `$ of $`(\mathrm{\Sigma },h,N)`$ is pointwise negative, and (c) The kernel of the inclusion map $`i_{}:\mathrm{\Pi }_1(\stackrel{~}{\mathrm{\Sigma }})\mathrm{\Pi }_1(\stackrel{~}{\mathrm{\Sigma }})`$ is generated by the $`S^1`$ factor. Then the spacetime $`M^{n+1}=\times \mathrm{\Sigma },g=N^2dt^2h.`$ determined by $`(\mathrm{\Sigma },h,N)`$ is isometric to the AdS soliton (1).

The proof uses complete null lines constructed using the convexity of the boundary together with a null splitting theorem of Galloway to prove that the universal cover of the conformal compactification $`(\stackrel{~}{\mathrm{\Sigma }},\stackrel{~}{h})`$ splits as a metric product $`^{n2}\times W`$ with the metric on $`W`$ determined uniquely by the field equations.

## 2 Dynamical systems in general relativity

Claes Uggla reviewed some applications of dynamical systems methods in general relativity. Dynamical systems formulations arise naturally in the study of Einstein’s field equations for special classes of spacetimes, such as hypersurface–homogenous and hypersurface-self-similar spacetimes. Further examples are provided by the geodesic equations of certain classes of spacetimes, where dynamical systems techniques can be used to e.g., calculate the null geodesics and the corresponding cosmic microwave background radiation pattern in a cosmological model, and the density perturbations of spatially homogenous models. The classical example of application of dynamical systems methods to the Einstein equations is the study of spatially homogenous (Bianchi) space–times, with or without matter. The Bianchi family of spacetimes contains the Friedman spacetimes as a special case, and gives a rich set of examples of the influence of matter and geometry on the early and late time evolution of space–times. It is a natural problem to analyze completely the behavior of the resulting family of dynamical systems.

The book on dynamical systems in cosmology edited by Ellis and Wainwright  gives a good overview of the state of the art of this approach, some years ago. As an outgrowth of that work, substantial progress has been made towards understanding the Bianchi systems both near the initial singularity and in the asymptotically expanding direction. Hans Ringström proved that Bianchi VIII and IX have oscillating singularity, which is an important step in understanding the full BKL picture for these models. An important step in the analysis was to prove that Bianchi II is an attractor for Bianchi IX. This sort of hierarchical behavior in the Bianchi models, is an important organizing principle also for further work on systems with lower dimensional isometry groups. Ringström has also, as he described in his talk, been able to analyze in detail the asymptotic behavior of the expanding vacuum Bianchi class A spacetimes (this excludes Bianchi IX). For the most complicated of these, Bianchi VIII, the result is that the geometry is asymptotically Taub, with precise control on the asymptotic behavior of the commutator variables. For Bianchi VII0, the system approaches (in terms of scale free variables), a flat VII0 spacetime. The Cauchy surfaces of these spacetimes are all Seifert fibered. The results show that the scale free geometry of the spatial sections in these spacetimes collapse along the fiber direction. The result is consistent with the picture for the limiting behavior in the asymptotic expanding direction developed by M. T. Anderson Analyzing the hierarchical structure of systems with small symmetry group, such as Gowdy models, or without symmetry, is likely to be a crucial step in understanding the global dynamics of the Einstein equations.

## 3 Well-posed forms of the 3+1 trace–decoupled Einstein equations

The standard ADM form of the Einstein evolution equations fails to be well posed both from the mathematical and numerical points of view. Due to the demand for stable numerical procedures for solving the Einstein evolution equations, a great deal of interest has recently been focussed on well posed hyperbolic formulations of the Einstein equations. This is of particular importance in solving boundary value problems. Oscar Reula discussed various well–posed reformulations of the Einstein equations. It is possible to derive a well–posed form of the standard 3+1 Einstein equations by reducing to first order form, densitizing the lapse function and adding combinations of the constraints to the evolution equations. An alternative approach, which recently has been shown to possess striking computational advantages over the standard form, is to decompose the fundamental variables in order to extract the trace of the extrinsic curvature and the determinant of the 3-metric with the purpose of evolving them independently from the remaining components. This trace–decoupled form of the Einstein equation turns out, oddly enough, to be ill–posed. Reula and Frittelli have used techniques similar to those used for the standard evolution equations, to obtain versions of the 3+1 Einstein equations which are both trace-decoupled and well posed. This well posed version propagates the constraints in a stable manner, which is relevant to unconstrained evolution. It uses essentially the same variables as the original system, but requires a gauge condition, namely that the lapse to be proportional to the determinant of the intrinsic geometry of the surfaces.

## 4 Critical phenomena

Critical phenomena in general relativity have been extensively studied since the discovery by Choptuik of universal behavior in the collapse of a self–gravitating scalar field. Since then a number of matter models, such as charged self–gravitating scalar fields, Yang–Mills and dilaton fields have been studied and a systematic understanding of these phenomena has been gained, in particular through the study of the stability of models with continuous or discrete self–similarity. We now have a good qualitative understanding of criticality in gravitational collapse: there is an intermediate (codimension-1) attractor of the evolution flow in Phase Space, sitting on the threshold of black hole formation. The symmetries of that attractive solution determine the type of criticality: ‘type I’ if the solution is static, ‘type II’ when it is self-similar.

The Einstein–Vlasov system describes the evolution of a statistical ensemble of non-interacting particles coupled to gravity through their average properties. Jose M. Martin–Garcia described joint work with Carsten Gundlach on critical phenomena in the Einstein–Vlasov system, which is important in Critical Phenomena Theory because it is the only one where type II critical phenomena have not yet been found.

There are two different ways of checking for the existence of criticality. On the one hand, it is possible to evolve numerically families of initial conditions and study the collapse threshold, looking for intermediate attractors and universality. On the other hand, one can construct static or self-similar solutions and then investigate whether they are codimension-1 linearly stable.

There are two relevant works in this problem, both using the method of evolution of initial conditions: Rendall, Rein and Schaeffer have not found criticality and suggest that, if there is any, it must be type I. Olabarrieta and Choptuik confirm the apparent nonexistence of type II phenomena and have found some signs of type I criticality, though not conclusive because there is not universality in the critical exponents. Gundlach and Martin–Garcia look for type II criticality in the massless Einstein-Vlasov system using the second method.

In his talk Martin–Garcia reported on the first step: the numerical construction of self-similar solutions, which is an interesting problem in itself. The situation resembles what is already known analytically for static solutions of the same system. There is a free function of two variables, representing different distributions of the energy-momentum among the particles of the system, but all giving rise to the same gravitational field. Therefore a candidate critical solution can not be isolated. It is likely that a similar result in the massive case could explain the observed non-universality.

Recent numerical investigations of near-critical scalar field collapse in 2+1 dimensional AdS spacetime show a continuous self-similar (CSS) behaviour near the central singularity and power law scaling for the black hole mass. Garfinkle found that this behaviour is well approximated, at intermediate times, by a member of a one-parameter family of exact CSS solutions to the equations with vanishing cosmological constant $`\mathrm{\Lambda }`$. However the extension of these solutions to $`\mathrm{\Lambda }0`$, necessary to explain black hole formation (as vacuum black holes exist only for $`\mathrm{\Lambda }<0`$) is a hard problem which has not been solved up to now.

Gérard Clément in his talk discussed joint work with Alessandro Fabbri on quasi-CSS 2+1 dimensional scalar field spacetimes. They have derived by a limiting process a new class of $`\mathrm{\Lambda }=0`$ CSS solutions

$`ds^2`$ $`=`$ $`dudv(u)^{2/(1+c^2)}d\theta ^2,`$
$`\varphi `$ $`=`$ $`{\displaystyle \frac{c}{1+c^2}}\mathrm{ln}(u)`$ (2)
These present for $`c^21`$ a point singularity ($`u=0`$, $`v=+\mathrm{}`$), and for $`c^2<1`$ a null line singularity $`u=0`$. They can be extended to exact solutions of the $`\mathrm{\Lambda }=l^2<0`$ equations by means of the ansatz

$$ds^2=e^{2\sigma (x)}dudv(u)^{\frac{2}{c^2+1}}\rho ^2(x)d\theta ^2,\varphi =\frac{c}{c^2+1}\mathrm{ln}|u|+\psi (x),$$
(3)
where $`x=uv`$, with the initial conditions $`\rho (0)=1,\sigma (0)=\psi (0)=0`$. These extended solutions inherit the CSS behaviour close to the central singularity, and have the correct AdS behaviour at spatial infinity.

To show that these quasi-CSS solutions are indeed threshold solutions for black hole formation, one studies their linear perturbations. The case $`c=0`$ is at the same time simple and interesting. In this case the scalar field decouples, $`\varphi =0`$, and the corresponding quasi-CSS solution is the vacuum BTZ solution. The linear perturbation equations with appropriate boundary conditions lead to a unique solution, with eigenvalue $`k=2`$. This solution is found to be an exact linearization in $`M`$ of the BTZ black hole. In the case of genuine scalar perturbations, one finds two possible modes $`k_a=c^2+3/2`$ and $`k_b=c^2+2`$. For the $`a`$ mode the singularity and the apparent horizon appear simultaneously on the null line $`u=0`$ at the time $`v=0`$ and evolve in the region $`v>0`$ in a physically meaningful way, while for the $`b`$ mode the singularity still appears for $`v=0`$, but the apparent horizon seems to be eternal, as in the case of the static BTZ black hole. The $`b`$ mode, which cannot describe actual gravitational collapse with regular initial conditions, appears to be unphysical. The critical exponent is related to the eigenvalue $`k`$ by $`\gamma =1/k`$. In the BTZ case $`c=0`$ Clément and Fabbri obtain $`\gamma =1/2`$. In the scalar field case, choosing the critical value $`c^2=1`$ they obtain for the $`a`$ mode $`\gamma =0.4`$. This value differs signicantly from those found in the numerical analysis, i.e. $`\gamma 1.2`$ and $`\gamma 0.81`$, showing that further work is needed in order to clarify this issue.

## 5 Isolated systems

It is an open problem to construct a nonflat vacuum spacetime, which has a regular conformal infinity in the sense of Penrose. The conformal compactification of a nonflat spacetime must have a singularity at spatial infinity $`\text{i}_0`$, due to the slow falloff of the Weyl tensor at $`\text{i}_0`$. Sergio Dain described joint work with Helmut Friedrich, which constructs a large class of asymptotically flat initial data with non-vanishing mass and angular momentum for which the metric and the extrinsic curvature have asymptotic expansions at space-like infinity in terms of powers of a radial coordinate. These asymptotic expansions are of the form $`\stackrel{~}{h}_{ij}(1+\frac{2m}{\stackrel{~}{r}})\delta _{ij}+_{k2}\frac{\stackrel{~}{h}_{ij}^k}{\stackrel{~}{r}^k},\stackrel{~}{\mathrm{\Psi }}_{ij}_{k2}\frac{\stackrel{~}{\mathrm{\Psi }}_{ij}^k}{\stackrel{~}{r}^k},`$ where $`\stackrel{~}{h}_{ij}^k`$ and $`\stackrel{~}{\mathrm{\Psi }}_{ij}^k`$ are smooth function on the unit 2-sphere (thought as being pulled back to the spheres $`\stackrel{~}{r}=const.`$ under the map $`\stackrel{~}{x}^j\stackrel{~}{x}^j/\stackrel{~}{r}`$). In earlier work, stationary data were shown to admit expansions of this type.

In the study of isolated system, an important and difficult problem is to give an unambiguous definition of angular momentum. The difficulty lies in the fact that without introducing extra structure, there is no unambiguously defined Poincaré subgroup of the asymptotic isometries of an isolated system. This is known as the supertranslation ambiguity. This ambiguity is known as supertranslation ambiguity. In order to provide with an unambiguous notion of angular momentum in radiative spacetimes, it is essential to make use of a reference frame system, which embodies the notion of rest frame. Osvaldo Moreschi presented a definition of angular momentum for radiative spacetimes which does not suffer from supertranslation ambiguity. The definition uses the construction of so called ‘nice sections’ at $`^+`$. The defining equations for these have recently been proved to have solutions in terms of a 4-parameter family of translations. The solutions have the expected physical properties of a rest frame. The nice section construction singles out a Poincaré structure from the infinite dimensional BMS group. In particular, given a fixed observational point $`p`$ at $`^+`$, there is precisely a 3-degree of freedom of spacelike translations which generate all the nice sections that contain $`p`$. In contrast, without this constructions there is an infinite dimensional family of general sections that contain $`p`$, one for each supertranslation. The condition that the spatial part of the Bondi 4–momentum is zero singles out a one parameter family of nonintersecting sections $`S_{\mathrm{𝚌𝚖}}`$ of $`^+`$, which define a center of mass frame and can be used to describe the detailed asymptotic structure of the spacetime. The angular momentum is defined in terms of a charge integral over the $`S_{\mathrm{𝚌𝚖}}`$.

Detectors of gravitational waves can be considered as observers at future null infinity. The construction of center of mass sections, provides for each point $`p`$ at scri a prescription that singles out the unique center of mass section containing $`p`$, together with the appropriate rest frame to calculate intrinsic quantities like angular momentum and multipole moments, which are important for the description of the gravitational waves that one wishes to detect.

John L. Friedman described described in his talk joint work with Kōji Uryū and Masaru Shibata. They consider compact binary systems, modeled as vacuum or perfect-fluid spacetimes with a helical Killing vector $`k^a`$, heuristically, the generator of time-translations in a corotating frame. Systems that are stationary in this sense are not asymptotically flat, but have asymptotic behavior corresponding to equal constant fluxes of ingoing and outgoing radiation. For black-hole binaries, a rigidity theorem implies that the Killing vector lies along the horizon’s generators, and from this one can deduce the zeroth law (constant surface gravity of the horizon). Remarkably, although the mass and angular momentum of such a system are not defined, a finite Noether charge is defined on any sphere enclosing the matter and black holes, its value the same on each such sphere. An exact first law relates the change in this charge to the changes in the vorticity, baryon mass and entropy of the fluid and in the area of black holes. Modeling a binary system by a spacetime with a helical Killing vector is accurate if the energy emitted in a dynamical time is small compared to the kinetic energy of the system.

Other approaches use a post–Newtonian approximation, or conformally flat spacelike slices that satisfy a truncated set of field equations. Like the post-Newtonian spacetimes, these conformally flat spacetimes are nonradiative and asymptotically flat. In the isolated horizon framework, for a horizon with a single Killing vector, one shows the existence of a charge $`E`$ defined on an isolated horizon for which $`\delta E=\kappa \delta A`$. The first law for black holes with helical symmetry, in contrast, relates this change in the black-hole charges to the changes in the Noether charge of a sphere surrounding all black holes and all matter and to the changes in the entropy, baryon number, and circulation of the fluid. The existence of such a first law depends precisely on what is not assumed in the isolated horizon framework: a globally defined Killing vector.

### 5.1 Isolated horizons

A closely related topic, isolated horizons (IH) was discussed by G. Date. Here again the work is motivated by the question if it is possible to have analogues of the usual laws of BH mechanics in spacetimes more general than the usual stationary asymptotically flat black holes.

G. Date stated the following results: A Killing horizon (KH) is always an IH, and an IH is a KH provided the IH admits a neighbourhood which is either stationary or has two commuting Killing vectors. Further, a vacuum, non-rotating IH necessarily admits an EH and for a rotating IH , an EH exists if $`+\pi \overline{\pi }2\pi <0`$. Otherwise, the IH is accessible from infinity. An EH could exist but the IH will be ‘outside’ of it.

### 5.2 Stars

The Einstein equations for stationary axisymmetric perfect fluid spacetimes reduces to an elliptic system, known as the Ernst equations. These spacetimes are commonly studied as models of isolated stars in general relativity.

Christian Klein presented joint work with Jörg Frauendiener on exact solutions of the stationary axisymmetric Einstein equations for perfect fluid matter. A boundary value problem for the Ernst equations is solved using Riemann–Hilbert techniques. As an application he presented the solution for a family of counter-rotating dust disks which contains the rigidly rotating dust disk as a limiting case.

H. Pfister considered existence and non–existence results proved using Banach space fixed point techniques.

## 6 Variational formulation

Naresh Dadhich and Niall Ó Murchadha gave talks relating to the variational formulation of general relativity. Dadhich considered the conditions under which one can deduce the field equations from the equations of motion for particles. In the context of electomagnetic field, this question was, according to Dyson , first considered by Feynman. By considering commutation relations between position and velocity, he could obtain the Bianchi set of the equations. This gave rise to considerable activity in this direction involving rederivations, generalizations and extensions. However the attention centred on the Bianchi set. Singh and Dadhich derived the full set of the Maxwell equations by demanding the law of motion to be linear in velocity and derivable from a Lagrangian. A similar procedure has been applied to the Einstein equations. In this case, the fact that the particle law in general relativity is quadratic in the 4–velocity leads to the Einstein equations. This is a novel way of deducing the field equations from the particle law of motion.

Niall Ó Murchadha described joint work with Julian Barbour and Brendan Foster, where a derivation of General Relativity was given which does not depend on the relativity principle. The fundemental idea is to choose superspace (the space of all riemannian three geometries) as the configuration space and to consider a parametrised curve in this configuration space, with parameter $`\lambda `$. Ó Murchadha et al. consider a generalized Jacobi-type action of the form $`𝑑\lambda 𝑑x^3\sqrt{gPT}`$ on such curves where one first integrates on each given three manifold and then integrates along the curve. ‘$`P`$’ is the potential which is a scalar on the three manifold but otherwise arbitrary, while ‘$`T`$’ is the generalized kinetic energy which means that it is a quadratic function of $`g_{ij}/\lambda `$. It turns out that in order for the constraints of the theory to be propagated, it is necessary to choose the potential energy term to be the three scalar curvature and the kinetic energy to be defined with the DeWitt supermetric. This recovers the Baierlein, Sharp, Wheeler parametrised action for G.R. as the unique self-consistent Jacobi-type action on superspace.

## 7 Singularities

There were several talks concerned with the structure of cosmological and black hole singularities.

An important development has been the application of Fuchsian techniques to construct families of asymptotically velocity term dominated (AVTD) cosmological singularities. M. Narita discussed his work , applying the Fuchsian method to Gowdy spacetimes with stringy matter, showing that these have AVTD singularities. Stringy gravity has receieved a great deal of attention recently due to the work of Damour and Henneaux, who showed that for stringy gravity in of no symmetries in spacetime dimension 11 and higher, one expects mixmaster type behavior.

Brien Nolan has studied the gravitational collapse of spherically symmetric thick shells admitting a homothetic Killing vector field under the assumption that the energy momentum tensor corresponds to the absence of a pure outgoing component. The stability of the Cauchy horizon against linear perturbations is investigated. The behaviour of a massless scalar field $`\varphi `$ propagating on a background space-time is studied. The following boundary conditions are imposed: (i) finiteness of the flux of $`\varphi `$ as measured by any time-like observer crossing $`N`$ and (ii) finiteness of the flux of $`\varphi `$ in the ingoing null direction measured at $`𝒥^{}`$. Imposing these boundary conditions, it is found that the flux of $`\varphi `$ as measured by any time-like (geodesic) observer crossing the Cauchy horizon is always finite. This is in contrast to the case of Cauchy horizons inside black holes. Those arising here always lie outside any possible event horizon. This indicates that the Cauchy horizon may be stable against linear and possibly non-linear gravitational perturbations, and so members of this class of space-times may provide physically realisable examples of naked singularities.

A. Beesham also discussed naked singularities. In joint work with S.G. Ghosh and R.V. Saraykar he has studied gravitational collapse of radiation shells in a non self-similar higher dimensional spherically symmetric spacetime. Strong curvature naked singularities were found to form for a highly inhomogeneous collapse, violating the cosmic censorship conjecture for this class of models.

Deborah Konkowski reviewed the status of two instability conjectures for Cauchy horizons.

## 8 Miscellaneous

### 8.1 Lanczos and Bel–Robinson

In his talk, Brian Edgar reviewed his work with F. Andersson on the Lanczos potential for the Weyl tensor. The existence proof for the Lanczos potential is now on firm foundation, and the construction of the potential in several situations with symmetry is understood. Unfortunately, a really interesting application of this intriguing object is still lacking.

Jose M. Senovilla discussed his work on Bel–Robinson squares. This work generalizes and unifies the earlier literature on Bel–Robinson and Bel tensors. The Bel–Robinsons tensor was introduced as a candidate for a “stress–energy” tensor for the gravitational field, and in 3+1 dimensional vacuum spacetimes, this symmetric, trace–free 4–tensor is divergence free and positive definite in a certain sense. The Bel tensor plays a similar role for Einstein equations with matter. The work of Senovilla gives a systematic procedure for constructing similar objects in higher dimensional spacetimes.

### 8.2 Distributional Geometry

James Vickers and Michael Kunzinger both in their talks discussed distributional geometry and its relation to the cosmic censorship hypothesis. Kunzinger concetrated on applications of ideas from Colombeau theory in distributional geometry, while Vickers concentrated on some examples related to cosmic censorship.

The study of singular spacetimes by distributional methods faces the fundamental obstacle of the inherent nonlinearity of the field equations. Staying strictly within the distributional (in particular: linear) regime, excludes a number of physically interesting examples such as strings and point particles.

In recent years, several authors have therefore employed nonlinear theories of generalized functions, in particular Colombeau’s theory of generalized functions in order to derive a suitable mathematical framework for a general “nonlinear distributional geometry” adapted to the needs of general relativity.

Under the influence of these applications in general relativity the nonlinear theory of generalized functions itself has undergone a rapid development lately, resulting in a diffeomorphism invariant global theory of nonlinear generalized functions on manifolds. In particular, a generalized pseudo-Riemannian geometry allowing for a rigorous treatment of generalized (distributional) spacetime metrics has been developed.

According to the Cosmic Censorship hypothesis realistic singularities should be hidden by an event horizon. However there are many examples of physically realistic spacetimes which are geodesically incomplete, and hence singular according to the usual definition, which are not inside an event horizon.

Many of these counterexamples to the cosmic censorship conjecture have a curvature tensor which is reasonably behaved (for example bounded or integrable) as one approaches the singularity. James Vickers gave examples of a class of weak singularities which may be described as having distributional curvature.

The propagation of test fields on spacetimes with weak singularities has also been investigated. A class of singularities which do not disrupt the Cauchy development of test fields and result in spacetimes which satisfy Clarke’s criterion of ‘generalised hyperbolicity’ has been found. Vickers argued that points which are well behaved in this way, and where Einstein’s equations make sense distributionally, should be regarded as interior points of the spacetime rather than counterexamples to cosmic censorship.

### 8.3 The a–boundary

Susan M. Scott discussed the application of the a–boundary construction to cosmic censorship.

The singularity theorems of Hawking and Penrose prove the existence of incomplete causal geodesics in space-times which satisfy quite general physical conditions. It has long been conjectured that these incomplete causal geodesics would terminate at curvature singularities, but this has never been proven. recently established a theorem which goes a considerable distance towards closing this gap left by the singularity theorems. By application of the abstract boundary (a-boundary) construction of Scott and Szekeres, to strongly causal space-times, Ashley and Scott have shown that strongly causal space-times satisfying the conditions of the singularity theorems will always have an a-boundary essential singularity.

### 8.4 The Ehlers group

Marc Mars discussed his work on the Ehlers group. The Ehlers group, which is a symmetry of the Einstein vacuum field equations for strictly stationary spacetimes, is considered by Mars from a spacetime perspective. In this setting, the Ehlers group becomes a subgroup of the infinite dimensional group of transformations that maps Lorentz metrics into Lorentz metrics. Mars has analyzed the global conditions which are required on the spacetime for the existence of the Ehlers group and has found the transformation law for the Weyl tensor under Ehlers transformations. This allows one to study where, and under which circumstances, curvature singularities in the transformed spacetime will arise. As an application, one finds a local characterization of the Kerr-NUT metric.

### 8.5 Composite spacetimes

Katsuhito Yasuno considered the dynamics of spatially compact ‘composite’ spacetimes. Motivated by Thurston’s geometrization conjecture, it is interesting to look for examples of spacetimes with more complicated topology than is allowed by the spatially homogenous models. The composite spacetimes are built from the spatially compact locally homogeneous vacuum spacetimes which have two commuting local Killing vector fields and are homeomorphic to torus bundles over the circle.


text:
# Radiative M1-decays of heavy-light mesons in the relativistic quark model

## Abstract

Radiative magnetic dipole decays of heavy-light vector mesons into pseudoscalar mesons $`VP\gamma `$ are considered within the relativistic quark model. The light quark is treated completely relativistically, while for the heavy quark the $`1/m_Q`$ expansion is used. It is found that relativistic effects result in a significant reduction of decay rates. Comparison with previous predictions and recent experimental data is presented.

HU-EP-02/14

, ,

In this letter we consider radiative magnetic dipole (M1) transitions of the ground state vector ($`V`$) heavy-light mesons to the pseudoscalar ($`P`$) ones, $`VP\gamma `$ (in quark model notations, $`1^3S_11^1S_0+\gamma `$). For this purpose, we use the relativistic quark model based on the quasipotential approach in quantum field theory. Recently, this model has been successfully applied for the description of different properties of the heavy-light mesons, such as their mass spectra  and rare radiative decays . Our analysis showed that the light quark in the heavy-light mesons should be treated completely relativistically, while for the heavy quark it is useful to apply the expansion in powers of the inverse heavy quark mass $`1/m_Q`$, which considerably simplifies calculations. The first and sometimes second order corrections in $`1/m_Q`$ are also important for the heavy-light meson decay description. Analogously, for the radiative decay calculations considered here, the $`1/m_Q`$ expansion is carried out up to the second order, and the light quark is treated completely relativistically (i. e. without the unjustified expansion in inverse powers of the light quark mass). It follows from the obtained results that the relativistic effects give substantial contributions to the calculated decay rates.

The radiative $`VP\gamma `$ decay rate is given by

$$\mathrm{\Gamma }=\frac{\omega ^3}{3\pi }\left|𝓜\right|^2,\mathrm{where}\omega =\frac{M_V^2M_P^2}{2M_V},$$
(1)
$`M_V`$ and $`M_P`$ are the vector and pseudoscalar meson masses. The matrix element of the magnetic moment $`𝓜`$ is defined by

$$𝓜=\frac{i}{2}\left[\frac{}{𝚫}\times P\left|𝐉(0)\right|V\right]_{𝚫=\mathrm{𝟎}},𝚫=𝐏𝐐,$$
(2)
where $`P\left|J_\mu (0)\right|V`$ is the matrix element of the electromagnetic current between initial vector ($`V`$) and final pseudoscalar ($`P`$) meson states with momenta $`𝐐`$ and $`𝐏`$ respectively.

We use the relativistic quark model for the calculation of the matrix element of the magnetic moment $`𝓜`$ (2). In our model a meson is described by the wave function of the bound quark-antiquark state, which satisfies the quasipotential equation of the Schrödinger type in the center-of-mass frame :

$$\left(\frac{b^2(M)}{2\mu _R}\frac{𝐩^2}{2\mu _R}\right)\mathrm{\Psi }_M(𝐩)=\frac{d^3q}{(2\pi )^3}V(𝐩,𝐪;M)\mathrm{\Psi }_M(𝐪),$$
(3)
where the relativistic reduced mass is

$$\mu _R=\frac{M^4(m_1^2m_2^2)^2}{4M^3},$$
(4)
and $`b^2(M)`$ denotes the on-mass-shell relative momentum squared

$$b^2(M)=\frac{[M^2(m_1+m_2)^2][M^2(m_1m_2)^2]}{4M^2}.$$
(5)
Here $`m_{1,2}`$ and $`M`$ are quark masses and a heavy-light meson mass, respectively.

The kernel $`V(𝐩,𝐪;M)`$ in Eq. (3) is the quasipotential operator of the quark-antiquark interaction. It is constructed with the help of the off-mass-shell scattering amplitude, projected onto the positive energy states. An important role in this construction is played by the Lorentz-structure of the confining quark-antiquark interaction in the meson. In constructing the quasipotential of the quark-antiquark interaction we have assumed that the effective interaction is the sum of the usual one-gluon exchange term and the mixture of vector and scalar linear confining potentials. The quasipotential is then defined by

$$V(𝐩,𝐪;M)=\overline{u}_1(p)\overline{u}_2(p)𝒱(𝐩,𝐪;M)u_1(q)u_2(q),$$
(6)
with

$$𝒱(𝐩,𝐪;M)=\frac{4}{3}\alpha _sD_{\mu \nu }(𝐤)\gamma _1^\mu \gamma _2^\nu +V_{\mathrm{conf}}^V(𝐤)\mathrm{\Gamma }_1^\mu \mathrm{\Gamma }_{2;\mu }+V_{\mathrm{conf}}^S(𝐤),$$

where $`\alpha _s`$ is the QCD coupling constant, $`D_{\mu \nu }`$ is the gluon propagator in the Coulomb gauge and $`𝐤=𝐩𝐪`$; $`\gamma _\mu `$ and $`u(p)`$ are the Dirac matrices and spinors. The effective long-range vector vertex is given by

$$\mathrm{\Gamma }_\mu (𝐤)=\gamma _\mu +\frac{i\kappa }{2m}\sigma _{\mu \nu }k^\nu ,k^\nu =(0,𝐤),$$
(7)
where $`\kappa `$ is the Pauli interaction constant characterizing the nonperturbative anomalous chromomagnetic moment of quarks. Vector and scalar confining potentials in the nonrelativistic limit reduce to

$$V_{\mathrm{conf}}^V(r)=(1\epsilon )(Ar+B),V_{\mathrm{conf}}^S(r)=\epsilon (Ar+B),$$
(8)
reproducing

$$V_{\mathrm{conf}}(r)=V_{\mathrm{conf}}^S(r)+V_{\mathrm{conf}}^V(r)=Ar+B,$$
(9)
where $`\epsilon `$ is the mixing coefficient.

The quasipotential for the heavy quarkonia, expanded in $`p^2/m^2`$, can be found in Ref.  and for heavy-light mesons in . All the parameters of our model, such as quark masses, parameters of the linear confining potential, mixing coefficient $`\epsilon `$ and anomalous chromomagnetic quark moment $`\kappa `$, were fixed from the analysis of heavy quarkonium spectra  and radiative decays . The quark masses $`m_b=4.88`$ GeV, $`m_c=1.55`$ GeV, $`m_s=0.50`$ GeV, $`m_{u,d}=0.33`$ GeV and the parameters of the linear potential $`A=0.18`$ GeV2 and $`B=0.30`$ GeV have the usual quark model values. In Ref.  we have considered the expansion of the matrix elements of weak heavy quark currents between pseudoscalar and vector meson ground states up to the second order in inverse powers of the heavy quark masses. It has been found that the general structure of the leading, first, and second order $`1/m_Q`$ corrections in our relativistic model is in accord with the predictions of HQET. The heavy quark symmetry and QCD impose rigid constraints on the parameters of the long-range potential in our model. The analysis of the first order corrections  fixes the value of the Pauli interaction constant $`\kappa =1`$. The same value of $`\kappa `$ was found previously from the fine splitting of heavy quarkonia $`{}_{}{}^{3}P_{J}^{}`$\- states . The value of the parameter characterizing the mixing of vector and scalar confining potentials, $`\epsilon =1`$, was found from the comparison of the second order ($`1/m_Q^2`$) corrections in our model  with the same order contributions in HQET. This value is very close to the one determined from considering radiative decays of heavy quarkonia , especially the M1-decays (e. g. the calculated decay rate of $`J/\mathrm{\Psi }\eta _c\gamma `$ can be brought in accord with the experiment only with the above value of $`\epsilon `$).

In the quasipotential approach, the matrix element of the electromagnetic current $`J_\mu `$ between the states of a vector $`V`$ meson and a pseudoscalar $`P`$ meson has the form

$$P|J_\mu (0)|V=\frac{d^3pd^3q}{(2\pi )^6}\overline{\mathrm{\Psi }}_{P𝐏}(𝐩)\mathrm{\Gamma }_\mu (𝐩,𝐪)\mathrm{\Psi }_{V𝐐}(𝐪),$$
(10)
where $`\mathrm{\Gamma }_\mu (𝐩,𝐪)`$ is the two-particle vertex function and $`\mathrm{\Psi }_{V,P}`$ are the meson wave functions projected onto the positive energy states of quarks and boosted to the moving reference frame. The contributions to $`\mathrm{\Gamma }`$ come from Figs. 1 and 2. The contribution $`\mathrm{\Gamma }^{(2)}`$ is the consequence of the projection onto the positive-energy states. Note that the form of the relativistic corrections resulting from the vertex function $`\mathrm{\Gamma }^{(2)}`$ explicitly depends on the Lorentz structure of the $`q\overline{q}`$-interaction. Thus the vertex function is given by

$$\mathrm{\Gamma }_\mu (𝐩,𝐪)=\mathrm{\Gamma }_\mu ^{(1)}(𝐩,𝐪)+\mathrm{\Gamma }_\mu ^{(2)}(𝐩,𝐪)+\mathrm{},$$
(11)
where

$$\mathrm{\Gamma }_\mu ^{(1)}(𝐩,𝐪)=e_1\overline{u}_1(p_1)\gamma _\mu u_1(q_1)(2\pi )^3\delta (𝐩_2𝐪_2)+(12),$$
(12)
and

$`\mathrm{\Gamma }_\mu ^{(2)}(𝐩,𝐪)`$ $`=`$ $`e_1\overline{u}_1(p_1)\overline{u}_2(p_2)\{𝒱(𝐩_2𝐪_2){\displaystyle \frac{\mathrm{\Lambda }_1^{()}(k_1^{})}{ϵ_1(k_1^{})+ϵ_1(q_1)}}\gamma _1^0\gamma _{1\mu }`$ (13)
$`+\gamma _{1\mu }{\displaystyle \frac{\mathrm{\Lambda }_1^{()}(k_1)}{ϵ_1(k_1)+ϵ_1(p_1)}}\gamma _1^0𝒱(𝐩_2𝐪_2)\}u_1(q_1)u_2(q_2)+(12).`$
Here $`e_{1,2}`$ are the quark charges, $`𝐤_1=𝐩_1𝚫;𝐤_1^{}=𝐪_1+𝚫;𝚫=𝐏𝐐`$;

$$\mathrm{\Lambda }^{()}(p)=\frac{ϵ(p)\left(m\gamma ^0+\gamma ^0(𝜸𝐩)\right)}{2ϵ(p)},ϵ(p)=\sqrt{p^2+m^2}.$$

It is important to note that the wave functions entering the current matrix element (10) cannot be both in the rest frame. In the initial $`V`$ meson rest frame, the final $`P`$ meson is moving with the recoil momentum $`𝚫`$. The wave function of the moving $`P`$ meson $`\mathrm{\Psi }_{P𝚫}`$ is connected with the wave function in the rest frame $`\mathrm{\Psi }_{P\mathbf{\hspace{0.17em}0}}\mathrm{\Psi }_P`$ by the transformation

$$\mathrm{\Psi }_{P𝚫}(𝐩)=D_1^{1/2}(R_{L_𝚫}^W)D_2^{1/2}(R_{L_𝚫}^W)\mathrm{\Psi }_{P\mathbf{\hspace{0.17em}0}}(𝐩),$$
(14)
where $`R^W`$ is the Wigner rotation, $`L_𝚫`$ is the Lorentz boost from the rest frame to a moving one, and $`D^{1/2}(R)`$ is the rotation matrix in the spinor representation.

We substitute the vertex functions $`\mathrm{\Gamma }^{(1)}`$ and $`\mathrm{\Gamma }^{(2)}`$ given by Eqs. (12) and (13) in the decay matrix element (10) and take into account the wave function transformation (14). To simplify calculations we note that the mass of the heavy-light mesons $`M_{V,P}`$ is large (due to presence of the heavy quark $`M_{V,P}m_Q`$) and carry out the expansion in inverse powers of this mass up to the second order. Then we calculate the matrix element of the magnetic moment operator (2) and get
(a) for the vector potential

$`𝓜_V`$ $`=`$ $`{\displaystyle }{\displaystyle \frac{d^3p}{(2\pi )^3}}\overline{\mathrm{\Psi }}_P(𝐩){\displaystyle \frac{e_1}{2ϵ_1(p)}}\{𝝈_1+{\displaystyle \frac{(1\epsilon )(1+2\kappa )[𝐩\times [𝝈_1\times 𝐩]]}{2ϵ_1(p)[ϵ_1(p)+m_1]}}`$
$`+{\displaystyle \frac{(1\epsilon )(1+\kappa )[𝐩\times [𝝈_2\times 𝐩]]}{ϵ_1(p)[ϵ_2(p)+m_2]}}`$
$`{\displaystyle \frac{ϵ_2(p)}{M_V}}\left(1+(1\epsilon ){\displaystyle \frac{M_Vϵ_1(p)ϵ_2(p)}{ϵ_1(p)}}\right)i\left[𝐩\times {\displaystyle \frac{}{𝐩}}\right]`$
$`+{\displaystyle \frac{1}{2M_V}}[𝐩\times [𝐩\times ({\displaystyle \frac{𝝈_1}{ϵ_1(p)+m_1}}{\displaystyle \frac{𝝈_2}{ϵ_2(p)+m_2}})]]\}\mathrm{\Psi }_V(𝐩)+(12),`$
(b) for the scalar potential

$`𝓜_S`$ $`=`$ $`{\displaystyle }{\displaystyle \frac{d^3p}{(2\pi )^3}}\overline{\mathrm{\Psi }}_P(𝐩){\displaystyle \frac{e_1}{2ϵ_1(p)}}\{(1+\epsilon {\displaystyle \frac{ϵ_1(p)+ϵ_2(p)M_V}{ϵ_1(p)}})`$
$`\times \left(𝝈_1{\displaystyle \frac{ϵ_2(p)}{M_V}}i\left[𝐩\times {\displaystyle \frac{}{𝐩}}\right]\right){\displaystyle \frac{\epsilon [𝐩\times [𝝈_1\times 𝐩]]}{2ϵ_1(p)[ϵ_1(p)+m_1]}}`$
$`+{\displaystyle \frac{1}{2M_V}}[𝐩\times [𝐩\times ({\displaystyle \frac{𝝈_1}{ϵ_1(p)+m_1}}{\displaystyle \frac{𝝈_2}{ϵ_2(p)+m_2}})]]\}\mathrm{\Psi }_V(𝐩)+(12).`$
Note that the last terms in Eqs. (Radiative M1-decays of heavy-light mesons in the relativistic quark model), (Radiative M1-decays of heavy-light mesons in the relativistic quark model) result from the wave function transformation (14) from the moving reference frame to the rest one. It is easy to see that in the limit $`p/m0`$ the usual nonrelativistic expression for the magnetic moment follows.

Since we are interested in radiative transitions of the vector mesons to the pseudoscalar mesons it is possible to evaluate spin matrix elements using the relation $`𝝈_1=𝝈_2`$. Then assuming one quark to be light $`q`$ and the other one $`Q`$ to be heavy and further expanding Eqs. (Radiative M1-decays of heavy-light mesons in the relativistic quark model), (Radiative M1-decays of heavy-light mesons in the relativistic quark model) in the inverse powers of the heavy quark mass $`m_Q`$ up to the second order corrections to the leading contribution we get
(a) for the purely vector potential ($`\epsilon =0`$)

$`_V`$ $`=`$ $`{\displaystyle \frac{e_q}{2m_q}}\{{\displaystyle \frac{m_q}{ϵ_q(p)}}{\displaystyle \frac{m_q}{3}}{\displaystyle \frac{𝐩^2}{ϵ_q(p)[ϵ_q(p)+m_q]}}({\displaystyle \frac{1}{ϵ_q(p)}}+{\displaystyle \frac{1}{M_V}})`$ (33)
$`+{\displaystyle \frac{(1+\kappa )m_q}{3}}{\displaystyle \frac{𝐩^2}{ϵ_q^2(p)}}({\displaystyle \frac{2}{ϵ_q(p)+m_q}}{\displaystyle \frac{1}{m_Q}}){\displaystyle \frac{m_q}{6M_Vm_Q}}{\displaystyle \frac{𝐩^2}{ϵ_q(p)}}\}`$
$`{\displaystyle \frac{e_Q}{2m_Q}}\{1{\displaystyle \frac{2𝐩^2}{3m_Q^2}}+{\displaystyle \frac{1+\kappa }{3}}{\displaystyle \frac{𝐩^2}{m_Q}}({\displaystyle \frac{1}{m_Q}}{\displaystyle \frac{2}{ϵ_q(p)+m_q}})`$
$`{\displaystyle \frac{𝐩^2}{6M_V}}({\displaystyle \frac{1}{m_Q}}+{\displaystyle \frac{2}{ϵ_q(p)+m_q}})\},`$
(b) for the purely scalar potential ($`\epsilon =1`$)

$`_S`$ $`=`$ $`{\displaystyle \frac{e_q}{2m_q}}\{2{\displaystyle \frac{m_q}{ϵ_q(p)}}{\displaystyle \frac{m_q(M_Vm_Q)}{ϵ_q^2(p)}}{\displaystyle \frac{m_q}{6M_Vm_Q}}{\displaystyle \frac{𝐩^2}{ϵ_q(p)}}`$ (40)
$`+{\displaystyle \frac{m_q}{2}}{\displaystyle \frac{𝐩^2}{ϵ_q^2(p)}}({\displaystyle \frac{1}{m_Q}}{\displaystyle \frac{2}{3[ϵ_q(p)+m_q]}}{\displaystyle \frac{2ϵ_q(p)}{3M_V[ϵ_q(p)+m_q]}})\}`$
$`{\displaystyle \frac{e_Q}{2m_Q}}\{2{\displaystyle \frac{M_Vϵ_q(p)}{m_Q}}{\displaystyle \frac{𝐩^2}{6m_Q}}({\displaystyle \frac{1}{m_Q}}+{\displaystyle \frac{1}{M_V}})`$
$`{\displaystyle \frac{1}{3M_V}}{\displaystyle \frac{𝐩^2}{ϵ_q(p)+m_q}}\}.`$
Here $`\mathrm{}`$ denotes the matrix element between radial meson wave functions. For these matrix element calculations we use the wave functions of heavy-light mesons obtained in Ref. . It is important to note that in this reference while calculating the heavy-light meson mass spectra only the heavy quark was treated using the $`1/m_Q`$ expansion but the light quark was treated completely relativistically.

The values of decay rates of mesons with open flavour calculated on the basis of Eqs. (1), (33), (40) are displayed in Table 1. In the second column ($`\mathrm{\Gamma }^{\mathrm{NR}}`$) we give predictions for decay rates obtained in the nonrelativistic approximation ($`p/m0`$) for both heavy and light quarks. In the third ($`\mathrm{\Gamma }^V`$) and fourth ($`\mathrm{\Gamma }^S`$) columns we show the results obtained for the purely vector and scalar confining potentials, respectively. And in the last column ($`\mathrm{\Gamma }`$) we present predictions for the mixture of vector and scalar confining potentials (8) with the mixing parameter $`\epsilon =1`$. As seen from this Table relativistic effects significantly influence the predictions. Their inclusion results in a significant reduction of decay rates ($`\mathrm{\Gamma }^{\mathrm{NR}}/\mathrm{\Gamma }=2÷4.5`$). Both relativistic corrections to the heavy quark and the relativistic treatment of the light quark play an important role. The dominant decay modes of $`D^{}`$ mesons are the strong decay $`D^{}D\pi `$, which is considerably suppressed by the phase space, and the electromagnetic decay $`D^{}D\gamma `$. The corresponding branching ratios are known already for a long time and listed in PDG tables . However, the total decay rates of $`D^{}`$ mesons were not measured until recently. In Ref.  CLEO collaboration reported the first measurement of the $`D^+`$ decay width $`\mathrm{\Gamma }(D^+)=96\pm 4\pm 22`$ keV. Combining this value with the measured $`BR(D^+D^+\gamma )=(1.6\pm 0.4)\%`$ , the following experimental value of the decay rate can be obtained: $`\mathrm{\Gamma }(D^+D^+\gamma )=(1.5\pm 0.6)`$ keV. Our model prediction is in agreement with this experimental value. However, the experimental errors are still large in order to discriminate the relativistic and nonrelativistic results. In the case of $`B`$ mesons the pion emission is kinematically forbidden, so the dominant decay mode is electromagnetic decay $`B^{}B\gamma `$. None of $`B^{}`$ widths has been measured yet. We also present our predictions for the radiative decay rate of $`B_c^{}`$ meson, which consists of two heavy quarks ($`b`$ and $`c`$). Therefore the expressions (33) and (40) can be further expanded in inverse powers of both quark masses up to the second order.

In Table 2 we compare our predictions for radiative decay rates of vector heavy-light mesons with other theoretical results. We show the predictions obtained in quark models , in the framework of heavy quark effective theory (HQET) combined with vector meson dominance (VMD) hypothesis  and in QCD sum rules . These predictions vary quite significantly from each other. Our predictions are in rough agreement with the quark model calculations of Ref. , with HQET+VMD results of Ref.  and with some of the predictions of the QCD sum rules. 111The QCD sum rule results  for the $`D^{}`$ decays were obtained using the $`1/m_c`$ expansion, which could be inaccurate due to the large value of the $`1/m_c`$ corrections. It is important to note that in our calculations we do not need to introduce the anomalous electromagnetic moment of the light quark as it is done in Ref. , where it was found that its value should be rather large ($`0.5`$) in order to get agreement with the experimental (CLEO) value for $`D^+D^+\gamma `$ decay rate. 222In Table 2 we show predictions of Ref.  for the value of anomalous electromagnetic quark moment equal to zero. The large value of the anomalous electromagnetic moment is not justified phenomenologically (see e. g. Refs. ). The other differences of our calculations from those of Ref.  are the Lorentz structure of the confining potential and a more comprehensive account of relativistic effects. In particular, the relativistic transformation of the meson wave function from the rest frame to a moving one given by Eq. (14) is missing in Ref. .

In Table 3 we present the comparison of our results with the predictions of different quark models  for the rates of the radiative M1-transitions ($`1^3S_11^1S_0+\gamma `$) in the heavy-heavy $`B_c`$ meson. There we also give the predicted values of the photon energy, which is determined by the mass splitting of the vector and pseudoscalar ground states. In previous calculations  the nonrelativistic expression for the matrix element of the magnetic moment was used. We see that even in the heavy-heavy $`B_c`$ meson inclusion of the relativistic effects results in considerable reduction of the radiative M1-decay rate. As can be seen from Table 1, this reduction is evoked by significant contributions of relativistic effects for the $`c`$ quark, since it is not heavy enough, as well as by the special choice of the mixture of vector and scalar confining potentials in our model (8).

In summary we calculated radiative M1-decay rates of mesons with open flavour in the framework of the relativistic quark model. In our analysis the light quark was treated relativistically, while for the heavy quark the $`1/m_Q`$ expansion was carried out up to the second order. Relativistic consideration of the light quark, relativistic heavy quark corrections as well as Lorentz-structure of the confining potential considerably influence the predictions. We find that only the mixture of vector and scalar confining potentials (8), with the mixing coefficient fixed previously from quarkonium radiative decays  and weak decays of heavy-light mesons , is in agreement with recent CLEO data for the $`D^+D^+\gamma `$ decay rate. More precise measurement of this decay rate and the measurement of radiative M1-decays of other heavy-light mesons will be crucial for testing the relativistic quark dynamics.

The authors express their gratitude to M. Müller-Preussker and V. Savrin for support and discussions. We are grateful to S. Narison for the useful comment on the QCD sum rule results. Two of us (R.N.F and V.O.G.) were supported in part by the Deutsche Forschungsgemeinschaft under contract Eb 139/2-1, Russian Foundation for Fundamental Research under Grant No. 00-02-17768 and Russian Ministry of Education under Grant No. E00-3.3-45.