**S. L E W A N O W I C Z (Wroc law)**

### A FAST ALGORITHM FOR THE CONSTRUCTION OF RECURRENCE RELATIONS FOR MODIFIED MOMENTS

### Abstract. A new approach is presented for constructing recurrence rela- tions for the modified moments of a function with respect to the Gegenbauer polynomials.

### 1. Introduction. Let w be a weight function on the interval (−1, 1).

### We call the integrals

### (1.1) m k [w] ≡ m ^{λ} _{k} [w] :=

### 1

## R

### −1

### w(x)C _{k} ^{λ} (x) dx (λ > −1/2) modified moments of w with respect to the Gegenbauer polynomials

### C _{k} ^{λ} (x) := (−1) ^{n} (2λ) k

### 2 ^{k} k!(λ + ^{1} _{2} ) k (1 − x ^{2} ) ^{λ−1/2} **D** ^{k} {(1 − x ^{2} ) ^{k+λ−1/2} },

**D**

**where D = d/dx (see, e.g., [2, Vol. 1, §10.9], or [9, Vol. 1, §8.3]), or the** Gegenbauer moments, for short. Here (a) m is the Pochhammer symbol given by

**where D = d/dx (see, e.g., [2, Vol. 1, §10.9], or [9, Vol. 1, §8.3]), or the**

### (a) 0 := 1, (a) m := a(a + 1) . . . (a + m − 1) (m = 1, 2, . . .).

### Modified moments, provided they are accurately computable, are used in the generation of nonstandard orthogonal polynomials (see [3, 4, 5] and the references given therein) which have applications in many areas (e.g.

### numerical quadrature, summation of series, approximation). The Chebyshev moments,

### (1.2) τ k [w] :=

### 1

## R

### −1

### w(x)T k (x) dx = m ^{0} _{0} [w] (k = 0),

### 1

### 2 km ^{0} _{k} [w] (k ≥ 1)

*1991 Mathematics Subject Classification: Primary 33A65; Secondary 65D30, 65Q05.*

*Key words and phrases: modified moments, Gegenbauer polynomials, recurrence re-* lations.

### Supported by Komitet Bada´ n Naukowych (Poland) under grant 2 1029 91 01.

### are needed in the numerical evaluation of certain difficult integrals by the so- called modified Clenshaw–Curtis method (see, e.g., [10, 11, 12]). Sometimes, by a stroke of luck, modified moments are explicitly known. More frequently, however, they are computed from a recurrence relation of the form

### (1.3) Lm k [w] = %(k),

### judiciously employed. Here L is a difference operator, L :=

### u+r

### X

### j=u

### λ j (k)E ^{j} ,

### where λ u , λ u+1 , . . . , λ u+r (λ u 6≡ 0, λ _{u+r} 6≡ 0) are known (rational) functions in k, u ∈ Z, r ∈ Z + is referred to as the order of L, and E ^{j} (j ∈ Z) is the jth shift operator, acting on the variable k:

### E ^{j} µ(k) = µ(k + j) for any function µ : Z → Z. (We write I for E ^{0} .)

### In the cited references, such recurrence relations are constructed by ad hoc methods. A systematic way of constructing a recurrence (1.3) under the assumption that the function w satisfies a linear differential equation (1.4)

### n

### X

### i=0

### p ni **(x)D** ^{i} w = q,

**(x)D**

### where p ni are polynomials, and q is a known function, was described in [6].

### The first algorithm given therein (which, however, is the most complex) leads to a recurrence of the lowest possible order. This algorithm seems to be of great theoretical value. For instance, it helped us to obtain some partial results on certain hypergeometric sums which were later generalized to the form given in [7]. Unfortunately, the degree of complexity of the algorithm grows quickly with n, so that the calculations may be very tedious.

### The aim of the present paper is to propose a new algorithm, which in the author’s belief is equivalent to the best algorithm of [6]. As a final result we obtain a pair P, L of difference operators and a function ψ(k) such that the recurrence (1.3) holds with %(k) = P m k [q] − ψ(k). The order of the recurrence equals

### ord(P ) + 2 max

### 0≤i≤n, p

ni### 6≡0 (deg p ni − i), where ord(P ) is the order of the operator P .

### In applications, the order n of the differential equation (1.4) is usually

### not greater than 3. We give the closed-form expressions for P, L and %(k)

### for the cases n = 1 and n = 2, and show how the case n = 3 may be treated

### with a little effort, using the results given for n ≤ 2.

### Let us remark that the algorithm can be easily extended to the case of arbitrary n > 3.

### It should be noticed that P, L and %(k) are given in terms of certain basic difference operators (see Section 2). A scalar form of the recurrence relation (1.3) may be obtained using a language for symbolic manipulation, as for instance Maple [1].

### In Section 2, we give some important properties of the Gegenbauer mo- ments (1.1). Section 3 contains the main result of the paper—formulae for P, L and %(k). In Section 4 we give an illustrative example.

### 2. Basic identities

### Lemma 2.1 [6]. The Gegenbauer moments (1.1) satisfy the identities m k [xw(x)] = Xm k [w],

### (2.1)

### Dm k **[Dw] = m** k [w] + Dϕ k [w], (2.2)

**[Dw] = m**

### where X and D are the second-order difference operators X := k + 2λ − 1

### 2k + 2λ E ^{−1} + k + 1 2k + 2λ E, (2.3)

### D := 1

### 2k + 2λ {E ^{−1} − E}, (2.4)

### and

### (2.5) ϕ k [w] ≡ ϕ ^{λ} _{k} [w] := [w(x)C _{k} ^{λ} (x)] ^{x=1} _{x=−1} . It is easy to generalize equation (2.1) to the form (2.6) m k [pw] = p(X)m k [w] (p a polynomial).

### In Lemmas 2.2–2.4 we give identities which may be considered as gener- alizations of (2.2). We shall need some notation.

### For i = 0, 1, . . . and σ = ±1 define the difference operator (2.7) A ^{(σ)} _{i} := I − σα i (k)E,

### where

### α i (k) := (2k + 2λ + 1) 2

### (2k + 2λ + i + 1) 2

### . Notice that

### A ^{(σ)} _{0} = I − σE, A ^{(σ)} _{1} = I − σ 2k + 2λ + 1 2k + 2λ + 3 E.

### Further, let

### S _{ij} ^{(σ)} := A ^{(σ)} _{i} A ^{(σ)} _{i−1} · · · A ^{(σ)} _{j} , (2.8)

### P _{j} ^{(σ)} := S _{j−1,0} ^{(σ)} (j = 0, 1, . . .).

### (2.9)

### We adopt the convention that S _{ij} ^{(σ)} = I for i < j.

### Also, we will use the notation

### κ(k) := (k + 1)(k + 2λ − 1), (2.10)

### µ i (k) := 2 ^{−b(i+1)/2c} (2k + 2λ + 1) i (i = 0, 1, . . .).

### (2.11)

### Finally, let us introduce the differential operators (2.12)

### (2.13)

**U := (x** ^{2} **− 1)D + (3 − 2λ)xI,** **G := U D,**

**U := (x**

**− 1)D + (3 − 2λ)xI,**

**G := U D,**

**V** σ **:= (x + σ)D + (3/2 − λ)I,** **H** σ **:= V** σ **D** (σ = ±1).

**V**

**:= (x + σ)D + (3/2 − λ)I,**

**H**

**:= V**

**D**

**Here I is the identity operator.**

**Here I is the identity operator.**

### Now we are able to prove the following.

**Lemma 2.2. Let Q** 1 be any of the following first-order differential oper- ators:

**Lemma 2.2. Let Q**

### (2.14) **Q** 1 :=

**Q**

###

###

###

**U** (case 1A),

**U**

**V** σ (σ = ±1) (case 1B),

**V**

**D** (case 1C).

**D**

### Then the identity

### (∗) Q 1 m k **[Q** 1 w] = M 1 m k [w] + τ _{k} ^{(1)} [w]

**[Q**

### holds with Q 1 :=

###

###

###

### I (case 1A), P _{1} ^{(σ)} (case 1B), D (case 1C),

### M 1 :=

###

###

###

### κ(k)D (case 1A), µ 1 (k)P _{1} ^{(−σ)} (case 1B),

### I (case 1C),

### τ _{k} ^{(1)} [w] :=

###

###

###

### 0 (case 1A),

### P _{1} ^{(σ)} ϕ k [(x + σ)w] (case 1B), Dϕ k [w] (case 1C).

### P r o o f. C a s e 1A. As m k [(x ^{2} **− 1)Df ] = Hm** k [f ], where H := (k + 2λ

**− 1)Df ] = Hm**

### − 2) _{2} E ^{−1} − (k + 1) _{2} E (cf. [6, Eq. (20)]), we have

### m k **[Q** 1 w] = m k **[U w] = {H + (3 − 2λ)X}m** k [w] = (k + 1)(k + 2λ − 1)Dm k [w].

**[Q**

**[U w] = {H + (3 − 2λ)X}m**

### Thus identity (∗) holds with Q 1 , M 1 and τ _{k} ^{(1)} [w] given in the lemma.

### C a s e 1B. It can be checked that the operator Q := (k + 2λ − 1)I + σ(k + 2)E satisfies the following equations:

### P _{1} ^{(σ)} (X + σI) = QD,

### Q + ( ^{3} _{2} − λ)P _{1} ^{(σ)} = ^{1} _{2} (2k + 2λ + 1)P _{1} ^{(−σ)} .

### Using these identities and (2.2), we obtain

### P _{1} ^{(σ)} m k **[V** σ w] = QDm k **[Dw] + (** ^{3} _{2} − λ)P _{1} ^{(σ)} m k [w]

**[V**

**[Dw] + (**

### = ^{1} _{2} (2k + 2λ + 1)P _{1} ^{(−σ)} m k [w] + QDϕ k [w]

### = µ 1 (k)P _{1} ^{(−σ)} m k [w] + P _{1} ^{(σ)} ϕ k [(x + σ)w].

### Thus, equation (∗) holds with Q 1 , M 1 and τ _{k} ^{(1)} [w] specified in the lemma.

### C a s e 1C. Equation (∗) is a disguised form of (2.2).

### The next two lemmas can be proved in a similar manner.

**Lemma 2.3. Let Q** 2 be any of the following second-order differential op- erators:

**Lemma 2.3. Let Q**

### (2.15) **Q** 2 :=

**Q**

###

###

###

**G** (case 2A),

**G**

**H** σ (σ = ±1) (case 2B),

**H**

**D** ^{2} (case 2C).

**D**

### Then the identity

### Q 2 m k **[Q** 2 w] = M 2 m k [w] + τ _{k} ^{(2)} [w]

**[Q**

### holds with Q 2 :=

###

###

###

### I (case 2A), P _{2} ^{(σ)} (case 2B), D ^{2} (case 2C),

### M 2 :=

###

###

###

### κ(k)I (case 2A), µ 2 (k)E (case 2B),

### I (case 2C),

### τ _{k} ^{(2)} [w] :=

###

###

###

### 0 (case 2A),

### µ 2 (k)EDϕ k [w] + P _{2} ^{(σ)} ϕ k **[(x + σ)Dw]** (case 2B), Dϕ k [w] + D ^{2} ϕ k **[Dw]** (case 2C).

**[(x + σ)Dw]**

**[Dw]**

**Lemma 2.4. Let Q** 3 be any of the following third-order differential oper- ators:

**Lemma 2.4. Let Q**

### (2.16) **Q** 3 :=

**Q**

###

###

###

###

###

###

###

###

###

**GU** (case 3A),

**GU**

**V** σ **G** (σ = ±1) (case 3B),

**V**

**G**

**DG** (case 3C),

**DG**

**H** σ **V** σ (σ = ±1) (case 3D), **DH** σ (σ = ±1) (case 3E),

**H**

**V**

**DH**

**D** ^{3} (case 3F).

**D**

### Then the identity

### Q 3 m k **[Q** 3 w] = M 3 m k [w] + τ _{k} ^{(3)} [w]

**[Q**

### holds with

### Q 3 :=

###

###

###

###

###

###

###

###

###

### I (case 3A), P _{1} ^{(σ)} (case 3B), D (case 3C), P _{3} ^{(σ)} (case 3D), P _{2} ^{(σ)} D (case 3E), D ^{3} (case 3F).

### M 3 :=

###

###

###

###

###

###

###

###

###

### [κ(k)] ^{2} D (case 3A), µ 1 (k)P _{1} ^{(−σ)} κ(k)I (case 3B),

### κ(k)I (case 3C),

### µ 3 (k)EP _{1} ^{(−σ)} (case 3D), µ 2 (k)E (case 3E),

### I (case 3F).

### τ _{k} ^{(3)} [w] :=

###

###

###

###

###

###

###

###

###

###

###

###

###

### 0 (case 3A),

### P _{1} ^{(σ)} ϕ k **[(x + σ)Gw]** (case 3B),

**[(x + σ)Gw]**

### Dϕ k **[Gw]** (case 3C),

**[Gw]**

### µ 2 (k)EP _{1} ^{(σ)} {ϕ _{k} [(x + σ)w] + Dϕ k **[V** σ w]}

**[V**

### + P _{3} ^{(σ)} ϕ k **[(x + σ)DV** σ w] (case 3D), µ 2 (k)EDϕ k [w]

**[(x + σ)DV**

### +P _{2} ^{(σ)} {ϕ _{k} **[(x + σ)Dw] + Dϕ** k **[H** σ w]} (case 3E), Dϕ k [w] + D ^{2} ϕ k **[Dw] + D** ^{3} ϕ k **[D** ^{2} w] (case 3F).

**[(x + σ)Dw] + Dϕ**

**[H**

**[Dw] + D**

**[D**

### Observe that the difference operators Q 1 , Q 2 , Q 3 given in the above lemmas are always of the form

### P _{d} ^{(σ)} D ^{e}

### with d, e ≥ 0 and σ ∈ {−1, +1}. In Section 3.0, we shall need an operator which is a common multiple of two operators of the above form. Such an operator is given in Lemma 2.6 below. We must introduce a new family of difference operators first.

### For m = 0, 1, . . . and σ ∈ {−1, 1} define the difference operator R ^{(σ)} _{m} := (2k + 2λ) ^{−1} E ^{−1} + σ% m (k)I,

### where

### % m (k) := 2k + 2λ + 2m + 1 (2k + 2λ + m) 2

### . Further, let

### T _{ij} ^{(σ)} := R ^{(σ)} _{i} R _{i+1} ^{(σ)} . . . R _{j} ^{(σ)} , U _{h} ^{(σ)} := T _{0,h−1} ^{(σ)} (h = 0, 1, . . .).

### By convention, T _{ij} ^{(σ)} = I for i > j.

### Lemma 2.5. The identity

### (2.17) P _{v} ^{(σ)} D ^{r} = T _{v,v+r−1} ^{(σ)} P _{v+r} ^{(σ)} holds for v, r = 0, 1, . . .

### P r o o f. It is easy to verify that

### R ^{(σ)} _{0} A ^{(σ)} _{0} = D, R ^{(σ)} _{m} A ^{(σ)} _{m} = A ^{(σ)} _{m−1} R ^{(σ)} _{m−1} (m = 1, 2, . . .).

### Hence,

### D ^{r} = R ^{(σ)} _{0} R ^{(σ)} _{1} . . . R ^{(σ)} _{r−1} A ^{(σ)} _{r−1} A ^{(σ)} _{r−2} . . . A ^{(σ)} _{0} = U _{r} ^{(σ)} P _{r} ^{(σ)} (r ≥ 0), P _{v} ^{(σ)} U _{r} ^{(σ)} = A ^{(σ)} v − 1 . . . A ^{(σ)} _{1} A ^{(σ)} _{0} R ^{(σ)} _{0} R ^{(σ)} _{1} . . . R ^{(σ)} _{r−1}

### = R ^{(σ)} _{v} R ^{(σ)} _{v+1} . . . R ^{(σ)} _{v+r−1} A ^{(σ)} _{v+r−1} . . . A ^{(σ)} _{r+1} A ^{(σ)} _{r}

### = T _{v,v+r−1} ^{(σ)} S _{v+r−1,r} ^{(σ)} (v, r ≥ 0)

### and the result follows in view of S _{v+r−1,r} ^{(σ)} P r ^{(σ)} = P _{v+r} ^{(σ)} (cf. (2.8), (2.9)).

### Lemma 2.6. Let

### (2.18) Q 1 := P _{v} ^{(σ)} D ^{r} , Q 2 := P _{u} ^{(τ )} D ^{s} ,

### where v, r, u, s ≥ 0, v + r ≥ u + s and σ, τ ∈ {−1, 1}. Set q := 0 when σ = τ , and q := u otherwise. Then the operator Q := P _{d} ^{(σ)} D ^{e} , where e := max{r, q + s} and d := v + r − e, is a common multiple of Q 1 and Q 2 , (2.19) Q = Y i Q i (i = 1, 2),

### where

### (2.20) Y 1 := T _{d,v−1} ^{(σ)} , Y 2 :=

### ( S _{d−1,h} ^{(σ)} T _{h,u−1} ^{(σ)} (σ = τ ), P _{d} ^{(σ)} D ^{−h} U u ^{(τ )} (σ = −τ ), and where h := u + s − e.

### P r o o f. Let σ = τ . Substitute the expressions for Y i and Q i , given in (2.18) and (2.20), into the right-hand side of (2.19), and use Lemma 2.5.

### We obtain

### Y 1 Q 1 = T _{d,v−1} ^{(σ)} P _{v} ^{(σ)} D ^{r} = P _{d} ^{(σ)} D ^{v−d} D ^{r} = P _{d} ^{(σ)} D ^{e} = Q,

### Y 2 Q 2 = S _{d−1,h} ^{(σ)} T _{h,u−1} ^{(σ)} P _{u} ^{(σ)} D ^{s} = S _{d−1,h} ^{(σ)} P _{h} ^{(σ)} D ^{u−h} D ^{s} = P _{d} ^{(σ)} D ^{e} = Q.

### The case τ = −σ can be treated in an analogous way.

### 3. Formulae

**3.0. Introduction. Let P** n be a differential operator of the form

**3.0. Introduction. Let P**

### (3.1) **P** n =

**P**

### n

### X

### i=0

### p ni **(x)D** ^{i}

**(x)D**

### of order n ≤ 3, with polynomial coefficients p ni (i = 0, 1, . . . , n). Let the **operator P** n−1 of order n − 1 be defined by

**operator P**

**P** n−1 **w := P** n **w − Q** n (q n w),

**P**

**w := P**

**w − Q**

**where Q** n is an nth-order operator, of the form given in formula (2.{13+n}),

**where Q**

### and q n **is a polynomial such that Q** n (q n w) = p nn **D** ^{n} w + . . . and that the

**is a polynomial such that Q**

**D**

### corresponding difference operator Q n (given in Lemma 2.{1 + n}) has the least order. Repeating this process with n replaced by n − 1, n − 2 etc. we obtain the representation

### (3.2) **P** n w =

**P**

### n

### X

### i=0

**Q** i (q i w),

**Q**

**where we set Q** 0 **= I (identity operator), for convenience.**

**where we set Q**

**= I (identity operator), for convenience.**

### Let w be a solution of the differential equation **P** n w = q (n ≤ 3),

**P**

**where the differential operator P** n is of the form given in (3.2). Using Lemmas 2.2–2.4 we obtain difference operators Q i , M i and functionals τ _{k} ^{(i)} [·]

**where the differential operator P**

### such that

### (3.3) Q i m k **[Q** i w] = M i m k [w] + τ _{k} ^{(i)} [w] (i = 1, . . . , n).

**[Q**

### Now, using Lemma 2.6, a common multiple of the operators Q i can be obtained in the form

### (3.4) P = P _{d} ^{(σ)} D ^{e} ,

### where σ ∈ {−1, +1} and d, e ≥ 0 are integers; let the difference operators Z i be such that

### (3.5) P = Z i Q i (i = 1, . . . , n).

### Multiplying both sides of the equation (3.3) on the left by Z i , and using Lem- ma 2.6, (3.5) and (2.6), we obtain the result summarized in the following theorem.

### Theorem 3.1. Let w be a solution of the differential equation **P** n w = q (n ≤ 3),

**P**

**where the differential operator P** n is of the form given in (3.2) and q is a known function, and suppose the moments m k [w ^{(i)} ] (i = 0, 1, . . . , n) and m k [q] exist. Then we have the recurrence relation

**where the differential operator P**

### (3.6) Lm k [w] = %(k)

### with

### L :=

### n

### X

### i=0

### Z i M i q i (X), (3.7)

### %(k) := P m k [q] −

### n

### X

### i=1

### Z i τ _{k} ^{(i)} [q i w].

### (3.8)

### The order of the recurrence (3.6) is ord(P ) + 2 max

### 0≤i≤n, p

ni### 6≡0 (deg p ni − i),

### where ord(P ) = d + 2e is the order of the difference operator P.

### The last part of the theorem follows from [6, Eq. (80)].

**Now, the form of the operators Q** i (i = 1, . . . , n) in the representa- tion (3.2) can be deduced from the coefficients p ni of the operator (3.1).

**Now, the form of the operators Q**

### Thus, we can actually obtain closed-form expressions for P , L and %, at least for small n. In the next subsections we give such formulae for n = 1 and n = 2, and describe the way the case n = 3 may be treated with a small effort.

### 3.1. First-order differential equation. Assume that w satisfies the equa- tion

**P** 1 w ≡ p 11 (x)w ^{0} (x) + p 10 (x)w(x) = q(x).

**P**

### C a s e 3.1.1. p 11 (±1) = 0:

### P := I, L := κ(k)Dq 1 (X) + q 0 (X), %(k) := m k [q], where

### q 1 (x) := p 11 (x)/(x ^{2} − 1), q 0 (x) := p 10 **(x) − U q** 1 (x).

**(x) − U q**

### C a s e 3.1.2. p 11 (σ) 6= 0, p 11 (−σ) = 0 for σ = −1 or σ = 1:

### P := P _{1} ^{(σ)} , L := µ 1 (k)P _{1} ^{(−σ)} q 1 (X) + P _{1} ^{(σ)} q 0 (X),

### %(k) := P {m k [q] − ϕ k [p 11 w]}, where

### q 1 (x) := p 11 (x)/(x + σ), q 0 (x) := p 10 **(x) − V** σ q 1 (x).

**(x) − V**

### C a s e 3.1.3. p 11 (±1) 6= 0:

### P := D, L := q 1 (X) + Dq 0 (X), %(k) := P {m k [q] − ϕ k [p 11 w]}, where

### q 1 := p 11 , q 0 := p 10 **− Dq** _{1} .

**− Dq**

### 3.2. Second-order differential equation. Assume that w satisfies the equation

**P** 2 w ≡ p 22 w ^{00} (x) + p 21 (x)w ^{0} (x) + p 20 (x)w(x) = q(x).

**P**

### Define

### (3.9) l σ := p 21 (σ) − ^{1} _{2} (3 − 2λ)p ^{0} _{22} (σ).

### C a s e 3.2.1. p 22 (±1) = 0 and l −1 = l 1 = 0:

### P := I, L := κ(k){q 2 + Dq 1 (X)} + q 0 (X), %(k) := m k [q],

### where

### q 2 (x) := p 22 (x)/(x ^{2} − 1),

### q 1 (x) := [p 21 (x) − (3 − 2λ)xq 2 (x)]/(x ^{2} **− 1) − 2Dq** _{2} (x), q 0 := p 20 **− Gq** _{2} **− U q** _{1} .

**− 1) − 2Dq**

**− Gq**

**− U q**

### C a s e 3.2.2. p 22 (±1) = 0, l −σ = 0, l σ 6= 0 for σ = −1 or σ = 1:

### P := P _{1} ^{(σ)} , L := P _{1} ^{(σ)} {κ(k)q _{2} (X) + q 0 (X)} + µ 1 (k)P _{1} ^{(−σ)} q 1 (X),

### %(k) := P {m k [q] − ϕ k [(x + σ)q 1 w]}, where

### q 2 (x) := p 22 (x)/(x ^{2} − 1),

### q 1 (x) := [p 21 (x) − (3 − 2λ)xq 2 **(x)]/(x + σ) − 2(x − σ)Dq** 2 (x), q 0 := p 20 **− V** _{σ} q 1 **− Gq** _{2} .

**(x)]/(x + σ) − 2(x − σ)Dq**

**− V**

**− Gq**

### C a s e 3.2.3. p 22 (±1) = 0, l −1 6= 0, l _{1} 6= 0:

### P := D, L := q 1 (X) + D{κ(k)q 2 (X) + q 0 (X)},

### %(k) := D{m k [q] − ϕ k [q 1 w]}, where

### q 2 (x) := p 22 (x)/(x ^{2} − 1),

### q 1 (x) := p 21 (x) − (3 − 2λ)xq 2 (x) − 2(x ^{2} **− 1)Dq** 2 (x), q 0 := p 20 **− Dq** _{1} **− Gq** _{2} .

**− 1)Dq**

**− Dq**

**− Gq**

### C a s e 3.2.4. p 22 (σ) 6= 0, p 22 (−σ) = 0, l −σ = 0 for σ = −1 or σ = 1:

### P := P _{2} ^{(σ)} , L := A ^{(σ)} _{1} {µ _{1} (k)P _{1} ^{(−σ)} q 1 (X) + P _{1} ^{(σ)} q 0 (X)} + µ 2 (k)Eq 2 (X),

### %(k) := P {m k [q] − ϕ k [(x + σ)((q 2 w) ^{0} + q 1 w)]} − µ 2 (k)EDϕ k [q 2 w], where

### q 2 (x) := p 22 (x)/(x + σ),

### q 1 (x) := [p 21 (x) − ^{1} _{2} (3 − 2λ)q 2 **(x)]/(x + σ) − 2Dq** 2 (x), q 0 := p 20 **− V** σ q 1 **− H** σ q 2 .

**(x)]/(x + σ) − 2Dq**

**− V**

**− H**

### C a s e 3.2.5. p 22 (σ) 6= 0, p 22 (−σ) = 0, l −σ 6= 0 for σ = −1 or σ = 1:

### P := P _{1} ^{(σ)} D, L := W q 2 (X) + P _{1} ^{(σ)} {q 1 (X) + Dq 0 (X)},

### %(k) := P {m k [q] − ϕ k [q 1 w + p 22 w ^{0} ]} − W Dϕ k [q 2 w], where W := R ^{(σ)} _{1} µ 2 (k)E = µ 1 (k − 1)I + σµ 1 (k + 1)E, and

### q 2 (x) := p 22 (x)/(x + σ),

### q 1 (x) := p 21 (x) − ^{1} _{2} (3 − 2λ)q 2 **(x) − 2(x + σ)Dq** 2 (x),

**(x) − 2(x + σ)Dq**

### q 0 := p 20 **− Dq** _{1} **− H** _{σ} q 2 .

**− Dq**

**− H**

### C a s e 3.2.6. p 22 (±1) 6= 0:

### P := D ^{2} , L := q 2 (X) + Dq 1 (X) + D ^{2} q 0 (X),

### %(k) := P {m k [q] − ϕ k **[D(q** 2 w) + q 1 w]} − Dϕ k [q 2 w], where

**[D(q**

### q 2 := p 22 , q 1 := p 21 **− 2Dq** _{2} , q 0 := p 20 **− Dq** _{1} **− D** ^{2} q 2 .

**− 2Dq**

**− Dq**

**− D**

### 3.3. Third-order differential equation. Assume that w satisfies the equa- tion

### (3.10) **P** 3 w ≡ p 33 w ^{000} (x) + p 32 w ^{00} (x) + p 31 (x)w ^{0} (x) + p 30 (x)w(x) = q(x).

**P**

### The list of explicit formulae for the difference operators P , L and the function % covers ten cases. Thus, it seems too long to be given here. How- ever, one may obtain a recurrence for the modified moments of w with a small effort, using the results of Section 3.2. To this end, represent the **left-hand side of (3.10) in the form P** 3 **w = Q** 3 (q 3 **w) + P** 2 **w, where Q** 3 is of the form (2.16) (cf. Section 3.0). Then apply (formally) the results of **Section 3.2 to the equation P** 2 w = q ^{∗} with q ^{∗} **:= q − Q** 3 (q 3 w), which will yield operators P , L and a function ψ such that

**left-hand side of (3.10) in the form P**

**w = Q**

**w) + P**

**w, where Q**

**Section 3.2 to the equation P**

**:= q − Q**

### (3.11) Lm k [w] = %(k)

### with %(k) = P m k [q ^{∗} ] − ψ(k). Next, using Lemma 2.4, obtain operators Q 3 , M 3 and a functional τ _{k} ^{(3)} [·] satisfying Q 3 m k **[Q** 3 w] = M 3 m k [w] + τ _{k} ^{(3)} [w].

**[Q**

### Using Lemma 2.6 obtain a common multiple P ^{∗} of P and Q 3 , i.e. find difference operators W and Z such that P ^{∗} = W P = ZQ 3 . Then multiply (3.11) from the left by W and observe that

### W P m k [q ^{∗} ] = ZQ 3 m k [q ^{∗} ] = P ^{∗} m k [q] − ZQ 3 m k **[Q** 3 (q 3 w)]

**[Q**

### = P ^{∗} m k [q] − Z{M 3 q 3 (X)m k [w] + τ _{k} ^{(3)} [q 3 w]}.

### Now, it is easy to see that the recurrence relation in question has the form L ^{∗} m k [w] = % ^{∗} (k)

### with

### L ^{∗} := W L + ZM 3 q 3 (X), % ^{∗} (k) := P ^{∗} m k [q] − W ψ(k) − Zτ _{k} ^{(3)} [q 3 w].

### 4. An example. Consider the numerical evaluation of the integral

### (4.1) I =

### 1

## R

### 0

### f (x)x ^{α} J p (ωx) dx,

### where J p is the Bessel function of the first kind and of order p, and where α > −p−1, ω > 0 are real numbers. We assume that f is a smooth function.

### If ω is large, the use of standard integration rules is not efficient in view of

### the highly oscillatory character of the integrand. Therefore, special methods should be used, as, e.g. the method given in [8] or the modified Clenshaw–

### Curtis method (see [10] where the case α = 0 is discussed). The latter method is based on the approximation of f by a polynomial

### p(x) =

### n

### X

### k=0

### 0 c k T _{k} ^{∗} (x) (0 ≤ x ≤ 1).

### Here the symbol P _{0}

### denotes the sum with the first term halved, and T _{k} ^{∗} is the kth shifted Chebyshev polynomial, T _{k} ^{∗} (x) = T k (2x − 1). Replacing f by p in (4.1), we obtain

### I ≈

### n

### X

### k=0 0 c k

### 1

## R

### 0

### x ^{α} J p (ωx)T _{k} ^{∗} (x) dx

### = 2 ^{−α−1}

### n

### X

### k=0 0 c k

### 1

## R

### −1

### (1 + x) ^{α} J p ( ^{1} _{2} ω(1 + x))T k (x) dx

### = 2 ^{−α−1}

### n

### X

### k=0

### 0 c k τ k [w], where

### w(x) = (1 + x) ^{α} J p (a(1 + x)) (a = ^{1} _{2} ω).

### We show that the Chebyshev moments τ k [w] obey a sixth-order recurrence relation.

### The second-order differential equation for the Bessel function J p implies (1 + x) ^{2} w ^{00} + (1 − 2α)(1 + x)w ^{0} + [a ^{2} (1 + x) ^{2} + α ^{2} − p ^{2} ]w = 0.

### We have

### l ε = −2(1 + α)(1 + ε) (ε = ±1)

### (cf. (3.9)). It is easy to see that we have here case 3.2.4 with σ = 1, P = P _{2} ^{(1)} , and

### L = P _{2} ^{(1)} q 0 (X) + µ 2 (k)E[q 1 D + q 2 (X)],

### %(k) = − P _{2} ^{(1)} ϕ k [(1 + x) ^{2} w ^{0} − (2α + ^{3} _{2} )(1 + x)w]

### − µ _{2} (k)EDϕ k [(1 + x)w], where

### q 2 (x) = 1 + x, q 1 (x) = −2α − ^{5} _{2} , q 0 (x) = a ^{2} (1 + x) ^{2} + (α + ^{3} _{2} ) ^{2} − p ^{2} . Using (2.3), (2.4) and (2.11), we get

### L = 1

### 4 a ^{2} k − 2

### k E ^{−2} + a ^{2} (k + 3)(k − 1)

### k(k + 2) E ^{−1}

### + 1 4

### 2[b + 2(k − 1)(k − 2α − 1)] − a ^{2} k(14k + 31) (k + 2)(2k + 3)

### I + 2(k + 1)

### 2k + 3

### 3(1 − 2α) − b + 2(k − 1)(k + 3) + a ^{2} 2k ^{2} + 4k + 3 k(k + 2)

### E

### + 1

### 4(2k + 3)

### 2(2k + 1)[b + 2(k + 3)(k + 2α + 3)]

### − a ^{2} (k + 2)(14k − 3) k

### E ^{2} + a ^{2} (k − 1)(k + 3)

### (k + 2)(2k + 3) E ^{3} + 1

### 4 a ^{2} (k + 4)(2k + 1) (k + 2)(2k + 3) E ^{4} , where

### b := 3a ^{2} + 2α ^{2} − 2p ^{2} . From (2.5) we obtain

### %(k) = − 2 ^{α+3}

### k(k + 2)(2k + 3) {[2(k + 1) ^{2} − 3α − 5]J p (2a) + 6aJ _{p} ^{0} (2a)}.

### Now, the Gegenbauer moments m ^{0} _{k} [w] obey the recurrence relation Lm ^{0} _{k} [w] = %(k).

### Substituting m ^{0} _{k} [w] = ^{2} _{k} τ k [w] (cf. (1.2)), replacing k by k−1 and multiplying the resulting equation by 2(k ^{2} − 1)(2k + 1), we obtain the desired recurrence relation

### 6

### X

### j=0

### C j (k)τ k−3+j [w] = 2 ^{α+4} {(5 + 3α − 2k ^{2} )J p (2a) − 6aJ _{p} ^{0} (2a)}, where

### C 0 (k) = C 6 (−k) = a ^{2} (k + 1)(2k + 1), C 1 (k) = C 5 (−k) = 4a ^{2} (k + 1)(k + 2), C 2 (k) = C 4 (−k)

### = 2(k + 1)(2k + 1)[b + 2(k − 2)(k − 2α − 2)] − a ^{2} (k − 1)(14k + 17), C 3 (k) = 8(k ^{2} − 1)[2(k ^{2} − 4) − b + 3(1 − 2α)] + 8a ^{2} (2k ^{2} + 1).

### Note that in [10] an eighth-order homogeneous recurrence relation is given for the special case α = 0.

**References**

### [1] *B. W. C h a r et al ., Maple V Language Reference Manual , Springer, New York,*

### 1991.

### [2] A. E r d ´ *e l y i (ed.), Higher Transcendental Functions, McGraw-Hill, New York, 1953.*

### [3] *W. G a u t s c h i, Orthogonal polynomials—Constructive theory and applications, J.*

### Comput. Appl. Math. 12&13 (1985), 61–76.

### [4] *—, On generating orthogonal polynomials, SIAM J. Sci. Statist. Comput. 3 (1982),* 289–317.

### [5] *—, On certain slowly convergent series occurring in plate contact problems, Math.*

### Comp. 57 (1991), 325–338.

### [6] *S. L e w a n o w i c z, Construction of a recurrence relation for modified moments, J.*

### Comput. Appl. Math. 5 (1979), 193–205.

### [7] *—, Recurrence relations for hypergeometric functions of unit argument , Math.*

### Comp. 45 (1985), 521–535; corr. ibid. 47 (1987), 853.

### [8] *—, Evaluation of Bessel function integrals with algebraic singularity , J. Comput.*

### Appl. Math. 37 (1991), 101–112.

### [9] *Y. L. L u k e, The Special Functions and their Approximations, Academic Press, New* York, 1969.

### [10] *R. P i e s s e n s and M. B r a n d e r s, Modified Clenshaw–Curtis method for the compu-* *tation of Bessel function integrals, BIT 23 (1983), 370–381.*

### [11] *—, —, On the computation of Fourier transforms of singular functions, J. Comput.*

### Appl. Math. 43 (1992), 159–169.

### [12] R. P i e s s e n s, E. d e D o n c k e r - K a p e n g a, C. W. ¨ U b e r h u b e r and D. K. K a - *h a n e r, QUADPACK. A Subroutine Package for Automatic Integration, Springer,* Berlin, 1983.

STANIS LAW LEWANOWICZ

INSTITUTE OF COMPUTER SCIENCE UNIVERSITY OF WROC LAW

UL. PRZESMYCKIEGO 20 51-151 WROC LAW, POLAND E-mail: SLE@II.UNI.WROC.PL