Update Getting started authored by Doriann Blain's avatar Doriann Blain
In this section we will describe how to quickly run an *Exo-REM* simulation. We will use the `inputs/example.nml` coming with all distributed versions of *Exo-REM* as a starting point. In this section we will describe how to quickly run an *Exo-REM* simulation. We will use the `inputs/example.nml` coming with all [distributed versions](https://gitlab.obspm.fr/dblain/exorem/-/tree/master/dist) of *Exo-REM* as a starting point.
# Prerequisites # Prerequisites
If you downloaded the archive in the [_dist_](https://gitlab.obspm.fr/dblain/exorem/-/tree/master/dist) directory, you should have everything you need, except a stellar spectrum (see point 4 below). Otherwise, check the following:
1. Verify if the collision-induced absorption (CIA) files are inside the _data/cia_ directory. If you want to use other CIA files, respect the same format and use "CIAname.cia.txt" as file name (e.g. "H2-He.cia.txt"). 1. Verify if the collision-induced absorption (CIA) files are inside the _data/cia_ directory. If you want to use other CIA files, respect the same format and use "CIAname.cia.txt" as file name (e.g. "H2-He.cia.txt").
2. If you intend to use clouds, verify if the cloud optical constants (OCST) files are inside the _data/cloud_optical_constants_ directory. If you want to use other OCST files, respect the same format and use "cloudName.ocst.txt" as file name (e.g. "H2O.ocst.txt"). 2. If you intend to use clouds, verify if the cloud optical constants (OCST) files are inside the _data/cloud_optical_constants_ directory. If you want to use other OCST files, respect the same format and use "cloudName.ocst.txt" as file name (e.g. "H2O.ocst.txt").
3. Put k-coefficients of desired absorbers inside the _data/k_coefficients_tables_ directory. You can use directly *CKC* outputs converted into HDF5 using *txt2h5*, or e.g. *[petitRADTRANS](https://petitradtrans.readthedocs.io/en/latest/content/available_opacities.html)* k-tables. If your absorber name is "absorberName", the lines file needs to be named "absorberName.ktable.exorem.h5" (e.g. H2O.ktable.exorem.h5). You can download compressed *Exo-REM* k-tables [here](https://gitlab.obspm.fr/dblain/exorem/-/tree/master/data/k_coefficients). 3. Put k-coefficients of desired absorbers inside the _data/k_coefficients_tables_ directory. You can use directly *CKC* outputs converted into HDF5 using *txt2h5*, or e.g. *[petitRADTRANS](https://petitradtrans.readthedocs.io/en/latest/content/available_opacities.html)* k-tables. If your absorber name is "absorberName", the lines file needs to be named "absorberName.ktable.exorem.h5" (e.g. H2O.ktable.exorem.h5). You can download compressed *Exo-REM* k-tables [here](https://gitlab.obspm.fr/dblain/exorem/-/tree/master/data/k_coefficients).
...@@ -9,8 +11,99 @@ In this section we will describe how to quickly run an *Exo-REM* simulation. We ...@@ -9,8 +11,99 @@ In this section we will describe how to quickly run an *Exo-REM* simulation. We
6. Put a temperature profile as a priori inside the _inputs/atmospheres/temperature_profiles_ directory. The *Exo-REM* data format must be respected. You should have received an example of such a file with your *Exo-REM* distribution. 6. Put a temperature profile as a priori inside the _inputs/atmospheres/temperature_profiles_ directory. The *Exo-REM* data format must be respected. You should have received an example of such a file with your *Exo-REM* distribution.
# Setup # Setup
7. Copy and edit the file _inputs/example.nml_ to suit your needs. **Be careful:** In this example, we will simulate the atmosphere of CoRoT-4 b, a well studied planet. A good source of planetary information can be found [here](https://exoplanetarchive.ipac.caltech.edu/). We will use the parameters from Moutou et al. 2008.
- `n_cia` must match the number of values of `cia_names`;
- `n_species` must match the number of values of `species_names` **and** `species_mode`; Copy and edit the file _inputs/example.nml_, rename it _corot-4b.nml_. An extended description of the input parameters is available [here]().
- `n_clouds` must match the number of values of all the cloud parameters **except** `cloud_mode` and `cloud_fraction`;
- `species_names` and `cia_names` values must be the names used for the _data/cia_ and _data/k_coefficients_tables_ files. 1. Open _inputs/_corot-4b.nml_ with any notepad editor.
\ No newline at end of file 2. Edit the suffix of your future output files:
```text
output_files_suffix = 'corot-4b' ! suffix of the output files
```
3. Edit the planetary ("target") parameters. We will assume an internal temperature (the effective temperature without stellar irradiation) of 130 K (a good way to estimate the internal temperature of a planet is to use a model such as this one from [Rogers and Seager 2010](https://iopscience.iop.org/article/10.1088/0004-637X/712/2/974)):
```text
target_mass = 1.37e27 ! (kg) mass of the target
target_equatorial_gravity = 0 ! (m.s-2) surface gravity of the target
target_equatorial_radius = 85100e3 ! (m) equatorial radius of the target
target_polar_radius = 0 ! (m) polar radius of the target
target_flattening = 0 ! flattening of the target
latitude = 0 ! (deg) latitude of observation on the target
target_internal_temperature = 130 ! (K) internal temperature of the target
```
4. It is always better to use a stellar spectrum rather than a blackbody spectrum. Here we will download a [BT-Settl](http://svo2.cab.inta-csic.es/theory/newov2/index.php) spectrum model. Take T_eff = 6200 K, Log(g) = 4.5, and a metallicity of 0. Put the file into the _data/stellar_spectra_ directory and rename it e.g. "spectrum_BTSettl_6200K_logg4.5_met0.dat" (mind the .dat extension). Replace the header of the file by the following :
```text
# wavelength spectral_radiosity ! effective_temperature = 6200 K
# angstrom erg.s-1.cm-2.a-1
```
5. Edit the stellar ("light source") parameters:
```text
add_light_source = True ! if True, add the light source
use_irradiation = False ! if True, use irradiation instead of range to calculate the light source spectrum
use_light_source_spectrum = True ! if True, use a spectrum for the light source instead of a black body
light_source_radius = 814e6 ! (m) radius of the light source
light_source_range = 13.5e9 ! (m) distance between the target and the light source
light_source_effective_temperature = 6190 ! (K) light source effective temperature
light_source_irradiation = 0 ! (W.m-2) light source irradiation
light_source_spectrum_file = 'spectrum_BTSettl_6200K_logg4.5_met0.dat' ! spectrum of the light sourcetarget
```
6. As a starting point, we will use a metallicity of 1 times the solar metallicity, no cloud, and a fixed eddy diffusion coefficient. Because of the effective temperature of the planet, TiO, VO and FeH are unlikely to have significant absorptions, so we will remove them in order to speed-up the calculations:
```text
metallicity = 1.0 ! (solar metallicity) atmospheric metallicity
n_levels = 81 ! number of atmospheric levels
n_species = 10 ! number of absorbing species ; an error can happen if this number doesn't match the size of species_names or species_mode
n_cia = 3 ! number of collision induced absorptions ; an error can happen if this number doesn't match the size of cia_names
n_clouds = 0 ! number of clouds in atmosphere ; an error can happen if this number doesn't match the size of the cloud arrays
eddy_mode = 'constant' ! eddy diffusion coefficient mode ('constant'|'Ackerman'|'AckermanConvective'|'infinity')
eddy_diffusion_coefficient = 1e8 ! (cm2.s-1) eddy diffusion coefficient
```
7. Edit the species parameters accordingly:
```text
species_names = 'CH4', 'CO', 'CO2', 'H2O', 'H2S', 'HCN', 'K', 'Na', 'NH3', 'PH3' ! absorbing species in atmosphere
species_at_equilibrium = False, False, False, False, False, False, False, False, False, False ! if True, the species is at thermochemical equilibrium, else it is out of equilibrium
```
8. The CoRot-4 is quite a hot star, so we will increase the wavenumber range in order to correctly take into account its spectrum.:
```text
wavenumber_max = 50130 ! (cm-1) last wavenumber
```
9. We will use the reference temperature profile provided with the Exo-REM distribution, so we are likely to be far from the "right" temperature profile. So, increase the number of iterations and reduce the retrieval tolerance to 0:
```text
n_iterations = 50 ! number of iterations
```
```text
retrieval_tolerance = 0 ! tolerance for the flux convergence (0 to use the iterations limits)
```
10. Now we should be ready to go !
# Running
1. Open a terminal.
2. `cd` yourself into the Exo-REM _bin_ directory.
3. Launch the calculations by executing:
```bash
./exorem.exe ../inputs/corot-4b.nml
```
During the run, keep a look in the terminal at how much the temperature vary (`dT`), at the ratio of the calculated internal radiosity over the theoretical radiosity (`J_int / (sigma * T_int_th**4)`), and at the variation of Chi^2 (`Chi2_var`). You should aim for a ratio close to 1, meaning that the current solution satisfy what you required. A low `dT` mean and `Chi2_var` mean that you are well converged (a negative `Chi2_var` mean that your current solution is better that the solution at the last iteration). You should have well converged by iteration 20, so in this case we overestimated the number of iterations we needed. This just mean that we waited longer than necessary to have the results. However, if the convergence was bad, we would like to re-run the code with a higher `n_iterations`, but this may not be sufficient. More on that [here]().
By default, at the end of the calculations, the results are stored in the _outputs/exorem_ directory.
# Plotting
Let's take a look at our results. Go back to the _exorem_ directory and execute:
```bash
python exorem_plot.py corot-4b
```
The figures should be in the _outputs/figures_ directory.
The temperature profile figure (temperature_profile_corot-4b) should look like this:
![temperature_profile_corot-4b](uploads/c74ea47ee170fd1d56327929950222ac/temperature_profile_corot-4b.png)
The profile crosses the condensation curves of several species. It seems that clouds are likely to form on this planets.
The transmission spectrum should look like this:
![transmission_spectrum_corot-4b](uploads/52ad367b11a16e6032d66c12a7152142/transmission_spectrum_corot-4b.png)
This is nice, but the resolution is quite low.
# More precision !
This time, our goal will be to have more precise results. We will use our calculated temperature profile as input. To keep it simple, we will consider only KCl and Na2S clouds.