The construction of orthonormal wavelets using symbolic methods and a matrix analytical approach for wavelets on the interval Frederic Chyzak Peter Pauley Otmar Scherzerz Burkhard Zimmermanny Armin Schoisswohlz March 27, 2000 Abstract
In this paper we discuss closed form representations of lter coecients of wavelets on the real line, half real line and on compact intervals. We show that computer algebra can be applied to achieve this task. Moreover, we present a matrix analytical approach that uni es constructions of wavelets on the interval.
1 Introduction Wavelets are one of the most popular tools in signal and image processing. These functions are widely used in many practical applications such as data compression [BBH93, DJL92, SSK98], or for the solution of partial di erential equations (see, e.g., [Jaf92]). Wavelets are special functions which often have a fractal character. This makes it relatively dicult to work with them explicitly; for example, point evaluation of a wavelet function may already be a computationally expensive task. To work with wavelets one uses the nice feature that they are de ned by a small number of parameters, the so-called lter coecients. In general, any algorithm relying on wavelets uses the lter coecients only and not the wavelet function itself. In this paper we review the basic equations for the lter coecients. We show that these equations can be solved using computer algebra. In particular we can construct closed form representations of the wavelet coecients (see section 3). The most popular wavelets form an orthonormal basis of the space of square integrable functions on R (Daubechies [Dau92]). In many practical applications one requires a basis on the half-line or on a compact interval. In section 4 INRIA Rocquencourt, F{78153 Le Chesnay, France. (
[email protected]) y Research Institute for Symbolic Computation, Johannes Kepler University, A{4040
Linz, Austria. (fPeter.Paule,
[email protected]) z Industrial Mathematics Institute, Johannes Kepler University, A{4040 Linz, Austria. (fscherzer,
[email protected])
1
we review several constructions of wavelets on the interval. We construct the lter coecients of wavelets on the interval using a matrix analytical approach which has the advantage to unify several constructions in the literature [Mey91, CDV93, DKU99]. Moreover, our construction reveals that there exists a closed form representation of the lter coecients of wavelets on the interval. Another big advantage of using computer algebra is that one can avoid instabilities which occur in numerical calculations of lter coecients. In sections 4.5 and 5 we summarize the algorithms for calculating closed form coecients of wavelets on the real line and on compact intervals and present some results.
2 Orthonormal Wavelets on R The construction of compactly supported orthonormal wavelet bases for L2 (R ) is well understood (see, e.g., [Dau88, Dau92, CD93, Dau93]) but is summarized here for the reader's convenience. The construction of wavelets is related with the construction of a scaling function such that for xed m 2 Z the functions m;k := 2?m=2 (2?m x ? k), k 2 Z, are orthonormal with respect to L2 (R ). Moreover, the spaces Vm := spanf m;k ; k 2 Z g constitute a multiresolution analysis for L2 (R ), i.e.,
Vm Vm?1 ; for m 2 Z, with
\ m2Z
Vm = f0g and
[ m2Z
Vm = L2 (R ):
The wavelet spaces Wm are the orthogonal complements of Vm in Vm?1 , i.e.,
Wm := Vm? \ Vm?1 : One de nes the wavelet such that the functions m;k := 2?m=2 (2?m x ? k), k 2 Z, form an orthonormal basis for Wm . Since both Vm and Wm are contained in Vm?1 the scaling function must satisfy the dilation equation (1)
(x) =
X
k2Z
hk (2x ? k);
where the sequence fhk g is known as the lter sequence and satis es constraints to be recalled below. Correspondingly, the wavelet satis es (2)
(x) =
X
k2Z
gk (2x ? k);
where gk = (?1)k h1?k . Daubechies [Dau88] established conditions on the lter sequence in order to ensure that the dilation equation (1) has a solution 2 L2 (R ), with supp = [ 1 ? N; N ] for a given integer N , and that for xed m the functions m;k are 2
orthogonal with the property that polynomials up to degree N ? 1 can be represented as linear combinations of m;k . Compact support of in [ ?N + 1; N ] is ensured by
hk = 0;
(3)
(k < 1 ? N or k > N ):
A requirement for the existence of a solution of (1) is (4)
N X k=1?N
hk = 2;
(l = 0; : : : ; N ? 1);
R is equivalent to (x) dx = 1. Orthonormality of the translates of , i.e., Rwhich (x) (x ? l) dx = , can be translated into 0;l
(5)
N X k=1?N
(l = 0; : : : ; N ? 1);
hk hk?2l = 20;l ;
and requirement that polynomials are representable by the m;k leads to R xl the (x) dx = 0 for l = 0; : : : ; N ? 1, which is equivalent to (6)
N X k=1?N
(?1)k h1?k kl = 0;
(l = 0; : : : ; N ? 1):
3 Closed Form Representation of Filter Coecients In this section we reconsider the calculation of the lter coecients hk from equations (3{6) by using methods of computer algebra. Below we give a brief and informal account on Grobner bases, which we exemplify by the calculation of the lter coecients for the special case N = 2. Afterwards we pass on to the cases N > 2 and present more computational details on our symbolic approach.
3.1 Grobner Bases First of all note that due to the conditions imposed explicitly on the summation bounds in the equations (4{6), we can restrict our attention to the task of solving only those; the conditions (3) can be satis ed separately by mere de nition. But instead of solving the equations (4{6) numerically and for a xed integer value N , we try to nd closed forms for the coecients hk , i.e., to approach the problem from the symbolic computation point of view. To solve systems of polynomial equations symbolically, the obvious tools to use are Grobner bases: after their original introduction by B. Buchberger [Buc65] to answer ideal-theoretic questions, the solving of algebraic systems was soon realized to be one of their natural domains of application [Buc70]. For further introductory information see, e.g., [Win96] or [VzGG99]; an excellent source for additional references and for the state of art is [BW98]. 3
The case N = 1 is trivial; h0 = h1 = 1 is the only solution. Hence we illustrate the Grobner bases method for N = 2. In this case we are interested in all common roots of the ve polynomials in the four variables x1 ; x2 ; x3 ; x4 : (7)
?2 + x1 + x2 + x3 + x4 ; ?2 + x21 + x22 + x23 + x24 ; x1 x3 + x2 x4 ; x1 ? x2 + x3 ? x4 ; 2x1 ? x2 + x4 ;
for the sake of simplicity we introduced the following renaming of variables:
x1 = h?1 ; x2 = h0 ; x3 = h1 ; and x4 = h2 : Let I be the ideal in the polynomial ring C [x1 ; x2 ; x3 ; x4 ] generated by the
polynomials from (7). Applying Buchberger's algorithm with respect to a certain order imposed on the monomials of C [x1 ; x2 ; x3 ; x4 ] (here: \lexicographic with x4 > x3 > x2 > x1 ") delivers an alternative description of the ideal I , namely by the generators: (8)
?1 ? 4x1 + 8x21 ; ?1 ? 2x1 + 2x2 ; ?1 + x1 + x3 ; ?1 + 2x1 + 2x4 :
The polynomials (8) again generate the ideal I , and in particular, share the same variety of common roots as the generators from (7). But additionally, they form a Grobner basis of I . Due to the choice of a lexicographic monomial order, they furthermore possess the following \elimination property": the rst polynomial in the Grobner basis is a univariate polynomial (here in x1 ), the second one a bivariate polynomial that involves only one further variable (here x1 and x2 ), the third one a polynomial in three variables (here x1 , x2 , and x3 ), and so on. In other words, the role of the Grobner basis algorithm in solving systems of algebraic equations is the same as that of Gaussian elimination in solving systems of linear equations, namely to triangularize the system or to carry out the elimination, respectively. Remarkably in our situation of solving lter coecient equations, an even nicer pattern emerges. Namely, given the rst univariate Grobner basis polynomial p1 (x1 ) in x1 only, the second Grobner basis polynomial is the sum of a univariate polynomial in x1 and a linear polynomial in x2 ; the third Grobner basis polynomial is the sum of a univariate polynomial in x1 and a linear polynomial in x3 , and so on. This means that all other lter coecients xi for i > 1 nd a representation of the form
xi = pi(x1 ); where each pi (x1 ) is a polynomial from C [x1 ], i.e., depending on x1 only. Con(9)
sequently, there are as many di erent solutions of a system of lter coecient equations as there are di erent roots of the rst univariate Grobner basis polynomial p1 (x1 ). So far we have observed the nice pattern of a univariate polynomial and relations like (9) for all values of N up to 6, so that we conjecture that this situation also holds for arbitrary N . For readers interested in ideal theory we state this in the form of the following conjecture. (For more ideal-theoretic background see, for instance, the \shape lemma" in [Win96].) Conjecture 1. Polynomial ideals corresponding to Daubechies lter coecient equations are 0-dimensional and radical. 4
Also note that by choosing di erent lexicographic orders on the xi , it would be possible to obtain alternative descriptions of the solutions, parameterized by another choice of xj . Also, since we expect a representation of the form of (9), other non-lexicographic orders could be used, with enhanced eciency. (Speci cally, any order which sorts x2 , x3 , and x4 by any non-lexicographic order, but sorts them lexicographically higher than x1 .) In the same vein, extensive calculations on which we cannot report here in details have shown that for N = 2; : : : ; 6, the ideal generated by the system (4{6) contains a nonzero univariate polynomial in hi of degree 2N ?1 for each i between 1 ? N and N , and none of smaller degree. (See related explicit results in section 5.2.) To conclude this informal discussion of the Grobnerpbases approach,pwe state the solution of the case N = 2 explicitly. Since (1 + 3)=4 and (1 ? 3)=4 are the roots of the rst Grobner basis polynomial ?1 ? 4x1 + 8x21 , we obtain two solutions for the lter coecients: p p p p! 3 3 3 3 1 + 3 + 3 ? 1 ? (x1 ; x2 ; x3 ; x4 ) = 4 ; 4 ; 4 ; 4 and p p p! p 3 3 3 3 : 3 ? 3 + 1 + 1 ? (x1 ; x2 ; x3 ; x4 ) = 4 ; 4 ; 4 ; 4
3.2 Reduction of the Filter Coecient Equations For xed N , the system (4{6) consists of 2N + 1 equations in 2N unknowns. For each N such that the system admits a nice triangular representation of the form (9), 2N is also the number of Grobner basis polynomials we nally have to solve explicitly. In this section, as an important preprocessing step to Grobner basis computation, we transform the system (4{6) into a more economic form. More precisely, this system will consist of only N equations in N unknowns; the corresponding Grobner bases will then consist of N polynomials for which one again observes the nice shape that was described above (see the discussion preceding Conjecture 1). Not only is this system more compact, but it also induce faster Grobner basis calculations in practice, as was suggested by the doubly exponential upper bound of the degree of polynomials in a Grobner basis in terms of the number of variables. In a rst reduction step we introduce a normalization via multiplication by a binomial coecient; namely, for any xed positive integer N we de ne ak by 2N ? 1 (10) a ; (k = 1 ? N; : : : ; N ): h = k
N ?k
N ?k
This implicitly installs conditions (3) and thus enables to relax the explicit statement of the summation bounds in (4{6). More importantly, a second change of variables will prove successful in the sequel: for xed positive integer N we restrict ourselves to consider ak as a polynomial in k of degree at most N ? 1. To this end, we write (11)
ak =
NX ?1
Pj kj ; j =0 5
where the Pj are the new unknowns we have to solve for. Note that we have in total N of those, instead of 2N in the original setting? (4{6). In addition, we shall see below why it is convenient to work with the kj as basis elements instead of the kj . With ansatz (10) and (11) in hand, we return to equations (3{6). It is not dicult to see that only two of them remain: (3) is guaranteed due to the presence of the binomial coecient in (10); also, equation (6) is satis ed for arbitrary l = 0; : : : ; N ? 1 because of the following elementary combinatorial lemma (see for instance [GKP94, (5.42)]). Lemma 1. For any nonnegative integer n and complex numbers ci : (12)
n X
(?1)k
n
k=0
n n k (c0 + c1 k + + cn k ) = (?1) n! cn :
Now, with ansatz (10), equation (6) for 0 l N ? 1 is rewritten as 2X N ?1
k=0
(?1)k?N +1 2N ? 1
k
ak (k ? N + 1)l = 0;
and both ak and (k ? N +1)l are polynomials in k with degree less than or equal to N ? 1. Hence, by Lemma 1, equation (6) is satis ed for all l in question. In order to state new versions of the remaining equations (4) and (5) in the form it is convenient to renormalize Pj by introducing Qj = ?2N ?1ofPpropositions, . The nal ansatz now becomes j j
hk =
(13)
NX ?1 j =0
2N ? j ? 1
Qj N + k ? 1 ;
after substituting (11) into (10) and using the elementary fact
2N ? 1k 2N ? 12N ? j ? 1 = : k
j
j
2N ? k ? 1
Proposition 2. Under the assumption (13), equation (4) is equivalent to (14)
NX ?1 j =0
22N ?j ?2 Qj = 1:
P? Proof. The proposition follows from l 2N ?l j ?1 = 22N ?j ?1 , a special instance of the binomial theorem. Proposition 3. Under the assumption (13), equation (5) is equivalent to (15)
NX ?1 4N ? i ? j ? 2
i;j =0 2N + 2l ? i ? 1
Qi Qj = 20;l ;
6
(l = 0; : : : ; N ? 1):
Proof. Substituting (13) twice in (5) yields
N 2N ? j ? 1 2N ? i ? 1 X Qi Qj = 20;l : i;j =0 k=1?N N + k ? 1 N + k ? 2l ? 1 NX ?1
The inner sum can be evaluated? as follows: ? after changing k into N ? k and applying the binomial symmetry mn = n?nm to the second binomial, it becomes 2X N ?1
k=0
2N ? j ? 1 2N ? k ? 1
2N ? i ? 1 4N ? i ? j ? 2 = ; k + 2l ? i
2N + 2l ? i ? 1
where the last identity is a variant of the standard Vandermonde summation (e.g., [GKP94, (5.22)]). The preceding derivation is clearly invertible. The number of equations can be reduced further. In order to prove this, we need another elementary combinatorial result. Lemma 4. For non-negative integers m and n such that m + n 1: X m + n m+n?1 =2 : l m + 2l
P?
P
?
Proof. We have l m+l n = 2m+n and l (?1)l m+l n = 0 as a result of the binomial theorem. Taking the sum and the di erence of both identities yields
X m + n X m + n m+n?1 = =2 : l
2l
2l + 1
l
Now, the sum in the claim is one of the two sums above, depending on the parity of m, whence the result. Remark. Let us mention that proofs of binomial summations like the Vandermonde formula or Lemma 4 can nowadays be carried out in a purely automatic fashion due to Zeilberger's summation machinery [PWZ96]; see, for instance, the Mathematica package [PS95]. We are now ready to carry out the last reduction step. As opposed to Propositions 2 and 3 that relate identities between the hi and the Qi , the following proposition states an reduction between the equations in the Qi only. Proposition 5. Under the simultaneous assumption of the cases l = 1; : : : ; N ? 1 in (15), the case l = 0 in equation (15) is equivalent to (16)
12 0N ?1 X @ 22N ?j?2Qj A = 1; j =0
and is therefore a consequence of equation (14).
7
Proof. By Ll we denote the double sum on the left hand side of (15). We assume the cases l = 1; : : : ; N ? 1 in (15). Because of the symmetry property Ll = L?l and after applying Lemma 4, we have
L0 =
NX ?1 l=1?N
Ll =
NX ?1 i;j =0
0 1 NX ?1 Q 2 j A: 24N ?i?j ?3 Qi Qj = 2 @22(N ?1) j j =0 2
This proves that the relation L0 = 2 is equivalent to equation (16). Finally we summarize our reduction of the 2N + 1 Daubechies equations in 2N unknowns to N algebraic equations in N unknowns in the form of the following proposition. Proposition 6. Any solution of the N algebraic equations NX ?1 Q j = 1 ; j 22N ?2 j =0 2
(17) and
(18)
NX ?1 4N ? i ? j ? 2
i;j =0 2N + 2l ? i ? 1
Qi Qj = 0;
(l = 1; : : : ; N ? 1);
gives rise to a solution of the Daubechies lter coecient equations (3{6) via
(19)
hk =
NX ?1 2N ? j ? 1 j =0
N + k ? 1 Qj ;
(k = 1 ? N; : : : ; N ):
Conversely for any solution h of (3{6), any solution Q of (19) is a solution of (17{18).
At this stage, the rst part of the proposition provides us with a means to obtain solutions of (3{6). However, this procedure may miss solutions h: for h satisfying (3{6), we have no proof yet that the system (19) is solvable as a system in Q. But this gap can also be closed, leading to the following more speci c \Equivalence Theorem," whose proof is postponed to the next section. Theorem 7 (\Equivalence Theorem"). Both systems (3{6) and (17{18) are equivalent descriptions of the same algebraic variety. More precisely, there exists an explicit linear isomorphism between C N and the N -dimensional linear subspace de ned by (6) in C 2N which realizes a linear change of coordinates between the solution set of the algebraic system (4{6) and the solution set of the algebraic system (17{18). In particular, this theorem provides a bijection between the solution sets of both systems; additionally the whole structures of these solution sets, including dimension and multiplicities, are the same. 8
Remark. Concerning the aspect of experimental mathematics it might be instructive to note that concrete Grobner bases computations have led us to conjecture the Equivalence Theorem. Namely, it turned out that the Grobner bases computed with respect to the system (17{18) have the same nice triangulation property as those computed with respect to the system (4{6). In addition we observed that in all instances the rst univariate Grobner basis polynomial is the same in both cases.
3.3 Proof of the Equivalence Theorem This section is devoted to the proof of Theorem 7. Viewing (4{5) as de ning an algebraic variety H in C 2N and (17{18) as de ning an algebraic variety Q in C N , we are about to show that (19) installs a linear isomorphism between the N -dimensional linear subspace K de ned by (6) in C 2N and C N . In this way H \ K and Q can be viewed as the same algebraic variety in C N ' K, but expressed in di erent linear bases. The proof consists of three steps: rst, we introduce an injective linear map from C N to C 2N , which embeds Q into H; next, we introduce a surjective linear map from C 2N to C N , whose kernel is K; nally, by showing = 0, we obtain the isomorphism between the algebraic varieties Q and H \ K.
The Injective Linear Map . Following the idea of embedding Q into H, we introduce the mapping
: C N ! C 2N which maps Q = (Q0 ; : : : ; QN ?1 )t to h = (h1?N ; : : : ; hN )t given by hi =
NX ?1 2N ? j ? 1 j =0
N + i ? 1 Qj :
Let us also denote = ( i;j ) the matrix of this linear map on the canonical bases of C N and C 2N , where the indices i and j range from 1 to 2N , and 1 to N , respectively. In view of the unusual indexing in h and Q, we have
hi?N = and by identi cation
i;j =
N X j =1
i;j Qj?1 ;
2N ? j i?1
for 1 i 2N and 1 j N . Therefore,
S = U
for two N N square blocks S and U . The antidiagonal of U is obtained when i ? N = N + 1 ? j , i.e., i ? 1 = 2N ? j , for which i;j = 1; the triangular 9
block under the antidiagonal of U is obtained when i ? N > N + 1 ? j , i.e., i ? 1 > 2N ? j , for which i;j = 0; the block U is therefore upper anti-triangular with \1"'s on the antidiagonal. It follows that has maximal rank, namely N . By the rank theorem we obtain that dim ker = 0 and thus that is injective. Now, let us consider Q 2 Q. In other words, Q satis es (17{18), so that by the rst part of Proposition 6, equations (4{6) are satis ed by h = Q. In particular, equations (4{5) are satis ed, so that h 2 H. Summarizing, the variety Q de ned by (17{18) is injectively mapped by into the variety H de ned by (4{5): (Q) H:
The Surjective Linear Map . Following the idea of viewing the variety de ned by (4{6) as an intersection by a suitable kernel K, we introduce the mapping
: C 2N ! C N which maps h = (h1?N ; : : : ; hN )t to z = (z0 ; : : : ; zN ?1 )t given in analogy with (6) by
zi =
N X j =1?N
(?1)j h1?j j i :
Again, we also denote = ( i;j ) the matrix of this linear map on the canonical bases of C 2N and C N , where the indices i and j range from 1 to N , and 1 to 2N , respectively. In view of the unusual indexing in z and h, we have
zi?1 = and by identi cation
2N X
j =1
i;j hj?N ;
i;j = (?1)N +1?j (N + 1 ? j )i?1 for 1 i N and 1 j 2N . Therefore, since the signs only depend on j , = (?1)N diag2N (1; ?1; : : : ; 1; ?1) where diagr (1 ; : : : ; r ) denotes the r r-diagonal matrix with entries i on the diagonal, and where is a matrix ?
= V W consisting of two N N square Vandermonde blocks V and W . For example, the rst block V is the Vandermonde matrix of the powers of N , on the rst column, to the powers of 1, on the last column. It follows from the non-nullity of the Vandermonde determinant that has maximal rank, namely N , and is thus surjective. Now, by the rank theorem the kernel K = ker has dimension N and by construction is precisely the linear subspace of C 2N de ned by (6). The algebraic variety de ned by (4{6) is therefore included in this kernel K; it is H\K, whose study is the topic of the next paragraph. 10
The Intersection H\K. Since is injective, it realizes a linear isomorphism between C N and the image (C N ). Proving that this image is the kernel K of now suces to obtain that Q is in bijection with the intersection H \ K by this linear isomorphism . This bijection is an isomorphism of algebraic varieties. To this end, let us compute z = Q for Q 2 C N : for 0 ` N ? 1,
z` =
N X
i=1?N 2X N ?1
(?1)i i` h1?i =
2X N ?1
i=0
(?1)i+1?N (i + 1 ? N )` hN ?i
NX ?1 2N ? j ? 1 N +1+ i ` Qj (?1) (i + 1 ? N ) = j =0 2N ? i ? 1 i=0 2N ? j ? 1 NX ?1 2X N ?1 N +1 i (i + 1 ? N )` (?1) = (?1) Qj i ? j j =0 i=0 2N ? j ? 1 2NX ?1?j NX ?1 i N +1+ j (i + j + 1 ? N )` : (?1) Qj = (?1) i i=0 j =0
We contend to say that the inner sum is zero, whence z = 0 and = 0. Indeed, this is a consequence of Lemma 1: the result is obtained when considering identity (12) for n = 2N ? j ? 1 and polynomials in i (taking the role of k in (12)), and in view of the fact that (i + j + 1 ? N )` is a polynomial in i of degree ` with ` < n since ` N ? 1 and j N ? 1, therefore such that cn = 0. We already know from the analysis of that (Q) H, and wish to obtain an equality. We have just obtained = 0, so that
(C N ) ker = K: Now, as already mentioned, by the injectivity of and the rank theorem applied to , the dimensions of (C N ) and ker are both n. Thus the equality (C N ) = ker = K: >From it and the known inclusion (Q) H therefore follows (Q) H \ K: For the converse inclusion, let h 2 H \ K. In other words, h satis es equations (4{6). By the second part of Proposition 6, equations (17{18) are satis ed by Q, so that Q 2 Q. Finally,
(Q) = H \ K: Remark. The inverse mapping of (from ker to C N ) is given by (20)
Qi =
N X j =1
(?1)N +i+j
N +j ?1
2N ? i ? 1 hj
for 0 i N ? 1, as can be obtained from elementary manipulations of binomials and the special case (1 + (?1))n of the binomial theorem. 11
3.4 Computational Aspects and a Conjecture We conclude this section with a few comments on various computational aspects of our approach, together with a conjecture that has been suggested by the observation of the lter coecients obtained for small values of N . First of all the reduction of the original system to N equations in N unknowns enables the computation of the corresponding Grobner bases up to N = 6. We used the computer algebra system Mathematica and the built-in procedure GroebnerBasis (see section 5). The case N = 6 takes about 20 seconds and uses less than 10 MB of memory on a computing platform equipped with a Pentium II processor; the case N = 7 is intractable within 128 MB. Another important aspect concerns the remark at the end of section 3.2 that the Grobner bases computed with respect to the N equations (17{18) from Proposition 6 have the same nice triangulation property as those computed with respect to (4{6). We make this link explicit in the form of a proposition, for which we need the following de nition: a polynomial explicit representation (with respect to x1 ) of a 0-dimensional ideal of C [x1 ; : : : ; xr ] is a generating system of the ideal of the form fp1 (x1 ); x2 ? p2 (x1 ); : : : ; xr ? pr (x1 )g for univariate polynomials pi . When an ideal has such a representation, the r generators are a Grobner basis of the ideal for any term order such that x1 < xi for i > 1, and the polynomial pi for i > 1 can be chosen of smaller degree than p1 . Proposition 8. Fix an integer N . Then, the algebraic system (4{6) admits a polynomial explicit representation with respect to hN if and only if the algebraic system (17{18) admits a polynomial explicit representation with respect to Q0 . Additionally, the univariate polynomials in both representations are equal (up to a renaming of indeterminates) when this property is satis ed. Note that in this proposition we consider univariate polynomials in hN , instead of univariate polynomials in h1?N = x1 as in section 3.1. Explicit values of the common univariate polynomial under consideration are provided for N between 2 and 6 in section 5.2. Proof. For a polynomial explicit representation (21) fp0 (hN ); hN ?1 ? p1 (hN ); : : : ; h1?N ? p2N ?1 (hN )g with respect to hN of the ideal generated by the system (4{6), we apply the mapping of section 3.3, which maps a solution h of (4{6) to a solution Q of (17{18). We let h be a solution of (4{6). First, setting i = 0 in (20) yields hN = Q0 , so that p0 (Q0 ) = p0 (hN ) = 0. Next, for any i > 0, equation (20) rewrites
Qi =
N X
(?1)N +i+j 2NN+?ji??11 pN ?j (hN ) = qi (Q0 ) j =1
for some polynomial qi . To prove that this makes fp0 (Q0 ); Qi ? q1 (Q0 ); : : : ; QN ?1 ? qN ?1 (Q0 )g 12
a polynomial explicit representation of the system (17{18) with respect to Q0 , there only remains to show that this set generates the same ideal as (17{18), or equivalently, that the latter ideal contains no univariate polynomial p~0 in Q0 of degree smaller than p0 . If it were so, by considering the mapping of section 3.3, we would obtain the polynomial p~0 (hN ) in the ideal generated by (4{6), contradicting that (21) is a polynomial explicit representation for this ideal. The converse claim is proved similarly, using (19) instead of (20). Since the Grobner bases for both systems seem to be of the same \triangular" shape with a common univariate polynomial, the degree of this polynomial is a bound on the number of solutions. In all the cases N = 1; : : : ; 6, it turns out to be of degree 2N ?1 . Besides, this number 2N ?1 is also the product of the degrees of the polynomials in system (17{18), the so-called Bezout bound of the system, which would bound the number of ane solutions if we could prove that this system has no solution at in nity (see [CLO98]). Conjecture 2. Both systems of algebraic equations (4{6) and (17{18) have at most 2N ?1 di erent solutions. In order to elaborate on the Bezout bound of the system, let us consider the system consisting of (18) and the variant of (17) obtained by setting its righthand side to 0. A proof that this system has no solution but the trivial solution Q0 = = QN ?1 = 0 would immediately imply by Bezout's theorem [CLO98] that both systems (3{6) and (17{18) have exactly 2N ?1 complex solutions counted with multiplicity. In other words, this would prove conjecture 2 and the zero-dimensionality in conjecture 1. To further corroborate conjecture 2 on the number of solutions, we also performed the calculations of Grobner bases for a total degree order up to N = 34. With this choice of an order, the bases computed do not have the elimination property and the triangular shape discussed in section 3, but suce to derive the dimension and degree of the algebraic systems. To this end we used the specialized software Gb by Jean-Charles Faugere [Fau94], with which each basis was obtained in a matter of seconds. (We could not go further due to a limitation in the size of integers.) The result is that in each case the variety has dimension 0 and degree 2N ?1 , which proves the conjecture up to N = 34. Also here more seems to be true. For instance, up to N = 6 the common univariate Grobner basis polynomial always has 2N ?1 di erent solutions. In particular, we have two real solutions in the cases N = 2 and N = 3; four real solutions in the cases N = 4 and N = 5; and eight real solutions if N = 6. In section 5.1 we give the Mathematica procedure we have used together with some Grobner bases output. Univariate polynomials in indeterminates other than hN or Q0 have been obtained by further extensive calculations for 2 N 6 (see section 5.2). The observation is that for each hi and each Qi but QN ?1 , the degree of the univariate polynomial is 2N ?1 , whereas for QN ?1 it is only 2. To end this section, we display two of the four solutions corresponding to the case N = 3. In order to obtain those, rst one has to nd all solutions of the univariate Grobner basis polynomial, which is: p1 (x1 ) = 9 ? 96x1 ? 1536x21 ? 4096x31 + 16384x41 : 13
>From this, one computes the two real solutions: a rst real value for the list (h?2 ; h?1 ; h0 ; h1 ; h2 ; h3 ) is
p
p p
p
p p
p
p p
p
p p
1+ 10+ 5+2 10 ; 5+ 10+3 5+2 10 ; 5? 10+ 5+2 10 ; 16 16 8
p
p p
p
p p
5? 10? 5+2 10 ; 5+ 10?3 5+2 10 ; 1+ 10? 5+2 10 8 16 16
;
and a second one is obtained by re ection. The four real solutions obtained for N = 4 are also expressible as explicit expressions in nested radicals, but are too large to be displayed here. The presence of 2N ?1 di erent solutions that can be expressed in terms of nested square roots for 2 N 4 suggests that this could hold for all N . However, Klappenecker seemingly proved by a Galois-theoretic result that the scaling coecients of the Daubechies wavelet cannot be expressed by radicals for all N between 6 and 100 (see [Kla97]).
4 Wavelets on the Interval 4.1 Meyer's Construction To our knowledge the rst construction of orthogonal wavelets on the interval was proposed by Meyer [Mey91]. His construction restricts compactly supported orthonormal wavelets on R (as considered in section 2) to the interval I := [ 0; 1 ] and manipulates the restricted functions in such a way that they form an orthonormal basis on I . Following his construction we develop a matrix analytical approach that allows to unify several constructions of wavelets on the interval in the same framework. To avoid notational diculties we restrict our attention to the construction of wavelets on R + := [ 0; 1). >From our presentation it will become evident how the construction can be generalized to obtain a wavelet basis for L2 (I ). We introduce the family of scaling functions restricted to R +
half (x) :=
(
0 if x < 0 m;k (x) if x 0
m;k
and the corresponding spaces
n
o
Vmhalf := span half m;k ; k 2 Z : The spaces Vmhalf form a multiresolution analysis for L2 (R + ). The corresponding wavelet spaces Wmhalf are given by
?
(22) Wmhalf := Vmhalf ? \ Vmhalf ?1 : As a consequence, we have the relation Vm?1 = Vm Wm , where the direct sum is orthogonal. By PWmhalf and PVmhalf we denote the orthogonal projection operators onto the spaces Wmhalf and Vmhalf , respectively. 14
Since the scaling function has support in [ 1 ? N; N ], half m;k = 0 for k ?N and half = for k N ? 1. m;k m;k Wavelets on R + can be constructed in the following way:
n
o
1. Orthonormalize the set of functions half ; k 1 ? N . The obtained m;k n o orthonormal basis of Vmhalf is denoted by edge ; k 1 ? N . m;k 2. Compute PWmhalf edge m?1;k and orthonormalize them to obtain an orthonoredge mal basis m;k of Wmhalf . half Orthonormalized scaling functions half m;k . The functions m;k are ortho-
normalized by making a basis transformation
half edge m = A m ;
(23) where
0 edge 1 m;?N +1 edge B C := edge m;?N +2 A @ m . ..
Using the notation
and
?
0 half 1 m;?N +1 half B m;?N +2 C := half A: @ m . ..
t := k ; l k;l1?N ; for any vector , we see that the orthonormality of the functions edge m is equivalent to the matrix equation edge t edge m m = I:
>From (23) it follows that half t t I = A half m m A :
If the matrix (24)
half t := half m m
of the inner products of the truncated scaling functions half m;k is known, then the matrix A in (23) can be obtained by the Cholesky factorization (25)
? ? = A?1 A?1 t ;
where A?1 is regular and of lower triangular form. Therefore the matrix A is also lower triangular. In view of (23) and of the supports of the functions half m;l , this in particular ensures staggered support of the functions edge , i.e., supp edge m;k m;k [ 0; 2m (N + k) ]. In the following we derive the re nement equations (similar to (1) and (2)) for edge edge m and the corresponding wavelets m . These equations are the basis for 15
the implementation of multiresolution cascade algorithms [Mal89], as they are used in data compression (see, e.g., [WA94]). By truncation of the dilation equations (1), the truncated scaling functions half m;k satisfy the dilation equation 1 X h half ; half r?2k m?1;r m;k = p2 r2Z which in matrix form rewrites half half m = H m?1 :
(26)
Here the dilation matrix H makes the transition from the truncated scaling functions at scale m ? 1 to those at the ner scale m, and is therefore called the re nement matrix. This is a 1 2 block-Toeplitz matrix:
Hk;l =
( hlp? k 2
0
2
if 1 ? N l ? 2k N , otherwise;
in other words, Hk;l is zero outside a band of slope 1/2. >From (23) and (26) it follows that half half ?1 edge edge m = A m = AH m?1 = AHA m?1 :
(27)
Thus the re nement matrix H edge for the dilation equation of the edge scaling functions edge m is
H edge = AHA?1 :
(28)
Note that H edge is no longer a 1 2 block-Toeplitz matrix (which is the case for wavelets on R ). This re ects the fact that the edge scaling functions cannot be obtained as shifts of a single function.
Projection onto Wmhalf . Now we construct the edge wavelets and derive their half re nement matrix. The projections of edge m?1;k onto Wm , which has been de ned by (22), are given by edge
edge
m;k := PWmhalf m?1;k = m?1;k ? half
or equivalently with matrices edge half m = m?1 ?
(29)
edge half m = m?1 ?
?
l
edge edge edge m?1;k ; m;l m;l ;
edge t edge : edge m m?1 m
In view of (27) and (28) it follows that (30)
X
edge t edge t H edge edge edge m?1 m?1 H m?1
half edge = I ? H edge t H edge edge m?1 =: G m?1 :
16
The matrix Ghalf does not have 1 2 lower block-triangular form, i.e., it does half not ful ll Ghalf k;l = 0 for l > N + 2k. Consequently the functions m;k do not have staggered support. In [CDV93] it is established that there exists a basis transformation U such that half stag m := U m
(31)
has staggered support. In section 4.5 we give a simple constructive algorithm for stag are orthonormalized to get the edge wavelets calculating U . The functions m;k edge
(32)
m
Orthonormalized edge wavelets
= B mstag :
. In the following we outline the orthonormalization procedure, i.e., the calculation of the matrix B . Let w be stag . Then from (30) and (31) it the matrix of inner products of the functions m;k follows that edge
m;k
w := mstag mstag t = UGhalf Ghalf t U t : On the other hand, we get from (32) and the orthonormality of the functions medge (33)
? ? w = B ?1 B ?1 t :
Thus the matrix B can be calculated from w by a Cholesky factorization and inversion. From (30{32) it follows that edge
m
= BUGhalf edge m?1 :
Thus the matrix Gedge for the re nement equation of the edge wavelets is given by
Gedge = BUGhalf :
(34)
To make the calculations complete we have to determine the matrix in (24). The matrix is independent of the scale m as one can see from the following argument. Since
p
half half m;k (2x) = 2 m?1;k (x)
it follows that (35)
half ; half = 2 half (2 ); half(2 ) = half ; half : m;k m;l m;k m;l m?1;k m?1;l
In view of (26) it follows that (36)
half t t = half m m = H H :
half Since k;l = half m;k ; m;l = k;l if k N ? 1 or l N ? 1, equation (36) above can be reduced to
(37)
t 0 = H0 ext 0 H0 ;
17
p
where H0 2 R (2N ?2)(4N ?4) with (H0 )k;l = hl?2k = 2 for 1 ? N k N ? 2 (2N ?2)(2N ?2) in and 1 ? N l 3N ? 4, and ext 0 extends the matrix 0 2 R the form 0 ext 0 = 00 I 2 R (4N ?4)(4N ?4) : Equation (37) is a non-homogeneous linear system in as many unknowns as equations. Whereas it gets numerically ill-conditioned, it can be solved symbolically. We have strong evidence that there exists a solution 0 of (37) but so far we have no proof: the existence of a solution is closely related to the eigenvalues of the matrix H1 2 R (2N ?2)(2N ?2) which is the restriction of H0 to the rst 2N ? 2 columns. If the absolute values of the eigenvalues of H1 are less than 1, then from (37) it follows that 0 =
1 X
n=0
? H1n H2 H2t H1t n ;
where H2 are the last 2N ? 2 columns of H0 . We mention a result from [Str96], which indicates that N eigenvalues of H1 are given by 2?k?1=2 ; (k = 0; : : : ; N ? 1): Since there is no estimate for the other N ? 1 eigenvalues available, this result does not give existence of a solution. However, p in all our considered examples the largest eigenvalue turned out to be 1= 2.
4.2 The Construction of Cohen, Daubechies and Vial The starting point is again a compactly supported orthogonal wavelet family on R . As in Meyer's approach the construction of Cohen, Daubechies and Vial [CDV93] retains the interior scaling functions and adds adapted edge scaling functions. In [CDJV93, CDV93] the family of transformed scaling functions restricted to R + is introduced as follows: (P ? N ?1?l half mod l N ?1?k m;l if 0 k N ? 1, m;k = half m;k if k N . The functions mod m;k generate all polynomials up to degree N ? 1 [CDV93, Proposition 4.1.]. In contrast to Meyer's construction this approach requires less edge scaling functions to ful ll this task. While in Meyer's construction the spaces Vmhalf are just the projections of Vm onto L2 (R + ), here the space
n
o
Vmhalf := span mod m;k ; k 2 N 0 = T (Vm ); where T = (Tk;l ) is a matrix with indices 0 k and 1 ? N l which satis es
Tk;l =
(? N ?1?l N ?1?k
k;l
if 0 k N ? 1, if k N . 18
Since Tk;l = 0 if l > k the family half mod m = T m
(38)
has staggered support. The spaces Vmhalf de ne a multiresolution analysis on L2 (R + ) and the corresponding wavelet spaces are given by
? Wmhalf := Vmhalf ? \ Vmhalf ?1 :
The functions mod m;k can be orthonormalized by a basis transformation mod edge m = A m : Again the orthonormalization matrix A is determined by the Cholesky decom-
(39)
position of (40)
mod t t ~ := mod m m = T T ;
where is as in (24), i.e.,
? ? ~ = A?1 A?1 t :
(41)
In the following we determine the lter matrix H edge ; once the lter matrix H edge is constructed, the re nement matrix Gedge of the edge wavelets can be calculated analogously to the construction presented in section 4.1. The lter matrix for the dilation equation satis es edge edge edge m = H m?1 :
(42) >From (38) and (26) we get (43)
half half mod m = T m = TH m?1 :
Suppose that there exists a dilation equation for mod m , i.e., (44)
mod mod mod m = H m?1 ;
then from (43) and (44) it follows that
TH = H mod T: Multiplication of this equation by a right inverse T y of T from the right gives (45)
H mod = THT y :
This yields the following condition on T , H and T y : (46)
THT y T = TH;
which is equivalent to (47)
N (T ) N (TH ); 19
where N denotes the nullspace. In particular this shows that the condition (46) is independent on the choice of the right inverse T y . >From (45) and (42) it follows that (48) H edge = A TH T yA?1 : Now the further procedure to construct Gedge is analogous as in section 4.1. For the reader's convenience we have summarized the calculation of the re nement matrices in section 4.5. This matrix analytical approach clearly reveals the similarity between the constructions proposed by Meyer and Cohen, Daubechies & Vial. In fact the only di erence in both constructions is that the construction of the lter matrix H edge incorporates the matrix T . Any right-invertible matrix T satisfying (47) can be used to construct wavelets on R + with di erent properties. The special form of the matrix T proposed in [CDV93] guarantees that the scaling functions have staggered support and that any polynomial up to degree N ? 1 can be represented as a linear combination of the scaling functions. Setting T = I gives the construction proposed by Meyer.
4.3 The Biorthogonal Case { The Constructions of Dahmen et al. In this section we show that the our matrix approach for the construction of wavelets on the interval can be generalized in a natural way to the construction of biorthogonal wavelets on the interval. This outlines the constructions proposed by Dahmen et al. [DKU99, DS98]. In the biorthogonal case one requires two scaling functions and ~ satisfying dilation equations
(x) =
N X k=1?N
hk (2x ? k)
and
~(x) =
N~ X k=1?N~
h~ k ~(2x ? k):
Both scaling functions satisfy (4), (6) and are biorthogonal, i.e., X ~ hk hk?2l = 20;l : k
The corresponding multiresolution analyses are given by Vm := span m;k ; k 2 Z ; V~m := span ~m;k ; k 2 Z : The wavelet spaces Wm and W~ m are then de ned by Wm = Vm?1 \ V~m? ; and W~ m = V~m?1 \ Vm? : For more background on biorthogonal wavelets we refer to [CDF92]. Following the notation of the previous sections we de ne the modi ed scaling functions on R + by half mod half (49) mod m = T m and ~m = T~ ~m ; 20
half where again half m and ~m are the restrictions to the positive real line. mod The two families mod m and ~m are biorthogonalized by two basis transforms ~ A and A, i.e., mod and ~edge := A~ ~mod (50) edge m := A m m m
satisfy edge t edge m ~m = I:
Analogously to (40), (41) the last equation is equivalent to (51)
?A?1?A~?1t = T T~t ;
half t where := half m ~m . For a given matrix T T~t the factorization into the matrices A?1 and A~?1 can be computed in several ways: one could for example use a factorization by means of a SVD, as suggested by Dahmen et al. [DKU99], an LU -decomposition, or simply set A = I and A~ = (T T~t )?1 . Each possible factorization results in di erent biorthogonal bases for the same multiresolution spaces Vmhalf and V~mhalf . For orthogonal wavelets we calculated the factorization by a Cholesky decomposition. The matrix can be calculated similarly to the orthogonal case (cf. (36)) as the solution of the following linear inhomogeneous system: = H H~ t :
The dilation matrices H edge and H~ edge are given by H edge = ATHT yA?1 and H~ edge = A~T~H~ T~y A~?1; where T y and T~y denote the right inverses of T and T~ satisfying (52) N (T ) N (TH ) and N (T~) N (T~H~ ): Note the similarity of the constructions of H edge in the orthogonal and biorthogonal case! The construction of the biorthogonal wavelet bases can be carried over from the orthogonal case. Since Vm Wm = Vm?1 and V~m W~ m = V~m?1 edge ~ m as we can write the projections of edge m?1 and ~m?1 onto Wm and W edge
edge
edge
half m := PWm m?1 = m?1 ? PVm m?1 ; edge edge ~mhalf := PW~ m ~edge m?1 = ~m?1 ? PV~m ~m?1 ;
and consequently half edge half m = G m?1 and
21
~mhalf = G~ half ~edge m?1 ;
where
Ghalf = I ? H~ edge t H edge and G~ half = I ? H edge t H~ edge :
In order to biorthogonalize the families of functions mhalf and ~mhalf we set edge half edge half (53) m := B m and ~m := B~ ~m ; where the matrices B and B~ satisfy
?B?1 ?B~ ?1 t = ; w
where
w := mhalf ~mhalf t = Ghalf G~ half t :
The construction presented above reveals that there is more freedom in generating biorthogonal wavelets on R + than for the the construction of orthogonal wavelets. The choice of the matrices T and T~ determines the properties of the multiresolution analyses. As in the orthogonal case, any T and T~, compatible with H and H~ in the sense of (52) can be used to construct biorthogonal wavelets on the interval. The choices of the biorthogonalizations (50) and (53) a ect the scaling functions and wavelets, but not the multiresolution and wavelet spaces. Dahmen et al. [DKU99] suggested transformations T and T~ for the construction of biorthogonal wavelet bases with certain polynomial exactness.
4.4 (Bi-)orthogonal Wavelets with Staggered Support The matrix analytical point of view of constructing wavelets on the half line clearly indicates how to impose additional properties on the wavelets and scaling functions. In the construction above we have not paid any attention to preserve staggered support of the scaling functions and wavelets. In the following we show how to construct (bi-)orthogonal wavelets and scaling functions with staggered support. To our knowledge biorthogonal wavelets on the interval with staggered support have not been considered in the literature so far. half The following lemma guarantees existence of a basis stag m of Vm with staggered support. Lemma 9. Let T0 2 R K (2N ?1) , K 2N ? 1 and let T0;1 be the K K submatrix consisting of the last K columns of T0 . If T0;1 is invertible, then there exists an invertible matrix S0 2 R K K such that S0 T0 is of lower triangular form.
? Proof. Since T0;1 is invertible, T0?;11 t exists and can be decomposed by an LU -factorization into ? P T0?;11 t = LU;
where L and U are lower and upper triangular matrices, respectively, and P is a permutation matrix. Thus
? Lt P t T0;1 = U ?1 t 22
? is a lower triangular matrix. (Note that U t is lower triangular, thus also U ?1 t , and Lt P t is invertible.) Let T0 = (T0;0 ; T0;1 ), then as a consequence ? Lt P t T0 = Lt P t (T0;0 ; T0;1 ) = Lt P t T0;0 ; (U ?1 )t is lower triangular. Thus the assertion is proved with S0 := Lt P t . T in (38) is of the form
T = T00 I0 :
Let S0 be de ned as in the lemma above, then T stag := S00 I0 T00 I0 = S00T0 I0 stag half is of lower triangular form and thus stag m := T m has staggered support. Proposition 4.3 in [CDV93] guarantees the existence of wavelets with staggered support. The above considerations can be easily carried over to biorthogonal wavelets. In order to get biorthogonal wavelets and scaling functions with staggered support an LU -factorization of (51) has to be performed, since only this factorization guarantees that the staggered support is preserved during the biorthogonalization procedure.
4.5 An Algorithm for the Calculation of the Re nement Matrices For the reader's convenience we summarize the computational steps for calculating the re nement matrices H edge and Gedge in the orthogonal case. The modi cations of this algorithm to calculate the re nement matrices in the biorthogonal case are obvious. Each step of the proposed algorithm can either be performed numerically or symbolically. Given the lter sequence hk of a compactly supported orthonormal wavelet family on R with hk = 0 if k ?N or k N +1 and the matrix T 2 R? K (2N?1) (T = I 2 R (2N ?1)(2N ?1) for Meyer's construction and (Tk;l ) = NN??kl??11 2 R N (2N ?1) for the construction of Cohen et al.). 1. De ne the lter matrix H := (Hk;l ) 2 R (2N ?1)(4N ?2) , with p Hk;l = hl?2k = 2 for 1 ? N k N ? 1 and 1 ? N l 3N ? 2. 2. Solve = H ext H t with 2 R (2N ?1)(2N ?1) and ext = 0 I0 2 R (4N ?2)(4N ?2) : 23
3. Compute the matrix of inner products ~ = T T t 2 R K K : 4. Compute A 2 R K K from the Cholesky decomposition ? ? ~ = A?1 A?1 t : 5. The dilation matrix for the edge scaling functions is then given by ? ? H edge = A TH T ext y Aext ?1 2 R K (K +2N ?1) ; where A 0 ext A = 0 I 2 R (K +2N ?1)(K +2N ?1) : and T ext = T0 I0 2 R (K +2N ?1)4N ?2 ; and (T ext )y is a right inverse of T ext . 6. Compute
C = I ? H edge t H edge 2 R (K +2N ?1)(K +2N ?1) ;
and de ne Ghalf 2 R N (K +2N ?1) as Ghalf = (Ck;l ) 0kN ?1
N ?K l2N ?2 7. Compute an upper triangular matrix U 2 R (N ?1)(N ?1) such that UGhalf is a lower triangular block matrix in the sense that (UGhalf )k;l = 0 for
l > N + 2k. This can be done using the following algorithm: (a) De ne the matrix C~ 2 R N N by C~k;l = Ghalf 0 k; l N ? 1: k;N +2l (b) Compute U by the unpivoted LU -decomposition C~ ?1 = LU:
8. Compute the matrix of inner products w = UGhalf Ghalf t U t 2 R N N : 9. Compute the Cholesky decomposition ? ? w = B ?1 B ?1 t : 10. The lter matrix Gedge is then given by Gedge = BUGhalf :
p
edge and Gedge for k N are given by H edge = hl?2k = 2 and The entries of H k;l p Gedge = g = 2, respectively. l?2k k;l
24
5 Results This section presents some closed form representations of the lter coecients for the Daubechies wavelets (section 5.1) and of the re nement matrices for the construction proposed by Cohen, Daubechies & Vial where N = 2 (section 5.3).
5.1 Mathematica Program and Output The following Mathematica program calculates the lter coecients of the Daubechies wavelets. (* Mathematica program *) Polys1[N_]:={Sum[Q[N,j]/2^j,{j,0, N-1}] - 1/2^(2N-2)} Polys2[N_]:=Table[Sum[Binomial[4N-i-j-2,2N+2l-i-1]* Q[N,i]*Q[N,j],{i,0,N-1},{j,0,N-1}], {l,1,N-1}] AllPolys[N_]:=Join[Polys1[N],Polys2[N]] Eqns[N_]:=Map[#==0&,AllPolys[N]] Unknowns[N_]:=Table[Q[N,j],{j,0,N-1}] CFSols[N_]:=Solve[Eqns[N],Unknowns[N]] GB[N_]:=GroebnerBasis[AllPolys[N],Reverse[Unknowns[N]]] cc[N_,k_]:=Binomial[2N-1,N-k] * Sum[Q[N,j] Binomial[N-k,j]/Binomial[2N-1,j], {j,0,N-1}] CoefficientTable[N_,rules_]:= Table[h[N,k]->Simplify[ cc[N,k] /.rules],{k,1-N,N}] (* Groebner basis in the case N=3 *) GB[3]
f9 ? 96 (3 0) ? 1536 (3 0)2 ? 4096 (3 0)3 + 16384 (3 0)4 21 (3 0) + 32 (3 0)2 ? 128 (3 0)3 + 3 (3 1) ? 3 ? 120 (3 0) ? 256 (3 0)2 + 1024 (3 0)3 + 12 (3 2)g Q
;
Q
Q
Q
;
;
;
Q
Q
Q
;
;
Q
Q
;
Q
;
;
;
;
Q
;
Q
;
;
(* A real solution for N=3 *) rules = Simplify[CFSols[3][[3]]] CoefficientTable[3, rules]
p p p p p p 3 5 + 2 10 f (3 ?2) ! 1 + 10 +16 5 + 2 10 (3 ?1) ! 5 + 10 + 16 p p p p p p 5 ? 10 + 5 + 2 10 10 ? 5 + 2 10 5 ? (3 1) ! (3 0) ! 8 p p 8p p p p (3 2) ! 5 + 10 ? 3 5 + 2 10 (3 3) ! 1 + 10 ? 5 + 2 10 g h
;
;h
h
h
;
;
16
25
;
;
;h
;
;h
;
;
16
(* Groebner basis in the case N=4 *) GB[4]
f625 + 16000 (4 0) ? 1433600 (4 0)2 + 22937600 (4 0)3 + 220200960 (4 0)4 ? 4697620480 (4 0)5 ? 60129542144 (4 0)6 ? 137438953472 (4 0)7 + 1099511627776 (4 0)8 125 + 389200 (4 0) ? 1469440 (4 0)2 ? 29245440 (4 0)3 + 124780544 (4 0)4 Q
;
Q
Q
;
Q
Q
Q
;
Q
Q
;
;
;
Q
;
Q
;
Q
;
;
;
;
Q
;
+ 2936012800 (4 0)5 + 7516192768 (4 0)6 ? 51539607552 (4 0)7 + 39200 (4 1) ? 1875 ? 6661200 (4 0) + 57164800 (4 0)2 + 775864320 (4 0)3 ? 9064939520 (4 0)4 ? 136113553408 (4 0)5 ? 323196289024 (4 0)6 + 2456721293312 (4 0)7 + 196000 (4 2) ? 11625 + 3553200 (4 0) ? 42470400 (4 0)2 ? 483409920 (4 0)3 + 7817134080 (4 0)4 + 106753425408 (4 0)5 + 248034361344 (4 0)6 ? 1941325217792 (4 0)7 + 98000 (4 3)g Q
;
Q
Q
Q
Q
;
;
Q
;
Q
;
;
Q
;
Q
Q
;
Q
;
;
Q
Q
Q
;
;
Q
;
;
;
Q
;
;
Q
Q
;
Q
;
Q
;
;
;
5.2 Additional Explicit Results In this section, we provide explicit values for the univariate polynomials that describe both systems (4{6) and (17{18), as described in sections 3.1 and 3.4. Speci cally, we choose to give the univariate polynomial satis ed by hN , or x1 in the notation of section 3.1, for any solution of (4{6), which is the same polynomial as the univariate polynomial satis ed by Q0 for any solution of (17{18). For xed N , univariate polynomials in another hi or Qi (except for Q0 ) would be of the same degree and with integer coecients of same typical size. All polynomials are normalized by enforcing integer coecients and no nontrivial integer content. Let us denote p = c0 + + cd X d the common univariate polynomial for given N and assume that d is its degree. To reduce the size of the coecients so as to display the polynomials, we remark that all coecients remain integers under the substitution of X with X=22N ?3 . The substituted polynomials are called p0 = c00 + + c0d X d . The degrees of p and p0 together with the sizes of their coecients are given in the following table, where jcj denotes the number of digits of the non-zero integer c.
N 2 3 4 5 6
d jcd j jc0 j maxi jci j d jc0d j jc00 j maxi jc0i j
2 4 8 16 32
1 5 13 37 87
1 1 3 13 29
1 5 13 37 87
2 4 8 16 32
The explicit values of the polynomials are:
p02 = 2X 2 ? 2 ? 1; 26
1 1 1 3 1
1 1 3 13 29
1 2 4 13 29
p03 = 4X 4 ? 8X 3 ? 24X 2 ? 12 + 9; p04 = X 8 ? 4X 7 ? 56X 6 ? 140X 5 + 210X 4 + 700X 3 ? 1400X 2 + 500 + 625; p05 = 256X 16 ? 2048X 15 ? 122880X 14 ? 1162240X 13 + 3672320X 12 + 82199040X 11 ? 239052800X 10 ? 2639571200X 9 + 21067452000X 8 ? 46192496000X 7 ? 73209920000X 6 + 440535480000X 5 + 344423450000X 4 ? 1907594500000X 3 ? 3529470000000X 2 ? 1029428750000 + 2251875390625; p06 = X 32 ? 16X 31 ? 3968X 30 ? 127120X 29 + 908488X 28 + 99001616X 27 ? 206896256X 26 ? 45046412656X 25 + 514227272860X 24 + 9384914783664X 23 ? 326335992812928X 22 + 3719423566862640X 21 ? 4725849211541640X 20 ? 321029601376721328X 19 +2420305333571518848X 18 + 16398235495598877648X 17 ? 211208519547389641914X 16 ? 1033088836222729291824X 15 + 9606191868945358307712X 14 + 80272488735445037902416X 13 ? 74446118321296204796040X 12 ? 3691291866649887797453520X 11 ? 20403669167515311931757952X 10 ? 36966997633084650250167888X 9 + 127608470131412725062780060X 8 + 704247243896852021529183888X 7 ? 203778389811329721233161344X 6 ? 6143110504249885426575754992X 5 + 3551450686163073382755632328X 4 + 31306969279922401804098069360X 3 ? 61565775711706432446473031552X 2 + 15639693023538327597289520112 + 61581291280182164914327485441: Finally, we provide univariate polynomials that are satis ed by QN ?1 in any solution of (17{18), for 2 N 8. They are all of degree 2 and prove that QN ?1 is the square root of a rational number.
p002 = 4X 2 ? 3; p003 = 8X 2 ? 5; p004 = 64X 2 ? 35; p005 = 128X 2 ? 63; p006 = 512X 2 ? 231; p007 = 1024X 2 ? 429; p008 = 16384X 2 ? 6435:
5.3 Re nement Matrices for the Construction Proposed by Cohen, Daubechies & Vial The lter coecients hk of the Daubechies wavelets for the case N = 2 are given by
p
p
p
p
h?1 = 1 +4 3 ; h0 = 3 +4 3 ; h1 = 3 ?4 3 ; h2 = 1 ?4 3 : The entries of the re nement matrices H edge and Gedge are the following: p p H edge = 2 (1137?119 3 ) 0;0
H0edge ;1 =
2182 p p p 2 (3969?2184 3 )(242883?140092 3 )(16589+9619 3 ) 14283372
q
27
q
p p 2 (3969?2184 3 )(123+85 3 ) H0edge = ? ;2 13092
H1edge ;0 H1edge ;1 H1edge ;2 H1edge ;4 H1edge ;4
q
p
p
p
2 3969?2184 3 )(242883?140092 3 )(32238331+17965009 3 ) = ( 501189240108 p2 999+238p3 ( ) = q 4364 p p 2 242883?140092 3 )(2963297+1744500 3 ) = ( 153128396 p3 q2 242883?140092p3 897+445p3 ( )( ) = 280712 q p p 2 242883?140092 3 )(897+445 3 ) =? ( 280712
Gedge 0;0 =
q
(2826138238+1021826769 p3)(3969?2184 p3)(6136686+2872499p3 )
12905163547593 p3 242883?140092p3 2102042+1389397p3 2826138238+1021826769 ( )( )( ) Gedge = ? 0;1 8603442365062 p p p 2826138238+1021826769 3 (80542?29121 3 ) Gedge = 0;2 7885831682
q
q
p
p
p
35089 (80542+29121 3 )(3969?2184 3 )(605486?295683 3 ) 8603442365062 q p3 3969?2184p3 1147827+503061p3 35089 80542+29121 ( )( )( ) Gedge 1;1 = 17206884730124 q p3 147657?115797p3 )( ) edge G1;2 = ? 35089 (80542+29121 15771663364 p3 q35089 80542+29121p3 ( ) edge G1;3 = ? 140356 q p3 ) edge G1;4 = 35089 (80542+29121 140356
Gedge 1;0 =
p
p
edge edge = h For k N the entries are given by Hk;l l?2k = 2 and Gk;l = gl?2k = 2.
Acknowledgment The work of F.C. and P.P. is partially supported by the SFB grant F1305 of the Austrian Science Foundation. The work of O.S. is partially supported by the SFB grant F1310 of the Austrian Science Foundation. The work of A.S. is partially supported by the Upper Austrian Government. The calculations of the total degree Grobner bases up to N = 34 have been performed on the machines of the UMS ME DICIS at E cole Polytechnique (Palaiseau, France). The authors thank Fabrizio Caruso, Philippe Dumas, and Zuhair Nashed for helpful comments and stimulating discussions.
28
References [BBH93] J.N. Bradley, C.M. Brislawn, and T. Hopper. The FBI wavelet/scalar quantisation standard for gray-scale ngerprint image compression. Technical report, Los Alamos National Laboratory, 1993. [Buc65] B. Buchberger. Ein Algorithmus zum Aunden der Basiselemente des Restklassenringes nach einem nulldimensionalen Polynomideal. PhD thesis, Philosophische Fakultat an der LeopoldFranzens-Universitat, Innsbruck, Austria, 1965. [Buc70] B. Buchberger. Ein algorithmisches Kriterium fur die Losbarkeit eines algebraischen Gleichungssystems. Aequationes Mathematicae, 4(3):271{272 and 374{383, 1970. English translation in [BW98]. [BW98] B. Buchberger and F. Winkler. Grobner bases and applications, volume 251 of London Math. Soc. Lecture Notes Series. Cambridge University Press, Cambridge, 1998. [CD93] A. Cohen and I. Daubechies. Orthonormal bases of compactly supported wavelets. III: Better frequency resolution. SIAM J. Math. Anal., 24(2):520{527, 1993. [CDF92] A. Cohen, I. Daubechies, and J.-C. Feauveau. Biorthogonal bases of compactly supported wavelets. Commun. Pure Appl. Math., 45(5):485{560, 1992. [CDJV93] A. Cohen, I. Daubechies, B. Jawerth, and P. Vial. Multiresolution analysis, wavelets and fast algorithms on an interval. C. R. Acad. Sci. Paris, Serie I 316(5):417{421, 1993. [CDV93] A. Cohen, I. Daubechies, and P. Vial. Wavelets on the interval and fast wavelet transforms. Appl. Comput. Harmon. Anal., 1(1):54{81, 1993. [CLO98] D. Cox, J. Little, and D. O'Shea. Using algebraic geometry. SpringerVerlag, New York, 1998. [Dau88] I. Daubechies. Orthonormal bases of compactly supported wavelets. Commun. Pure Appl. Math., 41(7):901{996, 1988. [Dau92] I. Daubechies. Ten Lectures on Wavelets, volume 61 of CBMS-NSF Regional Conference Series in Applied Mathematics. SIAM, Philadelphia, 1992. [Dau93] I. Daubechies. Orthonormal bases of compactly supported wavelets. II: Variations on a theme. SIAM J. Math. Anal., 24(2):499{519, 1993. [DJL92] R.A. DeVore, B. Jawerth, and B.J. Lucier. Image compression through wavelet transform coding. IEEE Trans. Inf. Theory, 38(2/II):719{746, 1992. 29
[DKU99] W. Dahmen, A. Kunoth, and K. Urban. Biorthogonal spline wavelets on the interval | stability and momentum conditions. Appl. Comput. Harmon. Anal., 6(2):132{196, 1999. [DS98] W. Dahmen and R. Schneider. Wavelets with complementary boundary conditions | functions spaces on the cube. Results in Mathematics, 34:255{293, 1998. [Fau94] J.C. Faugere. Resolution des systemes d'equations polynomiales. PhD thesis, Universite Paris VI | Pierre et Marie Curie, 1994. [GKP94] R.L. Graham, D.E. Knuth, and O. Patashnik. Concrete Mathematics. Addison-Wesley, 2nd edition, 1994. [Jaf92] S. Ja ard. Wavelet methods for fast resolution of elliptic problems. SIAM J. Numer. Anal., 29(4):965{986, 1992. [Kla97] A. Klappenecker. On algebraic properties of selfreciprocal polynomials and of Daubechies lters of low order. In Proc. of 1997 IEEE Int. Symp. on Inform. Theory, page 80, Ulm, Germany, 1997. [Mal89] S.G. Mallat. A theory for multiresolution signal decomposition: The wavelet representation. IEEE Trans. Pattern Anal. Mach. Intell., 11(7):674{693, 1989. [Mey91] Y. Meyer. Ondelettes sur l'intervalle. Rev. Mat. Iberoam., 7(2):115{ 133, 1991. [PS95] P. Paule and M. Schorn. A Mathematica version of Zeilberger's algorithm for proving binomial coecient identities. J. Symbolic Computation, 20:673{698, 1995. [PWZ96] M. Petkovsek, H.S. Wilf, and D. Zeilberger. A = B . A.K. Peters, 1996. [SSK98] O. Scherzer, A. Schoisswohl, and A. Kratochwil. Wavelet compression of 3D ultrasound data. Technical Report 7/1998, Industrial Mathematics Institute, Johannes Kepler University, Linz, 1998. [Str96] G. Strang. Eigenvalues of toeplitz matrices with 1 2 blocks. Z. Angew. Math. Mech., 76(2):37{39, 1996. [VzGG99] J. Von zur Gathen and J. Gerhard. Modern Computer Algebra. Cambridge University Press, Cambridge, 1999. [WA94] J.R. Williams and K. Amaratunga. Introduction to wavelets in engineering. Int. J. Numer. Methods Eng., 37(14):2365{2388, 1994. [Win96] F. Winkler. Polynomial Algorithms in Computer Algebra. RISC Texts and Monographs in Symbolic Computation. Springer, Vienna{ New York, 1996.
30
