diff --git a/Applications/ApplicationsLib/ProjectData.cpp b/Applications/ApplicationsLib/ProjectData.cpp
index 7189bc26a34b93e2df8434e99bfa1ddec3c444d6..a6f3693eebf2c5a6063319f285208cd00dd6c0b9 100644
--- a/Applications/ApplicationsLib/ProjectData.cpp
+++ b/Applications/ApplicationsLib/ProjectData.cpp
@@ -1043,8 +1043,8 @@ void ProjectData::parseChemicalSystem(
                 component_transport_process->getProcessIDToComponentNameMap();
 
             _chemical_system = ChemistryLib::createPhreeqcIO(
-                _mesh_vec[0]->getNumberOfBaseNodes(),
-                process_id_to_component_name_map, *config, output_directory);
+                *_mesh_vec[0], process_id_to_component_name_map, *config,
+                output_directory);
         }
         else
         {
diff --git a/ChemistryLib/CreatePhreeqcIO.cpp b/ChemistryLib/CreatePhreeqcIO.cpp
index 61b542d0f3e6dcf0ba41713b13a47cafca5905d5..1799565dbbdfca8c840d0b0f82d4cc3621e6766c 100644
--- a/ChemistryLib/CreatePhreeqcIO.cpp
+++ b/ChemistryLib/CreatePhreeqcIO.cpp
@@ -11,6 +11,7 @@
 #include "BaseLib/Error.h"
 #include "CreateOutput.h"
 #include "CreatePhreeqcIO.h"
+#include "MeshLib/Mesh.h"
 #include "PhreeqcIO.h"
 #include "PhreeqcIOData/AqueousSolution.h"
 #include "PhreeqcIOData/CreateAqueousSolution.h"
@@ -24,13 +25,13 @@
 namespace ChemistryLib
 {
 std::unique_ptr<PhreeqcIO> createPhreeqcIO(
-    std::size_t const num_nodes,
+    MeshLib::Mesh const& mesh,
     std::vector<std::pair<int, std::string>> const&
         process_id_to_component_name_map,
     BaseLib::ConfigTree const& config,
     std::string const& output_directory)
 {
-    auto const& num_chemical_systems = num_nodes;
+    auto const num_chemical_systems = mesh.getNumberOfBaseNodes();
 
     // database
     //! \ogs_file_param{prj__chemical_system__database}
@@ -88,18 +89,14 @@ std::unique_ptr<PhreeqcIO> createPhreeqcIO(
         num_chemical_systems, aqueous_solution_per_chem_sys);
 
     // equilibrium phases
-    auto const equilibrium_phases_per_chem_sys = createEquilibriumPhases(
+    auto equilibrium_phases = createEquilibriumPhases(
         //! \ogs_file_param{prj__chemical_system__equilibrium_phases}
-        config.getConfigSubtreeOptional("equilibrium_phases"));
-    std::vector<std::vector<EquilibriumPhase>> equilibrium_phases(
-        num_chemical_systems, equilibrium_phases_per_chem_sys);
+        config.getConfigSubtreeOptional("equilibrium_phases"), mesh);
 
     // kinetic reactants
-    auto const kinetic_reactants_per_chem_sys = createKineticReactants(
+    auto kinetic_reactants = createKineticReactants(
         //! \ogs_file_param{prj__chemical_system__kinetic_reactants}
-        config.getConfigSubtreeOptional("kinetic_reactants"));
-    std::vector<std::vector<KineticReactant>> kinetic_reactants(
-        num_chemical_systems, kinetic_reactants_per_chem_sys);
+        config.getConfigSubtreeOptional("kinetic_reactants"), mesh);
 
     // rates
     auto reaction_rates = createReactionRates(
diff --git a/ChemistryLib/CreatePhreeqcIO.h b/ChemistryLib/CreatePhreeqcIO.h
index 46652a3d1bf48296d30c58bcd5413c5dcb346eb1..addc2bfabef1468684a11d39060efbe99cea6823 100644
--- a/ChemistryLib/CreatePhreeqcIO.h
+++ b/ChemistryLib/CreatePhreeqcIO.h
@@ -13,6 +13,11 @@
 
 #include "PhreeqcIO.h"
 
+namespace MeshLib
+{
+class Mesh;
+}
+
 namespace BaseLib
 {
 class ConfigTree;
@@ -26,7 +31,7 @@ class Process;
 namespace ChemistryLib
 {
 std::unique_ptr<PhreeqcIO> createPhreeqcIO(
-    std::size_t const num_nodes,
+    MeshLib::Mesh const& mesh,
     std::vector<std::pair<int, std::string>> const&
         process_id_to_component_name_map,
     BaseLib::ConfigTree const& config,