Assigning probabilities to non-Lipschitz mechanical systems
AAssigning probabilities tonon-Lipschitz mechanical systems
Danny E. P. Vanpoucke ∗ Sylvia Wenmackers † January 29, 2020
Abstract
We present a method for assigning probabilities to the solutions ofinitial value problems that have a Lipschitz singularity. To illustratethe method, we focus on the following toy-example: ¨ r = r α , r ( t =0) = 0, and ˙ r | r ( t =0) = 0, where the dots indicate time derivatives and α ∈ ]0 , α = 1 / The nineteenth century mathematician and physicist Poisson [1] was thefirst to search for a mechanical interpretation of indeterministic Cauchyproblems. Later in the same century, Boussinesq [2] gave a gravitationalinterpretation of a broad class of such indeterministic Cauchy problems, byconsidering a mass placed at rest at the apex of a frictionless surface froma particular family of hill shapes. This work seems to have been largely for-gotten, but we want to alert physicists to a recent revival of this issue in thephilosophical literature: this question was raised again by a contemporaryphilosopher of science, Norton [3, 4], who focused on a particularly simplecase, now often referred to as
Norton’s dome (presented in section 1). Mala-ment [5] generalized Norton’s example to a family of problems that we willcall
Malament’s mounds (also in section 1). These examples involve initialvalue problems with a differential equation that exhibits a non-Lipschitz sin-gularity. Such non-Lipschitz Cauchy problems are prevalent in the context ∗ Maastricht University, The Netherlands and UHasselt, Faculty of Sciences, Belgium.E-mail: [email protected]. † KU Leuven, Centre for Logic and Philosophy of Science, Belgium. E-mail:[email protected]. a r X i v : . [ phy s i c s . c l a ss - ph ] J a n f physical applications such as turbulent flows and associated dispersion[6], shock waves [7], and collisions in Newtonian N -body problems [8]. Theyare also of interest for the foundations of physics, as a case study on deter-minism and causality, and may be suitable for didactic purposes. The maingoal of the present paper is to demonstrate a method for assigning proba-bilities to trajectories that are solutions to non-Lipschitz Cauchy problems.To illustrate the method, we focus on Norton’s and Malament’s problemsthroughout.Indeterministic theories can be supplemented by hidden variables to ar-rive at deterministic theories, that are empirically equivalent to the formerand that can be used to assign probabilities to the former (see, e.g. , [9]and [10]). As far as we know, this approach has not yet been applied tonon-Lipschitz Cauchy problems. One reason may be that this case requiresusing hidden variables that take on infinitesimal values in the sense of non-standard analysis.Our article is structured as follows. Section 1 presents the shape, initialvalue problem and standard solutions for Norton’s dome and Malament’sgeneralization. In section 2, we refine our research questions and specifyour working hypothesis. In section 3, we determine the finite differenceequation corresponding to the differential equation that represents Norton’sdome and Malament’s mounds. We present numerical results for initial valueproblems based on this equation. In section 4, we introduce concepts fromnon-standard analysis, allowing us to present an alternative model for Nor-ton’s dome and Malament’s mounds, in which we can consider infinitesimaldisplacements. This approach allows us to assign probabilities to the stan-dard solutions. We offer some discussion and review our main conclusionsin section 5. Norton’s problem represents a mass placed with zero velocity at the apex ofa particular dome in a uniform gravitational field. The shape of the domeis chosen such that Newton’s second law applied to the mass takes on aparticularly simple form, as we will see below. Its shape, also shown infig. 1, is: y ( x ) = − (cid:32) − (cid:18) − | x | (cid:19) (cid:33) , where x is the horizontal axis (orthogonal to the gravitational field), y isthe vertical height (anti-parallel to the gravitational field) and the apex is2t the point (0 , x -values.Figure 1: Cross-section of Norton’s dome.Define r ≥ r ( x ) = 1 − (cid:18) − | x | (cid:19) ,y ( r ) = − r . We assume that the gravitational field is constant with g = 1 and that themass is m = 1. Newton’s second law yields an ordinary differential equation(ODE) involving a non-Lipschitz continuous function. The relevant Cauchyproblem involves a second-order non-linear ODE: d r ( t ) dt = √ rr ( t = 0) = 0 dr ( t ) dt | r ( t =0) = 0 . (1)The solution of this problem is non-unique. Besides the trivial, singularsolution, r ( t ) = 0, there is a one-parameter family of regular solutions (see, e.g. , Theorem 2 in [11], due to Kneser [12]), which can be represented geo-metrically as a Peano broom (see fig. 2): r ( t ) = (cid:26) t ≤ T ( t − T ) if t ≥ T, (2)where T is a positive real number that can be interpreted as the time atwhich the mass starts sliding off the dome. The solution can be verified bysubstituting it in the ODE of (1).In the three-dimensional case, there is an additional continuum of possi-bilities regarding the direction of descent. Throughout this paper, we limitourselves to the two-dimensional case (as depicted in fig. 1), such that thisindeterminacy is reduced to two possible directions.3igure 2: Three regular solutions to Norton’s dome with T equal to 0 (black),5 (red), and 10 (green). We get a more general class of indeterministic Cauchy problems if we replacethe shape of Norton’s dome by: y ( x ) = − (cid:16) − (1 − ( α + 1) | x | ) α +1 (cid:17) α +1 α + 1 , with x and y as before and where α is any real number in ]0 , α = 1 / α as 1 / ( α + 1), which reduces to 2 / y as a function of the arc length r measured from theapex yields: y ( r ) = − r α α . The ODE for a mass moving in a gravitational field on a frictionless hillof this family can be obtained by replacing √ r in the ODE of (1) by a power α of r , with α ∈ ]0 , Malament’s mounds [5]. For each choice of α ∈ ]0 , d r ( t ) dt = r α r ( t = 0) = 0 dr ( t ) dt | r ( t =0) = 0 (3)has a trivial, singular solution (again, r ( t ) = 0) and a one-parameter familyof solutions: 4igure 3: Cross-section of five of Malament’s mounds: α = 1 / α = 1 / α = 1 / α = 2 /
3, and α = 9 /
10. The limiting cases with α → α → r ( t ) = t ≤ T (cid:16) (1 − α ) α ) (cid:17) − α ( t − T ) − α if t ≥ T, (4)where T is a positive real number, which can be interpreted as the timeat which the mass starts moving. The solution can again be verified bysubstitution. For α = 1 /
2, the general solution (4) reduces to the solutionof Norton’s special case (2).
Faced with indeterminism due to a lack of Lipschitz-continuity, some authorssearch for arguments that single out a unique solution, e.g. , by regulariza-tion or adding physical principles or heuristics not encoded by the Cauchyproblem itself. The assumption that there is one correct solution (moti-vated by additional physical constraints besides the mathematical equation)is widely—though perhaps not univocally—held in the field of fluid dynam-ics. For instance, the authors of [13] aim to regulate the solutions of Cauchyproblems with non-Lipschitz indeterminism to select a unique global solu-tion.Our current approach is slightly different: lacking a unique solution,we look for a probabilistic description for the trajectory of the mass. Thisapproach may be alien to Newtonian physics (the context in which the prob-lems of section 1 arose), but it is accepted in classical physics more gener-ally ( i.e. , in statistical physics). Hence, we start with the following researchquestion: • Given that there are multiple solutions to Norton’s Cauchy problem,is there a well-supported way to assign probabilities to them?5his research question leads us to two more specific questions: • Can we assign relative probabilities to the singular solution versus thefamily of regular solutions? • Can we assign relative probabilities to the various regular solutions?Our working hypothesis is that the intuition that the Cauchy problemought to have a unique solution, which some physicists express upon encoun-tering Norton’s dome, comes from experience with discrete models, heavilyused in physics. We will therefore aim to represent Norton’s dome in adiscrete, deterministic model.
In this section, we approximate the ODE for our problem (3) by an equationof finite differences and present the resulting numerical simulations.
By discretizing the time parameter ( t = n × ∆ t , with n ∈ N and ∆ t astrictly positive real number < r ( t ) can be transformed into a second-order difference equation forthe discrete sequence R n . This is standard praxis in numerical analysis,which is often used in physics.The first step is to approximate dt by ∆ t , a ‘sufficiently small’ strictlypositive real number. Euler’s method says that if t n +1 = t n + ∆ t , then dr ( t )is approximated by R n − R n − , where R n = R ( t n ) = R ( t n − + ∆ t ). Thisimplies: R n = R n − + R (cid:48) n − ∆ t , which has an error term of O (∆ t ) (see, e.g. , Ch. 16 of [14]). If the approximation is perfect, we have R n = r ( t n ).The smaller the finite difference ∆ t is chosen, the better the approximation.For practical computations, there are better numerical methods (more ac-curate, faster, and/or more stable), but for our current purposes we preferthe conceptually simpler Euler’s method.Since the ODE of Norton’s dome is of second order, we also need toconsider an approximation for the second-order derivative d r/dt . We useEuler’s method again, obtaining: dr/dt ≈ R n − R n − ∆ t ,d r/dt ≈ R n − R n − + R n − ∆ t . R n − R n − + R n − ∆ t = R αn − . Rewriting this as a recurrence relation (to make explicit that the solutioncan be obtained by iteration) and adding the initial conditions, we obtainthe following discrete initial value problem corresponding to Malament’smounds: R n = 2 R n − + ∆ t R αn − − R n − R = 0 R = 0 . (5)Note that we can continue R n for arbitrarily large values, i.e. beyond themaximal height of the mounds. However, we are interested only in the regionaround the apex.The discrete Cauchy problem (5) has a unique solution: it is the constantsequence R n = 0 for all n , which corresponds to the singular solution of thecontinuous Cauchy problem, r ( t ) = 0. This prompts the question: wherehave the other continuous solutions gone in the discrete picture?If we combine the recurrence relation of (5) with different initial condi-tions, it produces a different solution. So, we can regain the regular solu-tions by varying or perturbing the initial conditions: by replacing the zeroin R = 0 and/or R = 0 by some small, non-zero positive real number.Unfortunately, there is no analytic solution known for this non-linear dif-ference equation of second order. Nevertheless, we can continue numerically.In the general case, expressing R n as multiples of ∆ t / (1 − α ) is helpful, be-cause this constant factors out in all terms. (Consider: ∆ t × (cid:0) ∆ t / (1 − α ) (cid:1) α =∆ t / (1 − α ) .)For any values of α , R and R , and for large n , we can fit this numericalsolution to F α ( n ) = (cid:18) (1 − α ) α ) (cid:19) − α (cid:18) n − T ∆ t (cid:19) − α , for some value of T . For example, if we consider α = 1 /
2, ∆ t = 0 .
01, and R = 0 .
01 = R , the fitted T is about − . R = R thatare smaller, we find other regular solutions with a less negative or positive T -value. So, the smaller the perturbation, the later the onset of the motion,as expected. Moreover, the relation between the size of the perturbationand the T -value of the corresponding solution is highly non-linear. This isillustrated, for the case R = R , in fig. 4. The positive values are—in asense—a numerical artefact: after all, they represent cases where the massdid not start at the apex (albeit very close to it) and monotonically movedaway from it, so we should expect a very small negative value. The fact thatwe find positive values is due to the discrete approximation with a finite andnon-infinitesimal time-step ∆ t (in this case, equal to 0.01).7igure 4: T for the case where R = R . These T -values were obtainedby fitting continuous curves of the form F / ( n ) = 1 / n − T / ∆ t ) tosequences of the form R n = 2 R n − + ∆ t R αn − − R n − with α = 1 / t = 0 .
01. The inset gives the positive T -values on a logarithmic R -axis.In general, we need to vary both initial conditions R and R indepen-dently to study the dependence of the discrete perturbation on the solution.We have studied this dependence systematically, as discussed in the follow-ing section. Let us start with an illustration, which will make it easier to understandthe results of our systematic study below. To this end, we first show thenumerical results for two initial conditions that are close but not identical:figs. 5–7 are all based on the difference equation for Norton’s dome ( R n =2 R n − + ∆ t R αn − − R n − with α = 1 /
2) and use the time step ∆ t = 0 . R n − = 0 . R n = 0 . R n − = 0 . R n = 0 . ∼ . / t − T ) / ∆ t to both sequences,the first one (blue) results in T / ∆ t = − .
31 and the second one (orange)results in
T / ∆ t = 90 .
36. The negative T -value means that, at large n , thecurve looks as if the mass has started sliding off before t = 0, which isconsistent with the fact that it did not start at the apex and monotonicallymoved away from it. The positive T -value means that, at large n , the curvelooks as if the mass started sliding off after a delay, which is consistent withthe fact that it first moved closer to the apex before it started sliding off.8igure 5: Shown here are pairs of ( R n − , R n ) for R n = 2 R n − + ∆ t R αn − − R n − with α = 1 / t = 0 .
01. The blue trajectory passes through thepoint R n − = 0 .
01 = R n . The orange trajectory passes through the point R n − = 0 . R n = 0 . R n − , R n ). These discrete curves can be thoughtof as parametrized by time: subsequent dots are a temporal distance of∆ t = 0 .
01 apart. Fig. 6 shows the same two sequences R n as a functionof n ( ∼ t ), for large n : at this scale, the sequences on the one hand andthe continuous curves on the other hand coincide, allowing an excellent fitbetween them. We see that the orange curve reaches the distance of, e.g. ,100,000 at larger n ( i.e. , later in time) than the blue curve. So, the orangecurve is delayed compared to the blue one, which is consistent with thecurves’ T -values. Fig. 7 shows the two sequences R n as a function of n ,now for small n : at this scale, the sequences and the continuous curves arequalitatively different, although the fit between them is excellent for large n ,as we saw in fig. 6. The fitted T -values do not correspond to the minimumof the sequences: T / ∆ t is smaller than the n at which the correspondingsequence attains its minimum. We wrote a program in visual Pascal (Delphi), which allows us to studythe effect of initial conditions in the recurrence equation (5) on the fitted T -value. (If the paper is accepted, we will make the program and the codeavailable at https://github.com/DannyVanpoucke.) For example, for α =1 / t to the value 0.01, and considering the case9igure 6: Values for R n as a function of n ( ∼ t ) for large n for R n =2 R n − + ∆ t R αn − − R n − with α = 1 / t = 0 .
01. The blue trajectorypasses through the point R n − = 0 . R n . The orange trajectory passesthrough the point R n − = 0 . R n = 0 . n . Thedashed curves represent the continuous curves 1 / t − T ) / ∆ t that arefitted to the sequences. 10here R and R are both in the interval [0 , ∆ t ], we find a T -value for eachinitial condition. This is represented using a colour scale in our program (seefig. 8). In practice, the intervals start at a number slightly higher than zero(and much smaller than the upper bound): otherwise the singular solution(at R = 0 , R = 0) is in the field of view (at the left bottom corner),dominating the T -scale.Figure 8: The interface of our program with numerical results for R n =2 R n − + ∆ t R αn − − R n − with α = 1 / t = 0 . R ∈ [0 , ∆ t ] (horizontal, left to right) and R ∈ [0 , ∆ t ](vertical, bottom to top). The number of iterations (10,000) was chosensuch that the value for T was well-converged. We see that solutions with a positive T-value, visible as a narrow red band inthe figures occur for R smaller than R (below the main diagonal R = R ).This is to be expected: if the (sufficiently small) velocity is directed towardsthe top, the mass first moves towards it, before sliding off, thus increasingthe (apparent) T . However, for fixed R , R cannot be chosen arbitrarilysmall, otherwise we select a trajectory that goes over the top and slides offthe other side, resulting in a small positive or negative T .It is instructive to see the results of our program combined with par-ticular sequences or trajectories. This is shown in fig. 9. The part of thetrajectory below the main diagonal corresponds with a mass moving towardsthe apex and coincides with positive T -values.We also studied the dependence of T on R (keeping R fixed). Asan example, we considered α = 1 /
2, ∆ t = 0 .
01 and R = ∆ t = 10 − ;this corresponds to the righthand edge in fig. 9.d. When R is varied from0 to R , T monotonically increases to an asymptote (located near R =0 . × − ) and then monotonically decreases. This is shown infig 10. (Since we do not have an analytic solution to the recurrence equation,we cannot determine the position of the asymptote analytically either.)11igure 9: Numerical results for R n = 2 R n − +∆ t R αn − − R n − with α = 1 / R ∈ [0 , ∆ t m ] (horizontal) and R ∈ [0 , ∆ t m ] (vertical) with∆ t = 0 .
01 and m = 1 (panel a), m = 2 (panel b), m = 3 (panel c), and m = 4 (panel d). The overlay shows in grey one trajectory passing through( R n − , R n ) = (10 − , − ).Figure 10: T -dependence on R at constant R = 10 − for α = 1 / t = 0 .
01. 12ontinuing with the example, the initial condition R = R correspondsto a mass that is released from an arc distance of 10 − with velocity zero: itimmediately starts sliding off from the initial side. This leads to a negative T -value. Initial conditions with R between the asymptote and R corre-spond to a mass that is released from the same distance, but now with apositive velocity towards the top. As R is decreased in this interval, themass moves closer to the apex before sliding down, leading to a monotonicincrease in T -value towards the asymptote.The initial condition R = 0 corresponds to a mass that is released froma distance of the apex with a velocity directed towards the apex, such thatthe mass already reaches the top at n = 1 (or t = ∆ t ): this leads to a massthat rapidly slides off the dome at the other side, also corresponding with anegative T -value. As R is increased up to the aforementioned asymptote,the velocity decreases, leading to a slower descend from the other side, hencethe monotonically increasing T -values.The slope of the T -curve in the R -interval between 0 and the asymptoteis characterized by an infinite sequence of points of inflection. The first oneoccurs at R = 0, the second one at R = 0, the third one at R = 0,etc. In general, R = 2 R + ∆ t R α − R . Solving for R in the case where R = 10 − and R = 0 yields R = 0 . × − exactly. In principle, theposition of the other points of inflection can be computed from the generalexpression for R n , but this becomes impractical quickly. (For instance,the relevant equation for the inflection point corresponding with R = 0is R = 3 R + 2∆ t R α − R + ∆ t (2 R + ∆ t R α − R ) α .) Instead, we havedetermined numerically that the third and fourth inflection points occuraround R = 0 . × − and R = 0 . × − , respectively.The position of the asymptote can be regarded as the limit of this sequenceof inflection point positions, but this does not yield a more practical way ofcomputing it.So far, all numerical results we have shown were for Norton’s dome ( α =1 / α . As α increases, theregion with positive T -values becomes narrower and its slope approachesthe main diagonal. In other words, increasing α looks like zooming out anddecreasing α looks like zooming in as compared to the intermediate casewhere α = 1 /
2. Quantitatively, this scaling behaviour is consistent with thescale factor reported in section 3.1: ∆ t / (1 − α ) .These observations do not yet speak directly to the issue at hand: for anyunstable equilibrium, finite perturbations lead to a movement away from thesource. What is special about Norton’s dome and similar Cauchy problemsis that no perturbation is required—or at least no perturbation that anyreal number larger than zero can express. This brings us to the next stepin our analysis: studying infinitesimal variations in the initial conditions ofwhich the real-valued parts do not differ from zero.13igure 11: Numerical results for R n = 2 R n − + ∆ t R αn − − R n − with ∆ t =0 .
01 in the interval R ∈ [0 , ∆ t ] (horizontal) and R ∈ [0 , ∆ t ]. Panel (a): α = 1 /
10. Panel (b): α = 1 /
3. Panel (c): α = 2 /
3. Panel (d): α = 9 / Hyperfinite-time model for Norton’s dome
It has been suggested that physical praxis is related more closely to non-standard analysis than to standard calculus (see, e.g. , [15]). Following thissuggestion, we can give a hyperfinite description of Norton’s dome and showthat it yields a deterministic model for a mass on such a surface.In order to do this, some general definitions from non-standard analysisare needed. (For a more thorough introduction, see, e.g. , [16].) • In model theory (a branch of logic), a transfer principle relates sen-tences of a given language: it states that sentences are true for somestructure if and only if they are true of another structure. In the caseof non-standard analysis, the relevant sentences are elementary sen-tences ( i.e. , bounded quantifier formulas in first-order logic) applyingto real closed fields ( i.e. , the first-order axiomatization of the fieldof real numbers). The
Transfer Principle of non-standard analysisrelates the standard model of the real numbers with a non-standardmodel, also called the hyperreal numbers. • In what follows, we consider a fixed non-standard model of the naturalnumbers ∗ N ( hypernatural numbers ) and of the real numbers ∗ R ( hyperreal numbers ). • Objects in the non-standard model that can be obtained directly fromthe Transfer Principle are called internal ; objects that are not internalare called external . • An important example of an internal object is a hyperfinite set : it isa set of the form { , , , . . . , N } , where N is an infinite hypernaturalnumber ( N ∈ ∗ N \ N ). • The standard part map st is the function that maps any (limited)hyperreal number to the unique closest real number. A hyperrealnumber x is infinitesimal exactly when its standard part equals zero: st ( x ) = 0. So, taking the standard part can be thought of as ‘roundingoff’ the infinitesimal part. The standard map itself is an externalfunction. • Observe that the Transfer Principle does not apply to second-orderproperties: for instance, the sentence “Every non-empty subset of R which is bounded above has a least upper bound.” is true, but replac-ing R by ∗ R yields a false sentence. In particular, there is no ‘largestinfinitesimal’ in ∗ R . • In what follows, we use a hyperfinite grid to represent the time evo-lution. In particular, we consider a fixed infinite hypernatural number15 and a temporal interval of length 1, which we divide by N intoequal intervals of infinitesimal length ∆ t .Whereas the standard proof of Peano’s existence theorem does not yieldany insight in how to find additional solutions to an indeterministic Cauchyproblem, the non-standard proof does (see [15] p. 32). Our methodology inthe next section is similar to the approach in this proof. (See in particular[17].) Within the scope of standard analysis, the discrete approach using differenceequations is only an approximation to the continuous case described byODEs—an approach that improves as the time step ∆ t decreases. In thenon-standard approach, however, we choose an infinitesimal ∆ t , smallerthan any strictly positive real number. The discrete equations look thesame as before, except that we now assume T n and R n to be sequencesthat take values in the hyperreal numbers for all hypernatural numbers n inthe hyperfinite set { , , . . . , N } . In such a hyperfinite model, we considerinitial conditions such that R and R are infinitesimal: since the standardpart of these infinitesimals is zero, this is consistent with the real-valueddescription of the Cauchy problem. Since there is no largest infinitesimal,we will restrict R and R to the hyperreal interval [ − ∆ t, ∆ t ].In the numerical study of section 3, we observed the following: if wechange the value of ∆ t , we obtain the same visual result if we divide thevalues for R and R by ∆ t for the dome, or by ∆ t / (1 − α ) for Malament’smounds. Since this is true for any positive real ∆ t , it remains true forinfinitesimal ∆ t (by the Transfer Principle).If we take R = R = ∆ t , we observe the following: For finite n , R n is infinitesimal (its standard part is zero). For infinite n , st ( R n ) =1 / t − ( − ∆ t / )) . This corresponds with the general solution (2), with T = − ∆ t / —a negative infinitesimal. So, for all n , st ( R n ) = r ( t ), wherethe latter is one of the solutions to the continuous Cauchy problem. Inparticular, we find the undelayed solution: r ( t ) = 1 / t .The indeterminism in the usual model (continuous, without infinites-imals) can be interpreted as being due to rounding off the infinitesimalsfrom the hyperfinite model.Given that any physical measurement is only finitely precise, it is notpossible to distinguish experimentally between zero and infinitesimal quan-tities. This is similar for chaotic systems without any singularities: laterstates can be used to determine their initial positions beyond measurementprecision at the time. (But see [10], which suggested to reinterpret this asindeterminism after all.) Nevertheless, the resulting trajectories are all dif-ferent and their standard parts do agree with the family of solutions found16n the standard model. This allows us to assign probabilities to the standardsolutions.Even though we cannot do any numerical experiments with actual in-finitesimals on a physical computer, our numerical experiments from theprevious section are instructive: even if we change ∆ t to a hyperreal in-finitesimal, we obtain the same result if we divide the values for R and R by ∆ t / (1 − α ) . In other words, we do not need to compute with actual in-finitesimals to find R / ∆ t / (1 − α ) and R / ∆ t / (1 − α ) , and the results obtainedin the previous section are obtained here as well. We could propose to use the previous information as follows: first considerhyperreal initial conditions ( R , R ) that are randomly chosen from [0 , ∆ t ] and then compute the resulting probabilities for obtaining the singular so-lution and for the values of st ( T ) associated with the family of regularsolutions. The first step amounts to assuming a uniform prior on the phasespace. Admittedly, this approach may be contested, since we could haveused some other representation of the phase space via an arbitrary coordi-nate transformation. Our next question, then, is how to make a principledchoice. So far, we have used the Newtonian formalism with two generalized coordi-nates: the arc length and the arc velocity. It is well known that measures onthe phase space change due to coordinate transformations. A non-arbitrarychoice for the probability measure is a uniform measure that is invariantunder the dynamics [18]. In other words, we need phase space coordinatessuch that the density of states is conserved on the state space. For this, wehave to consider the Lebesgue measure on the canonical coordinates fromHamiltonian mechanics (such that Liouville’s theorem applies). For ourproblem, the canonical coordinate is the arc length and the conjugate mo-mentum equals mass times the arc velocity. Since we assumed m = 1, thephase space spanned by ( R , V ) is the proper starting point for a probabilis-tic analysis that allows us to start from a uniform probability distribution(maximum entropy).Hence, we use ( R , V ) to represent the initial conditions (where V n =( R n +1 − R n ) / ∆ t ), as shown in panel (b) of fig. 12. The vector field onthis phase plane is shown in fig. 13. The symmetry between the upperand lower half of the plane shows that the dynamics is time reversible. Ifwe reinterpret R as a position rather than an arc length and include theopposite side of the slope as negative R -values (not shown), the vector field17igure 12: Same numerical results represented in two ways. Numericalresults for R n = 2 R n − + ∆ t R αn − − R n − with α = 1 / t = 0 .
01 inthe interval R ∈ [0 , ∆ t ] (horizontal) and R ∈ [0 , ∆ t ]. Panel (a): ( R , R )-plane. Panel (b): ( R , V )-plane. The diagonal in panel (a) corresponds to R = R : that is the horizontal axis in panel (b), i.e. , V = 0.shows the signature of a saddle point at the singularity (0 ,
0) (indicating anunstable equilibrium).However, the prior need not be motivated by the dynamics, which wasthe reasoning behind privileging the uniform measure on the ( R , V )-space.Instead, we may consider a process that aims to place the mass as close tothe top as possible with a velocity as close to zero as possible. In that case,we may want to consider a Gaussian distribution around the origin in the( R , V )-plane instead of a uniform distribution. Fortunately, as we will seein section 4.2.2, our main results (in terms of real-valued probabilities for the T -values) turn out to be quite robust: they obtain for a uniform measure over( R , V ) as well as over ( R , R ) (which is related to the canonical choice viaa linear transformation) and Gaussian transformations of them. Moreover,our assignment of a real-valued probability to the singular solution is fullyrobust, in the sense that this result does not depend on the transformationat all. We already observed that exactly one combination of hyperreal-valued initialconditions leads to the equilibrium solution: R = R = 0. All other initialconditions are associated to a regular solution. This means that, as longas the prior does not assign a non-infinitesimal portion to the single point R = R = 0, the singular solution has infinitesimal probability—or zero ifwe take the standard part.Let us now focus on the vast majority of initial conditions that lead toa regular solution. We already observed that the relation between the in-18igure 13: Vector field on the phase space ( R , V ). Numerical results for R n = 2 R n − + ∆ t R αn − − R n − with α = 1 / t = 0 .
01 in the interval R ∈ [0 , t ] (horizontal) and V ∈ [ − . t, . t ].finitesimal initial conditions and the T -parameter in the standard part ofthe corresponding solution is strongly non-linear. This is shown explicitlyin fig. 4 for the case with R = R (for values corresponding to the diag-onal in panel (a) of fig. 12). These observations hold in general: only aninfinitesimal proportion of all infinitesimal initial conditions correspond tostandard solutions with st ( T ) >
0. Almost all infinitesimal initial conditionscorrespond to a standard solution with st ( T ) = 0.If we assume a uniform distribution of the infinitesimal initial conditions,we arrive at the following probabilities for the standard solutions: • The hyperreal-valued probability of the mass staying at the apex ofNorton’s dome or any of Malament’s mounds indefinitely (singularsolution) is infinitesimal; the standard part of this probability is zero. • The hyperreal-valued probability of the mass staying at the apex ofNorton’s dome or any of Malament’s mounds for some observable timeis infinitesimal; the standard part of this probability is zero. • The hyperreal-valued probability of the mass immediately starting toslide off the apex of Norton’s dome or any of Malament’s mounds isone minus an infinitesimal; the standard part of this probability isunity.In conclusion, a point mass placed at the apex of a frictionless Norton’sdome or any of Malament’s mounds in a gravitational field will immediatelystart sliding off the dome almost surely. This conclusion also applies if theuniform measure on the ( R , V )-plane is replaced by a Gaussian (cut-off19or non-infinitesimal values, since they contradict the conditions set by thestandard initial value problem).In addition, for a mass rolling to the apex, reaching it with a standardvelocity of exactly zero, it will roll off from the opposite side immediately,almost surely. In this case, the initial values are sampled from a differentpart of the phase plane, but it remains the case that those correspondingto any measurable delay are of infinitesimal measure as compared to thoseyielding a measurable delay. Moreover, whereas the original Cauchy problemgave no information regarding the direction of the infinitesimal perturbation,observing the initial direction does allow us to predict the final direction. First, we comment on two possible applications of this work. Then we drawgeneral conclusions.
So far, we focused on toy-examples. However, differential equations with anon-Lipschitz singularity are prevalent in the context of physical applica-tions such as shock formation and turbulent flows (in particular, in the highviscosity limit in the Kolmogorov theory of turbulence). In the case of weaksolutions to the inviscid Burgers equation, researchers can pick out the rightsolution from infinity. But this approach does not apply to all singulari-ties. After a trajectory encounters such a singularity in finite time (knownas ‘blowup’), the evolution is no longer deterministic: there are continuummany solutions. One way to understand this is that the initial state repre-sents a Dirac-delta-distribution of initial conditions: whereas this remainsa Dirac-delta-distribution for fully deterministic evolutions, non-Lipschitzsingularities make the delta-distribution spread out to a ‘spontaneous’ prob-ability distribution—a phenomenon known as ‘spontaneous stochasticy’.Using methods related to those in this article, we can re-interpret theDirac-delta distribution as a function from non-standard analysis (an in-finitely small Gaussian), non-singular points as finite dispersion (such thatinfinitesimal differences remain infinitesimal), and the singular point as aplace where there is infinite dispersion (such that infinitesimal differencesbecome non-infinitesimal). Viewed as such, the non-Lipschitz indeterminismof the continuous model can be connected to a case of deterministic chaosin a corresponding hyperfinite model.
Classical mechanics is often the first physics course in the curriculum forBachelor students of Physics, at which point these students may lack the20equired knowledge of calculus, and the implicit assumptions are rarely re-flected upon at later stages. Textbooks on standard calculus often mentionan example of a first-order non-linear ODE that has multiple solutions, onlyto move on quickly to those equations that are uniquely solvable. In text-books on physics, it is often assumed that the functions appearing in theequations are ‘sufficiently well-behaved’, the meaning of which is context-dependent, but which is typically made explicit for the theory of interestat least once in the text ( e.g. , infinitely differentiable potentials and wavefunctions in quantum mechanics). In the context of the Newtonian lawsof motion, however, it is usually left unspecified that one should assumeLipschitz continuity to arrive at a deterministic description.So, we suggest using Norton’s dome to alert students of the importance ofchecking the—often implicit—‘sufficiently well-behaved’ assumption. More-over, this case is relevant for showing the historical development of bothcalculus and classical mechanics (a topic which we plan to develop in aseparate article).
In this paper, we presented a method for assigning probabilities to the so-lutions of initial value problems that lack Lipschitz continuity. First, wereplaced the differential equations by finite difference equations, which aredeterministic and allow for systematic numerical studies. Second, startingfrom a uniform prior on the phase space spanned by canonical coordinates,we assigned probabilities to the solutions of the corresponding, continuousCauchy problem. Third and finally, to avoid the reliance on finite pertur-bations, we used infinitesimals in the sense of non-standard analysis; thisapproach is conceptually close to numerical methods from physical praxis.The current article focused on toy-examples (Norton’s dome and its gener-alization as Malament’s mounds). Although we set out to find a probabilitydistribution over the solutions, we find unit probability for one single solu-tion (the non-delayed, regular solution). So, although we did not assumeuniqueness at the outset, we do find ourselves in agreement with authorswho set out to find a unique continuation beyond the non-Lipschitz dis-continuity. Moreover, our method is applicable to study other initial valueproblems that lack Lipschitz continuity. Hence, we hope our approach willbe fruitful for application to the study of singularities in hydrodynamics andother blowup phenomena.
Acknowledgement
SW is grateful to Vieri Benci for helpful discussions on this topic. Partof SW’s work was supported by the Research Foundation Flanders (FondsWetenschappelijk Onderzoek, FWO), Grant No. G0B8616N.21 eferences [1] S.-D. Poisson,
M´emoire sur les solutions particuli`eres des ´equationsdiff´erentielles et des ´equations aux diff´erences , vol. 6. 1806.[2] J. Boussinesq, “Conciliation du v´eritable d´eterminisme m´ecanique avecl’existence de la vie et de la libert´e morale,”
M´emoires de la soci´et´e dessciences, de l’agriculture et des arts de Lille , vol. 6, pp. 35–251, 1879.[3] J. D. Norton, “Causation as folk science,”
Philosophers’ Imprint , vol. 3,2003.[4] J. D. Norton, “The dome: An unexpectedly simple failure of determin-ism,”
Philosophy of Science , vol. 75, pp. 786–798, 2008.[5] D. B. Malament, “Norton’s slippery slope,”
Philosophy of Science ,vol. 75, pp. 799–816, 2008.[6] D. Benveniste and T. D. Drivas, “Asymptotic results for backwardstwo-particle dispersion in a turbulent flow,”
Physical Review E , vol. 89,p. 041003, 2014.[7] D. Serre and A. F. Vasseur, “About the relative entropy method for hy-perbolic systems of conservation laws,” in
A panorama of mathematics:pure and applied (C. M. da Fonseca, D. Van Huynh, S. Kirkland, andV. K. Tuan, eds.), vol. 658 of
Contemporary mathematics , pp. 237–248,Providence, Rhode Island: American Mathematical Society, 2015.[8] L. F. Bakker, “Understanding the dynamics of collision and near-collision motions in the n-body problem,” in
Advances in Interdis-ciplinary Mathematical Research (B. Toni, ed.), vol. 37 of
SpringerProceedings in Mathematics & Statistics , pp. 99–115, New York, NY:Springer, 2013.[9] C. Werndl, “Are deterministic descriptions and indeterministic descrip-tions observationally equivalent?,”
Studies In History and Philosophyof Modern Physics , vol. 40, pp. 232–242, 2009.[10] N. Gisin, “Indeterminism in physics, classical chaos and Bohmian me-chanics: Are real numbers really real?.” Forthcoming in Erkenntnis;doi: 10.1007/s10670-019-00165-8, 2019.[11] G. Derks, “Existence and uniqueness of solutions of initial value prob-lems,” in
Encyclopedia of Complexity and Systems Science (R. Meyers,ed.), pp. 3255–3267, New York, NY: Springer, 2009.[12] H. Kneser, “ ¨Uber die L¨osungen eines Systems gew¨ohnlicher Differen-tialgleichungen das der Lipschitzschen Bedingung nicht gen¨ugt,”
S.-B.Preuss. Akad. Wiss. Phys.-Math. Kl. , pp. 171–174, 21923019.2213] T. D. Drivas and A. A. Mailybaev, ““Life after death” in ordinarydifferential equations with a non-Lipschitz singularity.” arXiv preprintarXiv:1806.09001, 2018.[14] W. H. Press, S. A. Teukolsky, W. T. Vetterling, and B. P. Flannery,
Numerical Recipees in C++ . New York, NY: Cambridge UniversityPress, 2005. second edition.[15] S. Albeverio, J. E. Fenstad, R. Hoegh-Krøhn, and T. Lindstrøm,
Non-Standard Methods in Stochastic Analysis and Mathematical Physics .Pure and Applied Mathematics, Orlando, FL: Academic Press, 1986.[16] V. Benci, M. Di Nasso, and M. Forti, “The eightfold path to nonstan-dard analysis,” in
Nonstandard Methods and Applications in Mathe-matics (N. J. Cutland, M. Di Nasso, and D. A. Ross, eds.), vol. 25 of
Lecture Notes in Logic , pp. 3–44, Wellesley, MA: Association for Sym-bolic Logic, AK Peters, 2006.[17] B. Birkeland and D. Normann, “A non-standard treatment of the equa-tion y (cid:48) = f ( y, t ),” Preprint series: Pure mathematics http://urn. nb.no/URN: NBN: no-8076 , 1980.[18] S. Goldstein, “Typicality and notions of probability in physics,” in