diff --git a/src/physical_constants.c b/src/physical_constants.c
index d00d63df1c1e6a4821ce4ba50dfef9c4e0def9d2..b726678ec1c0f2c85e409da3898dc3b6b842750f 100644
--- a/src/physical_constants.c
+++ b/src/physical_constants.c
@@ -40,65 +40,68 @@ void phys_const_init(struct UnitSystem* us, struct phys_const* internal_const) {
 
   const float dimension_G[5] = {-1, 3, -2, 0, 0};
   internal_const->const_newton_G =
-      const_newton_G_cgs / units_general_conversion_factor(us, dimension_G);
+      const_newton_G_cgs / units_general_cgs_conversion_factor(us, dimension_G);
 
   const float dimension_c[5] = {0, 1, -1, 0, 0};
   internal_const->const_speed_light_c =
       const_speed_light_c_cgs /
-      units_general_conversion_factor(us, dimension_c);
+      units_general_cgs_conversion_factor(us, dimension_c);
 
   const float dimension_h[5] = {1, -2, -1, 0, 0};
   internal_const->const_planck_h =
-      const_planck_h_cgs / units_general_conversion_factor(us, dimension_h);
+      const_planck_h_cgs / units_general_cgs_conversion_factor(us, dimension_h);
   internal_const->const_planck_hbar =
-      const_planck_hbar_cgs / units_general_conversion_factor(us, dimension_h);
+      const_planck_hbar_cgs /
+      units_general_cgs_conversion_factor(us, dimension_h);
 
   const float dimension_k[5] = {1, 2, -2, 0, -1};
   internal_const->const_boltzmann_k =
-      const_boltzmann_k_cgs / units_general_conversion_factor(us, dimension_k);
+      const_boltzmann_k_cgs /
+      units_general_cgs_conversion_factor(us, dimension_k);
 
   const float dimension_thomson[5] = {0, 2, 0, 0, 0};
   internal_const->const_thomson_cross_section =
       const_thomson_cross_section_cgs /
-      units_general_conversion_factor(us, dimension_thomson);
+      units_general_cgs_conversion_factor(us, dimension_thomson);
 
   const float dimension_ev[5] = {1, 2, -2, 0, 0};
   internal_const->const_electron_volt =
       const_electron_volt_cgs /
-      units_general_conversion_factor(us, dimension_ev);
+      units_general_cgs_conversion_factor(us, dimension_ev);
 
   const float dimension_charge[5] = {0, 0, -1, 1, 0};
   internal_const->const_electron_charge =
       const_electron_charge_cgs /
-      units_general_conversion_factor(us, dimension_charge);
+      units_general_cgs_conversion_factor(us, dimension_charge);
 
   const float dimension_mass[5] = {1, 0, 0, 0, 0};
   internal_const->const_electron_mass =
       const_electron_mass_cgs /
-      units_general_conversion_factor(us, dimension_mass);
+      units_general_cgs_conversion_factor(us, dimension_mass);
   internal_const->const_proton_mass =
       const_proton_mass_cgs /
-      units_general_conversion_factor(us, dimension_mass);
+      units_general_cgs_conversion_factor(us, dimension_mass);
   internal_const->const_solar_mass =
       const_solar_mass_cgs /
-      units_general_conversion_factor(us, dimension_mass);
+      units_general_cgs_conversion_factor(us, dimension_mass);
   internal_const->const_earth_mass =
       const_earth_mass_cgs /
-      units_general_conversion_factor(us, dimension_mass);
+      units_general_cgs_conversion_factor(us, dimension_mass);
 
   const float dimension_time[5] = {0, 0, 1, 0, 0};
   internal_const->const_year =
-      const_year_cgs / units_general_conversion_factor(us, dimension_time);
+      const_year_cgs / units_general_cgs_conversion_factor(us, dimension_time);
 
   const float dimension_length[5] = {0, 1, 0, 0, 0};
   internal_const->const_astronomical_unit =
       const_astronomical_unit_cgs /
-      units_general_conversion_factor(us, dimension_length);
+      units_general_cgs_conversion_factor(us, dimension_length);
   internal_const->const_parsec =
-      const_parsec_cgs / units_general_conversion_factor(us, dimension_length);
+      const_parsec_cgs /
+      units_general_cgs_conversion_factor(us, dimension_length);
   internal_const->const_light_year =
       const_light_year_cgs /
-      units_general_conversion_factor(us, dimension_length);
+      units_general_cgs_conversion_factor(us, dimension_length);
 }
 
 void phys_const_print(struct phys_const* internal_const) {
diff --git a/src/units.c b/src/units.c
index 8e0ff08f1d92f47afbf44366acef7038fc8675c5..123725a856e33f4e2ba7c75766794c189736b057 100644
--- a/src/units.c
+++ b/src/units.c
@@ -90,10 +90,10 @@ double units_get_base_unit(const struct UnitSystem* us,
 }
 
 /**
- * @brief Returns the base unit symbol
+ * @brief Returns the base unit symbol used internally
  * @param baseUnit The base unit
  */
-const char* units_get_base_unit_symbol(enum BaseUnits baseUnit) {
+const char* units_get_base_unit_internal_symbol(enum BaseUnits baseUnit) {
   switch (baseUnit) {
     case UNIT_MASS:
       return "U_M";
@@ -115,7 +115,7 @@ const char* units_get_base_unit_symbol(enum BaseUnits baseUnit) {
  * @brief Returns the base unit symbol in the cgs system
  * @param baseUnit The base unit
  */
-const char* units_get_base_unit_CGS_symbol(enum BaseUnits baseUnit) {
+const char* units_get_base_unit_cgs_symbol(enum BaseUnits baseUnit) {
   switch (baseUnit) {
     case UNIT_MASS:
       return "g";
@@ -275,13 +275,13 @@ void units_get_base_unit_exponants_array(float baseUnitsExp[5],
  * @param us The system of units in use
  * @param unit The unit to convert
  */
-double units_conversion_factor(const struct UnitSystem* us,
-                               enum UnitConversionFactor unit) {
+double units_cgs_conversion_factor(const struct UnitSystem* us,
+                                   enum UnitConversionFactor unit) {
   float baseUnitsExp[5] = {0.f};
 
   units_get_base_unit_exponants_array(baseUnitsExp, unit);
 
-  return units_general_conversion_factor(us, baseUnitsExp);
+  return units_general_cgs_conversion_factor(us, baseUnitsExp);
 }
 
 /**
@@ -316,13 +316,13 @@ float units_a_factor(const struct UnitSystem* us,
  * @brief Returns a string containing the exponents of the base units making up
  * the conversion factors
  */
-void units_conversion_string(char* buffer, const struct UnitSystem* us,
-                             enum UnitConversionFactor unit) {
+void units_cgs_conversion_string(char* buffer, const struct UnitSystem* us,
+                                 enum UnitConversionFactor unit) {
   float baseUnitsExp[5] = {0.f};
 
   units_get_base_unit_exponants_array(baseUnitsExp, unit);
 
-  units_general_conversion_string(buffer, us, baseUnitsExp);
+  units_general_cgs_conversion_string(buffer, us, baseUnitsExp);
 }
 
 /**
@@ -332,8 +332,8 @@ void units_conversion_string(char* buffer, const struct UnitSystem* us,
  * @param baseUnitsExponants The exponent of each base units required to form
  * the desired quantity. See conversionFactor() for a working example
  */
-double units_general_conversion_factor(const struct UnitSystem* us,
-                                       const float baseUnitsExponants[5]) {
+double units_general_cgs_conversion_factor(const struct UnitSystem* us,
+                                           const float baseUnitsExponants[5]) {
   double factor = 1.;
   int i;
 
@@ -386,8 +386,9 @@ float units_general_a_factor(const struct UnitSystem* us,
  * @param baseUnitsExponants The exponent of each base units required to form
  * the desired quantity. See conversionFactor() for a working example
  */
-void units_general_conversion_string(char* buffer, const struct UnitSystem* us,
-                                     const float baseUnitsExponants[5]) {
+void units_general_cgs_conversion_string(char* buffer,
+                                         const struct UnitSystem* us,
+                                         const float baseUnitsExponants[5]) {
   char temp[14];
   double a_exp = units_general_a_factor(us, baseUnitsExponants);
   double h_exp = units_general_h_factor(us, baseUnitsExponants);
@@ -430,12 +431,12 @@ void units_general_conversion_string(char* buffer, const struct UnitSystem* us,
       if (baseUnitsExponants[i] == 0.)
         sprintf(temp, " ");
       else if (baseUnitsExponants[i] == 1.)
-        sprintf(temp, "%s ", units_get_base_unit_symbol(i));
+        sprintf(temp, "%s ", units_get_base_unit_internal_symbol(i));
       else if (remainder(baseUnitsExponants[i], 1.) == 0)
-        sprintf(temp, "%s^%d ", units_get_base_unit_symbol(i),
+        sprintf(temp, "%s^%d ", units_get_base_unit_internal_symbol(i),
                 (int)baseUnitsExponants[i]);
       else
-        sprintf(temp, "%s^%7.4f ", units_get_base_unit_symbol(i),
+        sprintf(temp, "%s^%7.4f ", units_get_base_unit_internal_symbol(i),
                 baseUnitsExponants[i]);
       strncat(buffer, temp, 12);
     }
@@ -448,12 +449,12 @@ void units_general_conversion_string(char* buffer, const struct UnitSystem* us,
       if (baseUnitsExponants[i] == 0.)
         continue;
       else if (baseUnitsExponants[i] == 1.)
-        sprintf(temp, "%s ", units_get_base_unit_CGS_symbol(i));
+        sprintf(temp, "%s ", units_get_base_unit_cgs_symbol(i));
       else if (remainder(baseUnitsExponants[i], 1.) == 0)
-        sprintf(temp, "%s^%d ", units_get_base_unit_CGS_symbol(i),
+        sprintf(temp, "%s^%d ", units_get_base_unit_cgs_symbol(i),
                 (int)baseUnitsExponants[i]);
       else
-        sprintf(temp, "%s^%7.4f ", units_get_base_unit_CGS_symbol(i),
+        sprintf(temp, "%s^%7.4f ", units_get_base_unit_cgs_symbol(i),
                 baseUnitsExponants[i]);
       strncat(buffer, temp, 12);
     }
@@ -461,3 +462,59 @@ void units_general_conversion_string(char* buffer, const struct UnitSystem* us,
 
   strncat(buffer, "]", 2);
 }
+
+/**
+ * @brief Are the two unit systems equal ?
+ *
+ * @param a The First #UnitSystem
+ * @param b The second #UnitSystem
+ * @return 1 if the systems are the same, 0 otherwise
+ */
+int units_are_equal(const struct UnitSystem* a, const struct UnitSystem* b) {
+
+  if (a->UnitMass_in_cgs != b->UnitMass_in_cgs) return 0;
+  if (a->UnitLength_in_cgs != b->UnitLength_in_cgs) return 0;
+  if (a->UnitTime_in_cgs != b->UnitTime_in_cgs) return 0;
+  if (a->UnitCurrent_in_cgs != b->UnitCurrent_in_cgs) return 0;
+  if (a->UnitTemperature_in_cgs != b->UnitTemperature_in_cgs) return 0;
+
+  return 1;
+}
+
+/**
+ * @brief Return the unit conversion factor between two systems
+ *
+ * @param from The #UnitSystem we are converting from
+ * @param to The #UnitSystem we are converting to
+ * @param baseUnitsExponants The exponent of each base units required to form
+ * the desired quantity. See conversionFactor() for a working example
+ */
+double units_general_conversion_factor(const struct UnitSystem* from,
+                                       const struct UnitSystem* to,
+                                       const float baseUnitsExponants[5]) {
+
+  const double from_cgs =
+      units_general_cgs_conversion_factor(from, baseUnitsExponants);
+  const double to_cgs =
+      units_general_cgs_conversion_factor(to, baseUnitsExponants);
+
+  return from_cgs / to_cgs;
+}
+
+/**
+ * @brief Return the unit conversion factor between two systems
+ *
+ * @param from The #UnitSystem we are converting from
+ * @param to The #UnitSystem we are converting to
+ * @return unit The unit we are converting
+ */
+double units_conversion_factor(const struct UnitSystem* from,
+                               const struct UnitSystem* to,
+                               enum UnitConversionFactor unit) {
+
+  float baseUnitsExp[5] = {0.f};
+
+  units_get_base_unit_exponants_array(baseUnitsExp, unit);
+
+  return units_general_conversion_factor(from, to, baseUnitsExp);
+}
diff --git a/src/units.h b/src/units.h
index 24e37e177480d7f84e41df1b73e2036aa00b7220..29b4563163b4fe224c881b1c1055f2cbfbcc95a1 100644
--- a/src/units.h
+++ b/src/units.h
@@ -94,13 +94,14 @@ enum UnitConversionFactor {
 
 void units_init(struct UnitSystem*, const struct swift_params*,
                 const char* category);
+int units_are_equal(const struct UnitSystem* a, const struct UnitSystem* b);
+
+/* Base units */
 double units_get_base_unit(const struct UnitSystem*, enum BaseUnits);
-const char* units_get_base_unit_symbol(enum BaseUnits);
-const char* units_get_base_unit_CGS_symbol(enum BaseUnits);
-double units_general_conversion_factor(const struct UnitSystem* us,
-                                       const float baseUnitsExponants[5]);
-double units_conversion_factor(const struct UnitSystem* us,
-                               enum UnitConversionFactor unit);
+const char* units_get_base_unit_internal_symbol(enum BaseUnits);
+const char* units_get_base_unit_cgs_symbol(enum BaseUnits);
+
+/* Cosmology factors */
 float units_general_h_factor(const struct UnitSystem* us,
                              const float baseUnitsExponants[5]);
 float units_h_factor(const struct UnitSystem* us,
@@ -109,9 +110,24 @@ float units_general_a_factor(const struct UnitSystem* us,
                              const float baseUnitsExponants[5]);
 float units_a_factor(const struct UnitSystem* us,
                      enum UnitConversionFactor unit);
-void units_general_conversion_string(char* buffer, const struct UnitSystem* us,
-                                     const float baseUnitsExponants[5]);
-void units_conversion_string(char* buffer, const struct UnitSystem* us,
-                             enum UnitConversionFactor unit);
+
+/* Conversion to CGS */
+double units_general_cgs_conversion_factor(const struct UnitSystem* us,
+                                           const float baseUnitsExponants[5]);
+double units_cgs_conversion_factor(const struct UnitSystem* us,
+                                   enum UnitConversionFactor unit);
+void units_general_cgs_conversion_string(char* buffer,
+                                         const struct UnitSystem* us,
+                                         const float baseUnitsExponants[5]);
+void units_cgs_conversion_string(char* buffer, const struct UnitSystem* us,
+                                 enum UnitConversionFactor unit);
+
+/* Conversion between systems */
+double units_general_conversion_factor(const struct UnitSystem* from,
+                                       const struct UnitSystem* to,
+                                       const float baseUnitsExponants[5]);
+double units_conversion_factor(const struct UnitSystem* from,
+                               const struct UnitSystem* to,
+                               enum UnitConversionFactor unit);
 
 #endif /* SWIFT_UNITS_H */