Beyond integrability: Baryon-baryon backward scattering in the massive Gross-Neveu model

TL;DR

This paper investigates how adding a mass term to the Gross-Neveu model breaks its integrability, leading to backward and inelastic baryon-baryon scattering, analyzed through a combination of analytical and numerical methods near the non-relativistic limit.

Contribution

It provides the first detailed study of baryon-baryon scattering in the massive Gross-Neveu model, demonstrating the emergence of backward scattering due to broken integrability.

Findings

Backward elastic scattering occurs when integrability is broken.

The study confirms inelastic reactions in the massive model.

Analytical and numerical methods successfully analyze the scattering processes.

Abstract

Due to integrability, baryon-baryon scattering in the massless Gross-Neveu model at large N features only forward elastic scattering. A bare mass term breaks integrability and is therefore expected to induce backward elastic scattering as well as inelastic reactions. We confirm these expectations by a study of baryon-baryon scattering in the massive Gross-Neveu model near the non-relativistic limit. This restriction enables us to solve the time-dependent Hartree-Fock equations with controlled approximations, using a combination of analytical methods from an effective field theory and the numerical solution of partial differential equations.

Equations

$$\mathcal{L}=\sum _{i=1}^N{\bar{\psi }}^{(i)}\left(i{\gamma }^{\mu }{\partial }_{\mu }-m_0\right){\psi }^{(i)}+\frac{1}{2}{g}^2{\left(\sum _{i=1}^N{\bar{\psi }}^{(i)}{\psi }^{(i)}\right)}^2.$$

$$\frac{\pi }{N{g}^2}=\gamma +\ln \frac{\Lambda }{m},\gamma :=\frac{\pi }{N{g}^2}\frac{m_0}{m}.$$

$$\frac{m[m_0]}{m[0]}=e^{\gamma }.$$

$$\left(i{\gamma }^{\mu }{\partial }_{\mu }-S\right)\psi =0.$$

$$S=-g^2⟨\bar{\psi }\psi ⟩=-g^2\sum _{\alpha }^{\mathrm{occ}}{\bar{\psi }}_{\alpha }{\psi }_{\alpha },$$

$$i\partial_{t}\left(\begin{array}[]{c}\psi_{1}\\ \psi_{2}\end{array}\right)=\left(\begin{array}[]{cc}S-m&-i\partial_{x}\\ -i\partial_{x}&-S-m\end{array}\right)\left(\begin{array}[]{c}\psi_{1}\\ \psi_{2}\end{array}\right).$$

$$i{\partial }_t{\psi }_1=(S-m){\psi }_1-{\partial }_x\frac{1}{S+m+i{\partial }_t}{\partial }_x{\psi }_1,$$

$$\frac{1}{S+m+i{\partial }_t}\approx \frac{1}{2m}-\frac{(S-m+i{\partial }_t)}{4m^2}$$

$$i{\partial }_t{\psi }_1=\left(S-m-\frac{{\partial }_x^2}{2m}\right){\psi }_1+{\partial }_x\frac{\left(S-m+i{\partial }_t\right)}{4m^2}{\partial }_x{\psi }_1.$$

$$i{\partial }_t{\psi }_1$$

$$\sqrt{m^2+p^2}-m-\frac{p^2}{2m}\approx -\frac{p^4}{8m^3}.$$

$$Q=\int dx\left(|{\psi }_1{|}^2+|{\psi }_2{|}^2\right)\approx \int dx{\psi }_1^{\ast }\left(1-\frac{{\partial }_x^2}{4m^2}\right){\psi }_1.$$

$$H⟶\tilde{H}$$

$${\psi }_1⟶{\tilde{\psi }}_1$$

$$i{\partial }_t{\tilde{\psi }}_1$$

$$S-m$$

$$\sigma$$

$$S-m$$

$${\tilde{\psi }}_1$$

$${\tilde{\psi }}_2$$

$$i{\partial }_t{\tilde{\psi }}_{k,1}^{(0)}$$

$$i{\partial }_t{\tilde{\psi }}_{k,1}^{(2)}$$

$$\frac{{\nu }_k}{2}=\frac{\arcsin y_k}{\pi }+\frac{\gamma }{\pi }\frac{y_k}{\sqrt{1-y_k^2}}.$$

$${\nu }_k=ϵ{\nu }_k^{(1)}+{ϵ}^3{\nu }_k^{(3)}$$

$${\nu }_k^{(1)}$$

$${\nu }_k^{(3)}$$

$$\sigma$$

$${\sigma }^{(2)}$$

$${\sigma }^{(4)}$$

$${\tilde{\psi }}_{\ell ,2}^{(0)}=-\frac{i}{2m}{\partial }_x{\tilde{\psi }}_{\ell ,1}^{(0)}.$$

$$\bar{\psi }\psi =|{\psi }_1{|}^2-|{\psi }_2{|}^2.$$

$$S^{(2)}$$

$$S^{(4)}$$

$$i{\partial }_t{\psi }_k^{(0)}=\left(-\frac{{\partial }_x^2}{2m}+S^{(2)}\right){\psi }_k^{(0)},S^{(2)}=-\sum _{\ell }2y_{\ell }|{\psi }_{\ell }^{(0)}{|}^2(\mathrm{LO}),$$

$$i{\partial }_t{\psi }_k^{(2)}$$

Full text

I Introduction

The massive Gross-Neveu (GN) model 1 is the quantum field theory of $N$ flavors of Dirac fermions with a scalar-scalar four-fermion interaction and Lagrangian

[TABLE]

We will only consider the case of 1+1 dimensions, where $g^{2}$ is dimensionless and the theory is renormalizable. The bare mass $m_{0}$ breaks explicitly the discrete chiral symmetry of the original, massless GN model ($\psi\to\gamma_{5}\psi,\bar{\psi}\psi\to-\bar{\psi}\psi$) and renders the model non-integrable. A number of explicit, analytical results have been obtained in the ’t Hooft limit ($N\to\infty,Ng^{2}=$const.) of the massive GN model in the past, using semi-classical methods. Thus, the vacuum 2 , baryons 3 ; 4 , multi-baryon bound states 5 , cold and dense matter 6 and the phase diagram at finite chemical potential and temperature 2 ; 7 ; 8 are by now well understood. Somewhat surprisingly, explicit results for static properties have turned out to be equally accessible in the massive and massless GN models, despite the fact that only the massless one is integrable. In particular, scalar mean fields for baryons are transparent and those for inhomogeneous phases of dense matter are finite gap periodic potentials, irrespective of whether the bare mass $m_{0}$ is included or not.

The situation changes once we look at dynamical problems. In the case of the massless GN model, it has proven possible to solve time-dependent scattering problems of multi-fermion bound states and write down general results in closed analytical form 9 . The method used was based on a relativistic version of the time-dependent Hartree-Fock (TDHF) approximation, supposed to become exact in the large $N$ limit. Since these results show only elastic forward scattering and factorized transmission amplitudes, there is little doubt that integrability is at work here. By contrast, as pointed out in 10 , a similar ansatz method does not yield self-consistent mean fields for the massive GN model, despite the fact that individual static baryons can be correctly described.

The aim of this paper is to elaborate on the difference between integrable, massless and and non-integrable, massive GN models. Static properties apparently give no clue about this difference, at least at large $N$. One way of going beyond static properties is to head towards non-equilibrium thermodynamics, e.g., by studying the bulk viscosity 11 . In that reference one also finds a pedagogical discussion of the (non-)integrability of GN models in terms of Feynman diagrams for inelastic processes. Here we propose to follow another route. We generalize previous baryon-baryon scattering calculations of the massless GN model to the massive one, looking directly for non-forward elastic scattering and inelastic reactions. Clearly, we cannot hope to carry out such studies without numerical computations. Due to the Dirac sea, it would be very challenging to do a numerical TDHF calculation from scratch. Therefore, in this exploratory study, we set ourselves a more modest goal. We try to identify the leading order contribution to backward and inelastic scattering in the vicinity of the non-relativistic limit only. This will enable us to build on a previously developed effective low energy theory for the GN model 12 , while keeping numerical computations manageable. At the same time, by combining analytical and numerical tools, we hope to get more insight than with a purely numerical approach.

We finish this introduction with a reminder about regularization and renormalization of the massive GN model 2 ; 3 . The Lagrangian (1) has two bare parameters, $g^{2}$ and $m_{0}$. After the regularization/renormalization procedure at large $N$, all observables can be expressed through two physical parameters, $m$ and $\gamma$. The relation between the bare quantities, an UV cut-off $\Lambda$ and the physical parameters is given by the vacuum gap equation

[TABLE]

In the massless GN model ($m_{0}=0$), the dimensionless coupling constant $g^{2}$ gets traded for the dynamical fermion mass $m$, an example of dimensional transmutation. In the massive model, in addition, the bare mass $m_{0}$ gets replaced by a physical parameter $\gamma$. In condensed matter physics where the massless and massive Gross-Neveu models can be used for instance to model trans- and cis-polyacetylene 13 ; 14 , $\gamma$ is called “confinement parameter”. It can also be related to the ratio of the dynamical fermion mass at $m_{0}$ and $m_{0}=0$,

[TABLE]

The physical parameters ($m,\gamma$) are renormalization group invariant.

The plan of the present paper is as follows. In Sect. II, we derive the TDHF equations for the massive GN model in the vicinity of the non-relativistic limit. After introducing an appropriate expansion parameter, in Sect. III we simplify these equations further, exhibiting the leading order (LO) and next-to-leading order (NLO) equations in detail. Sect. IV explains our method of solving these equations, whereas numerical results will be presented in Sect. V. We close this paper with a summary and conclusions, Sect. VI.

II TDHF near the non-relativistic limit

Our starting point is the Dirac equation with scalar potential $S$ in 1+1 dimensions,

[TABLE]

In TDHF theory appropriate for the large $N$ limit of the GN model, $S$ is the self-consistent mean field,

[TABLE]

where the sum runs over all occupied states, i.e., the filled Dirac sea and positive energy valence states. The aim of the present section is to set up the first two terms of a systematic, non-relativistic approximation to the full TDHF problem.

Using the Dirac-Pauli representation ($\gamma^{0}=\sigma_{3},\gamma^{1}=i\sigma_{2},\gamma_{5}=\sigma_{1}$) and pulling out the fast factor $e^{-imt}$ from the spinor $\psi$, we cast Eq. (4) into the Hamiltonian form

[TABLE]

Next, we eliminate the “small” component $\psi_{2}$ formally from (6),

[TABLE]

use the non-relativistic expansion

[TABLE]

and arrive at a Schrödinger-type equation for the “large” component $\psi_{1}$,

[TABLE]

Replacing $i\partial_{t}\psi_{1}$ on the right hand side by the LO expression $(S-m-\partial_{x}^{2}/2m)\psi_{1}$, we find

[TABLE]

The first term on the right is of LO, the second term the NLO correction to the relativistic kinetic energy,

[TABLE]

The third term is non-hermitean, reflecting the fact that after elimination of the lower component $\psi_{2}$, the norm of the upper component $\psi_{1}$ is not conserved. To NLO, the conserved charge is

[TABLE]

Accordingly, we redefine the Hamiltonian and the wave functions as follows,

[TABLE]

This yields the amended version of the Schrödinger equation (10), now with a manifestly hermitean Hamiltonian,

[TABLE]

We have generated the analogue of the Darwin term for a scalar potential (for a vector potential, replace the anti-commutators by commutators.) The reader will have noticed that this calculation follows closely the textbook evaluation of the fine structure of the hydrogen atom, see e.g. 15 , except for the missing spin-orbit term in 1+1 dimensions.

What we have done so far is the non-relativistic reduction of the Dirac equation with a classical potential, including fine structure corrections. From the hydrogen atom, we already know that this is not the whole story: the Lamb shift is still missing. It arises from vertex corrections and vacuum polarization due to Dirac sea. Fortunately, we have the tools at our hands to include such quantum-field corrections systematically as well. To this end we turn to an “effective no-sea theory” of the massive GN model, derived in Ref. 12 by “integrating out” all negative energy states. Referring to 12 for the technical details and derivations, we immediately jump to the results. Keeping only the first two terms of a systematic expansion, one finds for the mean field of the massive GN model

[TABLE]

Let us briefly explain the various symbols and terms. To LO, $S=m+\sigma$. Here, $m$ is the dynamical fermion mass which would arises from the Dirac sea in the full HF calculation but has to be put in by hand in the effective theory. The scalar field $\sigma$ — the field of the $\sigma$ meson — has a similar self-consistent structure as $S$ in (5). But here, the sum over occupied states only extends over the positive energy valence levels. The coefficient $\nu_{\ell}=N_{\ell}/N$ denotes the occupation fraction of the $\ell$-th bound state, a continuous parameter in the limit $N\to\infty$. The bare coupling constant $g^{2}$ gets replaced by a renormalized coupling constant $g_{\rm eff}^{2}=-\frac{\pi}{N(1+\gamma)}$. This value can be derived from the tadpole contribution to the fermion self energy by projecting out the negative energy states, see Ref. 12 . The other two terms on the right hand side are NLO and arise from vacuum polarization effects. As explained in Ref.  12 , one needs to re-sum an infinite number of diagrams contributing to the $\sigma$-meson propagator. In this process, the bare coupling constant drops out, so that the result is truly non-perturbative in $g^{2}$. To arrive at the final form, one expands the inverse $\sigma$-propagator in powers of $k^{2}$, so that the truncation underlying Eq. (15) is indeed consistent with the non-relativistic reduction of the Dirac equation.

The task is now to solve Eq. (14) for all bound states, using as self-consistency condition, Eq. (15). Although these equations look more complicated than the original TDHF equations (4,5), they are much easier to solve. It is sufficient to determine positive energy bound states self-consistently, rather than the whole Dirac sea plus valence levels. Thus, the structure resembles that of a non-relativistic TDHF problem, albeit with a more complicated interaction and self-consistency condition. Relativistic corrections are included by means of the fine structure terms as well as the higher order corrections to the no-sea effective theory. There are no more divergences, as regularization and renormalization have already been performed when deriving (15). Since the scalar condensate $\bar{\psi}\psi$ in (15) still contains the original two-component spinors $\psi$ rather than $\tilde{\psi}_{1}$, the equations are not yet in a form ready to be solved. Further simplifications will appear once we introduce a small expansion parameter and truncate all the equations consistently, the goal of the following section.

III Expansion parameter, LO and NLO TDHF equations

In order to arrive at a tractable set of equations and to avoid the unnecessary computation of complicated higher order terms, we introduce a formal expansion parameter $\epsilon$. The regime we are interested in is characterized by $v\sim\epsilon,y\sim\epsilon$, with $v$ the baryon velocity and $y$ the (inverse) baryon size parameter 3 ; 4 . The first condition is self-evident for a non-relativistic expansion. The 2nd one induces a matching non-relativistic expansion for the internal structure of the baryon. The characteristic exponential in the single baryon is exp$[2y(x-vt)]$, hence we treat $\partial_{x}\sim\epsilon,\partial_{t}\sim\epsilon^{2}$. Guided by the single baryon results, we assume the following expansions for $S-m$ and the spinors,

[TABLE]

The small component $\tilde{\psi}_{2}$ will be needed later on for the condensate. Inserting these expressions into the Dirac equation (14) with $S-m$ from (15) and equating powers of $\epsilon$ then yields the following LO and NLO equations,

[TABLE]

We have added a subscript $k$ labeling the bound states, as required by the general TDHF problem.

We still need the relation between $S^{(2)},S^{(4)}$ and the spinors $\tilde{\psi}_{k,1}^{(0)},\tilde{\psi}_{k,1}^{(2)}$ arising from the self-consistency condition (15) and expansions (16). For each baryon, there is a size parameter $y_{k}$ and an occupation fraction $\nu_{k}$ related in a non-linear fashion, namely 3 ; 4

[TABLE]

Since $y_{k}\sim\epsilon$, $\nu_{k}$ can be expanded as

[TABLE]

with

[TABLE]

Inserting Eqs. (16,19,20) into Eq. (15) yields

[TABLE]

where

[TABLE]

The various terms in (21) can be understood as follows. The scalar density in our representation of the $\gamma$-matrices is

[TABLE]

Remembering the different powers of $\epsilon$ in (16), this explains the LO term $\sigma^{(2)}$. To this order, the difference between $\psi$ and $\tilde{\psi}$ does not matter. $\sigma^{(4)}$ contains all NLO terms, taking into account the fact that $\nu_{k}$ has the expansion (19). The last line in Eq. (21) is due to interference terms between $\tilde{\psi}_{\ell,1}^{(0)}$ and $\tilde{\psi}_{\ell,1}^{(2)}$, the only NLO contribution to $|\tilde{\psi}_{\ell,1}|^{2}$, and a derivative term coming from the transformation from $\psi$ to $\tilde{\psi}$, see (13). In the $|\psi_{2}|^{2}$ term in (23), only the lowest order is needed, so that we get away with expression (22).

According to Eqs. (15,16) and (21), the relationship between $S^{(2,4)}$ and $\sigma^{(2,4)}$ is

[TABLE]

Time derivatives in $\partial_{\mu}\partial^{\mu}\sigma$ are of higher order than what is needed here. Eqs. (17,21,24) together with the normalization conditions are a closed set of equations determining the positive energy bound state spinors in LO and NLO.

Before going on, it is useful to ease the notation. Since we can now express everything through the “large” components $\psi_{1}$, we drop the subscript 1 from all the spinors. We also omit the tilde on all wave functions, replace $\sigma^{(2)}$ everywhere by $S^{(2)}$ and insert $\nu_{k}^{(1,3)}$. The basic TDHF equations then read to LO

[TABLE]

and to NLO

[TABLE]

These equations are valid for arbitrary $\gamma$.

Let us first look at the LO equation (25) which has the form of the multi-component non-linear Schrödinger (NLS) equation 17 . We notice that it does not contain $\gamma$ at all, so that the same equation is valid independently of the bare fermion mass. How is this possible? The answer is the same as for static baryons: The wave functions and mean fields for a baryon in the massive theory are identical to those for a baryon in the massless theory, but for a different fermion number. Indeed, $\gamma$ still appears in Eq. (18) relating the occupation fraction and the size parameter. In the potential $S^{(2)}$, this is cancelled exactly by the $\gamma$-dependence of the effective coupling constant. Thus the TDHF equation of the massive GN model reduces to coupled NLS equations, independently of the bare mass. In this strict non-relativistic limit, we can solve baryon-baryon scattering problems by simply taking the solution of the massless model and expanding it in powers of $y,v$ to LO. At this level, there will be neither backward scattering nor any inelastic reaction. The model does not yet loose integrability.

We now turn to NLO, Eq. (26). The potential $S^{(4)}$ consists of several terms with different $\gamma$-dependencies. Somewhat surprisingly, one can eliminate the $\gamma$ dependence as follows. Let us denote the solutions of Eq. (26) at $\gamma=0$ (chiral limit) by hatted quantities $\hat{\psi}_{k}^{(2)},\hat{S}^{(4)}$ (we just saw that $\hat{\psi}_{k}^{(0)}=\psi_{k}^{(0)},\hat{S}^{(2)}=S^{(2)}$). We then make the following “scaling” ansatz for the NLO quantities,

[TABLE]

Inserting (27) into (26) and using the fact that the hatted spinors satisfy Eq. (26) at $\gamma=0$, we arrive at our final set of equations

[TABLE]

now independent of $\gamma$. They have the structure of a system of inhomogeneous, linear PDE’s. The only input needed from the massless GN model are the LO quantities $\psi_{k}^{(0)},S^{(2)}$. Since the solutions of Eqs. (26) in the chiral limit can be inferred from the analytically known exact solutions simply by Taylor expansion, there is no need to ever solve the complicated equations (26) at $\gamma=0$. Solving Eq. (28) yields the solutions for all values of $\gamma$ at once.

The inhomogeneity has now been isolated in the term ${\cal F}\psi_{k}^{(0)}$. At this point, before turning to the solution, we can perform a non-trivial consistency check of our formalism. Since we know that all corrections must vanish in the static case, we expect that ${\cal F}=0$ if all solitons are at rest. In the static case, the $x$ and $t$ variables can be separated by going to stationary states,

[TABLE]

where $\phi_{k}^{(0)}(x)$ is the eigenfunction of the LO Hamiltonian with eigenvalue $\epsilon_{k}$,

[TABLE]

$S^{(2)}$ is defined as above, Eq. (25). Choosing the $\phi_{k}^{(0)}$ to be real, it becomes

[TABLE]

A simple exercise with Computer Algebra is the following: Evaluate $\partial_{x}^{2}S^{(2)}$ and eliminate the 2nd derivatives of the time independent wave functions with the help of the NLS equation (30,31). Then differentiate the resulting ${\cal F}$ with respect to $x$. If one eliminates once again $\partial_{x}^{2}\phi_{k}^{(0)}$ using (30), one finds that $\partial_{x}{\cal F}=0$, i.e., ${\cal F}$ is constant in the static case. Since it obviously vanishes asymptotically, it has to be identically zero. For time dependent problems this proof fails and we cannot avoid solving the inhomogeneous, linear system of PDE’s, Eq. (28).

IV Method of solution

So far, things are valid for any number and type of scatterers. Let us first look at the trivial case of a single baryon, where we would expect no correction at all, since there is a frame in which it is static. In this case the inhomogeneous term in Eq. (28) must vanish, so that the system of equations admits the trivial solution $\chi_{1}^{(2)}=0$. As the single baryon of the massive GN model has the same structure as the original Dashen-Hasslacher-Neveu (DHN) baryon of the massless model 16 , we shall refer to the single baryon as DHN baryon. For the single DHN baryon with size parameter $y$ and velocity $v$ and using units where $m=1$ from now on, we have

[TABLE]

This yields indeed ${\cal F}=0$ when inserted into (28).

The first non-trivial example is scattering of two DHN baryons. In view of the exploratory character of our study, we simplify things as much as possible. We take two identical DHN baryons ($y_{1}=y_{2}=y$) in the center-of-mass (cm) frame ($v_{1}=-v_{2}=v$). Eq. (28) then becomes

[TABLE]

The solution of the LO problem is well known 9 ; 10 ,

[TABLE]

The prefactors $\lambda^{-1},\lambda$ in $U_{1,2}$ just shift the variables $x,t$ by a constant and have been chosen so as to slightly simplify the NLO terms in the chiral limit. The entries in Eq. (34) can be inferred from the known exact solution by transforming from the chiral representation of Dirac matrices used in Ref. 9 to the Dirac-Pauli representation, replacing $y\to\epsilon y,v\to\epsilon v,x\to\epsilon^{-1}x,t\to\epsilon^{-2}t$ and performing a Taylor expansion in $\epsilon$ to LO. The phases of the two bound state spinors are of course arbitrary, but must be chosen such that there are no linear terms in $\epsilon$ in Eq. (16). One can easily verify that the quantities in (34) satisfy the LO equations and the normalization condition. Using this input, ${\cal F}$ turns out to be significantly simpler than its three constituent terms,

[TABLE]

This is non-zero only when the two baryons overlap. Asymptotically, ${\cal F}$ vanishes and one has to deal with coupled, homogeneous, linear PDE’s for $\chi_{k}^{(2)},\chi_{k}^{(2)*}$.

Consider the homogeneous system first, since this already provides us with valuable information about the possible outcome of a baryon-baryon collision in the massive GN model. We discuss separately the situation before and after the collision. The incoming channel consists of two well separated DHN baryons. Here, we already know that there are no bare mass corrections, except for the relationship between $y$ and $\nu$. Hence we can assume the initial condition $\chi_{k}^{(2)}=0$ for $t\to-\infty$ for $k=1,2$. After the collision, the possible final states are determined by the non-trivial solutions of the homogeneous system

[TABLE]

Since the TDHF approach cannot give the complete information about individual reaction channels but treats them only in an average way (due to the assumption of a single Slater determinant), we expect a superposition of different final states. The weights of specific states can only be determined by solving the full, inhomogeneous system of PDE’s (33) numerically. The homogeneous system (36) can actually be solved analytically as follows. Suppose we can find a small deformation of the unperturbed solution $\psi_{k}^{(0)}$ of the NLS equation (25) such that the result is again a solution,

[TABLE]

Linearizing Eq. (37) in $\delta\psi_{k}^{(0)}$ yields a solution of Eq. (36), namely $\chi_{k}^{(2)}=\delta\psi_{k}^{(0)}$. If the deformation can be taken to be infinitesimal, i.e., if the solution $\psi_{k}^{(0)}+\delta\psi_{k}^{(0)}$ is continuously connected to the unperturbed solution $\psi_{k}^{(0)}$, this solution will be exact. Thus, in order to survey the possible asymptotic solutions $\chi_{k}^{(2)}$, all we have to do is list the solutions of the NLS equation that can be obtained by a continuous deformation of the unperturbed solution (35). This should already enable us to characterize the possible final states in a baryon-baryon collision of the massive GN model.

Consider the elastic channel first. A forward scattered baryon will emerge with a time delay and a phase factor. The first change corresponds to the solutions

[TABLE]

with real coefficients $A_{x,t}$. To confirm that these are indeed exact solutions of (36), start from the NLS equation for $\psi_{k}^{(0)}$ and just differentiate this equation with respect to $x$ or $t$. In order to describe the phase shift in an analogous way, multiply $\psi_{k}^{(0)}$ by $e^{i\delta_{k}}$ and insert it into the NLS equation. Differentiating with respect to $\delta_{k}$ then yields another solution of (36),

[TABLE]

again with real $A_{k}$. Clearly, the solutions (38,39) cannot account for elastic backward scattering, expected in the massive GN model. It is not hard to find the corresponding deformation. The multi-component NLS equation

[TABLE]

remains valid under unitary transformations of the $\psi_{k}$, here under the group U(2) since $k=1,2$ only. An infinitesimal U(2) transformation can be parameterized as

[TABLE]

The U(1) part and the $\tau_{3}$ rotation have already been accounted for by (39), so that the only new solution we get is

[TABLE]

with complex coefficient $C$. Since $\psi_{1}^{(0)}$ and $\psi_{2}^{(0)}$ are moving in opposite directions at the same speed, this is exactly what it takes to describe elastic backward scattering.

We now look for deformations of $\psi_{k}^{(0)}$ related to inelastic processes. The simplest possibility is that the baryon changes its velocity, i.e.,

[TABLE]

with real $A_{v}$. The baryon may also change its size parameter (and thereby fermion number). Due to the factor of $2y$ in the potential $S^{(2)}$, we cannot simply differentiate $\psi_{k}^{(0)}$ with respect to $y$ in this case. A simple calculation shows that the following modified expression is an exact solution of (36),

[TABLE]

Presumably, this does not yet exhaust all possibilities. It is known that the (multi-component) NLS equation possesses solutions with more than one soliton in each component 17 . We have not found a simple way of relating these multi-soliton solutions continuously to the standard solution, but cannot rule out such a possibility. Then these more complicated solutions of the NLS equation might also play some role in inelastic processes. In the present work, we shall not consider multi-soliton deformations any further, but see how far we can get with the above, simplest solutions.

Eventually, we must solve the inhomogeneous system of PDE’s (33) numerically. Since this problem amounts to solving an inhomogeneous, time dependent Schrödinger equation with time dependent Hamiltonian, a number of numerical methods are available in the literature. We follow the method described by Puzynin et al. 18 . It is a higher-order stable operator-difference scheme, generalizing the Crank-Nicolson scheme. It can be derived from the Magnus expansion of the evolution operator for one time step. We actually used the 2nd order variant described in detail in Ref. 18 . In order to solve the equation

[TABLE]

one divides the temporal interval $[0,T]$ into $K$ steps of length $\tau$ ($t_{k}=k\tau,k=0,1,...,K$). We also need the intermediate points $t_{k+1/2}=t_{k}+\tau/2$ and denote $\psi^{k}=\psi(t_{k}),\psi^{k+1/2}=\psi(t_{k+1/2})$. Define

[TABLE]

Then the step from $t_{k}$ to $t_{k+1}$ goes as follows,

[TABLE]

with $\alpha=1/\sqrt{3}-i$. If we choose a spatial grid with $M$ points, $H$ and $F$ will be $4M\times 4M$ matrices due to the coupling of $\chi_{1},\chi_{2},\chi_{1}^{*},\chi_{2}^{*}$. Similarly, $\psi$ and $Q_{1,2}$ are $4M$ component vectors. Consequently, Eqs. (47) are two systems of $4M$ linear, inhomogeneous algebraic equations which can be solved by standard methods. In the particular case at hand, due to the symmetry of the scatterers and of the initial conditions, it is actually possible to reduce the dimension by a factor of 2 by restricting $x$-space to a half line. The number of mesh points in time and space were chosen such that a further increase would not show up in the figures below. We found that 250 points on the time axis and 250 points on half the space axis were adequate. In this case, all numerical computations could still be done with Maple.

V Results

A necessary condition for being able to trust the results is that we work in a regime where the non-relativistic expansion converges well in the massless limit. Since we know the exact result in this case, this is something which can be checked. To this end, we have computed the scalar mean field $S$ and the bound state spinors $\psi_{1,2}$ for scattering of two DHN baryons in the cm frame analytically. We then introduce the parameter $\epsilon$ ($y\to\epsilon y,v\to\epsilon v,x\to x/\epsilon,t\to t/\epsilon^{2}$) and expand potential and spinors in $\epsilon$ to NLO. The resulting expressions are too long to be shown here, but can readily be generated by Computer Algebra. We have also checked that this indeed solves Eq. (26) at $\gamma=0$, a useful further test of the above formalism. We then compare the exact results with the LO and NLO approximations by looking at animations of the corresponding plots throughout the whole collision process. We identify the region in ($y,v$)-parameter space where the difference between LO and NLO is of the order of 10$\%$, but the difference between NLO and the exact result of the order of a few $\%$ only. One example of a snapshot of such a survey is shown in Figs. 1,2,3 for $y=0.5,v=0.3$. The result of this search is the region $0.4\leq y\leq 0.6$, $0.1\leq v\leq 0.6$. There we are in a situation where the NLO approximation is quantitatively reliable in the massless limit. The corrections which we compute for the non-integrable, massive model then also have a good chance of being trustworthy. We cover this preferred parameter region with 18 points in steps of $\Delta y=0.1,\Delta v=0.1$. Outside this region, the NLO corrections are either negligible or too big to truncate the non-relativistic expansion after two terms.

We now turn to the results of solving the inhomogeneous system of PDE’s (28) numerically. In the following subsections, we discuss backward and forward scattering separately.

V.1 Backward scattering

Backward scattering has the unique feature that the result for $\gamma=0$ is strictly zero, due to the integrability of the massless GN model. Hence the correction proportional to $\chi_{k}^{(2)}$ in Eq. (27) represents the full wave function $\psi_{k}$ in this region. Since this is of O($\epsilon^{2}$), the density of backward scattered fermions resulting from the $k$-th bound state,

[TABLE]

is correct, although it is of O($\epsilon^{4}$) and we have only been working to O($\epsilon^{2}$). Integration over the negative half-axis then yields the reflection coefficient $R$ and a “reduced” reflection coefficient $R^{\prime}$,

[TABLE]

These are also of O($\epsilon^{4}$), but can be trusted for the same reason.

In Fig. 4, we give an overview of our results for the backscattered, “reduced” density (i.e., expression (48) without the $\gamma$-dependent factor) for $y=0.5$, after the collision. At the lower velocities, one sees a clear peak, gradually decreasing towards the highest velocity. For $y=0.4$ and $y=0.6$, the results look similar and need not be shown here.

In order to understand the nature of the peak, we confront it with our theoretical expectation based on solutions of the homogeneous PDE’s (36). The only candidate for backward scattering which we have identified in Sect. IV is Eq. (42) describing elastic scattering. In Figs. 5 and 6 we compare the numerically calculated $|\chi_{1}^{(2)}|^{2}$ and $\chi_{1}^{(2)}$ at $y=0.5$ and the lowest velocity, $v=0.1$, to (42). We get perfect agreement between numerical calculation and analytical prediction for the value $C=0.110-0.022i$. The fact that backward scattering at the lowest energy is purely elastic is confirmed by comparable fits at the other two values of $y$, yielding $C=0.066-0.016i$ at $y=0.4$ and $C=0.167-0.027i$ at $y=0.6$.

If we move to higher velocities, a background is appearing below the elastic peak, see Fig. 7 for an example of the reduced density at $y=0.5,v=0.4$. At the same time, the fit to the elastically backscattered wave function deteriorates, see Fig. 8. The background can only be due to the onset of inelastic processes at higher energies, but so far we are lacking a way of analyzing it quantitatively. At $y=0.5$ and the highest velocity, $v=0.6$, the background becomes even more prominent, see Fig. 9, and a fit with the elastic wave function impossible due to interference with the large, unknown background. To illustrate the $y$ dependence, we show the corresponding plots at $y=0.4$, Fig. 10, and $y=0.6$, Fig. 11. Evidently, the ratio of elastic peak to inelastic background increases with decreasing $y$. This is what one would expect qualitatively, since the fermions are more loosely bound for smaller $y$, cf. Eq. (30).

Finally, Fig. 12 shows the reduced reflection coefficient $R^{\prime}$ for all calculated points. It differs from the reflection coefficient $R$ by a $\gamma$-dependent factor, see Eq. (49). It is impossible to neatly separate elastic scattering from the inelastic background, since everything is coherent in TDHF. To get at least a rough idea, we have fitted the smooth background by a 2nd order polynomial in $x$ and subtracted it, assuming incoherence. Fig. 12 shows both the total and elastic values of $R^{\prime}$ thus obtained.

V.2 Forward scattering

In the case of forward scattering it would not make sense to plot the density $|\chi_{k}^{(2)}|^{2}$. Since the LO contribution is non-vanishing, a quantity of O($\epsilon^{4}$) cannot be trusted here. The relevant physical observable of O($\epsilon^{2}$) is the NLO change in density due to the bare fermion mass,

[TABLE]

We first give an overview of our numerical results for $\delta\rho_{1}^{\prime}$ at $y=0.5$ in Fig. 13. The solution of the TDHF equation yields the whole time evolution of $\delta\rho_{1}^{\prime}$, but here we only show a snapshot taken after the collision. As $v$ increases, the curves transform from an antisymmetric into a symmetric shape. In order to highlight the transition between these two extrema, we have included one extra point along the $v$-axis where $\delta\rho_{1}^{\prime}$ is close to zero ($v=0.45$). The shapes of the extremal curves at $v=0.1$ and $v=0.6$ have a simple interpretation. Since $\delta\rho_{1}^{\prime}$ is the difference between two bell-shaped densities, at $v=0.1$ the curves indicate a spatial shift between the densities, as expected in purely elastic scattering. The shape at $v=0.6$ on the other hand is suggestive of a change in the width of the density, consistent with inelastic reactions where the size parameter $y$ (and therefore fermion number) changes.

Let us test this interpretation against the analytical, asymptotic predictions discussed in Sect. IV. If baryon-baryon scattering in the massive GN model is purely elastic, $\chi_{1}^{(2)}$ must be proportional to a linear combination of $\partial_{x}\psi_{k}^{(0)}$, $\partial_{t}\psi_{k}^{(0)}$ and $i\psi_{k}^{(0)}$ with real coefficients, see Eqs. (38,39). The first two solutions account for a change in time delay and cannot be distinguished here, the 3rd accounts for a change in the forward scattering phase shift. Figs. 14 and 15 show a corresponding fit at $y=0.5,v=0.1$,

[TABLE]

with $A=-0.130,B=-0.0707$. Thus at the lowest velocity $v=0.1$, everything is consistent with purely elastic scattering, both in forward and backward direction.

At higher values of $v$, it is difficult to find any satisfactory fit with the available basis functions. By way of example, consider $y=0.5,v=0.6$ where inelastic processes are expected to be important. In Fig. 16, we have tried to fit $\delta\rho_{1}^{\prime}$ with a linear combination of solutions (38) and (44), i.e.,

[TABLE]

Even the best fit values $A=0.00195,B=-0.00599$ are not satisfactory, although the shape of the curve is qualitatively reproduced. Taking into account the further inelastic solution (43) does not improve matters, so that our parametrization of the inelastic solutions is obviously incomplete. As mentioned in Sect. IV, one could perhaps generate further candidates by considering multi-soliton solutions of the single component NLS equation, but this is left for future work.

So far, all results for forward scattering shown refer to $y=0.5$. We have performed the same kind of calculations at $y=0.4$ and $0.6$ as well. The plots analogous to Fig. 13 look similar, except that the transition from “shift” to “broadening” type shapes happens around $v=0.35$ at $y=0.4$ and $v=0.55$ at $y=0.6$. Thus, elastic scattering prevails up to higher energies for larger values of $y$, as already noticed in backward scattering. For all 3 values of $y$, the fit (51) is excellent at $v=0.1$, with parameters $A=-0.0826,B=-0.0527$ at $y=0.4$ and $A=-0.150,B=-0.081$ at $y=0.6$, confirming elastic scattering dominance at the lowest velocity considered.

Recall that there are two independent corrections to the fermion density of O($\epsilon^{2}$). The first one arises from relativistic corrections to the NLS equation at $\gamma=0$ (“fine structure” and vacuum polarization type effects), discussed above in the context of Figs. 1,2,3. It is present in the massless GN model, does not destroy integrability and can be computed in closed analytical form. The second one is $\delta\rho_{1}$ from Eq. (50) originating from the bare mass term and to be computed numerically. It is interesting to compare these two corrections to the LO density and to each other. This is done in Fig. 17 for $y=0.5,v=0.1$ (purely elastic scattering) and in Fig. 18 for $y=0.5,v=0.6$ (important inelastic contribution). For better visibility, the two NLO corrections have been multiplied by 5. In the elastic scattering case, Fig. 17, the two corrections almost coincide, but this is most likely accidental. Both corrections go into the same direction of decreasing the time delay. In Fig. 18 where inelastic scattering plays a role, the relativistic correction is an order of magnitude larger than the mass correction and has opposite sign, but the shapes of the curves are again similar. The sign of the mass correction is such that the unperturbed density gets broadened, as one would expect from inelastic reactions where $y$ (and therefore fermion number) should decrease.

We finish with a remark about conservation of fermion number. The integral over the $x$-axis of both NLO contributions to the density in forward scattering vanishes, as we have checked numerically. The LO density is normalized to 1. Hence the total number of forward scattered fermions is the same in the massless and massive models. Nevertheless, we have computed a reflection coefficient of O($\epsilon^{4}$) in the massive GN model, whereas it vanishes in the massless limit. There is no contradiction between these findings. In forward direction, all we have shown is the vanishing of the O($\epsilon^{2}$) contribution. In order to check the expected reduction of the forward scattered fermions of O($\epsilon^{4}$), it would be necessary to do a full NNLO calculation. This is clearly impossible with the formalism developed so far.

VI Summary and conclusions

Recent progress on baryon-baryon scattering in the massless GN model owes much to the fact that this is an integrable quantum field theory. This enables us to compute analytically processes as complex as the collision between relativistic, composite particles. Unfortunately, the physical phenomena which can be studied in this manner are also severely limited by integrability. As is often the case, the same facts which make a system exactly solvable also render it somewhat unphysical. Only elastic forward scattering is allowed, no matter at what energy one collides which types of projectiles. This motivates us to turn our attention to non-integrable toy models promising a richer dynamics, like backscattering and inelastic reactions. Unlike real particle production, break-up processes should not be suppressed at large $N$, so that the massive GN model is a good candidate for exploring the transition from integrable to non-integrable systems with semi-classical methods.

Since a full numerical TDHF calculation seems to be exceedingly hard even in such a toy model, we had to compromise and focus on the vicinity of the non-relativistic regime. In this case, one can take advantage of the benefits offered by an effective no-sea theory applicable here. This enables us to formulate the TDHF problem in a manner close to non-relativistic many-body theory, based on a time dependent Schrödinger equation and systematically calculable relativistic and quantum field theoretic corrections. The result of such a calculation is a quantitative field theoretic calculation with controlled approximations and without adjustable parameters. This paper is actually the first successful application of the effective no-sea theory of Ref. 12 to a problem which has not been addressed by any other method before.

Our first important result was rather easy to get at, but is nonetheless quite interesting: In the non-relativistic limit, the multi-component NLS equation can be used to solve scattering problems both in the massless and massive GN models. To LO, integrability is maintained, even in time dependent problems. This extends the region of applications where the massive and massless GN models are equally tractable by analytical methods from static and thermodynamic problems to dynamical problems, at least at very low energies.

The second and more difficult part of the investigation concerns the NLO corrections. This is more interesting, because it shows the onset of physical processes absent in integrable models, like backscattering and inelastic reactions. Relativistic and quantum field theoretic corrections usually associated with spectroscopy (fine structure, Lamb shift) now break integrability and induce these forbidden processes. Thanks to our restriction to a certain parameter range, we could reduce the task to a system of coupled, inhomogeneous, linear PDE’s amenable to numerical solution by standard methods. The possible final states of a baryon-baryon collision can be obtained analytically by solving the corresponding homogeneous system. They are related to a certain class of solutions of the multi-channel NLS equation.

The picture emerging from backward scattering is perhaps the cleanest. We find a prominent elastic peak at low velocities above a smooth background. This background is steadily rising with increasing energy. Quantitative agreement with density and wave functions at the lowest velocity has been achieved by assuming purely elastic backward scattering. It is possible to determine the reflection coefficient and the phase of the backscattered wave function quantitatively. At higher velocities, a similar analysis is hampered by our inability to parametrize the inelastic background. A well known problem characteristic of the TDHF approach is the fact that different reaction channels are hard to disentangle, since they enter in a coherent, average way due to the assumption of a single Slater determinant. Interesting findings are the fact that the $\gamma$ dependence is somehow trivial, so that one does not have to repeat the calculation for different $\gamma$’s, and that a reflection coefficient of O($\epsilon^{4}$) can be computed reliably, even though the whole calculational scheme is truncated at O($\epsilon^{2}$), This is unique for a correction to an integrable model where backscattering vanishes to LO.

In forward scattering, since the LO term does not vanish, the NLO terms are always interference terms and harder to interpret. Again, the cleanest result is elastic scattering which exhausts what we see at $v=0.1$. We can compute the change in time delay and in scattering phase shift due to the bare mass, again with a factorized $\gamma$-dependence. The interference with inelastic processes on the other hand is impossible to analyze in detail with our methods. At the highest velocity studied, we see qualitatively that the density is broadened in $x$-space, corresponding to lower $y$ or loss of fermion number. The qualitative change in the density from shift to broadening with increasing energy and at all $y$-values is conspicuous. Thus there is no doubt that we have seen both backscattering and inelastic reactions in the NLO calculation, unlike at LO.