Quantum phenomena modelled by interactions between many classical worlds

We investigate whether quantum theory can be understood as the continuum limit of a mechanical theory, in which there is a huge, but finite, number of classical 'worlds', and quantum effects arise solely from a universal interaction between these worlds, without reference to any wave function. Here a `world' means an entire universe with well-defined properties, determined by the classical configuration of its particles and fields. In our approach each world evolves deterministically; probabilities arise due to ignorance as to which world a given observer occupies; and we argue that in the limit of infinitely many worlds the wave function can be recovered (as a secondary object) from the motion of these worlds. We introduce a simple model of such a 'many interacting worlds' approach and show that it can reproduce some generic quantum phenomena---such as Ehrenfest's theorem, wavepacket spreading, barrier tunneling and zero point energy---as a direct consequence of mutual repulsion between worlds. Finally, we perform numerical simulations using our approach. We demonstrate, first, that it can be used to calculate quantum ground states, and second, that it is capable of reproducing, at least qualitatively, the double-slit interference phenomenon.


I. INTRODUCTION
The role of the wave function differs markedly in various formulations of quantum mechanics. For example, in the Copenhagen interpretation it is a necessary tool for calculating statistical correlations between a priori classical preparation and registration devices [1]; in the de Broglie-Bohm interpretation it acts as a pilot wave that guides the world's classical configuration [2]; in the many-worlds interpretation it describes an everbranching tree of noninteracting quasiclassical worlds [3]; and in spontaneous collapse models its objective 'collapse' creates a single quasiclassical world [4].
In other formulations the wave function does not even play a primary role. For example, in Madelung's quantum hydrodynamics [5], Nelson's stochastic dynamics [6], and Hall and Reginatto's exact uncertainty approach [7], the fundamental equations of motion are formulated in terms of a configuration probability density P and a momentum potential S (or the gradient of the latter), with a purely formal translation to a wave function description via Ψ := P 1/2 exp[iS/ ]. These approaches can describe the evolution of any scalar wave function on configuration space, which includes any fixed number of spinless particles, plus bosonic fields. In this paper we will similarly treat spinless and bosonic degrees of freedom.
More recently, it has been observed by Holland [8] and by Poirier and coworkers [9][10][11] that the evolution of such quantum systems can be formulated without reference even to a momentum potential S. Instead, nonlinear Euler-Lagrange equations are used to define trajectories of a continuum of fluid elements, in an essentially hydrodynamical picture. The trajectories are labelled by a continuous parameter, such as the initial position of each element, and the equations involve partial derivatives of * Electronic address: H.Wiseman@Griffith.edu.au up to fourth order with respect to this parameter.
In the Holland-Poirier hydrodynamical approach the wave function plays no dynamical role. However, it may be recovered, in a nontrivial manner, by integrating the trajectories up to any given time [8]. This has proved a useful tool for making efficient and accurate numerical calculations in quantum chemistry [9,10]. Schiff and Poirier [11], while "drawing no definite conclusions", interpret their formulation as a "kind of "many worlds" theory", albeit they have a continuum of trajectories (i.e. flow lines), not a discrete set of worlds.
Here we take a different but related approach, with the aim of avoiding the ontological difficulty of a continuum of worlds. In particular, we explore the possibility of replacing the continuum of fluid-elements in the Holland-Poirier approach by a huge but finite number of interacting 'worlds'. Each world is classical in the sense that it has determinate properties which are functions of its configuration. In the absence of the interaction with other worlds, each world evolves according to classical Newtonian physics. All quantum effects arise from, and only from, the interaction between worlds. We therefore call this the 'many interacting worlds' (MIW) approach to quantum mechanics. A broadly similar idea has been independently suggested by Sebens [12], although without any explicit model being given.
The MIW approach can only become equivalent to standard quantum dynamics in the continuum limit, where the number of worlds becomes uncountably infinite. However, we will show that even in the case of just two interacting worlds it is a useful toy model for modeling and 'explaining' quantum phenomena, such as the spreading of wavepackets and tunneling through a potential barrier. Regarded as a fundamental physical theory in its own right, the MIW approach may also lead to new predictions arising from the restriction to a finite number of worlds. Moreover, it provides a natural discretisation of the Holland-Poirier approach which may be useful for numerical purposes.

A. Comparative ontology of the MIW approach
The natural ontology of the MIW approach provides a realist interpretation of quantum phenomena that has some advantages when compared to other realist interpretations, such as the de Broglie-Bohm (dBB) and many-worlds (MW) interpretations. To further motivate our work, we briefly discuss these here.
In the MIW approach each world is well-defined by a definite configuration that evolves deterministically. Probabilities arise only because observers are ignorant of which world they actually occupy, and so assign an equal weighting to all worlds compatible with the macroscopic state of affairs they perceive. In a typical quantum experiment, where the outcome is indeterminate in orthodox quantum mechanics, the final configurations of the worlds in the MIW approach can be grouped into different classes based on macroscopic properties corresponding to the different possible outcomes. The orthodox quantum probabilities will then be approximately proportional to the number of worlds in each class.
In contrast, the dBB interpretation postulates a single world, deterministically guided by a physical wave function. However, this world only occupies a single point of configuration space, and does not exert any back reaction on the guiding wave, which has no source but which occupies the entire configuration space. This makes it difficult to give an ontology for the wave function in parts of configuration space so remote from the 'real' configuration that it will never affect it (a nice analogy can be found in Feynman's criticism of classical electromagnetism [13]). Furthermore, this wave function also determines a probability density for the initial world configuration [2,14]. From a Bayesian perspective this dual role is not easy to reconcile [15].
In the Everett or MW interpretation, the 'worlds' are orthogonal components of a universal wave function [3]. The particular decomposition at any time, and the identity of worlds through time is argued to be defined (at least well-enough for practical purposes) by the quantum dynamics which generates essentially independent evolution of these quasiclassical worlds into the future (a phenomenon called effective decoherence). The inherent fuzziness of Everettian worlds is in contrast to the corresponding concepts in the MIW approach, of a well-defined group of deterministically-evolving configurations. In the MW interpretation it is meaningless to ask exactly how many worlds there are at a given time, or exactly when a branching event into subcomponents occurs, leading to criticisms that there is no precise ontology [16]. Another difficult issue is that worlds are not equally 'real' in the MW interpretation, but are 'weighted' by the modulus squared of the corresponding superposition coefficients. As noted above, in the MIW approach all worlds are equally weighted, so that Laplace's theory of probability is sufficient to account for our experience and expectations.
B. Outline of the paper In Sec. II A we introduce the MIW approach in detail. The approach allows flexibility in specifying the precise dynamics of the interacting worlds. Hence, we consider a whole class of MIW models which, as we discuss in Sec. II B, should agree with the predictions of orthodox quantum mechanics in the limit of infinitely many worlds. We present a particularly simple model in Sec. III and demonstrate in Sec. IV that it is capable of reproducing some well-known quantum phenomena: Ehrenfest's theorem, wavepacket spreading, barrier tunneling, and zero-point energy. We conclude with a numerical analysis of oscillator ground-states and double-slit interference phenomena in Secs. V and VI, respectively.

A. From dBB to MIW
For pedagogical reasons we introduce MIW from the perspective of the dBB interpretation of quantum mechanics, and hence, start with a brief review of the latter. We regard the MIW approach as fundamental, but to show that it should (in an appropriate limit) reproduce the predictions of orthodox quantum mechanics, it is convenient to build up to it via the equations of dBB mechanics. It is well known that dBB interpretation reproduces all the predictions of orthodox quantum mechanics, so far as the latter is well defined [14,17].
Consider then a universe comprising J scalar nonrelativistic distinguishable particles, in a D-dimensional space. (As mentioned in the introduction we could also include bosonic fields, but for simplicity we omit them here.) In dBB mechanics there is a universal (or 'world') wave function Ψ t (q) defined on configuration space. Here the argument of wavefunction is a vector q = q 1 , · · · , q K of length K = DJ. We order these variables so that for D = 3, the vector (q 3j−2 , q 3j−1 , q 3j ) describes the position of the jth particle. For convenience we associate a mass with each direction in configuration space, so that the mass of the jth particle appears here 3 times, as m 3j−2 = m 3j−1 = m 3j . Then we can write Schrödinger's equation for the world wave function Ψ t (q) as As well as the world wave function, there is another part of the dBB ontology, corresponding to the real positions of all the particles. We call this the world-particle, with position x(t) = x 1 (t), · · · , x K (t) , which we also call the world configuration. Reflecting "our ignorance of the initial conditions" [2], the initial position x(0) of the world-particle is a random variable distributed according to probability density P 0 (x), where The velocity of the world-particle,ẋ, is then defined by where This equation of motion guarantees that the probability density for the world-configuration x(t) at any time t is given by P t (x). This property is known as equivariance [14,18]. In Bohm's original formulation [2], the law of motion (3) is expressed equivalently by the second-order equation with Eq. (3) applied a constraint on the velocity at the initial time (t = 0). Here the force has been split into classical (f ) and quantum (r) contributions, the latter called r because of its locally repulsive nature, which will be shown later. These are defined by (with the kth component of ∇ being ∂/∂q k ). Here was called the quantum potential by Bohm [2], and vanishes for = 0. It can be shown that Eq. (5) reproduces Eq. (3) at all times. Although Eq. (6) looks like Newtonian mechanics there is no conserved energy for the world-particle alone, because Q t is time-dependent in general. Moreover, the wave function Ψ t (q) evolves in complete indifference to the world particle, so there is no transfer of energy there. Thus, these dynamics are quite unlike those familiar from classical mechanics.
Suppose now that instead of only one world-particle, as in the dBB interpretation, there were a huge number N of world-particles co-existing, with positions (worldconfigurations) x 1 , . . . , x n , . . . , x N . If each of the N initial world-configurations is chosen at random from P 0 (q) as described above, then by construction. The approximation is in the statistical sense that the averages of any sufficiently smooth function ϕ(q), calculated via either side, will be approximately equal, and it becomes arbitrarily good in the limit N → ∞. Clearly, a similar approximation also holds if the empirical density on the right hand side is replaced by a suitably smoothed version thereof. By equivariance the quality of of either approximation is then conserved for all times. One can thus approximate P t (q), and its derivatives, from a suitably smoothed version of the empirical density at time t. From this smoothed density one may also obtain a corresponding approximation of the Bohmian force (6) in terms of the list of world-configurations at time t. Note in fact that since only local properties of P t (q) are required for r t (q), the approximation r N (q; X t ) requires only worlds from the set X t which are in the K-dimensional neighbourhood of q. That is, the approximate force is local in K-dimensional configuration space. We now take the crucial step of replacing the Bohmian force (6), which acts on each world-particle x n (t) via (5), by the approximation r N (x n (t); X t ). Thus, the evolution of world-configuration x n (t) is directly determined by the other configurations in X t . This makes the wave function Ψ t (q), and the functions P t (q) and S t (q) derived from it, superfluous. What is left is a mechanical theory, referred to as MIW, which describes the motion of a 'multiverse' of N co-existing worlds x 1 (t), . . . x n (t), . . . x N (t), where each world-configuration x n (t) is a K-vector specifying the position of J = K/D particles.
While the MIW approach has been motivated above as an approximation to the dBB interpretation of quantum mechanics, we have the opposite in mind. We would like to regard MIW as the fundamental theory, from which under certain conditions dBB can be recovered as an effective theory provided N is sufficiently large; see Sec. II B below. Note that the MIW approach is conceptually and mathematically very different from dBB. Its fundamental dynamics are described by the system of N × J × D second-order differential equations In the absence of an interworld interaction, corresponding to the classical limit r k N (x n (t); X t ) = 0 in Eq. (11), the worlds evolve independently under purely Newtonian dynamics. Hence, all quantumlike effects arise from the existence of this interaction. It will be seen in Secs. III-VI that this nonclassical interaction corresponds to a repulsive force between worlds having close configurations, leading to a simple and intuitive picture for many typical quantum phenomena.
Note also that we use the term 'MIW approach' rather than 'MIW interpretation'. This is because while some predictions of quantum mechanics, such as the Ehrenfest theorem and rate of wavepacket spreading, will be seen to hold precisely for any number of worlds, other predictions can be accurately recovered only under certain conditions in the limit N → ∞. This has two immediate implications: the possibility of experimental predictions different from standard quantum mechanics, due to the finiteness of N , and the possiblity of using the MIW approach for approximating the dynamics of standard quantum systems in a controlled manner. In the latter case Eq. (11) must be supplemented by suitable initial conditions, corresponding to choosing the initial world-configurations randomly from P 0 , and the initial world-particle velocities from S 0 via Eq. (3) (see also below and Sec. VI).

B. Probabilities and the quantum limit
While each world evolves deterministically under Eq. (11), which of the N worlds we are actually living in, compatible with the perceived macroscopic state of affairs, is unknown. Hence, assertions about the configuration of the J particles in our world naturally become probabilistic, as it was the case in dBB. In particular, for a given function ϕ(x) of the world configuration, only an equally-weighted population mean over all the worlds compatible with observed macroscopic properties, can be predicted at any time.
We shall now show that under certain conditions the MIW expectation values (12) are expected to converge to the ones predicted by quantum theory when the number of worlds N tends to infinity. For this let us suppose we are provided a solution to the Schrödinger equation Ψ t on K-dimensional configuration space, and the initial configurations of the N worlds, x 1 (0), x 2 (0), . . . , x N (0) are approximately distributed according to the distribution P 0 (q) = |Ψ 0 (q)| 2 . Hence, at t = 0 one has, for any smooth function ϕ on configuration space, with the approximation becoming arbitrarily good as N → ∞ for ϕ sufficiently regular. Suppose further that Ψ t is such that at initial time t = 0 the velocities of the worlds fulfill Eq. (3). In consequence of (9), the trajectory of the nth world generated by (11) should stay close to the corresponding Bohmian trajectory generated by Eq. (5). Hence, by equivariance of the latter trajectories, the world configurations x 1 (t), . . . , x N (t) at time t will be approximately distributed according to the distribution as desired.
In summary, the configuration-space expectation values of the MIW approach should coincide with those for the quantum state Ψ t as N → ∞ provided that at some time τ the following conditions are met: 2. The velocities of the worlds fulfill Eq. (3) at t = τ .
Since configuration space expectation values are all that is required to establish empirical equivalence with orthodox quantum mechanics (just as in dBB [2]), this establishes the viability of the MIW approach.
It is a remarkable feature of the MIW approach that only a simple equal weighting of worlds, reflecting ignorance of which world an observer occupies, appears to be sufficient to reconstruct quantum statistics in an suitable limit. Similarly, in this limit, the observer should also see statistics as predicted by quantum mechanics when carrying out a sequence of experiments in his or her single world. In particular, from the typicality analysis by Dürr et al. for the dBB interpretation [18], it follows that Born's rule holds for dBB trajectories belonging for typical initial configurations, where in dBB typical stands for almost all with respect to the |Ψ t | 2 measure. Hence, since the MIW trajectories are expected to converge to dBB trajectories and the world configurations are |Ψ t | 2 distributed in the limit described above, Born's rule will hold for typical MIW worlds, where we emphasise again that in our MIW approach typicality simply means for the great majority of worlds, since each world is equally weighted.

C. Interworld interaction potential
The general MIW approach described in Sec. II A is only complete once the form of the force r N (x; X) between worlds, in Eq. (11), is specified. There are different possible ways of doing so, each leading to a different version of the approach.
It is natural to seek a formulation of the MIW approach in which the interaction force between worlds is guaranteed to be conservative. That is, a force of the form for some potential function U N (X) defined on the N world configurations X = (x 1 , . . . , x N ), where ∇ xn denotes the gradient vector with respect to x n . This imme-diately allows the equations of motion (11) to be rewritten in the equivalent Hamiltonian forṁ with all its attendant advantages. Here P = (p 1 , . . . , p n ) defines the momenta of the worlds, with components and the Hamiltonian is given by (17) We will refer to U N (X) as the interworld interaction potential.
To motivate the form of suitable interworld potentials U N , note that, to reproduce quantum mechanics in the limit N → ∞, the average energy per world, E N ≡ N −1 H N (X, P), should approach the quantum average energy in this limit. Note that we cannot associate a definite energy with each world, because of the interworld interaction. Note also that we have substituted the subscript N for the subscript X (or X, P as would be needed in this case), for ease of notation. Now, if the configurations x 1 , . . . , x N sample the distribution |Ψ(x)| 2 , the quantum average energy can be written, (for N sufficiently large). Moreover, the average energy per world is given via Eq. (17) as Assuming that the trajectories generated by the Hamiltonian H N are close to the dBB trajectories for sufficiently large N , then p k n = m kẋk n ≈ ∂S(x n )/∂x k n , and comparing the two averages shows that U N (X) should be chosen to be of the form where Here the left-hand side is to be understood as an approximation of the right-hand side, obtained via a suitable smoothing of the empirical density P (q) in Eq. (8), analogous to the approximation of the quantum force r N (q) by r N (x; X) in Eq. (9). It is important to note that a good approximation of the force (which is essential to obtain QM in the large N limit) is not guaranteed by a good approximation in Eq. (19). (This is easy to see, since for large N an approximation to U N (X) could use only a subset of the N worlds, but then the force (14) would be zero for all worlds not contained in this subset.) Thus Eq. (19) is a guide only, and Eq. (9) must also be checked. The interworld potential U N in Eq. (18) is positive by construction. This leads directly to the existence of a minimum energy and a corresponding stationary configuration for the N worlds, corresponding to the quantum groundstate energy and (real) groundstate wavefunction P −1/2 (q) in the limit N → ∞. Moreover, g k N (q; X)/m k in Eq. (19) approximates Nelson's osmotic velocity [6], suggesting that U N in Eq. (18) may be regarded as a sum of nonclassical kinetic energies, as is explored in Sec. IV. Further connections between U N (X), g N (q; X) and r N (x; X) are discussed in Appendix A.

III. SIMPLE EXAMPLE
A full specification of the MWI dynamics requires generation of a suitable force function r N (x n ; X) in Eq. (11) or a suitable function g k N (q; X) in Eq. (18). As discussed in the previous section, these may be viewed as corresponding to approximations of −∇Q and ∇P/P , respectively. One may obtain many candidates for interworld forces and potentials in this way, and it is a matter of future interest to determine what may be the most natural ones.
For the purposes of initial investigations, we will primarily restrict attention to the simplest case of onedimensional worlds (K = D = J = 1). In this section we give an example of the MIW approach for this onedimensional case, using the Hamiltonian formulation. In Sec. IV we show that this example provides a nice toy model for successfully describing various quantum phenomena.

A. One-dimensional toy model
For N one-dimensional worlds, it is convenient to distinguish them from the general case by denoting the configuration x n and momentum p n of each world by x n and p n , respectively. It is also convenient to label the configurations such that x 1 < x 2 < · · · < x N (it will be seen that this ordering is preserved by the repulsive nature of the interaction between worlds). Similarly, we will use X = (x 1 , . . . , x N ) in place of X, and P in place of P. Note that we can thus consider X and P as vectors. The aim here is to obtain a simple form for the interworld potential U N in Eq. (18). It is easiest to first approximate the empirical distribution of worlds by a smooth probability density, P X (q), and use this to obtain a suitable form for g k N (q; X) in Eq. (19). Now, any smooth interpolation P X (q), of the empirical density of the N worlds, must satisfy, for sufficiently slowly varying functions ϕ(x) and sufficiently large N . This suggests the following ansatz for the smoothed distribution P X , which relies on an assumption that the inter-world separation is slowly varying. The success of this ansatz requires that the dynamics we derive based upon it preserve this slow variation, which appears to be the case from the simulations in Sec. VI. Under this assumption, the two expressions in Eq. (20) have a relative difference O(1/N ). Further, for large N , we can approximate the derivative of P X (q) at q = x n by Hence, ∇P X /P X is approximated by It follows that one may take the interworld potential U N (X) in Eq. (18) to have the rather simple form To ensure the summation is well defined for n = 1 and n = N , we formally define x 0 = −∞ and x N +1 = ∞.

B. Basic properties of the toy model
The MIW toy model corresponding to Eq. (22) is defined by the equations of motion (15) and Hamiltonian (17). Since the total energy is conserved, two adjacent worlds x n and x n+1 cannot approach each other arbitrarily closely, as the potential U N would become unbounded. This implies a mutually repulsive force that acts between neighbouring worlds, which will be seen to be responsible for a number of generic quantum effects in Sec. IV. Note also that U N is invariant under translation of the configurations. This leads directly to an analog of the quantum Ehrenfest theorem, also shown in Sec. IV.
While of a simple form, the interworld potential U N in Eq. (22) is seen to be a sum of three-body terms, rather than the more usual two-body interactions found in physics. Further, the corresponding force, acting on configuration x n , may be evaluated as where σ n (X) actually only involves three worlds: (defining x n = −∞ for n < 1 and = ∞ for n > N ). It is shown in Appendix A that this force corresponds to a particular approximation of the Bohmian force r t (x n ), as per Eq. (9).

IV. QUANTUM PHENOMENA AS GENERIC MIW EFFECTS
Here we show the MIW approach is capable of reproducing some well-known quantum phenomena, including Ehrenfest's theorem, wavepacket spreading, barrier tunneling, zero-point energy, and a Heisenberg-like uncertainty relation. We primarily use the toy model of Sec. III for this purpose, although a number of the results, such as the Ehrenfest theorem and wavepacket spreading, are shown to hold for broad classes of MIW potentials.

A. Ehrenfest theorem: translation invariance
where ϕ(x, p) N := N −1 N n=1 ϕ(x n , p n ). Thus, the quantum Ehrenfest theorem corresponds to the special case N → ∞, with the more general result holding for any value of N . For example, when V (x) is no more than quadratic with respect to the components of x, the centroid of any phase space trajectory follows the classical equations of motion, irrespective of the number of worlds.
To prove the first part of Eq. (25), note from Eqs. (15) and (17) that To prove the second part, note one also has from Eqs. (15) and (17) that

Now, translation invariance is the condition
and taking the gradient thereof with respect to y at y = 0 yields N n=1 ∇ xn U N (X) = 0. The second part of Eq. (25) immediately follows.

B. Wavepacket spreading: inverse-square scaling
It is well known that the position variance of a free one-dimensional quantum particle of mass m, with initial wave function Ψ 0 , increases quadratically in time [19]: 2m .
Here E Ψ and p Ψ denote the average energy and momentum, respectively, and Cov Ψ (x, p) denotes the position and momentum covariance, Ψ|(xp + px)/2|Ψ − x Ψ p Ψ . This is often referred to as the spreading of the wavepacket [19].
Here we show that an equivalent result holds for N one-dimensional worlds, for any interworld potential U N satisfying translation invariance and the inverse-square scaling property such as the toy model in Eq. (22). In particular, the position variance for N onedimensional worlds at time t may be written via Eq. (12) as Differentiating with respect to time then gives, using Eq. (15) and (17) with V (x) ≡ 0 and noting that (d/dt) p N = 0 from the Ehrenfest theorem (25), Now, differentiating scaling condition (27) with respect to λ, and setting λ = 1, gives N n=1 x n (∂U N /∂x n ) = −2U N . Hence, recalling that the average energy per world is E N = N −1 H N , one obtains Finally, integration yields the variance at time t to be which is of precisely the same form as the quantum case above. Like the generalised Ehrenfest theorem (25), this result holds for any number of worlds N .
The spreading of variance per se is due to the repulsive interaction between worlds having close configurations (see previous section). The above result demonstrates that a spreading which is quadratic in time is a simple consequence of an inverse-square scaling property for the interworld potential. Such a scaling is generic, since the quantum potential in Eq. (7), that is replaced by U N in the MIW approach, also has this property. Note that it is straightforward to generalise this result to configuration spaces of arbitrary dimension K, with V N replaced by the K × K tensor X · X T N − X N · X T N . Here the transpose refers to the configuration space index k, while the dot-product refers to the world-index n as previously.

C. Barrier tunneling: mutual repulsion
In quantum mechanics, a wavepacket incident on a potential barrier will be partially reflected and partially transmitted, irrespective of the height of the barrier. The probabilities of reflection and transmission are dependent on the energy of the wavepacket relative to the energy of the barrier [19]. In the MIW formulation, the same qualitative behaviour arises as a generic consequence of mutual repulsion between worlds having close configurations. This is investigated here analytically for the toy model in Sec. III, for the simplest possible case of just two worlds, N = 2.

Nonclassical transmission
Consider two one-dimensional worlds, the configurations of which initially approach a potential barrier V (x) from the same side, with kinetic energies less than the height of the barrier (Fig. 1). In the absence of an interaction between the worlds, the configurations will undergo purely classical motion, and so will be unable to penetrate the barrier. However, in the MIW approach the mutual repulsion between worlds will boost the speed of the leading world. This boost can be sufficient for this world to pass through the barrier region, with the other world being reflected, in direct analogy to quantum tunneling.
To demonstrate this analytically, consider the toy model in Sec. III for N = 2. The Hamiltonian simplifies via Eqs. (17) and (22) to where V (x) = 0 outside the barrier and has a maximum value V 0 > 0 in the barrier region. For simplicity we will restrict to the case of equal initial velocities v 0 in the direction of the barrier, with 1 2 mv 2 0 < V 0 (Fig. 1). Thus, in the classical limit = 0, penetration of the barrier is impossible.
Defining relative and centre-of-mass coordinates by q := x 2 − x 1 ,q := (x 2 + x 1 )/2, with conjugate momenta p = (p 2 − p 1 )/2 andp = p 1 + p 2 , respectively, then while the configurations of each world are outside the barrier region their evolution is governed by the Hamiltonian with initial conditions p(0) = 0 andp(0) = M v 0 , where M = 2m and µ = m/2. Thus, the equations of motion for q andq decouple, and may be integrated to givẽ where the 0 subscript indicates initial values. These are valid up until one of the configurations reaches the barrier region. The second equation may be further integrated to give the exact solution Note that this increasing mutual separation is in agreement with Eq. (28), with V = q 2 /4, and is solely due to the inverse-square interaction between the worlds in Eq. (29). It follows from Eq.(30) that the maximum possible separation speed isq ∞ := /(mq 0 ), which can be arbitrarily closely reached if the particles start sufficiently far to the left of the barrier. In this case, the speed of the leading particle is well approximated bẏ at the time it reaches the barrier region. Hence, transmission through the barrier is always possible providing that the corresponding kinetic energy is greater than whrre v classical is the initial speed required for a classical particle to penetrate the barrier. Eq. (32) clearly demonstrates that one world configuration can pass through the barrier if the initial separation between the worlds, q 0 , is sufficiently small. In terms of energy conservation, such a small separation gives rise to a correspondingly large interworld potential energy in Eq. (29), which is converted into a kinetic energy sufficiently large for barrier transmission. Further, the second configuration will suffer a corresponding loss of kinetic energy, leading to its reflection by the barrier.
It is seen that even the simplest case of just two interacting worlds provides a toy model for the phenomenon of quantum tunneling, with an interpretation in terms of the energy of repulsion between close configurations. While this case can be treated analytically, the interactions between N > 2 worlds are more complex. It would therefore be of interest to investigate the general case numerically, including tunneling delay times.

Nonclassical reflection
The same toy model also captures the nonclassical phenomenon that a barrier can reflect a portion of a quantum wavepacket, even when the incident wavepacket has a large average kinetic energy [19]. In particular, consider the case where the initial kinetic energies of both worlds are larger than the barrier height for the 'toy' Hamiltonian in Eq. (29), i.e., mv 2 0 /2 > V 0 . Hence, in the classical limit = 0 both configurations will pass through the barrier region.
It follows from the above analysis that the leading configuration is always transmitted -the kinetic energy of this world is only further increased by the mutual repulsion energy. However, choosing the initial distance from the barrier to be large enough for the separation speed to approachq ∞ , it follows that the second configuration will not be transmitted if Indeed, it will not even reach the barrier region if the left hand side is less than zero. This result demonstrates the converse to Eq. (32): it is always possible for one configuration to be reflected, if the intial separation between the worlds is sufficiently small, as a consequence of being slowed by the mutual repulsion between worlds.

D. Zero point energy: mutual exclusion
The classical groundstate energy of a confined system corresponds to zero momentum and a configuration that minimises the classical potential V (x). However, the corresponding quantum groundstate energy is always greater, with the difference referred to as the quantum zero point energy.
For example, for a one-dimensional quantum system the Heisenberg uncertainty relation (∆x) Ψ (∆p) Ψ ≥ /2 implies for the average energy of any state Ψ (noting that p 2 Ψ ≥ Var Ψ p). Hence, for any state with finite average energy, the system cannot be confined to a single point, as this would imply Var Ψ x = 0. In particular, it cannot be confined to the classical groundstate position x min , corresponding to the classical groundstate energy V min = V (x min ), and therefore E Ψ > V min .
In the MIW formulation, a zero point energy similarly arises because no two worlds can be confined to the same position -and to the classical groundstate configuration x min in particular. This mutual exclusion of configurations is a consequence of the repulsion between worlds, which forces the interworld potential U N in Eq. (18) to be strictly positive. For example, as shown in Appendix B, the average energy per world for the toy model in Sec. III satisfies This is clearly very similar to the quantum bound in Eq. (34), with the latter being precisely recovered in the limit N → ∞.
Remarkably, both the quantum and the MIW lower bounds are saturated by the groundstate of a onedimensional oscillator. The quantum case is well known [19]. For the toy model the corresponding groundstate energy is as demonstrated in Appendix C. Note that this vanishes in the classical limit = 0 (or, alternatively, N = 1), and approaches the quantum groundstate energy 1 2 ω in the limit N → ∞, as expected.

E. Heisenberg-type uncertainty relation
It was noted in Sec. II C that the interworld potential in Eq. (18) has the form of a sum of nonclassical kinetic energies, where the corresponding 'nonclassical momentum' of the nth world is given by p nc n := g N (x n ; X), with the components of g N explicitly defined in Eq. (19). For the toy model in Eq. (22), this nonclassical momentum has the explicit form (defining x 0 = −∞ and x N +1 = ∞ as always). It follows immediately that the average of the nonclassical momentum vanishes, i.e., Hence, the energy bound for U N in Eq. (B1) can be rewritten in the rather suggestive form reminiscent of the Heisenberg uncertainty relation for the position and momentum of a quantum system. Indeed, a nonclassical momentum observable p nc Ψ may also be defined for quantum systems, which satisfies the uncertainty relation ∆ Ψ x ∆ Ψ p nc Ψ ≥ /2, and which implies the usual Heisenberg uncertainty relation [20]. Thus, Eq. (38) has a corresponding quantum analogue in the limit N → ∞. It would be of interest to further investigate the role of the nonclassical momentum in the MIW approach.

V. SIMULATING QUANTUM GROUNDSTATES
The MIW equations of motion (15) and (17)  for all n. In particular, the forces acting internally in any world are balanced by the forces due to the configurations of the other worlds. Unlike quantum systems, the number of stationary states in the MIW approach (with N finite) is typically finite. For example, for two one-dimensional worlds as per Eq. (29), with a classical potential V (x) symmetric about x = 0, Eq. (39) implies that the stationary configurations are given by x 1 = −a and x 2 = a, where a is any positive solution of V (a) = 2 /(16ma 3 ). For the particular case of a harmonic oscillator potential, V (x) = (1/2)mω 2 x 2 , there is just one stationary state, corresponding to a = 1 2 ( /mω) 1/2 . More generally, the number of stationary state configurations will increase with the number of worlds N .
The MIW formulation suggests a new approach for numerically approximating the groundstate wave function and corresponding groundstate energy for a given quantum system. In particular, the groundstate probability density is approximately determined by finding the global minimum of N n=1 V (x n ) + U N (X), for suitably large N . Below we will give and test a 'dynamical' algorithm for doing so. First, however, for the purposes of benchmarking this algorithm, we calculate the exact groundstate configurations for the toy model of Sec. III, for a harmonic oscillator potential.

A. Oscillator groundstates: exact MIW calculation
For the one-dimensional toy model defined in Sec. III, with harmonic oscillator potential V (x) = (1/2)mω 2 x 2 , it is convenient to define the dimensionless configuration coordinates As shown in Appendix C, the unique groundstate configuration for N worlds is then determined by the recurrence relation subject to the constraints These equations are straightforward to solve for any number of worlds. As discussed in Appendix C, for small values of N , the recurrence relation can be solved analytically, and in general it can be efficiently solved numerically for any given number of worlds N . For example, for N = 11 one obtains the configuration depicted in Fig. 2. The corresponding Gaussian probability density for the exact quantum groundstate is also plotted, as the smooth curve in Fig. 2. It is seen that the groundstate configuration of just 11 worlds provides a good approximation to the quantum groundstate for many purposes. Note also that the corresponding mean groundstate energy follows from Eq. (36) as (5/11) ω ≈ 0.45 ω, which is reasonably close to the quantum value of 0.5 ω. The approximation improves with increasing N .  20) is expected to approximate the quantum groundstate distribution for a one-dimensional oscillator, PΨ 0 (ξ) = (2π) −1/2 e −ξ 2 /2 , for large N . The latter distribution is plotted for comparison (smooth magenta curve). All quantities are dimensionless.

B. Oscillator groundstates: dynamical MIW algorithm
The net force acting on each world is zero for a stationary configuration of worlds, as per Eq. (39). We exploit this fact in the following algorithm to compute a stationary configuration.
Given an arbitrary configuration of N worlds X(0) we iterate the following two-step algorithm: 1) Set the ve-locitiesẊ(0) to zero; 2) Integrate Eq. (11) over a small time interval [0, ∆t], and replace the initial configuration X(0) by X(∆t). Note that each iteration tends to reduce the total energy from its initial value, as the worlds move away from maxima of the total potential energy, and the velocities are reset to zero in each iteration. The fixed-points of this iterative map are clearly the stationary states of the MIW equations. Hence, after a sufficient number of iterations, the configuration will converge to a stationary configuration. An alternative algorithm could, instead of setting velocities to zero, add a viscous term to the evolution, in analogy to the dBB-trajectory based algorithm of Maddox and Bittner [21].
We have tested this dynamical algorithm for the case of the one-dimensional harmonic oscillator. We slightly modified the interworld force in Eq. (23) by placing auxiliary worlds in fixed positions on either side of the configurations x 1 < x x · · · < x N , with two auxiliary worlds to the far left and two to the far right, rather than at −∞ and ∞. The auxiliary worlds have only have a tiny effect on the computations, but ensure the interworld force is well-defined for n = 1, 2 . . . , N . We took advantage of Convergence of the algorithm was found to be extremely rapid; see Figs 3 and 4. For example, on an laptop running a simple Mathematica implementation of the algorithm, it took only 30 seconds to converge to the groundstate configuration for the case of N = 11 worlds, with a corresponding groundstate energy accurate to one part in 10 10 in comparison to the exact groundstate energy (5/11) ω following from Eq. (36).
This dynamical algorithm can be employed to compute ground states for all Hamiltonians for which the groundstate wave function has a constant phase. The numerical advantage over contemporary techniques to compute groundstate wave functions is evident when considering higher dimensional problems, i.e., K = DJ > 1. Instead of solving the partial differential equation HΨ g = E g Ψ g that lives on a K-dimensional configuration space (which memorywise starts to become unfeasible already for three particles in three dimensions), one only needs to solve for the stationary state of K × N coupled ordinary differential equations. The quantum distribution |Ψ g | 2 and energy E g of the ground state can then be approximated from the stationary distribution of worlds.

VI. SIMULATING QUANTUM EVOLUTION
Given the apparent success of the MIW approach as a tool for simulating stationary quantum states (Sec. V), it is of considerable interest to also investigate whether this approach can provide a similarly useful controlled approximation of the Schrödinger time evolution and, in particular, whether it is capable of describing quantum interference, which is one of the most striking quantum phenomena. Here we offer evidence that it can, at least for the canonical 'double-slit' problem.
In its simplest form, the 'double-slit' scenario comprises the free evolution of a 1-dimensional wave function Ψ t for an initial value Ψ 0 given by a symmetric superposition of two identical separated wavepackets. We chose real, Gaussian wavepackets with spread σ = 1 and initial separation of 4 units; see Fig 5. As described in Sec. II B, the corresponding initial values for the MIW approach are N worlds, X(0) = (x 1 (0), . . . , x N (0)), distributed according to |Ψ 0 | 2 with zero initial velocities. In order to numerically integrate the MIW equations of motion Eq. 11 to find the corresponding world-particle trajectories X(t) we employed a standard Verlet integration scheme [22]. In our example we used N = 41. The computed trajectories X(t) are shown in the top two plots of Fig. 6 and Fig. 7  initial world configuration x i for i = 1, 2, . . . , N . These trajectories are shown in two bottom plots again using a linearly smoothed density and histogram plot. Note that due to equivariance, the Bohmian configurations x dBB i (t) are distributed according to |Ψ t | 2 for all times t as the x dBB i (0) are distributed according to |Ψ 0 | 2 . Therefore, the bottom density plot in Fig. 6 is a discretized plot of |Ψ t | 2 as shown in Fig. 5 (this discretization, corresponding to using just N = 41 Bohmian worlds and linear smoothing, is responsible for the blurriness relative to Fig. 5).
By comparison of the plots we conclude that the MIW approach is at least qualitatively able to reproduce quantum interference phenomena. Furthermore, it should be stressed that with respect to usual numerical PDE techniques on grids the MIW approach could have several substantial computational advantages: As in dBB, the worlds X(t) should remain approximately |Ψ t | 2 -distributed for all times t which implies that without increasing N , the density is naturally well-sampled in regions of high density even when the support of Ψ t becomes very large. The computational effort is then automatically focussed on high-density regions of configuration space, i.e., on the crucial regions where numerical errors have to be minimized to ensure convergence in the physically relevant L 2 norm. Furthermore, transmission and reflection coefficient can be computed simply by counting worlds. These advantageous also hold true for methods of integration the Schrödinger equation on co-moving grids as, e.g., proposed by Wyatt et al. [23]. However, there are reasons to expect that the MIW approach may be numerically more stable (see Ref. [24] for a discussion of such numerical instabilities in conventional approaches) and less performance intensive (as it is merely an ODE initial value problem).
There is much work to be done on the question of time evolution using our MIW approach. Even in the example above, the convergence as a function of N , the number of worlds, remains to be investigated. Also, there is the issue of wave functions with nodes (or even finite regions of zero values). These allow for different wave functions Ψ 0 (q) having identical P 0 (q), and identical ∇S 0 (q) everywhere that P 0 (q) = 0. From Sec. II B these different wave functions will all generate the same initial conditions for our ensemble of many worlds, for any finite value of N , and thus are indistinguishable in the MIW approach, contrary to quantum predictions (see also Sec. 10 of [12]). However it is plausible that restricting the moments of the quantum energy to be finite would typically eliminate all, or all but one, of such wavefunctions, by enforcing the continuity of ∇Ψ 0 (and thus also of the dBB velocity field) [25,26]. Finally, beyond the toy model, there is the question of whether the general approach of many interacting worlds allows successful numerical treatments of more general situations, with common potentials and K = DJ > 1.

VII. CONCLUSIONS
We have introduced an approach to quantum phenomena in which all quantum effects are due to interactions between a large but finite number N of worlds, and probabilities arise from assigning an equal weighting to each world. A number of generic quantum effects, including wavepacket spreading, tunnelling, zero point energies, and interference have been shown to be consequences of the mutual repulsion between worlds (Sec. IV). This alternative realistic interpretation of quantum phenomena is also of interest in not requiring the concept of a wave function. The wave function does not exist other than as an epiphenomenon in the notional limit that the initial distribution of worlds approaches |Ψ 0 (x)| 2 and the initial velocity of each world approaches Eq. (3) as N → ∞.
For finite N , our MIW approach can only ever give an approximation to quantum mechanics, but since quantum mechanics is such an accurate theory for our observations, we require this approximation to be very good. In this context, it is worth asking the following question: what restrictions must one place on the initial distribution of world configurations and velocities so that nearly quantum behaviour (that is, observations consistent with quantum mechanics for some wavefunction) will be experienced by macroscopic observers in almost all worlds? The answer may depend on the details of the inter-world potential. It may be necessary to impose the quantum velocity condition of Eqs. (3) and (4) using some smoothly varying (on the scale of the separation between nearby worlds) single-valued complex function Ψ 0 on configuration space. But it is conceivable that other (perhaps even most, by some measure of typicality) initial conditions would relax, under our interacting-world dynamics, to conditions approximating quantum behaviour, at least at the scale which can be probed by a macroscopic observer. This is one of the big questions which remains to be investigated.
Other matters for future investigation include how spin and entanglement are modelled in the MIW approach. The latter will require studying the case of worlds with configuration spaces of at least two dimensions (corresponding to two one-dimensional systems), and will also allow analysis of the quantum measurement problem (where one system acts as an 'apparatus' for the other). This may help clarify the ontology and epistemology of any fundamental new theory based upon the MIW approach to quantum mechanics.
In the context of entanglement, it is worth comparing our MIW approach with conventional many-worlds approaches. The latter are often motivated by the desire, first, to restrict reality to only the wave function, and second, to avoid the explicit nonlocality which arises from entanglement in other realist versions of quantum mechanics. Our approach is, by contrast, motivated by the desire to eliminate the wave function. It furthermore elevates the nonlocality of quantum mechanics to a kind of "super-nonlocality": particles in different worlds are nonlocally connected through the proposed MIW interaction, thus leading, indirectly, to nonlocal interactions between particles in the same world.
Turning from questions of foundations and interpetations to applied science, the MIW approach provides a promising controlled approximation for simulating quantum groundstates and the time-dependent Schrödinger equation, as discussed in Secs. V and VI. In particular, we have shown that it is able to reproduce quantum interference phenomena, at least qualitatively. Quantitative comparisons with different initial conditions, convergence as a function of the number of worlds N , and generalizations to higher dimensions, is a matter for immediate future work.
the Hamiltonian where Q(q) is the quantum potential corresponding to Ψ t (q) [identifying the momentum component p k with the right hand side of Eq. (3)]. Hence, the evolution of N Bohmian world-particles, each guided by Ψ t (q), are precisely described via the time-dependent Hamiltonian Thus, if Q t (q) (and its gradient) can be approximated by some functionQ(q; X t ), depending on the configurations X t of the N worlds (assumed to sample |Ψ t (q)| 2 ), then one immediately has a suitable corresponding Hamiltonian formulation of the MIW approach, with an interworld potential U N in Eq. (17) of the alternative form to that in Eq. (18).
To show how the above form is related to the positive definite form in Eq. (18), consider the relation dx |Ψ(q)| 2 Q(q) = dq |Ψ(q)| 2 K k=1 2 8m k P 2 ∂P ∂q k 2 following from integration by parts and Eq. (7) [2]. In particular, since the N world-particles have configurations sampled according to |Ψ(q)| 2 , it immediately follows that Hence, U N in Eq. (A1) may alternatively be replaced by the form in Eq. (18), corresponding to a suitable approximation of P and its derivatives in terms of the configurations of the N worlds. Here "suitable" means that we must check that the corresponding MIW force is given by r N (x n ; X) = −∇ xn U N (X) in the limit N → ∞.
For the toy MIW model of Sec. III, we can demonstrate directly that the force in Eq. (23) is an approximation of the Bohmian force, as per Eq. (9) of Sec. II. To so so, we first note that, for the one-dimensional case, the derivative of any regular function ϕ(x) can be approximated at x = x n by either of for sufficiently large N . Hence, using Eq. (20), one obtains where the difference operators D ± are defined by D + ϕ n := ϕ n+1 − ϕ n , D − ϕ n := ϕ n − ϕ n−1 for any sequence {ϕ n }. Now, as may be checked by explicit calculation from the definition of the quantum potential in Eq. (7), for a one-dimensional particle the Bohmian force in Eq. (6) can be written in the form i.e., where the derivative operators act on all terms to their right. Applying Eq. (A2) then gives as desired, where the second line uses Eq. (20), and the last line follows by expansion and direct comparison with Eq. (23).

Appendix B: Lower bound for the interworld potential
To obtain the bound in Eq. (35), first let f 1 , . . . , f N and g 0 , g 1 , . . . , g N be two sequences of real numbers such that g 0 = g N = 0. It follows that  with equality if and only if f n = α(g n − g n−1 ) for some α. Combining these results then yields The particular choices f n = x n − x N and g n = 1/(x n+1 − x n ) (with x 0 = −∞ and x N +1 = ∞) immediately yield, using the definition of U N in Eq. (22), the lower bound (B1) for the interworld potential energy, with equality if and only if the Schwarz inequality saturation condition is satisfied. Eq. (35) then follows via Eq. (17), as desired. Similar bounds can be obtained for other choices of the interworld potential, but will not be discussed here.
Appendix C: Groundstate energy and configuration of toy model oscillator To derive the groundstate energy in Eq. (36), for the 1D toy model oscillator, note first that for a harmonic oscillator potential V (x) = (1/2)mω 2 x 2 one has The second line follows via x 2 N ≥ V N and inequality (B1), and the last line via z + 1/z ≥ 2.
To show this inequality chain can be saturated, as per Eq. (36), the conditions for equality must be checked at each step. It follows that the lower bound is achievable if and only if (i) the momentum of each world vanishes, i.e., p 1 = · · · = p N = 0; and (ii) the configuration satisfies both and the Schwarz inequality saturation condition from Eq. (B2), for some constant α (using x N = 0). These conditions generate a 2nd-order recurrence relation with fixed boundary conditions, yielding a unique solution for any given number of worlds N . The corresponding unique groundstate configurations converge to the quantum Gaussian groundstate wave function for N → ∞, as discussed in Sec. V A.
To determine the groundstate configuration, for a given number of worlds N , note that summing each side of Eq. (C2) from n = 1 up to any n < N gives 1 x n+1 − x n = x 1 + · · · + x n α .