diff --git a/ChemistryLib/CreateChemicalSolverInterface.cpp b/ChemistryLib/CreateChemicalSolverInterface.cpp
index 998581e63582affc2b027484444d60fae45a9c17..83f3afaa1b499129a0e4e94c6865286c3ef820d3 100644
--- a/ChemistryLib/CreateChemicalSolverInterface.cpp
+++ b/ChemistryLib/CreateChemicalSolverInterface.cpp
@@ -104,9 +104,7 @@ createChemicalSolverInterface<ChemicalSolver::Phreeqc>(
         config.getConfigSubtreeOptional("rates"));
 
     // surface
-    auto surface = PhreeqcIOData::createSurface(
-        //! \ogs_file_param{prj__chemical_system__surface}
-        config.getConfigSubtreeOptional("surface"));
+    auto const& surface = chemical_system->surface;
 
     // exchange
     auto const& exchangers = chemical_system->exchangers;
@@ -147,8 +145,8 @@ createChemicalSolverInterface<ChemicalSolver::Phreeqc>(
     return std::make_unique<PhreeqcIOData::PhreeqcIO>(
         *linear_solver, std::move(project_file_name),
         std::move(path_to_database), std::move(chemical_system),
-        std::move(reaction_rates), std::move(surface), std::move(user_punch),
-        std::move(output), std::move(dump), std::move(knobs));
+        std::move(reaction_rates), std::move(user_punch), std::move(output),
+        std::move(dump), std::move(knobs));
 }
 
 template <>
diff --git a/ChemistryLib/PhreeqcIO.cpp b/ChemistryLib/PhreeqcIO.cpp
index c2e1e1e138714dd5cd7dde6bfed080203a3dccb9..66ac3d659beebc6bb2f14d06f4ebb8f90a6ea01a 100644
--- a/ChemistryLib/PhreeqcIO.cpp
+++ b/ChemistryLib/PhreeqcIO.cpp
@@ -253,7 +253,6 @@ PhreeqcIO::PhreeqcIO(GlobalLinearSolver& linear_solver,
                      std::string&& database,
                      std::unique_ptr<ChemicalSystem>&& chemical_system,
                      std::vector<ReactionRate>&& reaction_rates,
-                     std::vector<SurfaceSite>&& surface,
                      std::unique_ptr<UserPunch>&& user_punch,
                      std::unique_ptr<Output>&& output,
                      std::unique_ptr<Dump>&& dump,
@@ -263,7 +262,6 @@ PhreeqcIO::PhreeqcIO(GlobalLinearSolver& linear_solver,
       _database(std::move(database)),
       _chemical_system(std::move(chemical_system)),
       _reaction_rates(std::move(reaction_rates)),
-      _surface(std::move(surface)),
       _user_punch(std::move(user_punch)),
       _output(std::move(output)),
       _dump(std::move(dump)),
@@ -527,7 +525,7 @@ std::ostream& operator<<(std::ostream& os, PhreeqcIO const& phreeqc_io)
                << "\n";
         }
 
-        auto const& surface = phreeqc_io._surface;
+        auto const& surface = phreeqc_io._chemical_system->surface;
         if (!surface.empty())
         {
             os << "SURFACE " << chemical_system_id + 1 << "\n";
@@ -617,7 +615,7 @@ std::istream& operator>>(std::istream& in, PhreeqcIO& phreeqc_io)
     auto const& output = *phreeqc_io._output;
     auto const& dropped_item_ids = output.dropped_item_ids;
 
-    auto const& surface = phreeqc_io._surface;
+    auto const& surface = phreeqc_io._chemical_system->surface;
     auto const& exchangers = phreeqc_io._chemical_system->exchangers;
 
     if (!surface.empty() && !exchangers.empty())
diff --git a/ChemistryLib/PhreeqcIO.h b/ChemistryLib/PhreeqcIO.h
index bca62ca832b530ebd32401b5eaf6575b45e66f47..5b48e35d61d3f6391ae74b5fc46aaa3838a83a33 100644
--- a/ChemistryLib/PhreeqcIO.h
+++ b/ChemistryLib/PhreeqcIO.h
@@ -27,7 +27,6 @@ namespace PhreeqcIOData
 struct ChemicalSystem;
 struct ReactionRate;
 struct Output;
-struct SurfaceSite;
 struct Dump;
 struct UserPunch;
 
@@ -39,7 +38,6 @@ public:
               std::string&& database,
               std::unique_ptr<ChemicalSystem>&& chemical_system,
               std::vector<ReactionRate>&& reaction_rates,
-              std::vector<SurfaceSite>&& surface,
               std::unique_ptr<UserPunch>&& user_punch,
               std::unique_ptr<Output>&& output,
               std::unique_ptr<Dump>&& dump,
@@ -111,7 +109,6 @@ private:
     std::string const _database;
     std::unique_ptr<ChemicalSystem> _chemical_system;
     std::vector<ReactionRate> const _reaction_rates;
-    std::vector<SurfaceSite> const _surface;
     std::unique_ptr<UserPunch> _user_punch;
     std::unique_ptr<Output> const _output;
     std::unique_ptr<Dump> const _dump;
diff --git a/ChemistryLib/PhreeqcIOData/ChemicalSystem.h b/ChemistryLib/PhreeqcIOData/ChemicalSystem.h
index 9e50e9cb5ff2ba2160f78cb412f2a564f5a1fb57..22b209ca175badf3ddde1d84126dbb055b6eb19c 100644
--- a/ChemistryLib/PhreeqcIOData/ChemicalSystem.h
+++ b/ChemistryLib/PhreeqcIOData/ChemicalSystem.h
@@ -11,12 +11,14 @@
 #pragma once
 
 #include <memory>
+#include <variant>
 #include <vector>
 
 #include "AqueousSolution.h"
 #include "EquilibriumReactant.h"
 #include "Exchange.h"
 #include "KineticReactant.h"
+#include "Surface.h"
 
 namespace ChemistryLib
 {
@@ -27,11 +29,14 @@ struct ChemicalSystem
     ChemicalSystem(std::unique_ptr<AqueousSolution>&& aqueous_solution_,
                    std::vector<KineticReactant>&& kinetic_reactants_,
                    std::vector<EquilibriumReactant>&& equilibrium_reactants_,
-                   std::vector<ExchangeSite>&& exchangers_)
+                   std::vector<ExchangeSite>&& exchangers_,
+                   std::vector<std::variant<DensityBasedSurfaceSite,
+                                            MoleBasedSurfaceSite>>&& surface_)
         : aqueous_solution(std::move(aqueous_solution_)),
           kinetic_reactants(std::move(kinetic_reactants_)),
           equilibrium_reactants(std::move(equilibrium_reactants_)),
-          exchangers(std::move(exchangers_))
+          exchangers(std::move(exchangers_)),
+          surface(std::move(surface_))
     {
     }
 
@@ -41,6 +46,8 @@ struct ChemicalSystem
     std::vector<KineticReactant> kinetic_reactants;
     std::vector<EquilibriumReactant> equilibrium_reactants;
     std::vector<ExchangeSite> exchangers;
+    std::vector<std::variant<DensityBasedSurfaceSite, MoleBasedSurfaceSite>>
+        surface;
 };
 }  // namespace PhreeqcIOData
 }  // namespace ChemistryLib
diff --git a/ChemistryLib/PhreeqcIOData/CreateChemicalSystem.cpp b/ChemistryLib/PhreeqcIOData/CreateChemicalSystem.cpp
index 6319c899a41d0e290c14e2449a75f34e91618911..1aef5665d1ea563ab0d045d99539690318fdc29e 100644
--- a/ChemistryLib/PhreeqcIOData/CreateChemicalSystem.cpp
+++ b/ChemistryLib/PhreeqcIOData/CreateChemicalSystem.cpp
@@ -16,9 +16,11 @@
 #include "CreateEquilibriumReactants.h"
 #include "CreateExchange.h"
 #include "CreateKineticReactant.h"
+#include "CreateSurface.h"
 #include "EquilibriumReactant.h"
 #include "Exchange.h"
 #include "KineticReactant.h"
+#include "Surface.h"
 
 namespace ChemistryLib
 {
@@ -47,10 +49,16 @@ std::unique_ptr<ChemicalSystem> createChemicalSystem(
         //! \ogs_file_param{prj__chemical_system__exchangers}
         config.getConfigSubtreeOptional("exchangers"), mesh);
 
+    // surface
+    auto surface = createSurface(
+        //! \ogs_file_param{prj__chemical_system__surface}
+        config.getConfigSubtreeOptional("surface"), mesh);
+
     return std::make_unique<ChemicalSystem>(std::move(aqueous_solution),
                                             std::move(kinetic_reactants),
                                             std::move(equilibrium_reactants),
-                                            std::move(exchangers));
+                                            std::move(exchangers),
+                                            std::move(surface));
 }
 }  // namespace PhreeqcIOData
 }  // namespace ChemistryLib