text:
# MADPH-00-1164 hep-ph/0003248 March, 2000 Effective potential calculation of the MSSM lightest CP-even Higgs boson massContribution to PASCOS99: 7th International Symposium on Particles, Strings and Cosmology, Granlibakken, Tahoe City, California, 10-16 Dec 1999.

## Abstract

I summarize results of two-loop effective potential calculations of the lightest CP-even Higgs boson mass in the minimal supersymmetric standard model.

Computing the lightest CP-even Higgs boson mass is the most important loop calculation in the minimal supersymmetric standard model because of the paramount importance of a precise $`m_{h^0}`$ value to the Higgs boson experimental discovery. Tree-level supersymmetry relations require that the Higgs field quartic coupling be related to the electroweak gauge couplings; therefore they impose a strict upper bound $`m_{h^0}m_Z`$, which is already in conflict with the current lower limit from LEP 2.

It is well-known that this tree-level limit can be drastically changed by radiative corrections. One-loop calculations show that incomplete cancellations of the top and stop loops give positive corrections of the form

$$\mathrm{\Delta }m_{h^0}^2=\frac{3h_t^2m_t^2}{4\pi ^2}\mathrm{ln}\frac{m_{\stackrel{~}{t}}^2}{m_t^2},$$
(1)
where $`m_t`$ and $`m_{\stackrel{~}{t}}`$ are top and stop masses respectively. This formula, however, suffers from an ambiguity in the definition of $`m_t`$. Numerically, using running or on-shell top-quark mass can amount to about $`20\%`$ difference in $`\mathrm{\Delta }m_{h^0}^2`$. The problem can only be resolved by an explicit two-loop calculation.

Two-loop calculations in the existing literature have used two different approaches: (a) a renormalization group resummation approach, and (b) a two-loop diagrammatic approach. In the first approach, leading and next-to-leading logarithmic corrections are calculated by integrating one- and two-loop renormalization group equations. However, two-loop non-logarithmic finite corrections are not calculable in principle. The second approach was initiated by Hempfling and Hoang using an effective potential method; they restricted their calculation to specific choice of supersymmetry parameters: i.e. large $`\mathrm{tan}\beta \mathrm{}`$ and zero left-right stop mixing. Two-loop QCD corrections were later computed at more general cases in the effective potential approach. $`m_{h^0}`$ to the same two-loop QCD order was also computed using an explicit diagrammatic method. These calculations incorporate both two-loop logarithmic and non-logarithmic finite corrections. In the following, I shall concentrate on the effective potential approach.

The general way of calculating corrections to CP-even Higgs boson mass is to compute Higgs self-energy and tadpole diagrams to the required loop order. In an effective potential approach, these diagrams can be derived from a generating functional, i.e. the effective potential, by taking explicit derivatives with respect to the Higgs fields. These quantities then enter the MSSM CP-even Higgs boson mass-squared matrix as follows

$$_h^2=\left[\begin{array}{cc}m_Z^2c_\beta ^2+m_{A^0}^2s_\beta ^2+\mathrm{\Delta }_{11}^2& (m_Z^2+m_{A^0}^2)s_\beta c_\beta +\mathrm{\Delta }_{12}^2\\ (m_Z^2+m_{A^0}^2)s_\beta c_\beta +\mathrm{\Delta }_{21}^2& m_Z^2s_\beta +m_{A^0}^2c_\beta +\mathrm{\Delta }_{22}^2\end{array}\right],$$
(2)
where $`\mathrm{\Delta }_{ij}^2`$ represents radiative corrections to the $`ij`$-entry. We note that all these corrections are computed at the zero external momentum limit; sometimes it is necessary to calculate self-energy diagrams directly to account for corrections from non-zero external momenta.

The CP-even Higgs boson masses can be calculated by diagonalizing the above matrix in eq. (2). This computation is tedious but can be greatly simplified when one considers the case $`m_{A^0}m_Z`$, where $`m_{A^0}`$ is the mass of the pseudoscalar $`A^0`$. In this case, we find the corrections to $`m_{h^0}^2`$ is

$$\mathrm{\Delta }m_{h^0}^2=\frac{4m_t^4}{v^2}\left(\frac{d}{dm_t^2}\right)^2V\mathrm{Re}\mathrm{\Pi }_{hh}(m_{h^0}^2)+\mathrm{Re}\mathrm{\Pi }_{hh}(0).$$
(3)
where $`V`$ is the effective potential, $`v`$ the Higgs field VEV, and the last two terms account for non-zero external momentum corrections.

We have carried out this calculation procedure to the two-loop order including leading QCD and top Yukawa corrections. To illustrate our analysis, we present an approximation formula which is derived under the following assumptions: the soft masses for left and right stops, gluino, heavy Higgs bosons and Higgsinos have a common mass $`M_S`$, where $`M_S`$ can be identified as the supersymmetry scale. The two eigenvalues and mixing angle of stops are then accordingly $`m_{\stackrel{~}{t}_1}^2=m_{\stackrel{~}{t}}^2+m_tX_t`$, $`m_{\stackrel{~}{t}_2}^2=m_{\stackrel{~}{t}}^2m_tX_t`$ and $`s_t=c_t=\frac{1}{\sqrt{2}}`$, where the average top-squark mass $`m_{\stackrel{~}{t}}^2=M_S^2+m_t^2`$, and $`X_t=A_t+\mu /\mathrm{tan}\beta `$ is the left-right stop mixing parameter.

We find the approximation formula for two-loop QCD+top Yukawa corrections is (in terms of on-shell mass parameters)

$`\mathrm{\Delta }m_{h^0}^2={\displaystyle \frac{3m_t^4}{2\pi ^2v^2}}\left(\mathrm{ln}{\displaystyle \frac{m_{\stackrel{~}{t}}^2}{m_t^2}}+\widehat{X}_t^2{\displaystyle \frac{\widehat{X}_t^4}{12}}\right)`$ (4)
$`+`$ $`{\displaystyle \frac{\alpha _sm_t^4}{\pi ^3v^2}}\left(3\mathrm{ln}^2{\displaystyle \frac{m_{\stackrel{~}{t}}^2}{m_t^2}}6\mathrm{ln}{\displaystyle \frac{m_{\stackrel{~}{t}}^2}{m_t^2}}+6\widehat{X}_t3\widehat{X}_t^2\mathrm{ln}{\displaystyle \frac{m_{\stackrel{~}{t}}^2}{m_t^2}}{\displaystyle \frac{3\widehat{X}_t^4}{4}}\right)`$
$`+`$ $`{\displaystyle \frac{3\alpha _tm_t^4}{16\pi ^3v^2}}\{s_\beta ^2(3\mathrm{ln}^2{\displaystyle \frac{M_S^2}{m_t^2}}+13\mathrm{ln}{\displaystyle \frac{M_S^2}{m_t^2}})1{\displaystyle \frac{\pi ^2}{3}}+c_\beta ^2(60K+{\displaystyle \frac{13}{2}}+{\displaystyle \frac{4\pi ^2}{3}})`$
$`+`$ $`\left[3s_\beta ^2\mathrm{ln}{\displaystyle \frac{M_S^2}{m_t^2}}c_\beta ^2\left({\displaystyle \frac{69}{2}}+24K\right)+41\right]\widehat{X}_t^2\left(1+{\displaystyle \frac{61}{12}}s_\beta ^2\right)\widehat{X}_t^4+{\displaystyle \frac{s_\beta ^2}{2}}\widehat{X}_t^6`$
$`+`$ $`c_\beta ^2[(316K\sqrt{3}\pi )(4\widehat{X}_t\widehat{Y}_t+\widehat{Y}_t^2)+(16K+{\displaystyle \frac{2\pi }{\sqrt{3}}})\widehat{X}_t^3\widehat{Y}_t`$
$`+`$ $`({\displaystyle \frac{4}{3}}+24K+\sqrt{3}\pi )\widehat{X}_t^2\widehat{Y}_t^2({\displaystyle \frac{7}{12}}+8K+{\displaystyle \frac{\pi }{2\sqrt{3}}})\widehat{X}_t^4\widehat{Y}_t^2]\},`$
where the constant $`K0.195`$. We note that two-loop QCD corrections depend only on $`\widehat{X}_t=X_t/m_{\stackrel{~}{t}}`$ while the top Yukawa corrections depend on $`\widehat{Y}_t=(A_t\mu \mathrm{tan}\beta )/m_{\stackrel{~}{t}}`$ as well. This approximation formula is good to a level of 0.5 GeV for most of the parameter space.

Fig. 1 shows the Higgs boson mass $`m_{h^0}`$ vs. the stop mixing parameter $`\widehat{X}_t`$, for different choices of $`M_S`$, $`\mu `$ and $`\mathrm{tan}\beta `$. The two-loop QCD corrections agree well with other approaches. They generally decrease $`m_{h^0}`$ from their one-loop values by $`1020`$ GeV depending on the parameter choice. Two-loop Yukawa corrections are sizeable for large stop mixings, in particular, for $`\widehat{X}_t\pm 2`$ two-loop Yukawa corrections can increase $`m_{h^0}`$ by about $`5`$ GeV.

Another interesting feature observed in the literature is that two-loop corrections shift the maximal mixing peaks. At the one-loop level, these peaks are at $`\widehat{X}_t=\pm \sqrt{6}`$. It is easy to see from eq. (4) that the size of shifts is about $`10\%`$, i.e. the peaks move to $`\widehat{X}_t\pm 2`$. This is confirmed by Fig. 1.

Finally, renormalization group resummation technique can be used to derive a particularly nice mass correction formula which has clearer physical interpretations. We find eq. (4) can be transformed into the following form by using solutions to the renormalization group equations

$$\mathrm{\Delta }m_{h^0}^2=\frac{3\overline{m}_t^4(Q_t)}{2\pi ^2\overline{v}^2(Q_1^{})}\mathrm{ln}\frac{m_{\stackrel{~}{t}}^2(Q_{\mathrm{th}})}{\overline{m}_t^2(Q_t^{})}+\frac{3\overline{m}_t^4(Q_{\mathrm{th}})}{2\pi ^2\overline{v}^2(Q_2^{})}\left[\widehat{X}_t^2(Q_{\mathrm{th}})\frac{\widehat{X}_t^4(Q_{\mathrm{th}})}{12}\right]+\mathrm{\Delta }_{\mathrm{th}}^{(2)},$$
(5)
where $`Q_1^{}=e^{1/3}m_t`$, $`Q_2^{}=e^{1/3}m_t`$, $`Q_t=\sqrt{m_tm_{\stackrel{~}{t}}}`$, $`Q_t^{}=(m_tm_{\stackrel{~}{t}}^2)^{1/3}`$ and $`Q_{\mathrm{th}}=m_{\stackrel{~}{t}}`$, $`\overline{v}`$ and $`\overline{m}`$ are the Standard Model $`\overline{\mathrm{MS}}`$  parameters. These choices of scales for evaluating one-loop corrections automatically take into account two-loop leading and next-to-leading logarithmic effects. The leftover finite correction term $`\mathrm{\Delta }_{\mathrm{th}}^{(2)}`$ is understood as two-loop threshold corrections and numerically small; its detail form can be found in a forthcoming paper.

I thank J. R. Espinosa for collaborations. This work was supported in part by a DOE grant No. DE-FG02-95ER40896 and in part by the Wisconsin Alumni Research Foundation.


format: text

text:
# Electronic structure of self-assembled quantum dots: comparison between density functional theory and diffusion quantum Monte Carlo

## I Introduction

The electronic structure of quantum dots resembles that of atoms. Both systems display three dimensional electronic confinement leading to level degeneracies, shell structure, and spin correlation to only mention the most studied atomic properties of quantum dots.  Direct as well as exchange Coulomb interactions are also present in quantum dot systems, the importance of each depending on the dot dimensions.  That is why some of the theoretical methods used to calculate the electronic structure of quantum dots are the same as the ones used to study atoms and molecules. One of these methods, local spin density approximation (LSDA) density functional theory (DFT) within the effective mass approximation (EMA), has been widely used to calculate the electron-electron interaction in dots  because of its simple implementation and low computer demand. LSDA is an approximate theory, well known to predict incorrectly the physical properties of some molecules and solids, while performing well on many other systems. However, experience with LSDA on atoms and molecules may not carry over to quantum dots since the confining potentials and electron densities can be very different. As a general rule, density functional theory tends to work better for high effective density, and fails for low density systems, where correlation effects become important. There have been several investigations of LSDA and exact treatment of electron interactions in parabolic quantum dots, with particular emphasis on low density systems and effects of external magnetic fields. The electron interactions in parabolic dots have been treated exactly using diffusion Monte Carlo (DMC) as well as exact diagonalization methods. The parabolic potentials are popular models because they are computationally convenient and, with two adjustable parameters, have had relative success in explaining the influence of many-body effects on the electronic properties of dots.

In this paper we compare EMA results obtained with LSDA in self-assembled InAs/GaAs quantum dots against an exact treatment of the EMA many-body interactions obtained with diffusion quantum Monte Carlo.  These small dots, which can hold about six electrons, have complicated confining potentials that differ considerably from the parabolic dots considered in previous DMC and exact diagonalization calculations. The single particle states of realistically modeled pyramidal dots are significantly modified from the states in a parabolic dot,  so many-body interactions should also differ from those in parabolic models. The purpose of this comparison is to quantify errors due to LSDA for calculations on realistically modeled InAs/GaAs dots. Although error of LSDA is expected to be small in these systems, it is large enough to be a limitation to the comparison of model results to experiments. For the system reported here, we consider LSDA errors of up to 10 meV to be acceptable given the current precision of the models and experiments. Determination of whether LSDA meets this acceptability criterion requires careful comparison with many-body calculations on realistically modeled potentials.

In these calculations we consider the ground state energy as a function of the number of electrons in the dot, which is the quantity for which DFT should apply. Electron addition and removal energies are experimentally measurable. These may be rigorously defined as differences in the total energies of the ground states of the dot which differ by one electron. A common approach to calculating electron addition and removal energies is to use the eigenvalues of the DFT equations, which are well known not to be interpretable as electron addition and removal energies. We have addressed this issue in a previous paper  where we have shown that addition and removal energies are to a very good approximation given by the “Slater transition rule”,  which uses the eigenvalues at one-half occupation. Thus the present work, which compares total energies, is a rigorous test of the LSDA, and combined with our previous work  is a test of the accuracy of the addition and removal energies obtained using LSDA eigenvalues and the Slater transition rules.

The ground state energy is the lowest eigenvalue of the many-body Schrödinger equation, which we take to be

$$\left\{\underset{i=1}{\overset{N}{}}\frac{\mathrm{}^2}{2}_i\left[M^1_i\right]+V_{\mathrm{cb}}(R)+V_{\mathrm{ee}}(R)\right\}\mathrm{\Psi }(R)=E\mathrm{\Psi }(R),$$
(1)
where $`\mathrm{\Psi }(R)=\mathrm{\Psi }(𝐫_1,𝐫_2,\mathrm{}𝐫_N)`$ is the $`N`$-electron wave function, $`M`$ is the electron effective mass tensor and the potential energy terms $`V_{\mathrm{cb}}(R)`$ and $`V_{\mathrm{ee}}(R)`$ are given by

$`V_{\mathrm{cb}}(R)`$ $`=`$ $`{\displaystyle \underset{i=1}{\overset{N}{}}}v_{\mathrm{cb}}(𝐫_i)`$ (2)
$`V_{\mathrm{ee}}(R)`$ $`=`$ $`{\displaystyle \underset{i<j}{}}{\displaystyle \frac{e^2}{ϵ|𝐫_i𝐫_j|}}`$ (3)
where $`v_{\mathrm{cb}}(𝐫)`$ is the potential for a single electron due to the offset and strain potential of the conduction band,  and $`V_{\mathrm{ee}}(R)`$ is the Coulomb interaction of the conduction electrons with charge $`e`$ in a dielectric characterized by the dielectric constant $`ϵ`$. The offset of the conduction band edge causes a step potential at the edge of the dot, which is reduced by strain, resulting in a potential $`v_{\mathrm{cb}}`$ that is zero outside the dot and a roughly constant depth of several tenths of an eV in the interior of the dot, which typically has dimensions $`10\mathrm{nm}`$. Strain causes the electron effective mass to become anisotropic leading to a mass tensor given by $`\mathrm{diag}(M)=(m_{\mathrm{xx}}m_{\mathrm{yy}}m_{\mathrm{zz}})`$ and zero off-diagonal terms.  In the usual case of sample growth along the crystal direction (0 0 1), the electron masses along the plane perpendicular to the growth direction are equal, i.e. $`m_{\mathrm{xx}}=m_{\mathrm{yy}}`$.

## II Density Function Theory and the Local Spin Density Approximation

The DFT approach to obtaining the ground state energy is to replace the rather complex N-electron ground state wave function and the associated Schrödinger equation by the much simpler ground state electron density $`\rho (𝐫)`$ and the corresponding functional forms $`T[\rho ]`$ and $`V[\rho ]`$ of the kinetic and potential energy operators, respectively.  However, those functional forms are unknown, and approximations are necessary. In Kohn-Sham theory the density is given by a sum of densities of single particle orbitals,

$$\rho (𝐫)=\underset{\sigma }{}\underset{i=1}{\overset{N}{}}n_i|\psi _i(𝐫,\sigma )|^2.$$
(4)
The total energy can then written as

$$E=E_{\mathrm{kin}}+E_{\mathrm{cb}}[\rho ]+E_\mathrm{H}[\rho ]+E_\mathrm{x}[\rho ]+E_\mathrm{c}[\rho ],$$
(5)
where

$$E_{\mathrm{kin}}=\frac{\mathrm{}^2}{2}\underset{\sigma ,i=1}{\overset{N}{}}\psi _i^{}(𝐫)(M^1)\psi _i(𝐫,\sigma )𝑑𝐫$$
(6)
is the kinetic energy of the Kohn-Sham orbitals,

$$E_{\mathrm{cb}}[\rho ]=\rho (𝐫)v_{\mathrm{cb}}(𝐫)𝑑𝐫,$$
(7)
is the potential energy from the conduction band offset and strain,

$$E_\mathrm{H}[\rho ]=e\rho (𝐫)\varphi (𝐫)𝑑𝐫,$$
(8)
is the Hartree energy, and $`E_\mathrm{x}[\rho ]`$ and $`E_\mathrm{c}[\rho ]`$ are functionals of the electron density. The self-consistent electrostatic potential $`\varphi (𝐫)`$ is obtained as a solution to the Poisson equation

$$^2\varphi (𝐫)=\frac{e}{ϵ}\rho (𝐫)$$
(9)
with the boundary condition $`\varphi (𝐫)0`$ as $`|𝐫|\mathrm{}`$. Minimization of the total energy defined by Eqs. (4 \- 9) yields a one particle Schrödinger equation which may be solved self-consistently with the expression for the density and the Poisson equation to determine the exact many-body ground state energy. However, the nonlocal functionals $`E_\mathrm{x}[\rho ]`$ and $`E_\mathrm{c}[\rho ]`$ are in general not known, so an approximation must be made in order to proceed with the Kohn-Sham approach.

In the local spin density approximation (LSDA) the exchange and correlation functionals are taken to be functions of the local charge and spin density, which are then uniquely defined by the requirement that the approximation be exact for the homogeneous electron gas. These functions have been fit to the exchange and correlation of the uniform electron gas. The validity of LSDA has been a subject of much research in atomic, crystalline, molecular, and parabolic quantum dot systems,  but has not been checked for self-assembled dot potentials, which are qualitatively like finite square wells with rather high electron density. The density gradients and shell structures of self-assembled dots cause the exchange and correlation energy to differ from that of a uniform gas. The purpose of this paper is to quantify this difference as a function of the number of electrons in a realistically modeled dot.

## III Quantum Monte Carlo Methods

The simplest QMC method is variational Monte Carlo (VMC), in which expectation values for a trial wave function are evaluated using a Metropolis algorithm. Using the variational principle of the ground state energy, parameters in the wave function can be optimized to minimize the total energy, thus providing an approximation to the ground state energy and wave function. In our calculation we use VMC to optimize the guiding wave function for our diffusion Monte Carlo algorithm, described below.

We also use Monte Carlo integration to evaluate the exact exchange energy. In the Kohn-Sham formalism, the exchange energy is defined as the difference

$$E_\mathrm{x}[n]=\mathrm{\Psi }_{\mathrm{KS}}|V_{\mathrm{ee}}|\mathrm{\Psi }_{\mathrm{KS}}E_\mathrm{H}[\rho ],$$
(10)
where $`|\mathrm{\Psi }_{\mathrm{KS}}`$ is a Slater determinant of the true Kohn-Sham orbitals; those that give the exact density. We assume that LSDA is accurate enough to provide approximations to the Kohn-Sham orbitals, i. e. that the LSDA density is close to the exact density. We then evaluate the integral $`\mathrm{\Psi }_{\mathrm{KS}}|V_{\mathrm{ee}}|\mathrm{\Psi }_{\mathrm{KS}}`$ using a Slater determinant of our calculated LSDA wave functions.

Diffusion Monte Carlo is a stochastic method that is able to project the N-electron ground state wave function $`\mathrm{\Phi }_0`$ from a trial wave function $`\mathrm{\Psi }_\mathrm{T}`$. A position basis is used to described the state of the system. A good trial wave function is important because the variance of the Monte Carlo sampling decreases as the trial wave function approaches the true ground state wave function. We chose the trial wave function to be a product of Slater determinants for spin up and spin down electrons with a Jastrow factor that introduces correlation,

$$\mathrm{\Psi }_\mathrm{T}(R)=det[\theta _\mathrm{k}(𝐫_\mathrm{i})]\times \mathrm{exp}[\underset{\mathrm{i}<\mathrm{j}}{}u(r_{\mathrm{ij}})].$$
(11)
For convenience, we take the single particle wave functions $`\theta _\mathrm{k}(𝐫)`$ to be the non-interacting eigenstates for the external potential, and choose the Jastrow factor $`u(r)=ar/(1+br)`$ with $`a=m/ϵ`$ and $`b=1/2`$.

Following the method described in Ref. \[\], we used importance sampling from a function $`\mathrm{\Phi }(R)=\mathrm{\Psi }_\mathrm{T}(R)\mathrm{\Phi }_0(R)`$. The anisotropic mass is taken into account using the mass tensor, $`M`$. The projection of a state $`R`$ to the ground state distribution $`\mathrm{\Psi }_\mathrm{T}\mathrm{\Phi }_0`$ is accomplished by repeatedly applying the projection operator $`e^{H\tau }`$, each time sampling a new configuration. Here $`H`$ is the Hamiltonian and $`\tau `$ is a parameter chosen to be small so that the projection operator is approximated by

$`R|`$ $`\mathrm{\Psi }_\mathrm{T}(R)e^{H\tau }\mathrm{\Psi }_\mathrm{T}^1(R^{})|R^{}`$ (14)
$`={\displaystyle \frac{\mathrm{exp}\{[E_\mathrm{L}(R)E_0]\tau \}}{[(2\pi \tau )^3detM^1]^{\frac{n}{2}}}}\times `$
$`{\displaystyle \underset{i=1}{\overset{n}{}}}\mathrm{exp}\left[{\displaystyle \frac{(𝐫_\mathrm{i}𝐫_\mathrm{i}\tau 𝐯_{\mathrm{d},\mathrm{i}})M(𝐫_\mathrm{i}𝐫_\mathrm{i}\tau 𝐯_{\mathrm{d},\mathrm{i}})}{2\tau }}\right],`$
where $`E_0`$ is a constant parameter to control the population, $`E_\mathrm{L}(R)=\mathrm{\Psi }_\mathrm{T}(R)^1H\mathrm{\Psi }_\mathrm{T}(R)`$ is the local energy, and $`𝐯_\mathrm{d}=M^1\mathrm{log}\mathrm{\Psi }_\mathrm{T}`$. The algorithm is thus: (1) Start with an ensemble of configurations (walkers) distributed by $`|\mathrm{\Psi }_\mathrm{T}|^2`$, (2) propagate each configuration with a drift $`𝐯_\mathrm{d}`$ and Gaussian displacement with covariance matrix $`\sigma _{\mathrm{ij}}^2=M_{\mathrm{ij}}^1\tau `$, (3) reweight each configuration by a factor $`\mathrm{exp}\{(E_\mathrm{L}E_0)\tau \}`$, (4) repeat, collecting statistics once the steady state distribution is reached. The exact ground state energy is obtained by sampling the local energy. We use branching at each step to improve the efficiency of the process.

A complication known as the fermion sign problem arises when applying diffusion Monte Carlo to fermions. Electron exchange introduces negative signs into the projection operator, which decreases the efficiency of the Monte Carlo sampling. In the present discussion the short time approximation, Eq. (14), breaks down when walkers cross nodes of $`\mathrm{\Psi }_\mathrm{T}`$. We handle this problem with the fixed node approximation,  in which walkers are given a weight of zero when they cross a node of $`\mathrm{\Psi }_\mathrm{T}`$. This approximation has a lower bound property, so the true ground state energy can only be lower than the fixed node energy, and the exact choice of $`\mathrm{\Psi }_\mathrm{T}=\mathrm{\Phi }_0`$ would give the exact Fermionic ground state energy.

The correlation energy can be deduced from a knowledge of all other energies,

$$E_c[\rho ]=E_{\mathrm{tot}}E_{\mathrm{kin}}E_{\mathrm{cb}}E_\mathrm{H}E_\mathrm{x},$$
(15)
or equivalently

$$E_c[\rho ]=\mathrm{\Phi }_0|V_{\mathrm{ee}}|\mathrm{\Psi }_0\mathrm{\Psi }_{\mathrm{KS}}|V_{\mathrm{ee}}|\mathrm{\Psi }_{\mathrm{KS}}.$$
(16)
## IV Numerical Comparison

Using DFT and DMC we have calculated the total ground state energy of the system of electrons, defined for DFT by Eq. 5, and defined within DMC as the fixed node ground state energy, as defined in the previous section. Figure 1 shows the system we have considered, a pyramidal shaped InAs quantum dot in a GaAs matrix. However, our method of analysis is not restricted to this particular shape, and we expect our results to be valid for truncated pyramidal or lens shaped dots. We have introduced several simplifications in the model. Since our DMC scheme does not allow for a position dependent dielectric constant, we have used only one value, $`\overline{ϵ}=14ϵ_0`$, for both GaAs and for InAs. In the absence of strain the GaAs - InAs conduction band offset between the two materials is 770 meV. To partially take into account strain effects, we have assumed a constant shift of 400 meV in the InAs conduction band, and no strain effect in the GaAs conduction band, resulting in a conduction band offset of $`770\mathrm{m}\mathrm{e}\mathrm{V}400\mathrm{m}\mathrm{e}\mathrm{V}=370\mathrm{m}\mathrm{e}\mathrm{V}`$. We also take the electron effective masses to be constant in each material.

In order to study a realistic potential, we have used a nonuniform grid basis for the LSDA calculation and the solution of the Poisson equation. For the Schrödinger equation we set the wave functions to zero at the edge of the grid, which is reasonable due to the exponential decay of bound states, but the Poisson equation requires more care. We used a multipole expansion up to the quadrupole term to simulate the boundary conditions at infinity. This difficulty with the boundary conditions of the Poisson equation is one source of error in our comparison. We used cubic interpolation to map the gridded wave functions and conduction band potential to continuous functions for Monte Carlo evaluation of the exchange energy, for determination of the trial nodes in the DMC calculation, and for guiding wave functions in the DMC calculation. Although the cubic interpolation of the potential allows us to compare LSDA and DMC for similar external potentials $`v_{\mathrm{cb}}`$, errors in the kinetic energy operator are more difficult. The use of a finite difference approximation to the Laplacian in the computation of the LSDA solution creates an operator that has no simple counterpart in a continuous formulation of the problem. The error between the physical Laplacian operator and the artificial finite difference operator can be made arbitrarily small by the use of finer grids, but computer memory and CPU time constraints caused this to be a significant source of error in our comparison. The size of the errors are several meV, which are comparable to the rather small LSDA errors, especially in the correlation energy. Below we discuss a way we estimated this error so that we can compare correlation energy.

For the set of calculations we describe now, we have taken all masses to have the isotropic value of $`0.05m_e`$, for both the InAs and GaAs regions. The case of anisotropic masses will be discussed later. We have also accounted for numerical discrepancies between the two calculations due to grid interpolation error. To estimate this error we have compared eigenvalues of single electron states with both methods, and found that the DMC single particle eigenvalues typically lie about 0.25 meV above the eigenvalues computed by the grid method used in the DFT calculation. We attribute that difference to a systematic error in the finite difference kinetic energy operator and shift up our LSDA calculations to compensate for this.

Figure 2 shows the values of the different components of the total energy as a function of dot occupation. Differences between LSDA and QMC are not apparent on the scale of the figure. This figure clearly shows that the external potential energy and kinetic energy are much larger than the interaction energies. In other words the effects of interactions enter as a perturbative correction to the non-interacting system. The reason for this can be seen from the scaled electron density. The energy and length scales for the electron interaction are scaled by the dielectric constant and mass, so that the effective Bohr radius is $`a_0^{}=ϵ/m^{}a_0150\mathrm{\AA }`$ and the effective Hartree is $`\mathrm{Ha}^{}=ϵ^2m^{}\mathrm{Ha}7\mathrm{m}\mathrm{e}\mathrm{V}`$. If we estimate the density by assuming the $`N`$ electrons uniformly occupy the interior of the dot, we find an effective conduction electron density of $`r_s0.46N^{1/3}`$. To see the consequences of such a high effective density, consider the uniform electron gas at $`r_s=0.6`$. The electron gas has a ground state energy expansion for small $`r_s`$,

$$E=2.2099r_s^20.9163r_s^10.094+0.0622\mathrm{ln}(r_s)+\mathrm{},$$
(17)
where the first term is the kinetic energy, the next term is the exchange, and remaining terms are correlation energy. For the case of six electrons in the dot $`r_s0.25`$ and the expansion gives $`E_{\mathrm{kin}}1440\mathrm{m}\mathrm{e}\mathrm{V}`$, $`E_\mathrm{x}150\mathrm{m}\mathrm{e}\mathrm{V}`$, and $`E_\mathrm{c}0.6\mathrm{meV}`$. Although the comparison between these very different electronic systems cannot be pushed too far, this does show that our exchange and correlation energies are reasonable for this electron density. The leading effect of the interactions is the Hartree energy, with small corrections for exchange and very small correlation corrections.

The most important contribution of LSDA in this system is the exchange term, and we plot the comparison of the LSDA exchange energy in Figure 3. For finite systems, the definition for exchange includes the correction for self-interaction in the Hartree energy,  which we also show in the figure. Merely correcting for self-interaction will recover at least 75% of the exchange energy, and LSDA is able to do better, recovering about 90% for an error in the exchange energy of less than 5 meV. As is well known in many other systems , the error in the local approximation to the exchange is the largest error in LSDA.

By comparison the correlation energy is much smaller and its accuracy is difficult to assess due to grid errors. Correlation is defined as the difference in total energies, Eq. (16), and in this case the total energy difference is less than 1% . Thus our estimates of the exact correlation energy are fairly uncertain. In Figure 4 we show our best attempt at checking the LSDA correlation in the quantum dot. Since the LSDA wave functions are not the true Kohn-Sham orbitals (due mostly to grid errors) the correlation is overestimated due to relaxation of our interpolated LSDA states. As stated earlier, we correct for this by performing single particle calculations to estimate the relaxation of our LSDA states, and have removed this contribution from our calculated correlation energy. These results are shown as “QMC” and “QMC with corrections” in the figure. We find that LSDA gives a reasonable estimate of the magnitude of the correlation energy in this system.

Our comparison of the energy errors from the LSDA are summarized in Figure 5. The largest error is due to LSDA underestimating the magnitude of the exchange energy, $`E_\mathrm{x}`$, which is a 6 meV error for six electrons in the dot. There is a small error in the correlation energy, $`E_\mathrm{c}`$, which slightly compensates for the error in the exchange. The error in the Hartree energy $`E_{\mathrm{Ha}}`$ is due primarily to truncation of the multipole expansion at quadrupole terms for determination of the boundary conditions of the Poison equation. There are also errors in $`E_{\mathrm{kin}}`$ and $`E_{\mathrm{cb}}`$ (not shown in the figure) due to the grid, in particular the discretization of the kinetic energy operator discussed earlier. We therefore find that LSDA leads to an error in the total energy of up to 6 meV, which in our calculation is partially compensated for by grid errors in other energy terms, leading to an overall error in the total energy of no more than 2 meV.

Because of strain, the electron effective mass in the plane parallel to the base of the dot is usually different from the mass perpendicular to it.  The LSDA does not explicitly account for mass anisotropy in the exchange and correlation functionals, thus for the calculation of the exchange-correlation potential $`V_{\mathrm{xc}}`$ we have assumed a constant average mass given by

$$m_{\mathrm{xc}}=(m_{}^2m_{})^{1/3},$$
(18)
where $`m_{}`$ is the in-plane mass, $`m_{}=m_{\mathrm{xx}}=m_{\mathrm{yy}}`$, and $`m_{}`$ is the perpendicular mass, $`m_{}=m_{\mathrm{zz}}`$. The choice of the average mass $`m_{\mathrm{xc}}`$ is somewhat arbitrary. We have performed a test of this approximation by calculating the exchange and correlation energy for an anisotropic system using (i) LSDA with $`m_{\mathrm{xc}}`$ as given above, (ii) LSDA with $`m_{\mathrm{xc}}`$ optically averaged, $`m_{\mathrm{xc}}^1=(2/m_{}+1/m_{})/3`$, and (iii) DMC, which can explicitly treat anisotropic mass. Figure 6 shows the calculated exchange and correlation energies for several electron mass ratios, which we have again fixed to constant across the InAs and GaAs. We find that either choice is acceptable, as both give errors comparable to the isotropic mass case.

## V Conclusions

In conclusion, we have shown that DFT offers an accurate approximation for the ground state energy including many-body interactions in small self-assembled quantum dots, while providing a simple and fast means of modeling systems containing several electrons. This is in large part due to the small size of the self-assembled quantum dots in relation to the effective Bohr radius. In this regime the largest error is in the local approximation to the exchange, and we have verified that the errors due to this approximation are small. This is quantified in Figures 2 and 3 where we see that the exchange energy for six electrons in the dot is 72 meV with an LSDA error of 6 meV, compared to a total energy of 1490 meV. We expect this applicability of LSDA to hold in general for small quantum dots of various shapes, but we emphasize that LSDA is expected to become progressively worse if the dot becomes much larger than the effective Bohr radius. Also, the use of the EMA Schrödinger equation, Eq. (1), although reasonable for the chosen problem, is predicted to fail for very small dots.  Although our model assumes pure InAs dots surrounded by GaAs, in some systems indium and gallium may in fact intermix.  Again, we expect the relative magnetude of exchange and correlation energies and the acceptability of LSDA to be unchanged by intermixing or alloying of the dots.

Based on this comparison we conclude that the addition and removal energies found by EMA-LSDA calculations in realistic dots  are accurate to 1-2 meV per electron as far as the many-body interactions within the EMA are concerned. Thus, comparison of such calculated charging energies with experiment can be considered a direct test of the models and the uncertainties in the analysis of experiments.  Furthermore, we have carefully tested the grid errors and concluded that energies are accurate to $`\pm 10`$ meV including all errors.

We would like to point out that application of LSDA to a system of several coupled self-assembled dots may have difficulties. Systems that have weak intersite transitions along with intrasite interactions are a well-known many-body problem where LSDA is known to fail. The problem is closely related to the $`\mathrm{H}_2`$ molecule, where LSDA describes well the electronic energy at equilibrium distance, but gives an incorrect broken symmetry solution if the atoms are pulled apart beyond a critical distance. The problem of weakly coupled dots would be of particular importance for future studies.

## Acknowledgments

This work was supported by CRI from the University of Illinois and NSF Grants No. ECS 95-09751 and No. DMR 94-22496 and computer resources at NCSA.


text:
# Dynamics of Multiparticle Systems with non – Abelian Symmetry

## 1 Introduction and Overview

The consideration of the influence of internal symmetries on the final state of a many body system begun with the pioneering work of Bethe . Much of the subsequent interest in the subject arises from the realization that in the study of hadronic interactions and in particular in studies involving quark confinement, these constraints may be of decisive importance. An important progress in treating equilibrium systems was made employing group projection techniques. This allowed for a consistent treatment of abelian  and nonabelian  symmetries of compact groups and a consistent formulation of thermodynamics of many particle systems with internal symmetries taken into account . Application of these methods to specific processes demonstrated in which circumstances the presence of symmetry is of physical relevance .

However, it is not fully understood how the symmetry-modified properties of the equilibrium system arise from kinetic formulation of the dynamical evolution. When an internal symmetry is not at work, Boltzmann’s H theorem in principle assures that the statistical Bose/Fermi/Boltzmann distributions are the asymptotic (equilibrium) distributions, irrespective of the nature of microscopic interaction. However, in presence of exact symmetries the equilibrium distributions are modified, see e.g. . This implies that symmetry constraints introduce effective interactions of potentially far more complex nature than is the usual two body Boltzmann collision term. In fact it can be argued that quantum symmetry constraints are the heart of the nonlocality of quantum physics. However, in the limit of classical Boltzmann equation evolution these are implemented by a strictly local (though non-linear) consideration of Fermi blocking and Bose enhancement in phase space evolution. Our aim in this work is to make a step towards understanding how the microscopic nonabelian symmetry constrains operate within the kinetic master equation description of the time evolution, leading on to the symmetry modified (constrained) macroscopic many particle equilibrium state.

It is first important to convince oneself that an underlying symmetry of microscopic interactions does not lead in general to the desired symmetry properties of a (macroscopic) many body interacting system. To do this, we consider the high energy nuclear (heavy ion) collisions and specifically here two symmetry examples:

a) SU(2) Isospin symmetry: The initial state transforms under a given representation of the isospin $`SU(2)`$ group. All elementary high energy interactions are governed by the strong interaction, which preserves the isotopic symmetry. A final state results as a multiparticle state formed by many individual hadron – hadron collisions. In any of such microprocesses the isospin is conserved. However, proceeding ‘as usual’ without symmetry constrained treatment of local interactions does not assure that the final multiparticle state (macrostate) transforms under the same representation of the isospin group as the initial state, which is required for symmetry reasons.

b) SU(3) Colour symmetry: A similar situation appears in the context of the quark-gluon interactions, especially in case that local deconfinement occurs. The initial state is a colour singlet state, and quark-gluon and gluon-gluon interaction, although invariant under the colour $`SU(3)_c`$ symmetry group do not assure that during its evolution a (macroscopic) many particle state, once a singlet, always remains a singlet colour state, which, however, it must do because of exact colour symmetry of strong interactions.

As these examples show, in a dynamical (quantum) transport theory description of the approach to equilibrium there must exist a subsidiary condition which should be taken into account by corresponding kinetic equations governing the evolution. This condition is independent from other constraints related to dynamical gauged internal symmetries. For classical fluid dynamics a dynamical evolution equation addressing gauge symmetry has been proposed by Wong . In many current studies of the dynamics of classical non-abelian fields this proposal continues till today to attract considerable interest . However, these are constraints which have no relation to the intrinsic non-locality of the quantum system which we address here.

We first note that in the case of an abelian symmetry there are no additional constraints to consider. Quantum number conservation on a microscopic level is fully equivalent to preservation of all symmetry properties on the macroscopic level. This is easily seen considering the $`U(1)`$ symmetry related to microscopic particle-antiparticle formation: since the microscopic mechanisms produce equal number of particles and antiparticles (pair production), initial particle-antiparticle number difference is exactly preserved in the macroscopic many body state.

Thus only presence of a nonabelian symmetry poses a true challenge. A suitable mathematical method how to approach this problem is identified considering the previously treated statistical equilibrium case. Here one decomposes a general ‘macrostate’ consisting of many particles into possible irreducible representations of the symmetry group. Then a projection technique exploiting character function properties of the group is used to constrain the final state. We will here use this approach in order to describe multiparticle evolution applying microscopic kinetic theory scheme. We will show that the behaviour of particle phase - space distribution functions will depend not only on properties of basic interaction, but that it also depends on global properties of the macroscopic system. Those global properties provide subsidiary constraints needed, so that the asymptotic equilibrium state has properties consistent with the non-abelian constraints.

## 2 The Projection Method

Let $`G`$ be a compact internal symmetry group of our system consisting of particles (objects) transforming under irreducible representations of the symmetry group. These representation are denoted as $`\alpha _i`$ with corresponding dimensions $`d(\alpha _i)`$. One denotes $`f_{(\zeta )}^{(\alpha _i,\nu _i)}(\mathrm{\Gamma },\stackrel{}{r},t)`$ a distribution function of the particle which belongs to the multiplet $`\alpha _i`$ of the symmetry group. Members of this multiplet are numbered by indexes $`\nu _i(\nu _i=1,\mathrm{},d(\alpha _i))`$ which correspond to given values of charges related to the symmetry group. A subscript $`\zeta `$ denotes other quantum numbers characterizing different multiplets of the same representation $`\alpha `$. The variables $`(\mathrm{\Gamma },\stackrel{}{r})`$ denotes a set of the phase - space variables such as $`(\stackrel{}{p},\stackrel{}{r})`$ and $`t`$ is time.

The number of particles of the specie $`\{\alpha ,\nu _\alpha ,\zeta \}`$ is:

$$N_{\nu _\alpha ;(\zeta )}^{(\alpha )}(t)=𝑑V𝑑\mathrm{\Gamma }f_{(\zeta )}^{(\alpha ,\nu _\alpha )}(\mathrm{\Gamma },\stackrel{}{r},t);$$
(1)
We consider a system of $`\{N_{\alpha _1,\nu _{\alpha _1}}^{(\zeta _1)}(t),\mathrm{},N_{\alpha _n,\nu _{\alpha _n}}^{(\zeta _n)}(t)\}`$ particles at time $`t`$ . The distribution functions fulfill the generalized Vlasov \- Boltzmann kinetic equations, which can be written in the general form:

$`{\displaystyle \frac{f_{(\zeta _i)}^{(\alpha _i,\nu _i)}(\mathrm{\Gamma }_i,\stackrel{}{r},t)}{t}}`$ $`+`$ $`\stackrel{}{v}f_{(\zeta _i)}^{(\alpha _i,\nu _i)}(\mathrm{\Gamma }_i,\stackrel{}{r},t)`$ (2)
$`=`$ $`{\displaystyle \underset{\alpha _j,\alpha _k,\alpha _l}{}}{\displaystyle \underset{\nu _j,\nu _k,\nu _l}{}}{\displaystyle \underset{\zeta _j,\zeta _k,\zeta _l}{}}{\displaystyle 𝑑\mathrm{\Gamma }_j𝑑\mathrm{\Gamma }_k𝑑\mathrm{\Gamma }_l𝒲_{\nu _i\nu _j;\nu _k\nu _l}^{(\zeta _i,\zeta _j;\zeta _k,\zeta _l)}(\mathrm{\Gamma }_k,\mathrm{\Gamma }_l;\mathrm{\Gamma }_j,\mathrm{\Gamma }_i)}`$
$`[_{(\zeta _i)}^{(\alpha _i,\nu _i)}(\mathrm{\Gamma }_i,\stackrel{}{r},t)_{(\zeta _j)}^{(\alpha _j,\nu _j)}(\mathrm{\Gamma }_j,\stackrel{}{r},t)f_{(\zeta _k)}^{(\alpha _k,\nu _k)}(\mathrm{\Gamma }_k,\stackrel{}{r},t)f_{(\zeta _l)}^{(\alpha _l,\nu _l)}(\mathrm{\Gamma }_l,\stackrel{}{r},t)`$
$`_{(\zeta _k)}^{(\alpha _k,\nu _k)}(\mathrm{\Gamma }_k,\stackrel{}{r},t)_{(\zeta _l)}^{(\alpha _l,\nu _l)}(\mathrm{\Gamma }_l,\stackrel{}{r},t)f_{(\zeta _i)}^{(\alpha _i,\nu _i)}(\mathrm{\Gamma }_i,\stackrel{}{r},t)f_{(\zeta _j)}^{(\alpha _j,\nu _j)}(\mathrm{\Gamma }_j,\stackrel{}{r},t)];`$
Factors $`_{(\zeta )}^{(\alpha ,\nu )}(\mathrm{\Gamma },\stackrel{}{r},t)`$ are related to quantum statistics and they are equal to $`1`$ for classical particles, and equal to $`[1\pm f_{(\zeta )}^{(\alpha ,\nu )}(\mathrm{\Gamma },\stackrel{}{r},t)]`$ for bosons/fermions correspondingly. Since by assumption the whole system transforms under given representation $`\mathrm{\Lambda }`$ of an exact symmetry group, the system under consideration must preserve its transformations properties during its time evolution, provided that it is governed by a symmetry invariant interaction.

We now focus on the case of a quantum system and consider state vectors in particle number representation: $`|N_{\nu _{\alpha _1}}^{(\alpha _1)},\mathrm{},N_{\nu _{\alpha _n}}^{(\alpha _n)}.`$ These vectors describe symmetry properties of our systems and all other variables, related to phase-space properties of the system are suppressed here. They transform as a direct product representation of the symmetry group $`G`$. This representation is of the form:

$$\alpha _1^{N^{(\alpha _1)}}\alpha _2^{N^{(\alpha _2)}}\mathrm{}\alpha _n^{N^{(\alpha _n)}};$$
(3)
A multiplicity $`N^{(\alpha _j)}`$ of the representation $`\alpha _j`$ in this product is equal to a number of particles which transform under this representation:

$$N^{(\alpha _j)}=\underset{j}{}\left(\underset{\zeta _j}{}N_{\nu _{\alpha _j};(\zeta _j)}^{(\alpha _j)}\right)=\underset{j}{}N_{\nu _{\alpha _j}}^{(\alpha _j)};$$
(4)
The representation given by Eq. (3) can be decomposed into direct sum of irreducible representations $`\mathrm{\Lambda }_k`$. Corresponding states are denoted as $`|\mathrm{\Lambda }_k,\lambda _{\mathrm{\Lambda }_k};𝒩`$ where $`\lambda _{\mathrm{\Lambda }_k}`$ is an index numbering members of the representation $`\mathrm{\Lambda }`$ and $`𝒩`$ is a total number of particles

$$𝒩=\underset{k}{}N_{\nu _{\alpha _k}}^{(\alpha _k)};$$
(5)
Each physical state can be decomposed into irreducible representation base states with amplitudes depending on phase space variables $`\mathrm{\Gamma }`$:

$`|N_{\nu _{\alpha _1}}^{(\alpha _1)},\mathrm{},N_{\nu _{\alpha _n}}^{(\alpha _n)};\mathrm{\Gamma }={\displaystyle \underset{k}{\overset{}{}}}{\displaystyle \underset{\xi _{\mathrm{\Lambda }_k}}{\overset{}{}}}|\mathrm{\Lambda }_k,\lambda _{\mathrm{\Lambda }_k};𝒩;\xi _{\mathrm{\Lambda }_k}a_{\{N_{\nu _{\alpha _1}}^{(\alpha _1)},\mathrm{},N_{\nu _{\alpha _n}}^{(\alpha _n)}\}}^{\mathrm{\Lambda },\lambda _\mathrm{\Lambda }}(\xi _{\mathrm{\Lambda }_k};\mathrm{\Gamma });`$ (6)
Here appear new variables $`\xi _\mathrm{\Lambda }`$ which are degeneracy parameters required for the full description of a state in the ”symmetry space”.

Let us define an average weight

$$\overline{P_{\{N_{\nu _{\alpha _1}}^{(\alpha _1)},\mathrm{},N_{\nu _{\alpha _n}}^{(\alpha _n)}\}}^{\mathrm{\Lambda },\lambda _\mathrm{\Lambda }}}=\frac{\underset{\xi _\mathrm{\Lambda }}{}|a_{\{N_{\nu _{\alpha _1}}^{(\alpha _1)},\mathrm{},N_{\nu _{\alpha _n}}^{(\alpha _n)}\}}^{\mathrm{\Lambda },\lambda _\mathrm{\Lambda }}(\xi _\mathrm{\Lambda };\mathrm{\Gamma })|^2}{\underset{N_{\nu _{\alpha _1}}^{(\alpha _1)}+\mathrm{}+N_{\nu _{\alpha _n}}^{(\alpha _n)}=𝒩}{}\underset{\xi _\mathrm{\Lambda }}{}|a_{\{N_{\nu _{\alpha _1}}^{(\alpha _1)},\mathrm{},N_{\nu _{\alpha _n}}^{(\alpha _n)}\}}^{\mathrm{\Lambda },\lambda _\mathrm{\Lambda }}(\xi _\mathrm{\Lambda };\mathrm{\Gamma })|^2};$$
(7)
This expression gives the probability that $`N_{\nu _{\alpha _1}}^{(\alpha _1)},\mathrm{},N_{\nu _{\alpha _n}}^{(\alpha _n)}`$ particles transforming under the symmetry group representations $`\alpha _1,\mathrm{},\alpha _n`$ combine into $`𝒩`$ particle state transforming under representation $`\mathrm{\Lambda }`$ of the symmetry group.

We make the statistical hypothesis that average weights (7) do not depend on phase - space variables and can be calculated alone on basis of symmetry group consideration. This also can be proved under the stronger assumption that in Eq. (6) any state with fixed $`\mathrm{\Lambda },\lambda _\mathrm{\Lambda }`$ has the same weight (see e.g. )

Let us consider a projection operator $`𝒫^\mathrm{\Lambda }`$ on the subspace spanned by all states transforming under representation $`\mathrm{\Lambda }`$.

$$𝒫^\mathrm{\Lambda }|N_{\nu _{\alpha _1}}^{(\alpha _1)},\mathrm{},N_{\nu _{\alpha _n}}^{(\alpha _n)}=\underset{\xi _\mathrm{\Lambda }}{\overset{}{}}|\mathrm{\Lambda },\lambda _\mathrm{\Lambda };\xi _\mathrm{\Lambda }𝒞_{\{N_{\nu _{\alpha _1}}^{(\alpha _1)},\mathrm{},N_{\nu _{\alpha _n}}^{(\alpha _n)}\}}^{\mathrm{\Lambda },\lambda _\mathrm{\Lambda }}(\xi _\mathrm{\Lambda });$$
(8)
This operator has the generic form (see e.g. ):

$$𝒫^\mathrm{\Lambda }=d(\mathrm{\Lambda })\underset{G}{}𝑑\mu (g)\overline{\chi }^{(\mathrm{\Lambda })}(g)U(g);$$
(9)
Here $`\chi ^{(\mathrm{\Lambda })}`$ is the character of the representation $`\mathrm{\Lambda }`$, $`d\mu (g)`$ is the invariant Haar measure on the group, and $`U(g)`$ is an operator transforming a state under consideration. We will use the matrix representation:

$`U(g)|N_{\nu _{\alpha _1}}^{(\alpha _1)},\mathrm{},N_{\nu _{\alpha _n}}^{(\alpha _n)}`$ (10)
$`=`$ $`{\displaystyle \underset{\nu _1^{(1)},\mathrm{},\nu _n^{(N_{\nu _n})}}{}}D_{\nu _1^{(1)}\nu _1}^{(\alpha _1)}\mathrm{}D_{\nu _1^{(N_{\nu _1})}\nu _1}^{(\alpha _1)}\mathrm{}D_{\nu _n^{(1)}\nu _n}^{(\alpha _n)}\mathrm{}D_{\nu _n^{(N_{\nu _n})}\nu _n}^{(\alpha _n)}|N_{\nu _{\alpha _1}}^{(\alpha _1)},\mathrm{},N_{\nu _{\alpha _n}}^{(\alpha _n)};`$
$`D_{\nu ,\nu }^{(\alpha _n)}`$ is a matrix elements of the group element $`g`$ corresponding to the representation $`\alpha `$. Notation convention in Eq. (10) arises since there are $`N_{\nu _{\alpha _j}}^{(\alpha _j)}`$ states transforming under representation $`\alpha _j`$ and having quantum numbers of the $`\nu _{\alpha _j}`$-th member of a given multiplet.

The statistical hypothesis identifies the average weight $`\overline{P_{\{N_{\nu _{\alpha _1}}^{(\alpha _1)},\mathrm{},N_{\nu _{\alpha _n}}^{(\alpha _n)}\}}^{\mathrm{\Lambda },\lambda _\mathrm{\Lambda }}}`$ with a norm of the vector $`𝒫^\mathrm{\Lambda }|N_{\nu _{\alpha _1}}^{(\alpha _1)},\mathrm{},N_{\nu _{\alpha _n}}^{(\alpha _n)}`$. This norm can be written as

$`N_{\nu _{\alpha _1}}^{(\alpha _1)},\mathrm{},N_{\nu _{\alpha _n}}^{(\alpha _n)}\left|𝒫^\mathrm{\Lambda }\right|N_{\nu _{\alpha _1}}^{(\alpha _1)},\mathrm{},N_{\nu _{\alpha _n}}^{(\alpha _n)}={\displaystyle \underset{\xi _\mathrm{\Lambda }}{}}|𝒞_{\{N_{\nu _{\alpha _1}}^{(\alpha _1)},\mathrm{},N_{\nu _{\alpha _n}}^{(\alpha _n)}\}}^{\mathrm{\Lambda },\lambda _\mathrm{\Lambda }}(\xi _\mathrm{\Lambda })|^2;`$ (11)
where the relation $`(𝒫^\mathrm{\Lambda })^2=𝒫^\mathrm{\Lambda }`$ was used.

Left hand side of this equation can be calculated directly from Eqs.(9) and (10). One gets finally

$$\overline{P_{\{N_{\nu _{\alpha _1}}^{(\alpha _1)},\mathrm{},N_{\nu _{\alpha _n}}^{(\alpha _n)}\}}^{\mathrm{\Lambda },\lambda _\mathrm{\Lambda }}}=𝒜^{\{𝒩\}}d(\mathrm{\Lambda })\underset{G}{}𝑑\mu (g)\overline{\chi }^{(\mathrm{\Lambda })}(g)[D_{\nu _1\nu _1}^{(\alpha _1)}]^{N_{\nu _{\alpha _1}}^{(\alpha _1)}}\mathrm{}[D_{\nu _n\nu _n}^{(\alpha _n)}]^{N_{\nu _{\alpha _n}}^{(\alpha _n)}};$$
(12)
where $`𝒜^{\{𝒩\}}`$ is a permutation normalization factor. For particles of the kind $`\{\alpha ,\zeta \}`$ we included in Eq. (12) the permutation factor:

$$𝒜_{(\zeta )}^\alpha =\frac{𝒩_{(\zeta )}^{(\alpha )}!}{\underset{\nu _\alpha }{}𝒩_{\nu _\alpha ;(\zeta )}^{(\alpha )}!};$$
(13)
The permutation factor $`𝒜^{\{𝒩\}}`$ is a product of all ”partial” factors

$$𝒜^{\{𝒩\}}=\underset{j}{}\underset{\zeta _j}{}𝒜_{(\zeta _j)}^{\alpha _j};$$
(14)
The permutation factor assures the normalization of state vectors:

$$N_{\nu _{\alpha _1}}^{(\alpha _1)},\mathrm{},N_{\nu _{\alpha _n}}^{(\alpha _n)}|N_{\nu _{\alpha _1}}^{(\alpha _1)},\mathrm{},N_{\nu _{\alpha _n}}^{(\alpha _n)}=𝒜^{\{𝒩\}};$$
(15)
This normalization reflects an invariance of the state vector with respect to permutations which shuffle indistinguishable particles.

## 3 Incorporation of Symmetry

The expression Eq. (12) is a starting point for further considerations. It provides together with Eq. (1) and Eq. (2) subsidiary constraints on distribution functions $`f^{(\alpha _i,\nu _i)}`$. These conditions assure that in a dynamical evolution the symmetry of the system is preserved. When symmetry is conserved, then all weights in Eq. (12) are constant in time. In a case of strong interaction and colour symmetry, all weights, except for the weight corresponding to the singlet state, must remain zero.

We now convert the global constraint into a time evolution condition and consider:

$$\frac{d}{dt}\overline{P_{\{N_{\nu _{\alpha _1}}^{(\alpha _1)},\mathrm{},N_{\nu _{\alpha _n}}^{(\alpha _n)}\}}^{\mathrm{\Lambda },\lambda _\mathrm{\Lambda }}}=0;$$
(16)
Introducing here the result of Eq. (12) one obtains:

$`0`$ $`=`$ $`{\displaystyle \frac{d𝒜^{\{𝒩\}}}{dt}}d(\mathrm{\Lambda }){\displaystyle \underset{G}{}}𝑑\mu (g)\overline{\chi }^{(\mathrm{\Lambda })}(g)[D_{\nu _1\nu _1}^{(\alpha _1)}]^{N_{\nu _{\alpha _1}}^{(\alpha _1)}}\mathrm{}[D_{\nu _n\nu _n}^{(\alpha _n)}]^{N_{\nu _{\alpha _n}}^{(\alpha _n)}}`$ (17)
$`+`$ $`{\displaystyle \underset{j=1}{\overset{n}{}}}{\displaystyle \underset{\nu _{\alpha _j}}{}}{\displaystyle \frac{dN_{\nu _{\alpha _j}}^{(\alpha _j)}}{dt}}𝒜^{\{𝒩\}}d(\mathrm{\Lambda }){\displaystyle \underset{G}{}}𝑑\mu (g)\overline{\chi }^{(\mathrm{\Lambda })}(g)[D_{\nu _1\nu _1}^{(\alpha _1)}]^{N_{\nu _{\alpha _1}}^{(\alpha _1)}}\mathrm{}[D_{\nu _n\nu _n}^{(\alpha _n)}]^{N_{\nu _{\alpha _n}}^{(\alpha _n)}}\mathrm{log}[D_{\nu _j\nu _j}^{(\alpha _j)}];`$
All integrals which appear in Eq. (12) and Eq. (17) can be expressed explicitly in an analytic form for any compact symmetry group.

To write an expression for the time derivative of the normalization factor $`𝒜^{\{𝒩\}}`$ we perform analytic continuation from integer to continuous values of variables $`N_{\nu _{\alpha _n}}^{(\alpha _n)}.`$ Thus we replace all factorials by the $`\mathrm{\Gamma }`$–function of corresponding arguments. We encounter here also the digamma function $`\psi `$ :

$$\psi (x)=\frac{d\mathrm{log}\mathrm{\Gamma }(x)}{dx};$$
(18)
This allows to write for Eq. (17):

$$\frac{d𝒜^{\{𝒩\}}}{dt}=𝒜^{\{𝒩\}}\underset{j}{}\underset{\zeta _j}{}\left[\frac{d𝒩_{(\zeta _j)}^{(\alpha _j)}}{dt}\psi (𝒩_{(\zeta _j)}^{(\alpha _j)}+1)\underset{\nu _{\alpha _j}}{}\frac{d𝒩_{\nu _{\alpha _j};(\zeta _j)}^{(\alpha _j)}}{dt}\psi (𝒩_{\nu _\alpha ;(\zeta _j)}^{(\alpha _j)}+1)\right];$$
(19)
To get a consistent analytical continuation in the number of particles one should define the time derivatives $`d𝒩_{\nu _\alpha ;(\zeta )}^{(\alpha )}/dt`$. We define these rates of particle number change from the integrated Boltzmann kinetic equation, Eq. (2), explicitly

$`{\displaystyle \frac{dN_{\nu _{\alpha _i}}^{(\alpha _i)}}{dt}}`$ $`=`$ $`{\displaystyle \underset{\alpha _j,\alpha _k,\alpha _l}{}}{\displaystyle \underset{\nu _j,\nu _k,\nu _l}{}}{\displaystyle \underset{\zeta _j,\zeta _k,\zeta _l}{}}{\displaystyle 𝑑V𝑑\mathrm{\Gamma }_j𝑑\mathrm{\Gamma }_k𝑑\mathrm{\Gamma }_l𝑑\mathrm{\Gamma }_i𝒲_{\nu _i\nu _j;\nu _k\nu _l}^{(\zeta _i,\zeta _j;\zeta _k,\zeta _l)}(\mathrm{\Gamma }_k,\mathrm{\Gamma }_l;\mathrm{\Gamma }_j,\mathrm{\Gamma }_i)}`$ (20)
$`[_{(\zeta _i)}^{(\alpha _i,\nu _i)}(\mathrm{\Gamma }_i,\stackrel{}{r},t)_{(\zeta _j)}^{(\alpha _j,\nu _j)}(\mathrm{\Gamma }_j,\stackrel{}{r},t)f_{(\zeta _k)}^{(\alpha _k,\nu _k)}(\mathrm{\Gamma }_k,\stackrel{}{r},t)f_{(\zeta _l)}^{(\alpha _l,\nu _l)}(\mathrm{\Gamma }_l,\stackrel{}{r},t)`$
$`_{(\zeta _k)}^{(\alpha _k,\nu _k)}(\mathrm{\Gamma }_k,\stackrel{}{r},t)_{(\zeta _l)}^{(\alpha _l,\nu _l)}(\mathrm{\Gamma }_l,\stackrel{}{r},t)f_{(\zeta _i)}^{(\alpha _i,\nu _i)}(\mathrm{\Gamma }_i,\stackrel{}{r},t)f_{(\zeta _j)}^{(\alpha _j,\nu _j)}(\mathrm{\Gamma }_j,\stackrel{}{r},t)];`$
Contributions from gradient terms of Eq. (2) vanish due to Gauss law. These terms are transformed in surface integrals and beyond the volume occupied by the system all distribution functions are equal to zero.

Eqs (17,19,20) in fact constitute the global subsidiary condition which should be fulfilled by the microscopic kinetic equations Eq. (2). These are the necessary conditions for preserving the internal symmetry on the macroscopic level. Rates of change $`d𝒩_{\nu _\alpha ;(\zeta )}^{(\alpha )}/dt`$ are related to “macrocurrents”, which are counterparts of “microcurrents” related directly to a symmetry on a microscopic level via the Noether theorem. Eq. (20) can be considered as a set of conditions on macrocurrents to provide consistency with the overall symmetry of the system. Therefore we believe that this equation can also be used as a starting point for multicomponent hydrodynamic equations with internal symmetry properties taken into account. One should notice that this subsidiary condition takes into account also surface effects for the finite volume systems. This is due to the space variables integration which is performed in Eqs (1) and (20).

One easily sees that for the case of abelian symmetry the two constraints Eq. (12) and Eq. (17) do not lead to new results: first we recall that all irreducible representations of abelian group are one-dimensional. Next, let basic particles have “charges” $`q_1,\mathrm{},q_n`$ , and let the global charge be $`Q`$. Then the only consequence of Eq. (12) follows for nonvanishing weight $`Q=N_1q_1+\mathrm{}+N_nq_n`$, which is a rather obvious result. New results appear only for nonabelian symmetries.

## 4 Example: Isospin

We now consider as an example the case of the $`SU(2)`$ symmetry with basic particles transforming under spinor $`(\frac{\mathrm{𝟏}}{\mathrm{𝟐}})`$ (fundamental) and vector $`(\mathrm{𝟏})`$ (adjoint) representations. This example can be realized by a gas mixture of nucleons and pions. To describe all group elements the three group’s parameters $`\alpha ,\beta ,\gamma `$ are chosen in such a way that diagonal matrix elements have the well known form :
i.) for the fundamental representation $`(\frac{\mathrm{𝟏}}{\mathrm{𝟐}})`$:

$`D_{mm}^{(1/2)}(\alpha ,\beta ,\gamma )=e^{im(\alpha +\gamma )}\mathrm{cos}{\displaystyle \frac{\beta }{2}};m=\pm {\displaystyle \frac{1}{2}};`$ (21)
ii.) and for the adjoint representation $`(\mathrm{𝟏})`$:

$`D_{\pm 1,\pm 1}^{(1)}(\alpha ,\beta ,\gamma )`$ $`=`$ $`{\displaystyle \frac{1}{2}}e^{\pm im(\alpha +\gamma )}(1+\mathrm{cos}\beta );`$ (22)
$`D_{0,0}^{(1)}(\alpha ,\beta ,\gamma )`$ $`=`$ $`\mathrm{cos}\beta ;`$ (23)
The Haar measure for the $`SU(2)`$ group in this parametrization has the form

$`{\displaystyle 𝑑\mu (g)f[g]}={\displaystyle \frac{1}{8\pi ^2}}{\displaystyle \underset{0}{\overset{2\pi }{}}}𝑑\alpha {\displaystyle \underset{0}{\overset{2\pi }{}}}𝑑\gamma {\displaystyle \underset{0}{\overset{\pi }{}}}𝑑\beta \mathrm{sin}\beta f[g(\alpha ,\beta ,\gamma )];`$ (24)
Any ‘macrostate’ is made of an arbitrary number: $`n_n,n_p,n_{},n_0,n_+`$; where subscripts refer to members of the fundamental representation neutrons, protons; and members of adjoint representations, $`\pi ^{},\pi ^0,\pi ^+`$, correspondingly. Let us consider the special case when the macrostate is a $`SU(2)`$ singlet. The weight of the singlet state is according to Eq. (12):

$`\overline{P_{\{n_n,n_p,n_{},n_0,n_+\}}^{0,0}}=`$ (25)
$`𝒜^{\{𝒩\}}{\displaystyle \frac{1}{8\pi ^2}}{\displaystyle \underset{0}{\overset{2\pi }{}}}𝑑\alpha {\displaystyle \underset{0}{\overset{2\pi }{}}}𝑑\gamma {\displaystyle \underset{0}{\overset{\pi }{}}}𝑑\beta \mathrm{sin}\beta e^{\frac{i}{2}(n_nn_p+2n_{}2n_+)(\alpha +\gamma )}\mathrm{cos}^{}{\displaystyle \frac{\beta }{2}}\mathrm{cos}^{n_0}\beta `$
$`𝒜^{\{𝒩\}}\stackrel{~}{P}_{\{n_n,n_p,n_{},n_0,n_+\}}^{0,0};`$
where

$$=n_n+n_p+2n_{}+2n_+;$$
(26)
The permutation normalization factor is here:

$$𝒜^{\{𝒩\}}=\frac{(n_{}+n_0+n_+)!(n_n+n_p)!}{n_{}!n_0!n_+!n_n!n_p!};$$
(27)
The real nonzero values of the weight is obtained only when the argument of the exponent in Eq.(̇25) vanishes:

$$n_nn_p+2n_{}2n_+=0;$$
(28)
This is equivalent to the conservation of the third component of the isospin.

Novel behaviour is obtained only when one considers time evolution of the system. Presence of an exact symmetry means that the corresponding weight Eq. (12) is constant, here we consider the expression:

$$\frac{d}{dt}\overline{P_{\{n_n,n_p,n_{},n_0,n_+\}}^{0,0}}=0;$$
(29)
We note that an appropriate analytical continuation should be made. First we evaluate the integral appearing in Eq. (25):

$$\frac{1}{𝒜^{\{𝒩\}}}\overline{P_{\{n_n,n_p,n_{},n_0,n_+\}}^{0,0}}=(1)^{n_0}\underset{i=0}{\overset{n_0}{}}(2)^i\left(\genfrac{}{}{0pt}{}{n_0}{i}\right)\frac{1}{+1+i};$$
(30)
This discrete form is not allowing an analytic continuation which would allow for all necessary differentiations. However, we can write the integral also as the hypergeometric $`{}_{2}{}^{}F_{1}^{}`$ function :

$$\frac{1}{𝒜^{\{𝒩\}}}\overline{P_{\{n_n,n_p,n_{},n_0,n_+\}}^{0,0}}=(1)^{n_0}\frac{1}{+1}{}_{2}{}^{}F_{1}^{}(n_0,+1,+2;2);$$
(31)
where $``$ is as defined in Eq. (26). We so obtain:

$`\overline{P_{\{n_n,n_p,n_{},n_0,n_+\}}^{0,0}}`$
$`={\displaystyle \frac{\mathrm{\Gamma }(n_{}+n_0+n_++1)\mathrm{\Gamma }(n_n+n_p+1)}{\mathrm{\Gamma }(n_{}+1)\mathrm{\Gamma }(n_0+1)\mathrm{\Gamma }(n_++1)\mathrm{\Gamma }(n_n+1)\mathrm{\Gamma }(n_p+1)}}{\displaystyle \frac{\mathrm{cos}n_0\pi }{\mathrm{\Gamma }(n_0)}}{\displaystyle \underset{i=0}{\overset{\mathrm{}}{}}}{\displaystyle \frac{\mathrm{\Gamma }(n_0+i)}{+1+i}}{\displaystyle \frac{2^i}{i!}};`$ (32)
Eqs. (25), (29) and (4), together with the condition (28) result in:

$`0`$ $`=`$ $`2{\displaystyle \frac{\stackrel{~}{P}_{\{n_n,n_p,n_{},n_0,n_+\}}^{0,0}}{}}\left({\displaystyle \frac{dn_n}{dt}}+2{\displaystyle \frac{dn_{}}{dt}}\right)+{\displaystyle \frac{\stackrel{~}{P}_{\{n_n,n_p,n_{},n_0,n_+\}}^{0,0}}{n_0}}{\displaystyle \frac{dn_0}{dt}}`$ (33)
$`+`$ $`{\displaystyle \frac{d\mathrm{log}𝒜^{(𝒩)}}{dt}}\stackrel{~}{P}_{\{n_n,n_p,n_{},n_0,n_+\}}^{0,0};`$
The projection integrals determine the coefficients which are, explicitly:

$`{\displaystyle \frac{\stackrel{~}{P}_{\{n_n,n_p,n_{},n_0,n_+\}}^{0,0}}{n_0}}`$ $`=`$ $`\pi {\displaystyle \frac{\mathrm{sin}n_0\pi }{\mathrm{\Gamma }(n_0)}}{\displaystyle \underset{i=0}{\overset{\mathrm{}}{}}}{\displaystyle \frac{\mathrm{\Gamma }(n_0+i)}{+1+i}}{\displaystyle \frac{2^i}{i!}}`$
$`+{\displaystyle \frac{\mathrm{cos}n_0\pi }{\mathrm{\Gamma }(n_0)}}{\displaystyle \underset{i=0}{\overset{\mathrm{}}{}}}{\displaystyle \frac{\mathrm{\Gamma }(n_0+i)[\psi (n_0)\psi (n_0+i)]}{+1+i}}{\displaystyle \frac{2^i}{i!}}`$
$`=`$ $`(1)^{n_0}{\displaystyle \underset{i=1}{\overset{n_0}{}}}\left({\displaystyle \genfrac{}{}{0pt}{}{n_0}{i}}\right){\displaystyle \frac{[\psi (1+n_0)\psi (1+n_0i)]}{+1+i}}(2)^i;`$
$`{\displaystyle \frac{\stackrel{~}{P}_{\{n_n,n_p,n_{},n_0,n_+\}}^{0,0}}{}}`$ $`=`$ $`(1)^{n_0+1}{\displaystyle \underset{i=0}{\overset{n_0}{}}}(2)^i\left({\displaystyle \genfrac{}{}{0pt}{}{n_0}{i}}\right){\displaystyle \frac{1}{(+1+i)^2}};`$ (35)
and

$`{\displaystyle \frac{d\mathrm{log}𝒜^{(𝒩)}}{dt}}`$ $`=`$ $`{\displaystyle \frac{dn_n}{dt}}\left[\psi (n_N+1)\psi (n_n+1)\right]+{\displaystyle \frac{dn_p}{dt}}\left[\psi (n_N+1)\psi (n_p+1)\right]`$ (36)
$`+`$ $`{\displaystyle \frac{dn_{}}{dt}}\left[\psi (n_\pi +1)\psi (n_{}+1)\right]+{\displaystyle \frac{dn_0}{dt}}\left[\psi (n_\pi +1)\psi (n_0+1)\right]`$
$`+`$ $`{\displaystyle \frac{dn_+}{dt}}\left[\psi (n_\pi +1)\psi (n_++1)\right];`$
where $`n_N=n_n+n_p`$ is the total number of nucleons and $`n_\pi =n_{}+n_0+n_+`$ is the total number of pions. Eq. (36) can be also written in the form:

$`{\displaystyle \frac{d\mathrm{log}𝒜^{(𝒩)}}{dt}}`$ $`=`$ $`{\displaystyle \frac{dn_n}{dt}}{\displaystyle \underset{k=1+n_n}{\overset{n_N}{}}}{\displaystyle \frac{1}{k}}+{\displaystyle \frac{dn_p}{dt}}{\displaystyle \underset{k=1+n_p}{\overset{n_N}{}}}{\displaystyle \frac{1}{k}}`$ (37)
$`+`$ $`{\displaystyle \frac{dn_{}}{dt}}{\displaystyle \underset{k=1+n_{}}{\overset{n_\pi }{}}}{\displaystyle \frac{1}{k}}+{\displaystyle \frac{dn_0}{dt}}{\displaystyle \underset{k=1+n_0}{\overset{n_\pi }{}}}{\displaystyle \frac{1}{k}}+{\displaystyle \frac{dn_+}{dt}}{\displaystyle \underset{k=1+n_+}{\overset{n_\pi }{}}}{\displaystyle \frac{1}{k}};`$
Eqs. (3337) offer the final result for the $`SU(2)`$ case . Notably, they imply a relation for the number of neutral pions in a system. We thus see, that when the case of the non-abelian symmetry is carefully considered and not ignored, one can get relations determining also the ‘neutral’ members of multiplets. In the standard approach the multiplicity of neutral particles is obtained by introducing a subsidiary chemical potential which is related to the lack of chemical equilibrium of a system  or to the residual interaction with the environment .

## 5 Conclusions and Outlook

We have shown how constraints due to the preservation of the symmetry properties of a multiparticle macroscopic state define the path of evolution of the system. One should notice that results we presented are general and do not depend on the particular choice of the representation of the symmetry group. For different initial symmetry group representations one gets different paths, but they are all of a similar “shape”, considering the hypothesis that the global behaviour of a macrosystem (in a sense of statistical physics) should not be altered if a number of particles is changed by a very small (“microscopic”) amount.

We have explicitly presented the example how our constraint works in the simplest non - trivial case of $`SU(2)`$ symmetry group.

Although we have studied and implemented the discreet symmetry using quantum states, whenever we referred here to a dynamical equation we considered the limit of incoherent state evolution described by the Boltzmann equation. The dynamical evolution we consider thus is described in terms of diagonal density matrix. This is the appropriate approach given that our main objective is to arrive at a dynamical derivation of symmetry deformed statistical distribution.

We recall that quantum correlations (without symmetry) alone are responsible for the deformation of the Boltzmann distribution into Bose/Fermi distributions, and that the Boltzmann equation yields this result when we allow for Fermi blocking/Bose enhancement in the collision term. In that line of thought, the next step would be to show that it is possible to obtain now within a dynamical Boltzmann equation calculation the evolution of a many body system into symmetry-deformed statistical equilibrium distribution. We are also exploring the possibility that the methods here presented allow the formulation of a microscopic transport theory which would obey the long range correlations introduced by the macroscopic quantum and symmetry constraints.

## Acknowledgments

Work supported in part by a grant from the U.S. Department of Energy, DE-FG03-95ER40937 , and in part by the Polish Committee for Scientific Research under contract KBN - 2 P03B 030 18 .