cooling_uvfluc.c File Reference

#include <mpi.h>
#include <string.h>
#include <math.h>
#include "cooling_rates.h"
#include "physconst.h"
#include "bigfile.h"
#include "bigfile-mpi.h"
#include "utils/mymalloc.h"
#include "utils/interp.h"
#include "utils/endrun.h"

Include dependency graph for cooling_uvfluc.c:

Go to the source code of this file.

Functions
static double *	read_big_array (const char *filename, const char *dataset, int *Nread)
void	init_uvf_table (const char *UVFluctuationFile, const int UVFlucLen, const double BoxSize, const double UnitLength_in_cm)
struct UVBG	get_local_UVBG (double redshift, const struct UVBG *const GlobalUVBG, const double *const Pos, const double *const PosOffset)
void	InitMetalCooling (const char *MetalCoolFile)
double	TableMetalCoolingRate (double redshift, double temp, double nHcgs)

Variables
struct {
int   enabled
Interp   interp
double *   Table
ptrdiff_t   Nside
}	UVF
struct {
int   CoolingNoMetal
int   NRedshift_bins
double *   Redshift_bins
int   NHydrogenNumberDensity_bins
double *   HydrogenNumberDensity_bins
int   NTemperature_bins
double *   Temperature_bins
double *   Lmet_table
Interp   interp
}	MetalCool

Function Documentation

get_local_UVBG()

struct UVBG get_local_UVBG	(	double	redshift,
		const struct UVBG *const	GlobalUVBG,
		const double *const	Pos,
		const double *const	PosOffset
	)

Definition at line 91 of file cooling_uvfluc.c.

{
    if(!UVF.enabled) {
        /* directly use the TREECOOL table if UVF is disabled */
        return *GlobalUVBG;
    }
 
    struct UVBG uvbg = {0};
 
    uvbg.self_shield_dens = GlobalUVBG->self_shield_dens;
 
    double corrpos[3];
    int i;
    for(i = 0; i < 3; i++)
        corrpos[i] = Pos[i] - PosOffset[i];
    double zreion = interp_eval_periodic(&UVF.interp, corrpos, UVF.Table);
    if(zreion < redshift) {
        uvbg.zreion = zreion;
        return uvbg;
    }
    memcpy(&uvbg, GlobalUVBG, sizeof(struct UVBG));
    uvbg.zreion = zreion;
    return uvbg;
}

References CM_PER_MPC, endrun(), interp_init(), interp_init_dim(), message(), read_big_array(), UnitLength_in_cm, and UVF.

Referenced by cooling_direct(), get_helium_neutral_fraction_sfreff(), get_neutral_fraction_sfreff(), and get_sfr_eeqos().

Here is the call graph for this function:

Here is the caller graph for this function:

init_uvf_table()

void init_uvf_table	(	const char *	UVFluctuationFile,
		const int	UVFlucLen,
		const double	BoxSize,
		const double	UnitLength_in_cm
	)

Definition at line 91 of file cooling_uvfluc.c.

{
    if(strnlen(UVFluctuationFile, UVFlucLen) == 0) {
        UVF.enabled = 0;
        return;
    }
 
    /* Open and validate the UV fluctuation file*/
    BigFile bf;
    BigBlock bh;
    if(0 != big_file_mpi_open(&bf, UVFluctuationFile, MPI_COMM_WORLD)) {
        endrun(0, "Failed to open snapshot at %s:%s\n", UVFluctuationFile,
                    big_file_get_error_message());
    }
 
    if(0 != big_file_mpi_open_block(&bf, &bh, "Zreion_Table", MPI_COMM_WORLD)) {
        endrun(0, "Failed to create block at %s:%s\n", "Header",
                    big_file_get_error_message());
    }
    double TableBoxSize;
    double ReionRedshift;
    if ((0 != big_block_get_attr(&bh, "Nmesh", &UVF.Nside, "u8", 1)) ||
        (0 != big_block_get_attr(&bh, "BoxSize", &TableBoxSize, "f8", 1)) ||
        (0 != big_block_get_attr(&bh, "Redshift", &ReionRedshift, "f8", 1)) ||
        (0 != big_block_mpi_close(&bh, MPI_COMM_WORLD))) {
        endrun(0, "Failed to close block: %s\n",
                    big_file_get_error_message());
    }
    big_file_mpi_close(&bf, MPI_COMM_WORLD);
    double BoxMpc = BoxSize * UnitLength_in_cm / CM_PER_MPC;
    if(fabs(TableBoxSize - BoxMpc) > BoxMpc * 1e-5)
        endrun(0, "Wrong UV fluctuation file! %s is for box size %g Mpc/h, but current box is %g Mpc/h\n", UVFluctuationFile, TableBoxSize, BoxMpc);
 
    message(0, "Using NON-UNIFORM UV BG fluctuations from %s. Median reionization redshift is %g\n", UVFluctuationFile, ReionRedshift);
    UVF.enabled = 1;
 
    int size;
    UVF.Table = read_big_array(UVFluctuationFile, "Zreion_Table", &size);
 
    if(UVF.Nside * UVF.Nside * UVF.Nside != size)
        endrun(0, "Corrupt UV Fluctuation table: Nside = %ld, but table is %ld != %ld^3\n", UVF.Nside, size, UVF.Nside);
 
    int64_t dims[] = {UVF.Nside, UVF.Nside, UVF.Nside};
    interp_init(&UVF.interp, 3, dims);
    interp_init_dim(&UVF.interp, 0, 0, BoxSize);
    interp_init_dim(&UVF.interp, 1, 0, BoxSize);
    interp_init_dim(&UVF.interp, 2, 0, BoxSize);
 
    if(UVF.Table[0] < 0.01 || UVF.Table[0] > 100.0) {
        endrun(0, "UV Fluctuation out of range: %g\n", UVF.Table[0]);
    }
}

Referenced by init_cooling_and_star_formation().

Here is the caller graph for this function:

InitMetalCooling()

void InitMetalCooling

(

const char *

MetalCoolFile

)

Definition at line 193 of file cooling_uvfluc.c.

{
    /* now initialize the metal cooling table from cloudy; we got this file
     * from vogelsberger's Arepo simulations; it is supposed to be
     * cloudy + UVB - H and He; look so.
     * the table contains only 1 Z_sun values. Need to be scaled to the
     * metallicity.
     *
     * */
    /* let's see if the Metal Cool File is magic NoMetal */
    if(strlen(MetalCoolFile) == 0) {
        MetalCool.CoolingNoMetal = 1;
        return;
    } else {
        MetalCool.CoolingNoMetal = 0;
    }
 
    int size;
    //This is never used if MetalCoolFile == ""
    double * tabbedmet = read_big_array(MetalCoolFile, "MetallicityInSolar_bins", &size);
 
    if(size != 1 || tabbedmet[0] != 0.0) {
        endrun(123, "MetalCool file %s is wrongly tabulated\n", MetalCoolFile);
    }
    myfree(tabbedmet);
 
    MetalCool.Redshift_bins = read_big_array(MetalCoolFile, "Redshift_bins", &MetalCool.NRedshift_bins);
    MetalCool.HydrogenNumberDensity_bins = read_big_array(MetalCoolFile, "HydrogenNumberDensity_bins", &MetalCool.NHydrogenNumberDensity_bins);
    MetalCool.Temperature_bins = read_big_array(MetalCoolFile, "Temperature_bins", &MetalCool.NTemperature_bins);
    MetalCool.Lmet_table = read_big_array(MetalCoolFile, "NetCoolingRate", &size);
 
    int64_t dims[] = {MetalCool.NRedshift_bins, MetalCool.NHydrogenNumberDensity_bins, MetalCool.NTemperature_bins};
 
    interp_init(&MetalCool.interp, 3, dims);
    interp_init_dim(&MetalCool.interp, 0, MetalCool.Redshift_bins[0], MetalCool.Redshift_bins[MetalCool.NRedshift_bins - 1]);
    interp_init_dim(&MetalCool.interp, 1, MetalCool.HydrogenNumberDensity_bins[0],
                    MetalCool.HydrogenNumberDensity_bins[MetalCool.NHydrogenNumberDensity_bins - 1]);
    interp_init_dim(&MetalCool.interp, 2, MetalCool.Temperature_bins[0],
                    MetalCool.Temperature_bins[MetalCool.NTemperature_bins - 1]);
}

References endrun(), interp_init(), interp_init_dim(), MetalCool, myfree, and read_big_array().

Referenced by init_cooling_rates().

Here is the call graph for this function:

Here is the caller graph for this function:

read_big_array()

static double* read_big_array ( const char *  filename, const char *  dataset, int *  Nread  )

static

Definition at line 27 of file cooling_uvfluc.c.

{
    int N;
    double * buffer=NULL;
    int ThisTask;
    MPI_Comm_rank(MPI_COMM_WORLD, &ThisTask);
 
    if(ThisTask == 0) {
        BigFile bf[1];
        BigBlockPtr ptr;
        BigBlock bb[1];
        BigArray array[1];
        size_t dims[2];
        big_file_open(bf, filename);
        if(0 != big_file_open_block(bf, bb, dataset)) {
            endrun(1, "Cannot open %s %s: %s\n", filename, dataset, big_file_get_error_message());
        }
 
        N = bb->size;
 
        if(dtype_itemsize(bb->dtype) != sizeof(double))
            endrun(1, "UVflucatuation file %s should contain double-precision data, contains %s\n", filename, bb->dtype);
 
        buffer = (double *) mymalloc("cooling_data", N * dtype_itemsize(bb->dtype) * bb->nmemb);
        dims[0] = N;
        dims[1] = bb->nmemb;
 
        big_array_init(array, buffer, bb->dtype, 2, dims, NULL);
        if(0 != big_block_seek(bb, &ptr, 0))
            endrun(1, "Failed to seek block %s %s: %s\n", filename, dataset, big_file_get_error_message());
 
        if(0 != big_block_read(bb, &ptr, array))
            endrun(1, "Failed to read %s %s: %s", filename, dataset, big_file_get_error_message());
        /* steal the buffer */
        big_block_close(bb);
        big_file_close(bf);
    }
 
    MPI_Bcast(&N, 1, MPI_INT, 0, MPI_COMM_WORLD);
    if(ThisTask != 0)
        buffer = (double *) mymalloc("cooling_data",N * sizeof(double));
 
    MPI_Bcast(buffer, N, MPI_DOUBLE, 0, MPI_COMM_WORLD);
 
    *Nread = N;
    return buffer;
}

References endrun(), mymalloc, and ThisTask.

Referenced by get_local_UVBG(), and InitMetalCooling().

Here is the call graph for this function:

Here is the caller graph for this function:

TableMetalCoolingRate()

double TableMetalCoolingRate	(	double	redshift,
		double	temp,
		double	nHcgs
	)

Definition at line 235 of file cooling_uvfluc.c.

{
    if(MetalCool.CoolingNoMetal)
        return 0;
 
    double lognH = log10(nHcgs);
    double logT = log10(temp);
 
    double x[] = {redshift, lognH, logT};
    int status[3];
    double rate = interp_eval(&MetalCool.interp, x, MetalCool.Lmet_table, status);
    /* XXX: in case of very hot / very dense we just use whatever the table says at
     * the limit. should be OK. */
    return rate;
}

References interp_eval(), and MetalCool.

Referenced by get_heatingcooling_rate().

Here is the call graph for this function:

Here is the caller graph for this function:

Variable Documentation

CoolingNoMetal

int CoolingNoMetal

Definition at line 177 of file cooling_uvfluc.c.

enabled

int enabled

Definition at line 18 of file cooling_uvfluc.c.

HydrogenNumberDensity_bins

double* HydrogenNumberDensity_bins

Definition at line 182 of file cooling_uvfluc.c.

interp

Interp interp

Definition at line 19 of file cooling_uvfluc.c.

Referenced by mass_yield(), metal_return_copy(), metal_yield(), setup_metal_table_interp(), and test_yields().

Lmet_table

double* Lmet_table

Definition at line 187 of file cooling_uvfluc.c.

struct { ... } MetalCool

Referenced by InitMetalCooling(), TableMetalCoolingRate(), test_DoCooling(), test_heatingcooling_rate(), test_rate_network(), test_recomb_rates(), and test_uvbg_loader().

NHydrogenNumberDensity_bins

int NHydrogenNumberDensity_bins

Definition at line 181 of file cooling_uvfluc.c.

NRedshift_bins

int NRedshift_bins

Definition at line 178 of file cooling_uvfluc.c.

Nside

ptrdiff_t Nside

Definition at line 21 of file cooling_uvfluc.c.

NTemperature_bins

int NTemperature_bins

Definition at line 184 of file cooling_uvfluc.c.

Redshift_bins

double* Redshift_bins

Definition at line 179 of file cooling_uvfluc.c.

Table

double* Table

Definition at line 20 of file cooling_uvfluc.c.

Temperature_bins

double* Temperature_bins

Definition at line 185 of file cooling_uvfluc.c.

struct { ... } UVF

Referenced by get_local_UVBG().