{'fraction_non_alphanumeric': 0.05030342911457191, 'fraction_numerical': 0.01520801345324267, 'mean_word_length': 4.289056457849961, 'pattern_counts': {'":': 0, '<': 1, '<?xml version=': 0, '>': 1, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 1, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 6, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Addition theorems for Ck real functions and applications in ordinary differential equations

This work establishes the existence of addition theorems and double-angle formulas for Ck real scalar functions. Moreover, we determine necessary and sufficient conditions for a bivariate function to be an addition formula for a Ck real function. The double-angle formulas allow us to generate a duplication algorithm, which can be used as an alternative to the classical numerical methods to obtain an approximation for the solution of an ordinary differential equation. We demonstrate that this algorithm converges uniformly in any compact domain contained in the maximal domain of that solution. Finally, we carry out some numerical simulations showing a good performance of the duplication algorithm when compared with standard numerical methods


Introduction
We recall that a real function φ(t) has an addition theorem if φ(t + τ ) may be recovered from φ(t) and φ(τ ) through the following relation (1) G(φ(t), φ(τ ), φ(t + τ )) = 0, wherever the domains of G(x, y, z) and φ(t) are defined. Hereinafter, the function G(x, y, z) defining the addition theorem will be dubbed as the addition formula associated to the function φ(t). It is well known that trigonometric or exponential functions are examples of functions endowed with these features. For instance, the exponential function φ(t) = e at , with a ∈ R \ {0}, has the addition formula G(x, y, z) = z − xy = 0, because e a(t+τ ) − e at e aτ = 0 for all t, τ ∈ R. Traditionally, the theory of elliptic functions investigates the existence of addition formulas [8,12,13,15]. Indeed, quoting Weierstrass [19], we have that: "the problem of the theory of elliptic functions is to determine all functions of the complex argument for which there exists an algebraic addition theorem". In this regard, the following statement is proved in [14]. Every single-valued analytic function with an algebraic addition theorem is either an elliptic function or an algebraic function of x or e ax .
Equation (1) can be denominated as a scalar addition theorem because it depends only on one function. However, this equation can be generalized and consider an equation of the form G(φ 1 (t), φ 2 (τ ), φ 3 (t + τ )) = 0, where φ i (t) ∈ R N for i = 1, 2, 3, which are called vector addition theorems [2]. In [2,3] the authors related these kind of functional equations with the theory of one-dimensional integrable systems. Although in this paper we will restrict to the case of scalar addition theorems (1) we consider this investigation as the basis for the realization of a project that includes the case of vector addition theorems from the point of view of differential equations, which as we will show in this work could provide different perspectives for applications.
In the theory of addition theorems we distinguish between direct and inverse problems. More precisely, the problem of finding an addition theorem for a given function φ(t) will be called the direct problem of addition theorems. Conversely, one may also rise the question of whether a given function G could be an addition formula, which is referred to as inverse problem of addition theorems.
Although we have considered in (1) the general definition of an addition formula, we will focus throughout this paper on explicit addition theorems or formulas. That is to say, we call a bivariate function R(x, y) an explicit addition formula for φ(t) when relation (1) may be expressed as (2) φ(t + τ ) = R (φ(t), φ(τ )) , wherever the domains of the involved functions R(x, y) and φ(t) are defined. Therefore, we are interested in the direct and inverse problems for explicit addition formulas. The interest in the solutions and properties of explicit addition theorems is a classical problem. For instance, by fixing R(x, y) as a polynomial, rational or algebraic type, the existence and uniqueness of solutions for equations of the form (2) are discussed in [1,Sec. 2.2.4], and references therein, where a treatment of addition theorems is widely carried out from the point of view of functional equations.
The existence of addition theorems and double-angle formulas, which are the particular case of considering t = τ in (1) or (2), is of great practical interest. For instance, they are used in the construction of algebraic first integrals for differential equations. This fact is pointed out by Tsiganov [16], where the author recalls that the integrals of the Kepler system, given by the eccentricity and the angular momentum vectors, are particular cases of this procedure. Moreover, before the age of computers, double-angle formulas were traditionally used to compute tables of trigonometric, exponential or logarithm functions. Still today, this strategy is widely used in the fast and efficient computation of elliptic functions and integrals, see for example [4,5,6,7] and the references therein. Since these works relies on the use of double-angle formulas, their applicability could be extended to any function having a double-angle formula.
In this paper we deal with the problem of finding addition and double-angle theorems for C k real functions, with k ∈ N ∪ {∞, ω}, where C ω stands for analytic functions. In the applications we pay special attention to functions coming from solutions of autonomous ordinary differential equations. This problem was previously explored in [20]. In [17] the case of linear homogeneous ordinary differential equations with constant coefficients was analyzed, and also in [18], where the author places the binomial theorem and the addition theorems for exponential, trigonometric, and hyperbolic functions in the context of a single addition theorem generated by an initial value problem. Our approach allows to consider any autonomous ordinary differential equations of class C k .
As an application of the double-angle formulas, we propose an alternative numerical approximation of the solution of an initial value problem based on the duplication algorithm. This scheme requires the double-angle formula of the solution to be available, which restricts the applicability of the numerical method. Nevertheless, we will use the Taylor expansion of double-angle formulas associated to the solution of an initial value problem. Therefore, extending the applicability of the duplication algorithm and allowing to generate a polygonal approximation of the named solutions. In addition, we also establish the convergence of the proposed duplication algorithm.
This paper is organized as follows. In Section 2, we deal with results concerning the direct and inverse theory of addition theorems. Section 3 contains the Taylor expansion for double-angle formulas associated to the solution of an initial value problem, which allows to define a duplication algorithm in Section 3.2. Finally, in Section 4 we carry out numerical experiments comparing the duplication algorithm with the standard numerical methods provided by Wolfram Mathematica.

On the Theory of Addition Theorems
This section presents some results about addition theorems as existence or necessary and sufficient conditions for a function R(x, y) to be an explicit addition formula.
2.1. Direct Problem of Addition Theorems. The following theorem gives an affirmative answer to the question of whether exists addition theorems for a given function. This problem was tackled in [1], see page 256, as a functional equation of the form (2). Here we specialized to the case of C k functions, which allows for a shortened proof.
Previous to the existence theorem, we need some auxiliary notation and results. Let I ⊂ R is an open interval. If φ : I −→ R is a C k function with non vanishing derivative on I, then the inverse function theorem says that φ : I −→ ∆, where ∆ := φ(I) ⊂ R, is a diffeomorphism with a C k inverse function φ −1 : ∆ −→ I satisfying that φ −1 • φ = Id I and φ • φ −1 = Id ∆ . If we assume that the set {(t, τ ) ∈ I × I | t + τ ∈ I} is non empty, we can define the non empty set We note that the function is of class C k and that D I = χ −1 (I), which implies that D I is an open set. Therefore, in what follows D I will be named as the addition domain of φ(t).
Note that if a given function φ(t) is endowed with an explicit addition formula R(x, y), the above set D I is the maximal domain in which R(x, y) can be defined.
Proof. We define the functions Both functions are well-defined because of the hypothesis and they are of class C k . Hence, we have the following commutative diagram Then, for all t, τ, t + τ ∈ I, we have that φ(t), φ(τ ) ∈ D I . Thus, This proves that R(x, y) is a C k explicit addition theorem for φ(t). Finally, from the definition of Ψ and Φ we have The proof is completed.
Remark 2.2. Note that, for the case I = (a, b) with a, b > 0 and a ≥ b/2 the domain D I is empty. An analogous situation occurs for a, b < 0. Therefore, thinking in the applications, we will consider the case of I being a neighborhood of the origin, which guarantees D I = ∅. Remark 2.3. Theorem (2.1) gives the existence of a C k addition formula wherever it makes sense to add t, τ ∈ I. In this regard we claim that it is not local but a global result for the case of strictly monotone functions.

2.2.
Inverse Problem of Addition Theorems. Along this section we study necessary and sufficient conditions for a function R(x, y) to be an explicit addition formula for a C k function.
Theorem 2.4 (Necessary Conditions). Given an interval I containing the origin, a real C k function φ : I −→ R, with non empty addition domain D I and φ(0) = φ 0 , whose derivative vanishes only in a discrete subset of I, and its addition formula R(x, y) : D I −→ R, the following statements hold: Proof. By hypothesis we have the explicit addition theorem (2), that is, Hence, the symmetry of R(x, y) is obtained from the fact that φ(t + τ ) = φ(τ + t). Additionally, for all t, τ ∈ I we have Thus, statement (ii) follows by using φ(t) = φ(t + 0) = x and φ(τ ) = φ(0 + τ ) = y. We now differentiate (2) with respect to t and τ to obtain and equation (5) with τ = 0 becomeṡ From these two previous equations, and by considering t = τ in the last one, we obtain Therefore, by putting z = φ(τ ) we get statement (iii). Finally, statement (iv) is already given in [1]. We take t 1 , t 2 , t 3 ∈ I and such that x = φ(t 1 ), y = φ(t 2 ) and z = φ(t 3 ) are in D I . Then, the statement follows from the following equalities In [10], functions R(x, y) satisfying condition (iv) of previous theorem are analyzed, showing that they must also satisfy condition (i). Moreover, if such a function satisfies that R(0, 0) = 0, we have that it is given by one of the following possibilities where S 1 and S 2 are symmetric functions.
Next result follows directly from Theorem 2.4, particularly from (4) and (5). It says that a differentiable function φ(t) having an explicit addition formula is also the solution of a particular autonomous first order ordinary differential equation.
Moreover, φ(t) is the solution of the following initial value probleṁ Remark 2.6. The function R(x, y) may be found as the solution of the partial differential equation given by (6) y) is the explicit addition formula for the solution of the initial value problem Proof. For notational simplicity, we do not display the dependence of f (x; x 0 ) on x 0 . The function f (x) in (7) comes from condition (ii). We will prove the theorem by assuming that Let x(t) be the solution of (7). We define the following functions which are of class C k , with k ≥ 1 and k ≥ 2, respectively. We will show that φ(t) = ψ(t) for all t ∈ I, being I the maximal interval where x(t) is defined. As a consequence, we will obtain that x(tα + tβ) = R(x(tα), x(tβ)), and therefore, by setting α = 1 and τ = tβ, By using the above definition of φ(t), a straightforward computation shows that φ(t) is the unique solution of the following initial value problem Then, by applying condition (ii) we rewriteψ(t) aṡ Hence, by using (iii), (ii) and the definition of ψ(t) we geṫ Additionally, condition (i) implies that ψ(0) = R(x 0 , x 0 ) = x 0 Therefore, we have proved that ψ(t) satisfies the initial value problem From (8) and (9) is clear that φ(t) and ψ(t) are the solutions of the same initial value problem. Hence we have ψ(t) = φ(t). This completes the proof 3. Double-angle Formulas and applications.
In this section we will show how the results of previous section can be used to get an algorithm, based on the existence of double-angle formulas, to obtain accurate approximations of solutions of initial value problems in autonomous first order differential equations.
Recall that a double-angle formula for a real function φ(t) is a function R(x) satisfying the following relation φ(2t) = R(φ(t)). Under some assumptions the existence of a double-angle formula is established, as it showed by the next result, which is a straightforward consequence of Theorem 2.1. In what follows, we will denote both the explicit addition formula and the double-angle formula of a function φ(t) with the same symbol R, and the context will make clear which one are we referring to.
In general, the problem of finding the addition or double-angle formula of a given function is a difficult task. However, the previous results show the relationship between the direct and inverse problems of addition theorems for solutions of autonomous first order ordinary differential equations. More precisely, let x(t) be the solution of the initial value problem If f (x 0 ) = 0, then the constant function x(t) ≡ x 0 is the solution to the problem and trivially R(x, y) = x 0 is its explicit addition formula. If f (x 0 ) = 0, then Theorem 2.1 ensures the existence of R(x, y) for x(t). Recall that the exponential or trigonometric functions are solutions of this kind of initial value problems, which have explicit addition formulas.
On one hand, while the existence of R(x, y) is guaranteed, knowing its concrete expression is not easy and only a few examples are known in the literature. On the other hand, as far as we know, there is no systematic method to compute the expression of R(x, y) for a given function φ(t) or for the solution x(t) of (10), with f (x) an arbitrary function. Nevertheless, this problem could be tackled for the case when f (x) has a sort of simplicity. For instance, if n ∈ N, then for f (x) = x n+1 and f (x) = x 1/(n+1) we obtain that are the concrete explicit addition formulas for the solution x(t) of (10).
Analogously to the case of explicit addition formulas, the existence of a double-angle formula R(x), for the solution x(t) of (10), is guaranteed by Corollary 3.1 and its concrete expression of could be difficult to obtain. Nevertheless, we point out that it is possible take advantage of the existence of the double-angle formula R(x) to derive a procedure to compute the Taylor polynomial of any degree of R(x) at x 0 . Such a Taylor polynomial can be used to produce an algorithm that allow us to provide an approximation of x(t). More precisely, these Taylor polynomials can be used to produce a succession of functions converging uniformly to x(t) in a certain compact domain D, such a convergence and an iterative procedure of the Taylor polynomial of enough large degree allow us to give a accurate approximation to x(t) in D, which is contained in the maximal interval of the solution x(t). This mentioned numerical scheme will be dubbed as duplication algorithm which we will explain in the follows.
We will describe the procedure to compute the Taylor polynomial of R(x). We consider the non trivial case, that is, f (x 0 ) = 0 and by simplicity, we assume that f (x) in (10) is analytic, and so are x(t) and R(x). To obtain the Taylor polynomial of R(x), we compute the successive derivatives of R(x) at x 0 . Firstly, from (11) it is clear that We now take the derivative of (11) with respect to t: here R (·) means the derivative of R(x) with respect to x. Then, by using (10) in the leftand right-hand sides we get (13) 2f (x(2t)) = R (x(t))f (x(t)).
Then, the derivative of this last expression with respect to t is which, by using (10) and by setting t = 0, becomes Finally, by simplifying terms and by using (14) we arrive to The third order derivative R (x 0 ) is easily calculated by executing the same procedure Hence, by using this last equation and equations (12), (14), and (15), the 3th order Taylor polynomial of R(x) at x 0 is Following the same methodology, higher order derivatives of R(x) at x = x 0 can be computed. Thus, allowing to obtain the Taylor polynomials of arbitrary order for f (x) analytic and up to order n ≤ k for f (x) of class C k . We include the expression of the derivatives up to the 10th order in the Appendix A.
3.2. The duplication algorithm. In this section, we use double-angle formulas to provide an estimation for the solution x(t) of an initial value problem.
Our method is based on the use of the double-angle formula R(x) to generate local approximations for the function x(t) in a compact subset D containing t = 0 and contained in the maximal interval x(t). The following procedure will be dubbed as the duplication algorithm. Specifically, we compute x(t) in a small neighborhood V 0 of the origin; fort ∈ D far from the origin, we dividet by two as many times as it is needed to gett/2 n ∈ V 0 , and we compute x(t/2 n ); then, we can recover x(t) by iterating n-times the double-angle formula. That is, (16) x where R n denotes the n-th iteration of the double-angle formula R associated to x(t).
The previous process depends on the availability of the double-angle formula R, as well as the capacity of computing x(t) in V 0 . However, when x(t) is unknown, an approximation may be obtained by using a Taylor polynomial (17) x(t ) ≈ R n x m 1 t 2 n , being x m 1 (t) the m 1 -th order Taylor polynomial of x(t) at t = 0. Moreover, for the case in which we are not able to find the exact double-angle formula R, we will use its Taylor polynomials. For this purpose, the previous section gave an easy method to compute arbitrary derivatives of R at x = x 0 . This situation leads us to the following local approximation where R n m 2 is the n-th iteration of the m 2 -th order Taylor polynomial of R. In fact, straightforward computations prove that (19) x Keeping in mind the assumption that f (x) in (10) is analytic, we have that x(t) is analytic in its maximal domain I M = (α, β). Thus, the associated double-angle formula R(x) is analytic in D I M = x(α/2, β/2). Now we construct a polygonal function defined in a compact interval D = [a, b] ⊂ I M to estimate x(t), which is based on local approximations and their interpolation. More precisely, given the r/2 n -sized partition P(D) = {t 0 , t 1 , . . . , t 2 n }, with t 0 = a, t 2 n = b and r Remark 3.2. The duplication algorithm uses a Taylor approximation of x(t) in a small neighborhood of the origin V 0 = [−r 0 , r 0 ]. Note that the partition P(D) can be always be chosen in such a way that t i /2 n ∈ V 0 , for any t i ∈ P(D). In what follows we will assume that this condition is satisfied.
Remark 3.3. The above function R n m 2 admits an alternative definition as the m 2 -th order Taylor polynomial of the function R n , which can be obtained by combining the successive derivatives given in Appendix A and the formula of Faà di Bruno [9]. In that case, the numerical results could be improved. However, this subject is left for further study.
In the proof of the next theorem we need the following auxiliary polygonal function Proof. Let us consider t ∈ [t i , t i+1 ] ⊂ D. Then by the triangle inequality we have where P n (t) is defined in (21). By using (20), we get Hence, from (19) it follows that and we obtain the following upper bound by considering that |t j | ≤ r and that m ≥ 2, To finish the proof, we look for an upper bound of the second term in (22) where ∆x i := x(t i ) − x(t i+1 ) and ∆t i := t i − t i+1 . If now add and subtract t i+1 ∆t i x(t i+1 ), then we get The last equality follows from the fact that both factors in the previous step satisfy that The proof has been completed.
In previous proof is stated that We have showed that P m 1 ,m 2 n (t) − P n (t) depends on the minimum of the orders of the Taylor polynomials involved and P n (t) − x(t) depends on the size of the partition P(D).
The numerical experiments carried out in next section will show such dependence.

Numerical Simulations
In this section we analyze the efficiency and accuracy of the duplication algorithm as an approximation method for the solution of ordinary differential equations. We do so by comparing it with standard numerical methods. More precisely, our experiments have been carried out by using the software Wolfram Mathematica, version 12.1.1.0 Linux x86 (64-bit). This software is running on the platform Intel® Core™ i3-9100 CPU, 3.60GHz×4.
In all our simulations we consider the following notation: • x(t) is the exact solution of the initial value problem (10).
•x(t) is the approximate solution to (10) given by the duplication algorithm having the exact double-angle formula.
•x(t) is the approximate solution to (10) by using the duplication algorithm with the Taylor approximation of the double-angle formula.
•x(t) is the numeric approximation to (10) provided by NDSolve of Mathematica.

First Example. We considerẋ
The exact solution x(t) of this problem and its associated double-angle formula R(x) can be easily obtained after some algebraic manipulations. They have following expressions In order to test the accuracy and efficiency of the duplication algorithm in this toy example, we will compute the errorsx(t) − x(t),x(t) − x(t) andx(t) − x(t) as well as the computer time of these errors on a compact interval contained in the domain of x(t). To perform the numerical simulations we will consider x 0 = 1, thus x(t) is defined in (−∞, 1). It is easy to check that t = 1 is an asymptote for x(t). With the aim of studying the behaviour of the duplication algorithm in a generic interval and for t near the asymptote, we consider simulations in the intervals I 1 = [−0.5, 0.5] and I 2 = [−0.5, 0.99].
The approximationx(t) has been obtained by partitioning the integration intervals I 1 and I 2 in 240 and 10.000 points respectively. As we can see in Figure 1, the difference between the exact solution and the approximation by the duplication algorithm, with the exact double-angle formula R(x), produces an error of order of magnitude 10 −9 with a computer time of 0.0019 seconds in I 1 . Note that, as we approach the asymptote t = 1 in the interval I 2 , the accuracy of the simulation get worse and the computational effort increases.  In Figure 2 we show the error propagation between the exact solution x(t) andx(t), which we obtain by numerical simulation applying the duplication algorithm with an approximation of R(x) given by a 20-th order Taylor polynomial. The error in this simulation behaves as in the case ofx(t) with almost identical computer times.  Remarkably, the order of magnitude of the error in the previous procedures is less than the one provided by Wolfram Mathematica. Indeed, Figure 3 shows the difference of the exact solution versus the standard numerical integration method NDSolve Mathematica commad with options: Method = Automatic, AccuracyGoal= 20 and MachinePrecision which has error of order of magnitude 10 −8 and 10 −3 in I 1 and I 2 respectively, with computer times 0.0024 and 0.0020 seconds.
In the following experiments we analyze the relation between the number of interpolation point and the accuracy of the approximation. To this aim we consider again x 0 = 1, in Figure 4 we simulate x(t) in the interval [−0.5, 0.5] using 20th order Taylor polynomial for the double-angle formula. As we increase the number of points in the interpolation process, computation time and accuracy increase linearly. However, for more that 4000 interpolation points, the error does not get better than 10 −14 .

Second Example.
In the previous example we dealt with a simple equation and duplication algorithm. This was done with the aim of getting all the elements needed for the comparisons in a straightforward manner. Our next experiment addresses a more elaborated model, whose initial condition is close to the limit of the theory applicability. Precisely, the second example is described by the following initial value problem where | | 1. As the reader may notice, the case = 0 implies thatẋ(0) = 0. Thus, x(t) does not fulfill the hypothesis of Corollary 3.1, neither the initial value problem (23) satisfies the theorem of existence and uniqueness for ordinary differential equations. In this section we will consider an initial condition as close as possible to x 0 = −5/12 to check if the performance of our method is affected by the limits of its applicability. Equation (23) is written in the same way as the derivative of the ℘-Weierstrass elliptic function. However, the initial condition is evaluated in t = 0 preventing x(t) to have a pole in the origin as it is the case for ℘. Thus, we can consider the solution of the above equation as the restriction to the real domain of a ℘-Weierstrass elliptic function for which the lattice has been displaced from the origin in such a way that there is no poles in the real line. Alternatively, the explicit solution can be obtained in terms of the Jacobi elliptic functions [11]. Indeed, we have the following expression for the unique solution of (23) where δ and are related by the following formula By using the double-angle formula of the Jacobi dn-function, and after some cumbersome algebraic manipulations, we get the exact double-angle formula for (24) Therefore, as in the previous example, we can make comparisons between the exact solution x(t) of (23), the standard numerical approximation provided by Wolfram Mathematicax(t) and our duplication algorithm.  In Figure 5 we use the exact double-angle formula (25) to simulate the solution of (23) obtained with the duplication algorithm, then we compare it with the approximation given by the software package Wolfram Mathematica. We observe that the accuracy of the duplication algorithm is two orders of magnitude better. Moreover, the computer time is on the same order for both approaches, with a slight improvement in the duplication algorithm.
We also investigate how the order of the Taylor polynomials for the double-angle formula improves the approximations. In Figure 6, we show the numerical experiments for orders 10, 15, 20 and 30. Setting the number of interpolation points to 100, we obtain that the error remains fixed around 10 −10 for Taylor orders above 30.
Finally, we assess the role of the number of interpolation points in the accuracy of the approximation that we provide. More precisely, in Figure 7 we use a Taylor polynomials of order 30th to approximate the double-angle formula. Then we show the evolution of the error when comparing with the exact solution for 80, 160, 320, 640, 1280 and 2560 interpolation points. As we duplicate the number of points, the error magnitude is approximately multiplied by 10 −1 . However, our experiments shows that, for this fixed Taylor approximation of the double-angle formula, the precision stabilizes around 10 −14 as we increase the number of interpolation points. .