Continuous-Stage Runge-Kutta approximation to Differential Problems

In recent years, the efficient numerical solution of Hamiltonian problems has led to the definition of a class of energy-conserving Runge-Kutta methods named Hamiltonian Boundary Value Methods (HBVMs). Such methods admit an interesting interpretation in terms of continuous-stage Runge-Kutta methods, which is here recalled and revisited for general differential problems.


Introduction
Continuous-stage Runge-Kutta methods (csRK, hereafter), introduced by Butcher [34,35,36], have been used in recent years as a useful tool for studying energy-conserving methods for Hamiltonian problems (see, e.g., [22,25,38,41,42,44,43,45,46]). In particular, we shall at first consider methods, within this latter class, having Butcher tableau in the form: are suitable functions defining the method. Hereafter, we shall use the notation (2), in place of the more commonly used a cτ , b c , in order to make clear that these are functions of the respective arguments. For later use, we shall also denoteȧ (c, τ ) = d dc a(c, τ ).

Approximation of ODE-IVPs
Let us consider, at first, the initial value problem (4). As done in [2], for our analysis we shall use an expansion of the vector field along the orthonormal 1 In the sequel, we shall assume f to be analytical. 2 Clearly, because of consistenty, one requires that 1 0 b(c)dc = 1.
Legendre polynomial basis: where, as usual, Π i denotes the vector space of polynomials of degree i, and δ ij is the Kronecker symbol. Consequently, we can rewrite (4) aṡ with Integrating side-by-side, and imposing the initial condition, one then obtains that the solution of (10) is formally given by For c = 1, by considering that, by virtue of (9), 1 0 P j (x)dx = δ j0 , and taking into account (4), one obtains: i.e., the Fundamental Theorem of the Calculus. Interestingly, by setting one obtains that (10)- (11) and (12)- (13) can be rewritten, respectively, as: is the csRK "method" providing the exact solution to the problem.

Vector formulation
For later use, we now cast the formulation of (14) in vector form. For this purpose, let us introduce the infinite vectors, and the matrix, also recalling that, by virtue of (9), and due to the well-known relations between the Legendre polynomials and their integrals, with I the identity operator. Consequently, with reference to (14), one has and, in particular, one may regard the Butcher tableau, equivalent to (16) as the corresponding W -transformation [39] of the continuous problem. It is worth mentioning that, by defining the infinite vector (see (11)), then (12) can be rewritten as and, consequently, the vector γ satisfies the equation with (compare with (15)),

Polynomial approximation
In order to derive a polynomial approximation σ ∈ Π s of (15), it suffices to truncate the infinite series in (14) after s terms: so thatȧ and, therefore, (15) is approximated bẏ
then H(y 1 ) = H(y 0 ). In the case of Hamiltonian problems, H is the energy of the system. Consequently, the csRK method (28) is energy-conserving, as is shown in [25, Theorem 3]).
A corresponding vector formulation of the csRK method (28) can be derived by replacing the infinite vectors and matrix in (17)-(18) with . . .
and the matrices, such that with I r ∈ R r×r the identity matrix. Consequently, with reference to (25) and, in particular, one may regard the Butcher tableau, equivalent to (28) as the corresponding W -transformation of a HBVM(∞, s) method [2]. Further, according to [23], setting the vector with γ j (σ) defined as in (11), by formally replacing y by σ, one has that such vector satisfies the equation (compare with (23)) with (compare with (27) and (24)),

Approximation of special second-order ODE-IVPs
An interesting particular case is that of special second order problems, namely problems in the form which, in turn, are a special case of the more general problems Let us study, at first, the special problem (38), which is very important in many applications (as an example, separable Hamiltonian problems are in such a form), then discussing problem (39). Settingẏ(t) = p(t), and expanding the right-hand sides along the Legendre basis gives: Remark 3 It is worth mentioning that, with reference to (3) and (20), the two previous equations can be rewritten as: Integrating side by side (40), and imposing the initial condition then gives: with the vector γ formally still given by (22). Substitution of the second equation in (42) into the first one, taking into account that 4 In general, e i will denote the infinite vector whose jth entry is δ ij . and considering (19) and again (22), then gives, one has that (42) can be rewritten as: Remark 4 We observe that, from the second equation in (40) and (41), and considering (20), one derives: Moreover, it is worth mentioning that, see (44), (14), and (19), From (40) and (45) one obtains that the values at h will be given by (see (18)): y(h) ≡ p(h) =ẏ 0 +he 1 ⊤⊗I m γ =ẏ 0 +hγ 0 (y) ≡ẏ 0 +h 1 0 f (y(ch))dc, (48) and In other words, the exact solution of problem (38) is generated by the following csRKN "method":
(54) Consequently, we obtain the general csRKN method in place of (50). The bad news is that now the system of equations (53) has a doubled size, w.r.t. (45). On the other hand, the good news is that, upon modifying the definition of the vector γ in (22) as follows, one has that the equations (53) can be rewritten as Consequently, we obtain again a single equation for the vector γ defined in (56), Similarly, (48)-(49) respectively become: Remarkably enough, the equations (57)-(58) are very similar to (51)-(52).
An alternative formulation of the method (66) can be obtained by repeating arguments similar to those used in Section 3.1. In fact, by setting the vector (compare with (35)-(37)) . . .
by virtue of (60)-(61) such a vector satisfies the equation with the new approximations given by (compare with (52)): For the general problem (39), following similar steps as above, the generalized csRKN method (55) becomes, with reference to (33): Also in this case, we can derive a vector formulation, similar to (67)-(68). In fact, by defining the vector . . .
it satisfies the equation with the new approximations given by: Remark 7 By comparing the discrete problems (68) and (72), one realizes that they have the same dimension, independently of the fact that the latter one solves the general problem (39). This fact is even more striking since, as we are going to sketch in the next section, this will be the case for a general kth order ODE-IVP.

Approximation of general kth-order ODE-IVPs
Let us consider the case of a general kth order problem, y (k) (ch) = f (y(ch), y (1) (ch), . . . , y (k−1) (ch)), c ∈ [0, 1], Repeating similar steps as above, by defining the infinite vector 5 one has that it satisfies that equation with the values at t = h given by: , j = 0, . . . , k − 1 − i, being the (j + 1)st entry on the first row of the matrix A polynomial approximation of degree s (resulting, as usual, into an order 2s method) can be derived by formally substituting, in the equation (76), P ∞ (τ ) and I ∞ (τ ) with P s (τ ) and I s (τ ), respectively, 6 with the new approximations y (i) 1 ≈ y (i) (h), given by: satisfying the equation (compare with (36)) with (compare with (37)), 6 In so doing, the vector γ now belongs to R sm .
We observe that, by introducing the matrices (see (30)) one has, with reference to (31), Further, by also introducing the vector e = (1, . . . , 1) ⊤ ∈ R k , one obtains that (80) can be rewritten aŝ the vector of the stages of the corresponding Runge-Kutta method, from (84)-(85) one obtains the stage equation with the new approximation Summing all up, (84)-(87) define the k-stage Runge-Kutta method named HBVM(k, s) [15,16,22]. It is worth mentioning that the Butcher matrix of the method, with reference to the coefficients of the csRK (33), is given by: In particular, when k = s one obtains the s-stage Gauss collocation method. However, the use of values k > s (and even k ≫ s) is useful, in view of deriving energy-conserving methods for Hamiltonian systems [22,25,15,16].

Remark 8
As is clear, the formulation (80)-(81) is computationally more effective than (86)-(87), having the former (block) dimension s, independently of k, which is the dimension of the latter formulation. As observed in [23], this allows the use of relatively large values of k, without increasing too much the computational cost.
Similar arguments can be repeated in the case of the polynomial approximations for problems (38), (39), and (74): we here sketch only those concerning the csRKN (66), providing the correct implementation for the HBVM(k, s) method for the special second-order problem (38). By formally using the same approximate Fourier coefficients (79) in place of (67), one has that (68) is replaced by the following discrete counterpart, with the new approximations (compare with (69)) given bẏ Similarly as in the first order case, (89) can be rewritten aŝ where c is the vector of the abscissae. Again, the vector is the stage vector of a k-stage RKN method. In particular, from (91) and (92) one obtains and some algebra shows that the new approximations are given bẏ Consequently, we are speaking about the following RKN method: 7 The Butcher tableau (95) defines the RKN formulation of a HBVM(k, s) method [2,23,15]. Is is worth mentioning that, with reference to (60) and (65) 1 0 a s (c i , τ )a s (τ, c j )dτ ∈ R k×k . 7 Here, • denotes the Hadamard, i.e. componentwise, product.

Remark 9
As observed in Remark 8, also in this case, in consideration that k ≥ s, the formulation (90)-(91) of the method is much more efficient than the usual one, given by (93)-(94), also in view of the use of values of k ≫ s, as in the case of separable Hamiltonian problems.
We conclude this section by recalling that: • very efficient iterative procedures exist for solving the discrete problems (84) and (91) generated by a HBVM method (see, e.g., [23]); 8 • when HBVMs are used as spectral methods in time, i.e., choosing values of s and k > s so that no accuracy improvement can be obtained, for the considered finite-precision arithmetic and timestep h used [3,18,29], then there is no practical difference between the discrete methods and their continuous-stage counterparts.

Conclusions
In this paper we have recalled the basic facts, also reporting some new insight, on the energy-conserving class of Runge-Kutta methods named HBVMs. Such method have been here studied within the framework of continuous-stage Runge-Kutta methods. The extension to second-order problems has been also studied, providing a natural continuous-stage Runge-Kutta-Nyström formulation of the methods. Further, also the extension to general kth-order problems has been sketched. The relation with the fully discrete methods has been also recalled, thus showing the usefulness of using such a framework to study the fully discrete methods.