Time-Symmetric Rolling Tachyon Profile

We investigate the tachyon profile of a time-symmetric rolling tachyon solution to open string field theory. We algebraically construct the solution of [arXiv:0707.4472] at 6th order in the marginal parameter, and numerically evaluate the corresponding tachyon profile as well as the action and several correlation functions containing the equation of motion. We find that the marginal operator's singular self-OPE is properly regularized and all quantities we examine are finite. In contrast to the widely studied time-asymmetric case, the solution depends nontrivially on the strength of the deformation parameter. For example, we find that the number and period of oscillations of the tachyon field changes as the strength of the marginal deformation is increased. We use the recent renormalization scheme of [arXiv:1412.3466], which contains two free parameters. At finite deformation parameter the tachyon profile depends on these parameters, while when the deformation parameter is small, the solution becomes insensitive to them and behaves like previously studied time-asymmetric rolling tachyon solutions. We also show that convergence of perturbation series is not as straightforward as in the time-asymmetric case with regular OPE, and find evidence that it may depend on the renormalization constants.


Introduction and Conclusions
In a boundary CFT, the boundary condition can be deformed on any section of the boundary by exponentiating a marginal operator integrated along it, as in e λ dt V (t) . (1) The marginal parameter λ controls the strength of the deformation. In Open String Field Theory, allowed D-brane configurations are in one to one correspondence with classical solutions. The rolling tachyon is the time-dependent solution which corresponds to a decaying D-brane. There are two rolling tachyon solutions obtained by different marginal deformations of the D-brane CFT. The simpler case, the exponential rolling tachyon, involves the marginal operator V (t) = e X 0 (t) and represents a D-brane which exists in the infinite past and then decays at a finite time. This case has been studied using level truncation methods [3,4] as well as analytically [5,6,7], and is relatively simple because the OPE of the marginal operator with itself is finite. The more difficult case uses the time-symmetric marginal operator V = √ 2 cosh(X 0 ), which has the singular self-OPE V (0)V (x) ∼ 1 x 2 . This rolling tachyon corresponds to placing a D-brane at t = 0 and letting it decay at both t = −∞ and t = +∞.
In SFT, the tachyon profile is the tachyon component of the string field as a function of time. In the symmetric case it has the form where β (j) n are coefficients which can be calculated numerically and n/2 is the greatest integer less than or equal to n/2, so that n−2j ≥ 0. In this notation the deformation strength λ is taken to be negative for physical solutions [5]. In the exponential case, λ controls the time at which the D-brane decays, while for the time-symmetric case it determines the lifetime of the D-brane, with longer lifetimes corresponding to λ closer to zero. In the regular OPE case, instead of the double sum and time-symmetric cosh functions, energy conservation tells us that there is only a single sum of exponentials involving coefficients with β (0) n : While the β (j) n are in general gauge-dependent, for one choice of gauge it was proven that the sum in the regular OPE case converges for all λ, with the asymptotic behaviour β (0) n ∼ e −γn 2 [7]. Numerical data suggests that this is true in other gauges as well, as shown in figure 3b. The trouble with this is that the tachyon profile itself exhibits wild oscillations which grow exponentially in magnitude, while the vacuum without any D-branes is a well defined and finite point in string field space. How these two very different looking string fields are reconciled has been the source of much speculation (see, for example, [8]). Our results confirm that the tachyon profile has the same growing oscillatory behaviour in the time-symmetric case, and do not appear to exclude any of the current hypotheses.
When studying the marginal operator V = √ 2 cosh(X 0 ) leading to the time-symmetric rolling tachyon, we must be careful to avoid singularities arising from the operator's OPE. Analytic solutions for marginal deformations require the insertion of many copies of the marginal operator with separations that are integrated over, and there will be divergences when two operators approach each other. Fortunately there are several solutions which are intended to handle this issue [1,9,10,11]. The most recent work, by Erler and Maccaferri, does not apply to solutions which have a non-trivial time direction, so we cannot use it for the rolling tachyon. Fuchs, Kroyter and Potting's solution was designed with the photon marginal deformation in mind, but it is possible that it could describe the rolling tachyon as well. The solution of [10] is a generalization of [7] to operators with singular OPE, and it could be applied to the rolling tachyon. In fact it has been suggested that this solution could give the tachyon profile in the form (9), which would help settle the convergence issue.
Our focus, however, will be on the work of Kiermaier and Okawa. They proposed a general construction dependent on the existence of a suitable renormalization scheme [1], which was investigated and refined in [2]. A general renormalization scheme satisfying the necessary conditions was shown in [2] to have at least two free parameters, suggesting that the tachyon profile could have free parameters as well. Here we will perform the first explicit numerical calculations for this solution, and we will show that the tachyon profile is a finite function which does in fact depend on the free parameters.
When we implement the solution Ψ of [1] with the renormalization scheme of [2] we can find the tachyon profile for the symmetric rolling tachyon. This involves algebraically constructing the wedge states with insertions corresponding to the solution, taking expectation values, and then performing the required integrals numerically. Here this is done up to 6th order in λ. It will have the form (2), where now the function is symmetric in t and all the β (j) n are non-zero. The marginal operator √ 2 cosh(X 0 ) contains the operators e ±X 0 with both signs, and the coefficients β (j) n correspond to terms with n − j factors of one of the two operators and j factors of the other. Since renormalization has the effect of adding counterterms for collisions of operators with opposite sign, the j = 0 coefficients involve no counterterms and behave very similarly to exponential solutions. These show the same β (0) n ∼ e −γn 2 asymptotic behaviour, implying that the sum ∞ n=1 λ n β (0) n cosh(nt) converges absolutely for all λ. For |λ| 1, as is the case when the D-brane survives for a long time, the j > 0 coefficients are suppressed due to extra factors of λ. This results in a decay process which looks very much like the regular case, as the decay is well separated from the "anti-decay" by the D-brane's lifetime. Once this lifetime is long enough, further shrinking λ even has the same effect on the decay time as it would with the exponential rolling tachyon, simply shifting the time of the decay.
For the β (j) n coefficients with j > 0, each coefficient is calculated using a number of counterterms determined by j. The counterterms in turn are functions of the two parameters of the renormalization scheme, C 0 and C 1 . The bulk coefficients are therefore polynomial functions of C 0 and C 1 . Because these coefficients are not constants, patterns such as the asymptotic behaviour for j = 0 could depend on the choice of C 0 and C 1 . Considering only the j = 1 coefficients, with C 0 = C 1 = 0 they are quite a good fit to β In fact there is no choice of those constants for which the exponential quadratic behaviour β (1) n ∼ e −γ 1 (n−2) 2 fits as closely. This suggests that the sum ∞ n=2 λ n β (1) n cosh((n − 2)t) also converges, but there may still be some choices of C 0 and C 1 for which this is not the case, or for which the radius of convergence in λ is finite. For example, choosing the constants so that β (1) n are a best fit to the exponential cubic behaviour results in β (2) n which are increasing, at least for the three coefficients we can calculate with j = 2.
So how does the inclusion of all the β (j) n coefficients affect the shape of the tachyon profile? We show that for |λ| 1 these coefficients are negligible, but as the D-brane lifetime is decreased there comes a point where more coefficients must be considered. Some terms cease to dominate for any range of time, and the number of oscillations actually decreases. The missing oscillation means that the effective "period" is significantly decreased. What this means physically is not clear, since the period is a gauge dependent quantity related to the coefficient γ in the exponent of the asymptotic behaviour. The tachyon profile for small λ is very similar to that of [5], while for large λ it has features similar to the tachyon profile of [7], so perhaps the solution is interpolating between regular-OPE solutions in different gauges as the marginal deformation strength is changed. Understanding this phenomenon is left for future work.
This paper is organized as follows. Section 2 is focused on the details of the tachyon profile, and contains discussion and plots of results. In section 3 we briefly discuss how the numerical calculations were performed, and demonstrate that our results converge to appropriate values. There we also discuss sources of roundoff error which can influence the calculations.

The Tachyon Profile
The solution of [1] presents a promising framework for construction of a time-symmetric rolling tachyon solution, but it was not applied to any specific marginal deformation. Taking that approach and inserting the marginal deformation V = √ 2 cosh(X 0 ) we are able to numerically compute the tachyon profile up to 6th order in the deformation parameter λ. Since the tachyon profile has previously been calculated for several exponential rolling tachyon solutions, we can compare our results in order to determine what qualitative differences appear in the time-symmetric case. It is also useful to have explicit numerical evidence that the renormalization scheme studied in [2] is effective and the solution remains finite despite the singular OPE that the marginal operator has with itself.
The solution takes the form of a wedge state with insertions on the boundary. While one insertion will always be at a fixed location, the rest are integrated. The renormalization procedure replaces pairs of operators with appropriate counterterms under the integral. Each operator V contains two terms carrying ±1 unit of "momentum" in the time direction, but the counterterms are functions and carry no momentum. In (2) the coefficient β (j) n clearly contains the part of the tachyon profile with n factors of λ and k = n − 2j units of this momentum, so it follows that the coefficients β (0) n contain no counterterms. This is as it should be since operators e ±X 0 with the same sign have a regular OPE; the singular OPE of the cosh(X 0 ) marginal operator comes entirely from the collision of exponentials with opposite sign. The index j, which counts the momentum deficit, also has the effect of counting the maximum number of counterterm factors. For the coefficients β (j) n , table 1 shows their values as calculated by the Cuhre algorithm, and with the exception of two terms we will use those coefficients. For technical reasons explained in section 3.3, the two terms marked with asterisks will use values found by the Suave algorithm instead, and those are shown in table 2. Occasionally we will want to think of the tachyon profile in terms of these timelike momentum modes, and write This form is equivalent to (2) as long as we define β (j) n to vanish for non-integer j, as well as for n = j = 0.
The tachyon profile for several different solutions with regular OPE has been calculated before. It has the simpler form of T (t) = ∞ n=0 λ n β n √ 2 n e nt where the coefficients β n def = β (0) n are only non-zero for maximal momenta. In table 3 we compare the coefficients for those solutions to the ones we have found. We have changed the normalization of their coefficients by 2 − n 2 for better comparison with our coefficients, due to the relative normalizations of the marginal operators e X 0 and √ 2 cosh(X 0 ). Our coefficients show very similar falloff to [5] as n is increased, though we do not expect exact agreement between any of the sets of coefficients because the tachyon profile is a gauge dependent quantity. We believe that each of these lists is related to the others by such gauge transformations, but constructing them is beyond the scope of this work.
As in [5], we take λ to be negative in order to study physical solutions. With this assumption, the tachyon profile (2) can be rewritten as where in practice the sum over n only runs up to some cutoff N where the coefficients can be computed. When only the j = 0 coefficients and the first term in parentheses are considered, as in the regular OPE case, we can clearly see that a change of ln |λ| will only shift the time of the D-brane decay. For the singular case, however, the tachyon profile will have a different shape depending on the strength of the marginal deformation, controlled by ln |λ|. The renormalization scheme also contains the constants C 0 and C 1 , which will appear nontrivially in β Table 1: The non-zero coefficients β (j) n of the tachyon profile for the cosh rolling tachyon with singular self-OPE. Cuhre/QAG results shown. The constant C L is part of the renormalization scheme of [2], but it can not influence the solution, so we safely set it to zero in our analysis. C L was included in these numerical results only to demonstrate that it does not contribute to the solution at all. * These two coefficients found using the Cuhre algorithm appear to be unreliable, so the corresponding Suave results in table 2 will be used for analysis instead. 0.0506 7.732 · 10 −4 9.149 · 10 −4 4 6.54812 · 10 −7 4.18 · 10 −3 9.8145 · 10 −7 1.173 · 10 −6 5 4.93424 · 10 −11 1.54 · 10 −4 8.734 · 10 −11 1.275 · 10 −10 6 3.50136 · 10 −16 2.45 · 10 −6 7.903 · 10 −13 2.26 · 10 −15 7 2.41180 · 10 −22 1.64 · 10 −8 Table 3: The coefficients for the tachyon modes of a purely exponential rolling tachyon with regular self-OPE. Calculated from the solutions of [5], [7], and [3]. The n = k coefficients for our calculations based on [1] are included for comparison.

Small λ
We begin our analysis with the case |λ| 1, where only the coefficients β (0) n need to be considered. Following the notation of [5,7], we will refer to these coefficients as β n def = β (0) n . We will focus on (5), which receives significant contributions from the first term in parentheses when t > 0 and from the second term when t < 0. Knowing that T (t) is an even function, we will assume t > 0 and not need to consider the second term. Since we are considering −1 λ < 0, each term in the sum of (5) will be suppressed by the exponential until t is large compared to − ln |λ|. For a large fixed t, terms with j > 0 will be small relative to others, so only the j = 0 coefficients need to be considered. Since this is the case, the tachyon profile does not depend on the renormalization constants C 0 and C 1 at all. This had to be the case since there is no renormalization when all of the marginal operators have momentum in the same direction. We can then unambiguously plot the tachyon profile for small |λ|. In figure 1 we see ln |T (t)| for ln |λ| = −4. Each "peak" represents the range of t for which T (t) is dominated by a specific exponential in the sum. For different values of |λ| the shape of the oscillating part of the tachyon profile remains unchanged, and the whole half-plot shifts horizontally, with the size of the plateau in the middle changing as expected.
The size of the plateau which describes the time when the D-brane exists can be estimated by the time of the first zero of the tachyon profile. This time is plotted in figure 2a, and is linear for the region where |λ| 1 is valid. The slope is −1 as it had to be from (5) when t is significantly larger than 0. We can also examine the "period" of the oscillations. The oscillations result from each exponential overtaking the one before, so we can calculate an estimate of their spacing by setting adjacent terms to be equal.
β n e n(ln |λ|+tn) = β n+1 e (n+1)(ln |λ|+tn) (6a) n coefficients considered. This is plotted for ln |λ| = −4, where the approximation is valid. a) Black is for positive values of T (t) and red is for negative values. b) The tachyon profile with all Suave coefficients from table 2 is also shown in orange and grey, and where the Cuhre-only results deviate is indicated with a blue line. We can see that the Suave results are qualitatively equivalent to the ones we use.
While there is no reason to expect this a priori, let us suppose that ∆t is a constant. In this case we have which is a recursion relation with the solution The factor ρ n can always be removed by taking β n ρ n and simultaneously λ → λρ, which does not alter the tachyon profile. In one particular solution for the rolling tachyon with regular OPE, Kiermaier, Okawa, and Soler [7] found that their solution's coefficients had the asymptotic behaviour β n ∼ e −γn 2 +O(n ln n) , and in [6] it was shown that a solution equivalent to the one in [5] has coefficients which closely fit b n ∼ e −γn 2 without significant corrections. We have just shown that this same recurring pattern can be derived from the assumption of exponentially growing oscillations with constant period. In figure 3 we see the best fit lines for our j = 0 coefficients, as well as those of several other known solutions, to the form β n ∼ e −γn 2 . This was only predicted to be a fit for one solution at large n, but we see good agreement in all cases, even with n never rising past 6 or 7 for any of the solutions considered. In figure 3a the fit is to the deterministic results of table 1 for n ≤ 5 and the Suave result for n = 6, but the Suave results with smaller n are shown as the red points for reference. The Cuhre value for β (0) 6 cannot be shown on a logarithmic plot since it has the wrong sign. It is curious that the coefficients fall so close to the e −γn 2 lines without any correction, even such as choosing ρ = 1 in (7b). While this trend was derived in our case from a constant period of oscillation, if it holds at higher n it guarantees that the edge coefficients are a convergent series. The fact that all of the solutions appear to behave similarly suggests that they are also all convergent.

Large λ
Once we loosen the |λ| 1 restriction, we must consider all of the coefficients β (j) n and search for patterns there. Due to the small number of coefficients, and particularly the small number of rows with constant j, it is not possible to get a good understanding of any patterns or asymptotics for these coefficients, but we can speculate as to possible trends. The first thing we notice from table 1 is that the sign of the coefficients appears to alternate as (−1) j 2 . This is not strictly true even for the coefficients we have calculated, however, as choosing non-zero C 0 and especially C 1 will alter many of the coefficients and can affect their sign.
With only a small number of the coefficients known, we do not know whether the large n asymptotics are fixed or can be changed by a choice of the two free parameters. We can,  table 3. Square points are for [5], crosses are for [7], and circles for [3]. however, attempt to force a few patterns and see which appear more naturally.
As a first choice we pick C 0 = C 1 = 0 and notice that the j = 1 coefficients appear to be a good fit to β (1) n ∼ e −γ 1 k 3 with k = n − 2j, which is shown in figure 4a. While this can be made to fit even better by a choice of renormalization constants, this would lead to some of β (2) n being less than zero or to that row having increasing magnitudes. On the other hand, if we attempt to pick renormalization constants which are a fit to β (1) n ∼ e γ 1 k 2 , as shown in figure 4c, we do not find as good a fit. The same is true of β (1) n ∼ e γ 1 n 2 using n instead of k in the exponent. It appears that the j = 1 coefficients have a tendency towards the cubic exponential decay, while for j = 2 we lack enough points to reach any conclusions. The red points in figure 4 again represent the Suave coefficients, and we see that β (2) n have significantly different values once the renormalization constants are changed, but a look at table 2 suggests that this is mainly due to large errors in the Suave coefficients, so it is unlikely that the deterministic plots would change significantly if more sample points were used.
While five points is not a lot of data, the β (1) n coefficients suggest that each row with constant j may eventually be a convergent series for at least some choice of renormalization constants. Showing that the full tachyon profile converges when these rows are added together, however, remains impossible until much higher order calculations can be performed. In particular, the large dependence of coefficients such as β on C 1 is troubling since it suggests that if that constant is of order 1 then the sequence β ( n−k 2 ) n with constant k could have increasing magnitudes. Thinking of (4) as the effective coefficients β (k) eff (λ) would then be defined by series which do not converge for non-zero renormalization constants. Looking at table 1, even with vanishing renormalization constants, the magnitudes of the terms β j 2j do not drop off very fast, but with alternating sign the series may still converge. Now that we have seen what the bulk coefficients look like, we can begin examining the tachyon profile for larger values of |λ|. Of course as we do this we must be aware that we are missing all coefficients β (j) n with n ≥ 7, and as we increase the strength of the marginal deformation those coefficients will begin to play a larger role, but we can still get a qualitative idea of the impact of the bulk coefficients on the tachyon profile. Since none of the optimized fits in figure 4 were significantly better than simply setting the renormalization constants to zero, we will choose that from now on. Other reasonable choices will not give results that are qualitatively different. For large negative λ we see a tachyon profile in figure 5 with fewer oscillations than we had with just the edge coefficients from figure 1. As |λ| decreases, the additional oscillation appears at ln |λ| ≈ −1.948. Once ln |λ| −2.5 the profile has stabilized and the plateau continues growing as |λ| shrinks, just as we know it should from our discussion of the tachyon profile for small |λ|. The disappearance of this oscillation for n | plotted vs. functions of k = n − 2j with several different choices for C 0 and C 1 . Black points are Cuhre values while red points are from the Suave algorithm. In a) we set C 0 = C 1 = 0 and plot j = 1 coefficients on the left and j = 2 on the right. The j = 1 coefficients with C 0 and C 1 optimized for the best linear fit appear in b), and c) attempts the same linear fit assuming a k 2 horizontal axis rather than k 3 . d) attempts a linear fit to both the j = 1 and j = 2 sets of coefficients assuming a k 2 horizontal axis, with j = 1 coefficients on the left and j = 2 on the right. large |λ| is because the bulk coefficients cannot be neglected in this region, and they change the effective coefficients in (9).
The behaviour we see for these large values of |λ| is not unprecedented; in [7] the tachyon profile had coefficients (seen in table 3) which did not decrease as quickly as other timeasymmetric solutions. That tachyon profile had fewer oscillations for all λ because some of the exponentials did not dominate for any range of time. Aside from the obvious, that the singular OPE case is a symmetric function where the D-brane exists for a limited time while with regular OPE it exists until it decays at a finite time, the qualitative difference between the tachyon profiles seems to be that the period and number of oscillations can change this way. Because the strong deformation tachyon profile we have found is similar to the profile of [7], it suggests that changing the strength of the marginal deformation in the time-symmetric case is much like changing gauge in the time-asymmetric case. If the late time behaviour is equivalent to the tachyon vacuum under a time-dependent gauge transformation, as has been hypothesized [3], then in this case the gauge transformation should depend on both time and the marginal parameter in a non-trivial way. That our solution appears qualitatively like time-asymmetric ones for both weak and strong deformation parameter suggests that such gauge transformations remain a valid explanation of the oscillations in the time-symmetric case.

A Few Technical Details
The solution and necessary renormalized operators have been defined in [1] and [2], so it is possible to numerically compute many quantities up to a reasonable accuracy. We use Maple to algebraically manipulate wedge states with insertions and produce functions to be integrated. This means writing routines for the star product, BRST operator, and correlation functions, among other things. While the number of terms grows extremely fast, it is possible to algebraically construct the solution up to 7th order in λ. The integrands produced, however, are complicated enough that I was only able to compile them for integration up to 6th order. Although the algebraic construction of the wedge states with insertions is a lengthy process, it is an exact one which should not produce any numerical errors. Evaluating the resulting integrals, however, is done by sampling the integrand at many points and estimating a result, and this process unavoidably introduces errors.
The integration is handled by off-the-shelf C++ routines. The CUBA library appears to be a good choice in most cases [12]. 1 It is a collection of four algorithms for multidimensional numerical integration, three of which use pseudo-random sampling while the fourth is a deterministic algorithm. Since we are working at sixth order in λ and the solution has ghost number one (corresponding to the number of fixed moduli), there are never more than five integrated coordinates in a given integral. While Monte Carlo algorithms do scale better as the dimension rises, in five or less dimensions it appears that the deterministic algorithm, Cuhre, is slightly more efficient. As we will see, Cuhre is also more reliable in most cases. Unfortunately, it only integrates functions of more than one variable, so in the one-dimensional case we use the QAG routine from the GNU scientific library. 2 Each of the routines in the CUBA and GNU libraries provides its own error estimate, and the CUBA library routines also provide a chi-square estimate of the probability that the error is sufficient. A single quantity to be calculated numerically, such as an individual term in the tachyon profile of table 1 or in one of the consistency checks of tables 6 and 7, generally consists of a small number of integrals. Each integral contains all of the terms in the solution which are integrated over a given number of coordinates, or equivalently all of the terms with the same number of integrated operators. Most quantities are the sum of two integrals, but a few are only a single integral, and quantities in the action can consist of more than two.

Convergence
In order to get as much data as possible, the collection of integrals we look at here will include all of the ones used in calculating the tachyon profile as well as the action and several components of the equation of motion. The action and equation of motion will be discussed as checks of consistency in section 3.3.
It is difficult to study convergence of the one-dimensional integrals due to the fact that the QAG algorithm does not report the number of samples used. We can, however, take several full sets of data for these integrals and compare the different calculations of the same integrals. We find that they all agree with each other well within the error estimates. The only troubling one-dimensional case is that of a constant integrand. This can be seen in the kinetic energy of the solution at fourth order in λ, which is the particularly simple quantity 1 0 dt 3 2 − 3 2 . The integration algorithm performs operations on the constant causing tiny roundoff errors which, when the constant value is subtracted, causes the result to differ from zero. This would not be a problem except that error analysis in numerical integration is based on variation of the integrand, and as such gives an estimated error which is extremely small. While in principle error estimates should account for the roundoffs inherent in their algorithm, in practice this does not seem to be the case. The error estimate becomes small enough that the reported result can actually be incorrect, and even unstable with respect to changing the desired accuracy. Most integrands worth using a numerical algorithm to integrate undoubtedly vary enough that this is not normally an issue. This is the only instance of such a problem that we encounter, but since the integral is trivial we do not need to rely on numerical integration for its value.
We now turn our attention to the Cuhre algorithm, so we will only be considering integrals over two or more dimensions. The first question we will ask here is whether the error bars reported by the integration algorithms are sufficient. Because we do not know the correct  Figure 6: Plots showing the reliability of error estimates. The vertical axis is the difference of two calculations relative to the more precise of the two, and the horizontal axis is the relative error reported by the numerical integration. The Cuhre algorithm results are shown in a) and Suave results are in b). Points with differences greater than twice the reported errors are red, those with differences between one and two times the error are green, and those with error estimates large enough to cover their difference from the "best" value are blue.
results for most of the integrals, we evaluate each integral with at least three different choices for the sample size, N . Making the assumption that the calculation with the largest N is "correct", we can compare the difference between each computed integral and the most accurate one to the error estimate reported by that integral. This is shown in figure 6a, where blue points have sufficient error bars and green points are within twice the error bars. The few red points are the outliers which differ from their largest N partners by more than twice their error estimates. In a moment we will compare every computed integral with other calculations of the same integral, so the red points here will also have greater than 2σ difference in that comparison. Many of these points come from the same integrals, so there are actually not very many integrals which will need to be examined in detail. Our choice to prefer the Cuhre algorithm over the adaptive Monte Carlo algorithm Suave is justified by figure 6b, where we see that the Suave algorithm is as likely as not to underestimate the error. In its defence, Suave frequently reports a 100% χ 2 estimate that the reported error is insufficient, but this is not particularly helpful in finding accurate values.
The Suave algorithm can still be useful for comparison, however. Its error bars are not helpful, but if the same quantity computed in Cuhre and Suave differs by more than a few percent it is worth closer examination. This is how the two terms marked with asterisks in table 1 were identified. The integrals responsible for the troublesome behaviour of these two quantities are plotted in figures 7a and 7b. There were some other terms which were flagged by this test, but they converged reasonably well once the sample size was increased sufficiently.
With many integrals each calculated for several different values of N , we have a sizeable collection of data to examine. Among the 848 Cuhre integrals, there is only one instance of the error estimate increasing as the sample size was increased, so we can safely say that the reported error bars decrease monotonically as N → ∞. We can then examine the quantity for every pair of calculations of the same integral. We find that 88% of pairs are within σ of each other, while 94% are within 2σ. Those which disagree by more than 2σ can be studied individually, since they correspond to only 15 different integrals. Of those, four only disagree due to a single computation each with very low N (about 500 samples) that has a 50% χ 2 chance of being incorrect. One of the remaining 11 integrals is also the one responsible for the tachyon profile coefficient β The first two plots represent the parts of the tachyon profile for which we used the Suave results. In figure 7a the problem is not that the results are inconsistent, but that the errors are so large that the results are meaningless. It looks likely that as N is increased the results will continue to converge to something quite close to the Suave result. For the other integral, figure 7b shows that as N is made extremely large we finally find something like the much more consistent Suave results. The slope, however, does not yet appear to be significantly slowing down, so we cannot be sure that it is convergent. The rest of the plots in figures 7c-7e show the other integrals that contribute to the tachyon profile and have more than a 2σ variation between points. They all show signs that once N is sufficiently large they converge quite well. Only for smaller N do the error estimates appear to be insufficient.
Moving on to integrals which contribute to consistency checks, in figure 8 we see that the majority are fine. Only 8c does not appear to converge. As with the examples shown in figure  7, this may well be linear until some critical N where it begins to converge. In addition, while only the quantities composed of small sums of these integrals are supposed to vanish, in many cases each of the integrals will vanish independently. These seven examples all become closer to vanishing as N is increased, and only 8f is not getting very close to zero. Convergence of these results does not appear to be very much of an issue, and I expect that if we continued to increase N by another factor of ten they would all continue to approach zero.

Roundoff errors
Each integrand may contain a number of terms which are divergent either on the boundary of the region, t i ∈ {a, b}, or on a diagonal, t i = t j . While the renormalization is designed specifically so that these divergences will cancel, individual terms evaluated near these regions can be very large. Because we are limited to the double precision floating point datatype, each term has a relative precision of approximately 10 −16 . Any time an individual term is more than 10 16 times the theoretical value of the integrand evaluated at the same point, the machine uncertainty coming from that term can dominate the result. We would hope that since this only happens for a small subset of the points sampled the effect will be negligible as the number of points increases, but this is not the case. If we take a random sample of N points, as is done for Monte Carlo integration, we would expect the closest point to a given boundary (or other codimension 1 subspace) to be ∼ 1 N away. Individual terms, however, often have a 1 t 2 divergence from the OPE of the marginal operator, which would lead to ∼ N 2 divergence for the closest point. This grows faster than the denominator, N , so the roundoff error in the resulting integral should increase linearly with the number of points. The deterministic case is actually worse because some points are intentionally chosen near or even right on the boundary. To combat this, whenever a sample point is close to a boundary or a diagonal, we can replace it with a nearby point giving a decent approximation to the integrand. The integrand function is effectively replaced by one where the value is held constant on small strips. While this means that a perfect integration with no uncertainty would give an incorrect result, the errors introduced this way are less problematic than the roundoff errors when we sample many points without any regulation.
When we discussed the differences between the big G and little g renormalization schemes, we saw that the lack of a regulator was an advantage of the little g scheme. Here we have introduced another regulator, so we naturally ask why this is not a problem. The regulator in the big G scheme was required by the theory in that scheme, and we wanted the limit as it approached zero. This regulator is to prevent roundoff error, which is the unavoidable result of using a floating point datatype. Since we are regulating the integration region anyway, we might ask why (aside from issues regarding finiteness at higher orders) we did not use the big G scheme. By using the little g scheme, the integrand is independent of the regulator, which simplifies the integration process. There is not a different integrand for each value of the regulator, and instead of a limit, we only use a small value of the regulator, namely 3 · 10 −4 , for which the integrands always evaluate with negligible errors. The choice of regulator for these roundoff errors is arbitrary, but we can estimate the error we have introduced by using the same regulator to replace the value of the integrand with zero near boundaries and diagonals. Fortunately, the differences are minor compared with the statistical errors which are accounted for by the algorithms' reported uncertainties.

Consistency checks
Since the programs to construct wedge states with insertions and produce and evaluate integrals corresponding to the tachyon profile are quite complicated, it is worth using them to evaluate some known quantities. We will see that the numerical integration process gives results which are consistent with expectations the majority of the time, despite the presence of counterterms and the uncertain nature of numerical integration. An obvious choice for a quantity which we know is the equation of motion, which should vanish. The equation of motion, however, has ghost number two, which means that its expectation value by itself will trivially vanish because the ghosts are not saturated. In order to test that Q B Ψ + Ψ * Ψ vanishes we test that its overlaps with various other string fields all vanish. Because the equation of motion is stronger than just requiring that the equation of motion annihilates all states and actually tells us that it should vanish exactly, as long as the string field is constructed properly these correlation functions should work out to zero whether they themselves were computed correctly or not. In order to test that a non-trivial result also gives the correct answer, we look to the action. Because this is an exactly marginal solution, we expect the energy to vanish, and because the energy is proportional to the action, the action should vanish as well. The action has ghost number three and does not need any additional test states inserted. That this also vanishes is our first strong test that non-trivial expectation values are computed successfully. A summary of these test calculations and their results using the deterministic algorithms is found in tables 6 and 7 at the end of this section, and all of them are expected to vanish. The majority of the results are consistent with zero, but a few exceptions require detailed examination. These six examples are in table 4, which restates their values using the deterministic algorithms and then includes the corresponding results with Monte Carlo calculations and with deterministic calculations using a vanishing integrand near borders and diagonals where cancelling singularities may occur.
The values in tables 6 and 7 all use a border with width = 3 · 10 −4 . When the integrand is sampled within of a boundary or a diagonal, the closest point on the edge of this strip is used instead. In the last column of table 4, when the integrand was sampled at points within these strips, zero was returned instead. The difference between these results gives an estimate of how important the regulated region is to the final result of the integral, and we can see that in most cases it is small compared to the error estimates.
but consistent with zero. Perhaps requesting a higher accuracy causes the algorithm to notice the cusp where the regulated strips at the boundary are, and that increases the error estimate. The kinetic energy at fourth order is a peculiar case of rounding errors, and it was already discussed in the context of convergence. None of the integration algorithms is correct all of the time, and we should not expect them to be, but QAG is problematic because it gives less control over the sample size and does not report the total number of points it uses. We cannot show these results with different sample sizes, as we do for other quantities of interest. Because these are one-dimensional integrals, however, we can expect that the majority of them will be computed very accurately. For the other three quantities, the integrals are multidimensional so we can evaluate them using the Cuhre algorithm with several different sample sizes and see if they become closer to zero. This is shown in table 5. The first two of these both have error estimates that decrease as the sample size is increased, and values that decrease even faster. For the first, we see agreement once the sample size is large enough, and for the second we can suspect that the result will continue to tend towards zero. The most important integrals in these two quantities were shown in figures 8a and 8d respectively, where we can see the convergence to 0 as N is increased. In the case of the final quantity, however, the result clearly does not vanish. If we change the thickness of the regulating border, however, it causes a significant fluctuation, and for extremely thin borders both the result and the uncertainty become much larger. This suggests that it is a genuine case of the border region having a significant impact. Unfortunately, the only way to resolve this issue would be to use higher precision floating point datatypes.  Table 5: Three of the consistency checks in table 4 are shown with different sample sizes. The first two are computed as the sum of two integrals with different dimensions, while the third is a single integral. We expect to see results which tend towards zero as the sample size is increased.