diff --git a/src/parallel_io.c b/src/parallel_io.c index c7a0265ee4648e6ac13f8a188a6425d96793beda..489b2fd9207d44e5ee157d046008631edd12aca8 100644 --- a/src/parallel_io.c +++ b/src/parallel_io.c @@ -428,13 +428,28 @@ void prepareArray(struct engine* e, hid_t grp, char* fileName, FILE* xmfFile, /* Write unit conversion factors for this data set */ char buffer[FIELD_BUFFER_SIZE]; units_cgs_conversion_string(buffer, snapshot_units, props.units); - io_write_attribute_d( - h_data, "CGS conversion factor", - units_cgs_conversion_factor(snapshot_units, props.units)); + float baseUnitsExp[5]; + units_get_base_unit_exponents_array(baseUnitsExp, props.units); + const float a_factor_exp = units_a_factor(snapshot_units, props.units); + io_write_attribute_f(h_data, "U_M exponent", baseUnitsExp[UNIT_MASS]); + io_write_attribute_f(h_data, "U_L exponent", baseUnitsExp[UNIT_LENGTH]); + io_write_attribute_f(h_data, "U_t exponent", baseUnitsExp[UNIT_TIME]); + io_write_attribute_f(h_data, "U_I exponent", baseUnitsExp[UNIT_CURRENT]); + io_write_attribute_f(h_data, "U_T exponent", baseUnitsExp[UNIT_TEMPERATURE]); io_write_attribute_f(h_data, "h-scale exponent", 0); - io_write_attribute_f(h_data, "a-scale exponent", - units_a_factor(snapshot_units, props.units)); - io_write_attribute_s(h_data, "Conversion factor", buffer); + io_write_attribute_f(h_data, "a-scale exponent", a_factor_exp); + io_write_attribute_s(h_data, "Expression for physical CGS units", buffer); + + /* Write the actual number this conversion factor corresponds to */ + const float factor = units_cgs_conversion_factor(snapshot_units, props.units); + io_write_attribute_d( + h_data, + "Conversion factor to CGS (not including cosmological corrections)", + factor); + io_write_attribute_d( + h_data, + "Conversion factor to phyical CGS (including cosmological corrections)", + factor * pow(e->cosmology->a, a_factor_exp)); /* Add a line to the XMF */ if (xmfFile != NULL) diff --git a/src/serial_io.c b/src/serial_io.c index c3afef1283e4bc72eb7d75c1dd26c4517d2d97c3..8f2c011061cd9b2ce9b828c83f8d451ddb804f10 100644 --- a/src/serial_io.c +++ b/src/serial_io.c @@ -313,13 +313,28 @@ void prepareArray(const struct engine* e, hid_t grp, char* fileName, /* Write unit conversion factors for this data set */ char buffer[FIELD_BUFFER_SIZE]; units_cgs_conversion_string(buffer, snapshot_units, props.units); - io_write_attribute_d( - h_data, "CGS conversion factor", - units_cgs_conversion_factor(snapshot_units, props.units)); + float baseUnitsExp[5]; + units_get_base_unit_exponents_array(baseUnitsExp, props.units); + const float a_factor_exp = units_a_factor(snapshot_units, props.units); + io_write_attribute_f(h_data, "U_M exponent", baseUnitsExp[UNIT_MASS]); + io_write_attribute_f(h_data, "U_L exponent", baseUnitsExp[UNIT_LENGTH]); + io_write_attribute_f(h_data, "U_t exponent", baseUnitsExp[UNIT_TIME]); + io_write_attribute_f(h_data, "U_I exponent", baseUnitsExp[UNIT_CURRENT]); + io_write_attribute_f(h_data, "U_T exponent", baseUnitsExp[UNIT_TEMPERATURE]); io_write_attribute_f(h_data, "h-scale exponent", 0); - io_write_attribute_f(h_data, "a-scale exponent", - units_a_factor(snapshot_units, props.units)); - io_write_attribute_s(h_data, "Conversion factor", buffer); + io_write_attribute_f(h_data, "a-scale exponent", a_factor_exp); + io_write_attribute_s(h_data, "Expression for physical CGS units", buffer); + + /* Write the actual number this conversion factor corresponds to */ + const float factor = units_cgs_conversion_factor(snapshot_units, props.units); + io_write_attribute_d( + h_data, + "Conversion factor to CGS (not including cosmological corrections)", + factor); + io_write_attribute_d( + h_data, + "Conversion factor to phyical CGS (including cosmological corrections)", + factor * pow(e->cosmology->a, a_factor_exp)); /* Close everything */ H5Pclose(h_prop); diff --git a/src/single_io.c b/src/single_io.c index 4125b1ba91551a63af98e0286aa3d80be076716d..54761878887d5cf1d8173a188fed3a0ba72f9297 100644 --- a/src/single_io.c +++ b/src/single_io.c @@ -325,13 +325,28 @@ void writeArray(const struct engine* e, hid_t grp, char* fileName, /* Write unit conversion factors for this data set */ char buffer[FIELD_BUFFER_SIZE]; units_cgs_conversion_string(buffer, snapshot_units, props.units); - io_write_attribute_d( - h_data, "CGS conversion factor", - units_cgs_conversion_factor(snapshot_units, props.units)); + float baseUnitsExp[5]; + units_get_base_unit_exponents_array(baseUnitsExp, props.units); + const float a_factor_exp = units_a_factor(snapshot_units, props.units); + io_write_attribute_f(h_data, "U_M exponent", baseUnitsExp[UNIT_MASS]); + io_write_attribute_f(h_data, "U_L exponent", baseUnitsExp[UNIT_LENGTH]); + io_write_attribute_f(h_data, "U_t exponent", baseUnitsExp[UNIT_TIME]); + io_write_attribute_f(h_data, "U_I exponent", baseUnitsExp[UNIT_CURRENT]); + io_write_attribute_f(h_data, "U_T exponent", baseUnitsExp[UNIT_TEMPERATURE]); io_write_attribute_f(h_data, "h-scale exponent", 0); - io_write_attribute_f(h_data, "a-scale exponent", - units_a_factor(snapshot_units, props.units)); - io_write_attribute_s(h_data, "Conversion factor", buffer); + io_write_attribute_f(h_data, "a-scale exponent", a_factor_exp); + io_write_attribute_s(h_data, "Expression for physical CGS units", buffer); + + /* Write the actual number this conversion factor corresponds to */ + const float factor = units_cgs_conversion_factor(snapshot_units, props.units); + io_write_attribute_d( + h_data, + "Conversion factor to CGS (not including cosmological corrections)", + factor); + io_write_attribute_d( + h_data, + "Conversion factor to phyical CGS (including cosmological corrections)", + factor * pow(e->cosmology->a, a_factor_exp)); /* Free and close everything */ swift_free("writebuff", temp); diff --git a/src/units.c b/src/units.c index 55a1c9e1bdf62fce97b22189b4741be0fea2abb4..59596544fd6aa44abf6fafa28ea7ec71d9c703d6 100644 --- a/src/units.c +++ b/src/units.c @@ -225,8 +225,18 @@ const char* units_get_base_unit_cgs_symbol(enum base_units baseUnit) { return ""; } -void units_get_base_unit_exponants_array(float baseUnitsExp[5], +/** + * @brief Return an array of exponents of the base units required to construct a + * given unit. + * + * @param baseUnitsExp (return) The array of base-unit exponents. + * @param unit The unit we want the exponents for. + */ +void units_get_base_unit_exponents_array(float baseUnitsExp[5], enum unit_conversion_factor unit) { + /* Init everything */ + for (int k = 0; k < 5; ++k) baseUnitsExp[k] = 0.f; + switch (unit) { case UNIT_CONV_NO_UNITS: break; @@ -398,7 +408,7 @@ double units_cgs_conversion_factor(const struct unit_system* us, enum unit_conversion_factor unit) { float baseUnitsExp[5] = {0.f}; - units_get_base_unit_exponants_array(baseUnitsExp, unit); + units_get_base_unit_exponents_array(baseUnitsExp, unit); return units_general_cgs_conversion_factor(us, baseUnitsExp); } @@ -412,7 +422,7 @@ float units_h_factor(const struct unit_system* us, enum unit_conversion_factor unit) { float baseUnitsExp[5] = {0.f}; - units_get_base_unit_exponants_array(baseUnitsExp, unit); + units_get_base_unit_exponents_array(baseUnitsExp, unit); return units_general_h_factor(us, baseUnitsExp); } @@ -426,7 +436,7 @@ float units_a_factor(const struct unit_system* us, enum unit_conversion_factor unit) { float baseUnitsExp[5] = {0.f}; - units_get_base_unit_exponants_array(baseUnitsExp, unit); + units_get_base_unit_exponents_array(baseUnitsExp, unit); return units_general_a_factor(us, baseUnitsExp); } @@ -439,7 +449,7 @@ void units_cgs_conversion_string(char* buffer, const struct unit_system* us, enum unit_conversion_factor unit) { float baseUnitsExp[5] = {0.f}; - units_get_base_unit_exponants_array(baseUnitsExp, unit); + units_get_base_unit_exponents_array(baseUnitsExp, unit); units_general_cgs_conversion_string(buffer, us, baseUnitsExp); } @@ -452,13 +462,13 @@ void units_cgs_conversion_string(char* buffer, const struct unit_system* us, * the desired quantity. See conversionFactor() for a working example */ double units_general_cgs_conversion_factor(const struct unit_system* us, - const float baseUnitsExponants[5]) { + const float baseUnitsExponents[5]) { double factor = 1.; for (int i = 0; i < 5; ++i) - if (baseUnitsExponants[i] != 0) + if (baseUnitsExponents[i] != 0) factor *= pow(units_get_base_unit(us, (enum base_units)i), - baseUnitsExponants[i]); + baseUnitsExponents[i]); return factor; } @@ -470,12 +480,12 @@ double units_general_cgs_conversion_factor(const struct unit_system* us, * the desired quantity. See conversionFactor() for a working example */ float units_general_h_factor(const struct unit_system* us, - const float baseUnitsExponants[5]) { + const float baseUnitsExponents[5]) { float factor_exp = 0.f; - factor_exp += -baseUnitsExponants[UNIT_MASS]; - factor_exp += -baseUnitsExponants[UNIT_LENGTH]; - factor_exp += -baseUnitsExponants[UNIT_TIME]; + factor_exp += -baseUnitsExponents[UNIT_MASS]; + factor_exp += -baseUnitsExponents[UNIT_LENGTH]; + factor_exp += -baseUnitsExponents[UNIT_TIME]; return factor_exp; } @@ -488,10 +498,10 @@ float units_general_h_factor(const struct unit_system* us, * the desired quantity. See conversionFactor() for a working example */ float units_general_a_factor(const struct unit_system* us, - const float baseUnitsExponants[5]) { + const float baseUnitsExponents[5]) { float factor_exp = 0.f; - factor_exp += baseUnitsExponants[UNIT_LENGTH]; + factor_exp += baseUnitsExponents[UNIT_LENGTH]; return factor_exp; } @@ -511,58 +521,63 @@ float units_general_a_factor(const struct unit_system* us, */ void units_general_cgs_conversion_string(char* buffer, const struct unit_system* us, - const float baseUnitsExponants[5]) { - char temp[20]; - const double a_exp = units_general_a_factor(us, baseUnitsExponants); + const float baseUnitsExponents[5]) { + char temp[32]; + const double a_exp = units_general_a_factor(us, baseUnitsExponents); const double h_exp = 0.; /* There are no h-factors in SWIFT outputs. */ /* Check whether we are unitless or not */ - char isAllNonZero = 1; + char isAllZero = 1; for (int i = 0; i < 5; ++i) - if (baseUnitsExponants[i] != 0.) isAllNonZero = 0; + if (baseUnitsExponents[i] != 0.) isAllZero = 0; - if (isAllNonZero) { + if (isAllZero) { sprintf(buffer, "[ - ] "); return; } /* Add a-factor */ - if (a_exp == 0) - sprintf(buffer, " "); - else if (a_exp == 1) + if (a_exp == 0) { + /* Nothing to print */ + } else if (a_exp == 1) { sprintf(buffer, "a "); - else if (remainder(a_exp, 1.) == 0) + } else if (remainder(a_exp, 1.) == 0) { sprintf(buffer, "a^%d ", (int)a_exp); - else + } else { sprintf(buffer, "a^%7.4f ", a_exp); + } /* Add h-factor */ - if (h_exp == 0) - sprintf(temp, " "); - else if (h_exp == 1) + if (h_exp == 0) { + /* Nothing to print */ + } else if (h_exp == 1) { sprintf(temp, "h "); - else if (remainder(h_exp, 1.) == 0) + strcat(buffer, temp); + } else if (remainder(h_exp, 1.) == 0) { sprintf(temp, "h^%d ", (int)h_exp); - else + strcat(buffer, temp); + } else { sprintf(temp, "h^%7.4f ", h_exp); - strcat(buffer, temp); + strcat(buffer, temp); + } /* Add conversion units */ for (int i = 0; i < 5; ++i) - if (baseUnitsExponants[i] != 0) { - if (baseUnitsExponants[i] == 0.) + if (baseUnitsExponents[i] != 0) { + if (baseUnitsExponents[i] == 0.) { sprintf(temp, " "); - else if (baseUnitsExponants[i] == 1.) + } else if (baseUnitsExponents[i] == 1.) { sprintf(temp, "%s ", units_get_base_unit_internal_symbol((enum base_units)i)); - else if (remainder(baseUnitsExponants[i], 1.) == 0) + } else if (remainder(baseUnitsExponents[i], 1.) == 0) { sprintf(temp, "%s^%d ", units_get_base_unit_internal_symbol((enum base_units)i), - (int)baseUnitsExponants[i]); - else + (int)baseUnitsExponents[i]); + } else { sprintf(temp, "%s^%7.4f ", units_get_base_unit_internal_symbol((enum base_units)i), - baseUnitsExponants[i]); + baseUnitsExponents[i]); + } strcat(buffer, temp); } @@ -570,20 +585,21 @@ void units_general_cgs_conversion_string(char* buffer, strcat(buffer, " [ "); for (int i = 0; i < 5; ++i) { - if (baseUnitsExponants[i] != 0) { - if (baseUnitsExponants[i] == 0.) + if (baseUnitsExponents[i] != 0) { + if (baseUnitsExponents[i] == 0.) { continue; - else if (baseUnitsExponants[i] == 1.) + } else if (baseUnitsExponents[i] == 1.) { sprintf(temp, "%s ", units_get_base_unit_cgs_symbol((enum base_units)i)); - else if (remainder(baseUnitsExponants[i], 1.) == 0) + } else if (remainder(baseUnitsExponents[i], 1.) == 0) { sprintf(temp, "%s^%d ", units_get_base_unit_cgs_symbol((enum base_units)i), - (int)baseUnitsExponants[i]); - else + (int)baseUnitsExponents[i]); + } else { sprintf(temp, "%s^%7.4f ", units_get_base_unit_cgs_symbol((enum base_units)i), - baseUnitsExponants[i]); + baseUnitsExponents[i]); + } strcat(buffer, temp); } } @@ -619,12 +635,12 @@ int units_are_equal(const struct unit_system* a, const struct unit_system* b) { */ double units_general_conversion_factor(const struct unit_system* from, const struct unit_system* to, - const float baseUnitsExponants[5]) { + const float baseUnitsExponents[5]) { const double from_cgs = - units_general_cgs_conversion_factor(from, baseUnitsExponants); + units_general_cgs_conversion_factor(from, baseUnitsExponents); const double to_cgs = - units_general_cgs_conversion_factor(to, baseUnitsExponants); + units_general_cgs_conversion_factor(to, baseUnitsExponents); return from_cgs / to_cgs; } @@ -644,7 +660,7 @@ double units_conversion_factor(const struct unit_system* from, float baseUnitsExp[5] = {0.f}; - units_get_base_unit_exponants_array(baseUnitsExp, unit); + units_get_base_unit_exponents_array(baseUnitsExp, unit); return units_general_conversion_factor(from, to, baseUnitsExp); } diff --git a/src/units.h b/src/units.h index cc87696ccfdeb40b7cf4f23ccfcd853f8c8c1e56..b1eeeb33d63e3aee072d67b09768b7d3aa508f06 100644 --- a/src/units.h +++ b/src/units.h @@ -117,31 +117,34 @@ double units_get_base_unit(const struct unit_system*, enum base_units); const char* units_get_base_unit_internal_symbol(enum base_units); const char* units_get_base_unit_cgs_symbol(enum base_units); +void units_get_base_unit_exponents_array(float baseUnitsExp[5], + enum unit_conversion_factor unit); + /* Cosmology factors */ float units_general_h_factor(const struct unit_system* us, - const float baseUnitsExponants[5]); + const float baseUnitsExponents[5]); float units_h_factor(const struct unit_system* us, enum unit_conversion_factor unit); float units_general_a_factor(const struct unit_system* us, - const float baseUnitsExponants[5]); + const float baseUnitsExponents[5]); float units_a_factor(const struct unit_system* us, enum unit_conversion_factor unit); /* Conversion to CGS */ double units_general_cgs_conversion_factor(const struct unit_system* us, - const float baseUnitsExponants[5]); + const float baseUnitsExponents[5]); double units_cgs_conversion_factor(const struct unit_system* us, enum unit_conversion_factor unit); void units_general_cgs_conversion_string(char* buffer, const struct unit_system* us, - const float baseUnitsExponants[5]); + const float baseUnitsExponents[5]); void units_cgs_conversion_string(char* buffer, const struct unit_system* us, enum unit_conversion_factor unit); /* Conversion between systems */ double units_general_conversion_factor(const struct unit_system* from, const struct unit_system* to, - const float baseUnitsExponants[5]); + const float baseUnitsExponents[5]); double units_conversion_factor(const struct unit_system* from, const struct unit_system* to, enum unit_conversion_factor unit);