Computing Quark and Gluon Distribution Functions for Very Large Nuclei

TL;DR

This paper demonstrates that quark and gluon distribution functions for very large nuclei can be computed using weak coupling methods, recasting the problem as a many-body calculation with modified propagators and media-dependent vertices.

Contribution

It introduces a novel approach to compute distribution functions for large nuclei at small x using weak coupling and many-body techniques, including explicit lowest-order calculations.

Findings

Computed gluon distribution function at lowest order.

Proposed methods for higher-order calculations.

Argued the feasibility of weak coupling approaches for large nuclei.

Abstract

We argue that the distribution functions for quarks and gluons are computable at small x for sufficiently large nuclei, perhaps larger than can be physically realized. For such nuclei, we argue that weak coupling methods may be used. We show that the computation of the distribution functions can be recast as a many body problem with a modified propagator, a coupling constant which depends on the multiplicity of particles per unit rapidity per unit area, and for non-abelian gauge theories, some extra media dependent vertices. We explicitly compute the distribution function for gluons to lowest order, and argue how they may be computed in higher order.

Equations

$$\rho =\frac{1}{\pi R^2}\frac{dN}{dy}$$

$${\Lambda }_{QCD}^2<<{\alpha }_s\rho <<k_t^2<<\rho$$

$$\frac{1}{\pi R^2}\frac{dN}{dx{d}^2{k}_t}\sim \frac{{\mu }^2}{x{k}_t^2}$$

$$S=\sum _i{\overline{\psi }}_i(/P+M_i){\psi }_i+\frac{1}{4}{F}_{\mu \nu }{F}^{\mu \nu }-J\cdot A$$

$$F_a^{\mu \nu }={\partial }^{\mu }{A}_a^{\nu }-{\partial }^{\nu }{A}_a^{\mu }-g{f}_{abc}{A}_b^{\mu }{A}_c^{\nu },$$

$$P=\frac{1}{i}\partial -gA,$$

$$A^{\mu }=A_a^{\mu }{\tau }^a$$

$$a_{\pm }=(a_0\pm a_3)/\sqrt{2}$$

$${\alpha }^{\pm }=\frac{1}{\sqrt{2}}{\gamma }^0{\gamma }^{\pm }$$

$${\psi }_{\pm }={\alpha }^{\pm }\psi$$

$${\alpha }^{\pm }{\psi }_{\pm }={\psi }_{\pm }$$

$${\alpha }^{\pm }{\psi }_{\mp }=0$$

$$\psi ={\psi }_{+}+{\psi }_{-}$$

$$\sqrt{2}{P}_{-}{\psi }_{-}=-{\gamma }^0(/P_t+M){\psi }_{+}$$

$$A_{-}=-A^{+}=0$$

$${\psi }_{-}=\frac{1}{\sqrt{2}{P}^{+}}{\gamma }^0(/P_t+M){\psi }_{+}$$

$$S_F=-{\psi }_{+}^{†}{P}^{-}{\psi }_{+}+\frac{1}{2}{\psi }_{+}^{†}(M-/P_t)\frac{1}{P^{+}}(M+/P_t){\psi }_{+}$$

$$F^2=F_t^2-4F_{k+}{F}_{k-}+2F_{+-}{F}_{+-}$$

$$F_{+-}={\partial }_{+}{A}_{-}-{\partial }_{-}{A}_{+}-ig[A_{-},A_{+}]=-{\partial }_{-}{A}_{+},$$

$$F_{k+}={\partial }_k{A}_{+}-{\partial }_{+}{A}_k-ig[A_k,A_{+}],$$

$$F_{k-}=E_k=-{\partial }_{-}{A}_k$$

$$D_{\mu }{F}^{\mu \nu }=J^{\nu }$$

$$-{\partial }_{-}^2{A}^{-}=J_F^{+}+D_k{E}^k$$

$$J_{Fa}^{\mu }=\overline{\psi }{\gamma }^{\mu }{\tau }_a\psi$$

$$D_k{E}^k$$

$$\frac{1}{4}{F}^2=\frac{1}{4}{F}_t^2-D_k{E}^k{A}_{+}+\frac{1}{2}({\partial }_{-}{A}_{+}{)}^2-{\partial }_{-}{A}_k{\partial }_{+}{A}^k$$

$$S$$

$$P^{-}$$

$${\psi }_{\alpha }(x)={\int }_{k^{+}>0}\frac{d^{3}k}{(2\pi {)}^3}\left(b_{\alpha }(k)e^{ikx}+d_{\alpha }^{†}(k)e^{-ikx}\right)$$

$$A_i^a(x)={\int }_{k^{+}>0}\frac{d^{3}k}{(2\pi {)}^3\sqrt{2k^{+}}}\left(a_i^a(k)e^{ikx}+a_i^{a†}(k)e^{-ikx}\right)$$

Full text

TPI–MINN–93–44/T

NUC–MINN–93–24/T

HEP–UMN–TH–1220/93

August 1993

Computing Quark and Gluon Distribution Functions for Very Large Nuclei

Larry McLerran and Raju Venugopalan

School of Physics and Astronomy, University of Minnesota, Minneapolis, MN 55455

**Abstract

**

We argue that the distribution functions for quarks and gluons are computable at small x for sufficiently large nuclei, perhaps larger than can be physically realized. For such nuclei, we argue that weak coupling methods may be used. We show that the computation of the distribution functions can be recast as a many body problem with a modified propagator, a coupling constant which depends on the multiplicity of particles per unit rapidity per unit area, and for non-abelian gauge theories, some extra media dependent vertices. We explicitly compute the distribution function for gluons to lowest order, and argue how they may be computed in higher order.

1 Introduction

The problem of computing the distribution functions for quarks and gluons in hadrons is an old and difficult problem. Issues such as Bjorken scaling were greatly clarified by using light cone Hamiltonian methods [1] - [2]. There has been much progress recently in applying light cone Hamiltonian methods together with non-perturbative methods gleaned from lattice gauge theory technology to compute these distribution functions [3]. It has nevertheless been believed that the computation of these distribution functions is intrinsically non-perturbative.

On the other hand, for a very large nucleus or at very small Bjorken x, it is known that the density of quarks and gluons per unit area per unit rapidity

[TABLE]

is large. For nuclei, we expect that this density scales as $A^{1/3}$, or perhaps some larger power of $A$, so that for some sufficiently large $A$, $\rho>>\Lambda^{2}_{QCD}$. Even if $\rho$ is not sufficiently large for realistic nuclei, we certainly can imagine nuclei where this would be true. In any case, when $\rho>>\Lambda^{2}_{QCD}$, we expect that the coupling constant for strong interactions is weak, and weak coupling methods should be valid for computing the distribution functions. [4]

In this paper, we will show how to set up the problem of computing the distribution functions when the parton density $\rho$ is large. We will explicitly compute the lowest order gluon distribution function at small Bjorken $x$ for transverse momenta

[TABLE]

and find that in this kinematic region, the gluon distribution functions are of the Weiszacker-Williams form

[TABLE]

where $\mu$ is a parameter which we shall compute which behaves as $A^{1/3}$. In this theory, the dimensionful scale factor $\mu$ will set the scale of the coupling constant. All perturbation theory can be done in terms of $\alpha(\mu)$, and if $\alpha(\mu)<<1$ a weak coupling expansion is valid. This is equivalent to $\rho\sim\mu^{2}>>\Lambda_{QCD}$.

We will also argue that the quark distribution functions are computable in this kinematic region, and outline how to do the lowest order computation. It may be possible to extend the region of validity for computation of the distribution functions to smaller values of $k_{t}$ by including non-perturbative effects computable in weak coupling. We also argue that the power dependence of the distribution functions in Bjorken $x$ may be modified in higher orders of perturbation theory.

The outline of this paper is as follows: In the second section, we present a brief review of the light cone quantization method, and our notation. We also argue that the valence quarks inside the nucleus may be replaced by static quarks propagating along the light cone. We set up the formalism which allows us to compute the ground state properties of a large nucleus. In the third section, we do a simple analysis for QED and compute the infinite momentum frame wavefunction for an electron to lowest non-trivial order in weak coupling. We show that the photon distribution function is simply the Weiszacker-Williams distribution of Lorentz boosted Coulomb photons. In the fourth section, we do the problem analogous to QED for large nuclei. We argue that for large enough nuclei for $k_{t}^{2}<<\rho$ where $\rho$ was defined above, it is a good approximation to treat the local density of color charge classically. We argue that the local fluctuations in color charge may be integrated out of the problem for computing ground state properties of the system. After performing this integration, we generate the distribution functions for quarks and gluons. In lowest order in perturbation theory, we compute the distribution functions of gluons and show that they are of the Weiszacker-Williams form. We also show that there is an effective field theory with propagators modified from their vacuum form which allows for the systematic computation of corrections to the lowest order result. The problem of computing the distribution functions has therefore been converted into a many body problem which for sufficiently large nuclei is a weak coupling problem. In the summary, we outline the computation of the quark distribution functions. We also show how non-perturbative, although presumably weak coupling, effects modify the distribution functions at low $k_{t}$, $\Lambda^{2}_{QCD}<<k_{t}^{2}<<\alpha_{s}\rho$ We also argue that the $x_{Bj}$ dependence of the structure functions might be modified in higher orders in perturbation theory.

2 Notation and the Infinite Momentum Frame Hamiltonian

Before turning to a detailed computation, we first review the infinite momentum frame Hamiltonian method. We will also show how to include the effects of valence quarks which we will treat as static sources moving along the light cone.

We begin with the action for the Yang-Mills field coupled to fermions and also possibly some external current $J^{\mu}$,

[TABLE]

Here the Yang-Mills field strength is

[TABLE]

the covariant momentum is

[TABLE]

where we use the notation that

[TABLE]

We are using the metric where $a\cdot b=a^{i}b^{i}-a^{0}b^{0}$, but with the ordinary convention for gamma matrices.

This Hamiltonian may be re-expressed in terms of light cone variables by the identification of components of vectors as

[TABLE]

In these variables, the metric is $g_{+-}=-1=g_{-+}$, the transverse components of the metric are the unit matrix, and all other components vanish.

To quantize along the light cone, as usual we first express the fermion fields in terms of their light cone spinor components by introducing the projection operators

[TABLE]

These are Hermitian projection operators so that we can write

[TABLE]

where

[TABLE]

and

[TABLE]

so that

[TABLE]

We are going to be interested in the light cone Hamiltonian $P^{-}$ which generates displacements in $x^{+}$. When we write out the Euler-Lagrange equations for the fermion fields, one of the equations is a constraint equation for the fields on a surface of fixed $x^{+}$,

[TABLE]

This equation can be explicitly inverted in the light cone gauge

[TABLE]

We find

[TABLE]

The fermion contribution to the action is therefore

[TABLE]

where we have rescaled $\psi\rightarrow{1\over 2^{1/4}}\psi$. In terms of these variables, we see that $\psi_{+}^{\dagger}$ is the light cone momentum canonically conjugate to $\psi_{+}$.

To analyze the vector contribution to the action, we first write explicitly

[TABLE]

In light cone gauge, we have

[TABLE]

[TABLE]

and

[TABLE]

The equations of motion for the vector field are

[TABLE]

In particular, the equation for the $+$ component of the current is a constraint equation for $A^{-}$ on a fixed $x^{+}$ surface,

[TABLE]

where $J_{F}$ is the current generated by the fermion field

[TABLE]

plus whatever external currents there are in the system. The piece of the right hand side of this equation which is

[TABLE]

can be thought of as the bosonic contribution to the $+$ component of the charge density to which $A^{-}$ couples.

We can now write the bosonic part of the action as

[TABLE]

We see that the momentum canonically conjugate to the field $A_{k}$ is the momentum $\Pi_{k}=-E_{k}$ where the index $k$ runs only over transverse coordinates. We can then express the action in terms of momenta and coordinates as

[TABLE]

We therefore see that the generator of $x^{+}$ transformations is

[TABLE]

With the above identification of canonical momenta and coordinates, we see that we can quantize the fermion and vector fields as

[TABLE]

and

[TABLE]

The commutation relations for the operators $a,b$ and $d$ are

[TABLE]

and

[TABLE]

with all other (anti)commutators vanishing.

The above review of the light cone quantization procedure will standardize the notation for the following analysis. We have seen that it is possible to explicitly eliminate the constraints of gauge fixing in the light cone gauge and get a light cone Hamiltonian expressed in terms of the true dynamical degrees of freedom of the system.

Finally, we will discuss an essential part of our formalism: the reduction of the valence quarks for the infinite momentum frame wavefunction to static external sources of charge. To do this, we first make the assumption that in the end we will be applying our analysis to weakly coupled systems. In a theory like QED, this can be done automatically since in this theory the coupling constant is always weak, except possibly at very short distances. In QCD, we will have to restrict our attention to systems where the density of partons per unit area is very large.

If the coupling is weak, then the dominant mechanism for producing a cloud of partons around a fast moving valence particle is by bremstrahlung or chains of particles bremstrahlunging from bremstrahlung particles. In this case, for weak coupling, the typical momentum being transferred in the bremstrahlung process is soft, and the valence parton does not lose a large fraction of its momenta. If this is the case, then the valence parton typically moves with a velocity close to light velocity with only very small transverse components generated by coupling to multiparticle degrees of freedom. The valence parton is therefore a recoiless source of charge moving along the lightcone. It is well known that this approximation describes the soft photon dressing of the electron in QED, and has been used to study the infrared region of QED.

In QCD, we are claiming that the approximation should be valid for describing the parton distributions generated by the valence parton whenever the density of partons is sufficiently high that weak coupling methods can be used. This will however only work for describing the parton cloud generated at small $x_{F}$. In the fragmentation region, we are sensitive to details of the spatial distribution of quarks as generated by the nuclear wavefunction. For example, if we look at the distribution of nucleons in the rest frame of a nucleus, we see the extended nuclear structure. If we however look in the central region, the valence nucleons will appear to be Lorentz contracted to a pancake of dimensions of order $R_{nuc}/\gamma_{cm}$ where $R_{nuc}$ is the nuclear radius, and $\gamma_{cm}$ is the gamma factor for the nucleus as measured in the center of mass frame. In the central region, we expect that the parton distributions will be insensitive to the details of the distribution of partons in the fragmentation region, that is, the valence partons will act simply as a source of charge which is essentially a delta function along the light cone.

In addition, maintaining the constraint that the momentum transfer be small and the coupling weak is a little tricky in QCD. We will require that the momentum transfer be small compared to the total momentum of the valence partons, but large compared to the QCD scale. For this to be consistent, the low momentum range of integration must be cutoff for small momentum transfer for the range of dynamically important momenta. Presumably this happens due to media effects caused by the high density of partons.

Although the arguments given above are heuristic, we will find that the solution we generate for the quark and gluon distribution functions in the central region are self-consistent with the assumption that the valence partons may be treated as a delta function source of charge along the light cone. Perhaps it might be possible at some point to relax the assumptions about the fragmentation region distributions, but we do not at this point know how to do it.

We therefore will treat the valence parton distributions as recoiless sources of charge which are localized along the light cone, that is, we take

[TABLE]

where we have assumed the source is localized at $x^{-}=0$ and $\vec{x_{t}}=0$. (We could have placed the sources anywhere.) The charge $Q_{a}$ is an operator which has the charge algebra

[TABLE]

and is in some representation of the Lie group corresponding to the structure constant $f_{abc}$. If the representation of the Lie algebra is sufficiently large, we may be able to replace the charge operator $Q_{a}$ by a c-number classical source.

The $x^{+}$ dependence of the source $Q_{a}$ is determined by a covariant conservation law

[TABLE]

This equation is necessary so that the equations of motion for the vector fields are consistent, that is we must have

[TABLE]

This equation would be immediately true if we took the current for dynamical fermion fields. We can take the static limit for fermion fields and construct the current, and the above equation will therefore be valid.

We can solve the equation for the time dependence of the charge operator. If we let

[TABLE]

Further let the time ordered exponential be defined as

[TABLE]

where we are not writing out the dependence of $U$ upon $x^{-}$ and $x_{t}$. We therefore have that

[TABLE]

Note that $U$ is a unitary matrix and that

[TABLE]

so that the time evolution of the charge operator is just a rotation of the charge operator in charge space.

We emphasize here that the above equation for the time evolution of the charge operator is true for classical as well as quantum charge operators. For abelian theories, the charge is simply time independent. In non-abelian theories, the charge is specified by its initial value, but after this time rotates in charge space.

3 The Photon Distribution Function of an Electron

The discussion of the previous section has been highly formal, and largely a brief review of what is already known in the literature. The introduction of static sources along the light cone in the form presented above is perhaps less well known, and it is useful to see how this formalism works in the case of QED. Here we will compute the wavefunction for an electron moving along the lightcone, and the photon distribution function of the electron. It will not be too surprising that we find that the distribution function is precisely that for the Weiszacker-Williams distribution of photons generated by Lorentz boosting the Coulomb field.

We will work to lowest non-trivial order in $\alpha$ so that we can ignore $e^{+}e^{-}$ pair production. In this approximation, we ignore dynamical fermions in the light cone Hamiltonian. The source term for the valence electron is

[TABLE]

The light cone Hamiltonian is

[TABLE]

The ground state for this light cone Hamiltonian is a coherent state

[TABLE]

In this expression, the quantity $C$ is a normalization constant, the field $A^{op}$ is the operator value of the transverse component of the vector potential, and $E^{cl}$ is a c-number valued classical electric field.

The reason ground state is a coherent state is easily understood. In the Weiszacker-Williams approximation, the photon cloud is a classical distribution of electromagnetic fields generated by Lorentz boosting the Coulomb field. It is well known that coherent states describe classical fields.

The expectation value of the light cone Hamiltonian in this state is simply

[TABLE]

where we have used that the exponential in the coherent state simply shifts the electric field. Extremizing the Hamiltonian with respect to the classical electric field gives

[TABLE]

so that the classical electric field may be taken as purely longitudinal in the two dimensional transverse space. This is consistent with ignoring the contribution from $B^{2}$ which arises from transverse components. (A more careful treatment would have involved writing the Hamiltonian in terms of transverse and longitudinal creation and annihilation operators, constructing a coherent state in terms of transverse and longitudinal fields, and showing that the minimum energy configuration has zero transverse vector potential.)

An amusing feature of the QED analysis is that the eigenvalue of the light cone Hamiltonian for the ground state vanishes. This is because only a longitudinal field arises and hence there is no contribution from the transverse magnetic field. The light cone charge of the source is also identically canceled by the light cone charge generated by the transverse electric field.

We therefore have the classical electric field as

[TABLE]

and the vector potential as

[TABLE]

where our convention for Fourier transformation is

[TABLE]

Notice that the transverse electric field which we have generated is precisely the Weiszacker-Williams field generated by Lorentz boosting the Coulomb field. We could further reproduce the magnetic components of the field if we were to compute $A_{+}$, and its associated field strengths.

The photon distribution function is given by

[TABLE]

where $F(k^{+},k_{t})=dN/dk^{+}d^{2}k_{t}$ This can also be rewritten in terms of the Bjorken x variable as

[TABLE]

The re-expression in terms of Bjorken x variables may here seem a bit peculiar since nowhere has the longitudinal momentum of the source been introduced. It can be done because of the scaling property of the distribution function, peculiar to the $1/k^{+}$ behavior. If there was any other power of x, this could not be done. We might therefore ask what would happen in the QCD case or perhaps in QED if electron-positron pair creation was included. Presumably the distribution function becomes modified to be something like

[TABLE]

so that higher orders generate terms which involve logarithms of x. These logarithms of $k^{+}/P_{source}^{+}$ presumably arise in our formalism by a sensitivity to a high momentum cutoff necessary to regularize the $\delta$-function distributions for the external sources. Note that for a delta-function source, there is a contribution from arbitrarily large momenta. This is presumably cutoff at momenta of order that of the valence particle when recoil is included. If this is the case, then we expect that we will be able to compute the dependence of the distribution functions on $x$ and $k_{t}$, up to some overall constant. This constant cannot be computed without a better knowledge of the dynamics of the fragmentation region, that is understanding some of the relevant details of the valence particle recoil. In any case, it is remarkable that, to lowest order, the result is insensitive to such details.

4 Ground State Properties in QCD in the Presence

of Sources

In QCD, we shall be interested in computing the distribution functions for quarks and gluons generated by some valence distribution of quarks. In order to be able to use weak coupling techniques, we will have to require that the density of partons per unit area

[TABLE]

is large

[TABLE]

The simplest problem to consider is that of the distribution functions for a very large nucleus. Although our results could be generalized to a finite size nucleus, we will consider infinite nuclear matter in the transverse direction with a uniform transverse density distribution. We will assume that the nucleus is thick compared to a proton, that is that the local density of baryon number per square fermi is large. Roughly speaking, the density of baryon number per square fermi should be of order $A^{1/3}$ for an ordinary nucleus, and the transverse extent should also be of order $A^{1/3}$. We could also presumably take the distribution functions we compute here as a function of local baryon density per unit area and convolute them over a density distribution for a nucleus in order to determine realistic parton densities. For realistic nuclei, our weak coupling approximations are quite probably at the edge of being valid.

We shall therefore imagine that the baryons which generate the nuclear valence distribution are localized on the light cone but uniform in transverse space. Their charge will however in general be fluctuating around zero from one transverse position in space to the next. We expect that the average color charge associated with the valence distributions will be zero, and therefore the only way to generate a non-zero source for charge is by fluctuations. We will treat the nuclear valence quarks as static, that is recoiless, sources of charge.

To compute ground state properties of our ensemble of charges, we consider the ground state expectation value

[TABLE]

where $T$ is a parameter which will be taken to be large in the end. The state $\mid 0>$ is the ground state in the presence of a uniform density of sources of color charge.

If we take the limit that $T\rightarrow\infty$, we can generate ground state expectation values by

[TABLE]

The problem with evaluating the trace is that for quantized sources of color charge, the sum over different values of the color charges is difficult to evaluate. However for a large nucleus, so long as we resolve the system on a transverse size scale which is much larger than a typical transverse quark separation, we expect that many different quark charges will contribute to the valence charge density. A large number of charges corresponds to a high dimensional representation of the color algebra, and therefore the sum of the color charges of the quarks can be treated classically.

To be more specific, we look at a region of transverse extent where there are a large number of quarks. The total charge in this region on the average will be zero. There will however be fluctuations in the color charge. We would naively expect that the typical value of the fluctuating color charge would be of order $\sqrt{N}$ where $N$ is the number of quarks in this region [5]–[7]. If $N$ is large, the typical color charge in this region is large, and therefore may be treated classically.

This may be analyzed in a quantitative manner using our formula for $Z$ which generates ground state expectation values. The point is that in the trace, there is the sum over all possible color orientations of the external source corresponding to the valence quark distributions. Let us make a grid in the transverse space so that each grid has many valence quarks in it. This therefore restricts the validity of our analysis to spatial resolutions with

[TABLE]

or

[TABLE]

For realistic large nuclei, this corresponds to $q_{t}\sim 1-2$ GeV (recalling that in the center of the nucleus the effective value of $A^{1/3}$ is larger than average), and in this kinematic region, weak coupling methods are at best marginal. In principle we can imagine nuclei with very large $A$, so we can consider this as a theoretical laboratory. The situation might also be improved if the density of partons per unit area in the central region were to greatly exceed this and that the coupling constant were smaller than naively expected from the above consideration.

Now in our sum over states in the formula for $Z$, there is a large number of states being generated for each transverse area. We would like to determine the set of most probable configurations. To do this, we need to construct the density of states for such configurations. We know that the maximum will be centered around an average color charge of zero. Since it is a maximum, it must be true that it is quadratic in the charge density. Therefore the contribution to Z is of the form

[TABLE]

Here $\rho^{2}$ is the charge density squared measured in units of $g$, that is with the factor of $g$ extracted, and $\mu^{2}$ is the average charge density squared per unit area divided by $g$.

The above formula for the density of states is valid so long as the average charge is of the order of the $\sqrt{N_{0}}$. In this case, higher order corrections to the Gaussian approximation for the density of states go like $1/\sqrt{N_{0}}$. Notice that with the assumption that the number of particles in our box is of order $N_{0}$, we have that $d^{2}x\rho^{2}\sim 1$, and higher order terms are correspondingly smaller in powers of $N_{0}$ since $\rho\sim\rho_{Q}/\sqrt{N_{0}}$ where $\rho_{Q}$ is the density of valence quarks per unit area.

To see how such a density of charge works, we see that if we define

[TABLE]

then we have

[TABLE]

as we have assumed and

[TABLE]

We can determine the value of $\mu^{2}$ from elementary considerations. The number of quarks per unit area in a nucleus is

[TABLE]

and the average charge squared of a quark is

[TABLE]

so that

[TABLE]

Now that we understand the density of charge states associated with the external charge, we proceed to a path integral representation for Z. This is easily done in the standard way by introducing states which are momentum and coordinate basis states for the fields. The only tricky part is to recall that in integrating over the external charges, we must remember that the external charges must satisfy the extended current conservation law discussed in the second section. The result of some analysis is that in light cone gauge

[TABLE]

In the functional integral, we have that the path integral over the fields is at all space-time points. The path integral over $\rho$ is at only a fixed value of $x^{+}$ which we choose to be $x^{+}=0$, and we define $\rho$ at other values of $x^{+}$ by

[TABLE]

where $U$ is the time ordered exponential of $A_{+}$ as discussed in the second section.

To get an action expressed only in terms of the vector fields, we must perform the functional integral over fluctuations in the external charges. This is straightforward to do because the action is quadratic in the external charges. The part of the action involving $\rho$ is

[TABLE]

In this expression, the quantity $M^{ab}$ is

[TABLE]

Notice that in this expression there is non-locality in the time variable $x^{+}$, but locality in $x^{-}$ and $x_{t}$. The matrix $M$ is Hermitian since

[TABLE]

The integral over $\rho$ can be performed since the integral is Gaussian. The value of $\rho$ is given by

[TABLE]

Upon integrating out the sources, we obtain that the contribution of $S_{\rho}$ becomes

[TABLE]

In QED, the effect of integrating out the sources is simply to give

[TABLE]

that is, only the propagator is modified by the addition of an extra term, which is formally of order $\mu^{2}g^{2}$. In QCD, the problem is more complicated. In addition to the modification of the propagator in order $g^{2}$, there are in addition a non-polynomial and non-local set of vertices which contain arbitrary powers of the field $A_{+}$. These vertices will modify the Feynman rules in orders beyond the lowest order formally in $g^{2}$.

We must ask what is the expansion parameter of the theory. Later in this section, we will show that the distribution function for the number of gluons per unit area is given in lowest order by

[TABLE]

In general, when we are in the limit that the relevant momentum scales are

[TABLE]

we are in the low density small $g^{2}\mu^{2}$ limit so that the effects of the modified gluon propagator may be treated perturbatively. The other end of the limit, that $k^{2}<<\mu^{2}$ is just the statement that we are at sufficiently small momenta that the individual quarks which are the source of color charge cannot be individually resolved.

In the above kinematic limit, a treatment of our effective Lagrangian to lowest order in weak coupling and lowest order in $\mu^{2}$ is justified. The higher orders in $g^{2}$ should be small so long as $g^{2}$ is small due to the high density of partons. This may be only part of the story however. As is known in QCD at finite temperature, it may be true that at some order in perturbation theory for the computation of physical observables, the weak coupling expansion becomes essentially non-perturbative. This could happen here, and to know if it happens and what order it affects observables is beyond the scope of the present work.

In addition to the above problems, the constraints of the renormalization group equations must be satisfied. This will give a non-trivial dependence on coupling to the structure functions. Beyond this there may also be singularities which arise from approximating the source of color charge as recoiless. For example, if the structure functions end up having a power of Bjorken $x$ other than $1/x$, then a scale associated with the longitudinal momentum of the sources will have to come into the problem. This will arise from regularizing the recoiless charge distribution which for a delta function involves integrating over all longitudinal momentum. If we cut this off at some value of $P^{+}$, then if the corrections to the lowest order result for the distribution functions involve logarithms of this cutoff scale, they can sum up to something like

[TABLE]

Thus although there is only logarithmic dependence on the cutoff this can sum up to a power law behavior for the structure functions. If this is the case, then it might be that although the momentum dependence of the structure functions can be determined, their overall scale is not possible to evaluate. However, in weak coupling the effect of the overall scale factor which cannot be computed is of order $A^{\alpha}\sim 1+\alpha ln(A)$, and so long as all the kinematically singular terms have been included, the undetermined coefficient $A$ which is of order 1 and depends on the details of the cutoff only makes a weak coupling perturbative modification of the structure function.

Higher orders in $\mu^{2}$, the infrared structure of the theory is probed. There are several possibilities which might arise here. It might turn out that the infrared singularities of the theory are entirely screened. In this case, the infrared properties of all correlations functions can be computed self-consistently in weak coupling in the infrared. A more likely scenario is however that as in finite temperature QCD, one some of the screening lengths are essentially non-perturbative. If this is the case, we can compute structure functions only in specific kinematic regions, and the behavior of at least some of the structure functions in the infrared can be determined only non-perturbatively.

To determine the Feynman rules of the theory, we must evaluate the propagator. To zeroth order in $\mu^{2}$, the propagator is simply the light cone gauge propagator. This is

[TABLE]

where $n$ is the unit vector $n^{\mu}=\delta^{\mu-}$

To find the first order correction to the propagator arising from $\mu^{2}$ in our effective Lagrangian it is perhaps most easy to compute the expectation value of

[TABLE]

by computing with the source present and then integrating out the source. We see that this term is given in the presence of the sources to lowest order as

[TABLE]

Upon integrating over the sources, we have therefore that

[TABLE]

The change in the propagator has therefore been replaced by a product of two classical fields times an overall Kronecker delta function for color invariance and a two dimensional delta function for invariance under transverse translations. There is also no factor of $i$ in the change of the propagator. These are contributions to the dispersive part of the propagator. This is always the case in many body theory, where the dispersive part represents occupied states.

It is now straightforward to compute the distribution function for gluons. This may be done by for example computing the expectation value of the light cone current density. It may also be read off directly from the above equation in analogy to the result for QED. The result is that

[TABLE]

This is just the Weiszacker-Williams distribution. The strength of the gluon distribution is proportional to the typical charge squared per unit area, and should go like $A^{1/3}$

5 Summary and Conclusions

In the previous sections, we discussed the potential problems with computing the distribution functions in higher order in weak coupling. Much must be done to turn this hypothetical program into a reality. Some of the problems to be solved are associated with the ultraviolet. There are two type of modifications here. The first are corrections which generate the proper scaling behavior of the distribution functions in order to be consistent with the Altarelli-Parisi equations. This is presumably straightforward to do. Following Brodsky, we must introduce a transverse momentum cutoff, and then require that the theory have the correct renormalization group improved dependence on this cutoff. A second problem is more serious, and that is how to get any non-trivial behavior in $x_{Bj}$. This presumably arises from a dependence on the cutoff which regularizes the delta function source for the static light cone sources of charge. The resulting non-trivial scaling dependence of the structure functions must in some way be related to the Lipatov Pomeron [8]. Needless to say, establishing whether this may be done is entirely non-trivial and the success of this theory of the structure functions hinges crucially on being able to establish this fact. In the end, if the small x behavior is enhanced, the usefulness of weak coupling expansion might even be extended to small $A$ targets. This would make the above formalism much more attractive.

The properties of the infrared structure of this theory are also entirely non-trivial. Is the infrared structure computable in weak coupling, or is it like finite temperature field theory where some aspect of the infrared structure of the theory are intrinsically non-perturbative? Is there a hierarchy of different infrared scales as there is in finite temperature theory where the thermal wavelength is much less than an electric screening length which in turn is much less than a non-perturbative magnetic screening length? Are there non-trivial non-perturbative phenomena which occur at these non-perturbative length scales?

Finally, there remains the problem of computing the quark distribution functions. To compute them in lowest order demands a one loop computation of the fermion propagator. We intend to do this computation in a later paper. Of interest here is the ratio of the quark sea distribution to the small x gluon distribution. This should be computable in the same kinematic region where the gluon distribution is computed.

Finally, there is the issue of how the distribution functions for quarks and gluons described above are related to the distribution functions measured in deep inelastic scattering. Although at first sight this seems trivial, recall that we are in the limit where there is a high density of partons, where one might expect screening effects to be important. Moreover, the dependence of the structure functions on the $Q^{2}$ of the probe must be established through use of the Altarelli-Parisi equations.

In all of our analysis in QCD, we never explicitly computed a hadronic wavefunction. Is there anyway in QCD to establish what is this wavefunction, or is it as difficult and in the end unrewarding task as determining what the wavefunction is for a large system in the microcanonical ensemble? In this latter case, it is sufficient to consider density matrices to study properties of the system. The actual wavefunction is almost never useful.

There is also the question of the relationship between the distribution functions we have computed and the early time behavior in heavy ion collisions. How does one let the distributions thermalize and evolve through a quark-gluon plasma?

6 Acknowledgments

The authors especially thank Al Mueller for his comments and help in reading this manuscript prior to being issued as a preprint, and for guiding us through the current state of knowledge on this problem. We wish to thank Stan Brodsky for his help in understanding the light cone gauge formalism and providing the authors with extremely clear pedagogical lectures. We thank Joe Kapusta for encouragement, and Mikhail Voloshin for critical comments. We acknowledge support under DOE High Energy DE-AC02-83ER40105 and DOE Nuclear DE-FG02-87ER-40328. Larry McLerran wishes to thank the Aspen Center for Physics where part of this work was finished, and Klaus Kinder-Geiger who in his talk at Quark Matter 93 implied that to understand the possible formation of a quark-gluon plasma in heavy ion collisions, one must first understand the initial conditions.