Random Functions via Dyson Brownian Motion: Progress and Problems

We develope a computationally efficient extension of the Dyson Brownian Motion (DBM) algorithm to generate random function in C2 locally. We further explain that random functions generated via DBM show an unstable growth as the traversed distance increases. This feature restricts the use of such functions considerably if they are to be used to model globally defined ones. The latter is the case if one uses random functions to model landscapes in string theory. We provide a concrete example, based on a simple axionic potential often used in cosmology, to highlight this problem and also offer an ad hoc modification of DBM that suppresses this growth to some degree.


Introduction
Generating random functions locally via Dyson Brownian Motion (DBM) was first introduced in [1] following prior work by Dyson [2], see [3] for a textbook, and applied to potentials in cosmology in [1,4,5]. The algorithm starts with a Taylor expansion of the potential to second order and, after moving a set distance from the initial expansion point, adding random perturbations to the matrix of second derivatives 1 . Thus, a random function is generated locally along a trajectory. The algorithm as introduced in [1] generates a function in C 1 along a trajectory 2 with a Hessian in the Gaussian orthogonal ensemble. In [4], the method was generalized to generate functions in any differentiability class by relegating the perturbations to a higher order derivative tensor, with particular emphasis on functions in C 2 . Such functions are of particular interest to model the potential in inflationary cosmology, see [7,8] for recent reviews, while enabling numerical computation of the power-spectrum of cosmological fluctuations 3 without generating ringing or requiring a step size for DBM in line with the one used for the discretization of the required differential equations. However, the methods introduced in [4] are computationally more demanding than ordinary DBM, since a coordinate rotation is needed at each step.
In this brief technical note we improve upon the algorithm presented in [4] to generate a function in C 2 without the need to rotate the Hessian, see Sec. 2.2. As a consequence, field spaces with more than 100 fields, as needed to model certain landscapes in string theory/cosmology, can be handled. We also explain a severe shortcoming of functions generated via DBM or its generalization: while the elements of the Hessian have a stable distribution, the function itself performs a random walk and is thus unbounded. Simply using the appropriate distribution for a Hessian in a bounded potential, such as the axionic one investigated in [9], instead of one in the Gaussian orthogonal 1 Dyson Brownian motion [2] was developed primarily to introduce a "time" dependence into Random Matrix Theory, which enables the computation of transition probabilities. Analytic methods to compute such probabilities were derived in [6]. 2 The continuum limit of the algorithm yields a function in C 2 . If implemented numerically and the step length of the DBM function is kept in line with the discretization used elsewhere, the functions introduced in [1] are indistinguishable from functions in C 2 . However, using functions that are in C 2 by construction can be computationally advantageous. 3 In [5], semi-analytic methods were used to compute observables. ensemble does not alleviate the problem. This shortcoming renders DBM potentials ill suited to model many landscapes of interest in cosmology if one is interested in going beyond the first coherent patch. We explain this shortcoming in a concrete case-study in Sec. 3 and provide a minor modification of DBM that can suppress the unstable growth to some degree. We refer the interested reader to [1,4,5] for a more pedagogical introduction to DBM and cosmological implications, while focusing on technical aspects in this note.
2 Generating C 2 Random Potentials

Dyson Brownian Motion Potentials
To generate random functions, Dyson [2], see also [1], proposed to consider the Hessian to be a stochastic variable whose change after a step δs ≡ δφ a ≡ δφ a δφ a (2.1) along some path in a D-dimensional field space 4 with fields φ a , a = 1 . . . D is given by Here, Λ h is the horizontal correlation length of the Hessian, the δA ab are D(D + 1)/2 zero-mean stochastic variables and the term ∝ −H ab is a restoring force. Due to the latter, the Hessian is a (symmetric) matrix in the Gaussian Orthogonal Ensamble (GOE, Wigner matrix). The first two moments of the Hessian obey [1] where σ is the standard deviation of the GOE. To apply DBM and construct a random potential 5 , consider the Taylor expansion at some initial point p 0 where Λ v sets the vertical scale andφ a ≡ φ a /Λ h . Both, Λ v and Λ h , have units of mass. Truncating the series at second order and evaluating the potential at a nearby point p 1 a distance δφ a 1 from the initial point p 0 away, we get where δv ab | p 0 are elements of a random matrix with Repeated application along some path results in a locally defined, random potential V ∈ C 1 . See [1] for a more detailed discussion.

Extended Dyson Brownian Motion Potentials, V ∈ C 2
In [4], the DBM algorithm was extended to generate random potentials V ∈ C k with k = 2, 3, . . . , by relegating perturbations to higher order derivative tensors. Two methods were introduced in [4] for V ∈ C 2 , both computationally demanding, since they require a rotation of the coordinate system after each step. Thus, an application to large D is cumbersome. Here, we present an algorithm that doesn't require any coordinate transformations and is therefore efficient enough to generate random potentials with D = O(100). Let us start by expanding the potential to third order, The potential at a nearby point p 1 close to p 0 can thus be written as Note that the random variable enters now in the tensor of third derivatives, v abc . However, we still want the Hessian at well seperated points to act as a random variable with the same first two moments as in (2.7) and (2.8). To this end, we define the δv abc as random variables as in [4], leading to δv ab | pn = v abc | pn δφ c (2.14) Here, v abc and thus δv abc are symmetric under permutations of a, b and c. Note that the means and variances of δv abc appear in a sum, which lead to the use of rotations in [4] to either align one coordinate axis with the preceding step vector, or to diagonalise the Hessian and extract conditions on the means and variances of δv abc . Here, we make use of the fact that the above system of equations is under-determined to provide a more efficient algorithm.
For the D dimensional case, δv abc is a D × D × D totally symmetric tensor. We determine the mean and variance for each element separately, using the constraints given in (2.17) and (2.21) according to the following algorithm: 1. Set a = 1 and determine the first D × D elements in δv 1bc . Set the off-diagonal entries of this matrix to zero and obtain a unique value for the diagonal entries using (2.17) and (2.21). 4. If the step length in a certain direction is too small, instabilities can arise due to numerical errors. Besides increasing the computational precision, and thus the computational cost, one may also include a lower limit for the path in each coordinate direction, δφ limit . If δφ c < δφ limit , one may simply set δφ c ≡ δφ limit and avoid any computational problems. A plausible choice is, for instance, δφ limit = δs/(10Λ h √ D). This lower limit will alter the trajectory slightly, which, depending on the question one wishes to investigate, can be acceptable.

Demonstration of the Algorithm
Let us demonstrate the above procedure for D = 3.

1.
We set all the off-diagonal elements of δv 1bc to zero and choose the diagonal elements via The mean and variance for δv 122 and δv 133 are chosen analogously.

2.
For δv 2bc , we set the first row and the first column in line with their symmetric partners. We set the off-diagonal elements of the remaining 2 × 2 sub-matrix to zero and solve (2.17) as well as (2.21) for the diagonal ones . 3. We use the symmetry conditions to determine the elements of the first two rows and columns in δv 3bc . The only remaining element is δv 333 , for which we get The above procedure is summarized in Tab. 1 and easily implemented numerically for any D.

Generating Axionic Random Potentials
One possible application of DBM potentials is the modelling of globally defined potentials in cosmology/string theory by locally defined ones. Such an approach is motivated by the need for numerical experiments in the large D limit. To this end, the locally defined potential needs to have identical statistical properties as the global one. In this section, we explain at the example of a simple axionic potential used in prior work [9][10][11] some of the necessary steps as well as a severe problem (see [7,8] for reviews on inflationary and axionic cosmology): potentials generated via DBM are not bounded, rendering them of limited use without further modification if the aim is to model a bounded potential. To be concrete, while the original DBM potentials reviewed in Sec. 2.1 have a bounded Hessian, they show a runway behaviour for the gradient and the potential as the traversed distance in field space increases. The reason for this known, yet hardly ever mentioned, behaviour is simple: both perform, in essence, a random walk, which causes on average a growth of the potential's and gradient's variance (∝ N 3 steps for the potential, as we shall see below). The same growth is present in all cases investigated in this note.
The free parameters of DBM are the vertical scale Λ v , the horizontal correlation length Λ h of the Hessian and the concrete probability distribution function of the Hessian's elements (taken to be Gaussian in the original proposal expained above). The vertical scale is easily read of any globally defined potential and the potential's horizontal correlation length follows relatively straightforward as well. However, since the latter is somewhat less trivial to compute than the former and usually not identical with Λ h , we provide an examble in Sec. 3.1 for the axionic potential investigated in [9].

The Correlation Length
Consider an axionic potential defined via [9] V =ñ where D is the number of axions, n =ñ D is the number of shift symmetry breaking sources from non-perturbative effects, Λ J 1 ,...,J D sets the strength of the source terms, φ j denote the axions with decay constants f j that are included in the mixing matrix n ij which are set deterministically as n ij ≡ 3J j /ñ to enable comparison with analytic results in [9], and θ J 1 ,...,J D are the relative phases between different source terms. In line with [9], we let J j run from 1 toñ. Since results do not depend on the random phase factors, we set them to zero. For the overall scale of the potential, we choose a normal distribution for the Λ J 1 ,...,J D with mean and standard deviation Hence, the potential at the global minimum is V = 0 and at the maximum V = 2. To define a correlation function, let us introduce the normalized variable Since < V ( φ) >= 1, the above definition becomes f = V ( φ) − 1. A two point correlator can be defined as Since the potential and thus f is isotropic, C only depends on the distance r in field space. Note that C → 1 as r → 0 and C → 0 for large r. It is useful to introduce the Fourier space "power spectrum" We leave the normalization constants A and B unspecified, since they cancel out in the correlator C anyhow. Due to δ D ( k + k ) in (3.8), we get which can be simplified to Note that n ij = 3J j /ñ does not depended on the index i. Inserting the above expression into (3.9) and evaluating the integral yields n ij · r j ) > , (3.14) so that the correlator becomes The correlator may be expanded in a Taylor series if the correlation length is small, but it is simple to evaluate the sum numerically anyhow. Evidently, the potentials correlation length for the axionic case under consideration is independent of number of fields D and the number of shift symmetry breaking termsñ. Evaluating the correlation length numerically for the above mentioned distributions of parameters yields the ensemble average Λ h (V ) ∼ 0.75. However, Λ h (V ) is in general not identical with the correlation length of the Hessian, Λ h (V ab ). One may either adjust the latter iteratively in a concrete DBM potential to match the desired Λ h (V ) of the globally defined potential or compute Λ h (V ab ) semi-analytically, as we did with Λ h (V ) above, to match it directly. It should be noted that it is usually not possible to match the potential's and the Hessian's correlation length simultaneously, since only one free parameter can be adjusted.

DBM Potentials with Adjusted Mean of the Hessian's Entries
The entries of the Hessian for the axionic potential in (3.1) are, to a good approximation, Gaussian random variable with mean and standard deviation see [9] for details. Here, V c is the value of the potential at the point in field space at which the Hessian is to be evaluated. Note that the term 1 − V c makes maxima more likely for V > 1 and minima for V < 1. So, instead of the condition in (2.7) Enforcing (3.19) guarantees that the entries of the Hessian acquire the mean in (3.16). We do not impose the different standard deviations in (3.17). Random potentials generated this way show a common feature, namely, they oscillate with a growing amplitude, as shown in Fig. 1. After a few correlation lengths, the potential grows to an unreasonable height, far beyond the bounds of the globally defined potential. Let's see if we can pinpoint the origin of this feature: Our starting point for DBM is a small v ab | p 0 with δφ a /Λ h = 0.01 and Λ 4 v √ D = 1. Thus, the first term 9(1 − v 0 )/ñ 2 contributes much more than the second term, −v ab . The second term tries to compensate the effect of the last step while the first term changes its sign (and thus whether to push the potential upwards or downwards 6 ) after crossing v 0 = 1. Since the second term is smaller the first one, the DBM algorithm will always tend to push the potential downward(upwards), if v 0 is above(below) one. Thus, if we start with v 0 > 1, a large gradient results at v 0 = 1. After crossing 1 from above, the DBM algorithm tries to push the potential upwards, but the large gradient at V = 1 needs to be compensated in order to turn around, which in turn takes a longer distance the higher the magnitude of the gradient is. Therefore, once the potential reaches the next minimum, it is at a lower value than before which causes an even larger gradient at the next crossing of v 0 = 1. This mechanism leads to the observed oscillations with a growing amplitude.
The growth of V can be quantified by the variance of the potential as a function of N steps . For example, a numerical study to mimic the axionic landscape with D = 3,ñ = 5 via a DBM yields The dependence ∝ N 3 steps is generic. Thus, such DBM potentials do a poor job at modeling the axionic landscape: while the Hessian has the proper distribution, the potential and the gradient are completely dominated by these unstable oscillations. The same artifact is present for the Extended Dyson Brownian Motion algorithm of Sec. 2.2. Thus, such potentials can not be used to draw any conclusion with respect to the distribution of extrema or the likelihood of achieving cosmological inflation on axionic landscapes. It is possible to modify the DBM-algorithm to keep the potential under control for a while by increasing the magnitude of the second term and thus making the instantaneous back reaction stronger. However, the unstable growth of the potential's variance ∝ N 3 steps is still present, as in all DBM potentials considered in the literature, see the next section.

Modifying DBM, a Case Study
Let us consider the ad hoc introduction of a free parameter into the mean of δv ab , while keeping the variance unchanged, For appropriately chosen values of A, it is usually possible to reduce the numerical pre-factor in the growth of the variance considerably, but σ(v ab ) is altered as well; furthermore, the dependency on ∝ N 3 steps remains. To illustrate this point, let us try to mimic an axionic landscape with D = 3,ñ = 5.
In line with (3.19), we introduce the free parameter A into the mean of δv ab , δv ab | pn = 9 To compensate for the alteration of the Hessians variance, we allow for two additional free parameters, B and C, We first optimize A to match the desired standard deviation of the potential for a longer stretch in field space: starting from an initial value A = 20 we increase A in steps of one and generate 60 realizations for each case. All realizations run for 600 steps. We demand that v 0 lies in the interval [0, 2]. Once v 0 exceeds this range, we delete this realization and start a new run. We follow with the computation of the potential's standard deviation for the 600 steps, and calculate the mean value of the 60 realizations. We stop searching once a suitable value of A is found, for which we demand that the averaged standard deviation differs from the potential's desired one (see [9]) no more than 0.01, i.e. about 50% relative deviation. The smallest value of A found via this algorithm is A = 36 in our particular case.
Holding A fixed, we vary B and C at the same time. We start with B = 3, C = 2 and compute both, σ(V ) and σ(v ab ). For their deviations from the desired values we demand the difference to be smaller than 0.1 for σ(v ab ) and smaller than 0.05 for σ(V ). For a fixed B, we increase C in steps of one. If C grows over 5 and the two conditions haven't been satisfied, we increase B by one and continue with the same procedure. The set of parameter we found via this algorithm is A = 36 , B = 10 , C = 3 . (3.25) The above algorithm can of course be optimized, but we are primarily interested to provide an instructive example for which our method is sufficient. Three exemplary realizations are shown in Fig. 2. Evidently, the oscillations with increasing amplitude are not present any more. However, due to the still present unstable growth of the potentials variance with N 3 steps , one can not mimic the desired landscape for large step numbers, even with this modification. Furthermore, for low step numbers, the potential is dominated by the chosen initial conditions. To be concrete, the variance can be fitted well by Var(V ) ≈ −9 + 1.3 × 10 −9 N 3 steps , (3.26) see Fig. 3, showing clearly the dependence ∝ N 3 steps . Performing the same modification on the original DBM with mean and variance in (3.23) and (3.22) respectively, that is A = 36, B = 1 and C = 1, shows the same qualitative behaviour: while the prefactor is reduced by the choice of A, the proportionality with ∝ N 3 steps remains. A fit leads to Var(V ) ≈ −6 + 3 × 10 −9 N 3 steps , (3.27) see Fig. 4. The same unstable growth can be observed in C 2 -potentials generated via the algorithim introduced in Sec. 2. Similarly, the variance of the gradient is growing with the traversed distance in field space. Thus, due to this unstable growth, potentials generated via DBM have limited use if the goal is to mimic a globally defined one for an extended stretch in field space. For a local study, for example to study effects of a mild random component to the potential near a saddle point on which inflation takes place, DBM potentials can be used, see e.g. [5]. One needs to be cautious though, since the presence of a flat region suitable for inflation is due to ones initial conditions. Also, due to the growth of the gradient's variance, we expect that the end of inflation is strongly affected by the method by which the potential is generated, i.e. Dyson Brownian Motion.

Conclusions
We provided a computationally efficient extension of the Dyson Brownian Motion algorithm to generate random function in C 2 locally, as desirable for certain applications in cosmology. We further showed at the example of a simple globally defined potential that DBM potentials fail to recover basic features of the globally defined ones, due to the presence of an unstable growth of the gradient's and the potential's variance. We also showed that a minor ad hoc modification of the algorithm can weaken the unstable behaviour, such that a DBM potential can mimic a globally defined one reasonable well for some time. However, the unstable growth is still present. Thus, in their current form, DBM potentials can be used to model landscapes in String Theory or to draw conclusions in cosmological settings only for regions not much bigger than their coherence length.