CobayaSampler / cobaya

Code for Bayesian Analysis
http://cobaya.readthedocs.io/en/latest/
Other
128 stars 127 forks source link

how do i set fluid_equation_of_state correctly in the .yaml file? #220

Closed pasq-cat closed 2 years ago

pasq-cat commented 2 years ago

i have a local installation of class which has been modified to include a linear redshift model (LRS) and a squared redshift model (SRS). Everythings seems to work when i type ./class testlrs.ini but when i try to use jupyter+cobaya together with my definition of dark energy model i get an error. As an example, i used cobaya-cosmo-generator to generate a cosmo-g.yaml file for a clp model and then i modified it by including the line fluid_equation_of_state: LRS

The result looks like this


theory:
  classy:
    path: /home/pasquale/Documenti/class_public/
    extra_args:
      N_ncdm: 1
      N_ur: 2.0328
likelihood:
  sn.pantheon: null
params:
  logA:
    prior:
      min: 1.61
      max: 3.91
    ref:
      dist: norm
      loc: 3.05
      scale: 0.001
    proposal: 0.001
    latex: \log(10^{10} A_\mathrm{s})
    drop: true
  A_s:
    value: 'lambda logA: 1e-10*np.exp(logA)'
    latex: A_\mathrm{s}
  n_s:
    prior:
      min: 0.8
      max: 1.2
    ref:
      dist: norm
      loc: 0.965
      scale: 0.004
    proposal: 0.002
    latex: n_\mathrm{s}
  theta_s_1e2:
    prior:
      min: 0.5
      max: 10
    ref:
      dist: norm
      loc: 1.0416
      scale: 0.0004
    proposal: 0.0002
    latex: 100\theta_\mathrm{s}
    drop: true
  100*theta_s:
    value: 'lambda theta_s_1e2: theta_s_1e2'
    derived: false
  H0:
    latex: H_0
  omega_b:
    prior:
      min: 0.005
      max: 0.1
    ref:
      dist: norm
      loc: 0.0224
      scale: 0.0001
    proposal: 0.0001
    latex: \Omega_\mathrm{b} h^2
  omega_cdm:
    prior:
      min: 0.001
      max: 0.99
    ref:
      dist: norm
      loc: 0.12
      scale: 0.001
    proposal: 0.0005
    latex: \Omega_\mathrm{c} h^2
  Omega_m:
    latex: \Omega_\mathrm{m}
  omegamh2:
    derived: 'lambda Omega_m, H0: Omega_m*(H0/100)**2'
    latex: \Omega_\mathrm{m} h^2
  m_ncdm:
    value: 0.06
    renames: mnu
  w0_fld:
    prior:
      min: -3
      max: 1
    ref:
      dist: norm
      loc: -0.99
      scale: 0.02
    proposal: 0.02
    latex: w_{0,\mathrm{DE}}
  wa_fld:
    prior:
      min: -3
      max: 2
    ref:
      dist: norm
      loc: 0
      scale: 0.05
    proposal: 0.05
    latex: w_{a,\mathrm{DE}}
  Omega_Lambda: 0
  fluid_equation_of_state: LRS
  YHe:
    latex: Y_\mathrm{P}
  tau_reio:
    prior:
      min: 0.01
      max: 0.8
    ref:
      dist: norm
      loc: 0.055
      scale: 0.006
    proposal: 0.003
    latex: \tau_\mathrm{reio}
  z_reio:
    latex: z_\mathrm{re}
sampler:
  mcmc:
    covmat: auto

Initially i thought there was something wrong with my definition LRS but then i tested fluid_equation_of_state: CLP and it also fail:

[tools] ERROR Failed to load external function: 'NameError("name 'CLP' is not defined")'


NameError Traceback (most recent call last) ~/anaconda3/lib/python3.9/site-packages/cobaya/tools.py in get_external_function(string_or_function, name) 374 with PythonPath(os.curdir, when="import_module" in string_or_function): --> 375 function = eval(string_or_function, scope) 376 except Exception as e:

~/anaconda3/lib/python3.9/site-packages/cobaya/tools.py in

NameError: name 'CLP' is not defined

During handling of the above exception, another exception occurred:

LoggedError Traceback (most recent call last) /tmp/ipykernel_4970/93674339.py in 5 from cobaya.run import run 6 ----> 7 updated_info, sampler = run(info_from_yaml)

~/anaconda3/lib/python3.9/site-packages/cobaya/run.py in run(info_or_yaml_or_file, packages_path, output, debug, stop_at_error, resume, force, no_mpi, test, override) 132 out, info=updated_info["sampler"][sampler_name]) 133 # 4. Initialize the posterior and the sampler --> 134 with Model(updated_info["params"], updated_info["likelihood"], 135 updated_info.get("prior"), updated_info.get("theory"), 136 packages_path=info.get("packages_path"),

~/anaconda3/lib/python3.9/site-packages/cobaya/model.py in init(self, info_params, info_likelihood, info_prior, info_theory, packages_path, timing, allow_renames, stop_at_error, post, skip_unused_theories, dropped_theory_params) 141 if v not in (None, {}): 142 self._updated_info[k] = deepcopy_where_possible(v) # type: ignore --> 143 self.parameterization = Parameterization(self._updated_info["params"], 144 allow_renames=allow_renames, 145 ignore_unused_sampled=post)

~/anaconda3/lib/python3.9/site-packages/cobaya/parameterization.py in init(self, info_params, allow_renames, ignore_unused_sampled) 144 else: 145 self._input[p] = np.nan --> 146 self._input_funcs[p] = get_external_function(info["value"]) 147 self._input_args[p] = getfullargspec(self._input_funcs[p]).args 148 if is_sampled_param(info):

~/anaconda3/lib/python3.9/site-packages/cobaya/tools.py in get_external_function(string_or_function, name) 375 function = eval(string_or_function, scope) 376 except Exception as e: --> 377 raise LoggedError( 378 log, "Failed to load external function%s: '%r'", 379 " '%s'" % name if name else "", e)

LoggedError: Failed to load external function: 'NameError("name 'CLP' is not defined")'

so i think there is something wrong with the .yaml file itself. Any advice?

JesusTorrado commented 2 years ago

Hi @Rockdeldiablo,

Try passing it within extra_args instead. The parameter values in the params block must be floats or functions. Let us know whether it works!

pasq-cat commented 2 years ago

@JesusTorrado Hi , thanks for the advice, it actually worked for the CLP case, but i got errors in the LRS case that i didn't get with class:

[classy] Importing local CLASS from '/home/pasquale/Documenti/class_public/'. [classy] Initialized! [mcmc] Getting initial point... (this may take a few seconds) [classy] ERROR Serious error setting parameters or computing results. The parameters passed were {'A_s': 2.1096843541402204e-09, 'n_s': 0.9596304447630893, '100theta_s': 1.0419851109969605, 'omega_b': 0.02232934984579731, 'omega_cdm': 0.11932177393089169, 'm_ncdm': 0.06, 'w0_fld': -0.9697904479820014, 'wa_fld': -0.0439339879994118, 'Omega_Lambda': 0.0, 'tau_reio': 0.06032438609964339} and {'N_ncdm': 1, 'N_ur': 2.0328, 'fluid_equation_of_state': 'LRS', 'output': ''}. To see the original CLASS' error traceback, make 'debug: True'. [classy] ERROR Serious error setting parameters or computing results. The parameters passed were {'A_s': 2.115641683170898e-09, 'n_s': 0.9696622923862939, '100theta_s': 1.0412783700099566, 'omega_b': 0.022314289354474673, 'omega_cdm': 0.11901162440385196, 'm_ncdm': 0.06, 'w0_fld': -0.9643041922676895, 'wa_fld': 0.01067827263957695, 'Omega_Lambda': 0.0, 'tau_reio': 0.06588872216506923} and {'N_ncdm': 1, 'N_ur': 2.0328, 'fluid_equation_of_state': 'LRS', 'output': ''}. To see the original CLASS' error traceback, make 'debug: True'. [classy] ERROR Serious error setting parameters or computing results. The parameters passed were {'A_s': 2.1128119234359723e-09, 'n_s': 0.9665499244775799, '100*theta_s': 1.042025680113818, 'omega_b': 0.022244669827657305, 'omega_cdm': 0.11918836210183668, 'm_ncdm': 0.06, 'w0_fld': -0.9669670565125368, 'wa_fld': -0.09398703849893869, 'Omega_Lambda': 0.0, 'tau_reio': 0.05323826284155562} and {'N_ncdm': 1, 'N_ur': 2.0328, 'fluid_equation_of_state': 'LRS', 'output': ''}. To see the original CLASS' error traceback, make 'debug: True'.

etc etc

---------------------------------------------------------------------------
LoggedError                               Traceback (most recent call last)
/tmp/ipykernel_9097/93674339.py in <module>
      5 from cobaya.run import run
      6 
----> 7 updated_info, sampler = run(info_from_yaml)

~/anaconda3/lib/python3.9/site-packages/cobaya/run.py in run(info_or_yaml_or_file, packages_path, output, debug, stop_at_error, resume, force, no_mpi, test, override)
    141                 updated_info = recursive_update(updated_info, model.info())
    142                 out.check_and_dump_info(None, updated_info, check_compatible=False)
--> 143                 sampler = sampler_class(updated_info["sampler"][sampler_name],
    144                                         model, out, name=sampler_name,
    145                                         packages_path=info.get("packages_path"))

~/anaconda3/lib/python3.9/site-packages/cobaya/sampler.py in __init__(self, info_sampler, model, output, packages_path, name)
    266                 pass
    267         self._set_rng()
--> 268         self.initialize()
    269         model.set_cache_size(self._get_requested_cache_size())
    270         # Add to the updated info some values which are

~/anaconda3/lib/python3.9/site-packages/cobaya/samplers/mcmc/mcmc.py in initialize(self)
    154             self.log.info("Getting initial point... (this may take a few seconds)")
    155             initial_point, results = \
--> 156                 self.model.get_valid_point(max_tries=self.max_tries.value,
    157                                            random_state=self._rng)
    158             # If resuming but no existing chain, assume failed run and ignore blocking

~/anaconda3/lib/python3.9/site-packages/cobaya/model.py in get_valid_point(self, max_tries, ignore_fixed_ref, random_state)
    474                                             "likelihood. Set 'ref' to a different point "
    475                                             "or a pdf.")
--> 476             raise LoggedError(self.log, "Could not find random point giving finite "
    477                                         "posterior after %g tries", max_tries)
    478         return initial_point, results

LoggedError: Could not find random point giving finite posterior after 320 tries

i will post down here what i did for the LRS (linear redshift w= w0 +wa*z) model(i tried to follow the clp example since i don't know much about the c language):

1) in background.h i wrote: enum equation_of_state {CLP,EDE,LRS};

2) in input.c i wrote:

 if (pba->Omega0_fld != 0.) {
    /** 8.a.1) PPF approximation */
    /* Read */
    class_call(parser_read_string(pfc,"use_ppf",&string1,&flag1,errmsg),
               errmsg,
               errmsg);
    if (flag1 == _TRUE_){
      if(string_begins_with(string1,'y') || string_begins_with(string1,'Y')){
        pba->use_ppf = _TRUE_;
        class_read_double("c_gamma_over_c_fld",pba->c_gamma_over_c_fld);
      }
      else {
        pba->use_ppf = _FALSE_;
      }
    }

    /** 8.a.2) Equation of state */
    /* Read */
    class_call(parser_read_string(pfc,"fluid_equation_of_state",&string1,&flag1,errmsg),
               errmsg,
               errmsg);
    /* Complete set of parameters */
    if (flag1 == _TRUE_) {
      if ((strstr(string1,"CLP") != NULL) || (strstr(string1,"clp") != NULL )) {
        pba->fluid_equation_of_state = CLP;
      }
      else if ((strstr(string1,"LRS") != NULL) || (strstr(string1,"lrs") != NULL)) {
        pba->fluid_equation_of_state = LRS;
      }
      else if ((strstr(string1,"EDE") != NULL) || (strstr(string1,"ede") != NULL)) {
        pba->fluid_equation_of_state = EDE;
      }
      else {
        class_stop(errmsg,"incomprehensible input '%s' for the field 'fluid_equation_of_state'",string1);
      }
    }

    if (pba->fluid_equation_of_state == CLP) {
      /** 8.a.2.2) Equation of state of the fluid in 'CLP' case */
      /* Read */
      class_read_double("w0_fld",pba->w0_fld);
      class_read_double("wa_fld",pba->wa_fld);
      class_read_double("cs2_fld",pba->cs2_fld);
    }
            if (pba->fluid_equation_of_state == LRS) {
      /** 8.a.2.2) Equation of state of the fluid in 'LRS' case */
      /* Read */
      class_read_double("w0_fld",pba->w0_fld);
      class_read_double("wa_fld",pba->wa_fld);
      class_read_double("cs2_fld",pba->cs2_fld);
    }

    if (pba->fluid_equation_of_state == EDE) {
      /** 8.a.2.3) Equation of state of the fluid in 'EDE' case */
      /* Read */
      class_read_double("w0_fld",pba->w0_fld);
      class_read_double("Omega_EDE",pba->Omega_EDE);
      class_read_double("cs2_fld",pba->cs2_fld);
    }
  }

3) in background.c i added the 3 functions

    *w_fld = pba->w0_fld - pba->wa_fld/a  + pba->wa_fld  ;

    *dw_over_da_fld = pba->wa_fld/(a*a) ;

    *integral_fld = 3.*( pba->wa_fld -(1. +  pba->w0_fld  + pba->wa_fld ) *log(a) -   pba->wa_fld/a) ;
    break;

which i derived earlier. The last one should be the integral between a and a0 of 3*(1+w(a))/a da with a0 = 1, as the instructions in the comment said. I don't know if it's a class issue or cobaya, but with ./class testrls.ini i get no errors or warnings