Boltzmann-like kinetic models and Boltzmann maps

Vol. 25, No. 2, pp. 53-67, 1997 Copyright@1997 Elsevier Science Ltd Printed in Great Britain. All rights reserved 0895-7177/97 $17.00 + 0.00

Mathl. Comput.

Modelling

PII:SO8957177(97)00006-X

Boltzmann-Like Kinetic Models and Boltzmann Maps L.

RONDONI

Dipartimento di Matematica, Politecnico di Torino Corso Duca degli Abruzzi, 24, 10129 Torino, Italy rondoniQpolito.it (Received May 1996; accepted

July

1996)

Abstract-This

work shows how Boltzmann maps can be used in the study of a class of kinetic models formally similar to those recently proposed in population dynamics. Boltzmann maps are nonlinear, discrete, doubly-stochastic processes which have been introduced in statistical physics for the description of nonequilibrium phenomena. The key feature of Boltzmann maps is that they increase the entropy, when applied to nonstationary states, if their interaction operators have a spectral gap. This fact can be used to analyse the asymptotic behaviour of the relevant dynamics. In particular, in some circumstances we can prove that stationary states exist, are unique, and stable. Moreover, we can supplement these results with theorems on the global convergence of the dynamics to such stationary states. Keywords-Nonlinear

stochastic

processes, Entropy,

Stability.

1. INTRODUCTION Recently,

a model equation

for the evolution

of dominance

in populations

with anonymous

or-

In that work, dominance is assumed to be a scalar quantity, taking values in a finite interval, by which individuals of the same population can be discriminated. Technical generalizations and qualitative analysis of the initial value problem have been considered in [2,3]. In particular, [2] generalizes the original model to systems of equations, and

ganisms

has been proposed

in [l].

provides existence and uniqueness results for the initial value problem, under quite general assumptions on the initial data. The paper [3], on the other hand, deals with models with multiple collisions, and with the existence of equilibrium solutions, while leaving open the problem of their stability. Further developments of the models proposed in [2,3] make possible the modelling of cellular interactions in biological systems, as reported in [4], but we do not deal with such systems in this work. A common feature of the relevant equations is that their mathematical structure is reminiscent of classical kinetic theory, especially of the structure of the Boltzmann equation. In particular, the model of [l] was obtained under the assumption that interactions only involve pairs of individuals (binary interactions), thus producing quadratic nonlinearities in the interaction terms. In this paper, we consider a class of models, obtained from those mentioned above, with some restriction on the functions describing the interactions. Moreover, we introduce suitable approximations for such models, which can be viewed as Boltzmann maps [5,6]. The models of our interest may be derived under the following assumptions. This work has been supported in part by GNFM-CNR

(Italy),

and by the Australian

Ftesearch Council

Typ-t 53

by AA&-W

54

L.

F~OND~NI

ASSUMPTION 1.1. The system under investigation is constituted by the interacting individuals of p E N populations. The physical state (dominance) of the individuals within each population is described by a variable u E [-l,l].

ASSUMPTION 1.2. The probability distribution over the sample space [-1, l] at each instant of time t E [0, T], and for each population, is defined by the probability density fi :

[W] x

[-1, l] + ]w+,

i = 1,. . . ,p.

(1.1)

Then,

represents the probability of finding an individual of the ith population with state in [ul, uz] at time t.

ASSUMPTION 1.3. The number of interactions per unit time (encounter rate) involving one individual from population i (called individual i), with state v, and one individual k, with state 20, isgiven byqik(v,w), where~~~maps[-1,1]2into+IR+U{O}foralli,k=1,...,p. we mean the Cartesian product of I copies of [-1, 11.

Here, by [-l,l]“,

ASSUMPTION 1.4. The probability density for an individual i with state v, to change it to u in an encounter with an individual k with state w is denoted by $ik(v, w, u), where ?+!&maps [-1, 113 into R+ U (0) for all i, k = 1,. . . ,p.

ASSUMPTION 1.5. The evolution of the distribution of states is determined exclusively by the effects of the interactions. These only involve pairs of individuals (binary interactions). the total derivative of fi does not depend on u.

Moreover,

Using these hypotheses, a balance of losses and gains of individuals yields

(1.3) -,,t&J

k=l

%k(U,

[-1,11

w)fk(t,

W> dw,

i=l,...,p,

where the conditions

1

tiik(V,

i--1,11

w, u> du = 1,

forallv,wE[-l,l]

andall

i,k=l,...,

p

(1.4)

must be satisfied, so that (formally) the following holds:

J

I 1 1l fi(t,U)du - ,

= 1,

foralli=l,...,p.

Models of this form are called conservative.

REMARK 1.1. The gain term-i.e., the first term on the right-hand side of equation (1.3)implements the idea that one individual i may only result from interactions where at least one individual i is involved. If we relax the condition expressed by Remark 1.1, the relevant equations should be correspondingly modified. For instance, interactions could be allowed to also change the species of

Boltzmann-Like KineticModels individuals,

55

i.e., the populations to which they belong, similar to what happens in chemical re-

actions. In that case, equation (1.3) must be replaced by a different set of equations, propose the following:

and we

(1.5)

where &&(t&W, u) du = 1,

for all a, b = 1, . . . ,p

and all

w,w E [-l,l].

(1.6)

Here, the functions q,b have same meaning as for the models of equation (1.3), while $Ab(v, w, U) represents the probability density for an individual a with state TJ, and an individual b with state w, to produce one individual i with state u. We refer to the models of form equation (1.5), and satisfying equation (1.6) as quasi-conservative. Then, a quasi-conservative system is conservative if $& = &a?+!&,with 6ij the Kronecker symbol, for all i, a, b = 1, . . . ,p. To use this formalism in the framework of chemical kinetics, the language of population dynamics must be interpreted accordingly. In this case, individuals represent single molecules, populations identify chemical species, while the state of one individual may refer to difierent properties of the molecules. For instance, the state of a molecule could represent its degree of ionization, or its excitation level, and so on, according to the situation at hand. Section 2 of this paper describes the class of models which constitute the starting point of our investigation, and introduces for them suitable discrete approximations. In given circumstances, such approximations constitute examples of Boltzmann Maps (BM), which are introduced in Section 3. Section 4 contains our original results, which are proved in Section 5. In particular, we obtain existence, uniqueness, and stability of the stationary solutions both for the discrete, and for the continuous (with respect to time) models of our interest. Furthermore, we prove that such stationary states are global attractors for the relevant discrete dynamical systems, so that all initial probability distributions converge to them, in the long time limit. In order to get these results, Theorem

4.5 is proved first, which is also relevant to the theory of Boltzmann

maps.

Section 6 concludes the paper with some examples.

2. A CLASS

OF MODELS

AND

THEIR

APPROXIMATIONS

In this section, we consider a special class of models of the form given above, identified by certain conditions on the functions q and ?(I. We then introduce suitable approximations for them, whose qualitative analysis is reported in Section 4, and we observe that part of the relevant results holds for the time-continuous models as well. Note that our approximations are not very efficient from a computational point of view; their merit is that they can be thoroughly analyzed, so that global convergence theorems, based on the existence of strict Lyapunov functionals, can be obtained for them. Similar functionals are not easily (or at all) found in other systems, hence, the difficulty in more general settings to go beyond stability results. Let us introduce a finite partition {I,.}~!1 of the sample space I = [-l,l], such that 1,. = [u,., u,+l) for f = 1,. . . ,m - 1, and 1, = [um,um+l], where u, = -1 + 2(7- - 1)/m. We impose the following restriction. RESTRICTION 1. We restrict our attention to the cases where the functions rlab = q > 0 are constant for all a, b, while the functions $J:,(+, ., u) : [-1, 112 --+ lR+ U (0) are simple for every u E I: (2.1)

L.R,~ND~NI

56 and $&, E L1([-1,113) interval Ik. Furthermore,

for all i,a,b = l,... p. Here, xk is the characteristic we require that equation (1.6) be satisfied.

function

of the

The condition that $ib be simple is not a serious limitation, given the density of simple functions in L1, while the restriction implied by equation (1.6) has been adopted in the literature in various circumstances. restriction.

Therefore,

the new technical

RESTRICTION 2. For the functions

assumption

$&, obeying

we are going to make is the following

Restriction

1, we further

require

that

&,(v, w, .) dv dw = 2~.

Restrictions

1 and

consider.

Under

Equation

(2.3) suggests

it with respect

dfi x

A

dt

&$

-&2

m2 = 4

Equation associated

I xl xl k 1 I

approximations

we are going

to

(1.5) takes the form

of the original

model, formally

obtained

(2.3) by integrating

u over any Is:

?&,(v, w, u) a=1 k=l

J

whose

1, equation

k ~/&&k(t) - &z(t)

fb,l(t)

I=1

, I

for i = 1,. . . ,p, z = 1,. . . ,m.

(2.4)

dv dw du,

(2.4) constitutes a nonlinear system of ordinary differential equations (ODE) which, with initial conditions f” E Wmp, has a unique solution in a given time interval [0, T): f

Moreover,

the class of problems of Restriction

an approximation

to the variable

= 77

b=l

2 identify

the conditions

(2.2)

because

(t, f”) = {f&r (c (f&H 17

of equation

(1.6), the following

&rb(u) du = 1,

As shown in Sections

4 and 5, equation

holds:

$1 F i=l

fi,Z(t)

i=l

= const,

for all t E [0, T).

(2.5)

z=l

(2.5) implies

Ml = f:

fort E [O,T)*

that the L’-norm

of f,

2 ksl

z=l

is always 1,if f” is normalized and taken from the positive cone of lP”P. The model equation (2.4) may be further approximated, and replaced by a given mapping of IPP, by discretizing them with respect to time. Because we are interested in the qualitative analysis of the solutions of the initial value problems, it suffices to adopt Euler’s scheme, which is the simplest of all, and which can be shown to be a Boltzmann map (BM) in some instances. Assuming the normalization of f, Euler’s scheme for equation (2.4) reads

Boltzmann-Like KineticModels

57

where

and fi,= = f+ xbI , fk,r. In case p = 1, the superscripts (i, ab) may be removed, and we write

fx=

cpx,kl

fl,x, f@+‘) 5

=

fl:

fq;:;,

h,kl

=

&,::, (2.8)

t&k&)fi(n).

k,l=l

3. BOLTZMANN

MAPS

AND

NUMERICAL

SCHEMES

In this section, we introduce Boltzmann maps (BM) [5,6], and we show that they may be viewed as Euler’s schemes like equation (2.8). Therefore, given approximations of continuous models described in Section 2 happen to be BM, if some conditions are satisfied. To verify this fact, consider a population dynamics model with m species, and binary interactions. The species can be different populations, different states within a population, or both. We postulate that the dynamics of the system satisfies the following assumptions. ASSUMPTION 3.1.

Interactions only involve two individuals and preserve their number, i.e., the

encounter between the individuals i and j always results in another pair of individuals, k and 1 say The system can be represented by elements of the simplex Q c IV” of probability measures on the sample space Cl = (1, . . . , m}. One element f = { fi}& E Q represents a state of the system, and its component fi, the probability of finding one individual i in it.

ASSUMPTION 3.2.

The probability of finding individuals i is independent of the probability finding individuals j. Thus, the joint probability of finding both is the product of the two.

ASSUMPTION 3.3.

of

Assumption 3.3 is equivalent to Boltzmann’s hypothesis of molecular chaos. Let Q2 c RmP be the simplex of probability measures on R x Sz. A doubly stochastic matrix acting on Qz is an m2 x m2 nonnegative matrix, whose column sums and row

DEFINITION 3.1.

sums are 1. A doubly stochastic matrix transforms Q2 into itself, because it is stochastic [6]. The effect of interactions on the product of probability measures representing the state of the system is expressed by linear, doubly stochastic matrices.

ASSUMPTION 3.4.

DEFINITION 3.2. Given a probability measure f = {fl, . . . , fm} on fi = (1,. . . , m}, the map d : Q + Q2, defined by g(f) = f @ f, is called sampling map. The map IE : Qz -+ Q, acting expectation as lE(F) = {C3m,i Fij}g”=,, for all F = {{F&}~C”=l}~l E Q 2, is called conditional onto the first factor.

Assume the state of the system the state on Q is given by E(F).

ASSUMPTION 3.5.

{{F~~}~l}~l,

is known on Q2, and it takes the form F =

The system evolves exclusively in consequence of interactions of individuals. Such an evolution is observed at discrete times n E M.

ASSUMPTION 3.6.

With these assumptions, the evolution determined by a BM involves three steps. (I) From a probability measure for single individuals, produce a product probability. (II) Transform the product probability to account for interactions. (III) Reduce the resulting joint probability to a probability for single individuals.

L. RONDONI

58

Such steps constitute the foundations of the arguments which led Boltzmann to his famous equation, after limits to vanishing volumes and time intervals had been taken. Differently, BM evolve probability measures without resorting to limiting procedures, and they can be formally defined as follows. DEFINITION 3.3. Assume that T : &2 +

Q2 is a doubly stochastic matrix. The map r : & + Q

given by f A

f@f

Z+

for every f E &, is called a Boltzmann

T(f@f)

+% f’=Q),

(3.1)

map of order 2.

To express in explicit form the dynamics of the BM (3.1), let us say that the entries of T can be written as Tij,kl, where the indices i, j = 1, . . . , m label the rows of T, while k, 1 = 1,. . . , m label its columns. Then, we have for all i = 1,. . . ,m,

(3.2)

for the ith component of 7(f). Thus, equation (3.2) is formally identical to equation (2.8) if we assume that their coefficients obey the following relation: m 9,r,kl = c

Tij,kl,

for all i, k, 1 = 1,. . . , m.

j=l

(3.3)

Similarly, given equation (2.8), we can put them in the form of equation (3.2), by letting, for instance, % kl foralli,j,k,1=1,..., m. Tij,kl = A (3.4) m’ The basic properties of BM which we need in our study are that they preserve probability measures, and that, in given circumstances, they make the entropy of the system increase by finite amounts [6]. Given a discrete probability measure F = {FI, . . . , F,.} (e.g., for F E Qz, we have T = m2), its entropy is defined as

S(F) = -

2 Fi log Fi. i=l

(3.5)

For any function g : [0, co) 3 W, such that g(0) = 0, we define

with ai in an interval I c W, for every i. Then, the following theorem holds [7, pp. 27,281. THEOREM 3.1. Let a, b be in P, for I an interval in JR,and P the Cartesian product of m copies of I. Then a = Tb, with T a doubly stochastic matrix, if and only if C, (a) L C,(b) for all concave functions g : I -+ W.

with g(z) = 2 logz, and g is concave, the action of a doubly Because S(F) = -C,(F), stochastic operator on F does not decrease its entropy. Moreover, given f E &, we have that S(f @ f) = 2S(f), and that S(E(f @ f)) = S(f @ f)/2. Thus, the entropy of the system has not decreased if the action of the sampling map is followed immediately by the action of the conditional expectation. This is the case when f @ f is a fixed point of the matrix T, in the definition of a BM. If f @If is not fixed for T, T(f @ f) is not necessarily a product probability, and we have S(T(f@f)) 2 2S(f), b ut we cannot say that S(E(T(f @ f))) 2 S(f). For that to be the case, we need stronger properties, as in the following result [5].

Boltzmann-Like KineticModels PROPOSITION 3.2.

Consider a doubly stochastic

59

matrix T, whose elements obey

Tzj,kl = Tji,kr.

(3.6)

The associated BM of Definition 3.3 does not decrease the entropy This result can be further strengthened,

in some circumstances,

by the following [S, p. 5111.

Let T be a doubly stochastic matrix on Q2, such that equation (3.6) holds. Let T” be the adjoint of T, and assume that 1 is a simple eigenvalue of TT*, separated by a gap A > 0

LEMMA 3.3.

from the second largest eigenvalue. Then.

If F E Qg, define G E &-J through Gij = CT,_,

S(G) 2 S(F) +

Tij,klFkr.

+llF- e/l:,

where 11. 112is the L2 norm on Rm2, defined by the scalar product (F, G) = Ci,j FijGij, e = (1,. . . , 1}/m2 is the uniform distribution Therefore, S is a strict

(3.7) and

in Qz.

Lyapunov function for the cases satisfying the hypothesis of Lemmas 3.2

and 3.3.

4. RESULTS Here, we give the list of our original results. Their proofs are to be found in Section 5. Consider any approximate l-population model of form equation (2.81, on llP. be defined by equation (2.7), and &;E be defined by equation (2.4), where $J:, obeys

PROPOSITION 4.1.

Let @i$

Restrictions

1 and 2. The associated m x m2 matrix, M = { @i,kl} say, has the following properties.

(a) M is nonnegative. (b) cz1 (C) c&

@i,kl = 1. @Pi,kl= m.

Conversely, for every m x m2 matrix M satisfying conditions (a), (b), and (c), there is a function $J obeying Restrictions 1 and 2, which yields a model of form equation (2.8) on Rm, through relations

(2.4) and (2.7).

PROPOSITION 4.2. Any approximate p-population model of form equation (2.6), with m subdivisions of the interval I, is a given transformation of Iwmp which can be mapped into an equivalent l-population model of form equation (2.8), with mp subdivisions of I. The converse is also true.

It appears that an approximate,

conservative, p-population model can be viewed as a special

case of equation (2.8), with mp subdivisions of I, where all P,.,~~ vanish except (possibly) for m2p of them. The opposite is not true. Indeed, the matrix M = {@i,kl} of a conservative ppopulation model splits into a direct sum of p blocks of size m x m2p, while there are lpopulation models, with mp subdivisions of I with no vanishing coefficients. This is a consequence of the conservative condition for many populations, i.e., that production of individuals i occurs only through interactions involving other individuals i. This result does not depend on our technical assumptions, but only on the form of the approximate models equation (2.6) and equation (2.8). Consider an approximate l-population model defined as in Proposition 4.1, and denote by r such a mapping of W”. For 7 to be a BM, it is necessary and sufficient that 0 5 At I Ate = mini,r l/(n(l - pLi,ir)). PROPOSITION 4.3.

Clearly, Ate is not smaller than l/r], which is a positive quantity. Let Q c Wm be the simplex of probability measures on the sample space s1 = , m}. Consider the ODE system equation (2.4) with p = 1, where n > 0, and Qi, satisfies Restriction 1. The initial value problem equation (2.41, with initial data Ui(O>) = Vi01 E Q, *w25:z-c

PROPOSITION 4.4. {l,...

60

L. F&JND~NI

has a unique solution in Q for all times t L 0, f : [O,oo) + Q say Also, f is nonnegative and normalized, and Euler’s scheme equation (2.8) converge to it as At is made smaller and smaller. THEOREM

4.5.

Consider

(2.4) with p = 1, and assume that $i,

equation

obeys Restrictions

1

and 2. Denote by r the associated BM, equation (2.8), with time step At, E (0, Ate), and Ato defined in Proposition 4.3. T may be replaced by another BM, r’ say, which is an Euler’s scheme for the same T’ = {T&r}

ODE problem,

but has smaller or equal time step,

T&,kr = Tjli,kr = Tj+ Moreover,

and its interaction

operator

obeys

property

(4.1) is common

9

foralli,j,k,l=l,...,

to all maps obtained

m.

from r’ by reducing

(4.1) its time step.

Thus, the results of Lemma 3.3 can be used, when the conditions of Theorem 4.5 are satisfied. This is not always the case. Indeed, not all BM can be replaced by others with same dynamics and interaction operator obeying equation (4.1) (see Example 6.2). The BM we consider here, however, satisfies those conditions, their special feature being that they are derived as approximations of evolution equations. We now exploit this fact. COROLLARY 4.6. Consider an approximate l-population model of form equation (2.8), defined as in Proposition 4.1, and denote by r the relevant map of Rm. Assume that all terms @i,kr are positive, and that the dynamics of r can be obtained from a doubly stochastic matrix T say, with Tij,kl = Tji,kr for all i, j, k, 1 = 1,. . . , m. r is a BM with a unique fixed point in the simplex f = (l/m, l/m,. . . , l/m). Such a of probability measures & C IP, the uniform distribution stationary state is stable, and a global attractor in & under the dynamics of the map r. If equation (2.4) is considered independently of its connection with equation (1.3)) Restriction 1 may be relaxed to produce other models, as in the following. COROLLARY 4.7. Consider the ODE system equation (2.4) with p = 1. Assume that $J!, is a positive, integrable function $ : [-1, 113 -+ R+, such that equation (1.6), and Restriction 2 are satisfied. This model has a unique stationary state in &: the uniform distribution 1. Also, 3 is stable for the associated

dynamics.

However, Restriction 1 may still be retained because of the following results. PROPOSITION 4.8. The class of functions

of form (2.1) is dense in the set

{Ic, E Co (r-1,11”) according

: 1c,2 0)

to the uniform norm.

COROLLARY 4.9. Consider q > 0 for all a, b, and that initial value problem, with W+ u {0},for some T > 0, Let {Izc}y be the partition

a quasi-conservative model of form equation (1.3) such that nab = the functions $66 are continuous in [-1, 113. Assume the associated fi(O,u) = f!(u), h as a solution f = {fi); say, _f’i: [O,T) x [-1, l] 4 such that f(t, .) and $$ are integrable over [-1, l] for all t E [O,T). of [-1, l] introduced in Section 2, and

Fe(t) =

s I

2

f&7 u> d%

for all t E [O,T).

Take any E > 0. There is an m = m(T, E) E N such that II{F~(t)}

- {fG(t)}\I1

I &,

for allt E [O,T),

(4.2)

where {f$$} is the solution of the initial value problem equation (2.41, with fi,=(O) = F+,,(O), obtained from suitable functions $ib satisfying Restriction 1. Thus, if the functions $ib fall within the limits of Restriction 2, the corresponding ODE problem can be obtained as the limit for At 3 0 of a given BM. In this way, we have enslaved the evolution of our model, with continuous functions @ia, to the dynamics of a BM, for the time interval [0, T) .

Boltzmann-Like KineticModels

61

5. PROOFS Here we give the proofs of the results

announced

in Section

4.

PROOF OF PROPOSITION 4.1. Conditions (a) and (b) are satisfied by the coefficient matrix of equation (2.8), as seen from its definition equation (2.7), if equation (1.6) is satisfied. This is the case if Restriction 1 holds. Condition (c) holds as well. Indeed, for every u E [-1, 11,we have

s

$(v, w, u) dv dw =

IkXIl

and Restriction

2 implies

m

c

m b,kl

c

!m(u)du = m,

=

(5.2)

k&l

from which we get m

/m

c

@i,kl

=m+dt

\

(,,Pi,,,

-_!

(5.3)

.

k,l=l

It remains

to show that

given

M = {@i,kl}, for which conditions

a matrix

l-population model of form equation (2.8), satisfying Restrictions of it. To do that, choose 11= l/At (hence, @i&l = pi&l), let akZ(‘ZLi)

=

for all i,k,l

;@i.tl,

(a),(b),(c)

hold,

1 and 2 can be constructed

,m,

= 1,.

a

out

(5.4)

and define pkI(u)

=

(5.5)

&kl(ui)%(u). i=l

Now, the question

is whether

the function

ti(v,%u)

=

2

~kl(u)xk(v)xl(w)

(5.6)

k,l=l

satisfies integral

the necessary prerequisites. with respect to u yields

J:

$(v, w, u)

Indeed,

du =

2

it does. In the first place, II, is nonnegative,

xk(‘u)xl(w)

2

Et,=,

Xk(V)xr(w)

is identically

s, (-1

lla ti(v,w,u)dvdw

which is the condition

=

=

1, and property

(b) holds.

Furthermore,

eXi(U)y

@i,kl 1-1

i=l

required

(5.7)

($+,kL)

by Restriction

1,

i=l

k,l=l

because

and its

2

lea dvdw

=

we have

2,

115.8)

k&l

2. This completes

the proof.

I

PROOF OF PROPOSITION 4.2. The pairs of indices (i,z), (a, k), and (b, 1) in equation (2.6), which constitutes the rectangular array {{(h, q)}ph,l}!&, may be organized into a linear array, (1,. . . , mp}, with the substitutions (i, z) -+ r = (i - 1)m + 2, (a, k) --+ y = (u - 1)m + k, and (b, 1) -+ z = (b - l)m + 1. Then, equation (2.6) can be written as mp

w

(5.9)

L.

62

ROND~NI

where @r,yz = a:$. Thus, formally all our approximate models are equivalent to a l-population model with mp subdivisions of I. It remains to check that the corresponding matrices M = {ai,kl} have the correct properties. These are nonnegativity, the fact that columns sum up to 1, and that rows sum up to the number of subdivisions of I times the number of populations. Nonnegativity is obviously satisfied. The fact that columns sum up to 1 is a consequence of Restriction 1 and equations (2.5),(2.7), while Restrictions 1 and 2 combined with equation (2.7) ensure that the rows of M add up to mp. Hence, from Proposition 4.1, we conclude that a quasiconservative model with p populations and m subdivisions of I can be mapped into a l-population model, with mp subdivisions of I, and vice versa. I For the special case of conservative ppopulation as fjn+i)

am c

= f,‘“’ + g z=l

models, equation (5.9) can be explicitly written

,1

%r/tf$V)

[ y=(a-l)m+l

r=l,...,mp,

(5.10)

where the value of a is linked to that of T by a = [r/m]*, with r r * E m

if ?Z- E N, m

rn’

” {

integer part of

(

;

>

+l,

(5.11)

ifkeN.

PROOF OF PROPOSITION 4.3. Because of equation (2.4), and of the first identity in equation (2.5), we have for all i, 1 = 1, . . . , m,

hence,Ato=$n((1_~i,i1,9>

Zl>O.

Thus, for 0 I At I Ate, we have @i&l 2 0 and EE, @i,kl = 1. This means that the matrix T = (Tii,kl) defined by relation (3.4) is (column) stochastic. Furthermore, the same relations imply

Tij,ki = 1,

for all i,j = l,...,

m,

k,l=l

which is the row-stochasticity condition for T. Thus, T is doubly stochastic, and equation (2.8) defines a BM of order 2. However, taking At < 0 or At > Ate makes some coefficient @i&l < 0, i.e., some element of T necessarily negative. In that case, T is not doubly stochastic, and equation (2.8) is not a BM. I PROOF OF PROPOSITION4.4. The vector field of the ODE (2.4) is smooth, thus a local solution, f : [O,T) ---t Wm say, for some T > 0, exists. Moreover, such a solution can never cross the boundaries of a compact subset of Wm. Indeed, fi(t) > 0, for all i = 1,. . . , m, and CE”=,( g) = 0. Thus, existence and uniqueness of solutions extend to all times. Convergence of the Euler’s scheme is also guaranteed by the smoothness of the vector field. That is to say, that for every E > 0, every T E (0, oo), and all initial conditions in Q, there is a At’ > 0 such that Ilf(nAt) - f(n))(i < c, for all At 5 At’, and all n such that 0 5 nht 2 T. I PROOF OF THEOREM 4.5. Given a map of form (2.8), derived from equation (2.4), where $J:, satisfies Restrictions 1 and 2, we introduce a doubly stochastic matrix T(At) whose entries are defined by Tdi,kl (At) = @i&r/m. In this way, the dynamics of r can be cast into the form of equation (3.2). Then, we introduce a set of parameters {Kij,kl} through Tij,ii(At) = 1 + Kij,ijAt, Tij,kl(At) = Kii,klAtr

for all i,j = l,...,

m,

for all i, j, Ic,1 = 1,. . . , m,

such that (i, j) # (Ic, 1).

(5.12)

Boltzmann-Like

We have Kij,ij

< 0, for all i,j

The rows and columns of T(At)

Kinetic Models

63

and Kij,kl 2 0 whenever the pairs (i,j)

and (k,l)

are different.

all sum up to 1, which means that

m

c Kij,kl = c Kij,kl =

for all Ic,l = l,...,

0,

m,

&j-l

(5.13)

m

foralli,j=l,...,

0,

m.

k,l=l

If we fix the constants Kij,kl, and let At vary, we obtain a doubly stochastic matrix T(At) At 5 Ato. Furthermore, if At is small compared to 1, the diagonal terms of T(At) form Tij,ij) are of order o(l), while all the others are of order o(At/m). The contribution

to f!“+”

given by the product (fp’fi(n)),

for all

(i.e., those of

in equation (3.2), can be cumulated

with that given by (fr(nifrln)) 1and the two may be symmetrized by replacing Tij,kl and Tij,lk with
=
=

Tij,/cl+ Tij,Zk 2

.

Thus, the dynamics of equation (3.2) is identical to that of (5.14)

whose coefficients are symmetric under the exchange k H 1. We express this by saying that such is symmetric by columns. In particular, we have a matrix {
for all i = 1,. . . ,m,

Eij,ij = Jij,ji =

k(1 + (Kij,ij


-

=


=

1

2

(Kij,kl

+ Kij,ji)At)

+ Kij,lk)

At,

,

for all i # j, for all pairs (i,j)

# (k,E).

Thus, we assume that the constants Kij,kl have been fixed, and that T(At) satisfies equation (5.15). We show that T(At) can be replaced by a doubly stochastic T’(At) whose entries obey T’ij,kl = T’ji,kl (i.e., are symmetric by rows), and preserve the sums (3.3), provided that At is small enough. One way to proceed is to take one entry of T, Trs,yZ say, such that Trs,yz 5 Tij,kl for all i,j, k, I = 1,. . . , m. For At small enough, T,.S,yz can only be an off-diagonal element of T, such that (T, s) # (z, 7~). TO see this, take any At’ 5 Ato, fix the values of all terms Kij,kl, and observe whether the above statement is satisfied or not. If it is not, keep the same values Kij,kl,, and reduce At, until the largest term of form Ktj,klAt is smaller than the the smallest term of form l/2 + Kij,ijAt. This happens for some finite value of At. Moreover, the resulting interaction operator determines a map which is Euler’s scheme for the original ODE problem, where the time step has been reduced. Given TrS,yz, if T = s, there is no need to change the matrix T, because Trs,yZ = Tsr,yZ. Thus, we move to another element of T. Let this one not be greater than any other, except (possibly) the entry previously considered. If T # s, we have T,.s,yZ 5 Tsr,yz by definition. We replace Tssr,yz with TA.,yz = Trs,yZ, and we add the difference lrs,yz = Tsr,yZ - Trs,yZ to the entry Tss,yZ, i.e., we The resulting matrix T is nonnegative, and its column sums are 1, define TglSlyz= Tss,yz + &,z. because the sum (3.3) defining @s,yz has been preserved (which implies that the dynamics of r has not changed either). However, the sum of the elements of row (ST) is now 1 - IrS,yz, while the sum of the elements of row (ss) is 1 + Zrs,yz. In other words, the resulting T’ is not doubly stochastic,

unless Irs,yz = 0. To fix this, subtract Irs,yz from T,,,,,,

i.e., let T,‘, ss = T,,,,, - ‘!rs,yz,

64

L. RONDONI

and add the same amount to T,,,,,

i.e., let TLr,ss = Tsr,ss + lrs,yz. At this point, we have a

matrix whose entries Trs,yz and Tsr,yt are equal, and which is doubly stochastic if Tss,+ - Irs,yz is nonnegative. Now, Tss+ has form 1 + K,,,,,At, and lrs,yZ is proportional to At. Therefore, there is a finite value of At below which the operation above yields a doubly stochastic matrix. After this, we move to the next smallest element Tij,klwhich has not already been symmetrized, and we repeat the process above, taking care that nothing is done to the terms that have already been symmetrized. For instance, assume that we change Tji,kl in order to match it with Tij,kl. Then, to compensate in the sum (3.3), we add E.’ y,kl to Tjj,kl, and to compensate the row sums, we change Tjj,jj and Tji,jj accordingly. However, Tji,jj could have already been symmetrized, and we do not want to change it again. Thus, we adjust our scheme by adding r!ij,kl to a different Tjs,kl> instead of Tjj,kl. In turn, this requires US to change different terms, i.e., Tj,,j, and Tji,j,. At a certain point, we may have to consider one element of T such that there are no terms where changes could be made without involving other elements that have already been symmetrized. This means that Tij,kl is already equal to Tji,kl. In fact, columns and rows sum up to 1, thus fixing one of their elements, when all the others have been fixed. At every step, the range of values of At which allows us to perform all the symmetrization operations may be reduced. However, there are only a finite number of steps to be performed, thus, the range remains finite at the end of the process. The result of this is that the original map r has been replaced by another one, r’ say, defined by (5.16)

which is the Euler’s scheme for the same ODE problem, but has smaller time step. Moreover, the elements of the interaction operator of 7’ have the important symmetry property T’ij,kl = T’ji,kl = T’ji,lk.

(5.17)

Such a property remains if we further reduce the value of the time step. Indeed, the symmetrization process has merely added or subtracted from each Tij,kl a quantity proportional to At. Thus, the elements of T’ can be written as the tij,kl in equation (5.15) (with different Kij,kl), whose symmetry is not altered by changes in At. In particular, reducing At preserves the doubly stochssticity of the matrix, while increasing it above a certain value does not. I PROOF OF COROLLARY 4.6. That the uniform distribution f^ is a fixed point follows from the fact that f@f = (l/m2,. , . , l/m2) is a uniform distribution on Wm2. Thus, I@ f^ is a fixed point of T. Then, the associated BM acts on f as follows: (5.18) That the fixed point f^ is unique in Q follows from the the uniqueness of p@ f^ as a fixed point of T in Qs, and from the fact that T increases the entropy of any other input probability. Indeed, we note that @i&l > 0, for all i, k, 2 = 1,...,m, implies that all entries of T, constructed as in the previous proof, are positive. Then, TT* is a positive matrix, and Perron-Frobenius theory implies that TT* has a spectral gap A > 0. Thus, take f # f^. Then, f @ f is not fixed for T, and the entropies go like

Wf) = W @f) < S(T(f @f)) I WWW @f)>) = 2W’h

(5.19)

where the last inequality is a consequence of the property Tij,kl = Tji,kl (see Proposition 3.2 and Lemma 3.3). Because the entropy of f’ is different from that of f, we have f # f’, i.e., f is not fixed for the BM.

Boltzmann-Like Kinetic Models Stability

of the uniform

of the uniform neighbourhood

distribution,

under

65

the BM, is explained

as follows.

The

entropy

distribution _/ is the maximum value that can be attained in Q. Consider the of f^ made of all points f with S(f) > S(f) -,E and E > 0. By varying E, we can

make this set contained in any other fixed neighbourhood of f. Because entropy cannot decrease under the action of our BM, for every neighbourhood of f^, Uf say, there is another one, centered in f^, such that the iterations of its points cannot leave Uf. Thus, f^ is stable. To get the convergence part of our result, we summarise the proof of Theorem that the increase in entropy during In fact, the entropy is a bounded

uniqueness of f^ as a fixed point of our BM, the previous fact implies that The fixed point f is then unique, stable, and a global attractor in Q. 4.7. If $J > 0, all coefficients

PROOF OF COROLLARY

@+l are positive.

stationary state follows from the mere fact that a fixed point vector field of equation (2.4), and vice versa. Indeed, if

f

=

fn+l)

G(f),

:=

f(“)

-+ f^ as n -+ co. I

Then,

uniqueness

of the

of the BM (2.8) is a zero of the

f’“’ + G (f’.“‘)At

is an ODE system with associated Euler’s scheme, both the continuous are stationary for G = 0, while they are not if G # 0. Stability

3 in [5]. Observe

the evolution of any f E Q, f # f^ must go to zero as n -+ co. quantity in Q. Because of Lemma 3.3, and because of the

of the fixed point follows from the fact that the evolution

and the discrete

systems

of the ODE problem

decrease the entropy. In fact, if there is a decrease of entropy in the with given initial condition, there are two instants of time tz > tl > 0 and c > 0. But in the interval of time [tr, t2], the solution of equation as accurately as desired, by choosing a sufficiently small time step At

solution

cannot

of equation

(2.4),

such that S(t2) = S(tl) -c (2.4) can be approximated (Proposition 4.4). Thus, a

decrease of entropy would also be observed in the evolution of all the discrete maps with a time step smaller than a given At, > 0. This cannot be the case, as the discrete schemes are BM. 1 PROOF OF PROPOSITION 4.8. Given described in Section 2, and define

a $J E Ce([-1,

s

113), introduce

the same partition

$(w, w, u) dwdw.

of [- 1, l]

(5.20)

Ikxh

Then,

we can write

for all 21,W E Ik x Il,

(5.21)

for some ?& E Ik and 61 E 11. But the value E, = maXi=i,...,m(Ui+l - Ui) converges to zero as m gets larger and larger. Therefore, given the continuity of $J in [-1, 113, for every E > 0, we can find one m E N such that for all (w, w, U) E [-1, 113. This completes PROOF

(5.22)

the proof.

OF COROLLARY

I 4.9. Because

of our assumptions

w, w,

I

27p2%n, m

u) -

and of Proposition

q!&@)(w, w, u))fa(t,

4.8, we can write

'U)fb(t, w) dw dw du

(5.23)

L.F~~NDoNI

66

where gab i’(m) is defined as in the previous proof, through equations (5.20),(5.21), for all i,a,b, and e, can be taken smaller and smaller as m grows. Then IFz(t)

- fg(t)j

I $aB,t,

for all i,5,

(5.24)

and I 2P3??GlJ.

(5.25)

Thus, for given T and E, we only need to take m sufficiently large.

I

IIVGHt)

- U$Ht)(ll

6. EXAMPLES The following examples address, in the simplest possible way, the following three questions. Question 1: Are there chemical systems expressible in terms of quasi-conservative models? Answer 1: Yes. Question 2: Can all doubly stochastic matrices be symmetrized by rows? Answer 2: No. Question 3: What sort of functions $ obey Restrictions 1 and 2? Answer 3: Among the others, all functions which are probability densities in all their variables v, w, 21.

EXAMPLE 6.1.(A simple BM) Consider a model with two species only. The relevant sample space is R = {1,2},

and a probability measure on it has the form f = {fi, fi},

with f2 = 1 - fi.

The

probabilities of picking two independent individuals from the system are (fi)2, if the individuals are both of species 1, fif2 = Fiji, if they belong to different species, and (f~)~, if they are of species 2. In particular, let T be of the form

1-a

P

P

T=

0

0

1-P

2 1-P 2 0

2 1-P

0

it

1

--

P--

P

1-P

0

(6.1)

The associated BM is given by

fin+11=fi’“’- p (fp2 - fi’“‘fq

fp+l) =fp+p(fp)a_

,

(6.2)

fpfp).

If we let p = KAt, equation (6.2) is Euler’s scheme for the evolution of the concentrations fi and f2 of the chemicals A and B, in the reaction 2A = A+B, with reaction rate K. equation (6.2) is a BM if p I l/2, and they can be derived from a quasi-conservative model of form (1.15), if we have p = 2, m = 1, constant encounter rate v, and

l-u ti,111(vJ,u) = 2 ' 1+u ~:2(vw)

=

-

4

’

ti:1(wv) = f, 1cIP2(w+)

=

1-V 4'

-

with Y E (0,l). In this case, we have p= VVAt. More complex models of chemical reactions can be thought of and produced with this technique. In particular, any convex combination of maps like these yields another such BM.

Boltzmann-Like EXAMPLE

there

6.2.

(On the symmetry

is no symmetry

Kinetic Models

67

a case where T is doubly

by rows) Consider

stochastic,

T=

If f # (l/2,1/2),

but

by rows, e.g.,

the entropy

(6.4)

of the relevant

2W) = w

@f)

probabilities

obeys

< WYf 8 f)) <

s (f’) + s (9’)

9

where f’ and g’ are the marginals of T(f 8 f). Nonetheless, 7(f) = f, and S(T(~)) = S(f). Moreover, we find that not all doubly stochastic matrices can be symmetrized by rows, in such a way that

they remain

doubly

stochastic,

the only way to get the same dynamics, is to fill such rows with zeroes.

and that the associated and equality

BM are not changed.

of the second with the third

Hence, we see the importance

of Theorem

Indeed,

row of our T,

4.5 for the subsequent

results. EXAMPLE 6.3. Consider a l-population model whose function $ : [-1, 113 --) R+ U (0) is a probability distribution in all of its three arguments, i.e., &I(V) = $(., w, .) dw and dVz(w) = measures on [-1, 11. Then, if $ is simple, it satisfies Restrictions 1 $J(.>.? w) dw are probability

and 2, as $(w, W, U) du, = 1 The corresponding are probability

l(6.5)

s [-1,112

quasi-conservative

distributions

$(v, w, u) dv dw = 2.

and

s [-1311

models are those for which both Cz=,

&,

and CL-,

T&,

in all their variables.

REFERENCES 1. E. Jager and L. Segel, On the distribution of dominance in a population of interacting anonymous organisms, SIAM J. Appl. Math. 52, 1442-1468 (1992). 2. L. Arlotti and N. Bellomo, On a new model of population dynamics with stochastic interaction, nansp. Theory Stat. Phys. 24, 431-443 (1995). 3. L. Arlotti and N. Bellomo, Solution of a new class of nonlinear kinetic models of population dynamics, Appl. Math. Lett. 9 (2), 65-70 (1996). 4. L. Preziosi, From population dynamics to modelling the competition between tumors and immune system, Mathl. Comput. Modelling 23 (6), 135-152 (1996). 5. R.F. Streater, Convergence of the iterated Boltzmann map, Publ. Research Inst. Math. Sciences, Kyoto University 20, 913-927 (1984). 6. RF. Streater, Statistical Dynamics, Imperial College Press, London, (1995). 7. P.M. Alberti and A. Uhlmann, Stochasticity and Partial Order, D. Reidel, Dordrecht, (1982).