Extensions¶
Accuracy¶
The accuracy module is designed to alter the accuracy parameters in a dynmaic way. For cosmological runs it is strongly recommended to set dTheta and nReplicas in the following way:
from accuracy import classic_theta_switch,classic_replicas_switch
dTheta = classic_theta_switch()
nReplicas = classic_replicas_switch()
- class accuracy.classic_replicas_switch(theta_switch=0.52)¶
This is the classic replicas switch. Use this by setting the parameter nReplicas to an instance of this class.
- Parameters:
theta_switch (number) – theta threshold to switch between 1 and 2 replicas
- Returns:
number of replicas
- Return type:
integer
- class accuracy.classic_theta_switch(before20=0.4, between20and2=0.55, after2=0.7)¶
This is the classic theta switch. Use this by setting the parameter dTheta to an instance of this class.
- Parameters:
before20 (number) – theta to use before redshift 20
between20and2 (number) – theta to use between redshift 20 and 2
after2 (number) – theta to use after redshift 2
- Returns:
theta
- Return type:
number
PKDGRAV Interface¶
Instead of simulation mode, you can enter analysis mode by importing the PKDGRAV
module. Once imported you are then responsible for calling the correct methods
to achieve your analysis tasks. For example, the following script can be used
to measure the power spectrum of a periodic simulation output.
# simulation mode is skipped when PKDGRAV imported
import PKDGRAV as msr
from sys import argv,exit
# Name and grid size are required; box size is optional
if len(argv)<=2:
print(f"Usage: {argv[0]} <filename> <grid> [box-length]")
exit(1)
name=argv[1]
nGridPk = int(argv[2])
L = 1 if len(argv)<4 else int(argv[3])
# Load the file and setup the tree, then measure the power
time = msr.load(name,nGridPk=nGridPk)
msr.domain_decompose()
msr.build_tree()
(K,PK,N,*rest) = msr.measure_pk(nGridPk,L=L)
with open(f'{name}.pk','w') as f:
for k,pk,n in zip(K,PK,N):
if n>0: f.write(f'{k} {pk} {n}\n')
The complete list of methods is as follows.
- PKDGRAV.add_analysis(callable, memory=None)¶
Add an analysis function to the simulation.
- Parameters:
callable – the analysis function
memory – the ephemeral memory required by the analysis function
If memory is not specified then an attempt is made to use the ephemeral method of the callable to get the memory requirements if it exists. Otherwise it is assumed that the callable does not require ephemeral memory.
- PKDGRAV.add_linear_signal(seed, L, a, fixed=False, phase=0.0, grid_index=0)¶
Add a linear signal to the grid
- Parameters:
seed (integer) – random seed
L (number) – length unit of the box
a (number) – expansion factor
fixed (Boolean) – use fixed amplitude for the power spectrum
phase (number) – phase of the initial conditions (in units of \(\pi\) radians, normally 0 or 1)
grid_index (integer) – which grid number to use (default 0)
- PKDGRAV.assign_mass(order=3, grid_index=0, delta=0.0, fold=1)¶
Assign mass to the grid
- Parameters:
order (integer) – mass assignment order, 0=NGP, 1=CIC, 2=TSC, 3=PCS
grid_index (integer) – which grid number to use
delta (number) – grid shift (normally 0.0 or 0.5)
fold (number) – number of times to fold
- PKDGRAV.bispectrum_calculate(grid_index0, grid_index1, grid_index2)¶
Calculate the bispectrum
- Parameters:
grid_index0 (integer) – which grid number to use
grid_index1 (integer) – which grid number to use
grid_index2 (integer) – which grid number to use
- PKDGRAV.bispectrum_normalize(k_min, k_max, target_grid_index=0)¶
Normalize the bispectrum
- Parameters:
target_grid_index (integer) – which grid number to use
k_min (number) – minimum k value
k_max (number) – maximum k value
- PKDGRAV.bispectrum_select(k_min, k_max, target_grid_index=0, source_grid_index=1)¶
Select the bispectrum
- Parameters:
target_grid_index (integer) – which grid number to use
source_grid_index (integer) – which grid number to use
k_min (number) – minimum k value
k_max (number) – maximum k value
- PKDGRAV.build_tree(ewald=None)¶
Builds a tree in each domain
- Parameters:
ewald (Boolean) – construct a moment for Ewald if true
- PKDGRAV.density_contrast(grid_index=0, k=True)¶
Compute the density contrast
- Parameters:
grid_index (integer) – which grid number to use
- PKDGRAV.domain_decompose(rung=0)¶
Particles are ordered spatially across all nodes and cores using the Orthagonal Recursive Bisection (ORB) method.
- Parameters:
rung (integer) – balance work based on this rung
- PKDGRAV.fof(tau, minmembers=10, statistics=True)¶
Friends of friends (fof) group finding
- Parameters:
tau (number) – linking length
minmembers (integer) – minimum group size (in particles)
- PKDGRAV.generate_ic(cosmology: Cosmology, *, grid: int, seed: int, z: float, L: float, order: int = None, fixed_amplitude: bool = False, phase_pi: float = 0, **kwargs) float¶
Generate initial conditions for a cosmological simulation.
- Parameters:
cosmology (Cosmology) – cosmology
grid (integer) – grid size of the initial conditions
seed (integer) – random seed
z (number) – starting redshift
L (number) – length unit of the box
order (integer) – IC order, 1=Zeldovich, 2=2LPT, 3=3LPT
fixed_amplitude (Boolean) – use fixed amplitude for the power spectrum
phase_pi (number) – phase of the initial conditions (in units of \(\pi\) radians, normally 0 or 1)
- Returns:
time
- PKDGRAV.get_array(field=None, time=1.0, marked=False)¶
Retrieves an array with requested field. OBSOLETE: use get_arrays instead.
- Parameters:
field (number) – the field to retrieve. Values are:
FIELD_POSITION
FIELD_ACCELERATION
FIELD_VELOCITY
FIELD_POTENTIAL
FIELD_GROUP
FIELD_MASS
FIELD_SOFTENING
FIELD_DENSITY
FIELD_BALL
FIELD_PARTICLE_ID
FIELD_GLOBAL_GID
- Parameters:
time (number) – simulation time
marked (Boolean) – retrieve only marked particles
- PKDGRAV.get_arrays(*fields, time=1.0, marked=False)¶
Retrieves arrays for one or more particle fields.
- Parameters:
fields – field name(s) to retrieve, e.g.
'position','velocity'. Valid names are the keys of_FIELD_MAP.time (number) – simulation time
marked (Boolean) – retrieve only marked particles
- Returns:
a single array when one field is given, otherwise a tuple of arrays
- Return type:
numpy.ndarray or tuple of numpy.ndarray
- PKDGRAV.grid_bin_k(bins, grid_index)¶
Bin the k-space grid
- Parameters:
bins (integer) – number of bins
grid_index (integer) – which grid number to use
- PKDGRAV.grid_create(grid)¶
Create a grid for the mass assignment
- Parameters:
grid (integer) – grid size for the mass assignment
- PKDGRAV.grid_delete()¶
Delete the grid for the mass assignment
- PKDGRAV.grid_ephemeral(grid, count=1)¶
Return the Ephemeral memory required for the grid
- Parameters:
grid (integer) – grid size for the mass assignment
count (integer) – number of grids
- PKDGRAV.grid_interlace(target_grid_index=0, source_grid_index=0)¶
Interlace the grid
- Parameters:
target_grid_index (integer) – which grid number to use
source_grid_index (integer) – which grid number to use
- PKDGRAV.grid_write(filename, k=False, grid_index=0)¶
Write the grid to a file
- Parameters:
filename (str) – the name of the file
k (Boolean) – write the k-space grid
grid_index (integer) – grid index
parallel (Boolean) – number of parallel tasks
- PKDGRAV.group_stats()¶
Calculate group statistics
- PKDGRAV.load(filename, **kwargs)¶
Read particles from an input file.
- Parameters:
filename (str) – the name of the file or list of files
- Returns:
time
- Return type:
number
- PKDGRAV.load_checkpoint(filename, species=None, classes=None, step=None, steps=None, time=None, delta=None, E=None, U=None, Utime=None, **kwargs)¶
Read particles from a checkpoint file.
- Parameters:
filename (str) – the name of the file
- Returns:
time
- Return type:
number
- PKDGRAV.mark_box(center, apothem, set_if_true=1, clear_if_false=1)¶
Mark particles inside a given box
- Parameters:
center (number[3]) – center coordinate
apothem (number[3]) – distance from center to edge
set_if_true (integer) – mark the particle if it is inside
clear_if_flase (integer) – unmark the particle if it is outside
- Returns:
number of particles marked
- Return type:
integer
- PKDGRAV.mark_cylinder(point1, point2, radius, set_if_true=1, clear_if_false=1)¶
Mark particles inside a cylinder
- Parameters:
point1 (number[3]) – center of the first end of the cylinder
point2 (number[3]) – center of the second end of the cylinder
radius (number) – radius of the cylinder
set_if_true (integer) – mark the particle if it is inside
clear_if_flase (integer) – unmark the particle if it is outside
- Returns:
number of particles marked
- Return type:
integer
- PKDGRAV.mark_species(*, species, species_mask, set_if_true=1, clear_if_false=1)¶
Mark particles of a given species
- Parameters:
species (integer) – species number
species_mask (integer) – species mask (one or more species)
set_if_true (integer) – mark the particle if it is inside
clear_if_flase (integer) – unmark the particle if it is outside
- Returns:
number of particles marked
- Return type:
integer
- PKDGRAV.mark_sphere(center, radius, set_if_true=1, clear_if_false=1)¶
Mark particles inside a given sphere
- Parameters:
center (number[3]) – center coordinate
radius (number) – distance from center to edge
set_if_true (integer) – mark the particle if it is inside
clear_if_flase (integer) – unmark the particle if it is outside
- Returns:
number of particles marked
- Return type:
integer
- PKDGRAV.measure_pk(grid, bins=0, a=1.0, interlace=True, order=3, L=1.0)¶
Measure the Power spectrum P(k) for the box.
- Parameters:
grid (integer) – grid size for the mass assignment
bins (integer) – number of bins for P(k), defaults to half the grid size
a (number) – expansion factor
interlace (Boolean) – use interlacing to reduce grid aliasing
order (integer) – mass assignment order, 0=NGP, 1=CIC, 2=TSC, 3=PCS
L (number) – length unit of the box to convert k and P(k) to physical units
- Returns:
k, P(k), N(k), Pall(k)
- Return type:
tuple of numpy arrays
- PKDGRAV.reorder()¶
Returns particles to their original order (usually done for output). If particles are unordered then this function has no effect.
- PKDGRAV.restore(filename, species=None, classes=None, step=None, steps=None, time=None, delta=None, E=None, U=None, Utime=None, **kwargs)¶
Restore a simulation from a file.
- Parameters:
filename (str) – the name of the file
species – species counts (optional)
classes – particle classes (optional)
step (integer) – current step (optional)
steps (integer) – total steps (optional)
time (number) – simulation time (optional)
delta (number) – time step delta (optional)
E (number) – total energy (optional)
U (number) – potential energy (optional)
Utime (number) – time of potential energy calculation (optional)
The species and classes parameters are lists of tuples. Each tuple contains the species number, mass, softening, and imat values.
The kwargs are additional parameters to set.
- PKDGRAV.rs_extract(filename)¶
Extract particles to a file
- Parameters:
filename_template (str) – the name of the file
- PKDGRAV.rs_halo_load_ids(filename, append=False)¶
Load halo IDs from a file
- Parameters:
filename (str) – the name of the file
bAppend (Boolean) – append to the existing IDs
- PKDGRAV.rs_load_ids(filename, append=False)¶
Load particle IDs from a file
- Parameters:
filename (str) – the name of the file
bAppend (Boolean) – append to the existing IDs
- PKDGRAV.rs_reorder_ids()¶
Reorder the IDs
- PKDGRAV.rs_save_ids(filename)¶
Save particle IDs to a file
- Parameters:
filename (str) – the name of the file
- PKDGRAV.save(filename, time=1.0)¶
Save particles to a file.
- Parameters:
filename (str) – the name of the file
time (number) – simulation time
- PKDGRAV.set_arrays(*, time=1.0, marked=False, **fields)¶
Writes arrays into one or more particle fields.
- Parameters:
time (number) – simulation time
marked (Boolean) – write only marked particles (currently unsupported)
fields – keyword arguments mapping field names to arrays, e.g.
position=r, mass=m. Valid names are the keys of_FIELD_MAP.
- PKDGRAV.simulate(**kwargs)¶
Directly enter simulation mode. Normally simulation mode is disabled if PKDGRAV has been imported. Calling this function proceeds normally.
- PKDGRAV.smooth(type, n=32, time=1.0, delta=0.0, symmetric=False, resmooth=False)¶
Smooths the density field with a given kernel
- Parameters:
type (integer) – smoothing kernel type
n (integer) – smoothing kernel size
time (number) – simulation time
delta (number) – time step delta
symmetric (Boolean) – use symmetric smoothing
Values for the smoothing kernel type are:
SMOOTH_TYPE_DENSITY
SMOOTH_TYPE_F1
SMOOTH_TYPE_M3
SMOOTH_TYPE_GRADIENT_M3
SMOOTH_TYPE_HOP_LINK
SMOOTH_TYPE_BALL
SMOOTH_TYPE_PRINTNN
SMOOTH_TYPE_HYDRO_DENSITY
SMOOTH_TYPE_HYDRO_DENSITY_FINAL
SMOOTH_TYPE_HYDRO_GRADIENT
SMOOTH_TYPE_HYDRO_FLUX
SMOOTH_TYPE_HYDRO_STEP
SMOOTH_TYPE_HYDRO_FLUX_VEC
SMOOTH_TYPE_SN_FEEDBACK
SMOOTH_TYPE_BH_MERGER
SMOOTH_TYPE_BH_DRIFT
SMOOTH_TYPE_BH_STEP
SMOOTH_TYPE_CHEM_ENRICHMENT
- PKDGRAV.window_correction(grid_index=0, order=3)¶
Apply window correction to the grid
- Parameters:
grid_index (integer) – which grid number to use
order (integer) – mass assignment order, 0=NGP, 1=CIC, 2=TSC, 3=PCS
- PKDGRAV.write_array(filename, field)¶
Writes an array to a file.
- Parameters:
filename (str) – the name of the file
field (number) – the field to write. Values are:
OUT_DENSITY_ARRAY
OUT_POT_ARRAY
OUT_AMAG_ARRAY
OUT_IMASS_ARRAY
OUT_RUNG_ARRAY
OUT_DIVV_ARRAY
OUT_VELDISP2_ARRAY
OUT_VELDISP_ARRAY
OUT_PHASEDENS_ARRAY
OUT_SOFT_ARRAY
OUT_POS_VECTOR
OUT_VEL_VECTOR
OUT_ACCEL_VECTOR
OUT_MEANVEL_VECTOR
OUT_IORDER_ARRAY
OUT_C_ARRAY
OUT_HSPH_ARRAY
OUT_RUNGDEST_ARRAY
OUT_MARKED_ARRAY
OUT_CACHEFLUX_ARRAY
OUT_CACHECOLL_ARRAY
OUT_AVOIDEDFLUXES_ARRAY
OUT_COMPUTEDFLUXES_ARRAY
OUT_HOP_STATS
OUT_GROUP_ARRAY
OUT_GLOBALGID_ARRAY
OUT_BALL_ARRAY
OUT_PSGROUP_ARRAY
OUT_PSGROUP_STATS
- PKDGRAV.write_group_stats(filename, *, parallel=None)¶
Write group statistics to a file
- Parameters:
filename (str) – the name of the file
parallel (integer) – number of parallel files to write
Cosmology Interface¶
The classes SimpleCosmology and ClassCosmology expose the interface for determining cosmological values used in the code.
The complete list of classes and methods is as follows.
- class cosmology.ClassCosmology(file, L=1.0, As=0.0, ns=0.0, linear_species=[], power_species=[])¶
This is a cosmology class that uses CLASS. Use this by setting
- Parameters:
file (str) – File name of the CLASS/Concept HDF5 file
L (number) – box size in Mpc/h
As (number) – amplitude of the power spectrum \(A_s\)
ns (number) – spectral index \(n_s\)
linear_species (list) – list of linear species
power_species (list) – list of power species
- delta_linear(a, k)¶
The \(\delta\) perturbation of linear species.
- Parameters:
a (number) – expansion factor \(a\)
k (number) – wavenumber \(k\)
- Returns:
delta perturbations \(\delta(a,k)\) of linear species
- Return type:
number
- delta_matter(a, k)¶
The \(\delta\) perturbation of matter species.
- Parameters:
a (number) – expansion factor \(a\)
k (number) – wavenumber \(k\)
- Returns:
delta perturbations \(\delta(a,k)\) of matter species
- Return type:
number
- delta_power(a, k)¶
The \(\delta\) perturbation of linear species to include in the power spectrum.
- Parameters:
a (number) – expansion factor \(a\)
k (number) – wavenumber \(k\)
- Returns:
delta perturbations \(\delta(a,k)\) of linear species to include in the power spectrum
- Return type:
number
- delta_rho_linear(a, a_next, k)¶
Calculates \(\delta\bar{\rho}\) of linear species.
- Parameters:
a (number) – expansion factor \(a\)
a_next (number) – next expansion factor \(a'\)
k (number) – wavenumber \(k\)
- Returns:
\(\delta(a,k)\bar{\rho}(a)\) of linear species
- Return type:
number
- delta_rho_power(a, k)¶
Calculates \(\delta\bar{\rho}\) for the species included in the power spectrum.
- Parameters:
a (number) – expansion factor \(a\)
k (number) – wavenumber \(k\)
- Returns:
\(\delta(a,k)\bar{\rho}(a)\)
- Return type:
number
- rho_bar_linear(a)¶
Mean density \(\bar{\rho}\) of linear species.
- Parameters:
a (number) – expansion factor \(a\)
- Returns:
mean density \(\bar{\rho}(a)\) of linear species
- Return type:
number
- rho_bar_matter(a)¶
Mean density \(\bar{\rho}\) of matter.
- Parameters:
a (number) – expansion factor \(a\)
- Returns:
mean matter density \(\bar{\rho}(a)\)
- Return type:
number
- rho_bar_power(a)¶
Mean density \(\bar{\rho}\) of linear species to include in the power spectrum.
- Parameters:
a (number) – expansion factor \(a\)
- Returns:
mean density \(\bar{\rho}(a)\) of linear species to include in the power spectrum
- Return type:
number
- theta_matter(a, k)¶
The \(\theta\) perturbation of matter species.
- Parameters:
a (number) – expansion factor \(a\)
k (number) – wavenumber \(k\)
- Returns:
delta perturbations \(\theta(a,k)\) of matter species
- Return type:
number
- zeta(k)¶
- Parameters:
k (number) – wavenumber \(k\)
- Returns:
zeta \(\zeta(k)\)
- Return type:
number
- class cosmology.Cosmology¶
This is the base class for cosmology. It is not meant to be constructed directly, but rather through one of the subclasses
SimpleCosmologyorClassCosmology.- property As¶
- Returns:
amplitude of the power spectrum \(A_s\)
- Return type:
number
- property comoving¶
- Returns:
whether to use comoving or physical units
- Return type:
bool
- comoving_drift_factor(time, delta)¶
- Parameters:
time (number) – time \(t\)
delta (number) – delta time \(\Delta t\)
- Returns:
comoving drift factor
- Return type:
number
- comoving_growth(a)¶
- Parameters:
a (number) – expansion factor \(a\)
- Returns:
comoving growth
- Return type:
number
- comoving_kick_factor(time, delta)¶
- Parameters:
time (number) – time \(t\)
delta (number) – delta time \(\Delta t\)
- Returns:
comoving kick factor
- Return type:
number
- expansion_to_hubble(a)¶
- Parameters:
a (number) – expansion factor \(a\)
- Returns:
hubble parameter \(H(a)\)
- Return type:
number
- expansion_to_time(a)¶
- Parameters:
a (number) – expansion factor \(a\)
- Returns:
time \(t(a)\)
- Return type:
number
- property h¶
- Returns:
hubble parameter \(h\)
- Return type:
number
- property ns¶
- Returns:
spectral index \(n_s\)
- Return type:
number
- property omega_baryon¶
- Returns:
baryon density \(\Omega_b\)
- Return type:
number
- property omega_dark_energy¶
- Returns:
dark energy density \(\Omega_{DE}\)
- Return type:
number
- property omega_lambda¶
- Returns:
dark energy density \(\Omega_\Lambda\)
- Return type:
number
- property omega_matter¶
- Returns:
matter density \(\Omega_m\)
- Return type:
number
- property omega_radiation¶
- Returns:
radiation density \(\Omega_\gamma\)
- Return type:
number
- property pivot¶
- Returns:
pivot scale for the power spectrum \(k_p\)
- Return type:
number
- radiation_matter_equivalence()¶
- Returns:
radiation-matter equivalence \(a_{eq}\)
- Return type:
number
- property running¶
- Returns:
running of the spectral index \(dn_s/d\ln k\)
- Return type:
number
- property sigma8¶
- Returns:
sigma8 \(\sigma_8\)
- Return type:
number
- time_to_expansion(time)¶
- Parameters:
time (number) – time \(t\)
- Returns:
expansion factor \(a(t)\)
- Return type:
number
- time_to_hubble(time)¶
- Parameters:
time (number) – time \(t\)
- Returns:
hubble parameter \(H(t)\)
- Return type:
number
- property w0¶
- Returns:
dark energy equation of state parameter \(w_0\)
- Return type:
number
- property wa¶
- Returns:
dark energy equation of state parameter \(w_a\)
- Return type:
number
- class cosmology.SimpleCosmology(*, omega_matter, omega_lambda=None, omega_baryon=0.0, omega_radiation=0.0, omega_dark_energy=0.0, w0=-1.0, wa=0.0, running=0.0, pivot=0.05, comoving=True, sigma8=0.0, As=0.0, ns=0.0, h=1.0, k=None, Tk=None, transfer_file=None, usecols=(0, 1))¶
This is a simple cosmology class. Use this by setting
- Parameters:
omega_matter (number) – matter density \(\Omega_m\)
omega_lambda (number) – dark energy density \(\Omega_\Lambda\)
omega_baryon (number) – baryon density \(\Omega_b\)
omega_radiation (number) – radiation density \(\Omega_\gamma\)
omega_dark_energy (number) – dark energy density \(\Omega_{DE}\)
w0 (number) – dark energy equation of state parameter \(w_0\)
wa (number) – dark energy equation of state parameter \(w_a\)
running (number) – running of the spectral index \(dn_s/d\ln k\)
pivot (number) – pivot scale for the power spectrum \(k_p\)
comoving (bool) – whether to use comoving or physical units
sigma8 (number) – sigma8 \(\sigma_8\)
As (number) – amplitude of the power spectrum \(A_s\)
ns (number) – spectral index \(n_s\)
h (number) – hubble parameter \(h\)
k (number) – array of wavenumbers \(k\)
Tk (number) – array of transfer functions \(T(k)\)
transfer_file (str) – file name of the transfer function (input to numpy.loadtxt)
usecols (tuple) – columns to use in the transfer function file (passed to numpy.loadtxt)