From: Andreas Unterkircher Newsgroups: sci.math.symbolic Subject: Re: Recurrence Relations Date: Wed, 23 Jul 1997 10:54:46 +0200 Adam R Klivans wrote: > > Does anyone know if there is a symbolic comp program > that will solve recurrence relations which involve more than one variable. > > I.e., a(n,m) = a(n-1,m) + a(n,m-1) > > > > Any help is appreciated. > Thanks > Adam > > (PS-- The latest version of mathematica does not do it) If you mean by "solve" to compute an "analytical" expression for a(n,m) I cannot help you, but Groebner Bases techniques (available in Mathematica) help to solve systems of difference equations in the following way: Given a system of difference equations over $\N^n$ (IMPORTANT: your system above is not of this form: n=0, m=0: a(0,0) = a(-1,0) + a(0,-1), but so far Groebner Bases techniques can only be applied for polynomial rings and not for Laurent polynomial rings - I will come back to this point later), a solution is a map a: N^n -> K (K usually some field, I assume). Now you can identify a subset B of N^n such that all the values of a are uniquely determined by giving the values of a on B (this is sometimes called solving the Cauchy Problem). The values on B can be assigned arbitrarely, B can be an infinite set. Nevertheless the value of a on any point of N^n can be computed. I assume that you are familiar with Groebner Bases, if not you could have a look at the book "Ideals, Varietes and Algorithms" by Cox, Little, O'Shea (Springer), for example. Example: a(m+2,n) - a(m,n+1) = 0, a(m,n+2) - a(m,n) = 0 for all m,n >= 0 Now I want to know a(4,1). One could compute: m=2, n=1: a(4,1) - a(2,3) = 0, a(2,3) - a(2,1) = 0 thus a(4,1) = a(2,1) m=0, n=1: a(2,1) - a(0,1) = 0, a(0,3) - a(0,1) = 0 thus a(4,1) = a(0,1) = a(0,3) One sees that a(4,1) is determined by a(0,1), which may be assigned an arbitrary value. However this procedure is not practical, here is how to do it with Groebner Bases: The system of differnce equations has to be "translated" into polynomial form: a(m+2,n) - a(m,n+1) = 0 -> x^2y - y = 0, a(m,n+2) - a(m,n) = 0 -> y^2 - 1 = 0 Now using MMA: In[1]:= f = {x^2y - y, y^2 - 1}; In[2]:= gb = GroebnerBasis[f,{x,y}] Out[2]= {-1+y^2,-1+x^2} In[3]:= PolynomialReduce[x^4y,gb,{x,y}] Out[3]= {{0,y+x^2y},y} x^4y corresponds to a(4,1); y is the remainder of x^4y after division by gb and y corresponds to a(0,1) and so we see that a(4,1) = a(0,1). Furthermore: In[4]:= PolynomialReduce[y^3,gb,{x,y}] Out[4]= {{y,0},y} which means: a(0,3) = a(0,1), as we have already seen. The general procedure is as follows: 1. Translate the difference equations into polynomial form and obtain a set f of polynomials. 2. Compute a Groebner Bases gb of the ideal I = with respect to any term order. The polynomials of gb define the same difference equation as the ones in f. The equation is just rewritten in a "reduced" form. 3. Lt(I) is the set of all monomials which occur as leading term of some polynomial in I. The set {all monomials}\Lt(I) is the set were the values of a solution of the difference equation may be assigned arbitrarely. In the above example this is {1,x,y,xy} = {(0,0),(1,0),(0,1),(1,1)}. This set may be infinite. 4. If you want to know the solution a on a point p = (m_1,...,m_n) outside {all monomials}\Lt(I) you just have to compute remainder r of the monomial p = x_1^m_1...x_n^m_n after division by gb (this can be done with PolynomialReduce). The polynomial r has support in {all monomials}\Lt(I) (where you can assign the values as you want) and the value of the soulution on p is: a(p) = r. As I said this works only for difference equations over N^n, but it can be extendet to equations over Z^n as follows: Consider the algebra of Laurent polynomials K[x_1,...,x_n,x_1^-1,...,x_n^-1] as the factor algebra K[x_1,...,x_n,y_1,...,y_n]/ and compute the Groebner Bases of the inverse image of an ideal in the Laurent poly ring in K[x_1,...,x_n,x_1^-1,...,x_n^-1]. For more details to this approach consider: S. Zampieri: "A solution of the Cauchy Problem for multidimensional discrete linear shift-invariant systems", Linear Algebra and its Applications 202: 143-162,1994. Another approach is to try to extend Groebner Basis theory to the ring of Laurent Polynomial ring. There is a paper of F. Pauer and myself on this subject (available on request). For more information on Groebner Bases and difference equations you could have a look at: U. Oberst: "Multidimensional constant linear systems", Acta Applicandae Mathematicae, 20: 1-175,1990. Hope this helps, Andreas -- Andreas Unterkircher Institut fuer Umformtechnik ETH Zentrum / CLA, Tannenstr. 3,CH-8092 Zuerich, Switzerland Phone: ++43 1 632 6621 www.ifu.ethz.ch ============================================================================== From: Andreas Unterkircher Newsgroups: sci.math.symbolic Subject: Re: Recurrence Relations Date: Thu, 24 Jul 1997 11:21:55 +0200 Paul Abbott wrote: > > Andreas Unterkircher wrote: > > > The system of differnce equations has to be "translated" into polynomial > > form: > > > > a(m+2,n) - a(m,n+1) = 0 -> x^2y - y = 0, > > Should this be x^2 - y = 0? Thank you very much for the pointing out this mistake ! Sorry, the correct lines are: a(m+2,n) - a(m,n+1) = 0 -> x^2 - y = 0, a(m,n+2) - a(m,n) = 0 -> y^2 - 1 = 0 Now using MMA: In[1]:= f = {x^2 - y, y^2 - 1}; In[2]:= gb = GroebnerBasis[f,{x,y}] Out[2]= {-1+y^2,x^2-y} In[3]:= PolynomialReduce[x^4y,gb,{x,y}] Out[3]= {{x^2,1+x^2y},y} In[4]:= PolynomialReduce[y^3,gb,{x,y}] Out[4]= {{y,0},y} Luckily (or unfortunately) the result is the same (thus I did not recognize my initial error). > > As I said this works only for difference equations over N^n, but it can > > be extendet to equations over Z^n as follows: > > > > Consider the algebra of Laurent polynomials > > > > K[x_1,...,x_n,x_1^-1,...,x_n^-1] > > > > as the factor algebra > > > > K[x_1,...,x_n,y_1,...,y_n]/ > > > > and compute the Groebner Bases of the inverse image of an ideal in the > > Laurent poly ring in K[x_1,...,x_n,x_1^-1,...,x_n^-1]. For more details > > to this approach consider: > > S. Zampieri: "A solution of the Cauchy Problem for multidimensional > > discrete linear shift-invariant systems", Linear Algebra and its > > Applications 202: 143-162,1994. > > This is most interesting to me. At present I have an Honours student > working on a project of Factorizing Wavelet Transforms which involves > Euclidean Factorization over the ring of Laurent polynomials. The > interesting feature is that there is no unique factorization. There is another paper with the same results as Zampieri, but more algebraic: E. Zerz and U. Oberst: "The canonical Chauchy Problem for Linear Systems of Partial Difference Equations with Constant Coefficients over the complete r-Dimensional Integral lattice \bf Z^r", Acta Applicandae mathematicae 31: 249-273,1993. Best regards, Andreas -- Andreas Unterkircher Institut fuer Umformtechnik ETH Zentrum / CLA, Tannenstr. 3,CH-8092 Zuerich, Switzerland Phone: ++43 1 632 6621 www.ifu.ethz.ch