Skip to content

Commit

Permalink
Updated compliation unit and README to reflect addition of weak_rates…
Browse files Browse the repository at this point in the history 
… module.
  • Loading branch information
Chris Sullivan committed Sep 2, 2015
1 parent 33a3034 commit 9eba865
Show file tree
Hide file tree
Showing 17 changed files with 218 additions and 1,569 deletions.
154 changes: 119 additions & 35 deletions README
View file
Expand Up @@ -6,7 +6,7 @@
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!

DISCLAIMER:
DISCLAIMER:
-----------

Please note that the routines provided here come with absolutely no
Expand All @@ -22,7 +22,7 @@ please e-mail us at evanoc@tapir.caltech.edu or cott@tapir.caltech.edu
COPYRIGHT:
----------

While NuLib is open source, its copyright is held by Evan O’Connor and
While NuLib is open source, its copyright is held by Evan O'Connor and
Christian Ott. In the absence of suitable open scientific software
licenses, we release this version of NuLib to the community under the
Creative Commons attribution-noncommercial-share alike license:
Expand All @@ -37,16 +37,15 @@ and if so, only under the same license.
Introduction:
-------------

The goal of NULIB is to provide a basic standard set of neutrino
The goal of NuLib is to provide a basic standard set of neutrino
matter interaction routines that can be readily incorporated in
radiation-hydrodynamics codes for physics benchmarking.

NuLib v1.0 includes the basic neutrino emissivities and absorption
opacities (including pair processes) as well as neutrino-nucleon and
neutrino-nucleus scattering processes. Inelastic scattering (with
energy redistribution) is neglected, but will be included in future
revisions. Other inelastic processes will also be including in future
versions.
opacities (including pair processes) as well as neutrino-nucleon,
eutrino-nucleus elastic scattering processes and neutrino-electron
inelastic scattering. Other inelastic processes will also be
including in future versions.

If anyone would like to contribute to the development of NuLib please
let me know, I'm hosting this on GitHub with this in mind. I am
Expand All @@ -59,20 +58,21 @@ more than happy to make NuLib as accessible as it needs to be to allow
others to benefit from it.

NuLib, in its current form, is used by me to make tables of neutrino
emissivities and opacities (scattering and absorption). It is not yet
optimized for on-the-fly calculations of quantities. In fact, the
routines are coded in such a way as to be as clear and accurate as
possble, with little or no regard for computational speed. For
example, the weak magnitism correction for scattering processes is
small but many terms long, I do the full calculation. In the future I
hope NuLib will have routines capable of on the fly calculations, I
expect this is necessary if one whats fully differential cross
sections (in energy and angle, as a function of energy).
emissivities, opacities (scattering and absorption), and scattering
kernels. It is not yet optimized for on-the-fly calculations of
quantities. In fact, the routines are coded in such a way as to be as
clear and accurate as possble, with little or no regard for
computational speed. For example, the weak magnitism correction for
scattering processes is small but many terms long, I do the full
calculation. In the future I hope NuLib will have routines capable of
on the fly calculations, I expect this is necessary if one whats fully
differential cross sections (in energy and angle, as a function of
energy).

NuLib v1.0 Neutrino Interactions.
---------------------------------

Emissivities:
Emissivities:
-------------

1. electron-positron annhilation to \nu - \bar{\nu}. This is currently
Expand Down Expand Up @@ -134,23 +134,80 @@ blocking. Weak magnitism, phase space and recoil corrections
3. \nu_e absorption on heavy nuclei: This follows the simple treatment
of Bruenn 85 (among others), placing cuts on the cross section based
on the average A and average Z of the nucleus. Much better treatment
is desired and someday will be implemented.
is desired and someday will be implemented.

Neutrino-Electron Inelastic Scattering Kernels:
-----------------------------------------------

For a temperature and electron chemical potential, NuLib calculates
the first two terms in a Legendre expansion of the scattering kernels
for neutrinos on electrons. We essentially follow Bruenn 1985 and
references there in.


Electron-Capture Rates on Nuclei:
---------------------------------

Chris Sullivan and Evan O'Connor et al. have implemented a new module
for microphysical electron-capture rates on nuclei in NuLib. It utilizes the
formalism discussed in:
------------------------------------------------------------------------------------
| Sullivan, C., O'Connor, E., Zegers, R. G. T., Grubb, T., & Austin, S. M. (2015). |
| The Sensitivity of Core-Collapse Supernovae to Nuclear Electron Capture. |
| http://arxiv.org/abs/1508.07348 |
| Contact: Chris Sullivan <sullivan@nscl.msu.edu> |
------------------------------------------------------------------------------------
The primary calculation is electron-neutrino type emissivities from
electron captures on medium-heavy nuclei. These calculations rely on
a library ofelectron-capture rate tables that have been compiled and
are availableas a part of the weak_rates module (set in the parameters
file). In addition, number densities (abundances) and nuclear masses
are needed for a large set of nuclei. These are calculated via Matthias
Hempel's NSE mass distributions discussed below. To include emissivities
and opcaities for this interaction, set the corresponding flag in
requested_interactions.inc, as well as WEAK_RATES=1 in make.inc. If
these routines are utilized in a work, cite the above paper as well the
following publications from which the weak-rate tables derive:
--------------------------------------------------------------------------------
| Fuller, G. M., Fowler, W. A., & Newman, M. J. (1982). |
| Stellar weak interaction rates for intermediate-mass nuclei. |
| II - A = 21 to A = 60. The Astrophysical Journal, 252, 715. |
| http://doi.org/10.1086/159597 |
--------------------------------------------------------------------------------
| Oda, T., Hino, M., Muto, K., Takahara, M., & Sato, K. (1994). |
| Rate Tables for the Weak Processes of sd-Shell Nuclei in Stellar Matter. |
| Atomic Data and Nuclear Data Tables, 56(2), 231-403. |
| http://doi.org/10.1006/adnd.1994.1007 |
--------------------------------------------------------------------------------
| Langanke, K., & Mart\'{i}nez-Pinedo, G. (2000). |
| Shell-model calculations of stellar weak interaction rates: |
| II. Weak rates for nuclei in the mass range in supernovae environments. |
| Nuclear Physics A, 673(1-4), 481-508. |
| http://doi.org/10.1016/S0375-9474(00)00131-7 |
--------------------------------------------------------------------------------
| Langanke, K., & Mart\'{i}nez-Pinedo, G. (2003). |
| Electron capture rates on nuclei and implications for stellar core collapse. |
| Physical Review Letters 90, 241102. |
| http://prl.aps.org/abstract/PRL/v90/i24/e241102 |
--------------------------------------------------------------------------------


Sample Executables:
-------------------

make_table_example: by default this makes a horribly under resolved
10x10x10x24 (rho,temp,ye,energy) NuLib table in h5 format, the
1000*24*3*3 calculations takes abut 1 minute to generate. The table
boundaries, and number of data points are changable in the
make_table_example.F90 file. To get enough accuracy in the
interpolation I expect at least 10 points per decade in rho, 20 in
temperature and 1 for every 0.01 or 0.02 in ye. This makes a table
~1GB in size with 24 energy bins. The energy spacing is changable in
nulib.F90, right now it is a 4MeV bin, then a logarithmic spacing
starting at 1MeV going to ~300MeV, this may not be the best choice, if
you have a better suggestion, let me know, or code up a routine to
generate good energy spacing and send a pull request (?).
10x10x10x24 (rho,temp,ye,energy) + (10x10x24*24)
(temp,eta,energy_in,energy_out) NuLib table in h5 format. The
calculations takes abut 1 minute to generate. The table boundaries,
and number of data points are changable in the make_table_example.F90
file. To get enough accuracy in the interpolation I expect at least
10 points per decade in rho, 20 in temperature and 1 for every 0.01 or
0.02 in ye. This makes a table ~1GB in size with 24 energy bins. The
energy spacing is changable in nulib.F90, right now it is a 4MeV bin,
then a logarithmic spacing starting at 1MeV going to ~300MeV, this may
not be the best choice, if you have a better suggestion, let me know,
or code up a routine to generate good energy spacing and send a pull
request (?).

You must specify an equation of state, NuLib is set up to read in the
EOS tables on stellarcollapse.org, the filename is set in
Expand All @@ -169,11 +226,11 @@ me know, I want to make this as useful as possible.

! many people use different number of species, this is the possible
! summing scheme NuLib can currently do
!
!
! mytable_neutrino_scheme = 1 (three output species)
! species #1: electron neutrino #2 electron antineutrino
! #3: muon+tau neutrino+antineutrino
!
!
! neutrino_scheme = 2 (four output species)
! species #1: electron neutrino #2 electron antineutrino
! #3: muon+tau neutrino #4 mu and tau antineutrino
Expand All @@ -187,7 +244,10 @@ single_point_return_all appies Kirchoff's law to the emissivities and
absorption cross sections. This adds an contribution to the
emissivity from the absorption cross section (and vice-versa). This
is explained in BRT06 and explicitly showed in the
single_point_return_all routine.
single_point_return_all routine. There is a similar routine for the
inelastic scattering kernels, single_Ipoint_return_all. This routine
only calculates half of the terms, we use symmetry laws to calculate
the other half.

point_example: this program shows examples of how to call the NuLib
routines for a single point. Again, the energy spacing is changable in
Expand All @@ -205,7 +265,7 @@ Unlike single_point_return_all, the individual calls to emissivity
(e.g. return_emissivity_spectra_given_neutrino_scheme) or the cross
section routines
(e.g. return_absorption_opacity_spectra_given_neutrino_scheme) do not
apply Kirchoff's law.
apply Kirchoff's law.

nulibtable_driver: This routine is a driver routine for reading in a
NuLib table and using a trilinear interpolation (log rho, log temp,
Expand Down Expand Up @@ -239,7 +299,7 @@ configuring with your version of:

./configure --enable-fortran FC=ifort --prefix=/Users/evanoc/opt/hdf5-current-ifort12

and then
and then

make
make install
Expand All @@ -250,3 +310,27 @@ After this, a simple make should create three executables in the main
directory, a brief explanation of these is in the section
`executables' above.

Extras
------

There is support in NuLib for using Matthias Hempel's NSE mass
distributions available from
http://phys-merger.physik.unibas.ch/~hempel/eos.html. To enable these
use must download his code and tables from his website, place them in
the directory src/extra_code_and_tables/ and enable the preprocessor
flag NUCLEI_HEMPEL. We use the SFHo table as an example in the code,
you can change this by editting nuclei_distribution_helpers.F90
directly. Please see this file for more details.

A few small changes must be made to xxxx_xxxx_composition_module.f
as provided by M. Hempel if the weak_rates module (electron-capture
rates).Primarily a public (non-private) copy of the loaded nuclear
masses must be exposed. e.g. for the SFHo EOS, one must add the
following to the source file:
---------------------------------------------------
| double precision, dimension(kmax) :: sfho_mass |
|---------------------- & ------------------------|
| sfho_mass = mass |
---------------------------------------------------
We also found that an updated path is needed for the
composition binary included with the module.
34 changes: 30 additions & 4 deletions parameters
View file
Expand Up @@ -6,11 +6,37 @@ lmp_rates = "/projects/ceclub/gr1dnulib/GitHub/NuLib/src/extra_code_and_t
lmsh_rates = "/projects/ceclub/gr1dnulib/GitHub/NuLib/src/extra_code_and_tables/lmsh-extended.dat"
oda_rates = "/projects/ceclub/gr1dnulib/GitHub/NuLib/src/extra_code_and_tables/oda_hires.dat"
ffn_rates = "/projects/ceclub/gr1dnulib/GitHub/NuLib/src/extra_code_and_tables/ffn_hires.dat"
nuclear_masses = "/projects/ceclub/gr1dnulib/GitHub/NuLib/src/extra_code_and_tables/mass.mas12"

## Priority ##
ilmp = 1 # 1 indicates highest priority, if other tables have rates for a nuclei contained in this table, it will be ignored.
ilmsh = 2 # 0 indicates the table will not be used
## Rate table priority hierarchy ##
ilmp = 1
ilmsh = 2
ioda = 3
iffn = 0
iapprox = 4

# 0 indicates the table will not be used
# 1 indicates highest priority, if other tables have rates for a nuclei contained in this table, they will be ignored.
# CITATION NOTICE: By using the ilmp/ilmsh/iffn/ioda/iapprox flags you agree to cite the relevant publications in your work.
------------------------------------------------------------------------------------------
| iffn | Fuller, G. M., Fowler, W. A., & Newman, M. J. (1982). |
| | Stellar weak interaction rates for intermediate-mass nuclei. |
| | II - A = 21 to A = 60. The Astrophysical Journal, 252, 715. |
| | http://doi.org/10.1086/159597 |
------------------------------------------------------------------------------------------
| ioda | Oda, T., Hino, M., Muto, K., Takahara, M., & Sato, K. (1994). |
| | Rate Tables for the Weak Processes of sd-Shell Nuclei in Stellar Matter. |
| | Atomic Data and Nuclear Data Tables, 56(2), 231-403. |
| | http://doi.org/10.1006/adnd.1994.1007 |
------------------------------------------------------------------------------------------
| ilmp | Langanke, K., & Mart\'{i}nez-Pinedo, G. (2000). |
| | Shell-model calculations of stellar weak interaction rates: |
| | II. Weak rates for nuclei in the mass range in supernovae environments. |
| | Nuclear Physics A, 673(1-4), 481-508. |
| | http://doi.org/10.1016/S0375-9474(00)00131-7 |
------------------------------------------------------------------------------------------
| ilmsh | Langanke, K., & Mart\'{i}nez-Pinedo, G. (2003). |
| iapprox | Electron capture rates on nuclei and implications for stellar core collapse. |
| | Physical Review Letters 90, 241102. |
| | http://prl.aps.org/abstract/PRL/v90/i24/e241102 |
------------------------------------------------------------------------------------------

14 changes: 11 additions & 3 deletions src/Makefile
View file
@@ -1,7 +1,7 @@
include ../make.inc

SOURCES=nulib.F90 \
weak_rates.F90 \
SOURCES=nulib.F90\
weak_rates.F90\
fermi.F90 \
helpers.F90 \
gauss_laguerre_helpers.F90 \
Expand All @@ -20,11 +20,19 @@ NT_SOURCES=nulibtable.F90 \
nulibtable_reader.F90 \
linterp_many_mod.F90


#the weak_rates module requires M. Hempel's
#EOS dependent NSE distributions
ifeq ($(WEAK_RATES),1)
NUCLEI_HEMPEL=1
DEFS += -DWEAK_RATES
endif

ifeq ($(NUCLEI_HEMPEL),1)
F_SOURCES = extra_code_and_tables/sfho_frdm_comp/sfho_frdm_composition_module.f
F_OBJECTS=$(F_SOURCES:.f=.o )
MODINC += -I./extra_code_and_tables/sfho_frdm_comp
DEFS += -DNUCLEI_HEMPEL
DEFS += -DNUCLEI_HEMPEL
endif

#dual openmp+mpi jobs may work, but are not currently supported
Expand Down
4 changes: 3 additions & 1 deletion src/emissivities.F90
View file
Expand Up @@ -402,6 +402,7 @@ subroutine return_emissivity_spectra_given_neutrino_scheme(emissivity_spectra,eo
emissivity_spectra(ns,ng) = emissivity !ergs/cm^3/s/MeV/srad
enddo

#if WEAK_RATES
!eos composition modules (for NSE) require T>0.1MeV
if (eos_variables(tempindex).gt.1.0d-1) then

Expand All @@ -417,7 +418,8 @@ subroutine return_emissivity_spectra_given_neutrino_scheme(emissivity_spectra,eo
end if

end if

#endif

enddo

end subroutine return_emissivity_spectra_given_neutrino_scheme
8 changes: 6 additions & 2 deletions src/input_parser.F90
View file
Expand Up @@ -2,12 +2,15 @@
subroutine input_parser(fn)

use nulib
#if WEAK_RATES
use weak_rates

#endif

implicit none
character*(*) fn

call get_string_parameter(fn,'eos_table_name',eos_filename)
#if WEAK_RATES
call get_string_parameter(fn,'lmp_rates',files_to_load(1))
call get_string_parameter(fn,'lmsh_rates',files_to_load(2))
call get_string_parameter(fn,'oda_rates',files_to_load(3))
Expand All @@ -17,7 +20,8 @@ subroutine input_parser(fn)
call get_integer_parameter(fn,'ioda',file_priority(3))
call get_integer_parameter(fn,'iffn',file_priority(4))
call get_integer_parameter(fn,'iapprox',file_priority(5))

#endif

contains

subroutine get_string_parameter(fn,parname,par)
Expand Down

0 comments on commit 9eba865

Please sign in to comment.