Bayesian Analysis of Neutron Star Mass and Radius Observations

This document describes the opensource MPI implementation of a Bayesian analysis of mass and radius data to determine the mass versus radius curve and the equation of state of dense matter. This package will principally be useful for those physicists and astrophysicists who are already familiar with C++ and are interested in modifying this code for their own use.
This code was originally supported by Chandra grant TM112003X.
This is a beta version. Use at your own risk.
Currently, bamr is dualhosted as an SVN respostory at http://www.sourceforge.net/projects/bamr and a git repository at http://www.github.com/awsteiner/bamr . The corresponding entry at the Astrophysics Source Code Library (ASCL) is http://ascl.net/1408.020 and at the ADS is http://adsabs.harvard.edu/abs/2014ascl.soft08020S . This HTML documentation is hosted at http://bamr.sourceforge.net. While you are not required to do so, please consider citing the ASCL entry following the method described at http://ascl.net/wordpress/?page_id=351 , and the relevant references in the bibliography below.
If you are considering using this code for your research, I encourage you to contact me so that I can help you with the details and so that you can let me know if and how this code is useful to you. Nevertheless, you are not required to contact me and I will be improving documentation and updating this code as time permits.
The bamr executable requires the installation of GSL (versions 1.15 and later), HDF5 (versions 1.8.4 and later),
the current development version of
O2scl (in the
branches/dev folder of the
O2scl distribution),
and MPI (tested with openmpi1.4.2). After these four packages are successfully installed, you will need to edit the makefile
and then compile bamr before execution.
The most recent release version can be obtained from either of the following:
svn checkout svn://svn.code.sf.net/p/bamr/code/trunk bamr git clone https://github.com/awsteiner/bamr.git
The basic usage is something like
./bamr model twop run default.in mcmc run1
to perform a one day run with model twop
with the input file in default.in
.
There are several variables which can be modified with the set
command, e.g.
./bamr model twop set max_time 43200 run default.in mcmc run2
which runs for 12 hours instead of the default 24 hours.
An example of an MPI invocation is
mpirun np 4 ./bamr set model twop run default.in mcmc run3 &
which runs with four processors on the current machine.
Also try
./bamr help
which outputs some additional information on the relevant functions and parameters.
The input data files are HDF5 files (typically named with a .o2
extension) which contain one o2scl::table3d object giving the probability density of one neutron star observation as a slice in that table. The command adddata
, which adds a neutron star data set, takes four arguments (and a fifth optional argument):
It is assumed that the "x" grid in the data file refers to the radius and the "y" grid refers to the mass. The data need not be normalized. By default bm
renormalizes the data so that the maximum probability is 1. If the parameter norm_max
is set to false, then the data is renormalized to have a total integral of 1.
Output is stored in HDF files with a prefix given by the argument to the mcmc
command, one set of files for each processor. Output includes files with the following suffixes (where X is the processor index):
_X_out:
Main output file containing full Markov chain(s) and most of the parameter values_X_scr:
Running output of entire simulationIf the executable is run directly (without mpirun
) then X is always zero.
Representations of the EOS and the curve are stored on grids for each Monte Carlo point. The number of points in the grid is specified by the parameter grid_size
. The energy density grid is specified by the limits e_low
and e_high
, the baryon density grid is specified by the limits nb_low
and nb_high
, and the mass grid is specified by the limits m_low
and m_high
. A default run stores pressure as a function of energy density (in columns with prefix P_
), radius as a function of mass (in columns with prefix R_
), central pressure as a function of mass (in columns with prefix PM_
), pressure as a function of baryon density (in columns with prefix Pnb_
), and energy per baryon as a function of baryon density (in columns with prefix EoA_
).
The basic functionality is provided in the bamr::bamr_class and each Monte Carlo point is an object of type bamr::entry. All of the "models" (EOS parameterizations) are children of bamr::model class.
If the initial guess has no probability, then the code will fail. This is indicated by the line "Initial weight zero."
in the _scr
file. The easiest fix is just to change the initial guess.
In order to make the output more efficient, the table representing the full Markov chain is divided up into tables with about 10,000 rows each, named markov_chain0
, markov_chain1
, and so on. The total number of tables is stored in the integer n_chains
.
Different models have different optimal MC step sizes. The step size for each parameter is chosen to be the difference betwen the high and low limiting values divided by the value step_fac
. Increasing or decreasing this value may give better results.
Several EOS models are provided. New models (i.e. new children of the bamr::model class) must perform several tasks
ed
and the pressure in pr
with the correct units set for each column (currently only 1/fm^4
is supported)."nb"
, then bamr::model::compute_eos() should return one baryon density and energy density in bamr::model::baryon_density_point()."S"
and "L"
in the table (in ).
The toplevel functions in the call stack are:
"model"
: bamr::bamr_class::set_model()"adddata"
: bamr::bamr_class::add_data()"firstpoint"
: bamr::bamr_class::set_first_point()"mcmc"
: bamr::bamr_class::mcmc()"mcmc"
command.The operation of bamr::bamr_class::compute_weight() can be summarized with:
The process
program contains some code to analyze the bamr
output files. Principally, it attempts to remove autocorrelations from the data by separating the data into blocks, and using the fluctuations in the mean of each block to determine the uncertainty in the mean. For example,
creates a new file called out.o2
with an object of type o2scl::table which represents a probability distribution for from the bamr
output file named bamr_0_out
. The option set xscale 197.33
ensures that the new table converts from to . The five columns are reps
, avgs
, errs
, plus
and minus
, which give:
April 2015: Added process.cpp and created new functions bamr::model::setup_params and bamr::model::remove_params() . Added several new models.
I would like to thank Paulo Bedaque, Ed Brown, Farrukh Fattoyev, Stefano Gandolfi, Jim Lattimer, and Will Newton for their collaboration on these projects.
Some of the references which contain links should direct you to the work referred to using its DOI identifer.
Bedaque15sv: P. Bedaque and A.W. Steiner, Phys. Rev. Lett. 114 (2015).
Lattimer14co: J.M. Lattimer and A.W. Steiner, Eur. Phys. J. A 50 (2014) 40.
Lattimer14ns: J.M. Lattimer and A.W. Steiner, Astrophys. J. 784 (2014) 123.
Steiner10te: A.W. Steiner, J.M. Lattimer, E.F. Brown, Astrophys. J. 722 (2010) 33.
Steiner12cn: A.W. Steiner and S. Gandolfi, Phys. Rev. Lett. 108 (2012) 081102.
Steiner13tn: A.W. Steiner, J.M. Lattimer, E.F. Brown, Astrophys. J. Lett. 765 (2013) 5.
Steiner15un: A.W. Steiner, S. Gandolfi, F.J. Fattoyev, and W.G. Newton, Phys. Rev. C 91 (2015) 015804.
Documentation generated with Doxygen. Bamr documentation is under the GNU Free Documentation License.