philippkraft#67 - disabled Freundlich Adsorption. Broken Code is stil…

…l in extra files. Fixed doc warnings
jlu-ilr-hydro · Jan 19, 2024 · d266707 · d266707
1 parent f9c5129
commit d266707
Show file tree

Hide file tree

Showing 20 changed files with 3,139 additions and 3,956 deletions.
diff --git a/cmf/__init__.py b/cmf/__init__.py
@@ -22,7 +22,7 @@
 from .timetools import StopWatch, datetime_to_cmf, timerange
 
 
-__version__ = '2.0.0b8'
+__version__ = '2.0.0b10'
 
 from .cmf_core import connect_cells_with_flux as __ccwf
 

diff --git a/cmf/cmf_core.py b/cmf/cmf_core.py
diff --git a/cmf/cmf_core_src/CMakeLists.txt b/cmf/cmf_core_src/CMakeLists.txt
@@ -1,4 +1,4 @@
-cmake_minimum_required(VERSION 3.1)
+cmake_minimum_required(VERSION 3.15)
 project(cmf_core)
 
 # Super necessary! Without this cannot link to cmf_wrap.cxx to create a shared object!
@@ -38,6 +38,7 @@ add_library(cmf_core STATIC
         math/odesystem.cpp
         math/sparse_struct.cpp
         water/adsorption.cpp
+        # water/freundlich_adsorption.cpp
         water/boundary_condition.cpp
         water/collections.cpp
         water/flux_connection.cpp

diff --git a/cmf/cmf_core_src/cmf_wrap.cpp b/cmf/cmf_core_src/cmf_wrap.cpp
diff --git a/cmf/cmf_core_src/docstrings.i b/cmf/cmf_core_src/docstrings.i
@@ -1576,7 +1576,7 @@ and depth are ignored.
 'P' PipeReach, depth is ignored, width is the diameter of the pipe
 
 'S' SWATReachType, a trapezoid flow cross section, as used in the SWAT
-model, width (bank width) and depth are used the reach type
+model, width (bank width) and depth are used
 
 Parameters:
 -----------
@@ -4628,8 +4628,6 @@ cmf.fit_retention_curve.FitVanGenuchtenMualem.__init__";
 The connections in cmf hold the processes for the calculation of
 fluxes between water storages and model boundaries.
 
-Todo Elaborate on this
-
 Represents a connection between flux_nodes, where water fluxes occur.
 
 C++ includes: flux_connection.h ";
@@ -4965,44 +4963,7 @@ cmf::water::flux_connection::to_string";
 
 
 // File: classcmf_1_1water_1_1_freundlich_adsorbtion.xml
-%feature("docstring") cmf::water::FreundlichAdsorbtion "
-
-BROKEN: This class calculates the adsorption equilibrium between
-sorbat and sorbent using the Freundlich isotherme.
-
-Freundlich isotherme:
-
-
-
-.. math::
-
-    \\\\frac{x_{ad}}{m} = K c^n
-
-where  :math:`x_{ad} = x_{tot} - x_{free}` is the adsorbed tracer mass  :math:`x_{tot}` is the total tracer
-mass
-
-:math:`x_{free}` is the dissolved tracer mass
-
-:math:`m` is the mass of the sorbent in the same unit as the tracer mass
-
-:math:`K` is the Freundlich sorption coefficient
-
-:math:`c = \\\\frac{x_{free}}{V}` is the concentration of the tracer in
-tracer mass per m3
-
-:math:`n` is the Freundlich exponent
-
-CMF stores in a solute storage the total mass of a tracer and needs to
-calculate the free tracer mass. The eq. above can not be rearanged to
-get :math:`x_{free}` from :math:`x_{tot}`. Instead, the value is iterated
-usingregula falsi. If n is near to 1, using LinearAdsorption will
-speed up your calculations.
-
-The simplest physically based adsorption model by Langmuir (
-LangmuirAdsorption) has also a analytical solution and is hence
-calculated faster then Freundlich.
-
-C++ includes: adsorption.h ";
+%feature("docstring") cmf::water::FreundlichAdsorbtion "";
 
 %feature("docstring")
 cmf::water::FreundlichAdsorbtion::FreundlichAdsorbtion "FreundlichAdsorbtion(const FreundlichAdsorbtion &other)
@@ -5030,39 +4991,13 @@ cmf::water::FreundlichAdsorbtion::~FreundlichAdsorbtion "virtual
 cmf::water::FreundlichAdsorbtion::~FreundlichAdsorbtion";
 
 %feature("docstring")  cmf::water::FreundlichAdsorbtion::copy "FreundlichAdsorbtion * copy(real m=-1) const
-cmf::water::FreundlichAdsorbtion::copy returns a copy of the
-Adsorption object.
-
-If the adsorption is depending on the sorbent mass, you can give a
-positive value for the sorbent mass m. If the value is not given or
-negative, m is used from the original object. ";
+cmf::water::FreundlichAdsorbtion::copy";
 
 %feature("docstring")  cmf::water::FreundlichAdsorbtion::freesolute "real freesolute(real xt, real V) const
-cmf::water::FreundlichAdsorbtion::freesolute Returns the mass of
-dissolved tracer as a function of the total tracer mass in the solute
-storage and the water volume.
-
-Parameters:
------------
-
-xt:   :math:`x_t` the total tracer mass in the storage
-
-V:   :math:`V m^3` the water volume in the storage
-
-:math:`x_f` the dissolved mass of the tracer ";
+cmf::water::FreundlichAdsorbtion::freesolute";
 
 %feature("docstring")  cmf::water::FreundlichAdsorbtion::totalsolute "real totalsolute(real xf, real V) const
-cmf::water::FreundlichAdsorbtion::totalsolute Returns the total mass
-of the tracer from the dissolved concetration in tracer unit/m3.
-
-Parameters:
------------
-
-xf:   :math:`x_f` the dissolved tracer mass in the storage
-
-V:   :math:`V m^3` the water volume in the storage
-
-:math:`x_t` the total mass of the tracer ";
+cmf::water::FreundlichAdsorbtion::totalsolute";
 
 
 // File: classcmf_1_1upslope_1_1connections_1_1_gradient_macro_flow.xml
@@ -17994,6 +17929,9 @@ ymax=1) ";
 // File: flux__node_8h.xml
 
 
+// File: freundlich__adsorption_8h.xml
+
+
 // File: reaction_8h.xml
 
 
@@ -18075,18 +18013,6 @@ ymax=1) ";
 // File: overview_8md.xml
 
 
-// File: customization_8md.xml
-
-
-// File: extensions_8md.xml
-
-
-// File: tricks_8md.xml
-
-
-// File: _r_e_a_d_m_e_8md.xml
-
-
 // File: ems-2011-paper_8md.xml
 
 
@@ -18291,18 +18217,6 @@ ymax=1) ";
 // File: contrib_issues.xml
 
 
-// File: md_documentation_2doxygen-awesome-css_2docs_2customization.xml
-
-
-// File: md_documentation_2doxygen-awesome-css_2docs_2extensions.xml
-
-
-// File: md_documentation_2doxygen-awesome-css_2docs_2tricks.xml
-
-
-// File: md_documentation_2doxygen-awesome-css_2_r_e_a_d_m_e.xml
-
-
 // File: ems2011.xml
 
 
@@ -18471,15 +18385,9 @@ ymax=1) ";
 // File: dir_b5c4655dde8714643f73aa924cd61a50.xml
 
 
-// File: dir_410aa21f0837b9c4bec0934ed9f3214f.xml
-
-
 // File: dir_138aff360eb965c43b94267b8d1ce09e.xml
 
 
-// File: dir_f54148d295cc294d77e7dd4677dac0ad.xml
-
-
 // File: dir_e99589850f294dbf4b725494ab1c642e.xml
 
 

diff --git a/cmf/cmf_core_src/reach/ReachType.h b/cmf/cmf_core_src/reach/ReachType.h
@@ -376,7 +376,6 @@ namespace cmf {
 			/// - 'R' RectangularReach, depth is ignored
 			/// - 'P' PipeReach, depth is ignored, width is the diameter of the pipe
 			/// - 'S' SWATReachType, a trapezoid flow cross section, as used in the SWAT model, width (bank width) and depth are used
-			/// @returns the reach type
 			/// @param typecode Describes the geometry of the reach cross section.
 			/// @param length The length of the channel in m
 			/// @param width width of the reach cross section in m (ignored for typecode 'T')

diff --git a/cmf/cmf_core_src/water/adsorption.cpp b/cmf/cmf_core_src/water/adsorption.cpp
@@ -1,64 +1,11 @@
 #include "adsorption.h"
-#include "../math/root_finding.h"
 #include <cmath>
 #include <stdexcept>
 using namespace cmf::water;
 LinearAdsorption::LinearAdsorption( real _K,real _m )
 	: K(_K), m(_m)
 {
 
-}
-real cmf::water::FreundlichAdsorbtion::totalsolute( real xf, real V ) const
-{
-	double Ceq = xf / V;
-	double q = K * pow(Ceq, 1 / n);
-	return q * m + xf;
-}
-
-namespace cmf {
-	namespace water {
-		class FreundlichAdsorptionCalculator
-			: public cmf::math::root_finding::BrentsMethod
-		{
-		public:
-			const cmf::water::FreundlichAdsorbtion * owner;
-			FreundlichAdsorptionCalculator(const cmf::water::FreundlichAdsorbtion * fa)
-				: BrentsMethod(owner->epsilon, owner->maxiter), owner(fa)
-			{	}
-			virtual double f(double c) const {
-				return owner->totalsolute(c, 1.0);
-			}
-		};
-	}
-}
-
-
-real FreundlichAdsorbtion::freesolute(real xt, real V ) const
-{
-	//
-	// the Freundlich isotherm x_tot = x_ad + x_free = (m*K*c^n) + (c*V) cannot be rearranged for c if n!=1
-	// hence we have to iterate the solution using Brent's method
-	cmf::water::FreundlichAdsorptionCalculator fac(this);
-	double c_free = fac(0, 1, xt);
-	return c_free * V;
-}
-
-
-FreundlichAdsorbtion::FreundlichAdsorbtion( const FreundlichAdsorbtion& other ) 
-	: K(other.K), n(other.n), m(other.m), epsilon(other.epsilon), maxiter(other.maxiter)
-{}
-
-FreundlichAdsorbtion::FreundlichAdsorbtion( real _K,real _n, real _m, real _epsilon/*=1e-12*/, int _maxiter/*=100*/)
-	: K(_K), n(_n), m(_m), epsilon(_epsilon),maxiter(_maxiter)
-{
-    throw std::runtime_error("Freundlich Adsorption is currently broken and can not be used");
-}
-
-FreundlichAdsorbtion* FreundlichAdsorbtion::copy( real m/*=-1*/ ) const
-{
-	FreundlichAdsorbtion* res = new FreundlichAdsorbtion(*this);
-	if (m>=0) res->m = m;
-	return res;
 }
 real cmf::water::LangmuirAdsorption::totalsolute( real xf,real V ) const
 {

diff --git a/cmf/cmf_core_src/water/adsorption.h b/cmf/cmf_core_src/water/adsorption.h
@@ -98,53 +98,6 @@ namespace cmf {
 
 			}
 		};
-		/// @brief BROKEN: This class calculates the adsorption equilibrium between sorbat and sorbent using the Freundlich isotherme.
-		///
-		/// Freundlich isotherme:
-		///
-		/// \f[\frac{x_{ad}}{m} = K c^n\f]
-		/// where
-		///  - \f$x_{ad} = x_{tot} - x_{free}\f$ is the adsorbed tracer mass
-		///     - \f$x_{tot}\f$ is the total tracer mass
-		///     - \f$x_{free}\f$ is the dissolved tracer mass
-		///  - \f$m\f$ is the mass of the sorbent in the same unit as the tracer mass
-		///  - \f$K\f$ is the Freundlich sorption coefficient
-		///  - \f$c = \frac{x_{free}}{V}\f$ is the concentration of the tracer in tracer mass per m3
-		///  - \f$n\f$ is the Freundlich exponent
-		///
-		/// CMF stores in a solute storage the total mass of a tracer and needs to calculate the free tracer mass.
-		/// The eq. above can not be rearanged to get \f$x_{free}\f$ from \f$x_{tot}\f$. Instead, the value is iterated
-		/// using [regula falsi](http://en.wikipedia.org/wiki/False_position_method). If n is near to 1, 
-		/// using LinearAdsorption will speed up your calculations.
-		///
-		/// The simplest physically based adsorption model by Langmuir (LangmuirAdsorption) has also a analytical solution
-		/// and is hence calculated faster then Freundlich.
-		class FreundlichAdsorbtion: public Adsorption {
-		public:
-			/// Freundlich half saturation
-			real K;
-			/// Freundlich n
-			real n;
-			/// Mass of sorbent (CEC, clay mass etc.)
-			real m;
-			/// Tolerable error for Newton iteration
-			real epsilon;
-			/// Maximum number of iterations
-			int maxiter;
-			real freesolute(real xt,real V) const;
-			real totalsolute(real xf, real V) const;
-			/// @param K,n Freundlich coefficents
-			/// @param m Mass of sorbent in units of tracer
-			/// @param epsilon Tolerance of regula falsi iteration for the calculation of dissolved tracer from total trace, default = 1e-12
-			/// @param maxiter Maximum number of iterations, default = 100
-			FreundlichAdsorbtion(real K,real n, real m, real epsilon=1e-12, int maxiter=100);
-			FreundlichAdsorbtion(const FreundlichAdsorbtion& other);
-			FreundlichAdsorbtion* copy(real m=-1) const;
-			virtual ~FreundlichAdsorbtion()
-			{
-
-			}
-		};
 		/// @brief This class calculates the adsorption equilibrium between sorbat and sorbent using the Langmuir isotherme.
 		///
 		/// Langmuir Adsorption:

diff --git a/cmf/cmf_core_src/water/flux_connection.h b/cmf/cmf_core_src/water/flux_connection.h
@@ -37,8 +37,7 @@ namespace cmf {
 		/// @brief The connections between the nodes (boundary conditions, storages) of the water network
         ///
 		/// The connections in cmf hold the processes for the calculation of fluxes between water storages and model boundaries
-		/// @todo Elaborate on this
-
+		///
 		/// @ingroup connections
 		/// Represents a connection between flux_nodes, where water fluxes occur. 
 		class flux_connection 

diff --git a/cmf/cmf_core_src/water/freundlich_adsorption.cpp b/cmf/cmf_core_src/water/freundlich_adsorption.cpp
@@ -0,0 +1,59 @@
+#include "freundlich_adsorption.h"
+#include "../math/root_finding.h"
+#include <cmath>
+#include <stdexcept>
+
+
+real cmf::water::FreundlichAdsorbtion::totalsolute( real xf, real V ) const
+{
+	double Ceq = xf / V;
+	double q = K * pow(Ceq, 1 / n);
+	return q * m + xf;
+}
+
+namespace cmf {
+	namespace water {
+		class FreundlichAdsorptionCalculator
+			: public cmf::math::root_finding::Bisect
+		{
+		public:
+			const cmf::water::FreundlichAdsorbtion * owner;
+			FreundlichAdsorptionCalculator(const cmf::water::FreundlichAdsorbtion * fa)
+				: Bisect(owner->epsilon, owner->maxiter), owner(fa)
+			{	}
+			virtual double f(double c) const {
+				return owner->totalsolute(c, 1.0);
+			}
+		};
+	}
+}
+
+
+real cmf::water::FreundlichAdsorbtion::freesolute(real xt, real V ) const
+{
+	//
+	// the Freundlich isotherm x_tot = x_ad + x_free = (m*K*c^n) + (c*V) cannot be rearranged for c if n!=1
+	// hence we have to iterate the solution using Brent's method
+	cmf::water::FreundlichAdsorptionCalculator fac(this);
+
+	double c_free = fac(0, xt, xt);
+	return c_free * V;
+}
+
+
+cmf::water::FreundlichAdsorbtion::FreundlichAdsorbtion( const cmf::water::FreundlichAdsorbtion& other )
+	: K(other.K), n(other.n), m(other.m), epsilon(other.epsilon), maxiter(other.maxiter)
+{}
+
+cmf::water::FreundlichAdsorbtion::FreundlichAdsorbtion( real _K,real _n, real _m, real _epsilon/*=1e-12*/, int _maxiter/*=100*/)
+	: K(_K), n(_n), m(_m), epsilon(_epsilon),maxiter(_maxiter)
+{
+    // throw std::runtime_error("Freundlich Adsorption is currently broken and can not be used");
+}
+
+cmf::water::FreundlichAdsorbtion* cmf::water::FreundlichAdsorbtion::copy( real m/*=-1*/ ) const
+{
+	cmf::water::FreundlichAdsorbtion* res = new cmf::water::FreundlichAdsorbtion(*this);
+	if (m>=0) res->m = m;
+	return res;
+}