text:
# High strangeness dibaryons in the extended quark delocalization, color screening model

## I Introduction

One of the most remarkable achievements of theoretical physics in the past thirty years is the establishment and development of the fundamental theory of the strong interaction – quantum chromodynamics (QCD). Perturbative QCD has been verified by high energy experiments. However, the low energy physics of QCD, such as hadron structure, hadron interactions and the structure of exotic quark-gluon systems, are much harder to calculate directly from QCD. One needs effective theories and phenomenological models in these cases.

Borrowing the idea of quasiparticles from condensed matter and nuclear physics, one can approximately transform the complicated interactions between current quarks into dynamic properties of quasiparticles and what is left to be studied are the residual interactions between quasiparticles. One of the quasiparticles in QCD is the constituent quark. How to dress the current quark to be a constituent quark still poses a theoretical challenge; various effective theoriescahill ; tandy ; roberts  have been developed to derive constituent quarks from QCD. The common point of view is that the dynamical generation of the constituent quark mass is closely related to spontaneous chiral symmetry breaking initiated by the formation of a $`q\overline{q}`$ condensate in the QCD vacuum.

The constituent quark model has been quite successful in understanding hadron spectroscopy and hadron interactions even though we have not yet derived the constituent quark model directly from QCD. There are also various versions of the constituent quark model based on different effective degrees of freedom. De Rujula, Georgi and Glashowprd12  first put forward a quark-gluon coupling model based on constituent quark and gluon effective degrees of freedom. Isgur obtained a good description of hadron spectroscopy based on this modelisgur . However, extension of the model to baryon interactions does not reproduce the nucleon-nucleon intermediate and long range interaction.

One modification studied is the addition of Goldstone boson exchange on the quark levelfaesslerplb124 ; obukhoplb238 ; fernjpg19 ; fujiprl , in which the short-range part of the interaction is described by the quark-gluon degrees of freedom and the medium and long-range parts are attributed to meson-exchange. This quark gluon-meson exchange hybrid model achieves a quantitative fit of the nucleon-nucleon (NN) and nucleon-hyperon (NY) scattering data.

A different modification of the De Rujula-Georgi-Glashow-Isgur (RGGI) model, the quark delocalization and color screening model (QDCSM)prl69 , has also been suggested. It maintains the Isgur Hamiltonian for single hadrons but modifies it for baryon-baryon (BB) interactions with two new ingredients.

First, the two center single quark orbital wave function (WF) used in the quark cluster model is replaced by a delocalized quark orbital WF. The introduction of quark delocalization can be viewed either as taking into account the contribution of excited configurations or the distortion of each individual baryon due to their mutual interaction. This straightforward method enlarges the variational Hilbert space. Its advantage is that it permits the six-quark system to choose a more favorable configuration through its own dynamics, while maintaining a tolerable level of computational complexity.

Second, a different parametrization of the confinement interaction is assumed, in which the usual quadratic confinement and a color-screened quadratic confinement are used to parameterize the two body matrix elements of different quark orbits. The introduction of this parametrization is aimed at taking into account some of the nonlinear, nonperturbative properties of QCD which can not be described by two-body quark-quark interactions in multiquark systems, such as the formation of color flux tubes connecting many quarks, the three gluon interaction and the three body instanton interaction.

The QDCSM not only justifies the view of a nucleus as, approximately, a collection of nucleons rather than a single big bag with 3A quarks, but it also explains the long standing fact that the nuclear force and molecular forces are similar except for the obvious energy and length scale differences. It is also the model which requires the fewest adjustable parameters to fit the existing BB interaction dataprl69 ; prc53 ; npa673 ; prc62 ; prc65 ; prc65p ; lucpl03 .

There have been debates on which constituent quark model and which effective degrees of freedom are best to use for hadron structure and interaction studiesprc57 ; kfliu ; isgurprd6162 ; georgiprd59 ; valcarceprc616364 . In our phenomenological study of BB interactions with three different constituent quark models, we found that even though the QDCSM and the other two models appear to be quite different, they give similar BB interactions in 44 of the 64 lowest BB channels consisting of octet and decuplet baryonsprc65 . A preliminary analysis of the origin of this surprising similarity has been producednucl-th/0212012 . This result also implies that quark delocalization and color screening, working together, do provide the intermediate range attraction described by meson exchange in other models. On the other hand, the different models do give characteristically different results in some channels. For example, the QDCSM predicts a strong attraction in the I=0, J=3 channel, producing a dibaryon resonance, $`d^{}`$, in this channel, while the quark gluon-meson exchange hybrid model predicts a strong attraction in the I=0, J=0 $`\mathrm{\Omega }\mathrm{\Omega }`$ channel, implying the existence of a strong interaction stable di-$`\mathrm{\Omega }`$. One cannot expect scattering data to become available in these channels to test these model predictions, but the dibaryon states should be detectable to provide a check on whether these model predictions are realistic.

The study of dibaryon states not only checks constituent quark models but also searches for new hadronic matter. The $`H`$-particle has been assumed to be a six quark system from the very beginningjaffe  and has been both a theoretical and experimental topic for a long time. The S=0, $`J^p=0^{}`$ d’ dibaryon was assumed to be an NN$`\pi `$ system and was a hot topic in the 1990’sbilger ; prc62 . We showed that the $`S=0`$, $`I=0`$, $`J=3`$ $`d^{}`$prc39 ; prl69 ; prc51 ; mpla13 ; npa688 ; prc65 ; prc65p  is a tightly bound six quark system rather than a loosely bound nuclear-like system of two $`\mathrm{\Delta }`$’s. An $`S=3,I=1/2,J=2`$ N$`\mathrm{\Omega }`$ state was proposed as a high strangeness dibaryon candidategoldman . Kopeliovich predicted high strangeness dibaryons, such as the di-$`\mathrm{\Omega }`$ with $`S=6`$, using the flavor SU(3) Skyrmion modelkopelnpa639 , and Zhang suggested searching for the di-$`\mathrm{\Omega }`$ in ultrarelativistic heavy ion collisions065204 . Lomon predicted a deuteron-like dibaryon resonance using R-matrix theorylomon  and measurements at Saclay seem to offer experimental supportSacley  for its existence.

In this paper, we calculate promising dibaryons with strangeness $`S=2,3,4,5`$ using the extended quark-delocalization color-screening model to provide another reference spectrum for strange dibaryons and a possible further check of constituent quark models. The SIJ=003 $`d^{}`$ and -600 di-$`\mathrm{\Omega }`$ dibaryons, which we considered previously, are included in the discussion to demonstrate the differences between the predictions of the extended QDCSM and other quark models.

The extended QDCSM is briefly introduced in Section II. In Section III, we present our results. We discuss these results further in Sections IV and V, and conclude in Section VI.

## II Brief description of the extended QDCSM

The QDCSM was put forward in the early 90’s. Details can be found in Refs.prl69 ; prc51 ; npa688 . Although the intermediate range attraction of the $`NN`$ interaction is reproduced by the combination of quark delocalization and color screening, the effect of the long-range pion tail is missing in the QDCSM. Recently, the extended QDCSM was developedprc65p , which incorporates this long-range tail by adding $`\pi `$-exchange with a short-range cutoff. The extended QDCSM not only reproduces the properties of the deuteron well, but also improves agreement with NN scattering data as compared to previous worklucpl03 .

The Hamiltonian of the extended QDCSM, wave functions and the necessary equations used in the current calculation are given below. Here we do not take into account the effect of any tensor forces. The details of the resonating-group method (RGM) have been presented in Refs.npa688 ; Buchmann .

The Hamiltonian for the 3-quark system is the same as the well known quark potential model, the Isgur model. For the six-quark system, we assume

$`H_6`$ $`=`$ $`{\displaystyle \underset{i=1}{\overset{6}{}}}(m_i+{\displaystyle \frac{p_i^2}{2m_i}})T_{CM}+{\displaystyle \underset{i<j=1}{\overset{6}{}}}\left[V_{conf}(r_{ij})+V_G(r_{ij})+V_\pi (r_{ij})\right],`$
$`V_G(r_{ij})`$ $`=`$ $`\alpha _s{\displaystyle \frac{\stackrel{}{\lambda }_i\stackrel{}{\lambda }_j}{4}}\left[{\displaystyle \frac{1}{r_{ij}}}{\displaystyle \frac{\pi }{2}}\delta (\stackrel{}{r_{ij}})\left({\displaystyle \frac{1}{m_i^2}}+{\displaystyle \frac{1}{m_j^2}}+{\displaystyle \frac{4\stackrel{}{\sigma }_i\stackrel{}{\sigma }_j}{3m_im_j}}\right)\right],`$
$`V_\pi (r_{ij})`$ $`=`$ $`\theta (rr_0){\displaystyle \frac{g_8^2}{4\pi }}{\displaystyle \frac{m_\pi ^2}{12m_q^2}}{\displaystyle \frac{1}{r_{ij}}}e^{m_\pi r_{ij}}\stackrel{}{\sigma }_i\stackrel{}{\sigma }_j\stackrel{}{\tau }_i\stackrel{}{\tau }_j,`$ (1)
$`V_{conf}(r_{ij})`$ $`=`$ $`a_c\stackrel{}{\lambda }_i\stackrel{}{\lambda }_j\{\begin{array}{cc}r_{ij}^2\hfill & \text{if }i,j\text{ occur in the same baryon orbit},\hfill \\ \frac{1e^{\mu r_{ij}^2}}{\mu }\hfill & \text{if }i,j\text{ occur in different baryon orbits},\hfill \end{array}`$ (4)
$`\theta (r_{ij}r_0)`$ $`=`$ $`\{\begin{array}{cc}0\hfill & r_{ij}<r_0,\hfill \\ 1\hfill & \text{otherwise},\hfill \end{array}`$ (7)
where $`r_0`$ is the short range cutoff for pion exchange between quarks. All the symbols have their usual meanings, and the confinement potential $`V_{conf}(r_{ij})`$ has been discussed in Refs.prc65p ; npa688 .

The pion potential, $`V_\pi (r_{ij})`$, affects only the $`u`$ and $`d`$ quarks. We take these to have a common mass, $`m_q=m_d=m_u`$, ignoring isospin breaking effects as they are small on the scale of interest here.

The quark wave function in a given nucleon (orbit) relative to a reference center (defined by $`\stackrel{}{S}`$) is taken to have a Gaussian form characterized by a size parameter, $`b`$,

$$\varphi (\stackrel{}{r}\stackrel{}{S})=\left(\frac{1}{\pi b^2}\right)^{3/4}e^{\frac{1}{2b^2}(\stackrel{}{r}\stackrel{}{S})^2}.$$
(8)
The values of $`m_q`$, $`m_s`$, $`b`$, $`\alpha _s`$ and $`a_c`$ are determined by reproducing the $`\mathrm{\Delta }N`$ mass difference, the nucleon mass, a hyperon mass and by requiring a stability condition. The quark-pion coupling constant $`g_{qq\pi }`$ is obtained from the nucleon-pion coupling constant by a slight ($`<10`$%) correction to the classic symmetry relation, viz.,

$$\frac{g_{NN\pi }^2}{4\pi }=(M_N/m_q)^2\left(\frac{5}{3}\right)^2\frac{g_8^2}{4\pi }e^{m_\pi ^2b^2/2}$$
(9)
where $`M_N`$ is the nucleon mass and the last factor provides the correction due to the extent of the quark wavefunction in the nucleon. The color screening parameter, $`\mu `$, has been determined by matching our calculation to the mass of the deuteron. All of the parameters are listed in Table I.

Table I: Model Parameters

| $`m_q,m_s(MeV)`$ | $`b(fm)`$ | $`a_c(MeVfm^2)`$ | $`\alpha _s`$ | $`\frac{g_8^2}{4\pi }`$ | $`r_0(fm)`$ | $`\mu (fm^2)`$ |
| --- | --- | --- | --- | --- | --- | --- |
| $`313,634`$ | $`0.6015`$ | $`25.14`$ | $`1.5585`$ | $`0.5926`$ | $`0.8`$ | $`0.85`$ |

The model masses of all octet and decuplet baryons are listed in Table II.

Table II: Single Baryon Masses in Units of MeV

|  | N | $`\mathrm{\Sigma }`$ | $`\mathrm{\Lambda }`$ | $`\mathrm{\Xi }`$ | $`\mathrm{\Delta }`$ | $`\mathrm{\Sigma }^{}`$ | $`\mathrm{\Xi }^{}`$ | $`\mathrm{\Omega }`$ |
| --- | --- | --- | --- | --- | --- | --- | --- | --- |
| theor. | 939.0 | 1210.6 | 1113.6 | 1350.0 | 1232.0 | 1358.1 | 1497.5 | 1650.1 |
| expt. | 939 | 1193 | 1116 | 1318 | 1232 | 1385 | 1533 | 1672 |

We use the resonating group method to carry out a dynamical calculation. Introducing Gaussian functions with different reference centers $`S_i`$ i=1…n, which play the role of generating coordinates in this formalism, to expand the relative motion wave function of the two quark clusters, we have

$$\chi (\stackrel{}{R})=(\frac{3}{2\pi b^2})^{3/4}\underset{i}{}C_ie^{\frac{3}{4}(\stackrel{}{R}\stackrel{}{S}_i)^2/b^2}.$$

In principle, any set of base wave functions can be used to expand the relative motion wave function. The choice of a Gaussian with the same size parameter, $`b`$, as the single quark wave function given in Eq.(8), however, allows us to rewrite the resonating group wave function as a product of single quark wave functions; (see Eq.(10) below). This cluster wave function (physical basis) can be expressed in terms of the symmetry basis, classified by the symmetry properties, in a group chain which in turn allows the use of group theory methods to simplify the calculation of the matrix elements of the six quark Hamiltonianwpg . In our calculations, we typically use 12 Gaussian functions to expand the relative motion wave function over the range 0-8 fm. For a few channels, such as the deuteron and $`H`$ particle, 20 Gaussian functions are used to extend the range to 12 fm.

After including the wave function for the center-of-mass motion, the ansatz for the two-cluster wave function used in the RGM can be written as

$`\mathrm{\Psi }_{6q}`$ $`=`$ $`𝒜{\displaystyle \underset{k}{}}{\displaystyle \underset{i=1}{\overset{n}{}}}C_{k,i}{\displaystyle 𝑑\mathrm{\Omega }_{S_i}\underset{\alpha =1}{\overset{3}{}}\psi _R(\stackrel{}{r}_\alpha ,\stackrel{}{S_i},ϵ)\underset{\beta =4}{\overset{6}{}}\psi _L(\stackrel{}{r}_\beta ,\stackrel{}{S_i},ϵ)}`$ (10)
$`[\eta _{I_{1k}S_{1k}}(B_{1k})\eta _{I_{2k}S_{2k}}(B_{2k})]^{I,J=S}[\chi _c(B_1)\chi _c(B_2)]^{[\sigma ]},`$
where $`k`$ is the channel index. For example, for $`SIJ=2,0,0`$, we have $`k=1,2,3`$, corresponding to the channels $`\mathrm{\Lambda }\mathrm{\Lambda }`$, N$`\mathrm{\Xi }`$ and $`\mathrm{\Sigma }\mathrm{\Sigma }`$.

The delocalized orbital wavefunctions, $`\psi _R(\stackrel{}{r},\stackrel{}{S_i},ϵ)`$ and $`\psi _L(\stackrel{}{r},\stackrel{}{S_i},ϵ)`$, are given by

$`\psi _R(\stackrel{}{r},\stackrel{}{S_i},ϵ)`$ $`=`$ $`{\displaystyle \frac{1}{N(ϵ)}}\left(\varphi (\stackrel{}{r}{\displaystyle \frac{\stackrel{}{S_i}}{2}})+ϵ\varphi (\stackrel{}{r}+{\displaystyle \frac{\stackrel{}{S_i}}{2}})\right),`$
$`\psi _L(\stackrel{}{r},\stackrel{}{S_i},ϵ)`$ $`=`$ $`{\displaystyle \frac{1}{N(ϵ)}}\left(\varphi (\stackrel{}{r}+{\displaystyle \frac{\stackrel{}{S_i}}{2}})+ϵ\varphi (\stackrel{}{r}{\displaystyle \frac{\stackrel{}{S_i}}{2}})\right),`$ (11)
$`N(ϵ)`$ $`=`$ $`\sqrt{1+ϵ^2+2ϵe^{S_i^2/4b^2}},`$
where $`\varphi (\stackrel{}{r}\frac{\stackrel{}{S_i}}{2})`$ and $`\varphi (\stackrel{}{r}+\frac{\stackrel{}{S_i}}{2})`$ are the single-particle Gaussian quark wave functions referred to above in Eq.(8), with different reference centers $`\frac{S_i}{2}`$ and $`\frac{S_i}{2}`$, respectively. The delocalization parameter, $`ϵ`$, is determined by the dynamics of the quark system rather than being treated as an adjustable parameter.

The initial RGM equation is

$$H(\stackrel{}{R},\stackrel{}{R^{}})\chi (\stackrel{}{R^{}})𝑑\stackrel{}{R^{}}=EN(\stackrel{}{R},\stackrel{}{R^{}})\chi (\stackrel{}{R^{}})𝑑\stackrel{}{R^{}}.$$
(12)
With the above ansatz, the RGM Eq.(12) is converted into an algebraic eigenvalue equation,

$$\underset{j,k^{}}{}C_{j,k^{}}H_{i,j}^{k,k^{}}=E\underset{j}{}C_{j,k}N_{i,j}^k,$$
(13)
where $`N_{i,j}^k,H_{i,j}^{k,k^{}}`$ are the wave function overlaps and Hamiltonian matrix elements, respectively, obtained for the wave functions of Eq.(10).

## III Results and Discussion

Previously, we chose the di-$`\mathrm{\Omega }`$ as an example to study whether or not our model results were sensitive to the meson-exchange cut-off parameter, $`r_0`$, and the result demonstrates that they are notcommun38 . Hence, we consider it sufficient to calculate six-quark systems of different strangeness, $`S=2,3,4,5`$, with a representative cutoff value of $`r_0=0.8fm`$. Table III displays the masses (in MeV) calculated for the strange dibaryon states of interest here. The lowest (without taking tensor coupling into account) channel for each SIJ is identified by bold lettering; $`sc`$ and $`cc`$ denote single channel (the lowest one) and coupled channels, respectively. ’\*’ denotes those states without coupled channels. Additionally, it should be noted that in our calculation we assume the wavefunction to be zero at the boundary point, which is the usual boundary condition for bound states. If the state is unbound, we will not obtain a stable minimum eigenenergy in the course of extending the boundary point. The unbound states are denoted by ”-” in Table III.

Table III : Masses of Six-Quark Systems with Strangeness    S,I,J coupling channels Masssc Masscc    $`2,0,0`$ $`𝚲𝚲`$-$`N\mathrm{\Xi }\mathrm{\Sigma }\mathrm{\Sigma }`$  $`2225.5`$  $`3,1/2,2`$ $`𝐍𝛀\mathrm{\Sigma }\mathrm{\Xi }^{}\mathrm{\Xi }\mathrm{\Sigma }^{}\mathrm{\Xi }^{}\mathrm{\Lambda }\mathrm{\Xi }^{}\mathrm{\Sigma }^{}`$ $`2566.4`$ $`2549.1`$  $`3,1/2,1`$ $`𝚲𝚵N\mathrm{\Omega }\mathrm{\Lambda }\mathrm{\Xi }^{}\mathrm{\Xi }\mathrm{\Sigma }^{}\mathrm{\Sigma }^{}\mathrm{\Xi }^{}\mathrm{\Sigma }\mathrm{\Xi }^{}\mathrm{\Sigma }\mathrm{\Xi }`$    $`4,1,0`$ $`𝚵𝚵\mathrm{\Sigma }^{}\mathrm{\Omega }\mathrm{\Xi }^{}\mathrm{\Xi }^{}`$    $`4,0,1`$ $`𝚵𝚵\mathrm{\Lambda }\mathrm{\Omega }\mathrm{\Xi }\mathrm{\Xi }^{}\mathrm{\Xi }^{}\mathrm{\Xi }^{}`$    $`5,1/2,0`$ $`𝚵^{}𝛀`$ $`3145.0`$ *  $`5,1/2,1`$ $`𝚵𝛀\mathrm{\Xi }^{}\mathrm{\Omega }`$    $`6,0,0`$ $`𝛀𝛀`$ $`3298.2`$ *

In 1977, Jaffejaffe  studied the color-magnetic interaction of the one-gluon exchange potential in the multiquark system and found that the most attractive channel is the flavor singlet with quark content $`u^2d^2s^2`$. Moreover, the same symmetry analysis of the chiral boson exchange potential also leads to the very same conclusionnucl-th/0212012 .

However dibaryon physics can be very delicate. The deuteron channel is not a channel with strong attraction in any baryon interaction model. If the deuteron had not been found experimentally, it seems highly unlikely that any model would have been able to predict it to be a stable dibaryon. The $`H`$-particle ($`SIJ=200`$) is a six quark state consisting mainly of octet-baryons, similar to the deuteron, and we find only a weak attraction there in our model also. Hence, a qualitative analysis is insufficient to judge whether or not the $`H`$-particle is strong interaction stable. Systematically, we find that a strong attraction develops only in decuplet-decuplet channels and a mild attraction in octet-decuplet channels.

Moreover, in the $`H`$-particle case, the channel coupling effect may even be more important than the deuteron case. In fact, it is not bound without taking coupled channels into account. In our calculation, three channels have been taken into account. These are: $`\mathrm{\Sigma }\mathrm{\Sigma }`$, N$`\mathrm{\Xi }`$ and $`\mathrm{\Lambda }\mathrm{\Lambda }`$. The relative motion wave functions of each channel contribution to the $`H`$-particle are shown in Fig 1. We find that the $`\mathrm{\Lambda }\mathrm{\Lambda }`$ channel provides the largest contribution (67%), followed by the N$`\mathrm{\Xi }`$ channel (23%); the $`\mathrm{\Sigma }\mathrm{\Sigma }`$ channel contributes only (10%). However, according to the analysis by Jaffe, the biggest contribution is the $`N\mathrm{\Xi }`$ channel, and the $`\mathrm{\Lambda }\mathrm{\Lambda }`$ channel provides the smallest contribution.

The lowest mass we find for the $`u^2d^2s^2`$ system is $`2225.5MeV`$, which is $`6MeV`$ lower than the experimental threshold of $`\mathrm{\Lambda }\mathrm{\Lambda }`$ and $`1.7MeV`$ lower than our model threshold. These values are smaller than the rms uncertainty that may be inferred from our fit to the baryon octet and decuplet in Table II. Furthermore, these values are on the order of corrections one would expect from isospin violations which we have not included. Hence, we can draw no definite conclusion as to whether or not the $`H`$-particle is strong interaction stable in our model and we consider this to be consistent with recent experimental findingstak .

Besides the binding energy of the $`H`$, an interesting question regarding the $`H`$ is its compactness, i.e., whether the $`H`$ is a compact 6-quark object or a loosely bound $`\mathrm{\Lambda }\mathrm{\Lambda }`$ state. Fig.1 indicates that the maxima of the relative motion wave function of the dominant $`\mathrm{\Lambda }\mathrm{\Lambda }`$ channel occurs around 1.8 fm, and the delocalization parameter, $`ϵ`$, of the dominant channel at the maximum is $``$0.1; hence, the $`H`$ is a loosely bound system similar to the deuteron in our model. Such a similarity may well be more physically sensible than a compact 6-quark structure.

For systems with strangeness $`S=3`$, we have calculated the state N$`\mathrm{\Omega }`$ (SIJ=-3,1/2,2), which was shown to be mildly attractive, with energy below $`\mathrm{\Lambda }\mathrm{\Xi }\pi `$ thresholdgoldman . That conclusion was challenged by Okaokaprd38  and supported by Silvestre-Brac and Leandrislprd45 .

We have carried out a dynamical channel coupling calculation to examine this state further. The $`N\mathrm{\Omega },\mathrm{\Lambda }\mathrm{\Xi }^{},\mathrm{\Xi }\mathrm{\Sigma }^{},\mathrm{\Sigma }\mathrm{\Xi }^{},\mathrm{\Xi }^{}\mathrm{\Sigma }^{}`$ channels are all included. We find the eigenenergy to be 2549.1 MeV, $`24(54)MeV`$ lower than the $`\mathrm{\Lambda }\mathrm{\Xi }\pi `$ experimental(model) threshold. Relative motion cluster wave functions for the individual channels are shown in Fig. 2. The N$`\mathrm{\Omega }`$ channel is by far the dominant one (77%) and the maximum for its relative motion wave function occurs at around 0.8 fm. The value of $`ϵ`$ at this separation is 0.46. Mixing into the other channels is small. Hence, we find this to be a compact six quark state.

Note that the D-wave $`\mathrm{\Lambda }\mathrm{\Xi }`$ and $`\mathrm{\Sigma }\mathrm{\Xi }`$ channels have not been included since no tensor interaction has been included in this calculation. This coupling should be weak because there is no $`\pi `$exchange in these channels and K and $`\eta `$ exchanges are not important in our model (see next section). Its effect on the eigenenergy of the $`N\mathrm{\Omega }`$ (SIJ=-3,1/2,2) state should be small and the D-wave decay widths to $`\mathrm{\Lambda }\mathrm{\Xi }`$ and $`\mathrm{\Sigma }\mathrm{\Xi }`$ final states should be small also. Therefore we expect the $`SIJ=3,1/2,2`$ $`N\mathrm{\Omega }`$ should be a narrow dibaryon resonance. The tensor coupling calculation is currently in progress.

We have also calculated the state $`SIJ=3,1/2,1`$. In the course of extending the boundary point, the lowest eigenenergy closely approaches the model threshold but fails to come to a stable value within our limits of computation for the extension. This result is unchanged by taking into account all possible coupling channels. Due to the fact that there are only weak attractions in our model for the octet-octet channels and the size of the model uncertainty, it is difficult to conclude whether or not there are strong interaction stable states in these channels. We can only conclude that we have not found evidence for a strong interaction stable state with $`SIJ=3,1/2,1`$.

For systems with S=-4, we take the quantum numbers $`SIJ=4,1,0`$ as an example, because this case bears a number of similarities to the deuteron. For example, in both cases, the lowest mass channel is composed of two octet baryons from the same isodoublet. Also, the matrix element $`P_{36}^{sfc}`$ characterizing the symmetry property of the system is -1/81 for both cases. ($`P_{36}^{sfc}`$ is the permutation operator of the quarks between two clusters acting in spin, flavor and color space.) The result shows that the system with S=-4, I=1, J=0 is unbound, even when the $`\mathrm{\Xi }^{}\mathrm{\Xi }^{}`$ and $`\mathrm{\Sigma }^{}\mathrm{\Omega }`$ channel couplings are taken into account. Because the calculated energy is again very close to the $`\mathrm{\Xi }\mathrm{\Xi }`$ threshold this conclusion should be viewed as tentative.

For comparison, we have also calculated the $`SIJ=4,0,1`$ state as shown in Table III. The $`\mathrm{\Xi }\mathrm{\Xi }`$, $`\mathrm{\Xi }\mathrm{\Xi }^{}`$, $`\mathrm{\Lambda }\mathrm{\Omega }`$ and $`\mathrm{\Xi }^{}\mathrm{\Xi }^{}`$ coupling channels are included. Our model result is very similar to the $`SIJ=4,1,0`$, i.e., we do not find a bound state in this channel. The $`{}_{}{}^{3}S_{1}^{}^3D_1`$ tensor coupling should be small in the $`\mathrm{\Xi }\mathrm{\Xi }`$ channel so taking into account the tensor coupling will not change the unbound character.

For the system with S=-5, we take the $`SIJ=5,1/2,0`$ state as an example. This state is interesting due to its $`<P_{36}^{sfc}>`$ value, $`1/9`$, which makes it a Pauli principle favored state. If only two-baryon S-wave channels are taken into account, there is only one channel for this state. Our calculation shows that the contribution of the kinetic energy term, due to quark exchange and delocalization effects, contributes strongly towards the formation of a bound state. However, the one-gluon-exchange interaction largely compensates for this attraction and produces a mass of $`3145.0MeV`$, which is $`59MeV`$ lower than the experimental value of the $`\mathrm{\Xi }^{}\mathrm{\Omega }`$ threshold but only about $`2MeV`$ lower than the model threshold. (The one-gluon-exchange effect here is quite different from that in the $`d^{}`$ case, where large delocalization is favored for a wide range of cluster separations as well as there being a strong effective attraction due to the large reduction in the kinetic energy that accompanies significant delocalization.) We conclude that this state is not a good candidate for a dibaryon resonance search due to its small binding and its $`\mathrm{\Xi }^{}\mathrm{\Omega }`$ content.

Inclusion of the tensor interaction will mix the spin-2 D-wave $`\mathrm{\Xi }\mathrm{\Omega }`$ channel with the spin-0 S-wave $`\mathrm{\Xi }^{}\mathrm{\Omega }`$ channel. The tensor coupling should be weak also and we expect that its effect on the $`SIJ=5,1/2,0`$ state is small.

In the same strangeness sector, we also calculated the $`SIJ=5,1/2,1`$ state, since it includes the lowest channel, $`\mathrm{\Xi }\mathrm{\Omega }`$. The calculated energy is very close to the $`\mathrm{\Xi }\mathrm{\Omega }`$ threshold but a little higher so again there is no bound state, with $`SIJ=5,1/2,1`$, in our model.

To sum up, there are only a few high strangeness states worthy of experimental searches in our model. These are the $`H`$ particle, the N$`\mathrm{\Omega }`$ and the di-$`\mathrm{\Omega }`$. The di-$`\mathrm{\Omega }`$ was previously reported in Ref.(commun38 ) and the result is included in Table III.

The $`H`$-particle and the di-$`\mathrm{\Omega }`$ may be strong interaction stable. However in our model, the binding energies of both are small relative to the model uncertainty. The di-$`\mathrm{\Omega }`$ mass is about $`47MeV`$ lower than the experimental $`\mathrm{\Omega }\mathrm{\Omega }`$ threshold. However our model mass for the $`\mathrm{\Omega }`$ is $`1650MeV`$. If this model mass of $`\mathrm{\Omega }`$ were used to calculate the threshold, then the di-$`\mathrm{\Omega }`$ mass is no more than $`2MeV`$ below that threshold. Since our model mass for the single $`\mathrm{\Omega }`$ baryon deviates from the experimental value about $`22MeV`$, a reasonable estimate of the model uncertainty for the dibaryon would be at least that large. Therefore the di-$`\mathrm{\Omega }`$ should not be claimed as a strong interaction stable dibaryon within the model. Similarly, we cannot claim that the $`H`$-particle is strong interaction stable either. The $`SIJ=3,1/2,2`$ N$`\mathrm{\Omega }`$ case is certainly not strong interaction stable. However, the state is also certainly lower in mass than the N$`\mathrm{\Omega }`$ threshold, and quite possibly lower than the $`\mathrm{\Lambda }\mathrm{\Xi }\pi `$ threshold, as well. The tensor coupling to the $`\mathrm{\Lambda }\mathrm{\Xi }`$ and $`\mathrm{\Sigma }\mathrm{\Xi }`$ channels should be weak and the decay width should be small. This strongly suggests that it is a promising candidate for a narrow dibaryon resonance. This prediction can be tested by relativistic heavy ion reactions using RHIC detectors through the reconstruction of the vertex mass of the two body decay products, $`\mathrm{\Lambda }`$ and $`\mathrm{\Xi }`$.

## IV The $`K,\eta `$ meson exchange effect in the extended QDCSM

The effects of K and $`\eta `$ meson exchange have been studied in the di-$`\mathrm{\Omega }`$ channel and found to be negligible in our modelcommun38 . In this section, we carry out a further systematic study of the effect of ($`K,\eta `$) meson exchange on the masses of strange dibaryon candidates. A flavor-symmetric octet meson-quark coupling is assumed for all of the octet mesons ($`\pi ,K,\eta `$). The quark-meson exchange potential has the usual pseudoscalar meson exchange form except for a short range cut-off,

$$V_K(r_{ij})=\underset{a=4}{\overset{7}{}}\theta (rr_0)\frac{g_8^2}{4\pi }\frac{m_K^2}{12m_im_j}\frac{1}{r_{ij}}e^{m_Kr_{ij}}\stackrel{}{\sigma }_i\stackrel{}{\sigma }_j\stackrel{}{\lambda }_i^{f,a}\stackrel{}{\lambda }_j^{f,a},$$

$$V_\eta (r_{ij})=\theta (rr_0)\frac{g_8^2}{4\pi }\frac{m_\eta ^2}{12m_im_j}\frac{1}{r_{ij}}e^{m_\eta r_{ij}}\stackrel{}{\sigma }_i\stackrel{}{\sigma }_j\stackrel{}{\lambda }_i^{f,8}\stackrel{}{\lambda }_j^{f,8}.$$

A unified cutoff of $`r_0=0.8fm`$ is used. This is aimed at avoiding double counting since, in our model approachprc65 , short and intermediate range interactions have been accounted for by the combination of quark delocalization and color screening. The model parameters ($`b,a_c,\alpha _s`$) are refitted as before, and are listed in Table IV. The recalculated masses of octet and decuplet baryons are listed in Table V.

Table IV: Model Parameters

| $`m_q,m_s(MeV)`$ | $`b(fm)`$ | $`a_c(MeVfm^2)`$ | $`\alpha _s`$ | $`\frac{g_8^2}{4\pi }`$ | $`r_0(fm)`$ | $`\mu (fm^2)`$ |
| --- | --- | --- | --- | --- | --- | --- |
| $`313,634`$ | $`0.6022`$ | $`25.03`$ | $`1.5547`$ | $`0.5926`$ | $`0.8`$ | $`0.90`$ |

Table V: Single Baryon Masses in Units of MeV

|  | N | $`\mathrm{\Sigma }`$ | $`\mathrm{\Lambda }`$ | $`\mathrm{\Xi }`$ | $`\mathrm{\Delta }`$ | $`\mathrm{\Sigma }^{}`$ | $`\mathrm{\Xi }^{}`$ | $`\mathrm{\Omega }`$ |
| --- | --- | --- | --- | --- | --- | --- | --- | --- |
| theor. | 939.0 | 1217.47 | 1116.90 | 1357.56 | 1232.0 | 1359.61 | 1499.70 | 1652.27 |
| expt. | 939 | 1193 | 1116 | 1318 | 1232 | 1385 | 1533 | 1672 |

Comparing Tables I,II and IV,V, it is apparent that the addition of K and $`\eta `$ exchanges has modified the model parameters and single baryon masses only slightly. We also find that the properties of the deuteron are reproduced as well as before by a minor readjustment of the color-screening parameter $`\mu `$ commun38 .

With the full octet meson exchange, we recalculated every single channel case with different quantum numbers. For the $`H`$-particle, we calculated both the single channel and coupled channels, because it is very sensitive to channel coupling. The results are presented in Tables VI and VII. For comparison, the results with $`\pi `$ exchange only are also listed.

Table VI : Masses of Six-Quark Systems with $`\pi `$ exchange only and full octet pseudoscalar meson exchange. The masses are given in MeV.    S,I,J single channel  only $`\pi `$ exchange  $`\pi ,K,\eta `$ exchange    $`3,1/2,2`$ $`N\mathrm{\Omega }`$ $`2566.37`$ $`2556.95`$  $`3,1/2,1`$ $`\mathrm{\Lambda }\mathrm{\Xi }`$    $`4,1,0`$ $`\mathrm{\Xi }\mathrm{\Xi }`$    $`4,0,1`$ $`\mathrm{\Xi }\mathrm{\Xi }`$    $`5,1/2,1`$ $`\mathrm{\Xi }\mathrm{\Omega }`$    $`5,1/2,0`$ $`\mathrm{\Xi }^{}\mathrm{\Omega }`$ $`3145.01`$ $`3146.15`$  $`6,0,0`$ $`\mathrm{\Omega }\mathrm{\Omega }`$ $`3298.20`$ $`3300.00`$

Table VII : Mass of $`H`$-particle ($`SIJ=2,0,0`$) with $`\pi `$ exchange only and full octet pseudoscalar meson exchange. The masses are given in MeV.      only $`\pi `$ exchange  $`\pi ,K,\eta `$ exchange   $`\mathrm{\Sigma }\mathrm{\Sigma }`$ $`2280.40`$ $`2282.67`$  single channels $`N\mathrm{\Xi }`$ $`2263.95`$ $`2268.12`$   $`\mathrm{\Lambda }\mathrm{\Lambda }`$    coupled channels $`\mathrm{\Sigma }\mathrm{\Sigma }N\mathrm{\Xi }\mathrm{\Lambda }\mathrm{\Lambda }`$ $`2225.48`$ $`2230.28`$

Tables VI and VII show that those systems which are unbound in the extended QDCSM with $`\pi `$ exchange only, such as $`SIJ=3,1/2,1;4,1,0;4,0,1`$ and $`5,1/2,1`$, remain unbound after adding $`K`$ and $`\eta `$ exchanges, while those systems, which are bound in Table III, remain bound. Moreover, their masses are almost unaffected by the addition of $`K`$ and $`\eta `$ exchange. The largest difference is not more than $`10MeV`$. Especially for the $`H`$-particle, the $`\mathrm{\Lambda }\mathrm{\Lambda }`$ channel would be unbound without taking into account the $`N\mathrm{\Xi }`$ and $`\mathrm{\Sigma }\mathrm{\Sigma }`$ channel coupling; this character is not affected either by adding $`K`$ and $`\eta `$ exchange. Even the $`\mathrm{\Lambda }\mathrm{\Lambda }`$, $`N\mathrm{\Xi }`$ and $`\mathrm{\Sigma }\mathrm{\Sigma }`$ channel mixing fractions, $`56`$%, $`23`$% and $`21`$%, respectively, are very similar to the case with $`\pi `$ exchange only. All of these results confirm our expectation that heavier meson ($`K`$, $`\eta `$) exchange has already been mostly accounted for by the quark delocalization and color screening effects in our model approach. Hence, explicit inclusion of $`K`$ and $`\eta `$ exchanges beyond the cutoff scale is not important in our approach.

## V Further Discussion about the Mechanism of Intermediate Range Attraction

We have reported that the QDCSM gives very similar effective BB interactions to other modelsprc65 , in general. However a careful comparison found that our model results for the binding energy of high strangeness dibaryons are systematically smaller than those of the chiral quark modelzhangzy . The difference is mainly due to the different mechanism for the effective BB intermediate range attraction.

In the chiral quark modelfernjpg19 ; fujiprl ; zhangzy  the intermediate range attraction is attributed to $`\sigma `$ meson exchange. Because of its scalar-isoscalar character, $`\sigma `$ meson exchange provides a universal attraction independent of the flavor.

In the Bonn meson exchange modelmah1  the $`\sigma `$ meson is an effective description of correlated two $`\pi `$ exchange. This point has been confirmed by a chiral perturbation calculation of the NN interactionmah2 . Such an effective $`\sigma `$-N coupling can not be extended from the NN channel to other channels with strangeness by a universal $`\sigma `$-baryon coupling, because in the NN channel, there are NN, N$`\mathrm{\Delta }`$ and $`\mathrm{\Delta }\mathrm{\Delta }`$ intermediate states for two $`\pi `$ exchange, while for the $`\mathrm{\Lambda }\mathrm{\Lambda }`$ channel, there is only the $`\mathrm{\Sigma }\mathrm{\Sigma }`$ intermediate state, and for the $`\mathrm{\Omega }\mathrm{\Omega }`$ channel, no such intermediate state is possible. Therefore it is not justified to fix the parameters of $`\sigma `$ exchange in the NN channel and then directly extend that exchange to channels with strangeness.

There are arguments based on spontaneous chiral symmetry breaking for introduction of $`\sigma `$ meson quark coupling. In the SU(2) case, the non-linear realization of chiral symmetry can be linearized and in turn the $`\sigma `$ and $`\pi `$ meson coupling constants can be unifiedobukhoplb238 ; fernjpg19 . However one should note this universal coupling is restricted to the $`u`$ and $`d`$ quarks and the effective $`\sigma `$ is still due to even number multiple $`\pi `$ exchange. In the SU(3) case, one can introduce SU(3) chiral symmetry by neglecting the difference between the $`s`$ and light quark masses, followed by spontaneous chiral symmetry breaking. However, the nonlinear realization of SU(3) can not be linearized in the same way as in the SU(2) case to obtain a universal $`u`$, $`d`$ and $`s`$ quark $`\sigma `$ coupling. The internal quark structure of the $`\sigma `$ meson is a controversial issue, but an equal mix of $`u\overline{u}+d\overline{d}+s\overline{s}`$ is quite unexpected and so also is a universal $`u`$, $`d`$ and $`s`$ quark $`\sigma `$ coupling. Hence, a strong attraction in high strangeness channels arising from such a universal $`\sigma `$ quark coupling is quite questionable also.

In the QDCSM, quark delocalization and color screening work together to provide appropriate short-range repulsion and intermediate-range attraction for different channels. We illustrate this mechanism by showing contributions of the kinetic energy, confinement, color Coulomb, and color magnetic terms to the effective BB potential $`V_{BB}(S)`$, as well as the total sum, in curves a-e, respectively in Figs.3-6 for a few typical channels. The value of $`ϵ`$ varies with the separation $`S`$ and is also listed in Figs.3-6. The contribution of $`\pi ,K,\eta `$ exchange with a cut-off $`r_0=0.8fm`$ is small so we do not show it. These figures are ordered in terms of increasing values of $`<P_{36}^{sfc}>`$, from Pauli favored to Pauli forbidden.

The effective baryon-baryon interactions shown in these figures are obtained in the Born-Oppenheimer approximation,

$$V_{BB}(S_i)=\frac{\mathrm{\Psi }_{BB}(S_i)|H|\mathrm{\Psi }_{BB}(S_i)}{\mathrm{\Psi }_{BB}(S_i)|\mathrm{\Psi }_{BB}(S_i)}\frac{\mathrm{\Psi }_{BB}(S_i)|H|\mathrm{\Psi }_{BB}(S_i)}{\mathrm{\Psi }_{BB}(S_i)|\mathrm{\Psi }_{BB}(S_i)}|_{S_i\mathrm{},ϵ=0}.$$

Where the $`\mathrm{\Psi }_{BB}(S_i)`$ is the antisymmetric six quark cluster state at a specified separation, $`S_i`$, as given in Eq. (4) without summation over $`k`$ and $`i`$; $`H`$ is the six quark Hamiltonian (1). The contribution of each term in Eq. (1) is defined similarly, for example,

$$V_{BB}^{conf}(S_i)=\frac{\mathrm{\Psi }_{BB}(S_i)|V_{conf}|\mathrm{\Psi }_{BB}(S_i)}{\mathrm{\Psi }_{BB}(S_i)|\mathrm{\Psi }_{BB}(S_i)}\frac{\mathrm{\Psi }_{BB}(S_i)|V_{conf}|\mathrm{\Psi }_{BB}}{\mathrm{\Psi }_{BB}(S_i)|\mathrm{\Psi }_{BB}(S_i)}_{S_i\mathrm{},ϵ=0}.$$

Quark exchange alone induces a weak reduction in kinetic energy. Quark delocalization enhances this kinetic energy reduction. The kinetic energy reduction is also dependent on the strangeness of the channels due to the inverse quark mass dependence of the quark kinetic energy. The higher the strangeness, the smaller the contribution of quark kinetic energy and hence the smaller the reduction of the kinetic energy due to delocalization (see Figs.3-6(a)). This makes the intermediate range attraction weaker for the higher strangeness channels and is the main reason that our model gives less binding for the $`S=4,5,6`$ states than does the chiral quark modelzhangzy .

If the usual quadratic color confinement is used, it does not contribute to the effective BB interaction. With the introduction of quark delocalization, the usual color confinement contributes an effective repulsion, as shown in Fig.7. The color Coulomb has a similar behavior, as shown in Figs.3-6(c); it produces almost no contribution to the effective BB interaction without quark delocalization, and contributes an effective repulsion with quark delocalization. These two terms working together almost totally forbid quark delocalization. This implies that the internal structure of the baryon is unaffected by the mutual interaction. However, that is inconsistent with the observed difference between the nucleon structure function in a nucleus and in isolation as seen in deep inelastic lepton scattering (EMC effect).

On the other hand, there is no compelling reason to assume that the two body confinement potential is a good approximation for a multi-quark system. At the very least, the three gluon and the three body instanton interactions, which do not contribute to the $`q\overline{q}`$ meson and $`q^3`$ baryon but do contribute to the multi-quark system, have been omitted in every two body confinement model. The color flux scenario revealed in the lattice QCD calculation of two and three quark systems raises questions regarding an additive two body confinement approximation. To take these facts into account, the QDCSM reparametrizes the confinement by introducing color screeningprl69 ; prc53 ; npa673 ; prc65p ; prc51 ; npa688 . Figs.3-6(b) show that after the introduction of color screening, the confinement term contributes an additional attraction, which mainly reduces the repulsive core of the effective BB interactions (See Fig.7). There is a systematic uncertainty in our model related to this term which is not yet quantified.

The color magnetic term generally contributes to a repulsive core except for a few Pauli favored channels where it contributes an additional intermediate range attraction so that these channels develop a strong effective attraction. It is also dependent on the strangeness of the channels: The higher the channel strangeness, the weaker the color magnetic contribution, due to the inverse strange quark mass in the color magnetic term. Figs.3-6(d) display these results.

Altogether, these figures show that in the QDCSM the intermediate range attraction is mainly due to reduction of the quark kinetic energy, and the degree of reduction of the quark kinetic energy is connected with the degree of quark delocalization; the latter is determined by the competition between these four terms and that competition is different in different channels. We have mentioned before that in the $`d^{}`$ channel this competition produces a large quark delocalization for a wide range of separations between two $`\mathrm{\Delta }`$-like quark clusters so that the kinetic energy acquires a correspondingly large reduction, which in turn gives rise to a very strong attraction in this channel. (See Fig.3). The $`SIJ=6,0,0`$ has the same $`<P_{36}^{sfc}>`$ as that of the $`SIJ=0,0,3`$ channel and has a similar opportunity to develop large quark delocalizations and strong attraction. However, there the competition does not allow as large a quark delocalization to develop and so the resulting attraction is not as strong as in the $`d^{}`$ channel. (See Fig.4.)

## VI Conclusion

To sum up, we have carried out a dynamical calculation for the most promising dibaryon candidates with high strangeness by using the extended QDCSM. Only a few high strangeness states recommend themselves for experimental searches. These are the $`H`$ particle, the N$`\mathrm{\Omega }`$ and the di-$`\mathrm{\Omega }`$. We recommend searching for the $`(SIJ=3,1/2,2)`$ using $`N\mathrm{\Omega }\mathrm{\Lambda }\mathrm{\Xi }`$ D-wave decay vertex mass reconstruction with the RHIC detectors as the case with the best likelihood.

In the QDCSM, it is the quark dynamics that controls the competition among the four terms in the Hamiltonian (1): the kinetic energy, the confinement, the color Coulomb, and the color magnetic, and which determines the overall effective BB interactions. These are quite similar, though not identical, to those of the quark-gluon-meson hybrid modelprc65  for the majority of BB channels. The exchange of explicit $`K`$ and $`\eta `$ mesons is not important in our model for mass estimates of strange dibaryons because of our short range cut-off.

This model, which has the fewest parameters, describes the properties of the deuteron and the existing NN, N$`\mathrm{\Lambda }`$ and N$`\mathrm{\Sigma }`$ scattering data. Up to now it is the only model which gives an explanation of the long standing fact that the nuclear and molecular forces are similar in character despite the obvious length and energy scale differences and that nuclei are well described as a collection of A nucleons rather than 3A quarks. In view of the fact that the $`H`$-particle has not been observed experimentally, the BB interaction in the $`\mathrm{\Lambda }\mathrm{\Lambda }`$ channeltak  predicted by this model may be a good approximation of the real world. Further refinement is possible by including more channel couplings and spin-orbit and tensor interactions.

Of course, the QDCSM is only a model of QCD. The high strangeness dibaryon resonances may be a good venue for determining whether the QDCSM mechanism for the intermediate range attraction is more realistic than that of a universal $`u`$, $`d`$, $`s`$ quark-$`\sigma `$ meson coupling.

This work is supported by NSFC contract 90103018 and by the U.S. Department of Energy under contract W-7405-ENG-36. F. Wang would like to thank the ITP for their support through the visiting program.

FIGURE CAPTIONS

Fig.1 Relative motion wave functions of the coupled channels for the SIJ=-2,0,0 state with eigenenergy 2225.5 MeV.

Fig.2 Relative motion wavefunctions of the coupled channels for the SIJ=-3,1/2,2 state with eigenenergy 2549.1 MeV.

Fig.3 Contributions of kinetic energy, confinement, color Coulomb and color magnetic terms (a-d respectively) to the effective potential and the total (e) for the SIJ=0,0,3 $`\mathrm{\Delta }\mathrm{\Delta }`$ channel with $`P_{36}^{sfc}=1/9`$. In each subfigure, the dotted curve is for delocalization parameter $`ϵ`$=0.0, dashed curve for $`ϵ`$=0.5, dashed-dotted curve for $`ϵ=1.0`$; the solid curve is for the self-consistent value of $`ϵ`$ determined by the system dynamics.

Fig.4 The same as Fig.3 for the state $`SIJ=6,0,0`$ $`\mathrm{\Omega }\mathrm{\Omega }`$ with $`<P_{36}^{sfc}>=1/9`$.

Fig.5 The same as Fig.3 for the state $`SIJ=0,0,1`$ NN with $`<P_{36}^{sfc}>=1/81`$.

Fig.6 The same as Fig.3 for the state $`SIJ=2,0,0`$ $`\mathrm{\Lambda }\mathrm{\Lambda }`$ with $`<P_{36}^{sfc}>=0`$.

Fig.7 The contribution of conventional quadratic color confinement to the effective BB interaction. The solid curve is for $`ϵ`$=0.0, dotted curve for $`ϵ`$=0.5 and the dashed curve for $`ϵ`$=1.0.


format: text

text:
# Extinction by grey dust in the intergalactic medium

## 1 Introduction

During the last years two groups have investigated the behavior of distant supernovae with the result that their luminosity is lower than expected by standard cosmological models. In a number of papers (Perlmutter et al.), (Schmidt et al.; Riess et al.) both groups ascribe this luminosity lag to the action of a cosmological constant, that means to the existence of a medium with repelling gravity, which leads to an accelerated expansion of the universe. The authors argue that neither evolution effects nor dimming by intervening matter can be responsible for the weakening of the energy flux.

The main argument against extinction by cosmic dust, stressed by the two groups, is that dust would lead to an additional reddening of the spectra, which is not observed. In the reddening argument it is always assumed that the intervening dust grains have approximately the same size distribution as in our galaxy, thus producing the same strong wavelength dependence of extinction.

There has been some criticism on this argument by Aguirre \[1999a\], who describes a method, by which dust consisting of small carbon needles can be formed, showing less reddening than dust of galactic type. But as Riess et al. have shown in a recent paper, even the reddening by the carbon needles proposed by Aguirre (size $`0.1\mu \mathrm{m}`$) would be too strong to be consistent with observations. In a more recent paper Aguirre \[1999b\] argues that dust particles formed in galaxies can be ejected by radiation pressure. As this process favors the ejection of larger particles and disfavors their corrosion, the intergalactic dust produced in this way is less reddening than its galactic source. He shows that cutting the galactic grain size distribution at $`>0.1\mu \mathrm{m}`$ results in reddening compatible with observations.

For dust to be really ”grey” in the spectral range observed in the distant supernova search projects, the particle size should be in the order of $`1\mu \mathrm{m}`$ . To explain the observed extinction of 30% at z=0.46 (Riess et al.) this would require that a mass corresponding to about 0.2% of the critical density is contained in such dust grains. This can easily be estimated, assuming that the dust consists of spherical particles of similar composition as those in our galaxy (silicates) with a specific mass around $`2\mathrm{g}/\mathrm{cm}^3`$. If the density of ordinary matter in the universe is close to 0.2 of the critical density, this value corresponds to a ”metallicity” close to that observed in the solar system.

Though the existence of such grey dust cannot be ruled out by any observed effect, it has not yet been verified by other observations as well. By some straight forward considerations we will propose a mechanism in this paper, by which intergalactic dust can be formed, and which explains the characteristic difference in particle size compared to galactic dust. One even had to look for some new physics, if the extinction by grey dust should not exist. In contrast to the paper of Aguirre \[1999b\] in the mechanism proposed here the dust grains directly condense in an environment flowing at velocities exceeding the escape velocity at the site of formation.

## 2 Production of cosmic dust

It is generally believed that the chemical elements, from which cosmic dust is composed, are produced predominantly by explosions of stars, which have exhausted their nuclear fuel of hydrogen and helium. To understand the process of dust formation, we must look somewhat closer at the processes occurring during the collapse and subsequent explosion, observed as a supernova. While in type Ia supernovae (SNIa) it is assumed that a white dwarf exceeds its stability limit by accretion of mass and is completely destroyed in the explosion, other types of supernovae occur, when high mass stars collapse into a neutron star, expelling their outer shell in an explosive event. Common to both explosions is, that matter in the order of one solar mass consisting of elements heavier than He is ejected with velocities up to 10000 km/sec.

As has been verified by several observations, SNIa form a group with well established properties. If our understanding of SNIa is correct by principle, the mass involved is close to the Chandrasekar limit and is completely ejected in the explosion. The absence of hydrogen lines in the spectrum indicates that the first steps of nuclear fusion are completed. Also a large fraction of C, N and O should have been transformed into heavier elements. To simplify the model of dust formation we restrict our considerations to the processes in SNIa, to show that these processes can explain, why grey intergalactic dust should exist in a quantity as mentioned above.

We will consider the effects following the explosion of a ”template” supernova as described by Perlmutter , that means an SNIa with a peak magnitude $`M=19.5`$ and a decay time $`dM/dt`$ of 0.06m per day. The expanding gas cloud is assumed to have a mass close to the Chandrasekar limit ( $`M_\mathrm{C}=2.8\times 10^{33}\mathrm{g}`$), the element abundances being those of the solar system, but without the presence of H and He. The radial expansion velocity is assumed constant in the order of $`u`$ = 5000 km/sec. Only in the very first moments of the explosion energy is produced by fusion reactions and lost by emission of neutrinos, while during the period of decay radiation from the outer layers of the expanding cloud is the only loss mechanism. When the radius of the cloud is known, as long as the plasma is optically thick, the surface temperature can immediately be estimated from the measured luminosity by the relation

$$L=4\pi R^2\sigma T^4=10^{\frac{M_{\mathrm{peak}}M}{2.5}}\times L_{\mathrm{peak}}$$
(1)
($`M`$ is the absolute bolometric magnitude; $`\sigma `$ is the Stephan-Boltzmann constant).

Using the relation $`R=ut`$, for the data given above recombination into neutral atoms starts (at T=10000K) in the surface layer at about 20 days after maximum brightness and is nearly complete after 50 days (4000K). Though at this time the visible magnitude begins to deviate from the linear slope, as the outer layers of the plasma become optically thin, we can assume that the linear decrease of the bolometric magnitude continues for a longer period of time, as the radiation loss is shifted more and more to IR.

Using this extrapolation we can calculate the time dependence of the thermodynamic equilibrium composition of a gas of solar element abundances (without H and He; C, N and O abundances reduced to 50%). Assuming constant pressure in the expanding gas sphere, the element concentrations are determined by the total mass and the actual sphere radius. The results show that molecules begin to form after 75 days (SiO and CO at 2500K). 105 days from intensity maximum the temperature has dropped to 1500K. At this temperature the first condensed phases ($`\mathrm{MgSiO}_3`$ and $`\mathrm{SiO}_2`$) would occur under LTE conditions. ($`\mathrm{Fe}_2\mathrm{SiO}_4`$ would solidify at 1100K).

At that time the radius of the cloud is about $`4.5\times 10^{15}\mathrm{cm}`$, the particle density in the order of $`10^9\mathrm{cm}^3`$. Under these conditions homogeneous nucleation cannot be expected in the hot gas, but we must keep in mind that the cloud does not expand into vacuum, but into the interstellar medium (ISM), which already contains dust grains, which are swept up by the expanding gas and can serve as a nucleation seed.

Growth of the dust grains will preferably take place in the time period immediately following the onset of supersaturation, as the number of molecular collisions decreases with the third power of the cloud radius and in addition with the square root of the gas temperature. To estimate the final size, which a dust grain of initial radius $`r_0`$ can reach, we have to calculate the total number of collisions between the grain and the surrounding molecules. For simplicity we assume that there is only one type of molecules with molecular mass $`m`$ and that the particle density $`n`$ is constant in the gas cloud: $`n=M_C/(4/3\pi R^3)/m`$. Assuming hard core collision cross section, the mass of a spherical particle of density $`\rho _\mathrm{D}`$ changes at a rate

$$\frac{d}{dt}\left(\frac{4}{3}\pi r^3\rho _\mathrm{D}\right)=\pi r^2nmv_{\mathrm{th}}.$$
(2)
$`v_{\mathrm{th}}=\sqrt{3kT/m}`$ is the mean thermal velocity (the sticking coefficient has been set to 1). The time dependence of the temperature can be approximated by a power law $`T/T_0=(r/r_0)^\alpha `$ with $`\alpha `$ being close to one for the set of parameters used in the example. Using the relation $`R=ut`$, eq.2 can be integrated, starting from the time of beginning supersaturation to infinity with the result

$$r_{\mathrm{}}=r_0+\frac{3}{40\pi }\frac{M_\mathrm{C}v_{\mathrm{th}_0}}{\rho _\mathrm{D}uR_0^2},$$
(3)
The index 0 refers to the values at the time of beginning nucleation. Inserting the data of the example leads to a final radius of about $`2.8\mu \mathrm{m}`$. This value should be regarded only as an order of magnitude estimation, as the derivation contains very crude simplifications. It is assumed that all the material contained in the gas cloud is transformed into species, which can condense into solid particles near 1500K. Solid species occurring at 800K would reach only half the size, as the cloud has expanded further in the mean time. The second assumption, which has been made, is that the density of molecules in the cloud is not reduced by the condensation process itself.

Whether this condition is satisfied, depends on the number of particles, which form in the cloud. The order of magnitude should be equal to the number of seed grains swept up from the ISM. The density of the ISM and its dust content may vary strongly. Taking the orders of magnitude near the solar system as a guide line, the density is in the order of $`10^{23}\mathrm{g}/\mathrm{cm}^3`$ with a dust content in the order of 1%. Assuming a mean particle size of $`0.01\mu \mathrm{m}`$ this leads to a grain number of $`10^7\mathrm{cm}^3`$. Most of the grain growth in the expanding cloud will be completed, before the cloud has reached a volume of $`10^{48}\mathrm{cm}^3`$. This means that the number of grains would be in the order of $`10^{41}`$. If all the material contained in the SN would condense onto these grains, they would grow to a mass of about $`3\times 10^8\mathrm{g}`$ corresponding to a radius of $`15\mu \mathrm{m}`$. Thus we can conclude, that it is the collision rate and not the availability of material, which determines the final size of the dust particles. Grains, which start growing at the beginning of supersaturation can easily reach a size of more than $`1\mu \mathrm{m}`$. Seed grains swept up at later times will reach correspondingly smaller size.

## 3 Size selection

At the time of formation the radial velocity of the grains equals the expansion velocity of the gas cloud. If they were ejected into empty space, all grains would finally leave their parent galaxy, unless they are emitted in the direction towards the center of the galaxy. For an SNIa exploding at solar distance (7.7 kpc) from the center in an arm of a spiral galaxy of say $`10^{45}\mathrm{g}`$, the escape velocity is about 500 km/sec, much less then the velocity of the dust grains.

But in real galaxies the grains have first to traverse the ISM, before they reach free space. In every collision with the gas atoms they lose momentum. Even if the direction of emission is perpendicular to the disk plain, the column density of the gas will be in the order of $`n_\mathrm{c}=10^{21}\mathrm{cm}^2`$.

In every collision with a hydrogen atom the momentum loss of a grain is $`m_\mathrm{H}u`$, as the velocity of the gas atoms is small compared to the grain velocity. The momentum loss rate of a spherical grain of radius $`r`$ and density $`\rho _\mathrm{D}`$ thus is

$$\frac{4\pi }{3}\rho _\mathrm{D}r^3\frac{du}{ds}=nm_\mathrm{H}\pi r^2u,$$
(4)
or integrating along the path of the grain:

$$u_{\mathrm{}}=u_0\mathrm{exp}[3\mathrm{m}_\mathrm{H}\mathrm{n}_\mathrm{c}/(4\rho _\mathrm{D}\mathrm{r})].$$
(5)
If, as estimated above, the initial velocity is ten times the escape velocity, for a grain to leave the galaxy the exponential factor has to be $`>0.1`$. For a column density of $`10^{21}\mathrm{cm}^2`$ this means $`r>2.5\mu \mathrm{m}`$. Though this value can only be regarded as an order of magnitude estimation, it clearly shows that the ISM acts as an effective filter, which allows only the largest grains to escape, while the smaller particles are thermalized to form the seed for condensation in subsequent generations of SNe. The limit is just in the range, which has to be expected, when the dust grains escaping from the galaxy show ”grey” extinction, while the remaining dust is highly selective.

Aguirre and Haiman  have argued that too large dust grains will absorb too much of the visible and UV radiation from stars, which is then re-emitted in the far IR. This emission might be in conflict with observations. For graphite particles they estimate an upper size limit near $`1\mu \mathrm{m}`$. Due to the lower absorptivity of silicate grains, which are expected to form in supernova remnants, such grains will heat up much less, so that even with grains up to $`r=2.5\mu \mathrm{m}`$ the limit set by FIR emission is not exceeded.

## 4 Global distribution of dust

The low scatter in the luminosity data of SNIa requires that any dust causing the dimming must be uniformly distributed in space. Grains, which are formed in gas clouds moving with initial velocities in the order of 10000 km/sec, will reach free space retaining a considerable fraction of this velocity. Thus within the age of the universe they can traverse even the largest voids in space, so that a nearly uniform extinction is expected. This contrasts to the mechanism of acceleration by radiation pressure, proposed by Aguirre\[1999b\], where the final velocities are ten times lower.

By now it has been shown that the formation of large dust grains in SN remnants is possible and that their escape from the parent galaxy is highly selective with respect to grain size, but we have no good estimates, how effective this mechanism of dust formation is, how many stars undergo an SN explosion, and still more, what was the SN rate, when the galaxies were young.

Additional information can be obtained, however, from the x-ray spectra of galaxy clusters. In galaxy clusters dust grains, which leave their parent galaxies, are still confined in the gravitation potential of the cluster. A rich cluster may contain more than 1000 galaxies. But the mean distance of the galaxies from the cluster center is only 100 times larger than the distance of the SN, in which the dust is produced, from the center of the individual galaxy. Thus ten times higher momentum is required to leave the entire cluster, with the consequence that most of the dust is captured in the intracluster medium (ICM).

As is well known from x-ray observations (see e.g. Sarazin \[1988b\]), hot gas with particle densities of $`10^2`$ to $`10^3\mathrm{cm}^3`$ and temperatures between $`2\times 10^7`$ and $`10^8\mathrm{K}`$ fills the space between galaxies in clusters. By collisions with the high energy gas atoms the dust grains will be eroded and finally resolved in the gas. An order of magnitude estimation of the lifetime of grains can easily be obtained (similar to the growth rate in eq.2) from the collision rate and the efficiency of the erosion process. Detailed calculations of the sputter yield of various atoms colliding with silicate grains have been carried out by P.W.May et al. . Typical yields in the energy range considered here are in the order of $`\eta =10^2`$. For a spherical grain of radius $`r`$ the rate of mass loss is

$$\frac{d}{dt}\left(\frac{4}{3}\pi r^3\rho _\mathrm{D}\right)=\eta m_\mathrm{a}n\pi r^2v_{\mathrm{th}},$$
(6)
where $`m_\mathrm{a}`$ is the mass of the eroded atom and $`v_{\mathrm{th}}`$ is the thermal velocity in the gas, which is of the order $`10^8\mathrm{cm}/\mathrm{sec}`$. From eq.6 we find by integration $`\tau =r_0\rho _\mathrm{D}/(\eta m_\mathrm{a}nv_{\mathrm{th}})`$ for the lifetime of a grain with initial radius $`r_0`$. Typical data as discussed above result in lifetimes in the order of $`10^{15}`$ sec for $`1\mu \mathrm{m}`$ grains. That means that a large fraction of the grains is completely resolved, so that from the metallicity of the cluster gas we can immediately derive the amount of dust, which has been ejected during the existence of the galaxy cluster.

X-ray spectra of most clusters show metallicities around half that of the solar system (Sarazin \[1988a\]). Thus, if emission of dust from galaxies is the dominant mechanism for enrichment of heavy elements in the intergalactic medium, we can conclude, that the overall abundance of these elements in free space today should be of the same order as in the ICM. But in free space, where no erosion takes place, the heavy elements remain bound in large dust grains.

Of course there is the possibility that also outside of galaxies and clusters space is not empty, but filled with hot gas, which might for instance be ejected by the jets of quasars. But even if such an intergalactic plasma exists with a density close to the critical closing density of the universe, it would be 1000 times more dilute than the intracluster gas, so that the erosion time of dust grains would be in the order of $`10^{18}`$ sec, comparable to the age of the universe. Thus most grains, that have escaped into free space since the formation of galaxies, will still be there, causing dimming of the light from distant galaxies just at the level observed in the light curves of SNIa.

The study of distant supernovae has not proved the existence of a by now undetected medium with repulsive gravity. However, for the first time it has been proved that the space between galaxies is not as empty as thought before, but that it contains matter, which can considerably influence the observations of distant regions in space.


text:
# A Method for Solving Distributed Service Allocation Problems 1footnote 11footnote 1This work has been funded by Darpa under contract F30602-99-2-0513.

## 1 Introduction

The problem of dynamically allocating services to a changing set of consumers arises in many applications. For example, in an e-commerce system, the service providers are always trying to determine which service to provide to whom, and at what price ; in an automated manufacturing for mass customization scenario, agents must decide which services will be more popular/profitable ; and in a dynamic sensor allocation problem, a set of sensors in a field must decide which area to cover, if any, while preserving their resources.

While these problems might not seem related, they are instances of a more general service allocation problem in which a finite set of resources must be allocated by a set of autonomous agents so as to maximize some global measure of utility. A general approach to solving these types of problems has been used in many successful systems , such as    . The approach involves three general steps:

1. Assign each resource that needs to be preserved to an agent responsible for managing the resource.
2. Assign each goal of the problem domain to an agent responsible for achieving it. Achieving these goals requires the consumption of resources.
3. Have each agent take actions so as to maximize its own utility, but implement a coordination algorithm that encourages agents to take actions that also maximize the global utility.

In this paper we formalize this general approach by casting the problem as a search in a global *fitness landscape* which is defined as the sum of the agents’ utilities. We show how the choice of a coordination/communication protocol disseminates information, which in turn “smoothes” the global utility landscape. This smooth global utility landscape allows the agents to easily find the global optimum by simply making selfish decisions to maximize their individual utility.

We also present experiments that pinpoint the location of a *phase transition* in the time it takes for the agents to find the optimal allocation. The transition can be seen when the amount of communication allowed among agents is manipulated. It exists because communication allows the agents to align their individual landscapes with the global landscape. At some amount of communication, the alignment between these landscapes is good enough to allow the agents to find the global optimum, but less communication drives the agents into a random behavior from which the system cannot recuperate.

### 1.1 Task Allocation

The service allocation problem we discuss in this paper is a superset of the well known task allocation problem \[10, chapter 5.7\]. A task allocation problem is defined by a set of tasks that must be allocated among a set of agents. Each agent has a cost associated with each subset of tasks, which represents the cost the agent would incur if it had to perform those tasks. Coordination protocols are designed to allow agents to trade tasks so that the globally optimal allocation—the one that minimizes the sum of all the individual agent costs—is reached as soon as possible. It has been shown that this globally optimal allocation can reached if the agents use the contract-net protocol  with OCSM contracts . These OCSM contracts make it possible for the system to transition from any allocation to any other allocation in one step. As such, a simple hill-climbing search is guaranteed to eventually reach the global optimum.

In this paper we consider the service allocation problem, which is a superset of the task allocation because it allows for more than one agent to service a “task”. The service allocation problem we study also has the characteristic that every allocation cannot be reached from every other allocation in one step.

### 1.2 Service Allocation

In a service allocation problem there are a set of services, offered by service agents, and a set of consumers who use those services. A server can provide any one of a number of services and some consumers will benefit from that service without depleting it. A server agent incurs a cost when providing a service and can choose not to provide any service.

For example, a server could be an agent that sets up a website with information about cats. All the consumer agents with interests in cats will benefit from this service, but those with other interests will not benefit. Since each server can provide, at most, one service, the problem is to find the allocation of services that maximizes the sum of all the agents’ utilities, that is, an allocation that maximizes the global utility.

#### 1.2.1 Sensor Allocation

Another instance of the service allocation problem is the sensor allocation problem, which we will use as an example throughout this paper. In the sensor allocation problem we have a number of sensors placed in fixed positions in a two-dimensional space. Each sensor has a limited viewing angle and distance but can point in any one of a number of directions. For example, a sensor might have a viewing angle of 120 degrees, viewing distance of 3 feet, and be able to look in three directions, each one 120 degrees apart from the others. That is, it can “look” in any one of three directions. On each direction it can see everything that is in the 120 degree and 3 feet long view cone. Each time a sensor looks in a particular direction is uses energy.

There are also targets that move around in the field. The goal is for the sensors to detect and track all the targets in the field. However, in order to determine the location of a target, two or more sensors have to look at it at the same time. We also wish to minimize the amount of energy spent by the sensors.

We consider the sensor agents as being able to provide three services, one for each sector, but only one at a time. We consider the target agents as consuming the services of the sensors.

## 2 A Formal Model for Service Allocation

We define a service allocation problem $`SA`$ as a tuple $`SA=\{C,S\}`$ where $`C`$ is the set of consumer agents $`C=\{c_1,c_2,\mathrm{},c_{|C|}\}`$, and $`c_i`$ has only one possible state, $`c_i=0`$. The set of service agents is $`S=\{s_1,s_2,\mathrm{},s_{|S|}\}`$ and the value of $`s_i`$ is the value of that service. For the sensor domain in which a sensor can observe any one of three 120-degree sectors or be turned off, we have $`s_i\{0,1,2,\text{off}\}`$. An allocation is an assignment of states to the services (since the consumers have only one possible state we can ignore them). A particular allocation is denoted by $`a=\{s_1,s_2,\mathrm{},s_{|S|}\}`$, where the $`s_i`$ have some value taken from the domain of service states, and $`aA`$, where $`A`$ is the set of all possible allocations. That is, an allocation tells us the state of all agents (since consumers have only one state they can be omitted).

Each agent also has a utility function. The utility that an agent receives depends on the current allocation $`a`$, where we let $`a(s)`$ be the state of service agent $`s`$ under $`a`$. The agent’s utilities will depend on their state and the state of other agents. For example, in the sensor problem we define the utility of sensor $`s`$ as $`U_s(a)`$, where

$$U_s(a)=\{\begin{array}{cc}0\hfill & \text{if }a(s)=\text{ off}\hfill \\ K_1\hfill & \text{otherwise.}\hfill \end{array}$$
(1)
That is, a sensor receives no utility when it is off and must pay a penalty of $`K_1`$ when it is running.

The targets are the consumers, and each target’s utility is defined as

$$U_c(a)=\{\begin{array}{cc}0\hfill & \text{if }f_c(a)=0\hfill \\ K_2\hfill & \text{if }f_c(a)=1\hfill \\ K_2+n2\hfill & \text{if }f_c(a)=n\hfill \end{array}$$
(2)
where

$$f_c(a)=\text{number of sensors }s\text{ that see }c\text{ given their state }a(c)\text{.}$$
(3)
Finally, given the individual agent utilities, we define the global utility $`GU(a)`$ as the sum of the individual agents’ utilities:

$$GU(a)=\underset{cC}{}U_c(a)+\underset{sS}{}U_s(a).$$
(4)
The service allocation problem is to find the allocation $`a`$ that maximizes $`GU(a)`$. In the sensor problem, there are $`4^{|S|}`$ possible allocations, which would make a simple generate-and-test approach take exponential amounts of time. We wish to find the global optimum much faster than that.

### 2.1 Search Algorithms

Our goal is to design an interaction protocol whereby an allocation $`a`$ that maximizes the global utility $`GU(a)`$ is reached in a small number of steps. In each step of our protocol one of the agents will change its state or send a message to another agent. The messages might contain the state or utilities of other agents. We assume that the agents do not have direct access to the other agents’ states or utility values.

The simplest algorithm we can envision involves having each consumer, at each time, changing the state of a randomly chosen service agent so as to increase the consumer’s own utility. That is, a consumer $`c`$ will change the current allocation $`a`$ into $`a^{}`$ by changing the state of some sensor $`s`$ such that $`U_c(a^{})>U_c(a)`$. If the sensor’s state cannot be changed so as to increase the utility, then the consumer does nothing. In the sensor domain this amounts to a target picking a sensor and changing its state so that the sensor can see the target. We refer to this algorithm as individual hill-climbing.

The individual hill-climbing algorithm is simple to implement and the only communication needed is between the consumer and the chose server. This simple algorithm makes every consumer agent increase its individual utility at each turn. However, the new allocation $`a^{}`$ might result in a lower global utility, since $`a^{}`$ might reduce the utility of several other agents. Therefore, it does not guarantee that an optimal allocation will be eventually reached.

Another approach is for each agent to change state so as to increase the global utility. We call this a global hill-climbing algorithm. In order to implement this algorithm, an agent would need to know how the proposed state change affects the global utility as well as the states of all the other agents. That is, it would need to be able to determine $`GU(a^{})`$ which requires it to know the state of all the agents in $`a^{}`$ as well as the utility functions of every other agent, as per the definition of global utility (4). In order for an agent to know the state of others, it would need to somehow communicate with all other agents. If the system implements a global broadcasting method then we would need for each agent to broadcast its state at each time. If the system uses more specialized communications such as point-to-point, limited broadcasting, etc., then more messages will be needed.

Any protocol that implements the global hill-climbing algorithm will reach a locally optimal allocation in the global utility. This is because it is always true that, for a new allocation $`a^{}`$ and old allocation $`a`$, $`GU(a^{})GU(a)`$. Whether or not this local optimum is also a global optimum will depend on the ruggedness of the global utility landscape. That is, if it consists of one smooth peak then it is likely that any local optimum is the global optimum. On the other hand, if the landscape is very rugged then there are likely many local peaks. Studies in NK landscapes  tell us that smoother landscapes result when an agent’s utility depends on the state of smaller number of other agents.

Global hill-climbing is better than individual hill-climbing since it guarantees that we will find a local optima. However, it requires agents to know each others’ utility function and to constantly communicate their state. Such large amount of communication is often undesirable in multiagent systems. We need a better way to find the global optimum.

One way of correlating the individual landscapes to the global utility landscape is with the use of a bidding protocol in which each consumer agent tells each service the marginal utility the consumer would receive if the service switched its state to so as to maximize the consumer’s utility. The service agent can then choose to provide the service with the highest aggregate demand. Since the service is picking the value that maximizes the utility of everyone involved (all the consumers and the service) without decreasing the utility of anyone else (the other services) this protocol is guaranteed to never decrease the global utility. This bidding protocol is a simplified version of the contract-net  protocol in that it does not require contractors to send requests for bids.

However, in order for a consumer to determine the marginal utility it will receive from one sensor changing state, it still needs to know the state of all the other sensors. This means that a complete implementation of this protocol will still require a lot of communication (namely, the same amount as in global hill-climbing). We can reduce this number of messages by allowing agents to communicate with only a subset of the other agents and making their decisions based on only this subset of information. That is, instead of all services telling each consumer their state, a consumer could receive state information from only a subset of the services and make its decision based on this (assuming that the services chosen are representative of the whole). This strategy shows a lot of promise but its performance can only be evaluated on an instance-by-instance basis. We explore this strategy experimentally in Section 3 using the sensor domain.

#### 2.1.1 Theoretical Time Bounds of Global Hill-Climbing

Since we now know that global hill-climbing will always reach a local optimum, the next questions we must answer are:

1. How many local optima are there?
2. What is the probability that a local optimum is the global optimum?
3. How long does it take, on average, to reach a local optimum?

Let $`a`$ be the current allocation and $`a^{}`$ be a neighboring allocation. We know that $`a`$ is a local optimum if

$$_{a^{}N(a)}GU(a)>GU(a^{})$$
(5)
where

$$N(a)=\{x|x\text{ is a Neighbor of }a\}.$$
(6)
We define a Neighbor allocation as an allocation where one, and only one, agent has a different state.

The probability that some allocation $`a`$ is a local optimum is simply the probability that (5) is true. If the utility of all pairs of neighbors is not correlated, then this probability is

$$\mathrm{𝐏𝐫}[_{a^{}N(a)}GU(a)>GU(a^{})]=\mathrm{𝐏𝐫}[GU(a)>GU(a^{})]^b,$$
(7)
where $`b`$ is the branching factor. In the sensor problem $`b=3|S|`$ where $`S`$ is the set of all sensors. That is, since each sensor can be in any of four states it will have three neighbors from each state. In some systems it is safe to assume that the global utilities of $`a`$’s neighbors are independent. However, most systems show some degree of correlation.

Now we need to calculate the $`\mathrm{𝐏𝐫}[GU(a)>GU(a^{})]`$, that is, the probability that some allocation $`a`$ has a greater global utility that its neighbor $`a^{}`$, for all $`a`$ and $`a^{}`$. This could be calculated via an exhaustive enumeration of all possible allocations. However, often we can find the expected value of this probability.

For example, in the sensor problem each sensor has four possible states. If a sensor changes its state from sector $`x`$ to sector $`y`$ the utility of the target agents covered by $`x`$ will decrease while the utility of those in $`y`$ will increase. If we assume that, on average, the targets are evenly spaced on the field, then the global utilities for both of these are expected to be the same. That is, the expected probability that the global utility of one allocation is bigger than the other is $`1/2`$.

If, on the other hand, a sensor changes state from “off” to a sector, or from a sector to “off,” the global utility is expected to decrease and increase, respectively. However, there are an equal number of opportunities to go from “off” to “on” and vice-versa. Therefore, we can also expect that for these cases the probability that the global utility of one allocation is bigger than the other is $`1/2`$.

Based on these approximations, we can declare that for the sensor problem

$$\mathrm{𝐏𝐫}[_{a^{}N(a)}GU(a)>GU(a^{})]=\frac{1}{2^b}=\lambda .$$
(8)
If we assume an even distribution of local optima, the total number of local optima is simply the product of the total number of allocations times the probability that each one is a local optimum. That is,

$$\text{Total number of local optima}=\lambda |A|$$
(9)
For the sensor problem, $`\lambda =1/2^b`$, $`b=3|S|`$ and $`|A|=b^{|S|}`$, so the expected number of local optima is $`b^{|S|}/2^{3|S|}`$.

$$\mathrm{𝐏𝐫}[\text{a local optimum is also global}]=\frac{1}{\lambda |A|}=\frac{1}{2^b}.$$
(10)
We can find the expected time the algorithm will take to reach a local optimum by determining the maximum number of steps from every allocation to the nearest local optimum. This gives us an upper bound on the number of steps needed to reach the nearest local optimum using global hill-climbing. Notice that, under either individual hill-climbing or the bidding protocol it is possible that the local optimum is not reached, or is reached after more steps, since these algorithms can take steps that lower the global utility.

In order to find the expected number of steps to reach a local optimum, we start at any one of the local optima and then traverse all possible links at each depth $`d`$ until all possible allocations have been visited. This occurs when

$$\lambda |A|b^d>|A|.$$
(11)
Solving for $`d`$, and remembering that $`\lambda =1/2^b`$, we get

$$d>b\mathrm{log}_b2.$$
(12)
The expected worst-case distance from any point to the nearest local optimum is, therefore, $`b\mathrm{log}_b2`$ (this number only makes sense for $`b2`$ since smaller number of neighbors do not form a searchable space). That is, the number of steps to reach the nearest local optima in the sensor domain is proportional to the branching factor $`b`$, which is equal to $`3|S|`$. We can expect search time to increase linearly with the number of sensors in the field.

## 3 Simulations

While the theoretical results above give us some bounds on the number of iterations before the system is expected to converge to a local optimum, the bounds are rather loose and do not tell us much about the dynamics of the executing system. Also, we cannot show mathematically how changes in the amount of communication change the search. Therefore, we have implemented a service allocation simulator to answer these questions. It simulates the sensor allocation domain described in the introduction.

The simulator is written in Java and the source code is available upon request. It gathers and analyzes data from any desired number of runs. The program can analyze the behavior of any number of target and sensor agents on a two-dimensional space, and the agents can be given any desired utility function. The program is limited to static targets. That is, it only considers the one-shot service allocation problem. Each new allocation is completely independent of any previous one.

In the tests we performed, each run has seven sensors and seven targets, all of which are randomly placed on a two-dimensional grid. Each sensor can only point in one of three directions or sectors. These three sectors are the same for all sensors (specifically, the first sector is from 0 to 120 degrees, the second one from 120 to 240, and the third one from 240 to 360). All the sensors use the same utility function which is given by (1), while the targets use (2). After a sensor agent receives all the bids it chooses the sector that has the heighest aggregate demand, as described by the bidding protocol in Section 2.1.

During a run, each of the targets periodically sends a bid to a number of sensors asking them to turn to the sector that faces the target. We set the bid amount to a fixed number for these tests. Periodically, the sensors count the number of bids they have received for each sector and turn their detector (such as a radar) to face the sector with the highest aggregate demand. We assume that neither the targets nor the sensors can form coalitions.

We vary the number of sensors to which the targets send their bids in order to explore the quality of the solution that the system converges upon as the amount of communication changes. For example, at one extreme if the all the targets send their bids to all the sensors, then the sensors would always set their sector to be the one with the most targets. This particular service allocation should, usually, be the best. However, it might not always be the optimal solution. For example, if seven targets are clustered together and the eighth is on another part of the field, it would be better if six sensor agents pointed towards the cluster of targets while the remaining two sensor agents pointed towards the stray target rather than having all sensor agents point towards the cluster of targets. At the other extreme, if all the targets send their bids to only one sensor then they will minimize communications but then the sensors will point to the sector from which they received a message—an allocation which is likely to be suboptimal.

These simulations explore the ruggedness of the system’s global utility landscape and the dynamics of the agents’ exploration of this landscape. If the agents were to always converge on a local (non-global) optimum then we would deduce that this problem domain has a very rugged utility landscape. On the other hand, if they usually manage to reach the global optimum then we could deduce a smooth utility landscape.

## 4 Test Results

In each of our tests we set the number of agents that each target will send its bid to, that is, the number of neighbors, to a fixed number. Given this fixed number of neighbors, we then generated 100 random placements of agents on the field and ran our bidding algorithm 10 times on each of those placements. Finally, we plotted the average solution quality, over the 10 runs, as a function of time for each of the 100 different placements. The solution quality is given by the ratio

$$\alpha =\frac{\text{Current Utility}}{\text{Globally Optimal Utility}},$$
(13)
so if $`\alpha =1`$, then it means that the run has reached the global optimum. Since the number of agents is small, we were able to calculate the global optimum using a brute-force method. Specifically, there are $`3^7=2187`$ possible configurations times $`100`$ random placements leads to $`218700`$ combinations that we had to check for each run in order to find the global optimum using brute-force. Using more than $`7`$ sensors made the test take too long. Notice, however, that our algorithm is much faster than this brute-force search which we perform only to confirm that our search does find the global optimum.

In our tests there were always seven target agents and seven sensor agents. We varied the number of neighbors from 1 to 7. If the target can only communicate with one other sensor, the sensors will likely have very little information for making their decision, while if all targets communicate with all seven sensors, then each sensor will generally be able to point to the sector with the most targets. However, because these decisions are made in an asynchronous manner, it is possible that some sensor will sometimes not receive all the bids before it has to make a decision. The targets always send their bids to the sensors that are closest to them.

The results from our experiments are shown in Figure 1 where we can see that there is a transition in the system’s performance as the number of neighbors goes from three to five. That is, if the targets only send their bids to three sensors then it is almost certain that the system will stay in a configuration that has a very low global utility. However, if the targets send their bids to five sensors, then it is almost guaranteed ($`98\%`$ of the time) that the system will reach the globally optimal allocation. This is a huge difference in terms of the performance. We also notice in Figure 2 that for four neighbors there is a fairly even distribution in the utility of the final allocation.

## 5 Related Work

There is ongoing work in the field of complexity that attempts to study they dynamics of complex adaptive systems . Our approach is based on ideas borrowed from the use of NK landscapes for the analysis of co-evolving systems. As such, we are using some of the results from that field. However, complexity theory is more concerned with explaining the dynamic behavior of existing systems, while we are more concerned with the engineering of multiagent systems for distributed service allocation.

The Collective Intelligence (COIN) framework  shares many of the same goals of our research. They start with a global utility function from which they derive the rewards functions for each agent. The agents are assumed to use some form of reinforcement learning. They show that the global utility is maximized when using their prescribed reward functions. They do not, however, consider how agent communication might affect the individual agent’s utility landscape.

The task allocation problem has been studied in , but the service allocation problem we present in this paper has received very little attention. There is also work being done on the analysis of the dynamics of multiagent systems for other domains such as e-commerce  and automated manufacturing . It is possible that extensions to our approach will shed some light into the dynamics of these domains.

## 6 Conclusions

We have formalized the service allocation problem and examined a general approach to solving problems of this type. The approach involves the use of utility-maximizing agents that represent the resources and the services. A simple form of bidding is used for communication. An analysis of this approach reveals that it implements a form of distributed hill-climbing, where each agent climbs its own utility landscape and not the global utility landscape. However, we showed that increasing the amount of communication among the agents forces each individual agent’s landscape to become increasingly correlated to the global landscape.

These theoretical results were then verified in our implementation of a sensor allocation problem—an instance of a service allocation problem. Furthermore, the simulations allowed us to determine the location of a phase transition in the amount of communication needed for the system to consistently arrive at the globally optimal service allocation.

More generally, we have shown how a service allocation problem can be viewed as a distributed search by multiple agents over multiple landscapes. We also showed how the correlation between the global utility landscape and the individual agent’s utility landscape depends on the amount and type of inter-agent communication. Specifically, we have shown that increased communications leads to a higher correlation between the global and individual utility landscapes, which increases the probability that the global optimum will be reached. Of course, the success of the search still depends on the connectivity of the search space, which will vary from domain to domain. We expect that our general approach can be applied to the design of any multiagent systems whose desired behavior is given by a global utility function but whose agents must act selfishly.

Our future work includes the study of how the system will behave under perturbations. For example, as the target moves it perturbs the current allocation and the global optimum might change. We also hope to characterize the local to global utility function correlation for different service allocation problems and the expected time to find the global optimum under various amounts of communication.