CosmoMC: Adding additional parameter

[Vimal Simha]

Hi,

I am interested in adding an additional parameter to cosmomc - N_nu the number of neutrino species. This is already defined in CAMB.

1) I have modified the CMBParams to Params subroutine in params_CMB.f90 and the CMBtoCAMB sub-routine in CMB_Cls_simple.f90 for this?
I also changed the cmbtypes.f90 to define CMB%nu.

2)I use Param(13) as N_nu and avoid having to change settings.f90

3) What exactly does slow/fast do? How do I change this?

4) I am worried that CAMB does not accept values of N_nu higher than ~3.4 or so. Why is this and is there any way of getting around this?

5) Should I be concerned about the following message:

WARNING: computing tensors with num_cls=3 (BB=0)

I have got this despite setting compute tensors = F in params.ini



ALSO,

On the 1-D probability distribution plots, is it possible for me to plot \Delta \chi^2 (delta chisquared) on the vertical axis?

Thank You.

[Antony Lewis]

1) I have modified the CMBParams to Params subroutine in params_CMB.f90 and the CMBtoCAMB sub-routine in CMB_Cls_simple.f90 for this? I also changed the cmbtypes.f90 to define CMB%nu.

In general you need to change ParamsToCMBParams and CMBParamsToParams in params_CMB.f90, as well as the others you mention. Also RecomputeTransfers in CMB_Cls_simple.f90 if you are adding a slow parameter.

(As mentioned in the ReadMe faq it is often easier to re-interpret an exisiting parameter than adding a new one.)

2)I use Param(13) as N_nu and avoid having to change settings.f90
3) What exactly does slow/fast do? How do I change this?

By default params(13) is a fast parameter - not what you want for N_nu. See the CosmoMC paper for a discussion of fast and slow parameters. The numbers of fast and slow parameters is set in settings.f90.

4) I am worried that CAMB does not accept values of N_nu higher than ~3.4 or so. Why is this and is there any way of getting around this?

It hasn't really been tested for odd values - you can probably just edit CAMB_ValidateParams in camb.f90.

5) Should I be concerned about the following message:
WARNING: computing tensors with num_cls=3 (BB=0)

No, I accidentally missed a '.and. compute_tensors' in the May 2006 release. It will be fixed in the next update.

On the 1-D probability distribution plots, is it possible for me to plot \Delta \chi^2 (delta chisquared) on the vertical axis?

The plots are marginalized probability or mean likelihoods. Whether or not this is what you mean depends how you are definining your generalized \Delta \chi^2. There are no options to label the left axis with numbers (normalization is usually arbitrary).

[Anze Slosar]

I found it useful to define param 7 to be w if <0 and N_\nu if >0. This way you don't have to break the number of slow params plus you keep the old functionality of param 7 as w...

[Vimal Simha]

Thank You very much.

However, I am not clear about what I ought to do to the 'Recompute Transfers' subroutine in CMB_Cls.....

[Antony Lewis]

If N_nu has changed at the proposed parameter position, you need to re-compute the transfer functions to get the new C_l (unlike changing say n_s - a fast parameter - where the transfer functions are the same). So you need to check in RecomputeTransfers whether the old and new N_nu are the same or not.

[Vimal Simha]

But, am I correct in thinking that I do not need to do this if I set use_slow_fast = F in params.ini ?

[Antony Lewis]

No, use_fast_slow = F just stops it using fast subspace proposals. If you are starting new chains without a covmat you also need to edit RecomputeTransfers (unless you e.g. re-interpret w, in which case no need as that is already accounted for).

If you are using a covmat the chance of it making a proposal where no slow parameters change is negligible, so it probably doesn't matter.

[Vimal Simha]

Thank You very much.

By reinterpreting one of the unused parameters, I presume you mean that:

I change the CMBtoCAMB routine, so that

eg.

P%omegav=CMB%omv

is changed to

P%XYZ=CMB%omv

and then fix P%omegav.

[Antony Lewis]

Yes (see the w_is_w option for an example of re-interpreting w)

[Vimal Simha]

So, from the paper and from what you have said, am I correct in thinking that if I recompute the transfer functions for param 13, it would be equivalent to using it as a slow parameter.

[Antony Lewis]

If you also use use_fast_slow = F, then yes.

[Florian Beutler]

Hi everyone,
I found this discussion about adding the effective number of neutrinos as an additional parameter in cosmoMC. I tried to follow the steps discussed above, but it did not work out. Here is what I did...

1. I tried to reinterpret the parameter w. Hence I set

logical :: w_is_w = .false.

in CMB_Cls_simple.f90. I also made changes the subroutine CMBToCAMB to

...
if (w_is_w) then
w_lam = CMB%w
else
!P%InitialConditionVector(initial_iso_baryon) = CMB%w
P%Num_Nu_Massless = CMB%nnu - 3.046
w_lam = -1
end if
if (CMB%nnu < 3.04) call MpiStop('CMBToCAMB: nnu < 3.04, would give negative masless neutrinos')
!P%Num_Nu_Massless = CMB%nnu - 3.046 !assume three massive
...

where I replace "initial_iso_baryon" with the number of neutrinos... is Num_Nu_Massless the right parameter??

2. In the params.ini file I replaced

param[w] = -1 -1 -1 0 0

by

param[w] = 0 0 10 8.5 8.5

where the "start width" and "st. dev. estimate" are just guesses... as far as I understand the parameter Num_Nu_Massless is zero in LCDM?

Did I forget anything? Running this setup gives me the error message


$ ./cosmomc params.ini
Random seeds: 8194, 14740 rand_inst: 0
Computing tensors: F
Doing CMB lensing: T
lmax = 2100
Number of C_ls = 3
Varying 8 parameters ( 0 fast)
starting Monte-Carlo
Initialising BBN Helium data...
Warning: Rombint failed to converge;
integral, error, tol: 10241.672989512392 1.40870146103644345E-007 1.00000000000000009E-007
reading WMAP7 data
Initializing WMAP likelihood, version v4.1
Warning: Rombint failed to converge;
integral, error, tol: 10891.877962639144 -1.14582630930115670E-005 1.00000000000000009E-007
Warning: Rombint failed to converge;
integral, error, tol: 3239.4778743209367 -4.10066775891060331E-005 1.00000000000000009E-007
Warning: Rombint failed to converge;
integral, error, tol: 13090.461766344491 -9.09578591334891939E-007 1.00000000000000009E-007
Warning: Rombint failed to converge;
integral, error, tol: 4150.4179487180118 -2.22699991803132491E-005 1.00000000000000009E-007
Warning: Rombint failed to converge;
integral, error, tol: 11493.806585279612 -9.49991998844268437E-006 1.00000000000000009E-007

thanks a lot for your help
cheers
fl

[Florian Beutler]

Hi everyone,
ok I made some progress but it still doesn't seem to work, so I am still thankful for help and advice. The changes to the code I made are the following:

1. I set w_is_w to false
2. I changes in CMB_Cls_simple the following lines
...
if (w_is_w) then
w_lam = CMB%w
else
!P%InitialConditionVector(initial_iso_baryon) = CMB%w
P%Num_Nu_Massless = CMB%w - 3.046
w_lam = -1
end if
!if (CMB%nnu < 3.04) call MpiStop('CMBToCAMB: nnu < 3.04, would give negative masless neutrinos')
!P%Num_Nu_Massless = CMB%nnu - 3.046 !assume three massive
...

3. I changed in params_CMB the following lines

in SetForH()
...
!CMB%w = Params(7)
CMB%w = -1
CMB%nnu = Params(7)
...

and in CMBParamsToParams I changed

...
!Params(7) = CMB%w
Params(7) = CMB%nnu
...

the run works for a while, but the parameter distributions do not look correct... and at some point I get the error message


$ ./cosmomc params.ini
Random seeds: 19268, 13452 rand_inst: 0
Computing tensors: F
Doing CMB lensing: T
lmax = 2100
Number of C_ls = 3
Varying 8 parameters ( 0 fast)
starting Monte-Carlo
Initialising BBN Helium data...
reading WMAP7 data
Initializing WMAP likelihood, version v4.1
Input of bbn_splint out of interpolation range.
7.0000000000000000 8.3492631912231445 -3.0000000000000000

this error message appears to be dependent on the value of the upper prior on N_{nu} which I set to 10. If I set it to 15, the error message appears directly, while at 10 it takes a while...? Before the error message appears it gives an output, which I compare to the Nrel output on the NASA webpage (Komatsu et al.). For some parameter I find huge differences (see attachment)?? Since the constraint in \Omega_{c}h^{2} seems to suffer a huge amount because of the new parameter I guess I forgot something in the code?
Can somebody explain to me why \Omega_{c}h^{2} looks so bad in this 8 parameter model? It must have a huge degeneracy with N_{nu}? (red is Komatsu et al and blue is my result).

/cosmo_files.php?file=Florian_Beutler/c3.pdf

thanks a lot for help
cheers
fl

[Ali Vanderveld]

Florian -- In the end were you ever able to get this to work properly?
Thanks!

[Jan Hamann]

This error message is generated in the BBN module, which sets the primordial Helium fraction [tex]Y_p[/tex] as a function of [tex]\omega_{\rm b}[/tex] and [tex]N_{\rm eff}[/tex] by interpolating from a pre-computed grid supplied in data/helium.dat.
If your input values of [tex]\omega_{\rm b}[/tex] or [tex]N_{\rm eff}[/tex] are outside the parameter ranges covered by data/helium.dat, the code will stop with the error message you got.

The supported parameter ranges are [0.015,0.04] for [tex]\omega_{\rm b}[/tex] and [0.046,10.046] for [tex]N_{\rm eff}[/tex] .

You can fix this problem in two ways:
Either switch off the bbn_consistency option in params.ini, or, if you want to keep the BBN consistency relation, set an upper limit of 10.046 on [tex]N_{\rm eff}[/tex] in params.ini.

And by the way, [tex]\omega_{\rm c}[/tex] (or [tex]\omega_{\rm m}[/tex]) is indeed massively degenerate with [tex]N_{\rm eff}[/tex] for WMAP data. The reason is that the quantity actually constrained by WMAP data is the redshift of equality (see 0803.0547 for a more detailed discussion).