text:
# Untitled Document

Exponentiation and Euler measure

James Propp

University of Wisconsin

April 29, 2001

Revised May 29, 2002

dedicated to the memory of Gian-Carlo Rota

Abstract: Two of the pillars of combinatorics are the notion of choosing an arbitrary subset of a set with $`n`$ elements (which can be done in $`2^n`$ ways), and the notion of choosing a $`k`$-element subset of a set with $`n`$ elements (which can be done in $`\left(\genfrac{}{}{0pt}{}{n}{k}\right)`$ ways). In this article I sketch the beginnings of a theory that would import these notions into the category of hedral sets in the sense of Morelli and the category of polyhedral sets in the sense of Schanuel. Both of these theories can be viewed as extensions of the theory of finite sets and mappings between finite sets, with the concept of cardinality being replaced by the more general notion of Euler measure (sometimes called combinatorial Euler characteristic). I prove a “functoriality” theorem (Theorem 1) for subset-selection in the context of polyhedral sets, which provides quasi-combinatorial interpretations of assertions such as $`2^1=\frac{1}{2}`$ and $`\left(\genfrac{}{}{0pt}{}{1/2}{2}\right)=\frac{1}{8}`$. Furthermore, the operation of forming a power set can be viewed as a special case of the operation of forming the set of all mappings from one set to another; I conclude the article with a polyhedral analogue of the set of all mappings between two finite sets, and a restrictive but suggestive result (Theorem 2) that offers a hint of what a general theory of exponentiation in the polyhedral category might look like. (Other glimpses into the theory may be found in \[Pro\].)

1. Introduction

In the work of Hadwiger \[Ha1\] \[Ha2\], Lenz \[Len\], McMullen \[McM\], Rota \[KR\], Chen \[Ch1\] \[Ch2\], and others, the Euler characteristic of polytopes in $`𝐑^n`$ is extended to a valuation $`\chi `$ (the Euler measure) on the algebra of sets generated by these polytopes. Included in this algebra of sets (called “hedral sets” by Morelli \[Mor\] and “polyhedral sets” by Schanuel \[Sch\]) are the finite subsets of $`𝐑^n`$, for which the notion of Euler measure and the notion of cardinality coincide. This suggests the viewpoint (forcefully advocated by Schanuel) that we view polyhedral sets and their $`𝐙`$-valued Euler measures as a generalization of finite sets and their $`𝐍`$-valued cardinalities. Indeed, several of the key pillars on which combinatorics rests — concepts of union, Cartesian product and selection of a subset of fixed size — extend to the polyhedral setting straightforwardly, and this extension provides a natural interpretation for assertions like $`\left(\genfrac{}{}{0pt}{}{2}{3}\right)=4`$ (see section 2).

What else constitutes the foundation of combinatorics? Many combinatorialists would point to the notion of the set of all maps (or injections, or surjections, or bijections) from one finite set to another; see for instance the “twelve-fold way” of Rota (so dubbed by Spencer) as described in Stanley’s book \[Sta\]. As a special case, the set of all maps from a finite set to the set $`\{0,1\}`$ corresponds naturally with the set of all subsets of the given set.

I will show that a polyhedral analogue of this foundation can be built, using Morelli’s notion of hedral maps (piecewise-constant maps with finite range whose level sets are polyhedral sets) and Schanuel’s more inclusive notion of polyhedral maps (maps whose graphs are polyhedral sets). Since in neither case is the set of all such maps from $`A`$ to $`B`$ itself a polyhedral set (except when $`A`$ and $`B`$ are finite or $`A`$ or $`B`$ is empty), we cannot simply take Euler measure in the standard way; but, at least in the case where the polyhedral set $`A`$ is a subset of $`𝐑`$, we can divide this large set of maps into finite-dimensional strata by limiting the number of points of discontinuity/non-differentiability and use this stratification to define a power series in $`t`$ whose value at $`t=1`$ “should” be the Euler measure of the set of maps. In many cases the series possesses an analytic continuation with a well-defined value at $`t=1`$, which I call the regularized Euler measure.

For instance, let $`A=(0,1)`$ and $`B=\{0,1\}`$, with $`\chi (A)=1`$ and $`\chi (B)=2`$; then we will see in section 5 that $`B^A`$, defined as the set of all polyhedral maps from $`A`$ to $`B`$, has regularized Euler measure $`\frac{1}{2}`$ , and the set of all two-element subsets of $`B^A`$ has regularized Euler measure $`\frac{1}{8}`$. The preceding assertions correspond, in some way not yet formalized, to the algebraic assertions that $`\frac{1}{2}=2^1=\chi (B)^{\chi (A)}`$ and $`\frac{1}{8}=\left(\genfrac{}{}{0pt}{}{1/2}{2}\right)`$ (under the interpretation of $`\left(\genfrac{}{}{0pt}{}{x}{k}\right)`$ as a polynomial of degree $`k`$ in $`x𝐑`$).

A different theory of exponentiation of polyhedral sets arises from focusing on power sets, and making the restriction that one looks at only the finite subsets of a given polyhedral set. Since this theory is simpler and has been developed further, it is discussed first, in sections 3 and 4. Then, in section 5, I treat the set of maps from one polyhedral set to another, under both Morelli’s condition of hedrality and Schanuel’s weaker condition of polyhedrality.

I should warn the reader not to expect too much; the two theorems proved in this article are really more like examples than like theorems — and what is more, I do not give a clear statement of what sorts of definitions and theorems these examples are special cases of. It is my hope that at least some readers will see that there are issues here that are worthy of being pursued, and that there is some “spiritual” satisfaction in the prospect of a theory that might be a model of the rational numbers in the same way that ordinary combinatorics is a model of the non-negative integers.

I thank Scott Axelrod, John Baez, Beifang Chen, Timothy Chow, Ezra Getzler, Greg Kuperberg, Michael Larsen, Ayelet Lindenstrauss, Haynes Miller, Lauren Rose, and Dylan Thurston for helpful conversations. I also thank Joseph Kung for soliciting this article and putting this special issue together, and the anonymous referee for numerous helpful suggestions.

I acknowledge a special debt to Gian-Carlo Rota, who persistently encouraged the investigations described here, despite their inconclusive status both in 1996 (when I described my work to him) and now (after years of neglect). I also want to acknowledge the formative role he played in my career as the author of the beautiful book Finite Operator Calculus, which was part of what inspired my earliest explorations of combinatorics. His sagacity and warmth are missed by me and many others.

2. Background

Following Schanuel, we define a polyhedral set in $`𝐑^n`$ as one that can be defined by a Boolean predicate whose primitive terms are linear equations and inequalities in $`n`$ variables. Equivalently, the algebra of polyhedral sets is the algebra generated by the closed half-planes in $`𝐑^n`$ under intersection, union and complementation. Note that polyhedral sets can be unbounded; for instance, $`\{(x,y):|x|+|y|1\}\{(x,0):x>2\}\{(5,5)\}`$ is a polyhedral subset of $`𝐑^2`$.

There exists a unique function $`\chi ()`$ from the polyhedral algebra to the integers such that

$$\chi (AB)=\chi (A)+\chi (B)\chi (AB)$$
(1)
for all polyhedral sets $`A,B`$ and such that

$$\chi (A)=(1)^k$$

if the polyhedral set $`A`$ is homeomorphic to $`𝐑^k`$ or equivalently is homeomorphic to the interior of a $`k`$-dimensional disk (see e.g. \[Ch2\]). For instance, the polyhedral set described at the end of the preceding paragraph has $`\chi `$-value $`11+1=1`$. The function $`\chi `$, hereafter called simply Euler measure, is invariant under homeomorphism, and in the case of compact polyhedral sets, it coincides with the alternating sum of the Betti numbers. However, unlike the notion of Euler characteristic favored by topologists, it is not a homotopy invariant (e.g., $`\chi `$ of an open interval in $`𝐑`$ is $`1`$ whereas $`\chi `$ of a closed interval is $`+1`$). Property (1) extends naturally to an inclusion-exclusion principle:

$$\chi (A_1A_2\mathrm{}A_m)=$$

$$\underset{1im}{}\chi (A_i)\underset{1i<jm}{}\chi (A_iA_j)+\mathrm{}\pm \chi (A_1A_2\mathrm{}A_m).$$

For polyhedral sets $`A𝐑^m`$ and $`B𝐑^n`$, one obtains a Cartesian product $`A\times B`$ in $`𝐑^{m+n}`$ in the obvious way, and it can easily be shown that $`\chi (A\times B)=\chi (A)\chi (B)`$. Less trivially, given a polyhedral set $`A𝐑^n`$ and an integer $`k0`$, there is a way of defining a polyhedral set $`\left(\genfrac{}{}{0pt}{}{A}{k}\right)(𝐑^n)^k`$ that contains exactly one point in each free orbit in the $`k`$-fold Cartesian product $`A^k`$ under the permutation action of the symmetric group $`S_k`$ (where the free orbits are those of cardinality $`k!`$, i.e., the orbits of $`n`$-tuples of distinct real numbers); this construction is well-defined up to polyhedral isomorphism (defined below) and has the property that $`\chi (\left(\genfrac{}{}{0pt}{}{A}{k}\right))=\left(\genfrac{}{}{0pt}{}{\chi (A)}{k}\right)`$. Since the operation $`A\left(\genfrac{}{}{0pt}{}{A}{k}\right)`$ (mapping Schanuel’s category of polyhedral sets to itself) seems to be new, it is necessary to amplify on the preceding sentences in the following paragraphs.

Following Schanuel, we define a polyhedral map from a polyhedral set $`A𝐑^m`$ to a polyhedral set $`B𝐑^n`$ as a map $`AB`$ whose graph in $`𝐑^{m+n}`$ is itself a polyhedral set; that is, a polyhedral map is a piecewise-affine map that need not be continuous. A polyhedral isomorphism is a bijective polyhedral map, such as the map

$$f(x)=\{\begin{array}{cc}x1\hfill & \text{if }1<x<0\hfill \\ x+5\hfill & \text{if }x=0\hfill \\ 2x+1\hfill & \text{if }0<x<1\hfill \end{array}$$

between the interval $`(1,1)𝐑`$ and the set $`(2,1)(1,3)\{5\}𝐑`$. It can be shown that two polyhedrally isomorphic polyhedral sets have the same Euler measure. On the other hand, following Morelli, we define a hedral map from a polyhedral set $`A𝐑^m`$ to a polyhedral set $`B𝐑^n`$ to be a function $`f:AB`$ with finite range whose level-sets are all polyhedral.

In both cases, we identify a function $`f:𝐑^m𝐑^n`$ with its graph in $`𝐑^{m+n}`$.

Note that if $`A`$ is finite, the hedral and polyhedral maps from $`A`$ to $`B`$ are just the set-theoretic maps, so that, in the case where both $`A`$ and $`B`$ are finite, the set of hedral or polyhedral maps from $`A`$ to $`B`$ has cardinality $`|B|^{|A|}`$ (or, in other terms, Euler measure $`\chi (B)^{\chi (A)}`$). If $`A`$ is infinite, then most set-theoretic maps from $`A`$ to $`B`$ are neither hedral nor polyhedral.

Also note that when $`B`$ is finite, the the set of hedral maps from $`A`$ to $`B`$ coincides with the set of polyhedral maps from $`A`$ to $`B`$; when $`B`$ is infinite, every hedral map is polyhedral but not conversely.

We now define a polyhedral ordering on a polyhedral set $`A𝐑^n`$ as a (total) linear ordering whose graph in $`A\times A`$ is polyhedral; an example of such an ordering is the one induced by lexicographic ordering of $`n`$-tuples in $`𝐑^n`$. Henceforth, we restrict attention to lexicographic ordering. We form $`\left(\genfrac{}{}{0pt}{}{A}{k}\right)A^k`$ by selecting, from each free orbit in the $`k`$-fold Cartesian product $`A^k`$ under the permutation action of the symmetric group $`S_k`$, the one $`k`$-tuple that is lexicographically least in the induced ordering on $`A^k`$. That is, $`\left(\genfrac{}{}{0pt}{}{A}{k}\right)`$ consists of all $`k`$-tuples $`(v_1,v_2,\mathrm{},v_k)A^k`$ with $`v_1<v_2<\mathrm{}<v_k`$ in the lexicographic sense. (Note that any polyhedral set that contains exactly one point in each free orbit in $`A^k`$ under the action of $`S_k`$ is canonically equivalent to the one defined above, via a polyhedral isomorphism.)

Take for instance $`A=(0,1)(2,3)𝐑`$, with $`\chi (A)=2`$. Then we find that $`\left(\genfrac{}{}{0pt}{}{A}{3}\right)`$ is a union of four sets: $`\{(x,y,z):0<x<y<z<1\}`$, $`\{(x,y,z):0<x<y<1,2<z<3\}`$, $`\{(x,y,z):0<x<1,2<y<z<3\}`$ and $`\{(x,y,z):2<x<y<z<3\}`$. Each of these is homeomorphic to the interior of a 3-ball, so $`\chi (\left(\genfrac{}{}{0pt}{}{A}{3}\right))=(1)^3+(1)^3+(1)^3+(1)^3=4`$. Note that this is also the value of $`\left(\genfrac{}{}{0pt}{}{\chi (A)}{3}\right)`$ where for each fixed $`k0`$ and all $`x𝐑`$ we define $`\left(\genfrac{}{}{0pt}{}{x}{k}\right)`$ as $`x(x1)\mathrm{}(xk+1)/k!`$.

To show that $`\chi (\left(\genfrac{}{}{0pt}{}{A}{k}\right))=\left(\genfrac{}{}{0pt}{}{\chi (A)}{k}\right)`$ for all polyhedral sets $`A`$ and all integers $`k0`$, it suffices to show that the polyhedral set $`\{(v_1,v_2,\mathrm{},v_k)A^k:v_iv_j\text{for all }ij\}`$ has Euler measure $`\chi (A)(\chi (A)1)\mathrm{}(\chi (A)k+1)`$, since this larger polyhedral set is the disjoint union of $`k!`$ polyhedrally-isomorphic copies of $`\left(\genfrac{}{}{0pt}{}{A}{k}\right)`$. One way to prove this is by iterated application of a Fubini theorem (Theorem 3 in Morelli): If a polyhedral set $`AB\times C`$ has the property that every $`C`$-fiber has Euler measure $`m`$, then $`A`$ has Euler measure $`\chi (B)m`$. For example, in the case $`k=2`$, it is clear that $`\{(v_1,v_2)A^2:v_1v_2\}`$ projects to its first factor with fibers of Euler measure $`\chi (A)1`$, so by the Fubini principle its Euler measure must be $`\chi (A)(\chi (A)1)`$; the same approach works when applied iteratively in the case of larger values of $`k`$. This proof has the virtue of being entirely parallel to the standard way of counting ordered $`k`$-element subsets of a finite set, and it makes it clear where the factors $`\chi (A)i`$ come from. However, as an intimation of what is to come, it is more helpful to prove the result by poset-theoretic inclusion-exclusion, as we now do. (See Stanley’s book \[Sta\] for background on partially ordered sets.)

The partition lattice $`\mathrm{\Pi }_k`$ is the set of unordered partitions of a $`k`$-element set, ordered so that its top element is the 1-block partition $`\widehat{1}`$ and its bottom element is the $`k`$-block partition $`\widehat{0}`$. It is easy to see that for any finite set $`S`$, the product set $`S^n=\{(s_1,s_2,\mathrm{},s_n)\}`$ splits up into $`n`$ disjoint subsets, where the $`k`$th subset (with $`k`$ between 1 and $`n`$ inclusive) consists of those $`n`$-tuples that contain exactly $`k`$ distinct elements of $`S`$ and therefore has cardinality $`|S|^k`$. Ordinary inclusion-exclusion applied to the lattice $`\mathrm{\Pi }_k`$ lets us turn this around and conclude that the number of ordered $`k`$-element subsets of $`S`$ is

$$\underset{\pi \mathrm{\Pi }_k}{}\mu (\widehat{0},\pi )|S|^{N(\pi )},$$

where $`\mu (,)`$ is the Möbius function of $`\mathrm{\Pi }_k`$ and $`N(\pi )`$ is the number of blocks in the partition $`\pi `$. (See exercise 44 of \[Sta\].) Proceeding analogously in the case of polyhedral sets, we see that for any polyhedral set $`A`$, the product polyhedron $`A^n=\{(v_1,v_2,\mathrm{},v_n)\}`$ splits up into $`n`$ disjoint polyhedral sets, where the $`k`$th polyhedral set (with $`k`$ between 1 and $`n`$ inclusive) consists of those $`n`$-tuples that contain exactly $`k`$ distinct elements of $`A`$ and therefore has Euler measure $`\chi (A)^k`$. As before, inclusion-exclusion applied to the lattice $`\mathrm{\Pi }_k`$ tells us that the Euler measure of $`\left(\genfrac{}{}{0pt}{}{A}{k}\right)`$ is

$$\underset{\pi \mathrm{\Pi }_k}{}\mu (\widehat{0},\pi )\chi (A)^{N(\pi )}.$$

Since the two polynomials are the same (aside from the fact that one has $`\chi (A)`$ where the other has $`|S|`$), and since we know that the former is equal to $`|S|(|S|1)\mathrm{}(|S|k+1)`$, the desired result follows.

(The operations $`\lambda ^k:A\left(\genfrac{}{}{0pt}{}{A}{k}\right)`$ are akin to the $`\lambda `$-ring structures that Morelli defines on the abelian group of hedral functions, in Theorem 13 of his article and the Remark that follows it.)

As an aside, it is interesting to note that although finding a polyhedral analogue of $`n2^n`$ seem to require us to go beyond the polyhedral category and into an extension involving infinite-dimensional objects, finding a polyhedral analogue of the Fibonacci numbers does not. Given a polyhedral set $`P𝐑`$ of Euler measure $`n`$ (positive, negative, or zero), the set of all finite subsets $`SP`$ such that $`\chi ((PS)(t,t^{}))`$ is even for all $`t,t^{}S\{\mathrm{},+\mathrm{}\}`$ has Euler measure $`(\varphi ^{n+1}(1/\varphi )^{n+1})/\sqrt{5}`$ (with $`\varphi =(1+\sqrt{5})/2`$).

3. The Small Power Set

To motivate our definition of the Euler series, we begin by considering the “small power set” discussed earlier, that is, the set of all finite subsets of $`A`$, where $`A`$ is a polyhedral set. We denote the small power set of $`A`$ by $`2^A`$. $`2^A`$ is not itself a polyhedral set, but it can be realized as a formal union of disjoint polyhedral sets in spaces of increasing dimensionality. For instance, in the case where $`A`$ is an open (not necessarily bounded) interval, we can view $`2^A`$ as consisting of a point (the empty subset of $`A`$), an open interval (the set of singleton subsets of $`A`$), an open triangle (the set of 2-element subsets of $`A`$), an open 3-simplex (the set of 3-element subsets of $`A`$), and so on. These constituent pieces have Euler measure $`1,1,1,1,\mathrm{}`$ respectively. The terms do not yield a summable sequence, but there is no doubt what the incorrigible formalist Euler would have done: the series $`11+11+\mathrm{}`$ is after all a geometric series, so it can be assigned a notional sum $`1/(1(1))=1/2`$. More precisely, and more generally, given a polyhedral set $`A`$ one can define a formal power series $`E(t)`$ in which the coefficient of $`t^k`$ is the Euler measure of the set of $`k`$-element subsets of $`A`$. That is, $`E(t)=\left(\genfrac{}{}{0pt}{}{m}{0}\right)+\left(\genfrac{}{}{0pt}{}{m}{1}\right)t+\left(\genfrac{}{}{0pt}{}{m}{2}\right)t^2+\mathrm{}`$ where $`m=\chi (A)`$. This is just the binomial expansion of $`(1+t)^m`$; when the integer $`m`$ is non-negative, the series is a polynomial whose value at $`t=1`$ is $`2^m`$, while for $`m<0`$ the series converges in an open disk of radius 1 to a function that extends to the whole complex plane away from $`t=1`$, and whose value at $`t=1`$ is $`2^m`$. Either way, the value that one gets for $`\chi (2^A)`$ is $`2^{\chi (A)}`$. This “functoriality” is a sign that we are on the right track. (Here and elsewhere in the article, functoriality is said to hold when operations in the polyhedral category carry over, under the mapping $`\chi `$, to operations in the category of natural numbers, as described by Schanuel \[Sch\].)

Some less trivial corroboration of the validity of our approach comes from looking at $`\left(\genfrac{}{}{0pt}{}{B}{2}\right)`$, where $`B`$ is the power set of an open interval $`A`$. That is, we are looking at all the ways of choosing two distinct finite subsets $`S_1`$, $`S_2`$ of the interval $`A`$. Given $`\sigma =\{S_1,S_2\}\left(\genfrac{}{}{0pt}{}{B}{2}\right)`$, so that $`S_1S_2`$, define the support-set (or support) of $`\sigma `$ as $`S_1S_2`$, and define the rank of $`\sigma `$ as $`|S_1S_2|`$. We define an Euler series $`E(t)`$ in which the coefficient of $`t^k`$ is the Euler measure of the rank-$`k`$ part of $`\left(\genfrac{}{}{0pt}{}{B}{2}\right)`$. Given a $`k`$-element subset $`S`$ of $`A`$, there are $`3^{k1}+3^{k2}+\mathrm{}+3+1`$ ways to define $`\sigma \left(\genfrac{}{}{0pt}{}{B}{2}\right)`$ with support $`S`$. (Namely: for each $`i`$ between $`0`$ and $`k1`$ inclusive, one can assign the first $`i`$ elements of $`S`$ to both $`S_1`$ and $`S_2`$, the next element of $`S`$ to $`S_1`$ but not $`S_2`$, and each of the next $`ki1`$ elements of $`S`$ to $`S_1`$ only, $`S_2`$ only, or both $`S_1`$ and $`S_2`$, as one wishes. All $`3^{k1}+3^{k2}+\mathrm{}+3+1`$ such pairs $`S_1,S_2`$ are distinct from one another. On the other hand, every pair of distinct finite subsets of $`B`$ whose union has cardinality $`n`$ can be written in this form in a unique way, by labelling the two sets as $`S_1`$ and $`S_2`$ so that the smallest element in their symmetric difference belongs to $`S_1`$.) Thus the Euler measure of the rank-$`k`$ part of $`2^A`$ is $`\chi (\left(\genfrac{}{}{0pt}{}{A}{k}\right))(3^{k1}+3^{k2}+\mathrm{}+3+1)`$. The resulting Euler series $`t+(3+1)t^2(9+3+1)t^3+\mathrm{}`$ analytically extends to the rational function $`t/(1+t)(1+3t)`$, whose value at $`t=1`$ is $`1/8=\left(\genfrac{}{}{0pt}{}{1/2}{2}\right)`$.

Let us agree to limit ourselves to taking power sets of subsets of $`𝐑^2=𝐑`$. (Note that the Euler measure of a polyhedral subset of $`𝐑`$ can be any integer, so not much generality has been lost through this restriction, in terms of the sort of numbers that will arise as Euler measures of derived sets.) We are going to consider sets of the form

$$\left(\genfrac{}{}{0pt}{}{B}{k_1}\right),$$

$$\left(\genfrac{}{}{0pt}{}{\left(\genfrac{}{}{0pt}{}{B}{k_1}\right)}{k_2}\right),$$

$$\left(\genfrac{}{}{0pt}{}{\left(\genfrac{}{}{0pt}{}{\left(\genfrac{}{}{0pt}{}{B}{k_1}\right)}{k_2}\right)}{k_3}\right),$$

etc., where $`A`$ is a polyhedral subset of $`𝐑`$ and $`B=2^A`$; we will define them iteratively. We will call these sets gizmos (a more respectable term might be “selection schemes”). At each step, we will introduce a linear ordering on the gizmo, along with a support-set operation; together these allow us to proceed to the next step. Our notation for such gizmos will be $`G(2^A;k_1)`$, $`G(2^A;k_1,k_2)`$, $`G(2^A;k_1,k_2,k_3)`$, etc.

We start with $`2^A`$ itself. Given subsets $`S,T`$ of $`A`$ with $`|S|=|T|`$, we say that $`S`$ is less than $`T`$ if the smallest element in their symmetric difference belongs to $`S`$. (Equivalently: Writing $`S=\{s_1,s_2,\mathrm{},s_k\}`$ and $`T=\{t_1,t_2,\mathrm{},t_k\}`$ with $`s_1<s_2<\mathrm{}<s_k`$ and $`t_1<t_2<\mathrm{}<t_k`$, $`S`$ is less than $`T`$ if $`s_i<t_i`$ for the first location $`i`$ at which the two sequences differ.) The support-set operation from $`2^A`$ to itself is just the identity operation. Given a gizmo $`G`$ with a linear ordering and a support-set operation, and given an integer $`k0`$, we define the gizmo $`\left(\genfrac{}{}{0pt}{}{G}{k}\right)`$ as the set of all sequences $`(g_1,g_2,\mathrm{},g_k)G^k`$ with $`g_1<g_2<\mathrm{}<g_k`$ (with respect to the linear ordering on $`G`$). We order $`\left(\genfrac{}{}{0pt}{}{G}{k}\right)`$ lexicographically (so that $`(g_1,g_2,\mathrm{},g_k)`$ is less than $`(g_1^{},g_2^{},\mathrm{},g_k^{})`$ if and only if $`g_i<g_i^{}`$ at the first location $`i`$ at which the two sequences differ) and we define the support of $`(g_1,\mathrm{},g_k)`$ as the union of the supports of $`g_1`$ through $`g_k`$.

For any gizmo $`G`$, the number of elements of $`G`$ whose support is a particular $`k`$-element subset of $`A`$ does not depend on the particular subset of $`A`$, but is some number $`n_k`$ that depends only on $`k`$ (and $`G`$). We therefore define the Euler series of $`G`$ as

$$\underset{k=0}{\overset{\mathrm{}}{}}\chi (\left({\scriptscriptstyle \genfrac{}{}{0pt}{}{A}{k}}\right))n_kt^k.$$

4. Continued Binomial Coefficients

Theorem 1: The Euler series for $`G(2^A;k_1,k_2,\mathrm{},k_r)`$ converges in a neighborhood of $`t=0`$ to a rational function of $`t`$ whose value at $`t=1`$ is the iterated binomial coefficient $`C(2^{\chi (A)};k_1,k_2,\mathrm{},k_r)`$ (defined by the initial condition $`C(n)=n`$ and the recurrence relation $`C(n;k_1,k_2,\mathrm{},k_i)=\left(\genfrac{}{}{0pt}{}{C(n;k_1,k_2,\mathrm{},k_{i1})}{k_i}\right)`$).

Remark: While this result is not as general as one might like, it does back up my contention that my definition of the power set, or something like it, can interact in a functorially appropriate way with other more basic polyhedral operations, once the category of polyhedral sets and polyhedral maps has been extended in a suitable way.

Proof: As in the argument given at the end of Section 2, the proof of the main theorem uses inclusion-exclusion. We illustrate the proof for the concrete case $`G(2^{\chi (A)};2,2)`$ where $`A`$ varies over the polyhedral subsets of $`𝐑`$, since it contains all the key ideas. An element of this gizmo can be described as a 4-tuple $`(S_1,S_2,T_1,T_2)`$ of subsets of $`A`$ satisfying the following three conditions: (1) $`S_1<S_2`$; (2) $`T_1<T_2`$; (3) $`\{S_1,S_2\}<\{T_1,T_2\}`$. Note that condition (3) is equivalent to the condition that either (i) $`S_1<T_1`$ or else (ii) $`S_1=T_1`$ and $`S_2<T_2`$. We know that, in the case where $`A`$ is a finite set of points, the gizmo is itself just a finite set of points, and that its Euler series is a polynomial in $`t`$ whose value at $`t=1`$ is just $`C(2^{\chi (A)};2,2)`$. Our goal is to show that this same formula applies even when $`A`$ is not a finite set of points, essentially by showing that the two cases are governed by the same inclusion-exclusion calculation.

Before we do the comparison, it is helpful to switch to a re-statement of the problem more amenable to inclusion-exclusion analysis. We will the consider the “symmetrized gizmo” consisting of all 4-tuples obtained from the elements of $`G(2^{\chi (a)};2,2)`$ by the action of the 8-element group that allows one to switch $`S_1`$ with $`S_2`$, $`T_1`$ with $`T_2`$, or the pair $`S_1,S_2`$ with the pair $`T_1,T_2`$. Note that the elements of the symmetrized gizmo are exactly the 4-tuples $`(S_1,S_2,T_1,T_2)`$ of subsets of $`A`$ for which none of the following conditions holds: (a) $`S_1=S_2`$; (b) $`T_1=T_2`$; (c) $`S_1=T_1`$ and $`S_2=T_2`$; (d) $`S_1=T_2`$ and $`S_2=T_1`$. Regardless of the nature of $`A`$ (finite vs. infinite), it is clear that the Euler series of the new gizmo is 8 times the Euler series of the old, so it is enough to prove that the Euler series of the symmetrized gizmo is a rational function of $`t`$ whose value at $`t=1`$ is given by a polynomial in $`2^{\chi (A)}`$. Once we have shown that the value at $`t=1`$ depends only on $`\chi (A)`$, and not on $`A`$, our observations for the case where $`A`$ is finite will complete the proof.

To prove the required properties of the Euler series for symmetrized gizmos, the principle of inclusion-exclusion guarantees that it is enough to show that the claim holds when we look at the set of 4-tuples in which some specified (possibly empty) subset of the forbidden conditions (a),(b),(c),(d) does in fact hold. But the effect of imposing any or all of these four conditions is merely to require that one is not choosing four subsets of $`A`$ but three, two, or one. Specifically, if we let $`E_a`$ denote the Euler series for the set of 4-tuples of (finite) subsets of $`A`$ for which condition (a) holds, etc., and we let $`E(k)`$ denote the Euler series for the set of $`k`$-tuples of subsets of $`S`$, then the Euler series for $`G(2^A;2,2)`$ is equal to $`[E][E_a+E_b+E_c+E_d]+[E_{ab}+E_{ac}+E_{ad}+E_{bc}+E_{bd}+E_{cd}][E_{abc}+E_{abd}+E_{acd}+E_{bcd}]+[E_{abcd}]=[E_4][E_3+E_3+E_2+E_2]+[E_2+E_1+E_1+E_1+E_1+E_1][E_1+E_1+E_1+E_1]+[E_1]=E_42E_3E_2+2E_1`$. (Note that if we replace $`E_k`$ by $`x^k`$, the preceding expression reduces to

$$8\left(\genfrac{}{}{0pt}{}{\left(\genfrac{}{}{0pt}{}{x}{2}\right)}{2}\right),$$

as it ought to.) Here we see that the Möbius coefficients that enter into the inclusion-exclusion formula are independent of the set $`A`$, but depend only on $`k_1=2`$, $`k_2=2`$. This independence is true more generally. Hence the entire problem reduces to the case of $`(2^A)^j`$, with $`j`$ varying; that is, the gizmo consisting of all ways of choosing a $`j`$-tuple of (not necessarily distinct) subsets $`S_1,S_2,\mathrm{},S_j`$ of $`A`$.

We can now conclude the proof. When $`A`$ is a finite set of points, the Euler series for $`(2^A)^j`$ is $`(1+(2^j1)t)^{\chi (A)}`$, since for each of the $`\chi (A)`$ points of $`A`$ we can either choose to include it in none of the sets $`S_1,S_2,\mathrm{},S_j`$ or else choose to include it in some of the $`j`$ sets (which we can do in any of $`2^j1`$ ways). When $`A`$ is an open interval, the Euler series for $`(2^A)^j`$ is $`(1+(2^j1)t)^1=1(2^j1)t+(2^j1)^2t^2\mathrm{}`$, since for each subset of $`A`$ of size $`i`$, the number of ways of assigning each element of the subset to at least one (and possibly all) of the $`j`$ sets $`S_1,S_2,\mathrm{},S_j`$ is $`(2^j1)^i`$. In general, the Euler series for $`(2^A)^j`$ is $`(1+(2^j1)t)^{\chi (A)}`$, which when evaluated at $`t=1`$ gives $`(2^j)^{\chi (A)}=(2^{\chi (A)})^j`$. Note that this is a polynomial in $`2^{\chi (A)}`$, as claimed. Hence all the gizmos in question (the original gizmo $`G(2^A;k_1,k_2,\mathrm{},k_r)`$ and its symmetrized version) have Euler measure given by polynomials in $`\chi (A)`$, which we can identify by looking at the case where $`A`$ is a finite set. $`\mathrm{}`$

Theorem 1 in its current form is a bit eccentric in its focus. Ideally one would like to replace it by a result that says that the operation selecting a subset of size $`k`$ may be freely combined with the operations of disjoint union and Cartesian product (and perhaps other natural polyhedral operations as well) in such a fashion that all the resulting polyhedral sets satisfy formulas like $`\left(\genfrac{}{}{0pt}{}{\chi (A}{k)=\left(\genfrac{}{}{0pt}{}{\chi (A)}{k}\right)}\right)`$, $`\chi (AB)=\chi (A)+\chi (B)\chi (AB)`$, and $`\chi (A\times B)=\chi (A)\chi (B)`$.

5. Exponentiation

If $`A,B`$ are polyhedral sets ($`A𝐑^m`$, $`B𝐑^n`$), the set of affine maps from $`A`$ to $`B`$ has a well-defined Euler measure that can be associated with it in a number of ways. One simple way of defining it is to let $`v_1,v_2,\mathrm{},v_k`$ be points in $`𝐑^m`$ whose affine span contains $`A`$; then every map from $`A`$ to $`B`$ is determined by its values on $`v_1,v_2,\mathrm{},v_k`$, so that the set of affine maps from $`A`$ to $`B`$ becomes associated with a subset of $`B^k`$ (namely, the set of all $`k`$-tuples $`(f(v_1),f(v_2),\mathrm{},f(v_k))`$ arising from affine maps $`f:AB`$). One can show that this set in $`B^k`$ is itself polyhedral, and that its Euler measure does not depend on the spanning set $`v_1,v_2,\mathrm{},v_k`$ that was chosen.

For example, in the case $`A=(0,1)𝐑`$ it is natural to take $`v_1=0`$ and $`v_2=1`$. Assume that $`B`$ is compact. Then the set of affine maps from $`A`$ to $`B`$ corresponds to the set of pairs $`(w,w^{})𝐑^n\times 𝐑^n`$ such that $`tw+(1t)w^{}B`$ for all $`0<t<1`$. It is not hard to show that this set is polyhedral and that for each $`w𝐑^n`$, the set $`\{w^{}𝐑^n:tw+(1t)w^{}B\text{ for all }0<t<1\}`$ has Euler measure 1 or 0 according to whether or not $`wB`$ (this is where compactness gets used, via the fact that a starlike compact polyhedral set has Euler measure 1). Hence, by Fubini’s Theorem, the set of affine maps from $`(0,1)`$ to $`B`$ has Euler measure $`\chi (B)`$.

Let us now take $`A`$ to be a polyhedral subset of $`𝐑`$ with no isolated points. (It is clear that the case of general polyhedral sets $`A`$ can be handled by combining the case in which $`A`$ is finite and the case in which $`A`$ has no isolated points.) Then any polyhedral map $`f:AB`$ is locally affine at all but a finite number of points. We call a point at which $`f`$ fails to be locally affine (whether by virtue of discontinuity or by virtue of non-differentiability) a breakpoint of $`f`$, and we call the number of breakpoints the rank of $`f`$. (Note that this notion of rank, specialized to the case $`B=\{0,1\}`$, coincides with the notion introduced in section 3: since we are assuming that $`A`$ has no isolated points, the indicator function of a finite subset of $`A`$ with cardinality $`n`$ has exactly $`n`$ points of discontinuity, namely, the elements of $`A`$ themselves.)

Let $`S`$ be a finite subset of $`A`$. We wish to count (in the sense of Euler-measure) the polyhedral maps from $`A`$ to $`B`$ whose breakpoints are precisely the points of $`S`$. It is not hard to see that the set of maps from $`A`$ to $`B`$ that are locally affine on $`AS`$ is amenable to the same analysis as we applied in the first paragraph of this section; that is, it has a well-defined Euler measure. If a set is locally affine on $`AS`$, its breakpoints form a subset of $`S`$, and vice versa. Therefore the method of inclusion-exclusion gives us a way to calculate the Euler measure for the set of polyhedral maps from $`A`$ to $`B`$ with breakpoint-set $`S`$. It is easy to prove that this Euler measure has the same value for all $`S`$ of fixed size $`k`$. Call this value $`n_k`$. Then one can define the Euler measure of the set of polyhedral maps with exactly $`k`$ breakpoints as $`\chi (\left(\genfrac{}{}{0pt}{}{A}{k}\right))n_k=\left(\genfrac{}{}{0pt}{}{\chi (A)}{k}\right)n_k`$, and define $`E(t)=_{k=0}^{\mathrm{}}\left(\genfrac{}{}{0pt}{}{\chi (A)}{k}\right)n_kt^k`$.

For example, suppose $`A=(0,1)`$, $`B=\{0,1\}`$. A polyhedral map from $`A`$ to $`B`$ will have some number of breakpoints, say at $`0<x_1<x_2<\mathrm{}<x_k`$ ($`k0`$); in between these breakpoints, the map will be constant, since $`B`$ is discrete. I claim that the number of polyhedral maps from $`A`$ to $`B`$ with some fixed set of $`k`$ breakpoints is $`23^k`$. For, starting at the left, the map can take either the value 0 or the value 1 on the interval $`(0,x_1)`$. Thereafter, for each $`1ik`$, we can choose either value (0 or 1) for the function at $`x_i`$ and either value for the function on $`(x_i,x_{i+1})`$, as long as we do not have both values equal to the value on $`(x_{i1},x_i)`$ (for then there would not be a breakpoint at $`x_i`$); there are always 3 ways to implement this. Hence the Euler series is $`26t+18t^254t^3+\mathrm{}=2/(1+3t)`$, whose value at $`t=1`$ is $`2/4=2^1=\chi (B)^{\chi (A)}`$.

Continuing this example in the spirit of Theorem 1, we can look at the set of pairs of distinct maps from $`(0,1)`$ to $`\{0,1\}`$, where the rank of such a pair is the cardinality of the union of the sets of breakpoints of the two maps. Then one can check that the set of maps with $`k`$ specified breakpoints has Euler measure $`2(2^41)^k(2^21)^k=215^k3^k`$, so that $`E(t)=\frac{2}{1+15t}\frac{1}{1+3t}`$ and $`E(1)=\frac{2}{16}\frac{1}{4}=\frac{1}{8}=\left(\genfrac{}{}{0pt}{}{1/2}{2}\right)`$ as expected.

Let us now focus on the case where $`A=(0,1)`$ and $`B`$ is a bounded polyhedral set. There are two sorts of exponentiation we will consider: the “Morelli exponential” (the set of all hedral maps from $`A`$ to $`B`$) and the “Schanuel exponential” (the set of all polyhedral maps from $`A`$ to $`B`$); for definitions, see section 1. The two cases are similar, for in both cases, the set of affine maps from any sub-interval $`(x_i,x_{i+1})`$ to $`B`$ has Euler measure $`\chi (B)`$. (In the Morelli case, this is trivial, since one is merely looking at constant maps; in the Schanuel case, one appeals to the analysis given in the second paragraph of this section.) This allows us to prove:

Theorem 2: Let $`B`$ be any bounded polyhedral set, and let $``$ be the set of hedral maps from $`(0,1)`$ to $`B`$ or the set of polyhedral maps from $`(0,1)`$ to $`B`$. In either case, the regularized Euler measure of $``$ is 0 if $`\chi (B)=0`$ and $`1/\chi (B)=\chi (B)^1`$ otherwise.

Proof: For every set $`S`$ of $`k`$ points in $`(0,1)`$, the set of hedral or polyhedral maps from $`(0,1)`$ to $`B`$ whose breakpoints lie in $`S`$ has Euler measure $`\chi (B)^{2k+1}`$, with $`k`$ factors of $`\chi (B)`$ coming from the points in $`S`$ and $`k+1`$ factors coming from the intervals in between (by virtue of the fact mentioned in the paragraph preceding the statement of Theorem 2). Inclusion-exclusion on the lattice of subsets of $`S`$ implies that the set of maps from $`(0,1)`$ to $`B`$ whose breakpoint-set is $`S`$ has Euler measure $`\chi (B)(\chi (B)^21)^k`$. Hence $`E(t)=\chi (B)/(1+(\chi (B)^21)t)`$. If $`\chi (B)=0`$, all coefficients vanish; if $`\chi (B)0`$, $`E(1)=\chi (B)/(\chi (B)^2)=1/\chi (B)`$. $`\mathrm{}`$

6. Final thoughts

It is worth pointing out here that full functoriality is too much to expect; for instance, if $`A`$ has generalized Euler measure $`1/2`$ and $`B`$ has Euler measure $`1`$, then it is hard to believe that one could come up with a theory that would assign to $`C=B^A`$ a generalized Euler measure $`(1)^{\frac{1}{2}}=\pm i`$, and still harder to believe such a theory would know what to make of $`C^C`$ (given that $`i^i`$ takes on infinitely many values). Yet numerous examples demonstrate that as long as one does not push things too far, functoriality does hold.

I do not claim that my proposed definition of regularized Euler measure is the “right” one; if the reader finds its formulation artificial, I must to some extent agree. However, I think that the way the theory-in-progress hangs together indicates that there is probably a better-motivated definition that agrees with mine for the class of examples I consider. A homological formulation would be especially useful, and I hope some readers will pursue this line of investigation.

References

\[Ch1\] B. Chen, The Gram-Sommerville and Gauss-Bonnet theorems and combinatorial geometric measures for noncompact polyhedra, Advances in Mathematics 91 (1992), 269-291.

\[Ch2\] B. Chen, On the Euler measure of finite unions of convex sets, Discrete and Computational Geometry 10 (1993), 79-93.

\[Ha1\] H. Hadwiger, Eulers Characteristik und kombinatorische Geometrie, J. reine angew. Math. 194 (1955), 101–110.

\[Ha2\] H. Hadwiger, Erweiterter Polyedersatz und Euler-Shephardsche Additionstheoreme, Abhandlungen Mat. Seminar Hamburg 39 (1973), 120–129.

\[KR\] D. Klain and G.-C. Rota, Introduction to geometric probability. Cambridge University Press, 1997.

\[Len\] H. Lenz, Mengenalgebra und Eulersche Charakteristik, Abh. Math. Sem. Univ. Hamburg 34 (1970), 135–147.

\[McM\] P. McMullen, The polytope algebra, Advances in Mathematics 78 (1987), 76–130.

\[McS\] P. McMullen and R. Schneider, Valuations on convex bodies, Convexity and its Applications, P. Gruder and J. Wills, eds., Birkhäuser, Basel 1983, pp. 160–247.

\[Mor\] R. Morelli, A theory of polyhedra, Advances in Mathematics 97 (1993), 1–73.

\[Pro\] J. Propp, Euler measure as generalized cardinality; arXiv: math.CO/ 0203289.

\[Sch\] S. Schanuel, Negative sets have Euler characteristic and dimension, in: Lecture Note In Mathematics #1488, Springer, 1991.

\[Sta\] R. Stanley, Enumerative combinatorics I, Wadsworth, 1986.


format: text

text:
# The effect of mixing on metallicity gradients in the ICM

## 1 Introduction

The X-ray surface brightness of many clusters of galaxies shows a strong central peak which is generally interpreted as the signature of a cooling flow (Cowie & Binney 1977, Fabian & Nulsen 1977, Sarazin 1988, Fabian 1994). However, the simple cooling flow model conflicts with a growing number of observations that show that while the temperature is declining in the central region, gas with a temperature below $``$1 keV is significantly less abundant than predicted. Moreover, the level of star formation found in present-day cooling flows is insufficient to explain the mass deposition rates inferred from X-ray observations (Böhringer et al. 2001).

Many clusters of galaxies host a radio source at their centres. Around 71% of all cD galaxies in cooling flows show evidence for radio activity which implies that either the sources are active for longer than 1 Gyr or that they recur with high frequency (Burns 1990). These radio sources emit large quantities of hot relativistic plasma into the intracluster medium (ICM). The displacement of thermal plasma by the relativistic plasma can create ‘holes’ in the X-ray surface brightness such as observed in the Perseus cluster (Böhringer et al. 1993). Holes in the X-ray emission are also seen in the radio lobes of the Hydra A cluster (McNamara et al. 2000), Centaurus A (Saxton, Sutherland and Bicknell 2001) and other clusters (Fabian et al. 2000, Mazzotta et al. 2001). It has also been shown that bubbles are difficult to detect and that a substantial amount of energy may be hidden in this form in clusters (Brüggen et al. 2002).

The radio jets inflate pockets of low-density radio plasma that rise buoyantly through the cluster atmosphere. In the course of this rise the cavities become Rayleigh-Taylor unstable and form plumes and mushroom-like structures (Churazov et al. 2001, Brüggen & Kaiser 2001). It has been shown that the rising bubbles can lead to mixing and a significant restructuring of the inner regions of a cooling flow (Churazov et al. 2001, Brüggen & Kaiser 2001, Nulsen et al. 2001). This process can reduce the deposition of cool gas without directly introducing heat and thus can form an effective machanism to couple the mechanical power of the radio jets to the cooling gas. In some cases this may be sufficient to quench the cooling flow (Binney 2001, Binney & Tabor 1995).

However, pronounced metallicity gradients have been detected in the core of some clusters. Chandra observations of Hydra A have shown a factor 2 increase in the Fe and Si abundances within the central 100 kpc (David et al. 2001). Fukazawa et al. (2000) measured metal abundances in the ICM of 34 clusters using ASCA data and also found increments in Si and Fe abundances around cD galaxies. Similar results have been found by Ezawa et al. (1997), Ikebe et al. (1997), Finoguenov & Ponman (1999) and Finoguenov et al. (2000). Fukazawa et al. (2000) found that the central Fe abundances of cD type clusters are by up to a factor of 3 higher than their spatially averaged values, while non-cD type clusters do not show this Fe excess at the centre. This indicates that the excess metals are probably ejected from the central cD galaxy but the exact origin of these metals is still subject to debate. However, the metals have to be produced at the cluster centre because gravitational settling of the heavy ions takes much longer than the age of the Universe (Sarazin 1988). (The abundance gradients occur mostly around cD galaxies which are believed to sit still at the bottom of the cluster potential. Therefore, it seems unlikely that the metals have been stripped off by the ICM.)

The presence of steep metallicity gradients constrains the amount of mixing that can take place in a cluster. If stirring motions induced by buoyant bubbles erase the metallicity gradient faster than the metals can be replenished at the centre, it will be difficult to explain the observed metallicity gradients. In order to calculate the effect of rising bubbles on the metallicity gradients of a cluster, I performed three-dimensional hydrodynamical simulations and employed tracer particles to compute the metallicity gradients.

## 2 Method

The simulations were obtained using the parallel version of the ZEUS-3D code which was developed especially for problems in astrophysical hydrodynamics. The code uses finite differencing on an Eulerian grid and is fully explicit in time. For a detailed description of the algorithms and their numerical implementation see Stone & Norman (1992a, b). In order to study the mixing of the ICM, the ZEUS code was modified to follow the motion of tracer particles. These tracer particles are passively advected with the fluid and their positions are written to a file at regular intervals.

In our simulations we employed an ideal gas equation of state and we ignored the effects of magnetic fields and rotation. The cooling time in our cluster model is of the order of 500 Myrs which is sufficiently long compared to the time scales considered here so that we ignored it. Moreover, numerical experiments have shown that the neglect of cooling does not influence the large-scale dynamics of the fluid over the time-scales considered here. The simulations were computed on a Cartesian grid in three dimensions. The serial version of the code was run on a SUN ULTRA workstation and the parallel version on an SGI ORIGIN 3000. The computational domain spanned 30 kpc in height with a base area of 10 kpc $`\times `$ 10 kpc. This domain was covered by 450 $`\times `$ 150 $`\times `$ 150 equally spaced grid points. On the boundaries outflow conditions were chosen.

The initial mass and temperature distributions were modelled on the Virgo Cluster as given by Nulsen & Böhringer (1995). The gas density distribution was found by assuming hydrostatic equilibrium to maintain an initially static model. Here the results from two runs are reported: run 1 and run 2. At each timestep hot gas was injected into a spherical region which was situated at a distance of $`d=9`$ kpc from the gravitational centre. This spherical region had a radius of $`r_\mathrm{b}=0.7`$ kpc and $`r_\mathrm{b}=1.`$ kpc for run 1 and 2, respectively. The injected gas had zero initial velocity and the energy injection rates were $`L=4.4\times 10^{41}`$ erg s-1 (run 1) and $`L=1.3\times 10^{43}`$ erg s-1 (run 2). This particular set of parameters was chosen because they could reproduce a lot of observational features of late-stage radio activity. As shown elsewhere (Brüggen et al. 2002) the latter luminosity is close to the maximum rate of thermal energy injection that produces a buoyant plume. Even larger energy injection rates have the effect of evacuating the core before the bubble can start to rise. It was assumed that the material in the injection region is supplied by the jet of an Active Galactic Nucleus (AGN). The flow from the centre of the source to the injection region is taken to be ballistic while the injection region itself may be identified with the location of a strong internal shock that brings the jet into pressure equilibrium with the surrounding gas. This picture motivates the choice of parameters for the position of the injection region for which $`d=9`$ kpc is a reasonable choice. We have also explored slightly different parameters in Brüggen et al. (2002) and found that the resulting flow is not strongly dependent on the choice of these parameters.

At the start of the simulation 24,000 tracer particles were distributed uniformly throughout the computational domain (Fig. 1). The gas was treated as a single fluid and assumed to obey a polytropic equation of state with $`\gamma =5/3`$. For more details of the simulation the reader is referred to the paper by Brüggen et al. (2002).

Finally, I should address some issues related to the accuracy of these kinds of finite-difference hydrodynamical simulations. While the code can simulate large-scale mixing due to Rayleigh-Taylor and Kelvin-Helmholtz instabilities, it does not include real diffusion of particles. Any observed diffusion is therefore entirely numerical. Moreover, numerical viscosity is responsible for suppressing small-scale instabilities at the interface between the bubble and the cooler, surrounding, X-ray emitting gas. To assess the effects of numerical viscosity, I have repeated some simulations on grids with lower resolution. From these experiments I can conclude that ’global parameters’ such as the position and size of the ’mushrooms’ are relatively insensitive to the resolution. However, the detailed morphology on small scales does depend on the resolution and the choice of initial conditions.

## 3 Results

Fig. 1 shows the positions of the tracer particles at the start of the simulation and after a simulated time of 63 Myrs (run 2). The tracer particles are color coded according to their initial positions as shown in the left panel. One can see how the bubble displaces the material and lifts the tracer particles up from the cluster centre. In Fig. 2 the differences between the original iron abundance profiles and that at a time $`t`$ after the start of the energy injection are shown. I assumed an initial iron abundance given by $`A_{\mathrm{Fe}}=0.7030.0053(r/\mathrm{kpc})`$, where $`A`$ is measured in units of the solar abundance and $`r`$ denotes the distance from the centre of the cluster measured in kpc. This simple fit is based on the observations of Hydra A by David et al. (2000). The abundances have been averaged over horizontal slices through the computational volume. One should note that in computing the abundances I included only those metals of the ICM that already were present at the start of the simulation. I ignore those metals that may have been injected with the hot gas. Therefore, the changes in the abundances are caused solely by the displacement and mixing of the ICM. Thus, one obtains an upper limit on the abundance changes because the injected gas may contain metals, too. Furthermore, one can note that the bubbles might dredge up material from the central regions of the active galaxy which is likely to be enriched by stars of supersolar abundances.

Fig. 3 shows the same as Fig. 2 only that the abundances have been averaged over semispherical shells. This way the observations are simulated where the spectra are averaged over concentric rings about the cluster centre. Here I average over semi-spherical shells assuming that the jet is part of an ambipolar jet. Now, the abundance changes are much smaller at larger radii due to the larger volume to average over. The only noticeable decrease occurs within the inner 5 kpc. For run 1 the abundance changes found are smaller than the accuracy of the observed values which is of order of $`\pm 0.1`$ solar (David et al. 2001). For run 2, the metallicity at the inner edge of the core is decreased by nearly 0.25 solar but averaged over the inner 5 kpc the decrease is roughly 0.1 solar.

## 4 Conclusion

It was found that mixing by buoyant bubbles has a relatively weak impact on the metallicity gradients found in clusters. The simulations were stopped after 120 Myrs but during this time supernovae in the cD galaxy will have replenished some of these metals.

The fact that the Si/Fe abundance ratio of the central excess metals is close to the solar value suggests that typa Ia supernovae are the main suppliers for these metals. Over the course of our simulation the Fe abundance within the central 5 kpc was decreased on average by $`\mathrm{\Delta }A0.1`$ solar (run 2). Given that the mass within the same radius is $`M_{5\mathrm{k}\mathrm{p}\mathrm{c}}10^{10}M_{}`$ this corresponds to an Fe mass of $`\mathrm{\Delta }M_{\mathrm{Fe}}10^4M_{}`$. Adopting the supernovae yields that have been averaged over a Salpeter IMF as given by Finoguenov et al. (2000), $`y_{\mathrm{Fe}}=0.74M_{}`$ (type Ia), it requires less than $`N_{\mathrm{SN}}=\mathrm{\Delta }M_{\mathrm{Fe}}/y_{\mathrm{Fe}}2\times 10^4`$ supernovae to replenish the dredged-up mass of iron. Over the time of the simulation this would correspond to a supernova rate of $`2\times 10^4`$ yr-1 which is still less than the expected type Ia supernova rate for a typical central galaxy (van den Bergh 1987). Hence, the metals lost through the mixing by bubbles can be replenished by supernovae from the old stellar population in the central galaxy.

Therefore, one can conclude that the heating scenario in cooling flows is consistent with the measured metallicity gradients.

Fukazawa et al. (2000) found that the central iron excess is larger (by 0.2-0.4 solar) in cooler clusters than in hotter ones (by 0.1-0.2 solar). This observation could provide circumstantial evidence for the action of radio galaxies. In cooler clusters some time has passed since the last epoch of radio activity and heating has occurred. The hotter clusters, on the other hand, have recently experienced a period of radio activity or are in the middle of such a period. If the radio activity caused a redistribution of the metals, one would expect to find a correlation between the central metal excesses and the radio power of the cD galaxy. This conjecture could be checked observationally. Alternatively, the observations by Fukazawa et al. (2000) might merely reflect the fact that the hotter clusters may have been disturbed recently by a merger which has mixed the metals and reduced the central metallicity increment.

Some of the computations reported here were performed using the UK Astrophysical Fluids Facility (UKAFF). I thank Jim Pringle for helpful discussions and the referee for useful comments.


text:
# 1 Introduction

## 1 Introduction

In 1984 Beale, Kato and Maida  showed that sufficient and necessary condition for a smooth solution to the $`3D`$ Euler equations for ideal incompressible fluids on the time interval $`[0,t_0]`$ is the finiteness of the integral,

$$\underset{0}{\overset{t_0}{}}\text{sup}_r|\omega |𝑑t<\mathrm{}.$$
(1)
Here $`\omega `$ is the vorticity, connected with the velocity field $`𝐯`$ by the standard relation:

$$\omega =\text{curl}𝐯.$$

If the integral (1) is divergent, then the vorticity blows up (or collapses) in a finite time.

Proof of this criterion is based on the local existence theorem . According to this theorem a smooth solution to the 3D incompressible Euler equations exists if the initial conditions $`𝐯(0)`$ belong to the Sobolev space $`H^q(^3)`$ for $`q3`$ and the solution itself is from the class,

$$𝐯C([0,t_0];H^q)C^1([0,t_0];H^{q1}),$$

where the norm of $`H^q`$ is defined as follows,

$$𝐯_{H^q(^3)}=\left(\underset{\alpha q}{}\left(^\alpha 𝐯\right)^2𝑑𝐫\right)^{1/2},$$
(2)
with $`\alpha `$ being a multi-index. In accordance with , the violation of this property leads to divergence of the integral (1) and to collapse for vorticity.

However, in order to use the criterion (1), effectively one needs to have an explicit solution to the Euler equation, that is practically very difficult or impossible. Therefore, the main task is to find sufficient conditions which guarantee blow-up of the vorticity. The first step in this direction would be to construct such criteria for models more simple than the incompressible Euler equations which is one of the main aims of this paper.

We consider the simplest variant – the Euler equations for a compressible fluid without pressure (the hydrodynamics of dust) which are sometimes called the Hopf equations:

$$\frac{\rho }{t}+\text{div}\rho 𝐯=0,$$
(3)
$$\frac{𝐯}{t}+(𝐯)𝐯=0.$$
(4)
These equations can be successfully solved by means of the Lagrangian description. In terms of Lagrangian variables the equation (4) describes the motion of free fluid particles. For this reason breaking is possible in this model: it happens when trajectories of fluid particles intersect. In terms of the mapping, describing the change of variables from the Eulerian description to the Lagrangian one, the process of breaking corresponds to a vanishing Jacobian, $`J`$, of the mapping. This process can be considered as a collapse in this system: the density $`\rho `$ as well as the velocity gradient become infinite at a finite time in the points where the Jacobian vanishes. By virtue of this the breaking is sometimes called gradient catastrophe.

In this paper we construct sufficient integral criteria for breaking within the model (3),(4) and verify what the functional conditions of the theorem  correspond to those for the given system. Our main conjecture is the following: breaking in the 3D case within the system yields a violation of Sobolev space already for $`q=2`$ instead of $`q3`$ in the theorem .

To complete the introduction we would like to note that the model (3), (4) has a lot of astrophysical applications. According to the pioneering idea of Ya.B. Zeldovich , the formation of proto-galaxies due to breaking of stellar dust can be described by the system (3), (4) (see also the review ). It should be noticed also that recently  the mechanism of breaking was used to explain collapse in incompressible fluids when, instead of the breaking of fluid particles, the breaking of vortex lines can happen which results in a blow-up of vorticity. Therefore we believe that the results presented in this paper can also be useful from the point of view of collapse in incompressible fluids.

## 2 One-dimensional analysis

We begin with a one-dimensional calculations where Eq. (4) has the form

$$v_t+vv_x=0.$$
(5)
Consider the integrals

$$I_n=v_x^n𝑑x$$

with integer $`n`$, where integration is assumed to be from $`\mathrm{}`$ to $`+\mathrm{}`$.

It is easy to show that the evolution of these integrals is defined by the relations:

$$\frac{dI_n}{dt}=(n1)I_{n+1},$$
(6)
which can be obtained by means of (5) and integration by parts.

The relations for $`n=2,3`$ read:

$$\frac{d}{dt}v_x^2𝑑x=v_x^3𝑑x,$$
(7)
$$\frac{d}{dt}v_x^3𝑑x=2v_x^4𝑑x.$$
(8)
These differential relations allow one to find a closed differential inequality for the integral $`I_2=v_x^2𝑑x`$. Applying the Cauchy-Bunyakovsky inequality to (7) gives the estimation

$$\frac{dI_2}{dt}\left(v_x^4𝑑x\right)^{1/2}I_2^{1/2}.$$
(9)
Substituting,

$$I_4=v_x^4𝑑x=\frac{1}{2}\frac{d^2I_2}{dt^2},$$

into (9) we arrive at the closed differential inequality for $`I_2`$:

$$I_2\frac{d^2I_2}{dt^2}2\left(\frac{dI_2}{dt}\right)^20.$$
(10)
This inequality is solved by means of change of variables

$$I_2=X^\alpha >0,$$

with unknown exponent $`\alpha `$. We will define the exponent by requiring the absence in (10) of the terms $`X^{\alpha 2}X_t^2`$. Hence $`\alpha =1`$ and the inequality becomes,

$$X_{tt}0.$$
(11)
Using a mechanical interpretation of $`X`$ as the coordinate of a particle and taking into account that the particle acceleration is negative, we can immediately conclude that the particle can reach origin $`X=0`$ in a finite time if the initial particle velocity $`X_t(0)`$ is negative. For $`I_2`$ this means a blow-up:

$$I_2=\frac{1}{X}\mathrm{}.$$

Hence, by elementary calculation, we can estimate the collapse time $`t_0`$:

$$t_0<\frac{X(0)}{|X_t(0)|}\frac{I_2(0)}{|I_{2t}(0)|}.$$

After multiplying (11) by $`X_t`$ $`(<0)`$ and integrating the result over time, we arrive at the estimation from above for $`I_2`$ :

$$I_2(t)\frac{I_2^2(0)}{|I_{2t}(0)|(t_0t)}.$$
(12)
Thus, the blow-up of $`I_2`$ takes place if the initial velocity $`X_t(0)`$ is negative which is equivalent to the initial condition,

$$I_{2t}(0)>0$$

or,

$$I_3(0)=v_x^3𝑑x<0.$$

If the initial distribution of $`v_0(x)`$ is symmetric with respect to $`x`$, then $`I_3(0)0`$. However, as follows from (7), $`I_3`$ becomes negative already at $`t=+0`$ that results in a blow-up of the integral $`I_2`$.

To complete this section, let us compare the estimate (12) with the exact time dependence of $`I_2(t)`$ near the breaking point.

In order to define this dependence, we differentiate (5) with respect to $`x`$ that results in the equation,

$$\frac{dv_x}{dt}=v_x^2,\left(\frac{d}{dt}=\frac{}{t}+v\frac{}{x}\right).$$
(13)
After integration we have

$$v_x=\frac{v_0^{}(a)}{1+v_0^{}(a)t},$$
(14)
where $`a`$ is the initial coordinate of the fluid particle at $`t=0`$, $`v_0(a)`$ is the initial velocity, the prime denotes a derivative. The denominator in (14), $`1+v_a(a)t`$, represents the Jacobian $`J`$ of the mapping,

$$x=a+v(a)t.$$

The Jacobian tends to zero for the first time at $`t=t_0`$, defined by

$$t_0=\underset{a}{\mathrm{min}}[v_a^1(a)]>0.$$
(15)
Near the singular point, $`a=a_0`$, corresponding to the minimum, (15) $`J`$ can be approximated by the expression,

$$J\alpha (t_0t)+\gamma (aa_0)^2.$$

Here

$$\alpha =\frac{J}{t}|_{t_0,a_0}1/t_0,\mathrm{\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}2}\gamma =\frac{^2J}{a^2}|_{t_0,a_0}>0.$$

As a result, $`v_x`$ takes a singularity as $`\tau =t_0t0`$:

$$v_x=\frac{v_0^{}(a_0)}{\alpha \tau +\gamma (aa_0)^2}.$$
(16)
Hence, one can see that the contribution from the singularity to the integral $`I_2`$,

$$I_2\frac{v_0^2(a_0)}{\alpha \tau +\gamma (aa_0)^2}𝑑a\frac{1}{\tau ^{1/2}},$$
(17)
diverges as $`\tau 0`$ and satisfies the inequality (11).

## 3 Multi-dimensional breaking

In this section we generalize the above analysis to multi-dimensions.

In the multi-dimensional case, instead of $`v_x`$ it is convenient to introduce the matrix $`U`$ with matrix elements

$$U_{ij}=\frac{v_j}{x_i}.$$

The equations of motion for this matrix have a form analogous to (12):

$$\frac{dU}{dt}=U^2$$
(18)
where

$$\frac{d}{dt}=\frac{}{t}+(𝐯).$$

Our aim is now to find the inequality corresponding to (9). Consider two scalar characteristics of the matrix $`U`$: its trace, $`\text{tr}U\text{div}𝐯`$ and determinant, $`detU`$. From (18) we derive the following equations for these two quantities:

$$\frac{d}{dt}\text{tr}U=\text{tr}(U^2),$$
(19)
$$\frac{d}{dt}detU=\text{tr}UdetU.$$
(20)
Now we introduce the positive definite integral,

$$I=(detU)^2𝑑𝐫,$$

which, in the 1D case, coincides with $`I_2`$. Due to (20), we have

$$\frac{dI}{dt}=\text{tr}U(detU)^2𝑑𝐫.$$
(21)
This equation generalizes the equation (7) to the multi-dimensional case. The second derivative of $`I`$ will be given by the expression,

$$\frac{d^2I}{dt^2}=\left[(\text{tr}U)^2+\text{tr}(U^2)\right](detU)^2𝑑𝐫.$$
(22)
Applying the Cauchy-Bunyakovsky inequality to the r.h.s. of (20) yields

$$\frac{dI}{dt}I^{1/2}\left(\text{tr}(U^2)(detU)^2𝑑𝐫\right)^{1/2}.$$
(23)
From (22) we write the integral

$$(\text{tr}U)^2(detU)^2𝑑𝐫=\frac{d^2I}{dt^2}\text{tr}(U^2)(detU)^2𝑑𝐫.$$
(24)
Now, we shall estimate the integral on the r.h.s. of this equation. We shall assume that all eigenvalues of the matrix $`U`$ are real. Such an assumption means that the matrix $`U`$ is close to its symmetric part, $`S=1/2(U+U^T)`$, the so-called stress tensor (here $`T`$ denotes transpose). In this case the antisymmetric part of the matrix $`U`$, the vorticity tensor, $`\mathrm{\Omega }=1/2(UU^T)`$ 222The vorticity tensor $`\mathrm{\Omega }`$ is connected with the vorticity $`\omega `$ by the relation $`\mathrm{\Omega }_{ij}=\frac{1}{2}ϵ_{ijk}\omega _k.`$, is small compared to $`S`$. In particular, if $`\mathrm{\Omega }=0`$ the matrix $`U`$ coincides with $`S`$, representing the Hessian of the velocity potential $`\mathrm{\Phi }`$: $`U_{ij}=S_{ij}=^2\mathrm{\Phi }/x_ix_j`$.

Under this assumption the trace of the matrix $`U`$, $`\text{tr}(U^2)=\underset{i=1}{\overset{D}{}}\lambda _i^2`$, where $`\lambda _i`$ are eigenvalues of the $`U`$ and $`D`$ dimension, becomes positive. In this case the following relation between traces of $`U`$ and $`U^2`$ can be easily proven:

$$\underset{i=1}{\overset{D}{}}\lambda _i^2\frac{1}{D}\left(\underset{i=1}{\overset{D}{}}\lambda _i\right)^2.$$

This inequality generates the following estimate between integrals:

$$(detU)^2𝑑r\text{tr}(U^2)(detU)^2𝑑r\frac{1}{D}\left(\text{tr}U(detU)^2𝑑𝐫\right)^2\frac{1}{D}\left(\frac{dI}{dt}\right)^2.$$
(25)
Substitution (24) and (25) into (23) gives the desired differential inequality,

$$I\frac{d^2I}{dt^2}\left(1+\frac{1}{D}\right)\left(\frac{dI}{dt}\right)^20.$$
(26)
Its solution is sought in power form as before: $`I=x^\alpha `$. Excluding terms proportional to $`X^{\alpha 2}X_t^2`$ we find that $`\alpha =D`$ and

$$X=\frac{1}{I^D}.$$

For $`X`$ this results in the same inequality as (10): $`X_{tt}<0`$. The criterion for attaining the origin $`X=0`$ will be also analogous:

$$X_t(0)<0\text{or}I_t(0)>0.$$

Almost the same form will have the estimate for the collapse time

$$t_0<\frac{DI(0)}{I_t(0)}=\frac{1}{\lambda (0)},$$
(27)
where $`\lambda `$ is a mean eigenvalue defined in accordance with (20):

$$\lambda =\frac{1}{D}\underset{i}{}\overline{\lambda }_i=\frac{1}{DI}\text{tr}U(detU)^2𝑑𝐫.$$

For arbitrary $`D`$, instead of (11), the following estimation appears for $`I`$:

$$I(t)\frac{I^{D+1}(0)}{(DI_t(0)(t_0t))^D}.$$
(28)
## 4 Comparison with exact solution

In order to compare the estimation (28) with the exact dependence of $`I`$ we have to solve equation (18). This solution is

$$U=U_0(𝐚)(1+U_0(𝐚)t)^1.$$
(29)
Here $`𝐚`$ is the initial coordinates of a fluid particle and $`U_0(𝐚)`$ is the initial value of the matrix $`U`$. By introducing the projectors $`P^{(k)}`$ of the matrix $`U_0(𝐚)`$ ($`P^{(k)2}=P^{(k)}`$ corresponding to each of the eigenvalues $`\lambda _{0k}(𝐚)`$), this expression can be rewritten in the form of a spectral expansion:

$$U=\underset{k=1}{\overset{D}{}}\frac{\lambda _{0k}}{1+\lambda _{0k}t}P^{(k)}.$$
(30)
The projector $`P^{(k)}`$, being a matrix function of $`𝐚`$, is expressed through the eigenvectors for the direct ($`U_0(𝐚)\psi =\lambda _0\psi `$) and conjugated ($`\varphi U_0(𝐚)=\varphi \lambda _0`$) spectral problems for the matrix $`U_0(𝐚)`$:

$$P_{ij}^{(k)}=\psi _i^{(k)}\varphi _j^{(k)}.$$

where the vectors $`\psi ^{(n)}`$ and $`\varphi ^{(m)}`$ with different $`n`$ and $`m`$ are mutually orthogonal:

$$\psi _i^{(m)}\varphi _i^{(n)}=\delta _{mn}.$$

Hence, the determinant of the matrix $`U`$ is defined by the product,

$$detU=\underset{k=1}{\overset{D}{}}\frac{\lambda _{0k}}{1+\lambda _{0k}t}.$$

From (30) it follows also that singularity in $`U`$ first time appears at $`t=t_0`$, defined from the condition :

$$t_0=\underset{k,a}{\mathrm{min}}[1/\lambda _{0k}(𝐚)],$$
(31)
(compare this with (27)).

From (30), one can see that near the singular point only one term in the sum (30) survives,

$$U\frac{P^{(n)}}{\tau +\gamma _{\alpha \beta }\mathrm{\Delta }a_\alpha \mathrm{\Delta }a_\beta },$$
(32)
where the projector $`P^{(n)}`$ is evaluated at the point $`𝐚=𝐚_\mathrm{𝟎}`$ and $`k=n`$, corresponding to the minimum (31), $`\tau =t_0t`$, $`\mathrm{\Delta }𝐚=𝐚𝐚_\mathrm{𝟎}`$, and

$$2\gamma _{\alpha \beta }=\frac{^2\lambda _{0n}^1}{a_\alpha a_\beta }|_{a=a_0}$$

is a positive definite matrix.

The remarkable formula (32) demonstrates that i) the matrix $`U`$ tends to the degenerate one as $`tt_0`$ and ii) both parts of the matrix $`U`$ in this limit, i.e. the stress tensor $`S`$ and the vorticity tensor $`\mathrm{\Omega }`$, become simultaneously infinite (compare with ). It is interesting to note that at near singular time the ratio between both parts is fixed and governs by two relations following from the definition of the projector $`P`$:

$$P_S=P_S^2+P_A^2,P_A=P_SP_A+P_AP_S$$

where $`P_S`$ and $`P_A`$ are respectively symmetric (”potential”) and antisymmetric (vortical) parts of the projector $`P`$. In particular, the second relation provides the collapsing solution for the equation for vorticity

$$\frac{\omega }{t}=\text{curl}[𝐯\times \omega ]$$

which has the same form for both compressible and incompressible cases. It is also interesting to note that in the sense of the criterion (1), the collapsing solution (32) represents the marginal solution.

The asymptotic solution (32), far from the collapsing point should be matched with a ”regular” solution. The corresponding matching scale can be estimated as $`l_0\gamma ^{2/3}.`$ This scale $`l_0`$ alone can be taken as the size of the collapsing region for (32). This remark now allows one to calculate the contribution from the breaking area to the integral $`I`$.

Substituting (29) into $`I`$ and using a change of variables, from $`𝐫`$ to $`𝐚`$, one can get the expression for this contribution,

$$I=C_V\frac{d^Da}{\tau +\gamma _{\alpha \beta }a_\alpha a_\beta },$$
(33)
where the constant is

$$C=\lambda _{0n}detU_0\underset{kn}{}\frac{\lambda _{0k}}{1+\lambda _{0k}t_0}|_{a=a_0}.$$

The integral is taken over the spherical volume $`V`$ with coordinate center at $`𝐚=𝐚_\mathrm{𝟎}`$ and size $`l_0`$. Introducing a self-similar variable $`\xi =𝐚\tau ^{1/2}`$, one can see that the contribution depends significantly on the dimension $`D`$. At $`D=1`$ this integral behaves like $`\tau ^{1/2}`$ in full correspondence with (17). In this case, the integral over $`\xi `$ is convergent at large $`\xi `$ and is not sensitive to the cut-off size $`l_0`$. In the two-dimensional geometry, however, the integral (33) has a power dependence on $`\tau `$, but a logarithmic dependence on $`l_0`$ arising from integration on $`\xi `$:

$$I\mathrm{log}\frac{l_0}{\tau ^{1/2}}$$

that satisfies the inequality (28).

In the three-dimensional case, the integral (33) diverges at large scales as the first power of $`\xi `$, becomes proportional to the size of collapsing area $`l_0`$ as $`\tau 0`$:

$$Il_0.$$

This result for $`D=3`$ formally contradicts to the blow-up sufficient condition found above. This contradiction indicates only that in the three-dimensional case the blow-up of the integral $`I`$ has no universal behavior near the singular time which should be expected following to the universal asymptotics (32).

## 5 Concluding remarks

Thus, the initial condition,

$$\frac{dI}{dt}(0)>0,$$

represents a sufficient integral criterion for the collapse if the vorticity matrix $`U`$ is small in comparison to the stress tensor. Under this condition, the integral $`I=(detU)^2𝑑𝐫`$ becomes infinite in a finite time.

In turn, divergence of $`I`$ means that a solution can not be continued in the corresponding functional space. In the $`1D`$ case this is the Sobolev space $`H^1()`$. For $`D=2`$, from the inequality,

$$I|U|^4d^2r(|U|^2U_{ik}^2),$$

together with the embedding Sobolev inequality,

$$\left(|U|^4d^2𝐫\right)^{1/4}CU_{H^1(^2)}$$

it follows that $`U_{H^1(^2)}\mathrm{}`$ as $`I\mathrm{}`$. In terms of the velocity, this means that the solution is not continued in the Sobolev space $`H^2(^2)`$.

In the $`3D`$ case it is possible to write the following set of inequalities,

$$I^{1/6}\left(|U|^6d^3𝐫\right)^{1/6}CU_{H^1(^3)},$$

where the second inequality represents the partial case of the Sobolev embedding inequality .

Hence, one can see that $`I\mathrm{}`$ is equivalent to the divergence of the norm (2) for the Sobolev space $`H^2(^3)`$. Thus, the requirements for strong solutions in the hydrodynamic model (3), (4) are different from those for the 3D Euler equation for incompressible fluids. At the moment it is hard to say whether the results presented in this paper contradict to the theorem  (see also ) or not. In any case, it is a very interesting question. It should be added that collapse in incompressible fluids might happen through breaking of vortex lines (there are some arguments both analytical and numerical  in a favor of such point of view). In this case the corresponding norms would also blow up for the same Sobolev space $`H^2(^3)`$.

## 6 Acknowledgments

The author would like to thank to P.L. Sulem for helpful discussion and to R. Grauer who paid the author’ attention to the paper . The author is grateful to the Observatory of Nice, where this work was initiated within the Landau-CNRS agreement. This work was supported by RFBR (grant no. 00-01-00929) and by INTAS (grant no. 00-00292).


text:
# New Relativistic Effects in the Dynamics of Nonlinear Hydrodynamical Waves

\[

## Abstract

In Newtonian and relativistic hydrodynamics the Riemann problem consists of calculating the evolution of a fluid which is initially characterized by two states having different values of uniform rest-mass density, pressure and velocity. When the fluid is allowed to relax, one of three possible wave-patterns is produced, corresponding to the propagation in opposite directions of two nonlinear hydrodynamical waves. New effects emerge in a special relativistic Riemann problem when velocities tangential to the initial discontinuity surface are present. We show that a smooth transition from one wave-pattern to another can be produced by varying the initial tangential velocities while otherwise maintaining the initial states unmodified. These special relativistic effects are produced by the coupling through the relativistic Lorentz factors and do not have a Newtonian counterpart.

\]

More than a hundred years ago, Riemann addressed the problem of the one-dimensional time evolution of a perfect fluid which, at some given initial time $`t=0`$, is composed of two adjacent states characterized by different values of uniform velocity, pressure and density . Since then this problem has been referred to as the “Riemann problem” and is the prototype of the initial value problem for hyperbolic systems of partial differential equations with discontinuous initial conditions. The conclusion reached by Riemann in Newtonian hydrodynamics is that the one-dimensional flow that develops when the barrier separating the initial “left” (or 1) and “right” (or 2) states is removed, will allow for four different and distinct solutions. All of the solutions are composed of nonlinear waves, in the form of either shock waves or rarefaction waves, that propagate in opposite directions and join the two unperturbed left and right states. As a result, three different “wave-patterns” can be produced corresponding to two shock waves ($`2S`$), a shock and rarefaction wave ($`SR`$) and two rarefaction waves ($`2R`$), respectively. Schematically, the fluid solution for any $`t>0`$ can therefore be represented as

$$1𝒲_{}\mathrm{\hspace{0.33em}3}𝒞\mathrm{\hspace{0.33em}3}^{}𝒲_{}\mathrm{\hspace{0.33em}2},$$
(1)
where $`𝒲`$ denotes a shock or a rarefaction wave that propagates towards the left $`()`$ or towards the right $`()`$ of the initial discontinuity, $`1`$ and $`2`$ are the known initial states with “state-vectors” consisting of fixed and independent values of the pressure, rest-mass density and velocity $`𝐔_{1,2}(p,\rho ,v)_{1,2}`$. Similarly, $`3`$ and $`3^{}`$ represent the new hydrodynamic states that form behind the two waves propagating in opposite directions. The two nonlinear waves $`𝒲_{}`$ and $`𝒲_{}`$ are separated by a contact discontinuity $`𝒞`$ across which the rest-mass density is discontinuous, but not the velocity nor the pressure, i.e. $`\rho _3\rho _3^{},p_3=p_3^{},v_3=v_3^{}`$. These three quantities fully describe the states of the fluid in the flow regions $`3`$ and $`3^{}`$ behind the nonlinear waves. The formal solution of the Riemann problem consists then in determining the velocity, pressure and rest-mass density in the new states $`3`$ and $`3^{}`$ as well as calculating, at any time $`t>0`$, the positions of the waves separating the four states.

The solution of the Riemann problem attracted a wider interest when it was realized that its numerical solution could be implemented as building blocks in the construction of numerical Godunov-type methods for the accurate solution of the hydrodynamical equations . In such methods, the computational domain is discretized and each interface between two adjacent grid-zones is used to construct the initial left and right states of a “local” Riemann problem. The evolution of the hydrodynamical equations is then obtained through the solution of the sequence of local Riemann problems set up at the interfaces between successive grid-zones and basically consists of determining the fluid pressure in regions $`3`$ and $`3^{}`$ (see ref.  for a recent review).

The extension of the one-dimensional Riemann problem to relativistic fluid velocities has been performed by  for the “shock-tube” problem and by  for the general problem. Besides the obvious additional complications introduced by special relativity, the one-dimensional Riemann problem in relativistic hydrodynamics does not show qualitative differences from its Newtonian counterpart. In a recent paper , we have introduced a new approach to the solution of an exact relativistic one-dimensional Riemann problem which focuses on the relativistic invariant expression for the relative velocity between the two unperturbed initial states $`(v_{12}^x)_0[(v_1^xv_2^x)/(1v_1^xv_2^x)]_0`$, where we have assumed that the discontinuity is initially located on an $`x=\mathrm{const}.`$ surface and where $`v^x`$ represents the component of the 4-velocity normal to this surface. This new approach introduces great simplifications in the logical formulation of the Riemann problem and has been essential in revealing these new relativistic effects. For this reason we will briefly review it in what follows.

By construction, $`(v_{12}^x)_0`$ measures the relativistic jump of the fluid velocity normal to the discontinuity surface. Within this new approach, the solution to the relativistic Riemann problem is found after the pressure in the region between the two nonlinear waves has been calculated as the root of the nonlinear equation

$$v_{12}^x(p_3)(v_{12}^x)_0=0.$$
(2)
The notation used in equation (2) needs an additional comment. It is should be emphasized, in fact, that $`p_3`$ itself depends on the thermodynamical properties of the initial states so that $`v_{12}^x=v_{12}^x(p_1,\rho _1,p_2,\rho _2)`$. Note also that $`v_{12}^x(p_3)`$ does not depend on the value of the initial tangential velocities and has a functional form that is different for each of the three possible wave-patterns that might result from the decay of the initial discontinuity. In ref.  it was then shown that the wave-pattern produced in a one-dimensional relativistic Riemann problem by the decay of the discontinuity can be predicted from knowledge of the initial data. The possibility of predicting the wave-pattern from the initial conditions applies also to Newtonian hydrodynamics  and represents an important advantage since it allows one to deduce in advance which set of equations to use for the solution of the Riemann problem. From a mathematical point of view, the validity of the approach discussed in ref.  is based on the proof that the three branches of the function $`v_{12}^x=v_{12}^x(p_3)`$, corresponding to the three possible wave-patterns, are monotonically increasing with $`p_3`$. Furthermore, it was shown that the three different branches always join smoothly through specific values of $`v_{12}^x(p_3)`$ denoted respectively as $`(\stackrel{~}{v}_{}^{x}{}_{12}{}^{})_{_{2S}}`$ and $`(\stackrel{~}{v}_{}^{x}{}_{12}{}^{})_{_{SR}}`$. (In ref.  another limiting value, $`(\stackrel{~}{v}_{}^{x}{}_{12}{}^{})_{_{2R}}`$, was found but this is not relevant in the present discussion and will not be considered here.)

Exploiting these properties, an improved exact Riemann solver can be built in which the comparison of $`(v_{12}^x)_0`$ with the relevant limiting values $`(\stackrel{~}{v}_{}^{x}{}_{12}{}^{})_{_{2S}}`$ and $`(\stackrel{~}{v}_{}^{x}{}_{12}{}^{})_{_{SR}}`$ constructed from the initial conditions allows to determine, prior to solving equation (2), both the wave-pattern produced and the functional form of $`v_{12}^x(p_3)`$ to be used in equation (2. In practice, this approach leads to a straightforward numerical implementation and a more efficient numerical algorithm with computational costs reduced by up to 30%.

In generic multidimensional flows, the Riemann problem can still be cast in terms of a planar discontinuity surface across which a normal flow with velocity $`v^x`$ takes place and is responsible for the mass transfer. In general, however, the flow will also have components tangential to the discontinuity surface, which we will refer to as $`v^t[(v^y)^2+(v^z)^2]^{1/2}`$. The latter can also be discontinuous across $`𝒞`$ possibly leading to a Kelvin-Helmoltz instability there. In relativistic regimes, the tangential component cannot be removed by a single Lorentz boost and is constrained to satisfy $`|v^t|<[1(v^x)^2]^{1/2}`$.

A first study of a multidimensional relativistic Riemann problem involving two strong rarefactions was carried out in ref.  and has been extended to a more general treatment in ref.  where it was pointed out that the tangential velocities in multidimensional flows introduce differences between the Newtonian and the relativistic Riemann problems. In Newtonian hydrodynamics, in fact, the solution to the Riemann problem does not depend on the tangential component of the flow. In relativistic hydrodynamics, on the other hand, the different regions of the flow are coupled, both in the velocities and in the specific enthalpy, through the Lorentz factors $`W_i[1(v_i^x)^2(v_i^t)^2]^{1/2}`$, where $`i=1,2`$. This is an important difference whose consequences are much harder to assess.

To discuss the new effects that emerge in a multidimensional and relativistic Riemann problem we have plotted in Fig. 1 the function $`v_{12}^x`$ for different values of the pressure at the contact discontinuity $`p_3`$. The curves shown refer to the case of zero tangential velocities (continuous line) and to the cases where the tangential velocities are nonzero (dashed and dotted lines), respectively (Note that the curves do not depend on the sign of $`v^t`$.). In both cases the relative velocities are calculated for a perfect fluid with polytropic equation of state $`p=k\rho ^\gamma `$, with a standard set of the initial values of pressure and rest-mass density , i.e. $`p_1=1.0,p_2=0.1,\rho _1=1.0,\rho _2=0.125`$ (we have here chosen $`\gamma =5/3`$). The filled dots in each curve indicate the locations where the branches representing the three wave-patterns ($`2S,SR,`$ and $`2R`$) merge. The values of the relative normal-velocities where this happens are just the values of $`(\stackrel{~}{v}_{}^{x}{}_{12}{}^{})_{_{2S}}`$ and $`(\stackrel{~}{v}_{}^{x}{}_{12}{}^{})_{_{SR}}`$ mentioned above.

A rapid examination of Fig. 1 is sufficient to appreciate that the presence of tangential velocities can introduce new qualitative differences in a relativistic Riemann problem. When tangential velocities are present, in fact, the relative normal-velocity $`v_{12}^x`$ is a function of $`p_3`$ but, through the coupling introduced by the Lorentz factors, also of $`v_1^t`$ and $`v_2^t`$. This can introduce a new relativistic effect and shift the solution from one wave-pattern to another. To better appreciate this, let us restrict attention to a simpler situation in which only one of the two initial tangential velocities is allowed to vary.

Suppose now that the normal velocities are chosen to be $`v_1^x=0.5,v_2^x=0.0`$, and that there are no tangential velocities. In this case, $`(v_{12}^x)_0=0.5`$ and Fig. 1 shows that the solution to the Riemann problem falls in the $`SR`$ branch, hence producing a wave-pattern consisting of a shock and a rarefaction wave moving in opposite directions. This is presented in Fig. 2 which shows the solution of the Riemann problem at a time $`t>0`$ for the pressure, the rest-mass density and the velocity.

However, if we now maintain the same initial conditions but allow for nonzero tangential velocities in state 2, Fig. 1 also shows that the solution to the Riemann problem can fall in the $`2S`$ branch (cf. dotted line), hence producing a wave-pattern consisting of two shock waves moving in opposite directions. This is presented in Fig. 3 which shows the solution of the same Riemann problem as in Fig. 2 but with the initial tangential velocities being $`v_1^t=0.0`$ and $`v_2^t=0.9`$. Note that except for the tangential velocities, the solutions in Figs. 2 and 3 have the same initial state-vectors but different intermediate ones (i.e. $`p_3,\rho _3`$, and $`v_3^x`$).

The Riemann problems shown in Figs. 2 and 3 are only two possible examples but they are useful ones since they show that through the coupling among the different states introduced by the Lorentz factors, a difference in the tangential velocities can produce a smooth transition from one wave-pattern to another while maintaining the initial states unmodified. Interestingly, this transition does not need to always produce a solution consisting of two shock waves.

Suppose, in fact, that the normal velocities are now chosen to be $`v_1^x=0.0,v_2^x=0.5`$. We can then repeat the considerations made above and start examining the wave-pattern produced when there are zero tangential velocities. In this new setup, $`(v_{12}^x)_0=0.5`$ and Fig. 1 shows that the solution to the Riemann problem still falls in the $`SR`$ branch (cf. dashed line), with the corresponding solution at a time $`t>0`$ being presented in Fig. 4. Note that the wave-patterns in Fig. 2 and 4 both consist of a shock and a rarefaction wave, but have alternate initial normal velocities.

When nonzero tangential velocities are now considered in state 1, Fig.1 shows that $`(v_{12}^x)_0`$ can fall in the $`2R`$ branch, hence producing a wave-pattern consisting of two rarefaction waves moving in opposite directions. This is again shown in Fig. 5 where we have chosen initial tangential velocities $`v_1^t=0.9`$ and $`v_2^t=0.0`$. In this case too and except for the tangential velocities, the solutions in Figs. 4 and 5 have the same initial state-vectors but different intermediate ones.

As mentioned above, the appearance of these new relativistic effects is related to the behaviour of the function $`v_{12}^x`$ for different values of the initial tangential velocities and in particular to how the three branches composing the curve change under variation of $`v_{1,2}^t`$. As a result, the occurrence of these effects can be recast into the study of the functional behaviour of the values of $`v_{12}^x`$ delimiting the different branches when the tangential velocities are varied within the allowed ranges. A detailed mathematical treatment of this problem will be presented in a separate paper  but it is sufficient to point out in what follows the main results. Firstly, the dependence of the limiting values for the $`SR`$ and $`2S`$ branches on the tangential velocities can be expressed as (cf. Fig. 1)

$$(\stackrel{~}{v}_{}^{x}{}_{12}{}^{})_{_{SR}}=(\stackrel{~}{v}_{}^{x}{}_{12}{}^{})_{_{SR}}(v_1^t),(\stackrel{~}{v}_{}^{x}{}_{12}{}^{})_{_{2S}}=(\stackrel{~}{v}_{}^{x}{}_{12}{}^{})_{_{2S}}(v_2^t).$$
(3)
Secondly, these values converge to zero when the tangential velocities are allowed to reach their maximum values, i.e.

$$\underset{W_1\mathrm{}}{lim}(\stackrel{~}{v}_{}^{x}{}_{12}{}^{})_{_{SR}}=0=\underset{W_2\mathrm{}}{lim}(\stackrel{~}{v}_{}^{x}{}_{12}{}^{})_{_{2S}}.$$
(4)
Expressions (4) indicate that for tangential velocities assuming increasingly larger values, the $`SR`$ branch of the $`v_{12}^x`$ curve spans a progressively smaller interval of relative normal-velocities. When the tangential velocities reach their asymptotic values, the $`SR`$ branch reduces to a point. As a result, a wave-pattern consisting of a shock and a rarefaction wave is generically disfavoured by increasing tangential velocities.

The effects discussed so far have a purely special relativistic origin and might conflict with our physical intuition. However, a behaviour similar to the one reported here is found in the well-known relativistic transverse-Doppler effect, in which the wavelength of a photon received from a source moving at relativistic speeds changes also if the source has a velocity component orthogonal to the direction of emission of the photon . In this case too, a Lorentz factor including the transverse velocity is responsible for the effect.

A couple of comments are worth mention. Firstly, there exists a set of initial conditions for which these new relativistic effects will not occur. These initial conditions are those of the classic “shock-tube” problem, in which $`v_1^x=0=v_2^x`$. In this case, in fact, $`(v_{12}^x)_0=0`$ and, because of the limits (4), the solution of the Riemann problem will be given by a wave-pattern consisting of a shock and a rarefaction wave, independently of the values of the tangential velocities. Secondly, while the stability of shock fronts has been studied in the past (see  and references therein), little is still know about the stability properties of shock fronts across which the tangential velocities are discontinuous. It is plausible that a nonzero mass flux makes these fronts stable with the respect to the Kelvin-Helmoltz instability, but more work needs to be done in the investigation of the stability of generic shocks fronts.

The results reported in this letter can be used in astrophysical scenarios, such as those involving relativistic jets or $`\gamma `$ray bursts, in which nonlinear hydrodynamical waves with large Lorentz factors and complex multidimensional flows are expected . Furthermore, both the approach to the solution of the Riemann problem and the new relativistic effects discussed here could be relevant in the implementation of special relativistic exact Riemann solvers in multidimensional codes solving the hydrodynamical equations in curved spacetimes .

It is a pleasure to thank J.C. Miller and J.A. Pons for useful discussions. Support for this research has been provided by the MIUR and by the EU Network Programme (Research Training Network Contract HPRN-CT-2000-00137).