Analytical Techniques for a Numerical Solution of the Linear Volterra Integral Equation of the Second Kind

Abstract

In this work we use analytical tools—Schauder bases and Geometric Series theorem—in order to develop a new method for the numerical resolution of the linear Volterra integral equation of the second kind.

1. Introduction

Most mathematical models used in many problems of physics, biology, chemistry, engineering, and in other areas are based on integral equations like the linear Volterra integral equation of the second kind:

𝑥(𝑡)=𝑦₀(𝑡)+^𝑡_𝑎𝑘(𝑡,𝑠)𝑥(𝑠)𝑑𝑠(𝑡∈[𝑎,𝑏]), (1.1) where 𝑘∶[𝑎,𝑏]×[𝑎,𝑏]→ℝ and 𝑦₀∶[𝑎,𝑏]→ℝ are two known continuous functions and 𝑥∶[𝑎,𝑏]→ℝ is the unknown function to be determined.

Many authors have paid attention to the study of linear Volterra integral equation of the second kind from the viewpoint of their theoretical properties, numerical treatment, as well as its applications (see e.g., [1–4] and the references therein). Specifically, there are several numerical techniques to solve this equation such as the collocation method, finite element method, and spectral method.

In this paper a new technique for solving this linear Volterra integral equation is shown. The method is based on two classical analytical tools: the Geometric Series theorem and Schauder bases in a Banach space. Schauder bases in adequate Banach spaces have been used in other numerical methods for solving some integral, differential, or integrodifferential equations (see [5–10]), although in each problem the analytical techniques are quite different, as fixed point theorems, duality mapping in a Banach space, and generalized least-squares methods.

Among the main advantages that our method presents over the classical ones, as collocation or quadrature (see [11]), we can point out that it is not necessary to solve linear equations systems. In addition, the involved integrals in our method are immediate and therefore, we do not have to use any quadrature method to calculate them.

The paper is organized as follows: some basic facts and properties on the Volterra equation (1.1) and on Schauder bases in 𝐶([0,1]) and 𝐶([0,1]²) are reviewed in Section 2. In Section 3 we define the approximating functions and we study the error. Finally, in Section 4 two numerical examples taken from [12, 13] are given to illustrate the theoretical results.

2. Development of the Numerical Method: Classical Tools in Functional Analysis

In this section we show some analytical techniques and some related results, useful for us in order to give our numerical method.

Let 𝐶([𝑎,𝑏]) be the Banach space of all continuous and real-valued functions on [𝑎,𝑏], endowed with its usual sup-norm:

‖𝑓‖=sup_{𝑥∈[𝑎,𝑏]}|𝑓(𝑥)|,(𝑓∈𝐶([𝑎,𝑏])). (2.1) Let us write

𝐿(𝐶([𝑎,𝑏]))={𝑇∶𝐶([𝑎,𝑏])→𝐶([𝑎,𝑏])∶𝑇islinearandcontinuous}, (2.2) and let 𝑇∈𝐿(𝐶([𝑎,𝑏])) be the operator defined by

𝑇𝑥(𝑡)=^𝑡_𝑎𝑘(𝑡,𝑠)𝑥(𝑠)𝑑𝑠,(𝑎≤𝑡≤𝑏,𝑥∈𝐶([𝑎,𝑏])). (2.3) It can be shown by an induction argument and Fubini's theorem that

(𝑇^𝑛𝑥)(𝑡)=^𝑡_𝑎^𝑡1_𝑎^𝑡2_𝑎⋯^𝑡𝑛−1_𝑎𝑘(𝑡,𝑡₁)𝑘(𝑡₁,𝑡₂)⋯𝑘(𝑡_𝑛−1,𝑡_𝑛)𝑥(𝑡_𝑛)𝑑𝑡_𝑛𝑑𝑡_𝑛−1⋯𝑑𝑡₁ (2.4) for all 𝑛≥1,𝑥∈𝐶([𝑎,𝑏]) and 𝑡∈[𝑎,𝑏] (we adopt the convention 𝑡₀=𝑡). Hence,

|(𝑇^𝑛𝑥)(𝑡)|≤‖𝑘‖^𝑛‖𝑥‖^𝑡_𝑎^𝑡1_𝑎^𝑡2_𝑎⋯^𝑡𝑛−1_𝑎1𝑑𝑡_𝑛⋯𝑑𝑡₁=‖𝑘‖^𝑛‖𝑥‖(𝑡−𝑎)^𝑛𝑛!≤‖𝑘‖^𝑛‖𝑥‖(𝑏−𝑎)^𝑛,𝑛! (2.5) and therefore, 𝑇^𝑛≤‖𝑘‖^𝑛(𝑏−𝑎)^𝑛𝑛!, (2.6) where ‖𝑘‖=sup_{(𝑥,𝑦)∈[𝑎,𝑏]}²|𝑘(𝑥,𝑦)|. The convergence of the series of real numbers ∑_𝑛≥0(‖𝑘‖(𝑏−𝑎))^𝑛/𝑛! ensures that the series ∑_𝑛≥0𝑇^𝑛 is convergent in 𝐿(𝐶([𝑎,𝑏])).

This remark and the fact that (1.1) can be written equivalently as

(𝐼𝑑−𝑇)𝑥=𝑦₀ (2.7) lead us to consider the following result (see [14]).

Theorem 2.1 (geometric series theorem). Let 𝑋 be a Banach space and let 𝑇∶𝑋→𝑋 be a continuous and linear operator such that the series ∑_𝑛≥0𝑇^𝑛 converges. Then, 𝐼𝑑−𝑇 is a continuous, linear, and bijective operator and (𝐼𝑑−𝑇)⁻¹=∑_𝑛≥0𝑇^𝑛.

Then the unique solution 𝑥 of (2.7) is given by

𝑥=(𝐼𝑑−𝑇)⁻¹𝑦₀=_𝑛≥0𝑇^𝑛𝑦₀. (2.8) Thus, the sequence {𝑦_𝑛}_𝑛∈ℕ defined by 𝑦_𝑛=𝑦₀+𝑇𝑦_𝑛−1 (𝑛≥1) converges to the solution 𝑥 of (2.7).

By making use of an appropriate Schauder basis in the space 𝐶([𝑎,𝑏]²) we will replace each 𝑇𝑦_𝑗−1 (𝑗≥1) for a new function 𝑤_𝑗∈𝐶([𝑎,𝑏]), easier to calculate, and in such a way the error ‖𝑥−(𝑦₀+𝑤_𝑗)‖ is small enough.

Let us recall now that a sequence {𝑏_𝑛}_𝑛≥1 in a Banach space 𝐸 is said to be a Schauder basis provided that for all 𝑥∈𝐸 there exists a unique sequence of scalars {𝛽_𝑛}_𝑛≥1 in such a way that

𝑥=_𝑛≥1𝛽_𝑛𝑏_𝑛. (2.9) The associated sequence of (continuous and linear) biorthogonal functionals {𝑏^∗_𝑛}_𝑛≥1 in the topological dual 𝐸^∗ of 𝐸 is given by

𝑏^∗_𝑛(_𝑘≥1𝛽_𝑘𝑏_𝑘)=𝛽_𝑛, (2.10) and the sequence of (continuous and linear) projections {𝑃_𝑛}_𝑛≥1 is defined by the partial sums:

𝑃_𝑛(_𝑘≥1𝛽_𝑘𝑏_𝑘)=^𝑛_𝑘=1𝛽_𝑘𝑏_𝑘. (2.11) We now consider the usual Schauder basis for the space 𝐶([0,1]) (although the results in this paper work for equations stated over functions defined in [𝑎,𝑏], we shall just consider the Banach space 𝐶([0,1]) for the sake of simplicity): for a dense sequence of distinct points {𝑡_𝑖}_𝑖≥1 in [0,1], with 𝑡₁=0 and 𝑡₂=1, we define

𝑏₁(𝑡)∶=1,(𝑡∈[0,1]) (2.12) and for all 𝑛≥2 we stand 𝑏_𝑛 for the piecewise linear continuous function determined by the points {𝑡₁,…,𝑡_𝑛} satisfying

∀𝑘<𝑛,𝑏_𝑛(𝑡_𝑘𝑏)=0,_𝑛(𝑡_𝑛)=1. (2.13) It is easy to obtain the sequence of biorthogonal functionals (see [15]): if 𝑥∈𝐶([0,1]), then

𝑏^∗₁(𝑥)=𝑥(𝑡₁) (2.14) and for all 𝑛≥2,

𝑏^∗_𝑛(𝑥)=𝑥(𝑡_𝑛)−^𝑛−1_𝑘=1𝑏^∗_𝑘(𝑥)𝑏_𝑘(𝑡_𝑛). (2.15) In addition, the sequence of projections {𝑃_𝑛}_𝑛≥1 satisfies the following interpolation property: for all 𝑥∈𝐶([0,1]), for all 𝑛≥1, and for all 𝑖≤𝑛, we have that

𝑃_𝑛(𝑥)(𝑡_𝑖)=𝑥(𝑡_𝑖). (2.16)

In what follows, for a real number 𝑥, [𝑥] will denote its integer part and 𝜎=(𝜎₁,𝜎₂)∶ℕ→ℕ×ℕ the bijective mapping defined by

(√𝜎(𝑛)∶=𝑛,√𝑛),if[√𝑛]=√√𝑛,(𝑛−[𝑛]²,[√𝑛]+1),if0<𝑛−[√𝑛]²≤[√√𝑛],([𝑛]+1,𝑛−[√𝑛]²−[√𝑛]),if[√𝑛]<𝑛−[√𝑛]². (2.17) If {𝑏_𝑖}_𝑖≥1 and {𝑏̂_𝑗}_𝑗≥1 are Schauder bases for the space 𝐶([0,1]), then the sequence

𝐵_𝑛(𝑠,𝑡)∶=𝑏_𝑖(𝑠)𝑏̂_𝑗(𝑡),(𝑠,𝑡∈[0,1]) (2.18) with 𝜎(𝑛)=(𝑖,𝑗) is a Schauder basis for 𝐶([0,1]²) (the proof of this fact can be found in [15, 16]). Therefore, from now on, if {𝑡_𝑖∶𝑖≥1} is a dense subset of distinct points in [0,1], with 𝑡₁=0 and 𝑡₂=1, and {𝑏_𝑖}_𝑖≥1 is the associated usual Schauder basis, then we will write {𝐵_𝑛}_𝑛≥1 to denote the Schauder basis for 𝐶([0,1]²) obtained in this “natural” way. It is not difficult to check that this basis satisfies similar properties to the ones for the one-dimensional case: for all 𝑠,𝑡∈[0,1],

𝐵₁(𝑠,𝑡)=1 (2.19) and for 𝑛≥2

𝐵_𝑛(𝑡_𝑖,𝑡_𝑗)=1,if𝜎(𝑛)=(𝑖,𝑗),0,if𝜎⁻¹(𝑖,𝑗)<𝑛. (2.20) In particular, the sequence of biorthogonal functionals {𝐵^∗_𝑛}_𝑛≥1 can be easily obtained: if 𝑥∈𝐶([0,1]²), then

𝐵^∗₁(𝑥)=𝑥(𝑡₁,𝑡₁), (2.21) and for all 𝑛≥2, if 𝜎(𝑛)=(𝑖,𝑗), we have

𝐵^∗_𝑛(𝑥)=𝑥(𝑡_𝑖,𝑡_𝑗)−^𝑛−1_𝑘=1𝐵^∗_𝑘(𝑥)𝐵_𝑘(𝑡_𝑖,𝑡_𝑗). (2.22) As a consequence, the sequence of projections {𝑄_𝑛}_𝑛≥1 satisfies

𝑄_𝑛(𝑥)(𝑡_𝑖,𝑡_𝑗)=𝑥(𝑡_𝑖,𝑡_𝑗), (2.23) whenever 𝑛,𝑖,𝑗∈ℕ and 𝜎⁻¹(𝑖,𝑗)≤𝑛.

Under some weak condition, we can estimate the rate of the convergence of the sequence of projections in the bidimensional case. To this purpose, consider the dense subset {𝑡_𝑖}_𝑖≥1 of distinct points in [0,1] and let 𝑇_𝑛 be the set {𝑡₁,…,𝑡_𝑛} ordered in an increasing way for 𝑛≥2. Let Δ𝑇_𝑛 denotes the maximum distance between two consecutive points of 𝑇_𝑛.

The following result is derived easily from (2.23) and the Mean-Value theorem.

Proposition 2.2. Let 𝑥∈𝐶¹([0,1]²) and write 𝑀∶=max{‖𝜕𝑥𝜕𝑠‖,‖𝜕𝑥𝜕𝑡‖}. (2.24) If 𝑛≥2, then ‖𝑥−𝑄_𝑛2(𝑥)‖≤4𝑀Δ𝑇_𝑛. (2.25)

3. Numerical Study of the Linear Volterra Integral Equation of the Second Kind: Convergence and Error

We are now in a position to define the functions 𝑤_𝑗∈𝐶([0,1]) announced in the preceding section. Let 𝑚∈ℕ,𝑤₀=0 and, with the notation above, define inductively, for 𝑗∈{1,…,𝑚}, the functions

𝑤_𝑗(𝑡)=^𝑡₀𝑄_𝑛𝑗2(𝑘(𝑡,𝑠)(𝑦₀(𝑠)+𝑤_𝑗−1(𝑠)))𝑑𝑠, (3.1) where 𝑛_𝑗(for 𝑗∈{1,…,𝑚}) are natural numbers.

We obtain a first estimation of the error ‖𝑥−(𝑦₀+𝑤_𝑚)‖.

Proposition 3.1. Maintaining the notation, ‖𝑥−(𝑦₀+𝑤_𝑚)‖≤‖𝑇^𝑚+1‖‖(𝐼𝑑−𝑇)⁻¹‖‖𝑦₀‖+^𝑚−1_𝑛=0‖𝑇‖^{𝑚−1−𝑛}‖𝑇(𝑦₀+𝑤_𝑛)−𝑤_𝑛+1‖. (3.2)

Proof. The triangle inequality gives ‖𝑥−(𝑦₀+𝑤_𝑚)‖≤‖𝑥−𝑦_𝑚‖+‖𝑦_𝑚−(𝑦₀+𝑤_𝑚)‖. (3.3) For the first summand we have ‖𝑥−𝑦_𝑚‖=‖_𝑛≥0𝑇^𝑛𝑦₀−^𝑚_𝑛≥0𝑇^𝑛𝑦₀‖=‖_{𝑛≥𝑚+1}𝑇^𝑛𝑦₀‖≤‖𝑇^𝑚+1‖‖_𝑛≥0𝑇^𝑛𝑦₀‖=‖𝑇^𝑚+1‖‖(𝐼𝑑−𝑇)⁻¹𝑦₀‖≤‖𝑇^𝑚+1‖‖(𝐼𝑑−𝑇)⁻¹‖‖𝑦₀‖. (3.4) For the second one, by an induction argument we can show that ‖𝑦_𝑚−(𝑦₀+𝑤_𝑚)‖≤^𝑚−1_𝑛=0‖𝑇‖^{𝑚−1−𝑛}‖𝑇(𝑦₀+𝑤_𝑛)−𝑤_𝑛+1‖. (3.5) Indeed, for 𝑚=1, the result is clearly true. Suppose that it holds for 𝑝≤𝑚−1. Then, ‖𝑦_𝑚−(𝑦₀+𝑤_𝑚)‖=‖(𝑦₀+𝑇𝑦_𝑚−1)−(𝑦₀+𝑤_𝑚)‖=‖𝑇𝑦_𝑚−1−𝑤_𝑚‖≤‖𝑇𝑦_𝑚−1−𝑇(𝑦₀+𝑤_𝑚−1)‖+‖𝑇(𝑦₀+𝑤_𝑚−1)−𝑤_𝑚‖≤‖𝑇‖‖𝑦_𝑚−1−(𝑦₀+𝑤_𝑚−1)‖+‖𝑇(𝑦₀+𝑤_𝑚−1)−𝑤_𝑚‖≤‖𝑇‖(^𝑚−2_𝑛=0‖𝑇‖^{(𝑚−1)−1−𝑛}‖𝑇(𝑦₀+𝑤_𝑛)−𝑤_𝑛+1‖)+‖𝑇(𝑦₀+𝑤_𝑚−1)−𝑤_𝑚‖=^𝑚−1_𝑛=0‖𝑇‖^{𝑚−1−𝑛}‖𝑇(𝑦₀+𝑤_𝑛)−𝑤_𝑛+1‖. (3.6)

In order to control the sum in the right-hand term of the inequality stated in Proposition 3.1, let us assume that 𝑘,𝑦₀∈𝐶¹([0,1]²) and let us write

ℎ∶=max_{𝑗=1,…,𝑚}Δ𝑇_𝑛𝑗. (3.7) Then we derive from Proposition 2.2 that

|(𝑇(𝑦₀+𝑤_𝑗)−𝑤_𝑗+1)(𝑡)|=|^𝑡₀(𝑘(𝑡,𝑠)(𝑦₀(𝑠)+𝑤_𝑗(𝑠)))−𝑄_𝑛2𝑗+1(𝑘(𝑡,𝑠)(𝑦₀(𝑠)+𝑤_𝑗(𝑠)))𝑑𝑠|≤4𝑀_𝑗ℎ, (3.8) where

𝑀_𝑗=max{‖𝜕(𝑘(𝑦₀+𝑤_𝑗))𝜕𝑠‖,‖𝜕(𝑘(𝑦₀+𝑤_𝑗))𝜕𝑡‖}, (3.9) and hence, ‖𝑇(𝑦₀+𝑤_𝑗)−𝑤_𝑗+1‖≤4𝑀_𝑗ℎ. (3.10)

To arrive at the announced estimation we finally have the following.

Proposition 3.2. The sequences {‖𝑤_𝑗‖}_𝑗∈ℕ and {‖𝑤^′_𝑗‖}_𝑗∈ℕ are bounded and, as a consequence, the sequence {𝑀_𝑗}_𝑗∈ℕ is also bounded.

Proof. First we show that for all 𝑡∈[0,1], |𝑤_𝑗(𝑡)|≤‖𝑦₀‖(^𝑗_𝑛=1𝑡^𝑛𝑛!‖𝑘‖^𝑛), (3.11) and, as a consequence, ‖𝑤_𝑗‖≤‖𝑦₀‖(^𝑗_𝑛=1‖𝑘‖^𝑛𝑛!)≤‖𝑦₀‖𝑒^‖𝑘‖, (3.12) by using an inductive argument: since the Schauder basis {𝐵_𝑛}_𝑛≥1 is clearly monotone (norm-one projections), we have |𝑤₁(𝑡)|≤^𝑡₀‖𝑘‖‖𝑦₀‖𝑑𝑠≤‖𝑦₀‖‖𝑘‖𝑡. (3.13) Suppose that the result holds for 𝑝≤𝑗−1. Then |𝑤_𝑗(𝑡)|≤^𝑡₀‖𝑘‖‖𝑦₀‖𝑑𝑠+^𝑡₀‖𝑘‖|𝑤_𝑗−1(𝑠)|𝑑𝑠≤‖𝑦₀‖‖𝑘‖𝑡+‖𝑘‖‖𝑦₀‖^𝑡₀(^𝑗−1_𝑛=1𝑠^𝑛𝑛!‖𝑘‖^𝑛)𝑑𝑠≤‖𝑦₀‖‖𝑘‖𝑡+‖𝑦₀‖(^𝑗−1_𝑛=1𝑡^𝑛+1(𝑛+1)!‖𝑘‖^𝑛+1)≤‖𝑦₀‖(^𝑗_𝑛=1𝑡^𝑛𝑛!‖𝑘‖^𝑛). (3.14) On the other hand, 𝑤^′_𝑗(𝑡)=𝑄_𝑛𝑗2(𝑘(𝑡,𝑡)(𝑦₀(𝑡)+𝑤_𝑗−1(𝑡))), (3.15) and thus ‖𝑤^′_𝑗‖≤‖𝑘‖(‖𝑦₀‖+‖𝑤_𝑗−1‖). (3.16) Finally, the bounding of {𝑀_𝑗}_𝑗∈ℕ follows from ‖𝜕(𝑘(𝑦₀+𝑤_𝑗))𝜕𝑡‖≤‖𝜕𝑘𝜕𝑡‖(‖𝑦₀‖+‖𝑤_𝑗‖‖),𝜕(𝑘(𝑦₀+𝑤_𝑗))𝜕𝑠‖≤‖𝜕𝑘𝜕𝑠‖(‖𝑦₀‖+‖𝑤_𝑗‖)+‖𝑘‖(‖𝑦^′₀‖+‖𝑤^′_𝑗‖). (3.17)

Letting 𝛽=sup{𝑀_𝑗}_𝑗∈ℕ and 𝛿=4𝛽∑_𝑗≥0‖𝑇‖^𝑗, from (3.10) we have

^𝑚−1_𝑗=0‖𝑇‖^{𝑚−1−𝑗}‖𝑇(𝑦₀+𝑤_𝑗)−𝑤_𝑗+1‖≤𝛿ℎ. (3.18)

Then Propositions 3.1 and 3.2 imply the following result, which gives an upper bound of the error.

Theorem 3.3. Suppose 𝑘,𝑦₀∈𝐶¹([0,1]²) and let {𝑤₀,𝑤₁,…,𝑤_𝑚} be the functions defined in (3.1). Then, there exists 𝛿>0 such that ‖𝑥−(𝑦₀+𝑤_𝑚)‖≤‖𝑇^𝑚+1‖‖(𝐼𝑑−𝑇)⁻¹‖‖𝑦₀‖+𝛿ℎ. (3.19)

4. Numerical Examples

The behaviour of this method is illustrated by means of the following two examples. The computations associated with the numerical experiments were carried out using Mathematica 7.

The chosen dense subset of [0,1] is

1{0,1,2,14,341,…,2^𝑘,32^𝑘2,…,^𝑘−12^𝑘,…}. (4.1)

To construct the functions 𝑤_𝑗 introduced in Section 3, we have taken 𝑛_𝑗=𝑛 in the expression (3.1). For such a choice the value of ℎ in Section 3 is 1/(𝑛−1).

In both cases we exhibit, for 𝑛=9,17,33, and 65 the absolute errors committed in eight representative points of [0,1] when we approximate the exact solution 𝑥 by the function 𝑦₀+𝑤_𝑚, where 𝑚 has been determined in the following way: for each 𝑗∈ℕ we note

𝐸_𝑛,𝑗=max_{𝑖=1,…,𝑛}|(𝑦₀+𝑤_𝑗)(𝑡_𝑖)−𝑥(𝑡_𝑖)|, (4.2) and choose 𝑚 satisfying

𝐸_𝑛,𝑚𝐸_𝑛,𝑚+1<1+10⁻². (4.3)

Example 4.1. The equation (see [12]), 𝑥(𝑡)=𝑡⁵−𝑡⁸7+^𝑡₀𝑠𝑡𝑥(𝑠)𝑑𝑠,(𝑡∈[0,1]) (4.4) has exact solution 𝑥(𝑡)=𝑡⁵; see Table 1.
In Table 2 we compare the approximations 𝑦₀+𝑤_𝑚 obtained as described above with the collocation solution 𝑢 shown in [12] in terms of the norm ‖𝑣‖_ℎ,2∶=(ℎ^𝑛_𝑘=1𝑣(𝑡_𝑘)²)^1/2 (4.5) introduced in the cited article. We point out that ℎ=1/(𝑛−1) in both methods.

Table 1 

Numerical results for Example 4.1.

Table 2 

Compared results.

Example 4.2. Second example is taken from [13]. We consider the equation 𝑥(𝑡)=1+^𝑡₀(𝑠−𝑡)𝑥(𝑠)𝑑𝑠,(𝑡∈[0,1]), (4.6) whose exact solution is 𝑥(𝑡)=cos(𝑡); see Table 3.
The computed results by the suggested method for 𝑛=9 (ℎ=1/8) improve the obtained ones in [13] for ℎ=1/10.

Table 3 

Numerical results for Example 4.2.

Acknowledgments

This Research is Partially supported by M.E.C. (Spain) and FEDER project no. MTM2006-12533, and by Junta de Andalcía Grant FQM359.