diff options
| -rw-r--r-- | subsurface-core/dive.c | 61 | ||||
| -rw-r--r-- | subsurface-core/dive.h | 1 |
2 files changed, 50 insertions, 12 deletions
diff --git a/subsurface-core/dive.c b/subsurface-core/dive.c index f6fbc3261..8c69e9a6a 100644 --- a/subsurface-core/dive.c +++ b/subsurface-core/dive.c @@ -850,24 +850,61 @@ static void update_min_max_temperatures(struct dive *dive, temperature_t tempera * At high pressures air becomes less compressible, and * does not follow the ideal gas law any more. * - * This tries to correct for that, becoming the same - * as to_ATM() at lower pressures. - * - * THIS IS A ROUGH APPROXIMATION! The real numbers will - * depend on the exact gas mix and temperature. + * NOTE! The compressibility doesn't just depend on the + * gas, but on temperature too. However, this currently + * just follows the 300K curve for air, and ignores the + * gasmix. And the temperature we don't really even have + * a way to try to take into account. */ -static double surface_volume_multiplier(pressure_t pressure) -{ - double bar = pressure.mbar / 1000.0; +#define ARRAY_SIZE(array) (sizeof(array)/sizeof(array[0])) +double gas_compressibility_factor(struct gasmix *gas, double bar) +{ + static const struct z_factor { + double bar, z_factor; + } air_table[] = { + { 1, 0.9999 }, + { 5, 0.9987 }, + { 10, 0.9974 }, + { 20, 0.9950 }, + { 40, 0.9917 }, + { 60, 0.9901 }, + { 80, 0.9903 }, + { 100, 0.9930 }, + { 150, 1.0074 }, + { 200, 1.0326 }, + { 250, 1.0669 }, + { 300, 1.1089 }, + { 400, 1.2073 }, + { 500, 1.3163 } + }; + const struct z_factor *n; + int i; + + for (i = 0; i < ARRAY_SIZE(air_table); i++) { + const struct z_factor *n = air_table+i; + double frac; + + if (n->bar < bar) + continue; + if (!i) + return n->z_factor; + + /* How far from this one? */ + frac = (n->bar - bar) / (n->bar - n[-1].bar); + + /* Silly linear interpolation */ + return frac*n[-1].z_factor + (1.0-frac)*n->z_factor; + } - if (bar > 200) - bar = 0.00038 * bar * bar + 0.51629 * bar + 81.542; - return bar_to_atm(bar); + /* Over 500 bar? We make shit up */ + return air_table[ARRAY_SIZE(air_table)-1].z_factor; } int gas_volume(cylinder_t *cyl, pressure_t p) { - return cyl->type.size.mliter * surface_volume_multiplier(p); + double bar = p.mbar / 1000.0; + double z_factor = gas_compressibility_factor(&cyl->gasmix, bar); + return cyl->type.size.mliter * bar_to_atm(bar) / z_factor; } /* diff --git a/subsurface-core/dive.h b/subsurface-core/dive.h index 08100a270..0f5c43451 100644 --- a/subsurface-core/dive.h +++ b/subsurface-core/dive.h @@ -129,6 +129,7 @@ extern int units_to_sac(double volume); /* Volume in mliter of a cylinder at pressure 'p' */ extern int gas_volume(cylinder_t *cyl, pressure_t p); +extern double gas_compressibility_factor(struct gasmix *gas, double bar); static inline int get_o2(const struct gasmix *mix) |