## Abstract

We study the onset of sustained oscillations in a classical state-dependent delay (SDD) differential equation inspired by control theory. Owing to the large delays considered, the Hopf bifurcation is singular and the oscillations rapidly acquire a sawtooth profile past the instability threshold. Using asymptotic techniques, we explicitly capture the gradual change from nearly sinusoidal to sawtooth oscillations. The dependence of the delay on the solution can be either linear or nonlinear, with at least quadratic dependence. In the former case, an asymptotic connection is made with the Rayleigh oscillator. In the latter, van der Pol’s equation is derived for the small-amplitude oscillations. SDD differential equations are currently the subject of intense research in order to establish or amend general theorems valid for constant-delay differential equation, but explicit analytical construction of solutions are rare. This paper illustrates the use of singular perturbation techniques and the unusual way in which solvability conditions can arise for SDD problems with large delays.

## 1. Introduction

State-dependent delays (SDDs) are a common occurrence in natural or artificial systems subjected to feedback. In motion control, for instance, the delay of a signal reflected by an obstacle depends on the current position of the vehicle [1]. Another well-documented area of research where SDD problems are formulated concerns chatter instabilities. In cutting processes and in milling or rotary drilling systems used by oil industries, SDDs arise owing to the deformation of the cutting tool [2–4]. Furthermore, in population dynamics, the time required for the maturation level of a cell to achieve a certain threshold can be modelled as an SDD [5–7]. Finally, the electromagnetic two-body problem stands as a fundamental SDD problem in physics which has been regularly revisited over the years [8–11]. The mathematical analysis of SDD problems is still fairly recent compared to constant-delay differential equations and requires the development of new theoretical [12] and numerical [13] tools.

Some of these SDD models differ considerably from others, and most of them do not look simple (see the review by Hartung *et al.* [12]). In [2–5,14,15], the delay is not given explicitly as a function of the state variables but is defined implicitly by an integral equation that needs to be solved together with the remaining equations. The large diversity of SDD equations is such that a general theory is still lacking. This is why simple-looking equations are studied in detail as reference problems. In this work, we focus on one such problem, namely the scalar SDD equation
1.1where *k*>0 is our control parameter and *ε*≪1 is a small parameter. More precisely, we will study the first Hopf bifurcation of the trivial solution *x*=0 as *k* is varied. A general Hopf theorem, stipulating the conditions under which a periodic solution bifurcates from a steady state was only established very recently [16–18]. In this work, we construct explicitly the periodic solution of equation (1.1) in the weakly nonlinear limit of small *x*, taking advantage of the smallness of *ε*. Asymptotic treatments of this kind appear to be lacking in the literature of SDD differential equation and could thus be a useful complement to the rigorous theories that are being developed. Most of this paper is devoted to the simple SDD
1.2The more general case
1.3turns out to require a separate treatment and is treated at a second stage. Recently, a local bifurcation technique was applied to a weakly damped nonlinear oscillator, showing how the functional form of an SDD influences the direction of bifurcation [19]. Equation (1.1) with the linear SDD given by (1.2) has been intensively studied by John Mallet-Paret and Roger Nussbaum during the last 20 years [20–23] in order to develop analytical techniques that exploit the smallness of *ε*. Recently, an extension of this equation to the case of two SDDs was studied by Humphries *et al.* [24].

Taken with a constant delay, equation (1.1) is well documented in the control literature where the stabilization of the system
1.4is investigated with the control law *u*(*t*)=*Kx*(*t*) and for all possible values of *a* and *K* [25]. When the delay is state-dependent, the basic steady state *x*=0 and its stability properties remain unchanged but the problem is now violently nonlinear. Regarding (1.1) in the light of (1.4), *ε* can be interpreted as the ratio between the intrinsic dynamical time scale of the system (when *u*=0) and the delay of the control. Thus, *ε*≪1 stands as a large-delay limit.

Equation (1.1) displays a Hopf bifurcation at some value *k*=*k*_{H}(*ε*) but, as we shall demonstrate, this bifurcation is singular in the limit *ε*→0. Specifically, the small-amplitude expansion of the periodic solution in powers of (*k*−*k*_{H})^{1/2} becomes non-uniform as *ε*→0, making the Hopf bifurcation singular. Following the method of matched asymptotic expansions [26–28], a new asymptotic expansion of the solution where *k*−*k*_{H} is properly scaled with respect to *ε* is then considered as a separate perturbation problem. This approximation is connected to the first Hopf approximation through a process called ‘matching’ so that both solutions overlap. Specific singular Hopf bifurcation problems have been analysed in detail for neuronal and laser bifurcation problems [29–32] and for a constant-delay differential equation problem modelling high-speed machining [33]. A distinctive feature of the present problem from these earlier examples of singular Hopf bifurcation is that the second asymptotic expansion itself rapidly loses its validity with increasing *k*−*k*_{H}. Therefore, a third asymptotic expansion becomes necessary. Let us remark that Hopf bifurcations for slow–fast systems with one fast and two slow variables were studied in detail by Guckenheimer and co-workers [34,35]. In these references, the definition of a singular Hopf bifurcation is geometrical and is related to the observation of fold points of critical manifolds.

Mallet-Paret & Nussbaum [20] showed that the solution of equation (1.1) approaches a sawtooth limiting profile as *ε*→0 with *k*>1. Figure 1 shows a numerical example. The period of the oscillations is close to the value *k*+1, and the leading approximations of the extrema are
1.5A numerical bifurcation diagram of the periodic solutions is shown in figure 2 where the upper and lower straight lines are provided by (1.5). The parabolic curve emerging from *x*=0 in the vicinity of *k*=1 comes from a local Hopf calculation, which we provide here. The bifurcation diagram in figure 2 indicates that the branch of periodic solutions emerges from a Hopf bifurcation and that the extrema of the oscillations gradually increases until they approach the lines and corresponding to the fully developed sawtooth oscillations. A close examination of the oscillations near the Hopf bifurcation reveals that the transition from sinusoidal to sawtooth oscillations is very rapid as one gradually increases *k* from *k*_{H} (figure 3). In this paper, we show that the Hopf bifurcation becomes progressively more vertical near *k*=*k*_{H} as *ε*→0. The oscillations are nearly sinusoidal in the domain *k*−*k*_{H}=*O*(*ε*^{4}) but change from sinusoidal to sawtooth as soon as *k*−*k*_{H}=*O*(*ε*^{2}) for the linear SDD (1.2). On the other hand, sawtooth oscillations appear at *k*−*k*_{H}=*O*(*ε*) for the general SDD (1.3). By contrast to the previously mentioned singular Hopf bifurcation problems, the problem with SDD (1.2) requires the solution of two inner problems that result from two successive singularities in the perturbation expansions. Eventually, the transition from sinusoidal to sawtooth oscillations is captured by Rayleigh’s nonlinear ordinary differential equation. For the general case SDD (1.3), we find that the oscillations are described by van der Pol’s equation.

For nonlinear constant delay equations, the method of linearization is a standard tool in stability and bifurcation studies. But for SDD equations, this approach is delicate because the solution of the system is not differentiable with respect to the SDD [36]. Cooke & Huang [37] first addressed the question of the local stability of a constant equilibrium. Under reasonable assumptions, the formal linearization where the delay is evaluated at steady state does reflect the stability of the steady state of the SDD.

The rest of the paper is organized as follows. In §2, the linear stability analysis is briefly reviewed and some scalings are anticipated. Section 3 contains a first, textbook type, local analysis of the Hopf bifurcation point when the SDD is given by (1.2). Its validity is seen to break down already as *k*−*k*_{H} becomes *O*(*ε*^{4}). This is the limit studied in §4, where a corrected relationship between the oscillation amplitude and *k*−*k*_{H} is derived. However, a second breakdown of the analysis is found, which motivates the investigation of a new limit, *k*−*k*_{H}=*O*(*ε*^{2}) in §5. Section 6 relates the small amplitudes oscillations to those of Rayleigh’s oscillator. In §7, the general SDD (1.3) is studied, eventually leading to a description of oscillations with van der Pol’s equation. Finally, we comment on how different types of solvability conditions sequentially lead to informations on the solutions. We also discuss the interest of exploring singularly perturbed SDD problems.

## 2. Preliminaries

Equation (1.1) admits the steady state *x*=0 and is linearized simply by setting *τ*=1. In this linearized equation, solutions of the form , with real *ω*, exist at the Hopf bifurcation points *k*_{n}(*ε*) and Hopf bifurcation frequencies *ω*_{n}(*ε*) given by
2.1Assuming *k*>0, we have, as *ε*→0,
2.2We note that these bifurcation points become infinitely close as *ε*→0 (figure 4). Let us denote the smallest bifurcation point, *k*_{0}, and corresponding frequency, *ω*_{0}, by *k*_{H} and *ω*_{H}, respectively. They have the following asymptotic expression:
2.3and
2.4It is that bifurcation point which is responsible for the change of stability of the solution *x*=0.

In this paper, we perform a local analysis to study the amplitude and shape of the periodic solutions as a function of *k*−*k*_{H}(*ε*). From the linear stability analysis above, one may anticipate that the presence of the bifurcation points *k*_{1},*k*_{2},…, will be felt as soon as *k*−*k*_{H}=*O*(*ε*^{2}), so that the standard analysis has to be revised at this point. However, the standard analysis turns out to break down much earlier, yielding a first distinguished limit when *k*−*k*_{H}=*O*(*ε*^{4}). Then, a second distinguished limit is found for *k*−*k*_{H}=*O*(*ε*^{2}), where all the frequencies *ω*_{n} appear in the leading order solution. In this way, a scenario for the transition from harmonic to sawtooth oscillation is emerging.

## 3. Regular Hopf bifurcation (*k*−*k*_{0}≪*ε*^{4})

The construction of the periodic solutions of equations (1.1) and (1.2) bifurcating from the trivial solution *x*=0 initially follows a standard pattern [38]. We first introduce the time *s*≡*σt*, in which *x*(*s*) is 2*π*-periodic, and the small amplitude
3.1In the rescaled time, equations (1.1) and (1.2) become
3.2The solution as well as *σ* is then found by using the expansion
3.3
3.4
and
3.5where we have anticipated from Hopf bifurcation theorem that the expansion of *k* and *ω* admits no odd power corrections in *δ*. Inserting (3.3)–(3.5) into equation (1.1) and equating to zero the coefficients of each power of *δ*, we obtain a sequence of linear problems for the unknown functions *x*_{i}(*s*). They are of the form
3.6where the right-hand-side *H*_{m} is a function involving the solutions at the previous orders only. Given from the linear stability analysis that
3.7the solvability conditions of (3.6) is that the right-hand sides *H*_{m} contain no term proportional to . These conditions will provide the unknown coefficients *κ*_{2} and *σ*_{2}.

From the analysis of the three first problems, we find that the final solution has the form
3.8where *c*.*c*. means complex conjugate and *p*_{0}, *p*_{2}, *p*_{3} are constants depending on *ε*. In the limit *ε*→0, they are given by
3.9while *κ*_{2} and *σ*_{2} are given by
3.10Using (3.4) and the expression of *κ*_{2} in (3.10), we determine *δ* as a function of the deviation *k*−*k*_{H}. We obtain
3.11We next note from (3.9) that *p*_{3}∝*ε*^{−1}. As a result, expansion (3.3) loses its asymptotic character as soon as
3.12since the second- and third-order terms become of comparable size. Combining (3.12) and (3.11), the breakdown of the present analysis occurs when
3.13where *α*=*O*(1). Let us examine the limiting behaviour of expression (3.8) when (3.13) is assumed:
3.14This will be used as a matching condition when we develop a new expansion of the solution. The corresponding limiting expression for the frequency. Using (3.5), (3.10), (3.11) and (3.13), we find
3.15indicating that the first correction of the frequency coming from the nonlinear problem is proportional to *ε*^{3}.

In summary, we have found that the standard Hopf asymptotic approximation is valid only if *k*−*k*_{H}≪*ε*^{4}. If *k*−*k*_{H}∼*ε*^{4}*α*, a new expansion of the solution is needed where we use *ε* as the small parameter. As we shall see in §4, the new solution remains sinusoidal but its amplitude changes dramatically from a (*k*−*k*_{H})^{1/2} scaling law to an (*k*−*k*_{H})^{1/4} dependence.

## 4. Singular Hopf bifurcation (*k*−*k*_{H}=*ε*^{4}*α*)

Solution (3.14) suggests to seek a solution in powers of *ε*. We try the expansion
4.1
4.2
and
4.3where *s*≡*σt* and *x*(*s*,*ε*) is 2*π*-periodic in *s*. As suggested by (3.13) and (3.15), we assume that the first corrections of *k* and *σ* that contribute to the nonlinear problem are *ε*^{4} and *ε*^{3}, respectively. Inserting (4.1)–(4.3), we now obtain at each order in *ε* a problem of the form
4.4where *H*_{m}(*s*) is computed from the previous orders. The decomposition of the perturbation problem into a series of iterative maps is typical of delay differential equations with large delays [38–41]. Rewriting the above equation at time *s*−*π* and using the fact that *x*_{m}(*s*−2*π*)=*x*_{m}(*s*), we obtain
4.5Substracting (4.5) from (4.4) then yields the solvability condition for the right-hand side *H*_{m}
4.6The existence of such a solvability condition originates from the fact that the homogeneous problem *w*(*s*)+*w*(*s*−*π*)=0 admits non-trivial 2*π*-periodic solutions. Indeed, any combination of harmonic functions with *j* odd makes the left-hand side of (4.4) vanish. The general solution of the problem (4.4) under the condition (4.6) is given by
4.7Expressions (4.6) and (4.7) can be used to implement an automatic resolution of the perturbation analysis up to any order using a symbolic software such as Mathematica.^{1}

At first order, we have
4.8At second order, the right-hand side is
4.9Next, we obtain
4.10With *H*_{3}(*s*), a first non-trivial solvability condition is obtained, which yields the profile of *W*_{1}(*s*)
4.11At fourth order, a new solvability condition determines *W*_{2}(*s*)
4.12However, in order to prevent secular divergence of *W*_{2}, a secondary solvability condition, this time on equation (4.12), must be satisfied
4.13giving the first nonlinear correction to the frequency. Thus, , giving
4.14Finally, at fifth order, the solvability condition is a differential equation for *W*_{3}
4.15where
4.16As before, the solvability condition on equation (4.15) is that *q*_{1} vanishes, eventually providing the relationship between the solution amplitude and the deviation from the Hopf bifurcation point:
4.17

In summary, we found that
4.18where
4.19and
4.20Our solution is valid provided it matches (3.14) in the limit *A*→0. In this limit, we note from (4.19) that
4.21or equivalently with (3.13),
4.22With (4.22), we verify that expression (4.18) is exactly matching expression (3.14). We next examine the behaviour of solution (4.18) as the amplitude *A* increases. Specifically, we note that the expansion of the solution becomes itself non-uniform as the second term in the coefficient of *ε*^{2} becomes comparable with the leading term in *ε*. This occurs when
4.23Using (4.23), we note from (4.19) that the leading expression for *k*−*k*_{H} now is
4.24which implies that *k*−*k*_{H}=*O*(*ε*^{2}). Introducing
4.25the solution of (4.19) now has the approximation
4.26With (4.26), our expression (4.19) becomes
4.27It suggests to seek a new solution now in powers of *ε*^{1/2}. Finally, using (4.20) and (4.23), we determine the frequency as
4.28

In summary, the analysis of the solution in critical domain (3.13) indicates that the periodic solution remains nearly sinusoidal but that the amplitude of the oscillations, after the initial (*k*−*k*_{H})^{1/2} law, subsequently levels off as (*k*−*k*_{H})^{1/4} (figure 5).

## 5. Singular Hopf bifurcation (*k*−*k*_{H}=*ε*^{2}*β*)

Limiting expressions (4.25), (4.27) and (4.28) motivate the new expansion
5.1
5.2
and
5.3where *s*=*σt* and *x*(*s*) is 2*π*−periodic in *s*. Inserting (5.1)–(5.3) into equation (3.2) and equating to zero the coefficients of each power of *ε*^{1/2}, we obtain a new sequence of linear problems for the unknown functions *x*_{1}, *x*_{2}, *x*_{3},…. They are of the form (4.4) and their resolution follows from (4.7) as in §4. We obtain, successively, for the first four orders,
5.4
5.5
5.6
and
5.7The above right-hand sides are all even in *s* and therefore automatically satisfy solvability condition (4.6). The first non-trivial solvability condition arises at *O*(*ε*^{5/2}) and reads, with *x*_{1}(*s*)=*W*_{1}(*s*),
5.8Hence, the leading-order solution is no more sinusoidal but contains an arbitrary number of harmonics, as anticipated in the Introduction. We have to solve equation (5.8) with the condition
5.9This provides a necessary relationship between *σ*_{4} and *κ*_{4}. To be valid, *x*_{1}(*s*)=*W*_{1}(*s*) must match expression (4.25) in the small amplitude limit. To this end, we solve (5.8) in the small-amplitude limit. Introducing a new small parameter *η*≪1, we thus write
5.10and
5.11where we have omitted all zero contributions for clarity. From the analysis of the problems for *w*_{1}, *w*_{3} and *w*_{5}, and their solvability conditions, equation (5.8) eventually yields
5.12From the expression of *κ*_{4} in (5.12) and (5.2), we determine *η* as
5.13which we recognize as parameter *a* in (4.26). Hence, the small-amplitude limit (5.12) of approximate solution (5.1) matches the large-amplitude limit, (4.27), of previous expansion (4.1).

We have thus found that the oscillations are changing from nearly linear to nonlinear in the critical domain (4.25). In §6, we analyse these nonlinear oscillations.

## 6. Rayleigh’s equation

Equation (5.8) is nonlinear and is the final problem of our singular Hopf bifurcation analysis. This equation can be reformulated in a more meaningful form by introducing the new variables (*σ*_{4}<0)
6.1and
6.2Indeed, in terms of *y* and *T*, equation (5.8) can be recognized as Rayleigh’s equation [26,42]
6.3where prime means differentiation with respect to time *T*. The only parameter left, *Λ*, is defined by
6.4For any initial condition, the solution approaches a stable limit-cycle of period *P*. Note that the change of sign in *T* is necessary in order to have a stable limit-cycle of equation (6.3). Given *Λ*, we determine the period *P* numerically. The corresponding period in *s* must be 2*π* which implies, using (6.2), that
6.5From this, we extract *κ*_{4} as
6.6Combining (6.6) with (6.4) and (6.1), we further obtain
6.7Of particular interest is the asymptotic behaviour of the solution as , corresponding to in (6.3). The solution to that problem is well documented (see, for instance [26,42]) and gives, for the leading expressions of *P* and the extrema of *y*
6.8Using (6.6) and (6.7), the corresponding values for *κ*_{4} and *x*_{1} are
6.9The resulting analytical bifurcation diagram and time traces are compared to numerical simulations in figures 6 and 7. A good quantitative agreement is seen as long as *ε* is sufficiently small and *k* is close to 1. As *k* increases, however, the shape of the oscillations differ more and more between the original model and the asymptotic reduction. Indeed, the descending part of the oscillations become more and more vertical in the original model (1.1) (figure 1) which the Rayleigh oscillator cannot reproduce.

## 7. More general state-dependent delays

We now consider general SDDs. Expanding *τ*(*x*) in Taylor series, we may write, in all generality
7.1provided that *τ*′(0)>0. Given that *x* is small, higher order terms in (7.1) can be neglected in this analysis. It would seem at first that the quadratic term above should lead only to a small correction to the preceding results. Its influence on the dynamics, however, turns out to be predominant. Indeed, performing the same analysis as in §3, that is with an asymptotic development of the form (3.3)–(3.5), the expressions in (3.10) eventually become
7.2Hence, in the small-*ε* limit, the direction of bifurcation, given by the sign of *κ*_{2}, is entirely dictated by *b* (assuming that it is *O*(1)). In particular, *b* must be negative for the bifurcation to be supercritical. A similar observation was made by Mitchell & Carr [19] in their study of a weakly damped nonlinear oscillator with SDD. On the other hand, (3.11) is now
7.3while the solution has the asymptotic form
7.4Thus, the second- and third-order terms become of the same magnitude as soon as *δ*=*O*(*ε*^{2}). On the other hand, the third and fundamental harmonics balance when *δ*=*O*(*ε*). Let us focus directly on the latter case. The analysis of the singular Hopf bifurcation problem now requires an expansion of *x*(*s*) in power of *ε*
7.5where the *x*_{j}(*s*) are 2*π*−periodic functions of *s* defined by
7.6We also need to expand the bifurcation parameter as
7.7Substituting this development in (1.1) with *τ* given by (7.1), the first two problems arising at *O*(*ε*) and *O*(*ε*^{2}) are exactly as in §4. At the next order, the solvability condition reads
7.8Together with the periodicity condition *x*_{1}(*s*+2*π*)=*x*_{1}(*s*). Equation (7.8) is van der Pol’s equation. Indeed, with and , equation (7.8) becomes
7.9where prime means differentiation with respect to time *T*. Following the same reasoning as before, a solution of period *P* in *T* must be of period 2*π* in *s*, giving *κ*_{2}=(*P*^{2}−4*π*^{2})/8 and *σ*_{2}=*πPΛ*/4. A comparison between numerical simulations of equation (1.1) and its van der Pol asymptotic reduction is given in figure 8 for *ε*=0.05. The similarity of the shape of oscillations is striking. This time, however, and contrary to the comparison made in figure 7, the value of *k* in (1.1) needs some adjustment for a good fit. Better quantitative agreement is expected for smaller values of *ε*.

## 8. Discussion

In this paper, we have applied singular perturbation techniques to investigate how slow/fast oscillations are emerging from a Hopf bifurcation in a class of singularly perturbed SDD equations. In the same way as asymptotic methods have been tested for specific constant-delay differential equation problems [38], we concentrated here on a class of SDD equations and examined how these methods can be applied. In particular, we showed that they allow us to construct approximate solutions of physical significance. They are thus valuable complements to general theorems as well as precious guides for numerical simulations. In this work, we constructed the periodic solutions by adapting the Lindstedt–Poincaré method. Recently, the method of multiple scale was applied for the first time to an SDD turning problem [43].

When a perturbation problem is singular, a fundamental piece of information regarding the solution is usually lacking in the first few orders of the analysis. That piece of information is then deduced by consistency checks—solvability conditions—at higher orders. The problem we studied in this paper distinguishes itself by the form of the solvability conditions encountered. In §4, it is through a solvability condition, equation (4.11), that the sinusoidal shape of the leading-order oscillations is determined, whereas this would normally be known from the linear stability analysis already. The amplitude and frequency of the oscillation are subsequently derived from secondary solvability conditions, equations (4.13) and (4.17), in a rather unusual way. Further away from the Hopf bifurcation in §§5 and 7, the solvability conditions directly yields nonlinear differential equations, equations (6.3) and (7.9) for the amplitude and shape of the oscillations. These equations are the final objective of our analysis because they allow us to describe through a simple ordinary differential equation how the oscillations gradually change from sinusoidal to sawtooth in the vicinity of the Hopf bifurcation point. For the two cases we have investigated, we found Rayleigh’s equation for the linear SDD (1.2), and van der Pol’s equation for the more general SDD (1.3), containing a quadratic term. The latter is a surprising result. On account of the small amplitude of the solution, we were indeed expecting that the introduction of quadratic terms in (1.3) would only lead to a modification of Rayleigh’s equation. However, we found that this quadratic term controls the direction of bifurcation (supercritical versus subcritical) and the shape of the resulting oscillations.

A large delay is the reason why there is a small parameter *ε* multiplying the left-hand side of equation (1.1). It is this small parameter which is responsible for the singular character of the Hopf bifurcation. Already at the level of linearization, the smallness of *ε* (the largeness of the delay) is responsible for an accumulation of Hopf bifurcation points near the primary instability. This favours the emergence of higher harmonics very soon after the instability threshold, and hence, a rapid transition from sinusoidal to sawtooth oscillations.

Another case of a singular Hopf bifurcation, not treated here, arises when the frequency of the oscillations is small. We note from the linearized theory of equation (1.4), given in [25] that the Hopf bifurcation frequency tends to zero as . This situation was encountered by Stumpf [44] in his analysis of an SDD model for currency exchange fluctuations. For the constant delay case, there is a Hopf bifurcation exhibiting an infinite period and the local bifurcation diagram was captured by treating the problem as a singular perturbation problem [45]. We expect that the same technique could be used for the SDD case.

## Funding statement

This research was supported by the F.R.S.-FNRS (Belgium) and by the Interuniversity Attraction Poles program of the Belgian Science Policy Office under grant no. IAP P7-35.

## Acknowledgements

T.E. thanks John Mallet-Paret and Roger Nussbaum for many stimulating discussions.

## Footnotes

↵1 See Mathematica file in the electronic supplementary material.

- Received September 5, 2013.
- Accepted November 26, 2013.

- © 2013 The Author(s) Published by the Royal Society. All rights reserved.