A Lazarus

Discrete dynamical stabilization of a naturally diverging mass in a harmonically time-varying potential

Keywords: Dynamical systems, Stability analysis, Modal analysis, Floquet theory, Standing waves

We numerically investigate the stability and linear oscillatory behavior of a naturally diverging mass whose potential energy is harmonically modulated. It is known that in the Kapitza limit, i.e. when the period of modulation is much smaller than the diverging time, the collapsing mass can be dynamically stabilized and behave like an effective classic harmonic oscillator. We find that in the regime where the period of modulation is larger than the collapsing time of the mass, dynamical stabilization is still possible but in a discrete fashion. Only almost-periodic vibrational modes, or Floquet forms (FFs), are allowed that are located in independent stability stripes in the modulation parameter space. Reducing the FFs to their periodic eigenfunctions, one can transform the original equation of motion to a dimensionless Schrödinger stationary wave equation with a harmonic potential. This transformation allows for an analytical prediction of the stability stripes and the modal shapes of the vibrating mass. These results shed new light on the stability of linear dynamical systems, analytical solutions of Mathieu equations and on the relations between Initial and Boundary Value Problems.

Introduction

Modal analysis is a linear perturbation method that allows to characterize the local oscillatory and stability behavior of stationary states of dynamical Email address: arnaud.lazarus@sorbonne-universite.fr (A. Lazarus) systems [START_REF] Guckenheimer | Nonlinear oscillations, dynamical systems, and bifurcations of vector fields[END_REF][START_REF] Strogatz | Nonlinear dynamics and chaos: with applications to physics, biology and chemistry[END_REF]. It is used in various area of physics, from molecular vibrational frequencies [START_REF] Nakamoto | Infrared and Raman spectra of inorganic and coordination compounds Part A[END_REF] to the stability of engineered structures [START_REF] Nayfeh | Linear and nonlinear structural mechanics[END_REF]. Reduced to a single dimension in space, this concept is modeled by the archetypal example of a mass moving in a quadratic potential energy whose governing equation is a linear homogenous Ordinary Differential Equation (ODE) with initial conditions. In the case of a perturbed equilibrium, i.e. for a constant potential in time, two qualitative behaviors exist: the mass is either neutrally stable, harmonically oscillating in a potential well (this case is the classic harmonic oscillator) or unstable, exponentially diverging on a potential hill.

A less constrained situation eventually occurs when the potential energy of the perturbed stationary state is free to periodically vary with time [START_REF] Bergé | Order within chaos[END_REF][START_REF] Bolotin | The dynamic stability of elastic systems[END_REF], i.e. when the motion of the mass is mathematically governed by a linear homogenous ODE with periodically time-varying coefficients. This generalized framework explains a broader class of physical problems from parametric oscillators [START_REF] Turner | [END_REF]8] to the emergence of Faraday waves [9,10,11] or the motion of the lunar perigee [12,13]. Again, two different scenarios should be considered whether the mass is naturally diverging or not. In a periodically time-varying potential well, the quasi-periodically oscillating particle eventually destabilizes for certain regions in the modulation parameters space, a.k.a. Mathieu's tongues [14,15]. On a modulated potential hill, there is an asymptotic limit for which the modulation parameters allow to stabilize the naturally diverging mass [16,17]. This last model is generally called a Kapitza oscillator since Kapitza was the first to explain the dynamical stabilization of a collapsing inverted pendulum in which the pivot point vibrates in the vertical direction [18,19].

Although simplistic, the aforementioned 1D oscillators allow for the modeling and understanding of an outstanding number of physical systems and their complete theoretical description is therefore a necessity. Yet, to the best of our knowledge, we believe the linear vibratory analysis of a 1D Kapitzalike oscillator has only been partially studied, so far. Here, we present a numerical and theoretical study of the modal analysis of a naturally diverging mass on a harmonically modulated potential hill. Unlike focusing on the Kapitza asymptotic regime where the modulation periodicity is much smaller than the collapsing time of the particle, we are interested in the regime where the modulation periodicity is larger than the diverging time. In Section 2, we present the system under study and the numerical tools we use to perform its modal and stability analysis. Section 3 highlights the original mechanical behavior numerically observed when the modulation is slow as compared to the natural diverging time. Section 4 shows the theoretical reduction process that allows to analytically characterize the numerical results obtained in Section 3. We conclude and discuss those results in Section 5.

We find that, opposite to the Kapitza limit, the dynamical stabilization of the naturally diverging mass consists in almost-periodic vibrational modes, or Floquet Forms (FFs), that are located in independent stability stripes of the modulation parameter space. Based on simple numerical observations, we show it is possible to reduce the harmonically time-varying governing differential equation under study to a dimensionless Schrödinger stationary wave equation with a harmonic potential. The latter can be analytically solved and is used to predict the location of the stability stripes as well as the periodic modal shapes of the stable FFs of the modulated oscillator. These results should help improving our understanding of periodically time varying systems and the relation between initial and boundary value problems, two notions that are notably at the heart of the recently unveiled discrete time crystals [20,21,22,23]. This work could also benefit to the modeling of hydraulic quantum analogs that consist in a bouncing droplet "walking" on a harmonically vertically vibrated bath [24,25,26,27,28,[START_REF] Fort | Proc. Natl. Acad. Sci. USA[END_REF][START_REF] Perrard | [END_REF]. Although we believe the mathematics leading to Schrödinger's like equations are fundamentally different, it is interesting to notice that both dynamical systems involve a naturally collapsing or bouncing mass in a time-varying potential energy and lead to discrete mechanisms.

(a) (b)

EQUATION OF MOTION AND FLOQUET THEORY

Figs. 1(a),(b) show the 1D periodically time-varying linear dynamical system under study. A particle of mass m locally moves with a position r(τ ) on a harmonically modulated, quadratic potential hill, so that its total kinetic plus potential energy reads E(τ ) = 1/2m ṙ(τ ) 2 + 1/2k(t)r(τ ) 2 where ṙ(t) = dr(τ )/dτ is the velocity of the particle. The linear equation of motion of m associated with E(τ ) is modeled by the 1 degree-of-freedom mass-spring system shown in Fig. 1(b). In this periodically conservative system, the particle experiences a parametric excitation F (τ ) = -k(τ )r(τ ) that derives from the modulated potential energy and that can be represented by a spring with a T -periodic varying stiffness k(τ ) = k(τ + T ) = -k 0 + k 1 cos(2πτ /T ) where -k 0 and k 1 are the fundamental and modulated stiffness, respectively. The only difference with a classic parametric oscillator is that the fundamental stiffness -k 0 is negative: if k 1 = 0 N/m, the system is linearly unstable and the mass exponentially diverges following r(τ ) ∝ e 2πτ /T 0 where T 0 = 2π/ k 0 /m is the natural diverging time.

According to Newton's second law and Lagrangian mechanics, the dimensionless Initial Value Problem (IVP) governing r(τ ) reads

d 2 r(t) dt 2 -r(t) + α cos(2πt/T )r(t) = 0 ( 1
)
where t = 2πτ /T 0 is the dimensionless time. In this fundamental Mathieu equation [START_REF] Whittaker | A course of modern analysis[END_REF][START_REF] Magnus | Hill's equation[END_REF], the natural diverging time is now 2π and the dimensionless modulation period of the potential energy is T = 2πT /T 0 . The two relevant modulation parameters are the time scale ratio T /2π and the stiffness ratio α = k 1 /k 0 . For modal and stability analysis of the linear differential Eq.( 1), one can use Floquet theory [START_REF] Floquet | [END_REF][START_REF] Whittaker | A course of modern analysis[END_REF] to express r(t) as a linear combination of two almost-periodic vibrational modes, or Floquet forms,

r(t) = c 1 Ψ(t)e st + c 2 Ψ * (t)e -st (2)
where c 1 and c 2 are constants determined upon initial position r(0) and velocity ṙ(0). Replacing the Floquet form Ψ(t)e st in Eq.( 1) leads to an eigenvalue problem that can be numerically solved for each set of parameters (α, T ) [34,35,36]. According to Floquet theory, s is a complex eigenvalue called a Floquet exponent and the computed eigenfunction Ψ(t) and its complex conjugate Ψ * (t) are periodic with a dimensionless modulation period T . Because of the dimension and the linearity of Eq.(1), only three qualitative stability behaviors can be observed in the (α, T ) space:

• If (s) = 0 and 0 ≤ (s) ≤ π/T , the particle is neutrally stable and r(t) is an almost-periodic oscillation about r(t) = 0.

• If (s) > 0 and (s) = 0, the position r(t) = 0 is dynamically unstable and r(t) is a T -periodic motion that exponentially diverges from r(t) = 0 with a growth rate (s).

• If (s) > 0 and (s) = π/T , the position r(t) = 0 is dynamically unstable and r(t) is a 2T -periodic motion that exponentially diverges from r(t) = 0 with a growth rate (s).

MODAL AND STABILITY ANALYSIS (T > 2π)

The oscillating system represented by the governing equation Eq.( 1) has been well studied in the asymptotic limit where the modulation period T is much smaller than the natural diverging time T 0 , i.e. for T << 2π [START_REF] Magnus | Hill's equation[END_REF]. In this regime, first understood by Kapitza thanks to averaging techniques when studying the inverted pendulum whose pivot point is vertically vibrated [18,16], it is possible to stabilize the naturally diverging mass if (αT ) 2 /8π 2 > 1. Here, we will focus on the dynamical behavior of the modulated oscillator of Fig. 1 for T > 2π that is the overlooked regime we are interested in.

We perform a numerical stability analysis of Eq.( 1) by analyzing the growth rate (s) and fundamental frequency (s) of the computed Floquet forms Ψ(t)e st and Ψ * (t)e -st in the modulation parameter space (α, T ). Fig. 2 shows the linear stability chart of the particle for 0 ≤ α ≤ 10 and 0. T /2π ≤ 5. Following the stability classification given at the end of Section 2, we find an alternation of T -unstable (red dots in Fig. 2) and 2T -unstable regions (yellow dots) which is typical of classic Mathieu equations [START_REF] Whittaker | A course of modern analysis[END_REF][START_REF] Magnus | Hill's equation[END_REF]. If the system is not modulated enough, i.e. for α < 1 (the potential energy is never positive), the mass cannot be dynamically stabilized. For T < 2π, a stability region opens whose instability lower limit corresponds to the classic frontier of Kapitza, α = √ 8π/T , shown in black line in Fig. 2. For T > 2π, i.e. when the period of modulation is slower than the naturally diverging time, Kapitza's theory does not hold but the dynamical stabilization of the mass is still possible in the (α, T ) modulation space: our goal is to get some physical and mathematical insights into this stabilization process.

According to Floquet theory, the T and 2T instability regions cannot merge in the (α, T ) space so each alternation of colors indicates a tiny stability domain that we denote l = 0, 1, 2, . . . as shown in Fig. 2 [37]. Those regions of stability form independent "stripes", whose finite width drastically decreases as T increases (we barely reach inside the stability regions above T /2π > 5 because of machine epsilon of our computational software). Figs. 3(a) and (b) show the evolution of the spectrum of the two Floquet forms Ψ(t)e st and Ψ * (t)e -st as a function of T /2π for α = 2.25. According to Floquet theory, the eigenfunction Ψ(t) is T -periodic with a fundamental frequency 2π/T . Therefore, by Fourier expansion, the spectrum of the two Floquet forms read s + h i2πh/T and -s + h i2πh/T , where h is an integer. By the superposition principle of Eq.( 2), the solution r(t) contains the sum of both spectrum in it. Fig. 3(a) shows the location of the frequency spectrum of the two vibrational modes ± (s)+ h h as a function of T /2π for α = 2.25 where we introduced (s) = (s)T /2π (only the positive part of the spectrum is shown as the latter is symmetric with respect to the x-axis). For most modulation period T , the frequency spectra of the two modes are locked, alternatively in 0 + h h2π/T (T -periodic spectrum) or in π/T + h h2π/T (2T -periodic spectrum). As shown in Fig. 3(b), this lock-in is associated with a finite growth rate (s) so that one of the periodically oscillating Floquet form is exponentially diverging ( (s) > 0) and one is damped ( (s) < 0). Between those locked unstable regions, tiny zones exist where the two Floquet modes are unlocked and where the associated growth rate (s) is zero for both modes, a situation that leads to the stability stripes in Fig. 2 which correspond to neutrally stable almost periodic solutions r(t). The insets in Figs. 3(a) and (b) display a zoom on the first stability region l = 0. In stability regions, the fundamental frequency (s) evolves continuously from 0 to π/T for tiny variations of T . Note that the green and blue vertical lines are the limits of the stability zones corresponding to neutrally stable T and 2T -periodic solutions with (s) = 0, respectively, and located at modulation periods T which are characteristic values of the so-called periodic Mathieu functions [START_REF] Whittaker | A course of modern analysis[END_REF].

We now focus on the neutrally stable vibrational motions of the mass. Figs. 4(a),(b) and (c) display three typical examples of such vibrational modes over three periods T for modulation parameters in the stability regions l = 0, 1 and 2, respectively [37]. In each panel, we show the position r(t) = Ψ(t)e i (s)t , velocity ṙ(t) and mechanical energy E(r(t)) = 1/2 ṙ(t) 2 + 1/2(-1 + α cos(2πt/T ))r(t) 2 of the computed vibrational modes. The modes are normalized so that N = T Ψ * (t)Ψ(t)dt = 1. Full lines represent the computed almost-periodic motion when dotted lines show the periodic eigenfunctions Ψ(t), Ψ(t) and E(Ψ(t)). The stable almost-periodic vibrations of the mass described in Fig. 4 are synchronized with the modulated potential energy and can be decomposed in a successive repetition of similar motions that are scaled copies of their T -periodic eigenfunctions (the scaling factors can take all the values between -1 and 1 and will be discussed in next section). The eigenfunctions can themselves be decomposed in three parts over one period T : i) when the curvature of the quadratic potential energy in r(t) = 0, κ(t) = (-1 + α cos(2πt/T )), is negative, the mass exponentially diverges, ii) when κ(t) > 0 the mass oscillates, iii) when κ(t) < 0 the mass exponentially converges to a state very almost identical to the previous period, close to r(t) = 0 [37]. The qualitative difference between vibrational motions with ascending order of l from Fig. 4(a) to (c) is that in the l th stability region, the mass is able to do (l + 1) oscillations during the time κ(t) is positive. Like in a classical oscillator, the velocity is zero when the position is at a local maximum. Finally, the mechanical energy E(r(t)) illustrates the "diverging/vibrating" binary behavior of the studied oscillator: over one period T , it seems null when κ(t) < 0 and positive when κ(t) > 0.

REDUCTION OF THE DYNAMICS TO Ψ(t)

4.1. Reduction to the neutrally stable periodic solutions r l (t) = Ψ l (t)

As already mentioned and illustrated in Fig. 4, the almost-periodic neutrally stable modes of the l th stability region, r l (t) = Ψ l (t)e i (s)t , can be decomposed in a succession of cycles that are scaled copies of the periodic eigenfunctions Ψ l (t) = Ψ l (t + T ). To highlight this time-translational property, we superpose, on the primitive time cell -T /2 < t ≤ T /2, the 20 first consecutive periods of the position r(t) of some computed vibrational modes. The result is shown with T /2π = 5 and α = 1.14905152323442 for modes in the l = 0 stability region in Fig. 5(a) where the chronological order is indicated by 20 gradual colors from cold blue to hot red. The dotted line represents the periodic eigenfunction Ψ 0 (t). Around T /2π = 5, the width of the stability region has already decreased to ≈ 10 -13 in the (α, T ) parameter space. Throughout this width , the eigenvalue or fundamental frequency (s) varies in the primitive spectral cell 0 ≤ (s) ≤ π/T (see insets of Fig. 3(a) and (b) to see the evolution of s in a stability region); unlike the periodic eigenfunction Ψ(t) that remains unchanged. As a consequence, unlike its periodic eigenfunction Ψ(t), the position of the stable vibrational mode, r(t) = Ψ(t)e i (s)t , is sensitive to the width of the stability regions (the same is true for ṙ(t) and Ψ(t) or E(r(t)) and E(Ψ(t))). This sensitivity property is highlighted in the top and bottom panels of Fig. 5(a) that show vibrational modes in the l = 0 stability region for T /2π = 5 ± /10, corresponding to (s) ≈ 0.98 × π/T and (s) ≈ 0.62 × π/T , respectively. From one panel to another, Ψ 0 (t) as well as the various curves of r(t) have similar shapes, the only change lies in the chronological order in which r(t) appears. The same numerical observations on sensitivity can be made for vibrational modes in higher regions of stability as illustrated for l = 1 in Fig. 5(b), where the top and bottom panel correspond to a fundamental frequency (s) ≈ 0.25 × π/T and (s) ≈ 0.04 × π/T .

To go further, one needs to specify the scaling factors that relate the position r(t) = Ψ(t)e i (s)t to the eigenfunction Ψ(t) for each successive periods. According to Floquet theory [START_REF] Whittaker | A course of modern analysis[END_REF][START_REF] Magnus | Hill's equation[END_REF], the almost-periodic solution r(t) theoretically verifies the periodic mapping

r(nT ) = cos( (s)nT ) × r(0) = cos( (s)nT ) × Ψ(0) ( 3
)
where n is a positive integer and T is the dimensionless modulation period. What we observe in the numerical results illustrated in Fig. 5 is stronger than Eq.( 3) since it is the whole functions r([nT -T, nT ]), i.e. r(t) for t ∈ [nT -T, nT ], that are related to the unique T -periodic function Ψ([nT -T, nT ]) = Ψ([0, T ]). From numerical observations of r(t), we are able to generalize the quasi-periodicity property of Eq.( 3) from a single time value t to a full period [t, t + T ] so that

r([nT -T, nT ]) = ρ n × Ψ([0, T ]) ( 4
)
where ρ n = cos( (s)(nT -T )) is a real between -1 and 1. In Fig. 6(a r l ([nT -T, nT ]) all collapse on their respective eigenfunctions Ψ l ([0, T ]). If we increase the number of periods, the scaling factors ρ n take a very large number of different values between -1 and 1 and all the r l ([nT -T, nT ]) would collapse on Ψ l ([0, T ]).

The property of Eq.( 4), altogether with the aforementioned sensitivity property of thin regions of stability when T > 2π, offers a reduction opportunity. One could choose to neglect the width of the l th stability regions and reduce the latter to widthless stripes, i.e. to represent a neutrally stable mode of vibration r l (t) = Ψ l (t)e i (s)t by its sole T -periodic eigenfunction Ψ l (t) in the primitive periodic cell -T /2 < t ≤ T /2 by taking (s) = 0. By doing so, we loss informations on the vibrational modes r l (t) whose Floquet exponents actually span 0 ≤ (s) ≤ π/T , but we will gain analytical insights on the oscillator governed by Eq.( 1), notably on Ψ l (t) and their location in the (α, T ) modulation parameter space.

Reduction to a dimensionless Schrödinger stationary wave equation with

a harmonic potential Reducing the regions of stable vibrational modes to the limit branches of periodic Mathieu functions r l (t) = Ψ l (t) (that correspond to the stability frontiers s = 0 of the red regions in Fig. 2), Ψ l (t) are solutions of the initial value problem (IVP), Eq.( 1), so that we can write

d 2 Ψ(t) dt 2 -Ψ(t) + α cos(2πt/T )Ψ(t) = 0. (5)
Unlike the general Eq.(1), Eq.( 5) is periodic and can therefore be reduced to a single representative period T . The complete mathematical analysis of Eq.( 5) is beyond the scope of this paper but one could first resolve it thanks to numerical observations. By construction and for sufficiently thin stability regions with << 1, the normalized periodic solutions Ψ l (t) has a compact support over a period T and can be decomposed in two exponential functions, mirrored with respect to a central vertical axis, connected by a function with l + 1 extrema (see throughout Figs. 23456). Consequently, the intensity of the normalized Ψ l (t) is localized in the center of the period and almost zero elsewhere (this is especially pronounced as the width of the stability region decreases). So, multiplying Ψ l (t) by a T -periodic function F (t) would localized the intensity of F (t) at the center of the period. Choosing the primitive periodic range -T /2 < t ≤ T /2 as the representative cell for the whole solution Ψ(t), Ψ(t) cos(2πt/T ) would cancel the intensity of cos(2πt/T ) away from the origin t = 0 so that, assuming Ψ(t) is localized enough, the Taylor series approximation

Ψ(t) cos(2πt/T ) ≈ Ψ(t)(1 -2π 2 t 2 /T 2 ) (6)
would be legitimate. Replacing Eq.( 6) in Eq.( 5) and upon the change of variable χ = vt with v = (2π 2 α/T 2 ) 1/4 , Eq.( 5) governing Ψ l (t) can be rewritten in the form of a Boundary Value Problem (BVP) on the dual primitive cell 1) (green lines) are compared to the analytical results from the BVP Eq.( 7) (black lines). Each solution E l = l + 1/2 of Eq.( 8) corresponds to a branch of stability, s = 0. (b) In blue line, absolute square of the first four computed eigenfunctions Ψ l (t) of Eq.( 1) on the dual primitive cell -vT /2 < χ ≤ vT /2 for T /2π = 5 and v = (2π 2 α/T 2 ) 1/4 = 0.389351968371210, 0.416745399049464, 0.444673895477715, 0.472863276327510 from bottom to top. The eigenfunctions are normalized so that N = T |Ψ(t)| 2 dt = 1. Black lines are analytically obtained from Eq.( 9).

-vT /2 < χ ≤ vT /2

EΨ(χ) = 1 2 χ 2 Ψ(χ) - 1 2
d 2 Ψ(χ) dχ 2 (7)
where E = (α -1) / √ 2α × T /2π. The linear eigenvalue problem with variable coefficient in Eq.( 7) is well-known as it is the dimensionless form of a stationary Schrödinger equation with a harmonic potential, usually predicting the total energy E and wavefunction Ψ(χ) of a 1D quantum harmonic oscillator [START_REF] Messiah | Quantum Mechanics[END_REF]. For Ψ(-∞) = Ψ(+∞) = 0, the discrete set of eigenvalues E and eigenfunctions Ψ(χ) take the form

E l = (α -1) / √ 2α × T /2π = l + 1/2 (8)
and

Ψ l (χ) = H l (χ)e (-χ 2 /2) /(π 1/4 √ 2 l l!) (9)
where H l (χ) are Hermite polynomials and l = 0, 1, 2, . . . The analytical results of Eqs.( 8)-( 9) allow us to predict the shape of the T -periodic Mathieu functions Ψ l (t) as well as their location in the modulation parameter space (α, T ). Fig. 7(a) shows that each analytical E l of Eq.( 8) accurately predict the T -periodic limits of the l th numerical stability regions of Eq.(1) in the (α, T ) space. As expected from numerical observations, since the width of the l th stability region and the errors in the approximation Eq.( 6) decrease as T increases, the analytical prediction Eq.( 8) for each mode l becomes better with T as shown in Fig. 7(a). Fig. 7(b) illustrates, for T /2π = 5, that the analytical eigenfunctions Ψ l (χ) of Eq.( 9) are in good agreement with the periodic eigenfunctions Ψ l (t) that were numerically computed from Eq.( 1) and shown in Figs. 4(a), (b) and (c) for a α in the l = 0, l = 1 and l = 2 stability regions, respectively. Thus, in a certain asymptotic limit that will need to be rigorously defined in future work, the stationary wave equation Eq.( 7) allows to compute, in the (α, T ) modulation parameter space, the eigenfunctions Ψ l (t) of the stable vibrational motion of the oscillator. The actual position of the mass of the oscillator, r l (t) is not predicted by Eqs.( 8)-( 9), but it can be infer from the previous subsection that it will be a successive repetition of scaled periodic eigenfunctions Ψ l (t).

Conclusions and discussions

In summary, we presented a numerical and theoretical study of a fundamental vibratory mechanism: the discrete dynamical stabilization of a har-monically modulated 1D linear oscillator whose mass is naturally diverging. We have shown an intriguing asymptotic limit of a Kapitza like oscillator governed by a Mathieu equation in which the stable vibrational modes, or Floquet forms, get restricted to thin stability stripes of the modulation space. The multiple initial value problems (IVP) governing the motion of the neutrally stable particle in those modulation regions can be recast in a single boundary value problem (BVP) on a primitive periodic cell that is a dimensionless Schrödinger stationary wave equation with a harmonic potential. This stationary wave equation allows to analytically predict the modulation parameters leading to dynamical stabilization as well as the periodic modal shapes (which are the periodic Mathieu functions of the system) of the stable vibrational modes, albeit with a certain ignorance on the actual motion of the mass.

This paper is a first study, mainly numerical, that helps to shed new light on periodically time-varying oscillators governed by Floquet theory that are at the heart of many problems in physics. In future work, a more rigorous mathematical approach would be needed, notably to characterize the asymptotic limit where the original equation of motion can be approximated by a dimensionless Schrödinger's equation. Another mathematical task would be to determine whether the presented reduction from IVP to BVP is unique, or generalizable to a broader class of ordinary differential equations with periodic coefficients. Meanwhile, this theoretical work raises a practical question: is it possible to experimentally validate the described mechanism? A trivial answer would be to set up an experimental system whose linearized equation of motion is in the form of Eq.( 1). The simplest realization could be an inverted pendulum whose pivot point is vertically vibrated which is the archetypal example of a mass that naturally diverges under periodically modulated gravity. The major problem with such a basic experiment would lie in the theoretical width of the stability regions that we uncovered. The latter are intrinsically so small that any additional experimental noise would probably made the pendulum unstable. A possibility could be to add some damping to the oscillator that would decrease the growth rate (s) and widen the instability regions for any given modulation parameters. The drawback of his approach is that the eigenfunctions Ψ(t) as well as the actual positions r(t) would loose their symmetry in the primitive periodic cell. Another solution would be to start the experiment from modulation parameters for which the stability regions are reasonably large and strive to progressively explore the regions as it shrink. In any case, a study of the dependence of an additional noise on the stability of the motion in a nonlinear framework would be needed to get some physical insights in the transitions between stable regions.

Last but not least are the mathematical analogies of the presented harmonicallyvarying 1D oscillator with the 1D quantum harmonic oscillator [START_REF] Messiah | Quantum Mechanics[END_REF]. Quantum analogs have regained considerable attention due to the recent experimental discovery of a mechanism that consists in a bouncing fluid droplet "walking" on a vertically vibrated bath [24,25,26,27]. By varying the experimental setup, this new framework coined "dual walkers" has already allowed to describe quantum tunneling [28], quantization of classical orbits [START_REF] Fort | Proc. Natl. Acad. Sci. USA[END_REF] and the quantum harmonic oscillator [START_REF] Perrard | [END_REF] in an analogous fashion. Although powerful, the mathematical modeling underlying this hydraulic quantum analog is based on nonlinear interactions between a wave and a particle such as chaotic behavior and time-delay phenomena, which are rather elaborate ingredients when Schrödinger's equation could reduce to a 1D linear ODE in the case of the quantum harmonic oscillator [START_REF] Messiah | Quantum Mechanics[END_REF]. Since the existence of dual walkers is deeply related to Faraday instabilities that are parametric instabilities explained by Floquet theory, it is a possibility that the presented periodically time-varying oscillator framework could help in the physical understanding and modeling of quantization phenomena observed in macroscopic wave-particle interactions.

Figure 1 :

 1 Figure 1: The periodically time-varying oscillator under study. (a) Mass on a local quadratic potential hill, harmonically modulated over a period T . (b) Associated linear mass-spring system with a time-varying stiffness k(τ ) = -k 0 + k 1 cos(2πτ /T ).

Figure 2 :

 2 Figure 2: Numerical stability chart of the harmonically modulated collapsing mass for 0 ≤ α ≤ 10 and 0.2 ≤ T /2π ≤ 5. Red and yellow regions show T and 2T -periodic unstable solutions, respectively. Between those regions are tiny zones, denoted by l = 0, 1, 2, . . . in ascending order of T , where r(t) is a neutrally stable almost-periodic solution. The black line represents the Kapitza limit α = √ 8π/T .

2 ≤

 2

Figure 3 :

 3 Figure 3: Spectrum of the solution r(t) as a function of T /2π for α = 2.25. (a) Location of the frequency spectrum (s) + h. Green and blue vertical lines show the neutrally stable T and 2T -periodic solutions, respectively. Inset: Zoom on the first stability region l = 0 where the spectrum of the two Floquet forms are unlocked. (b) Evolution of the growth rate of the two vibrational modes Ψ(t)e st and Ψ * (t)e -st . Inset: Zoom on the first stability region l = 0 where the growth rate of both modes is (s) = 0.

Figure 4 :

 4 Figure 4: Time evolution position r(t), velocity ṙ(t) and energy E(r(t)) of neutrally stable vibrational modes Ψ(t)e i (s)t with N = T |Ψ(t)| 2 dt = 1. Periodic eigenfunctions Ψ(t), Ψ(t) and E(Ψ(t)) are shown in dotted lines. T -periodic cells are delimited by thin vertical lines. We fixed T /2π = 5. (a) Mode in the l = 0 stability regions for α = 1.14905152323442. (b) Mode in the l = 1 stability regions for α = 1.5081802859941. (c) Mode in the l = 2 stability regions for α = 1.9549540063305.

Figure 5 :

 5 Figure 5: Sensitivity analysis of the stable vibrational modes r(t) = Ψ(t)e i (s)t for T /2π = 5. (a) For a mode in the l = 0 stability region (α = 1.14905152323442), the 20 first consecutive periods of r(t), chronologically ordered from cold blue to hot red, are superposed on the primitive time cell -T /2 < t ≤ T /2. Dotted lines represent the periodic eigenfunction Ψ 0 (t). Top and bottom panel have slightly different T /2π corresponding to (s) ≈ 0.98 × π/T and (s) ≈ 0.62 × π/T , respectively. (b) Same as (a) but for a mode in the stable zone l = 1 (α = 1.9549540063305). Top and bottom panel have slightly different T /2π corresponding to (s) ≈ 0.25 × π/T and (s) ≈ 0.04 × π/T , respectively.

  ), we show for T /2π = 5 the evolution of cos( (s)(t -T )) as well as the scaling factors ρ n = cos( (s)(nT -T )) represented by color dots. From top to bottom, we varied (s) in the primitive spectral cell 0 ≤ (s) ≤ π/T by taking (s) ≈ {0.98, 0.62, 0.25, 0.04}×π/T . Those values of (s) are the ones we obtained in Fig.5from top to bottom. In Fig.6(b) and (c), we divide each position function r([nT -T, nT ]) of Fig. 5(a) and (b) by the corresponding scaling factor ρ n of Fig. 6(a): the 20 first periods of the various computed

Figure 6 :

 6 Figure 6: Periodic mapping and reduction of the neutrally stable vibrational modes to their periodic eigenfunctions for T /2π = 5. (a) Evolution of the n = 20 first scaling factors ρ n = cos( (s)(nT -T )) as a function of number of periods represented by color dots from cold blue to hot red. Black lines show the functions cos( (s)(t-T )). From top to bottom: (s) ≈ 0.98 × π/T , (s) ≈ 0.62 × π/T , (s) ≈ 0.25 × π/T and (s) ≈ 0.04 × π/T . (b) Reduction of the positions r([nT -T, nT ]) of Fig. 5(a) on the eigenfunction Ψ([0, T ]) by ρ n scaling. (c) Same as (b) but for the eigenfunction Ψ 1 (t) of Fig. 5(b).

Figure 7 :

 7 Figure 7: Analytical predictions from a dimensionless Schrödinger equation with a harmonic potential. (a) Limits of T -periodic instability regions of Fig. 2(a) in the (α, T ) space for 1 ≤ α ≤ 2 and 1 ≤ T /2π ≤ 20. Numerical outcomes from the IVP Eq.(1) (green lines) are compared to the analytical results from the BVP Eq.(7) (black lines). Each solution E l = l + 1/2 of Eq.(8) corresponds to a branch of stability, s = 0. (b) In blue line, absolute square of the first four computed eigenfunctions Ψ l (t) of Eq.(1) on the dual primitive cell -vT /2 < χ ≤ vT /2 for T /2π = 5 and v = (2π 2 α/T 2 ) 1/4 = 0.389351968371210, 0.416745399049464, 0.444673895477715, 0.472863276327510 from bottom to top. The eigenfunctions are normalized so that N = T |Ψ(t)| 2 dt = 1. Black lines are analytically obtained from Eq.(9).