Implementing modified gravity models

[sthagstotz]

Hello everyone, I would like to implement modified gravity models following the mu / gamma parametrisation also used e.g. in MGCAMB (http://arxiv.org/abs/1106.4543), where mu modifies the Poisson equation and gamma is the ratio between the psi and phi potentials (in absence of anisotropic stress). The implementation of gamma seems straightforward, but I'm a bit lost in the perturbation module considering mu.. Did anyone try something similar before?

Thanks a lot! Steffen

[lesgourg]

Hi Steffen,

I agree that it is not as simple as one may think for the following reason:

1) if you run CLASS in newtonian gauge, the code uses only 2 out of the 4 perturbed scalar Einstein equations, namely the one involving the shear to get psi (since phi is known at each step, it is integrated by the differential equation solver), and the one involving the velocity divergence (i.e. (0i)) to get phi_prime. So the equation (00) is not used, and unfortunately the Poisson equation relates to that one. In principle, using the Bianchi identities, I think that you should be able to find analytically the only modification of the 2 used equations that is consistent with a given modification of the Poisson equation. But this might not be so simple. I never tried myself. Maybe it is written somewhere in the literature?

2) if you run CLASS in the synchronous gauge, it might be easier. In that gauge CLASS uses the equation (00) involving the density to get h_prime (since eta is known at each step, it is integrated by the differential solver), and the equation (0i) involving theta to get eta_prime. The other two equations, the ones involving pressure and shear, also appear, but they are not so important, they are used only within some approximation schemes. The equation (00) has a term [k^2 eta], and of course the Poisson equation is related to that. However, be careful, there are non-trivial gauge transformations relating the Poisson equation to the einstein equations in synchronous gauge. [k^2 eta] is not just equal to [k^2 phi] or [k^2 psi]. But you can start from the Poisson equation, carefully perform the gauge transformation to the synchronous gauge using Ma & Bertschinger formulas, and then write your modifications in the synchronous Einstein equations

I don't know which line from 1) or 2) is easiest, maybe 2). Also the paper you are referring to probably used 2). It would be interesting to know if you (or somebody else) find the solution.

Best Julien

[sthagstotz]

Hi Julien, thanks a lot! I'm working on option 2) now and I'm trying to get a toy model to test everything. I implemented the evolution equations in synchronous gauge for alpha, eta_prime and h_prime from http://arxiv.org/abs/1106.4543 in perturb_einstein, but to simplify matters I only switch to the modified equations for z<30 (I just assume the shear terms are zero then, sigma and sigma_prime terms complicate things slightly). Then I get an error:

Error in perturb_init =>perturb_init(L:318) :error in perturb_solve(ppr, pba, pth, ppt, index_md, index_ic, index_k, pppw[thread]); =>perturb_solve(L:1908) :error in generic_evolver(perturb_derivs, interval_limit[index_interval], interval_limit[index_interval+1], ppw->pv->y, ppw->pv->used_in_sources, ppw->pv->pt_size, &ppaw, ppr->tol_perturb_integration, ppr->smallest_allowed_variation, perturb_timescale, ppr->perturb_integration_stepsize, ppt->tau_sampling, tau_actual_size, perturb_sources, perhaps_print_variables, ppt->error_message); =>evolver_ndf15(L:461) :condition (absh <= hmin) is true; Step size too small: step:2.22045e-16, minimum:2.22045e-16, in interval: [246.106:14159.9]

Is there anything obviously wrong with my approach?

Cheers, Steffen

[lesgourg]

Hi Steffen, Yes there is an obvious problem: any high-performance differential integration solver (with adaptative step) needs continuous equations. You are putting by hand a discontinuity at z=30. Evn if it is a small discontinuity, the integrator will find that around z=30, the time step must be reduced to zero, to avoid infinite derivatives... this causes the issue. Possible ways out:

1) code things in such a way that the code will integrated perturbations till z=30, then stop, then restart the integration (such that a time step will never cross z=30, it will just end or start at z=30. Doable, but a bit difficult. I would not recommend it to people not having a strong experience with the code.

2) use the modified equations at anytime. Why not? It would probably be the simplest thing?

3) use a continuous function transforming continuously the standard equations into the modified equations. For instance, suppose that you define the function s(z)=[tanh((z-30)/3)+1]/2 going from 0 to 1 around z=30, with a width of delta_z=3 that you can change. Then you can write your Einstein equations as [normal terms] + s(z) * [additional terms] This will work from the point of view of the code. From the point of view of physics, I am not sure: the Einstein equations need to satisfy some consistency conditions, maybe this is breaking them.

[sthagstotz]

Hi Julien, Again, thanks a lot for your help and comments. I am implementing the full equations valid at all times now and came across some questions (and a minor bug): 1) In perturb_sources there are two implementations to calculate the temperature source terms. The standard (efficient) way relies on alpha_prime, but I can't use the Einstein equation related to the shear anymore to get it easily. I switched to the other one, but there seems to be something wrong with the equations there - the TT spectrum is far off (I also checked that this is the case for the standard 2.4.2 version of class). Are they derived or written down somewhere so I can check the equations? As far as I can tell the other uses of alpha_prime are optional.

2) To implement shear_prime terms (the modified equations depend on them), I guess I would have to add entries to the perturbation vector y[...shear_prime] and derivatives dy[...shear_prime]. The calculation of y[...shear_prime] is already included (as dy[...shear]), but would it work to implement dy[...shear_prime] just once via the boltzmann hierarchy? Or do I have to worry about different expressions during tight coupling?

3) In perturb_einstein, the the order of calculations is: h_prime-> theta(h_prime) -> eta_prime(theta) -> alpha(eta_prime, h_prime) -> shear_g(alpha) the modified equations don't decouple in such a way anymore, all quantities depend on a combination of theta, rho, shear and shear_prime. I'm not sure how to handle the updates of theta and the shear consistently now... but calculating h', eta', alpha once in the beginning, then updating theta / shear with those values and finally updating h', eta' and alpha again seems to go wrong.

[ThomasTram]

Hi, I can give some hints regarding your question (1). The simple temperature source term which is commented out is in fact correct. The error you see in the TT spectrum is just numerical, and that is the reason that we are using a more complicated expression.

The simple expression can be cross checked against Hu&White 1998 for instance, Sec. 2.8 of http://arxiv.org/pdf/1305.3261.pdf, and Eq. 3.5 of http://arxiv.org/pdf/1312.2697.pdf. Our choice of source function splitting is given in Appendix A of 1312.2697.

I am not sure we ever tried to converge the simple expressions numerically, but it might be possible. Try to decrease the following two parameters: perturb_sampling_stepsize = 0.10 tol_perturb_integration = 1e-5 and possibly increase ppr->hyper_sampling_flat = 8.0

Cheers, Thomas

[sthagstotz]

Hi Thomas, thanks a lot! That does indeed improve the TT spectrum. However the increase in calculation time is quite high, so I guess eventually deriving an analytic expression for alpha_prime might be the better way to go.

Cheers, Steffen

[hasi03]

Hi, I'm looking for mu / gamma parametrisation in CLASS like in MGCAMB. Is there a patch already ? Has some one done it? Thank you

[MBilicki]

Hi, I'm looking for mu / gamma parametrisation in CLASS like in MGCAMB. Is there a patch already ? Has some one done it? Thank you

Hi about hi_class? http://miguelzuma.github.io/hi_class_public/

[hasi03]

Hi, Thank you, but I couldn't find a mu,gamma parameterize (for a model independent way) as in MGCAMB. It seemed it only has Horndeski implementation similar to EFTCAMB. Do you know it has already mu,gamma or mu,Sigma parameterization? Sorry I couldn't find it there. Thanks again