From bb6ceba4ac0310c504f931674a71e9e77c0afa1c Mon Sep 17 00:00:00 2001
From: "Robert C. Helling" <helling@atdotde.de>
Date: Fri, 12 May 2017 15:36:24 +0200
Subject: [PATCH] Compute and display gas density

This appears to be critical for work of breathing so it might be
worthwhile to compute. So far only in infobox.

For background, see

https://www.youtube.com/watch?v=QBajM3xmOtc

Signed-off-by: Robert C. Helling <helling@atdotde.de>
---
 core/dive.h      | 2 ++
 core/gas-model.c | 8 ++++++++
 core/profile.c   | 6 +++++-
 core/profile.h   | 1 +
 core/units.h     | 2 +-
 5 files changed, 17 insertions(+), 2 deletions(-)

diff --git a/core/dive.h b/core/dive.h
index fa89edc44..c2ae1e36d 100644
--- a/core/dive.h
+++ b/core/dive.h
@@ -139,6 +139,8 @@ extern int units_to_sac(double volume);
 extern int gas_volume(cylinder_t *cyl, pressure_t p);
 extern double gas_compressibility_factor(struct gasmix *gas, double bar);
 extern double isothermal_pressure(struct gasmix *gas, double p1, int volume1, int volume2);
+extern double gas_density(struct gasmix *gas, int pressure);
+
 
 
 static inline int get_o2(const struct gasmix *mix)
diff --git a/core/gas-model.c b/core/gas-model.c
index dd84e9b9d..49c677967 100644
--- a/core/gas-model.c
+++ b/core/gas-model.c
@@ -9,6 +9,8 @@
 #define virial_m1(C, x1, x2, x3) (C[0]*x1+C[1]*x2+C[2]*x3)
 
 /*
+ * Z = pV/nRT
+ *
  * Cubic virial least-square coefficients for O2/N2/He based on data from
  *
  *    PERRY’S CHEMICAL ENGINEERS’ HANDBOOK SEVENTH EDITION
@@ -73,3 +75,9 @@ double isothermal_pressure(struct gasmix *gas, double p1, int volume1, int volum
 
 	return p_ideal * gas_compressibility_factor(gas, p_ideal);
 }
+
+inline double gas_density(struct gasmix *gas, int pressure) {
+	int density =  gas->he.permille * HE_DENSITY + gas->o2.permille * O2_DENSITY + (1000 - gas->he.permille - gas->o2.permille) * N2_DENSITY;
+
+	return density * (double) pressure / gas_compressibility_factor(gas, pressure / 1000.0) / SURFACE_PRESSURE / 1000000.0;
+}
diff --git a/core/profile.c b/core/profile.c
index 86bc974b8..d46904cbe 100644
--- a/core/profile.c
+++ b/core/profile.c
@@ -1151,6 +1151,7 @@ static void calculate_gas_information_new(struct dive *dive, struct plot_info *p
 				       entry->pressures.n2 / amb_pressure * N2_DENSITY +
 				       entry->pressures.he / amb_pressure * HE_DENSITY) /
 				      (O2_IN_AIR * O2_DENSITY + N2_IN_AIR * N2_DENSITY) * 1000 - 10000;
+		entry->density = gas_density(&dive->cylinder[cylinderindex].gasmix, depth_to_mbar(entry->depth, dive));
 		if (entry->mod < 0)
 			entry->mod = 0;
 		if (entry->ead < 0)
@@ -1292,7 +1293,7 @@ static void plot_string(struct plot_info *pi, struct plot_data *entry, struct me
 {
 	int pressurevalue, mod, ead, end, eadd;
 	const char *depth_unit, *pressure_unit, *temp_unit, *vertical_speed_unit;
-	double depthvalue, tempvalue, speedvalue, sacvalue;
+	double depthvalue, tempvalue, speedvalue, sacvalue, density;
 	int decimals;
 	const char *unit;
 
@@ -1327,15 +1328,18 @@ static void plot_string(struct plot_info *pi, struct plot_data *entry, struct me
 		put_format(b, translate("gettextFromC", "MOD: %d%s\n"), mod, depth_unit);
 	}
 	eadd = lrint(get_depth_units(lrint(entry->eadd), NULL, &depth_unit));
+
 	if (prefs.ead) {
 		switch (pi->dive_type) {
 		case NITROX:
 			ead = lrint(get_depth_units(lrint(entry->ead), NULL, &depth_unit));
 			put_format(b, translate("gettextFromC", "EAD: %d%s\nEADD: %d%s\n"), ead, depth_unit, eadd, depth_unit);
+			put_format(b, translate("gettextFromC", "density: %.1fg/l\n"), entry->density);
 			break;
 		case TRIMIX:
 			end = lrint(get_depth_units(lrint(entry->end), NULL, &depth_unit));
 			put_format(b, translate("gettextFromC", "END: %d%s\nEADD: %d%s\n"), end, depth_unit, eadd, depth_unit);
+			put_format(b, translate("gettextFromC", "density: %.1fg/l\n"), entry->density);
 			break;
 		case AIR:
 		case FREEDIVING:
diff --git a/core/profile.h b/core/profile.h
index 0fd284d47..079aaa1db 100644
--- a/core/profile.h
+++ b/core/profile.h
@@ -64,6 +64,7 @@ struct plot_data {
 	int bearing;
 	double ambpressure;
 	double gfline;
+	double density;
 };
 
 struct ev_select {
diff --git a/core/units.h b/core/units.h
index 7e4c2e7d2..df63dba49 100644
--- a/core/units.h
+++ b/core/units.h
@@ -14,7 +14,7 @@ extern "C" {
 #define O2_IN_AIR 209 // permille
 #define N2_IN_AIR 781
 #define O2_DENSITY 1429 // mg/Liter
-#define N2_DENSITY 1251
+#define N2_DENSITY 1165
 #define HE_DENSITY 179
 #define SURFACE_PRESSURE 1013 // mbar
 #define SURFACE_PRESSURE_STRING "1013"