From 94b6c2988a32da57246b2adc4b264bec27573fb3 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Sebastian=20M=C3=BCller?= <mueller.seb@posteo.de>
Date: Fri, 26 Jul 2024 23:11:40 +0200
Subject: [PATCH 01/58] plot: handle len_scale of 0

---
 src/gstools/covmodel/plot.py | 75 ++++++++++++++++++++++++++++--------
 1 file changed, 60 insertions(+), 15 deletions(-)

diff --git a/src/gstools/covmodel/plot.py b/src/gstools/covmodel/plot.py
index 32148c14a..f7fa58db5 100644
--- a/src/gstools/covmodel/plot.py
+++ b/src/gstools/covmodel/plot.py
@@ -68,7 +68,10 @@ def plot_vario_spatial(
 ):  # pragma: no cover
     """Plot spatial variogram of a given CovModel."""
     if x_max is None:
-        x_max = 3 * model.len_scale
+        if np.isclose(model.len_scale, 0):
+            x_max = 1.0
+        else:
+            x_max = 3 * model.len_scale
     x_s = np.linspace(-x_max, x_max) + x_min
     pos = [x_s] * model.dim
     shp = tuple(len(p) for p in pos)
@@ -83,7 +86,10 @@ def plot_cov_spatial(
 ):  # pragma: no cover
     """Plot spatial covariance of a given CovModel."""
     if x_max is None:
-        x_max = 3 * model.len_scale
+        if np.isclose(model.len_scale, 0):
+            x_max = 1.0
+        else:
+            x_max = 3 * model.len_scale
     x_s = np.linspace(-x_max, x_max) + x_min
     pos = [x_s] * model.dim
     shp = tuple(len(p) for p in pos)
@@ -98,7 +104,10 @@ def plot_cor_spatial(
 ):  # pragma: no cover
     """Plot spatial correlation of a given CovModel."""
     if x_max is None:
-        x_max = 3 * model.len_scale
+        if np.isclose(model.len_scale, 0):
+            x_max = 1.0
+        else:
+            x_max = 3 * model.len_scale
     x_s = np.linspace(-x_max, x_max) + x_min
     pos = [x_s] * model.dim
     shp = tuple(len(p) for p in pos)
@@ -114,7 +123,10 @@ def plot_variogram(
     """Plot variogram of a given CovModel."""
     fig, ax = get_fig_ax(fig, ax)
     if x_max is None:
-        x_max = 3 * model.len_scale
+        if np.isclose(model.len_scale, 0):
+            x_max = 1.0
+        else:
+            x_max = 3 * model.len_scale
     x_s = np.linspace(x_min, x_max)
     kwargs.setdefault("label", f"{model.name} variogram")
     ax.plot(x_s, model.variogram(x_s), **kwargs)
@@ -129,7 +141,10 @@ def plot_covariance(
     """Plot covariance of a given CovModel."""
     fig, ax = get_fig_ax(fig, ax)
     if x_max is None:
-        x_max = 3 * model.len_scale
+        if np.isclose(model.len_scale, 0):
+            x_max = 1.0
+        else:
+            x_max = 3 * model.len_scale
     x_s = np.linspace(x_min, x_max)
     kwargs.setdefault("label", f"{model.name} covariance")
     ax.plot(x_s, model.covariance(x_s), **kwargs)
@@ -144,7 +159,10 @@ def plot_correlation(
     """Plot correlation function of a given CovModel."""
     fig, ax = get_fig_ax(fig, ax)
     if x_max is None:
-        x_max = 3 * model.len_scale
+        if np.isclose(model.len_scale, 0):
+            x_max = 1.0
+        else:
+            x_max = 3 * model.len_scale
     x_s = np.linspace(x_min, x_max)
     kwargs.setdefault("label", f"{model.name} correlation")
     ax.plot(x_s, model.correlation(x_s), **kwargs)
@@ -159,7 +177,10 @@ def plot_vario_yadrenko(
     """Plot Yadrenko variogram of a given CovModel."""
     fig, ax = get_fig_ax(fig, ax)
     if x_max is None:
-        x_max = min(3 * model.len_scale, model.geo_scale * np.pi)
+        if np.isclose(model.len_scale, 0):
+            x_max = 1.0
+        else:
+            x_max = min(3 * model.len_scale, model.geo_scale * np.pi)
     x_s = np.linspace(x_min, x_max)
     kwargs.setdefault("label", f"{model.name} Yadrenko variogram")
     ax.plot(x_s, model.vario_yadrenko(x_s), **kwargs)
@@ -174,7 +195,10 @@ def plot_cov_yadrenko(
     """Plot Yadrenko covariance of a given CovModel."""
     fig, ax = get_fig_ax(fig, ax)
     if x_max is None:
-        x_max = min(3 * model.len_scale, model.geo_scale * np.pi)
+        if np.isclose(model.len_scale, 0):
+            x_max = 1.0
+        else:
+            x_max = min(3 * model.len_scale, model.geo_scale * np.pi)
     x_s = np.linspace(x_min, x_max)
     kwargs.setdefault("label", f"{model.name} Yadrenko covariance")
     ax.plot(x_s, model.cov_yadrenko(x_s), **kwargs)
@@ -189,7 +213,10 @@ def plot_cor_yadrenko(
     """Plot Yadrenko correlation function of a given CovModel."""
     fig, ax = get_fig_ax(fig, ax)
     if x_max is None:
-        x_max = min(3 * model.len_scale, model.geo_scale * np.pi)
+        if np.isclose(model.len_scale, 0):
+            x_max = 1.0
+        else:
+            x_max = min(3 * model.len_scale, model.geo_scale * np.pi)
     x_s = np.linspace(x_min, x_max)
     kwargs.setdefault("label", f"{model.name} Yadrenko correlation")
     ax.plot(x_s, model.cor_yadrenko(x_s), **kwargs)
@@ -204,7 +231,10 @@ def plot_vario_axis(
     """Plot variogram of a given CovModel."""
     fig, ax = get_fig_ax(fig, ax)
     if x_max is None:
-        x_max = 3 * model.len_scale
+        if np.isclose(model.len_scale, 0):
+            x_max = 1.0
+        else:
+            x_max = 3 * model.len_scale
     x_s = np.linspace(x_min, x_max)
     kwargs.setdefault("label", f"{model.name} variogram on axis {axis}")
     ax.plot(x_s, model.vario_axis(x_s, axis), **kwargs)
@@ -219,7 +249,10 @@ def plot_cov_axis(
     """Plot variogram of a given CovModel."""
     fig, ax = get_fig_ax(fig, ax)
     if x_max is None:
-        x_max = 3 * model.len_scale
+        if np.isclose(model.len_scale, 0):
+            x_max = 1.0
+        else:
+            x_max = 3 * model.len_scale
     x_s = np.linspace(x_min, x_max)
     kwargs.setdefault("label", f"{model.name} covariance on axis {axis}")
     ax.plot(x_s, model.cov_axis(x_s, axis), **kwargs)
@@ -234,7 +267,10 @@ def plot_cor_axis(
     """Plot variogram of a given CovModel."""
     fig, ax = get_fig_ax(fig, ax)
     if x_max is None:
-        x_max = 3 * model.len_scale
+        if np.isclose(model.len_scale, 0):
+            x_max = 1.0
+        else:
+            x_max = 3 * model.len_scale
     x_s = np.linspace(x_min, x_max)
     kwargs.setdefault("label", f"{model.name} correlation on axis {axis}")
     ax.plot(x_s, model.cor_axis(x_s, axis), **kwargs)
@@ -249,7 +285,10 @@ def plot_spectrum(
     """Plot spectrum of a given CovModel."""
     fig, ax = get_fig_ax(fig, ax)
     if x_max is None:
-        x_max = 3 / model.len_scale
+        if np.isclose(model.len_scale, 0):
+            x_max = 1.0
+        else:
+            x_max = 3 / model.len_scale
     x_s = np.linspace(x_min, x_max)
     kwargs.setdefault("label", f"{model.name} {model.dim}D spectrum")
     ax.plot(x_s, model.spectrum(x_s), **kwargs)
@@ -264,7 +303,10 @@ def plot_spectral_density(
     """Plot spectral density of a given CovModel."""
     fig, ax = get_fig_ax(fig, ax)
     if x_max is None:
-        x_max = 3 / model.len_scale
+        if np.isclose(model.len_scale, 0):
+            x_max = 1.0
+        else:
+            x_max = 3 / model.len_scale
     x_s = np.linspace(x_min, x_max)
     kwargs.setdefault("label", f"{model.name} {model.dim}D spectral-density")
     ax.plot(x_s, model.spectral_density(x_s), **kwargs)
@@ -279,7 +321,10 @@ def plot_spectral_rad_pdf(
     """Plot radial spectral pdf of a given CovModel."""
     fig, ax = get_fig_ax(fig, ax)
     if x_max is None:
-        x_max = 3 / model.len_scale
+        if np.isclose(model.len_scale, 0):
+            x_max = 1.0
+        else:
+            x_max = 3 / model.len_scale
     x_s = np.linspace(x_min, x_max)
     kwargs.setdefault("label", f"{model.name} {model.dim}D spectral-rad-pdf")
     ax.plot(x_s, model.spectral_rad_pdf(x_s), **kwargs)

From abbb062cffff08391bf64e6c3ed96ffee5996a44 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Sebastian=20M=C3=BCller?= <mueller.seb@posteo.de>
Date: Fri, 26 Jul 2024 23:12:34 +0200
Subject: [PATCH 02/58] covmodel: init len_scale with integral_scale if given;
 remove unused _integral_scale

---
 src/gstools/covmodel/base.py | 11 +++++------
 1 file changed, 5 insertions(+), 6 deletions(-)

diff --git a/src/gstools/covmodel/base.py b/src/gstools/covmodel/base.py
index 23e198812..f5def43ca 100644
--- a/src/gstools/covmodel/base.py
+++ b/src/gstools/covmodel/base.py
@@ -208,6 +208,8 @@ def __init__(
         self._nugget = float(nugget)
 
         # set anisotropy and len_scale, disable anisotropy for latlon models
+        if integral_scale is not None:
+            len_scale = integral_scale
         self._len_scale, self._anis = set_len_anis(
             self.dim, len_scale, anis, self.latlon
         )
@@ -221,7 +223,6 @@ def __init__(
             self.var = var
         else:
             self._var = float(var_raw)
-        self._integral_scale = None
         self.integral_scale = integral_scale
         # set var again, if int_scale affects var_factor
         if var_raw is None:
@@ -486,8 +487,7 @@ def default_rescale(self):
 
     def calc_integral_scale(self):
         """Calculate the integral scale of the isotrope model."""
-        self._integral_scale = integral(self.correlation, 0, np.inf)[0]
-        return self._integral_scale
+        return integral(self.correlation, 0, np.inf)[0]
 
     def percentile_scale(self, per=0.9):
         """Calculate the percentile scale of the isotrope model.
@@ -1043,8 +1043,7 @@ def integral_scale(self):
         ValueError
             If integral scale is not setable.
         """
-        self._integral_scale = self.calc_integral_scale()
-        return self._integral_scale
+        return self.calc_integral_scale()
 
     @integral_scale.setter
     def integral_scale(self, integral_scale):
@@ -1054,7 +1053,7 @@ def integral_scale(self, integral_scale):
             integral_scale = self.len_scale
             # reset len_scale
             self.len_scale = 1.0
-            int_tmp = self.calc_integral_scale()
+            int_tmp = self.integral_scale
             self.len_scale = integral_scale / int_tmp
             if not np.isclose(self.integral_scale, integral_scale, rtol=1e-3):
                 raise ValueError(

From c39866a90096a71db11092b746ca979e78c991a4 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Sebastian=20M=C3=BCller?= <mueller.seb@posteo.de>
Date: Fri, 26 Jul 2024 23:13:14 +0200
Subject: [PATCH 03/58] bounds: allow interval with equal bounds if closed

---
 src/gstools/covmodel/tools.py | 7 +++++--
 1 file changed, 5 insertions(+), 2 deletions(-)

diff --git a/src/gstools/covmodel/tools.py b/src/gstools/covmodel/tools.py
index dddeb4413..cf310b94f 100644
--- a/src/gstools/covmodel/tools.py
+++ b/src/gstools/covmodel/tools.py
@@ -297,11 +297,14 @@ def check_bounds(bounds):
                 * "oo" : open - open
                 * "oc" : open - close
                 * "co" : close - open
-                * "cc" : close - close
+                * "cc" : close - close (default)
     """
+    typ = bounds[2] if len(bounds) == 3 else "cc"
     if len(bounds) not in (2, 3):
         return False
-    if bounds[1] <= bounds[0]:
+    if (typ == "cc" and bounds[1] < bounds[0]) or (
+        typ != "cc" and bounds[1] <= bounds[0]
+    ):
         return False
     if len(bounds) == 3 and bounds[2] not in ("oo", "oc", "co", "cc"):
         return False

From df9f2baa949201a6f54f41a2473f074909365f38 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Sebastian=20M=C3=BCller?= <mueller.seb@posteo.de>
Date: Sat, 27 Jul 2024 01:04:30 +0200
Subject: [PATCH 04/58] CovModel: add force_values logic

---
 src/gstools/covmodel/base.py  | 65 ++++++++++++++++++++++++-----------
 src/gstools/covmodel/fit.py   |  4 +++
 src/gstools/covmodel/tools.py | 10 +++---
 3 files changed, 55 insertions(+), 24 deletions(-)

diff --git a/src/gstools/covmodel/base.py b/src/gstools/covmodel/base.py
index f5def43ca..3a61e2294 100644
--- a/src/gstools/covmodel/base.py
+++ b/src/gstools/covmodel/base.py
@@ -139,6 +139,8 @@ class CovModel:
         used for the spectrum calculation. Use with caution (Better: Don't!).
         ``None`` is equivalent to ``{"a": -1, "b": 1, "N": 1000, "h": 0.001}``.
         Default: :any:`None`
+    fixed: :class:`set` or :any:`None`, optional
+        Names of fixed arguments. Default: :any:`None`
     **opt_arg
         Optional arguments are covered by these keyword arguments.
         If present, they are described in the section `Other Parameters`.
@@ -161,6 +163,7 @@ def __init__(
         spatial_dim=None,
         var_raw=None,
         hankel_kw=None,
+        fixed=None,
         **opt_arg,
     ):
         # assert, that we use a subclass
@@ -169,6 +172,16 @@ def __init__(
         if not hasattr(self, "variogram"):
             raise TypeError("Don't instantiate 'CovModel' directly!")
 
+        # indicate that arguments are fixed (True after __init__)
+        self._fix = False
+        # force input values
+        forced = self.force_values()
+        var = forced.get("var", var)
+        len_scale = forced.get("len_scale", len_scale)
+        nugget = forced.get("nugget", nugget)
+        integral_scale = forced.get("integral_scale", integral_scale)
+        for opt in opt_arg:
+            opt_arg[opt] = forced.get(opt, opt_arg[opt])
         # prepare dim setting
         self._dim = None
         self._hankel_kw = None
@@ -198,6 +211,16 @@ def __init__(
         # optional arguments for the variogram-model
         set_opt_args(self, opt_arg)
 
+        # check fixed
+        fixed = set(fixed) if fixed is not None else set()
+        fixed = fixed | self.always_fixed()
+        valid_fixed = set(self.iso_arg)  # | set(["anis"])
+        if not fixed <= valid_fixed:
+            raise ValueError(
+                f"CovModel: unknown names in 'fixed': {fixed - valid_fixed}"
+            )
+        self._fixed = fixed
+
         # set standard boundaries for variance, len_scale, nugget and opt_arg
         bounds = self.default_arg_bounds()
         bounds.update(self.default_opt_arg_bounds())
@@ -234,6 +257,11 @@ def __init__(
         self.check_arg_bounds()
         # additional checks for the optional arguments (provided by user)
         self.check_opt_arg()
+        # set fixed bounds after checking original bounds
+        for arg in self.fixed:
+            val = getattr(self, arg)
+            self.set_arg_bounds(check_args=False, **{arg: (val, val, "cc")})
+        self._fix = True
         # precision for printing
         self._prec = 3
 
@@ -438,6 +466,14 @@ def pykrige_kwargs(self):
 
     # methods for optional/default arguments (can be overridden)
 
+    def always_fixed(self):
+        """Provide set of fixed arguments."""
+        return set()
+
+    def force_values(self):
+        """:class:`dict`: Forced values for arguments."""
+        return {}
+
     def default_opt_arg(self):
         """Provide default optional arguments by the user.
 
@@ -795,6 +831,11 @@ def check_arg_bounds(self):
 
     # bounds properties
 
+    @property
+    def fixed(self):
+        """:class:`set`: Set with names of fixed arguments."""
+        return self._fixed
+
     @property
     def var_bounds(self):
         """:class:`list`: Bounds for the variance.
@@ -809,11 +850,7 @@ def var_bounds(self):
 
     @var_bounds.setter
     def var_bounds(self, bounds):
-        if not check_bounds(bounds):
-            raise ValueError(
-                f"Given bounds for 'var' are not valid, got: {bounds}"
-            )
-        self._var_bounds = bounds
+        self.set_arg_bounds(var=bounds)
 
     @property
     def len_scale_bounds(self):
@@ -829,11 +866,7 @@ def len_scale_bounds(self):
 
     @len_scale_bounds.setter
     def len_scale_bounds(self, bounds):
-        if not check_bounds(bounds):
-            raise ValueError(
-                f"Given bounds for 'len_scale' are not valid, got: {bounds}"
-            )
-        self._len_scale_bounds = bounds
+        self.set_arg_bounds(len_scale=bounds)
 
     @property
     def nugget_bounds(self):
@@ -849,11 +882,7 @@ def nugget_bounds(self):
 
     @nugget_bounds.setter
     def nugget_bounds(self, bounds):
-        if not check_bounds(bounds):
-            raise ValueError(
-                f"Given bounds for 'nugget' are not valid, got: {bounds}"
-            )
-        self._nugget_bounds = bounds
+        self.set_arg_bounds(nugget=bounds)
 
     @property
     def anis_bounds(self):
@@ -869,11 +898,7 @@ def anis_bounds(self):
 
     @anis_bounds.setter
     def anis_bounds(self, bounds):
-        if not check_bounds(bounds):
-            raise ValueError(
-                f"Given bounds for 'anis' are not valid, got: {bounds}"
-            )
-        self._anis_bounds = bounds
+        self.set_arg_bounds(anis=bounds)
 
     @property
     def opt_arg_bounds(self):
diff --git a/src/gstools/covmodel/fit.py b/src/gstools/covmodel/fit.py
index 8b19f4977..f19eaac35 100755
--- a/src/gstools/covmodel/fit.py
+++ b/src/gstools/covmodel/fit.py
@@ -165,10 +165,14 @@ def fit_variogram(
 
     The fitted parameters will be instantly set in the model.
     """
+    for para in model.fixed:
+        para_select[para] = False
     # preprocess selected parameters
     para, sill, constrain_sill, anis = _pre_para(
         model, para_select, sill, anis
     )
+    if not any(para.values()):
+        raise ValueError("fit: no parameters selected for fitting.")
     # check curve_fit kwargs
     curve_fit_kwargs = {} if curve_fit_kwargs is None else curve_fit_kwargs
     # check method
diff --git a/src/gstools/covmodel/tools.py b/src/gstools/covmodel/tools.py
index cf310b94f..ad731aeb2 100644
--- a/src/gstools/covmodel/tools.py
+++ b/src/gstools/covmodel/tools.py
@@ -458,6 +458,8 @@ def set_arg_bounds(model, check_args=True, **kwargs):
     # if variance needs to be resetted, do this at last
     var_bnds = []
     for arg, bounds in kwargs.items():
+        if model._fix and arg in model.fixed:
+            raise ValueError(f"Can't set bounds for fixed argument '{arg}'")
         if not check_bounds(bounds):
             raise ValueError(
                 f"Given bounds for '{arg}' are not valid, got: {bounds}"
@@ -468,11 +470,11 @@ def set_arg_bounds(model, check_args=True, **kwargs):
             var_bnds = bounds
             continue
         elif arg == "len_scale":
-            model.len_scale_bounds = bounds
+            model._len_scale_bounds = bounds
         elif arg == "nugget":
-            model.nugget_bounds = bounds
+            model._nugget_bounds = bounds
         elif arg == "anis":
-            model.anis_bounds = bounds
+            model._anis_bounds = bounds
         else:
             raise ValueError(f"set_arg_bounds: unknown argument '{arg}'")
         if check_args and check_arg_in_bounds(model, arg) > 0:
@@ -483,7 +485,7 @@ def set_arg_bounds(model, check_args=True, **kwargs):
                 setattr(model, arg, def_arg)
     # set var last like always
     if var_bnds:
-        model.var_bounds = var_bnds
+        model._var_bounds = var_bnds
         if check_args and check_arg_in_bounds(model, "var") > 0:
             model.var = default_arg_from_bounds(var_bnds)
 

From a3f1ec6d897ba161f93abd133db7480ee70fa5d5 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Sebastian=20M=C3=BCller?= <mueller.seb@posteo.de>
Date: Sat, 27 Jul 2024 23:05:10 +0200
Subject: [PATCH 05/58] generator: shortcut for 0 variance

---
 src/gstools/field/generator.py | 75 ++++++++++++++++++++++------------
 1 file changed, 50 insertions(+), 25 deletions(-)

diff --git a/src/gstools/field/generator.py b/src/gstools/field/generator.py
index d4cb0dea5..89914785a 100755
--- a/src/gstools/field/generator.py
+++ b/src/gstools/field/generator.py
@@ -152,6 +152,16 @@ def __call__(self, pos, add_nugget=True):
             the random modes
         """
 
+    @property
+    @abstractmethod
+    def model(self):
+        """:any:`CovModel`: Covariance model of the spatial random field."""
+
+    @property
+    def zero_var(self):
+        """:class:`bool`: Whether Covariance model has zero variance."""
+        return np.isclose(self.model.var, 0)
+
     @property
     @abstractmethod
     def value_type(self):
@@ -260,6 +270,10 @@ def __call__(self, pos, add_nugget=True):
             the random modes
         """
         pos = np.asarray(pos, dtype=np.double)
+        if self.zero_var:
+            shp = pos.shape[1:]
+            return self.get_nugget(shp) if add_nugget else np.full(shp, 0.0)
+        # generate if var is not 0
         summed_modes = _summate(
             self._cov_sample, self._z_1, self._z_2, pos, config.NUM_THREADS
         )
@@ -285,7 +299,7 @@ def get_nugget(self, shape):
                 size=shape
             )
         else:
-            nugget = 0.0
+            nugget = np.full(shape, 0.0) if self.zero_var else 0.0
         return nugget
 
     def update(self, model=None, seed=np.nan):
@@ -362,23 +376,26 @@ def reset_seed(self, seed=np.nan):
         self._z_1 = self._rng.random.normal(size=self._mode_no)
         self._z_2 = self._rng.random.normal(size=self._mode_no)
         # sample uniform on a sphere
-        sphere_coord = self._rng.sample_sphere(self.model.dim, self._mode_no)
-        # sample radii according to radial spectral density of the model
-        if self.sampling == "inversion" or (
-            self.sampling == "auto" and self.model.has_ppf
-        ):
-            pdf, cdf, ppf = self.model.dist_func
-            rad = self._rng.sample_dist(
-                size=self._mode_no, pdf=pdf, cdf=cdf, ppf=ppf, a=0
-            )
+        if self.zero_var:
+            self._cov_sample = np.full((self.model.dim, self._mode_no), 0.0)
         else:
-            rad = self._rng.sample_ln_pdf(
-                ln_pdf=self.model.ln_spectral_rad_pdf,
-                size=self._mode_no,
-                sample_around=1.0 / self.model.len_rescaled,
-            )
-        # get fully spatial samples by multiplying sphere samples and radii
-        self._cov_sample = rad * sphere_coord
+            sph_crd = self._rng.sample_sphere(self.model.dim, self._mode_no)
+            # sample radii according to radial spectral density of the model
+            if self.sampling == "inversion" or (
+                self.sampling == "auto" and self.model.has_ppf
+            ):
+                pdf, cdf, ppf = self.model.dist_func
+                rad = self._rng.sample_dist(
+                    size=self._mode_no, pdf=pdf, cdf=cdf, ppf=ppf, a=0
+                )
+            else:
+                rad = self._rng.sample_ln_pdf(
+                    ln_pdf=self.model.ln_spectral_rad_pdf,
+                    size=self._mode_no,
+                    sample_around=1.0 / self.model.len_rescaled,
+                )
+            # get spatial samples by multiplying sphere samples and radii
+            self._cov_sample = rad * sph_crd
 
     @property
     def sampling(self):
@@ -541,6 +558,10 @@ def __call__(self, pos, add_nugget=True):
             the random modes
         """
         pos = np.asarray(pos, dtype=np.double)
+        nugget = self.get_nugget(pos.shape) if add_nugget else 0.0
+        e1 = self._create_unit_vector(pos.shape)
+        if self.zero_var:
+            return self.mean_u * e1 + nugget
         summed_modes = _summate_incompr(
             self._cov_sample,
             self._z_1,
@@ -548,8 +569,6 @@ def __call__(self, pos, add_nugget=True):
             pos,
             config.NUM_THREADS,
         )
-        nugget = self.get_nugget(summed_modes.shape) if add_nugget else 0.0
-        e1 = self._create_unit_vector(summed_modes.shape)
         return (
             self.mean_u * e1
             + self.mean_u
@@ -667,7 +686,10 @@ def __call__(self, pos, add_nugget=True):
             the random modes
         """
         pos = np.asarray(pos, dtype=np.double)
-
+        if self.zero_var:
+            shp = pos.shape[1:]
+            return self.get_nugget(shp) if add_nugget else np.full(shp, 0.0)
+        # generate if var is not 0
         summed_modes = _summate_fourier(
             self._spectrum_factor,
             self._modes,
@@ -698,7 +720,7 @@ def get_nugget(self, shape):
                 size=shape
             )
         else:
-            nugget = 0.0
+            nugget = np.full(shape, 0.0) if self.zero_var else 0.0
         return nugget
 
     def update(self, model=None, seed=np.nan, period=None, mode_no=None):
@@ -800,10 +822,13 @@ def reset_seed(self, seed=np.nan):
         self._z_2 = self._rng.random.normal(size=np.prod(self._mode_no))
         # pre calc. the spectrum for all wave numbers they are handed over to
         # Cython, which doesn't have access to the CovModel
-        k_norm = np.linalg.norm(self._modes, axis=0)
-        self._spectrum_factor = np.sqrt(
-            self._model.spectrum(k_norm) * np.prod(self._delta_k)
-        )
+        if self.zero_var:
+            self._spectrum_factor = np.full(np.prod(self._mode_no), 0.0)
+        else:
+            k_norm = np.linalg.norm(self._modes, axis=0)
+            self._spectrum_factor = np.sqrt(
+                self._model.spectrum(k_norm) * np.prod(self._delta_k)
+            )
 
     def _fill_to_dim(
         self, values, dim, dtype=float, default_value=None

From 8b3ad819e3cf4431fa66f37bd11ace80d850db52 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Sebastian=20M=C3=BCller?= <mueller.seb@posteo.de>
Date: Mon, 29 Jul 2024 00:42:25 +0200
Subject: [PATCH 06/58] CovModel: safer usage of privat attr; anis set bug fix;
 simpler int scale

---
 src/gstools/covmodel/base.py | 51 +++++++++++++++++-------------------
 1 file changed, 24 insertions(+), 27 deletions(-)

diff --git a/src/gstools/covmodel/base.py b/src/gstools/covmodel/base.py
index 3a61e2294..dcbffe57e 100644
--- a/src/gstools/covmodel/base.py
+++ b/src/gstools/covmodel/base.py
@@ -173,7 +173,7 @@ def __init__(
             raise TypeError("Don't instantiate 'CovModel' directly!")
 
         # indicate that arguments are fixed (True after __init__)
-        self._fix = False
+        self._init = False
         # force input values
         forced = self.force_values()
         var = forced.get("var", var)
@@ -187,6 +187,7 @@ def __init__(
         self._hankel_kw = None
         self._sft = None
         # prepare parameters (they are checked in dim setting)
+        self._fixed = None
         self._rescale = None
         self._len_scale = None
         self._anis = None
@@ -246,7 +247,8 @@ def __init__(
             self.var = var
         else:
             self._var = float(var_raw)
-        self.integral_scale = integral_scale
+        if integral_scale is not None:
+            self.integral_scale = integral_scale
         # set var again, if int_scale affects var_factor
         if var_raw is None:
             self._var = None
@@ -261,7 +263,7 @@ def __init__(
         for arg in self.fixed:
             val = getattr(self, arg)
             self.set_arg_bounds(check_args=False, **{arg: (val, val, "cc")})
-        self._fix = True
+        self._init = True
         # precision for printing
         self._prec = 3
 
@@ -801,8 +803,8 @@ def default_arg_bounds(self):
         Given as a dictionary.
         """
         res = {
-            "var": (0.0, np.inf, "oo"),
-            "len_scale": (0.0, np.inf, "oo"),
+            "var": (0.0, np.inf, "co"),
+            "len_scale": (0.0, np.inf, "co"),
             "nugget": (0.0, np.inf, "co"),
             "anis": (0.0, np.inf, "oo"),
         }
@@ -1014,10 +1016,7 @@ def len_scale(self, len_scale):
         self._len_scale, anis = set_len_anis(
             self.dim, len_scale, self.anis, self.latlon
         )
-        if self.latlon:
-            self._anis = np.array((self.dim - 1) * [1], dtype=np.double)
-        else:
-            self._anis = anis
+        self._anis = anis
         self.check_arg_bounds()
 
     @property
@@ -1033,7 +1032,7 @@ def rescale(self, rescale):
     @property
     def len_rescaled(self):
         """:class:`float`: The rescaled main length scale of the model."""
-        return self._len_scale / self._rescale
+        return self.len_scale / self.rescale
 
     @property
     def anis(self):
@@ -1072,19 +1071,17 @@ def integral_scale(self):
 
     @integral_scale.setter
     def integral_scale(self, integral_scale):
-        if integral_scale is not None:
-            # format int-scale right
-            self.len_scale = integral_scale
-            integral_scale = self.len_scale
-            # reset len_scale
-            self.len_scale = 1.0
-            int_tmp = self.integral_scale
-            self.len_scale = integral_scale / int_tmp
-            if not np.isclose(self.integral_scale, integral_scale, rtol=1e-3):
-                raise ValueError(
-                    f"{self.name}: Integral scale could not be set correctly! "
-                    "Please just provide a 'len_scale'!"
-                )
+        int_scale, anis = set_len_anis(
+            self.dim, integral_scale, self.anis, self.latlon
+        )
+        old_scale = self.integral_scale
+        self.anis = anis
+        self.len_scale = self.len_scale * int_scale / old_scale
+        if not np.isclose(self.integral_scale, integral_scale, rtol=1e-3):
+            raise ValueError(
+                f"{self.name}: Integral scale could not be set correctly! "
+                "Please just provide a 'len_scale'!"
+            )
 
     @property
     def hankel_kw(self):
@@ -1139,7 +1136,7 @@ def sill(self):
     @property
     def arg(self):
         """:class:`list` of :class:`str`: Names of all arguments."""
-        return ["var", "len_scale", "nugget", "anis", "angles"] + self._opt_arg
+        return ["var", "len_scale", "nugget", "anis", "angles"] + self.opt_arg
 
     @property
     def arg_list(self):
@@ -1152,7 +1149,7 @@ def arg_list(self):
     @property
     def iso_arg(self):
         """:class:`list` of :class:`str`: Names of isotropic arguments."""
-        return ["var", "len_scale", "nugget"] + self._opt_arg
+        return ["var", "len_scale", "nugget"] + self.opt_arg
 
     @property
     def iso_arg_list(self):
@@ -1180,7 +1177,7 @@ def len_scale_vec(self):
         """
         res = np.zeros(self.dim, dtype=np.double)
         res[0] = self.len_scale
-        for i in range(1, self._dim):
+        for i in range(1, self.dim):
             res[i] = self.len_scale * self.anis[i - 1]
         return res
 
@@ -1226,7 +1223,7 @@ def __setattr__(self, name, value):
         """Set an attribute."""
         super().__setattr__(name, value)
         # if an optional variogram argument was given, check bounds
-        if hasattr(self, "_opt_arg") and name in self._opt_arg:
+        if getattr(self, "_init", False) and name in self.opt_arg:
             self.check_arg_bounds()
 
     def __repr__(self):

From a60b3d7e9e6066f986ea4baea60e0f78c7e1feb5 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Sebastian=20M=C3=BCller?= <mueller.seb@posteo.de>
Date: Mon, 29 Jul 2024 00:43:12 +0200
Subject: [PATCH 07/58] CovModel: check _init in bounds setter; also compare
 geo_scale; simpler repr

---
 src/gstools/covmodel/tools.py | 50 +++++++++++++----------------------
 1 file changed, 18 insertions(+), 32 deletions(-)

diff --git a/src/gstools/covmodel/tools.py b/src/gstools/covmodel/tools.py
index ad731aeb2..9d53de583 100644
--- a/src/gstools/covmodel/tools.py
+++ b/src/gstools/covmodel/tools.py
@@ -458,7 +458,7 @@ def set_arg_bounds(model, check_args=True, **kwargs):
     # if variance needs to be resetted, do this at last
     var_bnds = []
     for arg, bounds in kwargs.items():
-        if model._fix and arg in model.fixed:
+        if model._init and arg in model.fixed:
             raise ValueError(f"Can't set bounds for fixed argument '{arg}'")
         if not check_bounds(bounds):
             raise ValueError(
@@ -598,6 +598,7 @@ def compare(this, that):
     equal &= np.all(np.isclose(this.anis, that.anis))
     equal &= np.all(np.isclose(this.angles, that.angles))
     equal &= np.isclose(this.rescale, that.rescale)
+    equal &= np.isclose(this.geo_scale, that.geo_scale)
     equal &= this.latlon == that.latlon
     equal &= this.temporal == that.temporal
     for opt in this.opt_arg:
@@ -614,39 +615,24 @@ def model_repr(model):  # pragma: no cover
     model : :any:`CovModel`
         The covariance model in use.
     """
-    m = model
-    p = model._prec
-    opt_str = ""
+    m, p = model, model._prec
+    ani_str, ang_str, o_str, r_str = "", "", "", ""
     t_str = ", temporal=True" if m.temporal else ""
+    d_str = f"latlon={m.latlon}" if m.latlon else f"dim={m.spatial_dim}"
+    p_str = f", var={m.var:.{p}}, len_scale={m.len_scale:.{p}}"
+    p_str += "" if np.isclose(m.nugget, 0) else f", nugget={m.nugget:.{p}}"
     if not np.isclose(m.rescale, m.default_rescale()):
-        opt_str += f", rescale={m.rescale:.{p}}"
+        o_str += f", rescale={m.rescale:.{p}}"
     for opt in m.opt_arg:
-        opt_str += f", {opt}={getattr(m, opt):.{p}}"
+        o_str += f", {opt}={getattr(m, opt):.{p}}"
     if m.latlon:
-        ani_str = (
-            ""
-            if m.is_isotropic or not m.temporal
-            else f", anis={m.anis[-1]:.{p}}"
-        )
-        r_str = (
-            ""
-            if np.isclose(m.geo_scale, 1)
-            else f", geo_scale={m.geo_scale:.{p}}"
-        )
-        repr_str = (
-            f"{m.name}(latlon={m.latlon}{t_str}, var={m.var:.{p}}, "
-            f"len_scale={m.len_scale:.{p}}, nugget={m.nugget:.{p}}"
-            f"{ani_str}{r_str}{opt_str})"
-        )
+        if not m.is_isotropic and m.temporal:
+            ani_str = f", anis={m.anis[-1]:.{p}}"
+        if not np.isclose(m.geo_scale, 1):
+            r_str = f", geo_scale={m.geo_scale:.{p}}"
     else:
-        # only print anis and angles if model is anisotropic or rotated
-        ani_str = "" if m.is_isotropic else f", anis={list_format(m.anis, p)}"
-        ang_str = (
-            f", angles={list_format(m.angles, p)}" if m.do_rotation else ""
-        )
-        repr_str = (
-            f"{m.name}(dim={m.spatial_dim}{t_str}, var={m.var:.{p}}, "
-            f"len_scale={m.len_scale:.{p}}, nugget={m.nugget:.{p}}"
-            f"{ani_str}{ang_str}{opt_str})"
-        )
-    return repr_str
+        if not m.is_isotropic:
+            ani_str = f", anis={list_format(m.anis, p)}"
+        if m.do_rotation:
+            ang_str = f", angles={list_format(m.angles, p)}"
+    return f"{m.name}({d_str}{t_str}{p_str}{ani_str}{ang_str}{r_str}{o_str})"

From 92210adabaae3c5052830d027d712bd71d1fce26 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Sebastian=20M=C3=BCller?= <mueller.seb@posteo.de>
Date: Mon, 29 Jul 2024 00:45:59 +0200
Subject: [PATCH 08/58] CovModel: remove force arguments mechanic

---
 src/gstools/covmodel/base.py | 12 ------------
 1 file changed, 12 deletions(-)

diff --git a/src/gstools/covmodel/base.py b/src/gstools/covmodel/base.py
index dcbffe57e..d1d1cf0bb 100644
--- a/src/gstools/covmodel/base.py
+++ b/src/gstools/covmodel/base.py
@@ -174,14 +174,6 @@ def __init__(
 
         # indicate that arguments are fixed (True after __init__)
         self._init = False
-        # force input values
-        forced = self.force_values()
-        var = forced.get("var", var)
-        len_scale = forced.get("len_scale", len_scale)
-        nugget = forced.get("nugget", nugget)
-        integral_scale = forced.get("integral_scale", integral_scale)
-        for opt in opt_arg:
-            opt_arg[opt] = forced.get(opt, opt_arg[opt])
         # prepare dim setting
         self._dim = None
         self._hankel_kw = None
@@ -472,10 +464,6 @@ def always_fixed(self):
         """Provide set of fixed arguments."""
         return set()
 
-    def force_values(self):
-        """:class:`dict`: Forced values for arguments."""
-        return {}
-
     def default_opt_arg(self):
         """Provide default optional arguments by the user.
 

From 8999bf98088446a91d728c4be376db474f7d06b8 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Sebastian=20M=C3=BCller?= <mueller.seb@posteo.de>
Date: Mon, 29 Jul 2024 00:50:43 +0200
Subject: [PATCH 09/58] pylint fixes

---
 pyproject.toml               | 2 +-
 src/gstools/covmodel/base.py | 1 -
 2 files changed, 1 insertion(+), 2 deletions(-)

diff --git a/pyproject.toml b/pyproject.toml
index 090c7c635..88185bf30 100644
--- a/pyproject.toml
+++ b/pyproject.toml
@@ -155,7 +155,7 @@ target-version = [
     max-locals = 50
     max-branches = 30
     max-statements = 85
-    max-attributes = 25
+    max-attributes = 30
     max-public-methods = 80
 
 [tool.cibuildwheel]
diff --git a/src/gstools/covmodel/base.py b/src/gstools/covmodel/base.py
index d1d1cf0bb..56f776af0 100644
--- a/src/gstools/covmodel/base.py
+++ b/src/gstools/covmodel/base.py
@@ -21,7 +21,6 @@
 from gstools.covmodel.tools import (
     _init_subclass,
     check_arg_bounds,
-    check_bounds,
     compare,
     default_arg_from_bounds,
     model_repr,

From d643b450b86ff82f29f15c092c7d8fa271a1466d Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Sebastian=20M=C3=BCller?= <mueller.seb@posteo.de>
Date: Sun, 11 Aug 2024 12:16:13 +0200
Subject: [PATCH 10/58] TPL: remove var_raw and var_factor logic and only add
 intensity as property to TPL models

---
 src/gstools/covmodel/base.py       | 44 ++++--------------------------
 src/gstools/covmodel/fit.py        | 32 ++++------------------
 src/gstools/covmodel/tools.py      |  1 -
 src/gstools/covmodel/tpl_models.py | 26 +++++++++++-------
 4 files changed, 28 insertions(+), 75 deletions(-)

diff --git a/src/gstools/covmodel/base.py b/src/gstools/covmodel/base.py
index 56f776af0..fce1e64b9 100644
--- a/src/gstools/covmodel/base.py
+++ b/src/gstools/covmodel/base.py
@@ -126,12 +126,6 @@ class CovModel:
         If given, the model dimension will be determined from this spatial dimension
         and the possible temporal dimension if temporal is ture.
         Default: None
-    var_raw : :class:`float` or :any:`None`, optional
-        raw variance of the model which will be multiplied with
-        :any:`CovModel.var_factor` to result in the actual variance.
-        If given, ``var`` will be ignored.
-        (This is just for models that override :any:`CovModel.var_factor`)
-        Default: :any:`None`
     hankel_kw: :class:`dict` or :any:`None`, optional
         Modify the init-arguments of
         :any:`hankel.SymmetricFourierTransform`
@@ -160,7 +154,6 @@ def __init__(
         geo_scale=RADIAN_SCALE,
         temporal=False,
         spatial_dim=None,
-        var_raw=None,
         hankel_kw=None,
         fixed=None,
         **opt_arg,
@@ -232,20 +225,12 @@ def __init__(
             self.dim, angles, self.latlon, self.temporal
         )
 
-        # set var at last, because of the var_factor (to be right initialized)
-        if var_raw is None:
-            self._var = None
-            self.var = var
-        else:
-            self._var = float(var_raw)
+        self._var = None
+        self.var = var
+
         if integral_scale is not None:
             self.integral_scale = integral_scale
-        # set var again, if int_scale affects var_factor
-        if var_raw is None:
-            self._var = None
-            self.var = var
-        else:
-            self._var = float(var_raw)
+
         # final check for parameter bounds
         self.check_arg_bounds()
         # additional checks for the optional arguments (provided by user)
@@ -500,10 +485,6 @@ def fix_dim(self):
         """Set a fix dimension for the model."""
         return None
 
-    def var_factor(self):
-        """Factor for the variance."""
-        return 1.0
-
     def default_rescale(self):
         """Provide default rescaling factor."""
         return 1.0
@@ -963,24 +944,11 @@ def dim(self, dim):
     @property
     def var(self):
         """:class:`float`: The variance of the model."""
-        return self._var * self.var_factor()
+        return self._var
 
     @var.setter
     def var(self, var):
-        self._var = float(var) / self.var_factor()
-        self.check_arg_bounds()
-
-    @property
-    def var_raw(self):
-        """:class:`float`: The raw variance of the model without factor.
-
-        (See. CovModel.var_factor)
-        """
-        return self._var
-
-    @var_raw.setter
-    def var_raw(self, var_raw):
-        self._var = float(var_raw)
+        self._var = float(var)
         self.check_arg_bounds()
 
     @property
diff --git a/src/gstools/covmodel/fit.py b/src/gstools/covmodel/fit.py
index f19eaac35..e6d6b03d5 100755
--- a/src/gstools/covmodel/fit.py
+++ b/src/gstools/covmodel/fit.py
@@ -19,6 +19,7 @@
 __all__ = ["fit_variogram"]
 
 
+# this should be: set(mod.iso_arg) - set(mod.opt_arg)
 DEFAULT_PARA = ["var", "len_scale", "nugget"]
 
 
@@ -217,21 +218,13 @@ def fit_variogram(
 
 def _pre_para(model, para_select, sill, anis):
     """Preprocess selected parameters."""
-    var_last = False
-    var_tmp = 0.0  # init value
     for par in para_select:
         if par not in model.arg_bounds:
             raise ValueError(f"fit: unknown parameter in selection: {par}")
         if not isinstance(para_select[par], bool):
-            if par == "var":
-                var_last = True
-                var_tmp = float(para_select[par])
-            else:
-                setattr(model, par, float(para_select[par]))
+            setattr(model, par, float(para_select[par]))
             para_select[par] = False
-    # set variance last due to possible recalculations
-    if var_last:
-        model.var = var_tmp
+
     # remove those that were set to True
     para_select = {k: v for k, v in para_select.items() if not v}
     # handling the sill
@@ -426,9 +419,9 @@ def curve(x, arg1, *args):
         para_skip = 0
         opt_skip = 0
         if para["var"]:
-            var_tmp = args[para_skip]
+            model.var = args[para_skip]
             if constrain_sill:
-                nugget_tmp = sill - var_tmp
+                nugget_tmp = sill - model.var
                 # punishment, if resulting nugget out of range for fixed sill
                 if check_arg_in_bounds(model, "nugget", nugget_tmp) > 0:
                     return np.full_like(x, np.inf)
@@ -445,12 +438,6 @@ def curve(x, arg1, *args):
             if para[opt]:
                 setattr(model, opt, args[para_skip + opt_skip])
                 opt_skip += 1
-        # set var at last because of var_factor (other parameter needed)
-        if para["var"]:
-            model.var = var_tmp
-        # needs to be reset for TPL models when len_scale was changed
-        else:
-            model.var = var_save
         if is_dir_vario:
             if anis:
                 model.anis = args[1 - model.dim :]
@@ -469,13 +456,9 @@ def _post_fitting(model, para, popt, anis, is_dir_vario):
     fit_para = {}
     para_skip = 0
     opt_skip = 0
-    var_tmp = 0.0  # init value
     for par in DEFAULT_PARA:
         if para[par]:
-            if par == "var":  # set variance last
-                var_tmp = popt[para_skip]
-            else:
-                setattr(model, par, popt[para_skip])
+            setattr(model, par, popt[para_skip])
             fit_para[par] = popt[para_skip]
             para_skip += 1
         else:
@@ -491,9 +474,6 @@ def _post_fitting(model, para, popt, anis, is_dir_vario):
         if anis:
             model.anis = popt[1 - model.dim :]
         fit_para["anis"] = model.anis
-    # set var at last because of var_factor (other parameter needed)
-    if para["var"]:
-        model.var = var_tmp
     return fit_para
 
 
diff --git a/src/gstools/covmodel/tools.py b/src/gstools/covmodel/tools.py
index 9d53de583..1c3c9ca0c 100644
--- a/src/gstools/covmodel/tools.py
+++ b/src/gstools/covmodel/tools.py
@@ -592,7 +592,6 @@ def compare(this, that):
     equal = True
     equal &= this.name == that.name
     equal &= np.isclose(this.var, that.var)
-    equal &= np.isclose(this.var_raw, that.var_raw)  # ?! needless?
     equal &= np.isclose(this.nugget, that.nugget)
     equal &= np.isclose(this.len_scale, that.len_scale)
     equal &= np.all(np.isclose(this.anis, that.anis))
diff --git a/src/gstools/covmodel/tpl_models.py b/src/gstools/covmodel/tpl_models.py
index b728e7b98..40ed4ec8e 100644
--- a/src/gstools/covmodel/tpl_models.py
+++ b/src/gstools/covmodel/tpl_models.py
@@ -31,6 +31,19 @@
 class TPLCovModel(CovModel):
     """Truncated-Power-Law Covariance Model base class for super-position."""
 
+    @property
+    def intensity(self):
+        """:class:`float`: Intensity of variation."""
+        return self.var / self.intensity_scaling
+
+    @property
+    def intensity_scaling(self):
+        """:class:`float`: Scaling of Intensity to result in variance."""
+        return (
+            self.len_up_rescaled ** (2 * self.hurst)
+            - self.len_low_rescaled ** (2 * self.hurst)
+        ) / (2 * self.hurst)
+
     @property
     def len_up(self):
         """:class:`float`: Upper length scale truncation of the model.
@@ -55,13 +68,6 @@ def len_low_rescaled(self):
         """
         return self.len_low / self.rescale
 
-    def var_factor(self):
-        """Factor for C (intensity of variation) to result in variance."""
-        return (
-            self.len_up_rescaled ** (2 * self.hurst)
-            - self.len_low_rescaled ** (2 * self.hurst)
-        ) / (2 * self.hurst)
-
     def cor(self, h):
         """TPL - normalized correlation function."""
 
@@ -125,7 +131,7 @@ class TPLGaussian(TPLCovModel):
         * :math:`C>0` :
           scaling factor from the Power-Law (intensity of variation)
           This parameter will be calculated internally by the given variance.
-          You can access C directly by ``model.var_raw``
+          You can access C by ``model.intensity``
         * :math:`0<H<1` : hurst coefficient (``model.hurst``)
         * :math:`\ell_{\mathrm{low}}\geq 0` : lower length scale truncation
           of the model (``model.len_low``)
@@ -258,7 +264,7 @@ class TPLExponential(TPLCovModel):
         * :math:`C>0` :
           scaling factor from the Power-Law (intensity of variation)
           This parameter will be calculated internally by the given variance.
-          You can access C directly by ``model.var_raw``
+          You can access C by ``model.intensity``
         * :math:`0<H<\frac{1}{2}` : hurst coefficient (``model.hurst``)
         * :math:`\ell_{\mathrm{low}}\geq 0` : lower length scale truncation
           of the model (``model.len_low``)
@@ -398,7 +404,7 @@ class TPLStable(TPLCovModel):
         * :math:`C>0` :
           scaling factor from the Power-Law (intensity of variation)
           This parameter will be calculated internally by the given variance.
-          You can access C directly by ``model.var_raw``
+          You can access C by ``model.intensity``
         * :math:`0<H<\frac{\alpha}{2}` : hurst coefficient (``model.hurst``)
         * :math:`\ell_{\mathrm{low}}\geq 0` : lower length scale truncation
           of the model (``model.len_low``)

From 81e5421ff3f0e9f1a5b81f0f78d193091be655f0 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Sebastian=20M=C3=BCller?= <mueller.seb@posteo.de>
Date: Sun, 11 Aug 2024 12:16:38 +0200
Subject: [PATCH 11/58] tests: remove tests for var_raw

---
 tests/test_covmodel.py | 8 +-------
 1 file changed, 1 insertion(+), 7 deletions(-)

diff --git a/tests/test_covmodel.py b/tests/test_covmodel.py
index a2729dd68..c42e47746 100644
--- a/tests/test_covmodel.py
+++ b/tests/test_covmodel.py
@@ -201,15 +201,9 @@ def test_tpl_models(self):
                     self.gamma_x, self.gamma_y, sill=1.1, nugget=False
                 )
                 self.assertAlmostEqual(model.var, 1.1, delta=1e-5)
-                # check var_raw handling
-                model = Model(var_raw=1, len_low=0, integral_scale=10)
-                var_save = model.var
-                model.var_raw = 1.1
-                self.assertAlmostEqual(model.var, var_save * 1.1)
-                self.assertAlmostEqual(model.integral_scale, 10)
                 # integral scale is not setable when len_low is not 0
                 with self.assertRaises(ValueError):
-                    Model(var_raw=1, len_low=5, integral_scale=10)
+                    Model(len_low=5, integral_scale=10)
 
     def test_fitting(self):
         for Model in self.std_cov_models:

From 669e9a203a5d40d7b35f11896e64674f3c57e0b2 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Sebastian=20M=C3=BCller?= <mueller.seb@posteo.de>
Date: Sun, 11 Aug 2024 16:58:40 +0200
Subject: [PATCH 12/58] CovModel: remove mechanism for fixed args again

---
 src/gstools/covmodel/base.py  | 36 ++++-------------------------------
 src/gstools/covmodel/fit.py   |  2 --
 src/gstools/covmodel/tools.py |  2 --
 3 files changed, 4 insertions(+), 36 deletions(-)

diff --git a/src/gstools/covmodel/base.py b/src/gstools/covmodel/base.py
index fce1e64b9..4a1ee76c2 100644
--- a/src/gstools/covmodel/base.py
+++ b/src/gstools/covmodel/base.py
@@ -132,8 +132,6 @@ class CovModel:
         used for the spectrum calculation. Use with caution (Better: Don't!).
         ``None`` is equivalent to ``{"a": -1, "b": 1, "N": 1000, "h": 0.001}``.
         Default: :any:`None`
-    fixed: :class:`set` or :any:`None`, optional
-        Names of fixed arguments. Default: :any:`None`
     **opt_arg
         Optional arguments are covered by these keyword arguments.
         If present, they are described in the section `Other Parameters`.
@@ -155,7 +153,6 @@ def __init__(
         temporal=False,
         spatial_dim=None,
         hankel_kw=None,
-        fixed=None,
         **opt_arg,
     ):
         # assert, that we use a subclass
@@ -164,14 +161,13 @@ def __init__(
         if not hasattr(self, "variogram"):
             raise TypeError("Don't instantiate 'CovModel' directly!")
 
-        # indicate that arguments are fixed (True after __init__)
+        # indicator for initialization status (True after __init__)
         self._init = False
         # prepare dim setting
         self._dim = None
         self._hankel_kw = None
         self._sft = None
         # prepare parameters (they are checked in dim setting)
-        self._fixed = None
         self._rescale = None
         self._len_scale = None
         self._anis = None
@@ -196,16 +192,6 @@ def __init__(
         # optional arguments for the variogram-model
         set_opt_args(self, opt_arg)
 
-        # check fixed
-        fixed = set(fixed) if fixed is not None else set()
-        fixed = fixed | self.always_fixed()
-        valid_fixed = set(self.iso_arg)  # | set(["anis"])
-        if not fixed <= valid_fixed:
-            raise ValueError(
-                f"CovModel: unknown names in 'fixed': {fixed - valid_fixed}"
-            )
-        self._fixed = fixed
-
         # set standard boundaries for variance, len_scale, nugget and opt_arg
         bounds = self.default_arg_bounds()
         bounds.update(self.default_opt_arg_bounds())
@@ -213,6 +199,7 @@ def __init__(
 
         # set parameters
         self.rescale = rescale
+        self._var = float(var)
         self._nugget = float(nugget)
 
         # set anisotropy and len_scale, disable anisotropy for latlon models
@@ -225,9 +212,6 @@ def __init__(
             self.dim, angles, self.latlon, self.temporal
         )
 
-        self._var = None
-        self.var = var
-
         if integral_scale is not None:
             self.integral_scale = integral_scale
 
@@ -235,13 +219,10 @@ def __init__(
         self.check_arg_bounds()
         # additional checks for the optional arguments (provided by user)
         self.check_opt_arg()
-        # set fixed bounds after checking original bounds
-        for arg in self.fixed:
-            val = getattr(self, arg)
-            self.set_arg_bounds(check_args=False, **{arg: (val, val, "cc")})
-        self._init = True
         # precision for printing
         self._prec = 3
+        # initialized
+        self._init = True
 
     # one of these functions needs to be overridden
     def __init_subclass__(cls):
@@ -444,10 +425,6 @@ def pykrige_kwargs(self):
 
     # methods for optional/default arguments (can be overridden)
 
-    def always_fixed(self):
-        """Provide set of fixed arguments."""
-        return set()
-
     def default_opt_arg(self):
         """Provide default optional arguments by the user.
 
@@ -801,11 +778,6 @@ def check_arg_bounds(self):
 
     # bounds properties
 
-    @property
-    def fixed(self):
-        """:class:`set`: Set with names of fixed arguments."""
-        return self._fixed
-
     @property
     def var_bounds(self):
         """:class:`list`: Bounds for the variance.
diff --git a/src/gstools/covmodel/fit.py b/src/gstools/covmodel/fit.py
index e6d6b03d5..b40fc8457 100755
--- a/src/gstools/covmodel/fit.py
+++ b/src/gstools/covmodel/fit.py
@@ -166,8 +166,6 @@ def fit_variogram(
 
     The fitted parameters will be instantly set in the model.
     """
-    for para in model.fixed:
-        para_select[para] = False
     # preprocess selected parameters
     para, sill, constrain_sill, anis = _pre_para(
         model, para_select, sill, anis
diff --git a/src/gstools/covmodel/tools.py b/src/gstools/covmodel/tools.py
index 1c3c9ca0c..78a90a7bd 100644
--- a/src/gstools/covmodel/tools.py
+++ b/src/gstools/covmodel/tools.py
@@ -458,8 +458,6 @@ def set_arg_bounds(model, check_args=True, **kwargs):
     # if variance needs to be resetted, do this at last
     var_bnds = []
     for arg, bounds in kwargs.items():
-        if model._init and arg in model.fixed:
-            raise ValueError(f"Can't set bounds for fixed argument '{arg}'")
         if not check_bounds(bounds):
             raise ValueError(
                 f"Given bounds for '{arg}' are not valid, got: {bounds}"

From 01d1d8233ae44fe63ad84e170a0e68d1b296592b Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Sebastian=20M=C3=BCller?= <mueller.seb@posteo.de>
Date: Mon, 12 Aug 2024 00:34:57 +0200
Subject: [PATCH 13/58] CovModel: simplify fitting routines

---
 src/gstools/covmodel/fit.py | 62 +++++++++----------------------------
 1 file changed, 14 insertions(+), 48 deletions(-)

diff --git a/src/gstools/covmodel/fit.py b/src/gstools/covmodel/fit.py
index b40fc8457..b5836c0d3 100755
--- a/src/gstools/covmodel/fit.py
+++ b/src/gstools/covmodel/fit.py
@@ -19,10 +19,6 @@
 __all__ = ["fit_variogram"]
 
 
-# this should be: set(mod.iso_arg) - set(mod.opt_arg)
-DEFAULT_PARA = ["var", "len_scale", "nugget"]
-
-
 def fit_variogram(
     model,
     x_data,
@@ -263,8 +259,7 @@ def _pre_para(model, para_select, sill, anis):
     else:
         constrain_sill = False
     # select all parameters to be fitted
-    para = {par: True for par in DEFAULT_PARA}
-    para.update({opt: True for opt in model.opt_arg})
+    para = {par: True for par in model.iso_arg}
     # now deselect unwanted parameters
     para.update(para_select)
     # check if anisotropy should be fitted or set
@@ -363,7 +358,7 @@ def _init_curve_fit_para(model, para, init_guess, constrain_sill, sill, anis):
     low_bounds = []
     top_bounds = []
     init_guess_list = []
-    for par in DEFAULT_PARA:
+    for par in model.iso_arg:
         if para[par]:
             low_bounds.append(model.arg_bounds[par][0])
             if par == "var" and constrain_sill:  # var <= sill in this case
@@ -376,16 +371,6 @@ def _init_curve_fit_para(model, para, init_guess, constrain_sill, sill, anis):
                     default=init_guess[par],
                 )
             )
-    for opt in model.opt_arg:
-        if para[opt]:
-            low_bounds.append(model.arg_bounds[opt][0])
-            top_bounds.append(model.arg_bounds[opt][1])
-            init_guess_list.append(
-                _init_guess(
-                    bounds=[low_bounds[-1], top_bounds[-1]],
-                    default=init_guess[opt],
-                )
-            )
     if anis:
         for i in range(model.dim - 1):
             low_bounds.append(model.anis_bounds[0])
@@ -408,34 +393,23 @@ def _init_guess(bounds, default):
 
 def _get_curve(model, para, constrain_sill, sill, anis, is_dir_vario):
     """Create the curve for scipys curve_fit."""
-    var_save = model.var
 
     # we need arg1, otherwise curve_fit throws an error (bug?!)
     def curve(x, arg1, *args):
         """Adapted Variogram function."""
         args = (arg1,) + args
         para_skip = 0
-        opt_skip = 0
-        if para["var"]:
-            model.var = args[para_skip]
-            if constrain_sill:
-                nugget_tmp = sill - model.var
-                # punishment, if resulting nugget out of range for fixed sill
-                if check_arg_in_bounds(model, "nugget", nugget_tmp) > 0:
-                    return np.full_like(x, np.inf)
-                # nugget estimation deselected in this case
-                model.nugget = nugget_tmp
-            para_skip += 1
-        if para["len_scale"]:
-            model.len_scale = args[para_skip]
-            para_skip += 1
-        if para["nugget"]:
-            model.nugget = args[para_skip]
-            para_skip += 1
-        for opt in model.opt_arg:
-            if para[opt]:
-                setattr(model, opt, args[para_skip + opt_skip])
-                opt_skip += 1
+        for par in model.iso_arg:
+            if para[par]:
+                setattr(model, par, args[para_skip])
+                para_skip += 1
+        if constrain_sill:
+            nugget_tmp = sill - model.var
+            # punishment, if resulting nugget out of range for fixed sill
+            if check_arg_in_bounds(model, "nugget", nugget_tmp) > 0:
+                return np.full_like(x, np.inf)
+            # nugget estimation deselected in this case
+            model.nugget = nugget_tmp
         if is_dir_vario:
             if anis:
                 model.anis = args[1 - model.dim :]
@@ -453,21 +427,13 @@ def _post_fitting(model, para, popt, anis, is_dir_vario):
     """Postprocess fitting results and application to model."""
     fit_para = {}
     para_skip = 0
-    opt_skip = 0
-    for par in DEFAULT_PARA:
+    for par in model.iso_arg:
         if para[par]:
             setattr(model, par, popt[para_skip])
             fit_para[par] = popt[para_skip]
             para_skip += 1
         else:
             fit_para[par] = getattr(model, par)
-    for opt in model.opt_arg:
-        if para[opt]:
-            setattr(model, opt, popt[para_skip + opt_skip])
-            fit_para[opt] = popt[para_skip + opt_skip]
-            opt_skip += 1
-        else:
-            fit_para[opt] = getattr(model, opt)
     if is_dir_vario:
         if anis:
             model.anis = popt[1 - model.dim :]

From 418038f521a4f566b8016990db67547be1f39b69 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Sebastian=20M=C3=BCller?= <mueller.seb@posteo.de>
Date: Mon, 12 Aug 2024 11:34:51 +0200
Subject: [PATCH 14/58] CovModel: fix setting integral scale as list of values

---
 src/gstools/covmodel/base.py | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/src/gstools/covmodel/base.py b/src/gstools/covmodel/base.py
index 4a1ee76c2..a02a46516 100644
--- a/src/gstools/covmodel/base.py
+++ b/src/gstools/covmodel/base.py
@@ -1004,7 +1004,7 @@ def integral_scale(self, integral_scale):
         old_scale = self.integral_scale
         self.anis = anis
         self.len_scale = self.len_scale * int_scale / old_scale
-        if not np.isclose(self.integral_scale, integral_scale, rtol=1e-3):
+        if not np.isclose(self.integral_scale, int_scale, rtol=1e-3):
             raise ValueError(
                 f"{self.name}: Integral scale could not be set correctly! "
                 "Please just provide a 'len_scale'!"

From d8b037fe994c422fc29c147d3b01ef152e927df5 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Sebastian=20M=C3=BCller?= <mueller.seb@posteo.de>
Date: Mon, 12 Aug 2024 12:38:40 +0200
Subject: [PATCH 15/58] no need to set var last anymore

---
 src/gstools/covmodel/tools.py | 15 +++------------
 1 file changed, 3 insertions(+), 12 deletions(-)

diff --git a/src/gstools/covmodel/tools.py b/src/gstools/covmodel/tools.py
index 78a90a7bd..448dc6c11 100644
--- a/src/gstools/covmodel/tools.py
+++ b/src/gstools/covmodel/tools.py
@@ -455,18 +455,14 @@ def set_arg_bounds(model, check_args=True, **kwargs):
         <type> is one of ``"oo"``, ``"cc"``, ``"oc"`` or ``"co"``
         to define if the bounds are open ("o") or closed ("c").
     """
-    # if variance needs to be resetted, do this at last
-    var_bnds = []
     for arg, bounds in kwargs.items():
         if not check_bounds(bounds):
-            raise ValueError(
-                f"Given bounds for '{arg}' are not valid, got: {bounds}"
-            )
+            msg = f"Given bounds for '{arg}' are not valid, got: {bounds}"
+            raise ValueError(msg)
         if arg in model.opt_arg:
             model._opt_arg_bounds[arg] = bounds
         elif arg == "var":
-            var_bnds = bounds
-            continue
+            model._var_bounds = bounds
         elif arg == "len_scale":
             model._len_scale_bounds = bounds
         elif arg == "nugget":
@@ -481,11 +477,6 @@ def set_arg_bounds(model, check_args=True, **kwargs):
                 setattr(model, arg, [def_arg] * (model.dim - 1))
             else:
                 setattr(model, arg, def_arg)
-    # set var last like always
-    if var_bnds:
-        model._var_bounds = var_bnds
-        if check_args and check_arg_in_bounds(model, "var") > 0:
-            model.var = default_arg_from_bounds(var_bnds)
 
 
 def check_arg_bounds(model):

From 0b632253efcb6a1b69f76bf72b393ff2b1e94fe1 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Sebastian=20M=C3=BCller?= <mueller.seb@posteo.de>
Date: Mon, 12 Aug 2024 12:38:57 +0200
Subject: [PATCH 16/58] fit: add comment

---
 src/gstools/covmodel/fit.py | 4 ++--
 1 file changed, 2 insertions(+), 2 deletions(-)

diff --git a/src/gstools/covmodel/fit.py b/src/gstools/covmodel/fit.py
index b5836c0d3..5ba398ef9 100755
--- a/src/gstools/covmodel/fit.py
+++ b/src/gstools/covmodel/fit.py
@@ -169,7 +169,7 @@ def fit_variogram(
     if not any(para.values()):
         raise ValueError("fit: no parameters selected for fitting.")
     # check curve_fit kwargs
-    curve_fit_kwargs = {} if curve_fit_kwargs is None else curve_fit_kwargs
+    curve_fit_kwargs = curve_fit_kwargs or {}
     # check method
     if method not in ["trf", "dogbox"]:
         raise ValueError("fit: method needs to be either 'trf' or 'dogbox'")
@@ -212,13 +212,13 @@ def fit_variogram(
 
 def _pre_para(model, para_select, sill, anis):
     """Preprocess selected parameters."""
+    # if values given, set them in the model, afterwards all entries are bool
     for par in para_select:
         if par not in model.arg_bounds:
             raise ValueError(f"fit: unknown parameter in selection: {par}")
         if not isinstance(para_select[par], bool):
             setattr(model, par, float(para_select[par]))
             para_select[par] = False
-
     # remove those that were set to True
     para_select = {k: v for k, v in para_select.items() if not v}
     # handling the sill

From 21d0e4a27d8b3dce26a8c2956ea398c0db601c50 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Sebastian=20M=C3=BCller?= <mueller.seb@posteo.de>
Date: Wed, 14 Aug 2024 12:12:45 +0200
Subject: [PATCH 17/58] add ratio error class; better arg checking routines

---
 src/gstools/covmodel/tools.py | 30 ++++++++++++++++++------------
 1 file changed, 18 insertions(+), 12 deletions(-)

diff --git a/src/gstools/covmodel/tools.py b/src/gstools/covmodel/tools.py
index 448dc6c11..50ef04dd5 100644
--- a/src/gstools/covmodel/tools.py
+++ b/src/gstools/covmodel/tools.py
@@ -6,6 +6,7 @@
 The following classes and functions are provided
 
 .. autosummary::
+   RatioError
    AttributeWarning
    rad_fac
    set_opt_args
@@ -34,6 +35,7 @@
 from gstools.tools.misc import list_format
 
 __all__ = [
+    "RatioError",
     "AttributeWarning",
     "rad_fac",
     "set_opt_args",
@@ -52,6 +54,10 @@
 ]
 
 
+class RatioError(Exception):
+    """Error for invalid ratios in SumModel."""
+
+
 class AttributeWarning(UserWarning):
     """Attribute warning for CovModel class."""
 
@@ -311,7 +317,7 @@ def check_bounds(bounds):
     return True
 
 
-def check_arg_in_bounds(model, arg, val=None):
+def check_arg_in_bounds(model, arg, val=None, error=False):
     """Check if given argument value is in bounds of the given model."""
     if arg not in model.arg_bounds:
         raise ValueError(f"check bounds: unknown argument: {arg}")
@@ -320,6 +326,7 @@ def check_arg_in_bounds(model, arg, val=None):
     val = np.asarray(val)
     error_case = 0
     if len(bnd) == 2:
+        bnd = list(bnd)
         bnd.append("cc")  # use closed intervals by default
     if bnd[2][0] == "c":
         if np.any(val < bnd[0]):
@@ -333,6 +340,15 @@ def check_arg_in_bounds(model, arg, val=None):
     else:
         if np.any(val >= bnd[1]):
             error_case = 4
+    if error:
+        if error_case == 1:
+            raise ValueError(f"{arg} needs to be >= {bnd[0]}, got: {val}")
+        if error_case == 2:
+            raise ValueError(f"{arg} needs to be > {bnd[0]}, got: {val}")
+        if error_case == 3:
+            raise ValueError(f"{arg} needs to be <= {bnd[1]}, got: {val}")
+        if error_case == 4:
+            raise ValueError(f"{arg} needs to be < {bnd[1]}, got: {val}")
     return error_case
 
 
@@ -497,17 +513,7 @@ def check_arg_bounds(model):
     for arg in model.arg_bounds:
         if not model.arg_bounds[arg]:
             continue  # no bounds given during init (called from self.dim)
-        bnd = list(model.arg_bounds[arg])
-        val = getattr(model, arg)
-        error_case = check_arg_in_bounds(model, arg)
-        if error_case == 1:
-            raise ValueError(f"{arg} needs to be >= {bnd[0]}, got: {val}")
-        if error_case == 2:
-            raise ValueError(f"{arg} needs to be > {bnd[0]}, got: {val}")
-        if error_case == 3:
-            raise ValueError(f"{arg} needs to be <= {bnd[1]}, got: {val}")
-        if error_case == 4:
-            raise ValueError(f"{arg} needs to be < {bnd[1]}, got: {val}")
+        check_arg_in_bounds(model, arg, error=True)
 
 
 def set_dim(model, dim):

From e6374db261076351920c3fb1857690c558a57b90 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Sebastian=20M=C3=BCller?= <mueller.seb@posteo.de>
Date: Thu, 15 Aug 2024 23:26:21 +0200
Subject: [PATCH 18/58] CovModel: better dim setting

---
 src/gstools/covmodel/tools.py | 17 +++++++----------
 1 file changed, 7 insertions(+), 10 deletions(-)

diff --git a/src/gstools/covmodel/tools.py b/src/gstools/covmodel/tools.py
index 50ef04dd5..a70c8118a 100644
--- a/src/gstools/covmodel/tools.py
+++ b/src/gstools/covmodel/tools.py
@@ -511,8 +511,6 @@ def check_arg_bounds(model):
     """
     # check var, len_scale, nugget and optional-arguments
     for arg in model.arg_bounds:
-        if not model.arg_bounds[arg]:
-            continue  # no bounds given during init (called from self.dim)
         check_arg_in_bounds(model, arg, error=True)
 
 
@@ -557,16 +555,15 @@ def set_dim(model, dim):
     model._dim = int(dim)
     # create fourier transform just once (recreate for dim change)
     model._sft = SFT(ndim=model.dim, **model.hankel_kw)
-    # recalculate dimension related parameters
-    if model._anis is not None:
-        model._len_scale, model._anis = set_len_anis(
-            model.dim, model._len_scale, model._anis
+    # recalculate dimension related parameters (if model initialized)
+    if model._init:
+        model.len_scale, model.anis = set_len_anis(
+            model.dim, model.len_scale, model.anis
         )
-    if model._angles is not None:
-        model._angles = set_model_angles(
-            model.dim, model._angles, model.latlon, model.temporal
+        model.angles = set_model_angles(
+            model.dim, model.angles, model.latlon, model.temporal
         )
-    model.check_arg_bounds()
+        model.check_arg_bounds()
 
 
 def compare(this, that):

From dcb003bffec9fa2f4e0c1ebf0e174620a0183bd7 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Sebastian=20M=C3=BCller?= <mueller.seb@posteo.de>
Date: Fri, 16 Aug 2024 00:04:00 +0200
Subject: [PATCH 19/58] add sum_tools submodule

---
 src/gstools/covmodel/sum_tools.py | 245 ++++++++++++++++++++++++++++++
 src/gstools/covmodel/tools.py     |   5 -
 2 files changed, 245 insertions(+), 5 deletions(-)
 create mode 100644 src/gstools/covmodel/sum_tools.py

diff --git a/src/gstools/covmodel/sum_tools.py b/src/gstools/covmodel/sum_tools.py
new file mode 100644
index 000000000..a4a1647a3
--- /dev/null
+++ b/src/gstools/covmodel/sum_tools.py
@@ -0,0 +1,245 @@
+"""
+GStools subpackage providing tools for sum-models.
+
+.. currentmodule:: gstools.covmodel.sum_tools
+
+The following classes and functions are provided
+
+.. autosummary::
+   RatioError
+   ARG_DEF
+   default_mod_kwargs
+   sum_check
+   sum_default_arg_bounds
+   sum_default_opt_arg_bounds
+   sum_set_norm_var_ratios
+   sum_set_norm_len_ratios
+   sum_model_repr
+"""
+
+import numpy as np
+
+from gstools.covmodel.tools import check_arg_in_bounds
+from gstools.tools import RADIAN_SCALE
+from gstools.tools.misc import list_format
+
+__all__ = [
+    "RatioError",
+    "ARG_DEF",
+    "default_mod_kwargs",
+    "sum_check",
+    "sum_default_arg_bounds",
+    "sum_default_opt_arg_bounds",
+    "sum_set_norm_var_ratios",
+    "sum_set_norm_len_ratios",
+    "sum_model_repr",
+]
+
+
+class RatioError(Exception):
+    """Error for invalid ratios in SumModel."""
+
+
+ARG_DEF = {
+    "dim": 3,
+    "latlon": False,
+    "temporal": False,
+    "geo_scale": RADIAN_SCALE,
+    "spatial_dim": None,
+    "hankel_kw": None,
+}
+"""dict: default model arguments"""
+
+
+def default_mod_kwargs(kwargs):
+    """Generate default model keyword arguments."""
+    mod_kw = {}
+    for arg, default in ARG_DEF.items():
+        mod_kw[arg] = kwargs.get(arg, default)
+    return mod_kw
+
+
+def sum_check(summod):
+    """Check consistency of contained models."""
+    # prevent dim error in anis and angles
+    if any(mod.dim != summod.dim for mod in summod):
+        msg = "SumModel: models need to have same dimension."
+        raise ValueError(msg)
+    if any(mod.latlon != summod.latlon for mod in summod):
+        msg = "SumModel: models need to have same latlon config."
+        raise ValueError(msg)
+    if any(mod.temporal != summod.temporal for mod in summod):
+        msg = "SumModel: models need to have same temporal config."
+        raise ValueError(msg)
+    if not all(np.isclose(mod.nugget, 0) for mod in summod):
+        msg = "SumModel: models need to have 0 nugget."
+        raise ValueError(msg)
+    if not np.allclose([mod.geo_scale for mod in summod], summod.geo_scale):
+        msg = "SumModel: models need to have same geo_scale."
+        raise ValueError(msg)
+    if not all(np.allclose(mod.anis, summod.anis) for mod in summod):
+        msg = "SumModel: models need to have same anisotropy ratios."
+        raise ValueError(msg)
+    if not all(np.allclose(mod.angles, summod.angles) for mod in summod):
+        msg = "SumModel: models need to have same rotation angles."
+        raise ValueError(msg)
+
+
+def sum_default_arg_bounds(summod):
+    """Default boundaries for arguments as dict."""
+    var_bnds = [mod.var_bounds for mod in summod.models]
+    len_bnds = [mod.len_scale_bounds for mod in summod.models]
+    var_lo = sum([bnd[0] for bnd in var_bnds], start=0.0)
+    var_hi = sum([bnd[1] for bnd in var_bnds], start=0.0)
+    len_lo = min([bnd[0] for bnd in len_bnds], default=0.0)
+    len_hi = max([bnd[1] for bnd in len_bnds], default=0.0)
+    res = {
+        "var": (var_lo, var_hi),
+        "len_scale": (len_lo, len_hi),
+        "nugget": (0.0, np.inf, "co"),
+        "anis": (0.0, np.inf, "oo"),
+    }
+    return res
+
+
+def sum_default_opt_arg_bounds(summod):
+    """Defaults boundaries for optional arguments as dict."""
+    bounds = {}
+    for i, mod in enumerate(summod.models):
+        bounds.update(
+            {f"{opt}_{i}": bnd for opt, bnd in mod.opt_arg_bounds.items()}
+        )
+    return bounds
+
+
+def sum_set_norm_var_ratios(summod, ratios, skip=None, var=None):
+    """
+    Set variances of contained models by normalized ratios in [0, 1].
+
+    Ratios are given as normalized ratios in [0, 1] as relative ratio of
+    variance to remaining difference to total variance of the Sum-Model.
+
+    Parameters
+    ----------
+    ratios : iterable
+        Ratios to set. Should have a length of len(models) - len(exclude) - 1
+    skip : iterable, optional
+        Model indices to skip. Should have compatible lenth, by default None
+    var : float, optional
+        Desired variance, by default current variance
+
+    Raises
+    ------
+    ValueError
+        If number of ratios is not matching.
+    """
+    skip = skip or set()
+    if len(summod) != len(ratios) + len(skip) + 1:
+        msg = "SumModel.set_norm_ratios: number of ratios not matching."
+        raise ValueError(msg)
+    ids = range(len(summod))
+    if fail := set(skip) - set(ids):
+        msg = f"SumModel.set_norm_var_ratios: skip ids not valid: {fail}"
+        raise ValueError(msg)
+    var = summod.var if var is None else float(var)
+    check_arg_in_bounds(summod, "var", var, error=True)
+    var_sum = sum(summod.models[i].var for i in skip)
+    if var_sum > var:
+        msg = "SumModel.set_norm_var_ratios: skiped variances to big."
+        raise RatioError(msg)
+    j = 0
+    for i in ids:
+        if i in skip:
+            continue
+        var_diff = var - var_sum
+        # last model gets remaining diff
+        var_set = var_diff * ratios[j] if j < len(ratios) else var_diff
+        summod[i].var = var_set
+        var_sum += var_set
+        j += 1
+
+
+def sum_set_norm_len_ratios(summod, ratios, skip=None, len_scale=None):
+    """
+    Set length scales of contained models by normalized ratios in [0, 1].
+
+    Ratios are given as normalized ratios in [0, 1] as relative ratio of
+    len_scale * var / total_var to remaining difference to
+    total len_scale of the Sum-Model.
+
+    Parameters
+    ----------
+    ratios : iterable
+        Ratios to set. Should have a length of len(models) - len(exclude) - 1
+    skip : iterable, optional
+        Model indices to skip. Should have compatible lenth, by default None
+    len_scale : float, optional
+        Desired len_scale, by default current len_scale
+
+    Raises
+    ------
+    ValueError
+        If number of ratios is not matching.
+    """
+    skip = skip or set()
+    if len(summod) != len(ratios) + len(skip) + 1:
+        msg = "SumModel.set_norm_len_ratios: number of ratios not matching."
+        raise ValueError(msg)
+    ids = range(len(summod))
+    if fail := set(skip) - set(ids):
+        msg = f"SumModel.set_norm_len_ratios: skip ids not valid: {fail}"
+        raise ValueError(msg)
+    len_scale = summod.len_scale if len_scale is None else float(len_scale)
+    check_arg_in_bounds(summod, "len_scale", len_scale, error=True)
+    len_sum = sum(summod[i].len_scale * summod.ratios[i] for i in skip)
+    if len_sum > len_scale:
+        msg = "SumModel.set_norm_len_ratios: skiped length scales to big."
+        raise RatioError(msg)
+    j = 0
+    for i in ids:
+        if i in skip:
+            continue
+        len_diff = len_scale - len_sum
+        # last model gets remaining diff
+        len_set = len_diff * ratios[j] if j < len(ratios) else len_diff
+        summod[i].len_scale = (
+            0.0
+            if np.isclose(summod.ratios[j], 0)
+            else len_set / summod.ratios[j]
+        )
+        len_sum += len_set
+        j += 1
+
+
+def sum_model_repr(summod):  # pragma: no cover
+    """
+    Generate the sum-model string representation.
+
+    Parameters
+    ----------
+    model : :any:`SumModel`
+        The sum-model in use.
+    """
+    m, p = summod, summod._prec
+    ani_str, ang_str, o_str, r_str, p_str = "", "", "", "", ""
+    m_str = ", ".join([mod.name for mod in m.models])
+    t_str = ", temporal=True" if m.temporal else ""
+    d_str = f"latlon={m.latlon}" if m.latlon else f"dim={m.spatial_dim}"
+    if len(m) > 0:
+        m_str += ", "
+        p_str += f", vars={list_format(m.vars, p)}"
+        p_str += f", len_scales={list_format(m.len_scales, p)}"
+    p_str += "" if np.isclose(m.nugget, 0) else f", nugget={m.nugget:.{p}}"
+    for opt in m.opt_arg:
+        o_str += f", {opt}={getattr(m, opt):.{p}}"
+    if m.latlon:
+        if not m.is_isotropic and m.temporal:
+            ani_str = f", anis={m.anis[-1]:.{p}}"
+        if not np.isclose(m.geo_scale, 1):
+            r_str = f", geo_scale={m.geo_scale:.{p}}"
+    else:
+        if not m.is_isotropic:
+            ani_str = f", anis={list_format(m.anis, p)}"
+        if m.do_rotation:
+            ang_str = f", angles={list_format(m.angles, p)}"
+    return f"{m.name}({m_str}{d_str}{t_str}{p_str}{ani_str}{ang_str}{r_str}{o_str})"
diff --git a/src/gstools/covmodel/tools.py b/src/gstools/covmodel/tools.py
index a70c8118a..52a6e64c8 100644
--- a/src/gstools/covmodel/tools.py
+++ b/src/gstools/covmodel/tools.py
@@ -6,7 +6,6 @@
 The following classes and functions are provided
 
 .. autosummary::
-   RatioError
    AttributeWarning
    rad_fac
    set_opt_args
@@ -54,10 +53,6 @@
 ]
 
 
-class RatioError(Exception):
-    """Error for invalid ratios in SumModel."""
-
-
 class AttributeWarning(UserWarning):
     """Attribute warning for CovModel class."""
 

From 057c8884791059fdf06de95108f8c12d7c3ebd8a Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Sebastian=20M=C3=BCller?= <mueller.seb@posteo.de>
Date: Fri, 16 Aug 2024 00:04:21 +0200
Subject: [PATCH 20/58] add SumModel class

---
 src/gstools/covmodel/base.py | 503 +++++++++++++++++++++++++++++++++++
 1 file changed, 503 insertions(+)

diff --git a/src/gstools/covmodel/base.py b/src/gstools/covmodel/base.py
index a02a46516..8e694b26c 100644
--- a/src/gstools/covmodel/base.py
+++ b/src/gstools/covmodel/base.py
@@ -18,9 +18,19 @@
 
 from gstools.covmodel import plot
 from gstools.covmodel.fit import fit_variogram
+from gstools.covmodel.sum_tools import (
+    default_mod_kwargs,
+    sum_check,
+    sum_default_arg_bounds,
+    sum_default_opt_arg_bounds,
+    sum_model_repr,
+    sum_set_norm_len_ratios,
+    sum_set_norm_var_ratios,
+)
 from gstools.covmodel.tools import (
     _init_subclass,
     check_arg_bounds,
+    check_arg_in_bounds,
     compare,
     default_arg_from_bounds,
     model_repr,
@@ -1156,3 +1166,496 @@ def __setattr__(self, name, value):
     def __repr__(self):
         """Return String representation."""
         return model_repr(self)
+
+
+class SumModel(CovModel):
+    r"""Class for sums of covariance models.
+
+    Parameters
+    ----------
+    *models
+        tuple of :any:`CovModel` instances of subclasses to sum.
+    dim : :class:`int`, optional
+        dimension of the model.
+        Includes the temporal dimension if temporal is true.
+        To specify only the spatial dimension in that case, use `spatial_dim`.
+        Default: ``3``
+    vars : :class:`list` of :class:`float`, optional
+        variances of the models. Will reset present variances.
+    len_scales : :class:`list` of :class:`float`, optional
+        length scale of the models. Will reset present length scales.
+    nugget : :class:`float`, optional
+        nugget of the sum-model. All summed models will a nugget of zero.
+        Default: ``0.0``
+    anis : :class:`float` or :class:`list`, optional
+        anisotropy ratios in the transversal directions [e_y, e_z].
+
+            * e_y = l_y / l_x
+            * e_z = l_z / l_x
+
+        If only one value is given in 3D, e_y will be set to 1.
+        This value will be ignored, if multiple len_scales are given.
+        Default: ``1.0``
+    angles : :class:`float` or :class:`list`, optional
+        angles of rotation (given in rad):
+
+            * in 2D: given as rotation around z-axis
+            * in 3D: given by yaw, pitch, and roll (known as Tait–Bryan angles)
+
+         Default: ``0.0``
+    integral_scale : :class:`float` or :class:`list` or :any:`None`, optional
+        If given, ``len_scale`` will be ignored and recalculated,
+        so that the integral scale of the model matches the given one.
+        Default: :any:`None`
+    rescale : :class:`float` or :any:`None`, optional
+        Optional rescaling factor to divide the length scale with.
+        This could be used for unit conversion or rescaling the length scale
+        to coincide with e.g. the integral scale.
+        Will be set by each model individually.
+        Default: :any:`None`
+    latlon : :class:`bool`, optional
+        Whether the model is describing 2D fields on earths surface described
+        by latitude and longitude. When using this, the model will internally
+        use the associated 'Yadrenko' model to represent a valid model.
+        This means, the spatial distance :math:`r` will be replaced by
+        :math:`2\sin(\alpha/2)`, where :math:`\alpha` is the great-circle
+        distance, which is equal to the spatial distance of two points in 3D.
+        As a consequence, `dim` will be set to `3` and anisotropy will be
+        disabled. `geo_scale` can be set to e.g. earth's radius,
+        to have a meaningful `len_scale` parameter.
+        Default: False
+    geo_scale : :class:`float`, optional
+        Geographic unit scaling in case of latlon coordinates to get a
+        meaningful length scale unit.
+        By default, len_scale is assumed to be in radians with latlon=True.
+        Can be set to :any:`KM_SCALE` to have len_scale in km or
+        :any:`DEGREE_SCALE` to have len_scale in degrees.
+        Default: :any:`RADIAN_SCALE`
+    temporal : :class:`bool`, optional
+        Create a metric spatio-temporal covariance model.
+        Setting this to true will increase `dim` and `field_dim` by 1.
+        `spatial_dim` will be `field_dim - 1`.
+        The time-dimension is appended, meaning the pos tuple is (x,y,z,...,t).
+        Default: False
+    spatial_dim : :class:`int`, optional
+        spatial dimension of the model.
+        If given, the model dimension will be determined from this spatial dimension
+        and the possible temporal dimension if temporal is ture.
+        Default: None
+    var : :class:`float`, optional
+        Total variance of the sum-model. Will evenly scale present variances.
+    len_scale : :class:`float` or :class:`list`, optional
+        Total length scale of the sum-model. Will evenly scale present length scales.
+    **opt_arg
+        Optional arguments of the sub-models will have and added index of the sub-model.
+        Also covers ``var_<i>`` and ``length_scale_<i>``.
+    """
+
+    def __init__(self, *models, **kwargs):
+        self._init = False
+        self._models = []
+        add_nug = 0.0
+        for mod in models:
+            if isinstance(mod, type) and issubclass(mod, SumModel):
+                continue  # treat un-init sum-model as nugget model with 0 nugget
+            if isinstance(mod, SumModel):
+                self._models += mod.models
+                add_nug += mod.nugget
+            elif isinstance(mod, CovModel):
+                self._models.append(mod)
+            elif isinstance(mod, type) and issubclass(mod, CovModel):
+                self._models.append(mod(**default_mod_kwargs(kwargs)))
+            else:
+                msg = "SumModel: models need to be instances or subclasses of CovModel."
+                raise ValueError(msg)
+        # handle parameters when only Nugget models given
+        if models and not self.models:
+            for mod in models:
+                if not isinstance(mod, type):
+                    kwargs.setdefault("dim", mod.dim)
+                    kwargs.setdefault("latlon", mod.latlon)
+                    kwargs.setdefault("temporal", mod.temporal)
+                    kwargs.setdefault("geo_scale", mod.geo_scale)
+                    kwargs.setdefault("anis", mod.anis)
+                    kwargs.setdefault("angles", mod.angles)
+                    break
+        # pop entries that can't be re-set
+        self._latlon = bool(
+            kwargs.pop(
+                "latlon", self.models[0].latlon if self.models else False
+            )
+        )
+        self._temporal = bool(
+            kwargs.pop(
+                "temporal", self.models[0].temporal if self.models else False
+            )
+        )
+        self._geo_scale = float(
+            kwargs.pop(
+                "geo_scale",
+                self.models[0].geo_scale if self.models else RADIAN_SCALE,
+            )
+        )
+        var_set = kwargs.pop("var", None)
+        len_set = kwargs.pop("len_scale", None)
+        # convert nugget
+        self._nugget = float(
+            kwargs.pop(
+                "nugget",
+                sum((mod.nugget for mod in self.models), 0.0) + add_nug,
+            )
+        )
+        for mod in self.models:
+            mod._nugget = 0.0
+        # prepare dim setting
+        if "spatial_dim" in kwargs:
+            spatial_dim = kwargs.pop("spatial_dim")
+            if spatial_dim is not None:
+                kwargs["dim"] = spatial_dim + int(self.temporal)
+        self._dim = None
+        self._hankel_kw = {}
+        self._sft = None
+        self.dim = kwargs.get("dim", self.models[0].dim if self.models else 3)
+        # prepare parameters (they are checked in dim setting)
+        self._rescale = 1.0
+        anis = kwargs.get("anis", self.models[0].anis if self.models else 1)
+        angles = kwargs.get(
+            "angles", self.models[0].angles if self.models else 0
+        )
+        _, self._anis = set_len_anis(self.dim, 1.0, anis, self.latlon)
+        self._angles = set_model_angles(
+            self.dim, angles, self.latlon, self.temporal
+        )
+        # prepare parameters boundaries
+        self._var_bounds = None
+        self._len_scale_bounds = None
+        self._nugget_bounds = None
+        self._anis_bounds = None
+        self._opt_arg_bounds = {}
+        bounds = self.default_arg_bounds()
+        bounds.update(self.default_opt_arg_bounds())
+        self.set_arg_bounds(check_args=False, **bounds)
+        # opt arg determining
+        self._prec = 3
+        self._init = True
+        # set remaining args
+        for arg, val in kwargs.items():
+            setattr(self, arg, val)
+        if var_set is not None:
+            self.var = var_set
+        if len_set is not None:
+            self.len_scale = len_set
+        # check for consistency
+        self.check()
+
+    def __iter__(self):
+        return iter(self.models)
+
+    def __len__(self):
+        return len(self.models)
+
+    def __contains__(self, item):
+        return item in self.models
+
+    def __getitem__(self, key):
+        return self.models[key]
+
+    def check(self):
+        """Check consistency of contained models."""
+        sum_check(self)
+
+    def default_arg_bounds(self):
+        """Default boundaries for arguments as dict."""
+        return sum_default_arg_bounds(self)
+
+    def default_opt_arg_bounds(self):
+        """Defaults boundaries for optional arguments as dict."""
+        return sum_default_opt_arg_bounds(self)
+
+    def set_norm_var_ratios(self, ratios, skip=None, var=None):
+        """
+        Set variances of contained models by normalized ratios in [0, 1].
+
+        Ratios are given as normalized ratios in [0, 1] as relative ratio of
+        variance to remaining difference to total variance of the Sum-Model.
+
+        Parameters
+        ----------
+        ratios : iterable
+            Ratios to set. Should have a length of len(models) - len(exclude) - 1
+        skip : iterable, optional
+            Model indices to skip. Should have compatible lenth, by default None
+        var : float, optional
+            Desired variance, by default current variance
+
+        Raises
+        ------
+        ValueError
+            If number of ratios is not matching.
+        """
+        sum_set_norm_var_ratios(self, ratios, skip, var)
+
+    def set_norm_len_ratios(self, ratios, skip=None, len_scale=None):
+        """
+        Set length scales of contained models by normalized ratios in [0, 1].
+
+        Ratios are given as normalized ratios in [0, 1] as relative ratio of
+        len_scale * var / total_var to remaining difference to
+        total len_scale of the Sum-Model.
+
+        Parameters
+        ----------
+        ratios : iterable
+            Ratios to set. Should have a length of len(models) - len(exclude) - 1
+        skip : iterable, optional
+            Model indices to skip. Should have compatible lenth, by default None
+        len_scale : float, optional
+            Desired len_scale, by default current len_scale
+
+        Raises
+        ------
+        ValueError
+            If number of ratios is not matching.
+        """
+        sum_set_norm_len_ratios(self, ratios, skip, len_scale)
+
+    @property
+    def hankel_kw(self):
+        """:class:`dict`: :any:`hankel.SymmetricFourierTransform` kwargs."""
+        return self._hankel_kw
+
+    @property
+    def models(self):
+        """:class:`tuple`: The summed models."""
+        return self._models
+
+    @property
+    def geo_scale(self):
+        """:class:`float`: Geographic scaling for geographical coords."""
+        return self._geo_scale
+
+    @geo_scale.setter
+    def geo_scale(self, geo_scale):
+        self._geo_scale = abs(float(geo_scale))
+        for mod in self.models:
+            mod.geo_scale = geo_scale
+
+    @property
+    def dim(self):
+        """:class:`int`: The dimension of the model."""
+        return self._dim
+
+    @dim.setter
+    def dim(self, dim):
+        set_dim(self, dim)
+        for mod in self.models:
+            mod.dim = dim
+
+    @property
+    def var(self):
+        """:class:`float`: The variance of the model."""
+        return sum((var for var in self.vars), 0.0)
+
+    @var.setter
+    def var(self, var):
+        for mod, rat in zip(self.models, self.ratios):
+            mod.var = rat * var
+        if not self.models and not np.isclose(var, 0):
+            msg = f"{self.name}: variance needs to be 0."
+            raise ValueError(msg)
+        check_arg_in_bounds(self, "var", error=True)
+        check_arg_in_bounds(self, "len_scale", error=True)
+
+    @property
+    def vars(self):
+        """:class:`list`: The variances of the models."""
+        return [mod.var for mod in self.models]
+
+    @vars.setter
+    def vars(self, vars):
+        if len(vars) != len(self):
+            msg = "SumModel: number of given variances not matching"
+            raise ValueError(msg)
+        for mod, var in zip(self.models, vars):
+            mod.var = var
+        check_arg_in_bounds(self, "var", error=True)
+        check_arg_in_bounds(self, "len_scale", error=True)
+
+    @property
+    def len_scale(self):
+        """:class:`float`: The main length scale of the model."""
+        return sum(
+            (
+                mod.len_scale * rat
+                for mod, rat in zip(self.models, self.ratios)
+            ),
+            0.0,
+        )
+
+    @len_scale.setter
+    def len_scale(self, len_scale):
+        len_scale, anis = set_len_anis(
+            self.dim, len_scale, self.anis, self.latlon
+        )
+        old_scale = self.len_scale
+        self.anis = anis
+        for mod in self.models:
+            mod.len_scale = mod.len_scale * len_scale / old_scale
+        if not self.models and not np.isclose(len_scale, 0):
+            msg = f"{self.name}: length scale needs to be 0."
+            raise ValueError(msg)
+        check_arg_in_bounds(self, "len_scale", error=True)
+
+    @property
+    def len_scales(self):
+        """:class:`list`: The main length scales of the models."""
+        return [mod.len_scale for mod in self.models]
+
+    @len_scales.setter
+    def len_scales(self, len_scales):
+        if len(len_scales) != len(self):
+            msg = "SumModel: number of given length scales not matching"
+            raise ValueError(msg)
+        for mod, len_scale in zip(self.models, len_scales):
+            mod.len_scale = len_scale
+        check_arg_in_bounds(self, "len_scale", error=True)
+
+    @property
+    def anis(self):
+        """:class:`numpy.ndarray`: The anisotropy factors of the model."""
+        return self._anis
+
+    @anis.setter
+    def anis(self, anis):
+        _, self._anis = set_len_anis(self.dim, 1.0, anis, self.latlon)
+        for mod in self.models:
+            mod.anis = anis
+        check_arg_in_bounds(self, "anis", error=True)
+
+    @property
+    def angles(self):
+        """:class:`numpy.ndarray`: Rotation angles (in rad) of the model."""
+        return self._angles
+
+    @angles.setter
+    def angles(self, angles):
+        self._angles = set_model_angles(
+            self.dim, angles, self.latlon, self.temporal
+        )
+        for mod in self.models:
+            mod.angles = angles
+
+    @property
+    def ratios(self):
+        var = self.var
+        if np.isclose(var, 0) and len(self) > 0:
+            return np.full(len(self), 1 / len(self))
+        return np.array([mod.var / var for mod in self.models])
+
+    @ratios.setter
+    def ratios(self, ratios):
+        if len(ratios) != len(self):
+            msg = "SumModel.ratios: wrong number of given ratios."
+            raise ValueError(msg)
+        if ratios and not np.isclose(np.sum(ratios), 1):
+            msg = "SumModel.ratios: given ratios do not sum to 1."
+            raise ValueError(msg)
+        var = self.var
+        for mod, rat in zip(self.models, ratios):
+            mod.var = var * rat
+        check_arg_in_bounds(self, "var", error=True)
+        check_arg_in_bounds(self, "len_scale", error=True)
+
+    def calc_integral_scale(self):
+        return sum(
+            (
+                mod.integral_scale * rat
+                for mod, rat in zip(self.models, self.ratios)
+            ),
+            0.0,
+        )
+
+    def spectral_density(self, k):
+        return sum(
+            (
+                mod.spectral_density(k) * rat
+                for mod, rat in zip(self.models, self.ratios)
+            ),
+            np.zeros_like(k, dtype=np.double),
+        )
+
+    def correlation(self, r):
+        return sum(
+            (
+                mod.correlation(r) * rat
+                for mod, rat in zip(self.models, self.ratios)
+            ),
+            np.zeros_like(r, dtype=np.double),
+        )
+
+    @property
+    def opt_arg(self):
+        """:class:`list` of :class:`str`: Names of the optional arguments."""
+        return sum(
+            [
+                [f"{opt}_{i}" for opt in mod.opt_arg]
+                for i, mod in enumerate(self.models)
+            ],
+            [],
+        )
+
+    @property
+    def sub_arg(self):
+        """:class:`list` of :class:`str`: Names of the sub-arguments for var and len_scale."""
+        return sum(
+            [
+                [f"{o}_{i}" for o in ["var", "len_scale"]]
+                for i, mod in enumerate(self.models)
+            ],
+            [],
+        )
+
+    def __setattr__(self, name, value):
+        """Set an attribute."""
+        sub_arg = False
+        if getattr(self, "_init", False):
+            for i, mod in enumerate(self.models):
+                if name == f"var_{i}":
+                    mod.var = value
+                    sub_arg = True
+                if name == f"len_scale_{i}":
+                    mod.len_scale = value
+                    sub_arg = True
+                if name == f"integral_scale_{i}":
+                    mod.integral_scale = value
+                    sub_arg = True
+                for opt in mod.opt_arg:
+                    if name == f"{opt}_{i}":
+                        setattr(mod, opt, value)
+                        sub_arg = True
+                if sub_arg:
+                    break
+        if sub_arg:
+            self.check_arg_bounds()
+        else:
+            super().__setattr__(name, value)
+
+    def __getattr__(self, name):
+        """Get an attribute."""
+        # __getattr__ is only called when an attribute is not found in the usual places
+        if name != "_init" and getattr(self, "_init", False):
+            for i, mod in enumerate(self.models):
+                if name == f"var_{i}":
+                    return mod.var
+                if name == f"len_scale_{i}":
+                    return mod.len_scale
+                if name == f"integral_scale_{i}":
+                    return mod.integral_scale
+                for opt in mod.opt_arg:
+                    if name == f"{opt}_{i}":
+                        return getattr(mod, opt)
+        raise AttributeError(f"'{self.name}' object has no attribute '{name}'")
+
+    def __repr__(self):
+        """Return String representation."""
+        return sum_model_repr(self)

From 656125a4c37345b7bd4cea12639e52bb538f31ad Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Sebastian=20M=C3=BCller?= <mueller.seb@posteo.de>
Date: Fri, 16 Aug 2024 00:09:43 +0200
Subject: [PATCH 21/58] add pure Nugget model

---
 src/gstools/__init__.py          |  6 +++
 src/gstools/covmodel/__init__.py |  7 +++-
 src/gstools/covmodel/models.py   | 67 +++++++++++++++++++++++++++++++-
 3 files changed, 78 insertions(+), 2 deletions(-)

diff --git a/src/gstools/__init__.py b/src/gstools/__init__.py
index 11e63a2b3..7488dc9ce 100644
--- a/src/gstools/__init__.py
+++ b/src/gstools/__init__.py
@@ -54,6 +54,7 @@
 
 .. autosummary::
    CovModel
+   SumModel
 
 Covariance Models
 ^^^^^^^^^^^^^^^^^
@@ -62,6 +63,7 @@
 ~~~~~~~~~~~~~~~~~~~~~~~~~~
 
 .. autosummary::
+   Nugget
    Gaussian
    Exponential
    Matern
@@ -153,9 +155,11 @@
     JBessel,
     Linear,
     Matern,
+    Nugget,
     Rational,
     Spherical,
     Stable,
+    SumModel,
     SuperSpherical,
     TPLExponential,
     TPLGaussian,
@@ -198,6 +202,8 @@
 __all__ += ["transform", "normalizer", "config"]
 __all__ += [
     "CovModel",
+    "SumModel",
+    "Nugget",
     "Gaussian",
     "Exponential",
     "Matern",
diff --git a/src/gstools/covmodel/__init__.py b/src/gstools/covmodel/__init__.py
index 28ab81f21..76920c61e 100644
--- a/src/gstools/covmodel/__init__.py
+++ b/src/gstools/covmodel/__init__.py
@@ -19,6 +19,7 @@
    :toctree:
 
    CovModel
+   SumModel
 
 Covariance Models
 ^^^^^^^^^^^^^^^^^
@@ -27,6 +28,7 @@
 .. autosummary::
    :toctree:
 
+   Nugget
    Gaussian
    Exponential
    Matern
@@ -53,7 +55,7 @@
    TPLSimple
 """
 
-from gstools.covmodel.base import CovModel
+from gstools.covmodel.base import CovModel, SumModel
 from gstools.covmodel.models import (
     Circular,
     Cubic,
@@ -64,6 +66,7 @@
     JBessel,
     Linear,
     Matern,
+    Nugget,
     Rational,
     Spherical,
     Stable,
@@ -78,6 +81,8 @@
 
 __all__ = [
     "CovModel",
+    "SumModel",
+    "Nugget",
     "Gaussian",
     "Exponential",
     "Matern",
diff --git a/src/gstools/covmodel/models.py b/src/gstools/covmodel/models.py
index b1a9d68ec..80dbd705a 100644
--- a/src/gstools/covmodel/models.py
+++ b/src/gstools/covmodel/models.py
@@ -6,6 +6,7 @@
 The following classes are provided
 
 .. autosummary::
+   Nugget
    Gaussian
    Exponential
    Matern
@@ -27,11 +28,12 @@
 import numpy as np
 from scipy import special as sps
 
-from gstools.covmodel.base import CovModel
+from gstools.covmodel.base import CovModel, SumModel
 from gstools.covmodel.tools import AttributeWarning
 from gstools.tools.special import exp_int, inc_gamma_low
 
 __all__ = [
+    "Nugget",
     "Gaussian",
     "Exponential",
     "Matern",
@@ -48,6 +50,69 @@
 ]
 
 
+class Nugget(SumModel):
+    r"""Pure nugget model.
+
+    Parameters
+    ----------
+    dim : :class:`int`, optional
+        dimension of the model.
+        Includes the temporal dimension if temporal is true.
+        To specify only the spatial dimension in that case, use `spatial_dim`.
+        Default: ``3``
+    nugget : :class:`float`, optional
+        nugget of the model. Default: ``0.0``
+    anis : :class:`float` or :class:`list`, optional
+        anisotropy ratios in the transversal directions [e_y, e_z].
+
+            * e_y = l_y / l_x
+            * e_z = l_z / l_x
+
+        If only one value is given in 3D, e_y will be set to 1.
+        This value will be ignored, if multiple len_scales are given.
+        Default: ``1.0``
+    angles : :class:`float` or :class:`list`, optional
+        angles of rotation (given in rad):
+
+            * in 2D: given as rotation around z-axis
+            * in 3D: given by yaw, pitch, and roll (known as Tait–Bryan angles)
+
+         Default: ``0.0``
+    latlon : :class:`bool`, optional
+        Whether the model is describing 2D fields on earths surface described
+        by latitude and longitude. When using this, the model will internally
+        use the associated 'Yadrenko' model to represent a valid model.
+        This means, the spatial distance :math:`r` will be replaced by
+        :math:`2\sin(\alpha/2)`, where :math:`\alpha` is the great-circle
+        distance, which is equal to the spatial distance of two points in 3D.
+        As a consequence, `dim` will be set to `3` and anisotropy will be
+        disabled. `geo_scale` can be set to e.g. earth's radius,
+        to have a meaningful `len_scale` parameter.
+        Default: False
+    geo_scale : :class:`float`, optional
+        Geographic unit scaling in case of latlon coordinates to get a
+        meaningful length scale unit.
+        By default, len_scale is assumed to be in radians with latlon=True.
+        Can be set to :any:`KM_SCALE` to have len_scale in km or
+        :any:`DEGREE_SCALE` to have len_scale in degrees.
+        Default: :any:`RADIAN_SCALE`
+    temporal : :class:`bool`, optional
+        Create a metric spatio-temporal covariance model.
+        Setting this to true will increase `dim` and `field_dim` by 1.
+        `spatial_dim` will be `field_dim - 1`.
+        The time-dimension is appended, meaning the pos tuple is (x,y,z,...,t).
+        Default: False
+    spatial_dim : :class:`int`, optional
+        spatial dimension of the model.
+        If given, the model dimension will be determined from this spatial dimension
+        and the possible temporal dimension if temporal is ture.
+        Default: None
+    """
+
+    def __init__(self, **kwargs):
+        super().__init__(**kwargs)
+
+
 class Gaussian(CovModel):
     r"""The Gaussian covariance model.
 

From dbdfa27fabcf5eb072ac0c43085c77d6bd008d74 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Sebastian=20M=C3=BCller?= <mueller.seb@posteo.de>
Date: Fri, 16 Aug 2024 00:26:32 +0200
Subject: [PATCH 22/58] CovModel: let the sum magic happen

---
 src/gstools/covmodel/base.py | 8 ++++++++
 1 file changed, 8 insertions(+)

diff --git a/src/gstools/covmodel/base.py b/src/gstools/covmodel/base.py
index 8e694b26c..68af58d02 100644
--- a/src/gstools/covmodel/base.py
+++ b/src/gstools/covmodel/base.py
@@ -1163,6 +1163,14 @@ def __setattr__(self, name, value):
         if getattr(self, "_init", False) and name in self.opt_arg:
             self.check_arg_bounds()
 
+    def __add__(self, other):
+        """Add two covariance models and return a SumModel."""
+        return SumModel(self, other)
+
+    def __radd__(self, other):
+        """Add two covariance models and return a SumModel."""
+        return SumModel(self, other)
+
     def __repr__(self):
         """Return String representation."""
         return model_repr(self)

From 2858f063665ffb2b1d7b5a921ac25609f47c7835 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Sebastian=20M=C3=BCller?= <mueller.seb@posteo.de>
Date: Fri, 16 Aug 2024 00:35:53 +0200
Subject: [PATCH 23/58] pylint fixes

---
 src/gstools/covmodel/base.py      | 10 ++++++----
 src/gstools/covmodel/models.py    |  2 +-
 src/gstools/covmodel/sum_tools.py |  9 +++++----
 src/gstools/covmodel/tools.py     |  1 -
 4 files changed, 12 insertions(+), 10 deletions(-)

diff --git a/src/gstools/covmodel/base.py b/src/gstools/covmodel/base.py
index 68af58d02..a495d80b5 100644
--- a/src/gstools/covmodel/base.py
+++ b/src/gstools/covmodel/base.py
@@ -9,7 +9,7 @@
    CovModel
 """
 
-# pylint: disable=C0103, R0201, E1101, C0302, W0613
+# pylint: disable=C0103, R0201, E1101, C0302, W0613, W0231
 import copy
 
 import numpy as np
@@ -1480,11 +1480,11 @@ def vars(self):
         return [mod.var for mod in self.models]
 
     @vars.setter
-    def vars(self, vars):
-        if len(vars) != len(self):
+    def vars(self, variances):
+        if len(variances) != len(self):
             msg = "SumModel: number of given variances not matching"
             raise ValueError(msg)
-        for mod, var in zip(self.models, vars):
+        for mod, var in zip(self.models, variances):
             mod.var = var
         check_arg_in_bounds(self, "var", error=True)
         check_arg_in_bounds(self, "len_scale", error=True)
@@ -1555,6 +1555,7 @@ def angles(self, angles):
 
     @property
     def ratios(self):
+        """:class:`numpy.ndarray`: Variance ratios of the sub-models."""
         var = self.var
         if np.isclose(var, 0) and len(self) > 0:
             return np.full(len(self), 1 / len(self))
@@ -1593,6 +1594,7 @@ def spectral_density(self, k):
         )
 
     def correlation(self, r):
+        """SumModel correlation function."""
         return sum(
             (
                 mod.correlation(r) * rat
diff --git a/src/gstools/covmodel/models.py b/src/gstools/covmodel/models.py
index 80dbd705a..7c1cd836d 100644
--- a/src/gstools/covmodel/models.py
+++ b/src/gstools/covmodel/models.py
@@ -22,7 +22,7 @@
    JBessel
 """
 
-# pylint: disable=C0103, E1101, R0201
+# pylint: disable=C0103, C0302, E1101, R0201
 import warnings
 
 import numpy as np
diff --git a/src/gstools/covmodel/sum_tools.py b/src/gstools/covmodel/sum_tools.py
index a4a1647a3..c4b9a90ee 100644
--- a/src/gstools/covmodel/sum_tools.py
+++ b/src/gstools/covmodel/sum_tools.py
@@ -17,6 +17,7 @@
    sum_model_repr
 """
 
+# pylint: disable=W0212
 import numpy as np
 
 from gstools.covmodel.tools import check_arg_in_bounds
@@ -89,10 +90,10 @@ def sum_default_arg_bounds(summod):
     """Default boundaries for arguments as dict."""
     var_bnds = [mod.var_bounds for mod in summod.models]
     len_bnds = [mod.len_scale_bounds for mod in summod.models]
-    var_lo = sum([bnd[0] for bnd in var_bnds], start=0.0)
-    var_hi = sum([bnd[1] for bnd in var_bnds], start=0.0)
-    len_lo = min([bnd[0] for bnd in len_bnds], default=0.0)
-    len_hi = max([bnd[1] for bnd in len_bnds], default=0.0)
+    var_lo = sum((bnd[0] for bnd in var_bnds), start=0.0)
+    var_hi = sum((bnd[1] for bnd in var_bnds), start=0.0)
+    len_lo = min((bnd[0] for bnd in len_bnds), default=0.0)
+    len_hi = max((bnd[1] for bnd in len_bnds), default=0.0)
     res = {
         "var": (var_lo, var_hi),
         "len_scale": (len_lo, len_hi),
diff --git a/src/gstools/covmodel/tools.py b/src/gstools/covmodel/tools.py
index 52a6e64c8..69bf963b1 100644
--- a/src/gstools/covmodel/tools.py
+++ b/src/gstools/covmodel/tools.py
@@ -34,7 +34,6 @@
 from gstools.tools.misc import list_format
 
 __all__ = [
-    "RatioError",
     "AttributeWarning",
     "rad_fac",
     "set_opt_args",

From 947619a0428e716e7b3ee78952b631be3c24824d Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Sebastian=20M=C3=BCller?= <mueller.seb@posteo.de>
Date: Fri, 16 Aug 2024 00:40:47 +0200
Subject: [PATCH 24/58] fix doc-string of SumModel for docs

---
 src/gstools/covmodel/base.py | 12 ++++++++++++
 1 file changed, 12 insertions(+)

diff --git a/src/gstools/covmodel/base.py b/src/gstools/covmodel/base.py
index a495d80b5..252bdeb30 100644
--- a/src/gstools/covmodel/base.py
+++ b/src/gstools/covmodel/base.py
@@ -1356,6 +1356,18 @@ def __init__(self, *models, **kwargs):
         # check for consistency
         self.check()
 
+    def __init_subclass__(cls):
+        """Initialize gstools sum-model."""
+        _init_subclass(cls)
+        # overridden functions get standard doc if no new doc was created
+        ign = ["__", "variogram", "covariance", "cor"]
+        for att, attr_cls in cls.__dict__.items():
+            if any(att.startswith(i) for i in ign) or att not in dir(CovModel):
+                continue
+            attr_doc = getattr(CovModel, att).__doc__
+            if attr_cls.__doc__ is None:
+                attr_cls.__doc__ = attr_doc
+
     def __iter__(self):
         return iter(self.models)
 

From 480fa32e8d6d227ac9e8039e4dabc222951adea0 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Sebastian=20M=C3=BCller?= <mueller.seb@posteo.de>
Date: Fri, 16 Aug 2024 00:57:54 +0200
Subject: [PATCH 25/58] CovModel: add switch to add doc string

---
 src/gstools/covmodel/base.py | 19 ++++++-------------
 1 file changed, 6 insertions(+), 13 deletions(-)

diff --git a/src/gstools/covmodel/base.py b/src/gstools/covmodel/base.py
index 252bdeb30..e080ba20b 100644
--- a/src/gstools/covmodel/base.py
+++ b/src/gstools/covmodel/base.py
@@ -147,6 +147,8 @@ class CovModel:
         If present, they are described in the section `Other Parameters`.
     """
 
+    _add_doc = True
+
     def __init__(
         self,
         dim=3,
@@ -242,7 +244,8 @@ def __init_subclass__(cls):
         # modify the docstrings: class docstring gets attributes added
         if cls.__doc__ is None:
             cls.__doc__ = "User defined GSTools Covariance-Model."
-        cls.__doc__ += CovModel.__doc__[45:]
+        if cls._add_doc:
+            cls.__doc__ += CovModel.__doc__[45:]
         # overridden functions get standard doc if no new doc was created
         ign = ["__", "variogram", "covariance", "cor"]
         for att, attr_cls in cls.__dict__.items():
@@ -1259,6 +1262,8 @@ class SumModel(CovModel):
         Also covers ``var_<i>`` and ``length_scale_<i>``.
     """
 
+    _add_doc = False
+
     def __init__(self, *models, **kwargs):
         self._init = False
         self._models = []
@@ -1356,18 +1361,6 @@ def __init__(self, *models, **kwargs):
         # check for consistency
         self.check()
 
-    def __init_subclass__(cls):
-        """Initialize gstools sum-model."""
-        _init_subclass(cls)
-        # overridden functions get standard doc if no new doc was created
-        ign = ["__", "variogram", "covariance", "cor"]
-        for att, attr_cls in cls.__dict__.items():
-            if any(att.startswith(i) for i in ign) or att not in dir(CovModel):
-                continue
-            attr_doc = getattr(CovModel, att).__doc__
-            if attr_cls.__doc__ is None:
-                attr_cls.__doc__ = attr_doc
-
     def __iter__(self):
         return iter(self.models)
 

From 171efbad0cbb1f55523be315b6db053ba8414a44 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Sebastian=20M=C3=BCller?= <mueller.seb@posteo.de>
Date: Fri, 16 Aug 2024 00:58:09 +0200
Subject: [PATCH 26/58] Fourier: fix doc

---
 src/gstools/field/generator.py | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/src/gstools/field/generator.py b/src/gstools/field/generator.py
index 89914785a..83d138165 100755
--- a/src/gstools/field/generator.py
+++ b/src/gstools/field/generator.py
@@ -736,7 +736,7 @@ def update(self, model=None, seed=np.nan, period=None, mode_no=None):
             the seed of the random number generator.
             If :any:`None`, a random seed is used. If :any:`numpy.nan`,
             the actual seed will be kept. Default: :any:`numpy.nan`
-        period : :class:`list` or :any:`None, optional
+        period : :class:`list` or :any:`None`, optional
             The spatial periodicity of the field, is often the domain size.
         mode_no : :class:`list` or :any:`None`, optional
             Number of Fourier modes per dimension.

From e7838cf023ceff961054426e7b6eb0a169d32601 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Sebastian=20M=C3=BCller?= <mueller.seb@posteo.de>
Date: Mon, 19 Aug 2024 11:30:24 +0200
Subject: [PATCH 27/58] SumModel: models need to either be all instances or all
 subclasses of CovModel

---
 src/gstools/covmodel/base.py | 57 +++++++++++++++++++++++++++++-------
 1 file changed, 46 insertions(+), 11 deletions(-)

diff --git a/src/gstools/covmodel/base.py b/src/gstools/covmodel/base.py
index e080ba20b..dd7ab2ee0 100644
--- a/src/gstools/covmodel/base.py
+++ b/src/gstools/covmodel/base.py
@@ -1182,15 +1182,32 @@ def __repr__(self):
 class SumModel(CovModel):
     r"""Class for sums of covariance models.
 
+    This class represents sums of covariance models. The nugget of
+    each contained model will be set to zero and the sum model will
+    have its own nugget.
+    The variance of the sum model is the sum of the sub model variances
+    and the length scale of the sum model is the variance-weighted sum
+    of the length scales of the sub models. This is motivated by the fact,
+    that the integral scale of the sum model is equal to the variance-weighted
+    sum of the integral scales of the sub models.
+    An empty sum represents a pure Nugget model.
+    Resetting the total variance or the total length scale will evenly
+    scale the variances or length scales of the sub models.
+
     Parameters
     ----------
     *models
         tuple of :any:`CovModel` instances of subclasses to sum.
+        All models will get a nugget of zero and the nugget will
+        be set in the SumModel directly.
+        Models need to have matching temporal setting,
+        latlon setting, anis, angles and geo_scale.
+        The model dimension will be specified by the first given model.
     dim : :class:`int`, optional
         dimension of the model.
         Includes the temporal dimension if temporal is true.
         To specify only the spatial dimension in that case, use `spatial_dim`.
-        Default: ``3``
+        Default: ``3`` or dimension of the first given model (if instance).
     vars : :class:`list` of :class:`float`, optional
         variances of the models. Will reset present variances.
     len_scales : :class:`list` of :class:`float`, optional
@@ -1206,14 +1223,14 @@ class SumModel(CovModel):
 
         If only one value is given in 3D, e_y will be set to 1.
         This value will be ignored, if multiple len_scales are given.
-        Default: ``1.0``
+        Default: ``1.0`` or anis of the first given model (if instance).
     angles : :class:`float` or :class:`list`, optional
         angles of rotation (given in rad):
 
             * in 2D: given as rotation around z-axis
             * in 3D: given by yaw, pitch, and roll (known as Tait–Bryan angles)
 
-         Default: ``0.0``
+         Default: ``0.0`` or angles of the first given model (if instance).
     integral_scale : :class:`float` or :class:`list` or :any:`None`, optional
         If given, ``len_scale`` will be ignored and recalculated,
         so that the integral scale of the model matches the given one.
@@ -1234,20 +1251,20 @@ class SumModel(CovModel):
         As a consequence, `dim` will be set to `3` and anisotropy will be
         disabled. `geo_scale` can be set to e.g. earth's radius,
         to have a meaningful `len_scale` parameter.
-        Default: False
+        Default: False or latlon config of the first given model (if instance).
     geo_scale : :class:`float`, optional
         Geographic unit scaling in case of latlon coordinates to get a
         meaningful length scale unit.
         By default, len_scale is assumed to be in radians with latlon=True.
         Can be set to :any:`KM_SCALE` to have len_scale in km or
         :any:`DEGREE_SCALE` to have len_scale in degrees.
-        Default: :any:`RADIAN_SCALE`
+        Default: :any:`RADIAN_SCALE` or geo_scale of the first given model (if instance).
     temporal : :class:`bool`, optional
         Create a metric spatio-temporal covariance model.
         Setting this to true will increase `dim` and `field_dim` by 1.
         `spatial_dim` will be `field_dim - 1`.
         The time-dimension is appended, meaning the pos tuple is (x,y,z,...,t).
-        Default: False
+        Default: False or temporal config of the first given model (if instance).
     spatial_dim : :class:`int`, optional
         spatial dimension of the model.
         If given, the model dimension will be determined from this spatial dimension
@@ -1259,7 +1276,8 @@ class SumModel(CovModel):
         Total length scale of the sum-model. Will evenly scale present length scales.
     **opt_arg
         Optional arguments of the sub-models will have and added index of the sub-model.
-        Also covers ``var_<i>`` and ``length_scale_<i>``.
+        Also covers ``var_<i>`` and ``length_scale_<i>`` but they should preferably be
+        set by ``vars`` and ``length_scales``.
     """
 
     _add_doc = False
@@ -1268,15 +1286,31 @@ def __init__(self, *models, **kwargs):
         self._init = False
         self._models = []
         add_nug = 0.0
+        to_init = None
+        imsg = (
+            "SumModel: either all models are CovModel instances or subclasses."
+        )
         for mod in models:
             if isinstance(mod, type) and issubclass(mod, SumModel):
+                if to_init is not None and not to_init:
+                    raise ValueError(imsg)
+                to_init = True
                 continue  # treat un-init sum-model as nugget model with 0 nugget
             if isinstance(mod, SumModel):
+                if to_init is not None and to_init:
+                    raise ValueError(imsg)
+                to_init = False
                 self._models += mod.models
                 add_nug += mod.nugget
             elif isinstance(mod, CovModel):
+                if to_init is not None and to_init:
+                    raise ValueError(imsg)
+                to_init = False
                 self._models.append(mod)
             elif isinstance(mod, type) and issubclass(mod, CovModel):
+                if to_init is not None and not to_init:
+                    raise ValueError(imsg)
+                to_init = True
                 self._models.append(mod(**default_mod_kwargs(kwargs)))
             else:
                 msg = "SumModel: models need to be instances or subclasses of CovModel."
@@ -1315,7 +1349,7 @@ def __init__(self, *models, **kwargs):
         self._nugget = float(
             kwargs.pop(
                 "nugget",
-                sum((mod.nugget for mod in self.models), 0.0) + add_nug,
+                sum((mod.nugget for mod in self.models), 0) + add_nug,
             )
         )
         for mod in self.models:
@@ -1339,7 +1373,7 @@ def __init__(self, *models, **kwargs):
         self._angles = set_model_angles(
             self.dim, angles, self.latlon, self.temporal
         )
-        # prepare parameters boundaries
+        # prepare parameter boundaries
         self._var_bounds = None
         self._len_scale_bounds = None
         self._nugget_bounds = None
@@ -1348,17 +1382,18 @@ def __init__(self, *models, **kwargs):
         bounds = self.default_arg_bounds()
         bounds.update(self.default_opt_arg_bounds())
         self.set_arg_bounds(check_args=False, **bounds)
-        # opt arg determining
+        # finalize init
         self._prec = 3
         self._init = True
         # set remaining args
         for arg, val in kwargs.items():
             setattr(self, arg, val)
+        # reset total variance and length scale last
         if var_set is not None:
             self.var = var_set
         if len_set is not None:
             self.len_scale = len_set
-        # check for consistency
+        # check consistency of sub models
         self.check()
 
     def __iter__(self):

From 0e35a222ef110af52f59ec2048b0be0d4c3cf2de Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Sebastian=20M=C3=BCller?= <mueller.seb@posteo.de>
Date: Mon, 19 Aug 2024 11:40:08 +0200
Subject: [PATCH 28/58] SumModel: sum models have a constant rescale factor of
 1 since they are derived from correlation not cor

---
 src/gstools/covmodel/base.py | 13 ++++++-------
 1 file changed, 6 insertions(+), 7 deletions(-)

diff --git a/src/gstools/covmodel/base.py b/src/gstools/covmodel/base.py
index dd7ab2ee0..38a474912 100644
--- a/src/gstools/covmodel/base.py
+++ b/src/gstools/covmodel/base.py
@@ -1193,6 +1193,7 @@ class SumModel(CovModel):
     An empty sum represents a pure Nugget model.
     Resetting the total variance or the total length scale will evenly
     scale the variances or length scales of the sub models.
+    Sum models will have a constant rescale factor of one.
 
     Parameters
     ----------
@@ -1235,12 +1236,6 @@ class SumModel(CovModel):
         If given, ``len_scale`` will be ignored and recalculated,
         so that the integral scale of the model matches the given one.
         Default: :any:`None`
-    rescale : :class:`float` or :any:`None`, optional
-        Optional rescaling factor to divide the length scale with.
-        This could be used for unit conversion or rescaling the length scale
-        to coincide with e.g. the integral scale.
-        Will be set by each model individually.
-        Default: :any:`None`
     latlon : :class:`bool`, optional
         Whether the model is describing 2D fields on earths surface described
         by latitude and longitude. When using this, the model will internally
@@ -1364,7 +1359,6 @@ def __init__(self, *models, **kwargs):
         self._sft = None
         self.dim = kwargs.get("dim", self.models[0].dim if self.models else 3)
         # prepare parameters (they are checked in dim setting)
-        self._rescale = 1.0
         anis = kwargs.get("anis", self.models[0].anis if self.models else 1)
         angles = kwargs.get(
             "angles", self.models[0].angles if self.models else 0
@@ -1477,6 +1471,11 @@ def models(self):
         """:class:`tuple`: The summed models."""
         return self._models
 
+    @property
+    def rescale(self):
+        """:class:`float`: SumModel has a constant rescale factor of one."""
+        return 1.0
+
     @property
     def geo_scale(self):
         """:class:`float`: Geographic scaling for geographical coords."""

From 2ded2e0c47a4b107047caf2ad6ecbf88cf269938 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Sebastian=20M=C3=BCller?= <mueller.seb@posteo.de>
Date: Wed, 1 Jan 2025 12:48:45 +0100
Subject: [PATCH 29/58] CovModel: add 'needs_fourier_transform' attribute

---
 src/gstools/covmodel/base.py  | 35 ++++++++++++++++++++++-------------
 src/gstools/covmodel/tools.py |  3 ++-
 2 files changed, 24 insertions(+), 14 deletions(-)

diff --git a/src/gstools/covmodel/base.py b/src/gstools/covmodel/base.py
index 38a474912..be927ca03 100644
--- a/src/gstools/covmodel/base.py
+++ b/src/gstools/covmodel/base.py
@@ -7,6 +7,7 @@
 
 .. autosummary::
    CovModel
+   SumModel
 """
 
 # pylint: disable=C0103, R0201, E1101, C0302, W0613, W0231
@@ -495,6 +496,13 @@ def percentile_scale(self, per=0.9):
 
     # spectrum methods (can be overridden for speedup)
 
+    @property
+    def needs_fourier_transform(self):
+        """Whether the model doesn't implement an analytical spectral density."""
+        base_method = getattr(CovModel, "spectral_density")
+        instance_method = getattr(self.__class__, "spectral_density")
+        return base_method == instance_method
+
     def spectrum(self, k):
         r"""
         Spectrum of the covariance model.
@@ -1030,12 +1038,13 @@ def hankel_kw(self):
 
     @hankel_kw.setter
     def hankel_kw(self, hankel_kw):
-        if self._hankel_kw is None or hankel_kw is None:
-            self._hankel_kw = copy.copy(HANKEL_DEFAULT)
-        if hankel_kw is not None:
-            self._hankel_kw.update(hankel_kw)
-        if self.dim is not None:
-            self._sft = SFT(ndim=self.dim, **self.hankel_kw)
+        if self.needs_fourier_transform:
+            if self._hankel_kw is None or hankel_kw is None:
+                self._hankel_kw = copy.copy(HANKEL_DEFAULT)
+            if hankel_kw is not None:
+                self._hankel_kw.update(hankel_kw)
+            if self.dim is not None:
+                self._sft = SFT(ndim=self.dim, **self.hankel_kw)
 
     @property
     def dist_func(self):
@@ -1214,7 +1223,7 @@ class SumModel(CovModel):
     len_scales : :class:`list` of :class:`float`, optional
         length scale of the models. Will reset present length scales.
     nugget : :class:`float`, optional
-        nugget of the sum-model. All summed models will a nugget of zero.
+        nugget of the sum-model. All summed models will have a nugget of zero.
         Default: ``0.0``
     anis : :class:`float` or :class:`list`, optional
         anisotropy ratios in the transversal directions [e_y, e_z].
@@ -1355,7 +1364,7 @@ def __init__(self, *models, **kwargs):
             if spatial_dim is not None:
                 kwargs["dim"] = spatial_dim + int(self.temporal)
         self._dim = None
-        self._hankel_kw = {}
+        self._hankel_kw = None
         self._sft = None
         self.dim = kwargs.get("dim", self.models[0].dim if self.models else 3)
         # prepare parameters (they are checked in dim setting)
@@ -1461,11 +1470,6 @@ def set_norm_len_ratios(self, ratios, skip=None, len_scale=None):
         """
         sum_set_norm_len_ratios(self, ratios, skip, len_scale)
 
-    @property
-    def hankel_kw(self):
-        """:class:`dict`: :any:`hankel.SymmetricFourierTransform` kwargs."""
-        return self._hankel_kw
-
     @property
     def models(self):
         """:class:`tuple`: The summed models."""
@@ -1623,6 +1627,11 @@ def calc_integral_scale(self):
             0.0,
         )
 
+    @property
+    def needs_fourier_transform(self):
+        """Whether the model doesn't implement an analytical spectral density."""
+        return False
+
     def spectral_density(self, k):
         return sum(
             (
diff --git a/src/gstools/covmodel/tools.py b/src/gstools/covmodel/tools.py
index 69bf963b1..5fbd179fd 100644
--- a/src/gstools/covmodel/tools.py
+++ b/src/gstools/covmodel/tools.py
@@ -548,7 +548,8 @@ def set_dim(model, dim):
         )
     model._dim = int(dim)
     # create fourier transform just once (recreate for dim change)
-    model._sft = SFT(ndim=model.dim, **model.hankel_kw)
+    if model.needs_fourier_transform:
+        model._sft = SFT(ndim=model.dim, **model.hankel_kw)
     # recalculate dimension related parameters (if model initialized)
     if model._init:
         model.len_scale, model.anis = set_len_anis(

From aee1985f756df073b52b0adf86156eec9af75eed Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Sebastian=20M=C3=BCller?= <mueller.seb@posteo.de>
Date: Mon, 6 Jan 2025 13:17:49 +0100
Subject: [PATCH 30/58] covmodel.tools: allow setting bounds of sub-args

---
 src/gstools/covmodel/fit.py   |  8 +++-----
 src/gstools/covmodel/tools.py | 25 ++++++++++++++++---------
 2 files changed, 19 insertions(+), 14 deletions(-)

diff --git a/src/gstools/covmodel/fit.py b/src/gstools/covmodel/fit.py
index 5ba398ef9..e3a1b014d 100755
--- a/src/gstools/covmodel/fit.py
+++ b/src/gstools/covmodel/fit.py
@@ -231,11 +231,9 @@ def _pre_para(model, para_select, sill, anis):
         if not sill_low <= sill <= sill_up:
             raise ValueError("fit: sill out of bounds.")
         if "var" in para_select and "nugget" in para_select:
-            if model.var > sill:
-                model.nugget = model.arg_bounds["nugget"][0]
-                model.var = sill - model.nugget
-            else:
-                model.nugget = sill - model.var
+            if not np.isclose(model.var + model.nugget, sill):
+                msg = "fit: if sill, var and nugget are fixed, var + nugget should match the given sill"
+                raise ValueError(msg)
         elif "var" in para_select:
             if model.var > sill:
                 raise ValueError(
diff --git a/src/gstools/covmodel/tools.py b/src/gstools/covmodel/tools.py
index 5fbd179fd..cfba8259a 100644
--- a/src/gstools/covmodel/tools.py
+++ b/src/gstools/covmodel/tools.py
@@ -469,16 +469,14 @@ def set_arg_bounds(model, check_args=True, **kwargs):
         if not check_bounds(bounds):
             msg = f"Given bounds for '{arg}' are not valid, got: {bounds}"
             raise ValueError(msg)
-        if arg in model.opt_arg:
+        if arg in getattr(model, "sub_arg", []):
+            # var_<i> and len_scale_<i>
+            name, i = _split_sub(arg)
+            setattr(model[i], f"_{name}_bounds", bounds)
+        elif arg in model.opt_arg:
             model._opt_arg_bounds[arg] = bounds
-        elif arg == "var":
-            model._var_bounds = bounds
-        elif arg == "len_scale":
-            model._len_scale_bounds = bounds
-        elif arg == "nugget":
-            model._nugget_bounds = bounds
-        elif arg == "anis":
-            model._anis_bounds = bounds
+        elif arg in model.arg_bounds:
+            setattr(model, f"_{arg}_bounds", bounds)
         else:
             raise ValueError(f"set_arg_bounds: unknown argument '{arg}'")
         if check_args and check_arg_in_bounds(model, arg) > 0:
@@ -489,6 +487,15 @@ def set_arg_bounds(model, check_args=True, **kwargs):
                 setattr(model, arg, def_arg)
 
 
+def _split_sub(name):
+    if name.startswith("var_"):
+        return "var", int(name[4:])
+    if name.startswith("len_scale_"):
+        return "len_scale", int(name[10:])
+    msg = f"Unknown sub variable: {name}"
+    raise ValueError(msg)
+
+
 def check_arg_bounds(model):
     """
     Check arguments to be within their given bounds.

From f30ebb739662e59aa4b1e181c979e05c3996f441 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Sebastian=20M=C3=BCller?= <mueller.seb@posteo.de>
Date: Mon, 6 Jan 2025 13:18:50 +0100
Subject: [PATCH 31/58] sum-tools: prepare setting sub-args by weights

---
 src/gstools/covmodel/sum_tools.py | 97 +++++++++++++++----------------
 1 file changed, 48 insertions(+), 49 deletions(-)

diff --git a/src/gstools/covmodel/sum_tools.py b/src/gstools/covmodel/sum_tools.py
index c4b9a90ee..11c31efbd 100644
--- a/src/gstools/covmodel/sum_tools.py
+++ b/src/gstools/covmodel/sum_tools.py
@@ -31,8 +31,8 @@
     "sum_check",
     "sum_default_arg_bounds",
     "sum_default_opt_arg_bounds",
-    "sum_set_norm_var_ratios",
-    "sum_set_norm_len_ratios",
+    "sum_set_var_weights",
+    "sum_set_len_weights",
     "sum_model_repr",
 ]
 
@@ -113,103 +113,102 @@ def sum_default_opt_arg_bounds(summod):
     return bounds
 
 
-def sum_set_norm_var_ratios(summod, ratios, skip=None, var=None):
+def sum_set_var_weights(summod, weights, skip=None, var=None):
     """
-    Set variances of contained models by normalized ratios in [0, 1].
-
-    Ratios are given as normalized ratios in [0, 1] as relative ratio of
-    variance to remaining difference to total variance of the Sum-Model.
+    Set variances of contained models by weights.
 
     Parameters
     ----------
-    ratios : iterable
-        Ratios to set. Should have a length of len(models) - len(exclude) - 1
+    weights : iterable
+        Weights to set. Should have a length of len(models) - len(skip)
     skip : iterable, optional
-        Model indices to skip. Should have compatible lenth, by default None
+        Model indices to skip. Should have compatible length, by default None
     var : float, optional
         Desired variance, by default current variance
 
     Raises
     ------
     ValueError
-        If number of ratios is not matching.
+        If number of weights is not matching.
     """
-    skip = skip or set()
-    if len(summod) != len(ratios) + len(skip) + 1:
-        msg = "SumModel.set_norm_ratios: number of ratios not matching."
+    skip = set() if skip is None else set(skip)
+    if len(summod) != len(weights) + len(skip):
+        msg = "SumModel.set_var_weights: number of ratios not matching."
         raise ValueError(msg)
     ids = range(len(summod))
     if fail := set(skip) - set(ids):
-        msg = f"SumModel.set_norm_var_ratios: skip ids not valid: {fail}"
+        msg = f"SumModel.set_var_weights: skip ids not valid: {fail}"
         raise ValueError(msg)
     var = summod.var if var is None else float(var)
-    check_arg_in_bounds(summod, "var", var, error=True)
     var_sum = sum(summod.models[i].var for i in skip)
-    if var_sum > var:
-        msg = "SumModel.set_norm_var_ratios: skiped variances to big."
+    var_diff = var - var_sum
+    if var_diff < 0:
+        msg = "SumModel.set_var_weights: skipped variances to big."
         raise RatioError(msg)
+    weights_sum = sum(weights)
+    vars = summod.vars
     j = 0
     for i in ids:
         if i in skip:
             continue
-        var_diff = var - var_sum
-        # last model gets remaining diff
-        var_set = var_diff * ratios[j] if j < len(ratios) else var_diff
-        summod[i].var = var_set
-        var_sum += var_set
+        vars[i] = var_diff * weights[j] / weights_sum
         j += 1
+    summod.vars = vars
 
 
-def sum_set_norm_len_ratios(summod, ratios, skip=None, len_scale=None):
+def sum_set_len_weights(summod, weights, skip=None, len_scale=None):
     """
-    Set length scales of contained models by normalized ratios in [0, 1].
-
-    Ratios are given as normalized ratios in [0, 1] as relative ratio of
-    len_scale * var / total_var to remaining difference to
-    total len_scale of the Sum-Model.
+    Set length scales of contained models by weights.
 
     Parameters
     ----------
-    ratios : iterable
-        Ratios to set. Should have a length of len(models) - len(exclude) - 1
+    weights : iterable
+        Weights to set. Should have a length of len(models) - len(skip)
     skip : iterable, optional
-        Model indices to skip. Should have compatible lenth, by default None
+        Model indices to skip. Should have compatible length, by default None
     len_scale : float, optional
         Desired len_scale, by default current len_scale
 
     Raises
     ------
     ValueError
-        If number of ratios is not matching.
+        If number of weights is not matching.
     """
-    skip = skip or set()
-    if len(summod) != len(ratios) + len(skip) + 1:
-        msg = "SumModel.set_norm_len_ratios: number of ratios not matching."
+    skip = set() if skip is None else set(skip)
+    if len(summod) != len(weights) + len(skip):
+        msg = "SumModel.set_len_weights: number of weights not matching."
         raise ValueError(msg)
     ids = range(len(summod))
     if fail := set(skip) - set(ids):
-        msg = f"SumModel.set_norm_len_ratios: skip ids not valid: {fail}"
+        msg = f"SumModel.set_len_weights: skip ids not valid: {fail}"
         raise ValueError(msg)
     len_scale = summod.len_scale if len_scale is None else float(len_scale)
-    check_arg_in_bounds(summod, "len_scale", len_scale, error=True)
+    # also skip models with no variance (not contributing to total len scale)
+    j = 0
+    wei = []
+    for i in ids:
+        if i in skip:
+            continue
+        if np.isclose(summod.ratios[i], 0):
+            skip.add(i)
+        else:
+            wei.append(weights[j])
+        j += 1
+    weights = wei
     len_sum = sum(summod[i].len_scale * summod.ratios[i] for i in skip)
-    if len_sum > len_scale:
-        msg = "SumModel.set_norm_len_ratios: skiped length scales to big."
+    len_diff = len_scale - len_sum
+    if len_diff < 0:
+        msg = "SumModel.set_len_weights: skipped length scales to big."
         raise RatioError(msg)
+    weights_sum = sum(weights)
+    len_scales = summod.len_scales
     j = 0
     for i in ids:
         if i in skip:
             continue
-        len_diff = len_scale - len_sum
-        # last model gets remaining diff
-        len_set = len_diff * ratios[j] if j < len(ratios) else len_diff
-        summod[i].len_scale = (
-            0.0
-            if np.isclose(summod.ratios[j], 0)
-            else len_set / summod.ratios[j]
-        )
-        len_sum += len_set
+        len_scales[i] = len_diff * weights[j] / weights_sum / summod.ratios[j]
         j += 1
+    summod.len_scales = len_scales
 
 
 def sum_model_repr(summod):  # pragma: no cover

From a8ef32f457b891e9e1b20ecdd8c8924e9c7f5371 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Sebastian=20M=C3=BCller?= <mueller.seb@posteo.de>
Date: Mon, 6 Jan 2025 13:19:24 +0100
Subject: [PATCH 32/58] SumModel: prepare sum models for fitting

---
 src/gstools/covmodel/base.py | 86 +++++++++++++++++++++++-------------
 1 file changed, 55 insertions(+), 31 deletions(-)

diff --git a/src/gstools/covmodel/base.py b/src/gstools/covmodel/base.py
index be927ca03..f633ff10b 100644
--- a/src/gstools/covmodel/base.py
+++ b/src/gstools/covmodel/base.py
@@ -25,8 +25,8 @@
     sum_default_arg_bounds,
     sum_default_opt_arg_bounds,
     sum_model_repr,
-    sum_set_norm_len_ratios,
-    sum_set_norm_var_ratios,
+    sum_set_len_weights,
+    sum_set_var_weights,
 )
 from gstools.covmodel.tools import (
     _init_subclass,
@@ -1403,7 +1403,7 @@ def __iter__(self):
         return iter(self.models)
 
     def __len__(self):
-        return len(self.models)
+        return self.size
 
     def __contains__(self, item):
         return item in self.models
@@ -1423,58 +1423,56 @@ def default_opt_arg_bounds(self):
         """Defaults boundaries for optional arguments as dict."""
         return sum_default_opt_arg_bounds(self)
 
-    def set_norm_var_ratios(self, ratios, skip=None, var=None):
+    def set_var_weights(self, weights, skip=None, var=None):
         """
-        Set variances of contained models by normalized ratios in [0, 1].
-
-        Ratios are given as normalized ratios in [0, 1] as relative ratio of
-        variance to remaining difference to total variance of the Sum-Model.
+        Set variances of contained models by weights.
 
         Parameters
         ----------
-        ratios : iterable
-            Ratios to set. Should have a length of len(models) - len(exclude) - 1
+        weights : iterable
+            Weights to set. Should have a length of len(models) - len(skip)
         skip : iterable, optional
-            Model indices to skip. Should have compatible lenth, by default None
+            Model indices to skip. Should have compatible length, by default None
         var : float, optional
             Desired variance, by default current variance
 
         Raises
         ------
         ValueError
-            If number of ratios is not matching.
+            If number of weights is not matching.
         """
-        sum_set_norm_var_ratios(self, ratios, skip, var)
+        sum_set_var_weights(self, weights, skip, var)
 
-    def set_norm_len_ratios(self, ratios, skip=None, len_scale=None):
+    def set_len_weights(self, weights, skip=None, len_scale=None):
         """
-        Set length scales of contained models by normalized ratios in [0, 1].
-
-        Ratios are given as normalized ratios in [0, 1] as relative ratio of
-        len_scale * var / total_var to remaining difference to
-        total len_scale of the Sum-Model.
+        Set length scales of contained models by weights.
 
         Parameters
         ----------
-        ratios : iterable
-            Ratios to set. Should have a length of len(models) - len(exclude) - 1
+        weights : iterable
+            Weights to set. Should have a length of len(models) - len(skip)
         skip : iterable, optional
-            Model indices to skip. Should have compatible lenth, by default None
+            Model indices to skip. Should have compatible length, by default None
         len_scale : float, optional
             Desired len_scale, by default current len_scale
 
         Raises
         ------
         ValueError
-            If number of ratios is not matching.
+            If number of weights is not matching.
         """
-        sum_set_norm_len_ratios(self, ratios, skip, len_scale)
+        sum_set_len_weights(self, weights, skip, len_scale)
 
     @property
     def models(self):
         """:class:`tuple`: The summed models."""
         return self._models
 
+    @property
+    def size(self):
+        """:class:`int`: Number of summed models."""
+        return len(self._models)
+
     @property
     def rescale(self):
         """:class:`float`: SumModel has a constant rescale factor of one."""
@@ -1665,13 +1663,39 @@ def opt_arg(self):
     @property
     def sub_arg(self):
         """:class:`list` of :class:`str`: Names of the sub-arguments for var and len_scale."""
-        return sum(
-            [
-                [f"{o}_{i}" for o in ["var", "len_scale"]]
-                for i, mod in enumerate(self.models)
-            ],
-            [],
-        )
+        return [
+            f"{o}_{i}" for o in ["var", "len_scale"] for i in range(self.size)
+        ]
+
+    @property
+    def sub_arg_bounds(self):
+        """:class:`dict`: Names of the sub-arguments for var and len_scale."""
+        return {
+            f"{o}_{i}": mod.arg_bounds[o]
+            for o in ["var", "len_scale"]
+            for (i, mod) in enumerate(self.models)
+        }
+
+    @property
+    def arg_bounds(self):
+        """:class:`dict`: Bounds for all parameters.
+
+        Notes
+        -----
+        Keys are the arg names and values are lists of 2 or 3 values:
+        ``[a, b]`` or ``[a, b, <type>]`` where
+        <type> is one of ``"oo"``, ``"cc"``, ``"oc"`` or ``"co"``
+        to define if the bounds are open ("o") or closed ("c").
+        """
+        res = {
+            "var": self.var_bounds,
+            "len_scale": self.len_scale_bounds,
+            "nugget": self.nugget_bounds,
+            "anis": self.anis_bounds,
+        }
+        res.update(self.opt_arg_bounds)
+        res.update(self.sub_arg_bounds)
+        return res
 
     def __setattr__(self, name, value):
         """Set an attribute."""

From 050cc3776421153ea0f10e797bad16a84d21dbdd Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Sebastian=20M=C3=BCller?= <mueller.seb@posteo.de>
Date: Mon, 6 Jan 2025 13:19:53 +0100
Subject: [PATCH 33/58] fit: enable fitting of sum-models; refactor fitting in
 general

---
 src/gstools/covmodel/fit.py | 383 ++++++++++++++++++++++++++++++------
 1 file changed, 325 insertions(+), 58 deletions(-)

diff --git a/src/gstools/covmodel/fit.py b/src/gstools/covmodel/fit.py
index e3a1b014d..455b93b7b 100755
--- a/src/gstools/covmodel/fit.py
+++ b/src/gstools/covmodel/fit.py
@@ -163,11 +163,7 @@ def fit_variogram(
     The fitted parameters will be instantly set in the model.
     """
     # preprocess selected parameters
-    para, sill, constrain_sill, anis = _pre_para(
-        model, para_select, sill, anis
-    )
-    if not any(para.values()):
-        raise ValueError("fit: no parameters selected for fitting.")
+    para, sill, anis, sum_cfg = _pre_para(model, para_select, sill, anis)
     # check curve_fit kwargs
     curve_fit_kwargs = curve_fit_kwargs or {}
     # check method
@@ -186,11 +182,11 @@ def fit_variogram(
     _set_weights(model, weights, x_data, curve_fit_kwargs, is_dir_vario)
     # set the lower/upper boundaries for the variogram-parameters
     bounds, init_guess_list = _init_curve_fit_para(
-        model, para, init_guess, constrain_sill, sill, anis
+        model, para, init_guess, sill, anis, sum_cfg
     )
     # create the fitting curve
     curve_fit_kwargs["f"] = _get_curve(
-        model, para, constrain_sill, sill, anis, is_dir_vario
+        model, para, sill, anis, is_dir_vario, sum_cfg
     )
     # set the remaining kwargs for curve_fit
     curve_fit_kwargs["bounds"] = bounds
@@ -203,7 +199,9 @@ def fit_variogram(
     # fit the variogram
     popt, pcov = curve_fit(**curve_fit_kwargs)
     # convert the results
-    fit_para = _post_fitting(model, para, popt, anis, is_dir_vario)
+    fit_para = _post_fitting(
+        model, para, popt, sill, anis, is_dir_vario, sum_cfg
+    )
     # calculate the r2 score if wanted
     if return_r2:
         return fit_para, pcov, _r2_score(model, x_data, y_data, is_dir_vario)
@@ -212,59 +210,188 @@ def fit_variogram(
 
 def _pre_para(model, para_select, sill, anis):
     """Preprocess selected parameters."""
+    sum_cfg = {"fix": {}}
+    is_sum = hasattr(model, "sub_arg")
+    sub_args = getattr(model, "sub_arg", [])
+    valid_args = model.iso_arg + sub_args
+    bnd = model.arg_bounds
     # if values given, set them in the model, afterwards all entries are bool
-    for par in para_select:
-        if par not in model.arg_bounds:
+    for par, val in para_select.items():
+        if par not in valid_args:
             raise ValueError(f"fit: unknown parameter in selection: {par}")
-        if not isinstance(para_select[par], bool):
-            setattr(model, par, float(para_select[par]))
+        # if parameters given with value, set it and deselect from fitting
+        if not isinstance(val, bool):
+            # don't set sub-args, var or len_scale in sum model yet
+            if is_sum and par in sub_args + ["var", "len_scale"]:
+                sum_cfg["fix"][par] = float(val)
+            else:
+                setattr(model, par, float(val))
             para_select[par] = False
     # remove those that were set to True
     para_select = {k: v for k, v in para_select.items() if not v}
+    # handling sum models
+    if is_sum:
+        _check_sum(model, para_select, sum_cfg, bnd)
     # handling the sill
     sill = None if (isinstance(sill, bool) and sill) else sill
     if sill is not None:
         sill = model.sill if isinstance(sill, bool) else float(sill)
-        constrain_sill = True
-        sill_low = model.arg_bounds["var"][0] + model.arg_bounds["nugget"][0]
-        sill_up = model.arg_bounds["var"][1] + model.arg_bounds["nugget"][1]
-        if not sill_low <= sill <= sill_up:
-            raise ValueError("fit: sill out of bounds.")
-        if "var" in para_select and "nugget" in para_select:
-            if not np.isclose(model.var + model.nugget, sill):
-                msg = "fit: if sill, var and nugget are fixed, var + nugget should match the given sill"
-                raise ValueError(msg)
-        elif "var" in para_select:
-            if model.var > sill:
-                raise ValueError(
-                    "fit: if sill is fixed and variance deselected, "
-                    "the set variance should be less than the given sill."
-                )
-            para_select["nugget"] = False
-            model.nugget = sill - model.var
-        elif "nugget" in para_select:
-            if model.nugget > sill:
-                raise ValueError(
-                    "fit: if sill is fixed and nugget deselected, "
-                    "the set nugget should be less than the given sill."
-                )
-            para_select["var"] = False
-            model.var = sill - model.nugget
-        else:
-            # deselect the nugget, to recalculate it accordingly
-            # nugget = sill - var
-            para_select["nugget"] = False
-    else:
-        constrain_sill = False
-    # select all parameters to be fitted
-    para = {par: True for par in model.iso_arg}
+        _check_sill(model, para_select, sill, bnd, sum_cfg)
+    # select all parameters to be fitted (if bounds do not indicate fixed parameter)
+    para = {par: bnd[par][0] < bnd[par][1] for par in valid_args}
     # now deselect unwanted parameters
     para.update(para_select)
     # check if anisotropy should be fitted or set
     if not isinstance(anis, bool):
         model.anis = anis
         anis = False
-    return para, sill, constrain_sill, anis
+    sub_fitting = sum_cfg.get("var_size", 0) + sum_cfg.get("len_size", 0) > 0
+    if not (any(para.values()) or anis or sub_fitting):
+        raise ValueError("fit: no parameters selected for fitting.")
+    return para, sill, anis, sum_cfg
+
+
+def _check_sill(model, para_select, sill, bnd, sum_cfg):
+    is_sum = hasattr(model, "sub_arg")
+    sill_low = bnd["var"][0] + bnd["nugget"][0]
+    sill_up = bnd["var"][1] + bnd["nugget"][1]
+    if not sill_low <= sill <= sill_up:
+        raise ValueError("fit: sill out of bounds.")
+    if is_sum:
+        var_fixed = sum_cfg["var_fix"]
+    else:
+        var_fixed = "var" in para_select
+    if var_fixed and "nugget" in para_select:
+        if not np.isclose(model.var + model.nugget, sill):
+            msg = "fit: if sill, var and nugget are fixed, var + nugget should match the given sill"
+            raise ValueError(msg)
+    elif var_fixed:
+        if model.var > sill:
+            raise ValueError(
+                "fit: if sill is fixed and variance deselected, "
+                "the set variance should be less than the given sill."
+            )
+        para_select["nugget"] = False
+        # this also works for a pure nugget model
+        model.nugget = sill - model.var
+    elif "nugget" in para_select:
+        if model.nugget > sill:
+            raise ValueError(
+                "fit: if sill is fixed and nugget deselected, "
+                "the set nugget should be less than the given sill."
+            )
+        para_select["var"] = False
+        var = sill - model.nugget
+        if is_sum:
+            model.set_var_weights(
+                np.ones_like(sum_cfg["var_fit"]), sum_cfg["var_skip"], var
+            )
+            # if only 1 sub-var was to fit, this is now also fixed
+            if len(sum_cfg["var_fit"]) == 1:
+                i = sum_cfg["var_fit"][0]
+                sum_cfg["fix"].setdefault(f"var_{i}", model.vars[i])
+                sum_cfg["var_fit"] = []
+                sum_cfg["var_skip"] = list(range(model.size))
+            sum_cfg["fix"].setdefault("var", var)
+            sum_cfg["var_fix"] = True
+            # number or sub-var parameters
+            sum_cfg["var_size"] = max(
+                model.size - len(sum_cfg["var_skip"]) - 1, 0
+            )
+        else:
+            # in case of a nugget model, this should raise an error
+            model.var = var
+    else:
+        # deselect the nugget, to recalculate it accordingly
+        # nugget = sill - var
+        para_select["nugget"] = False
+
+
+def _check_sum(model, para_select, sum_cfg, bnd):
+    """Check for consistent parameter selection in case of a SumModel."""
+    # check len_scale
+    if "len_scale" in para_select:
+        for par in para_select:
+            if par.startswith("len_scale_"):
+                msg = (
+                    "fit: for sum-models you can only fix "
+                    f"'len_scale' or the sub-arguments 'len_scale_<i>', not both."
+                )
+                raise ValueError(msg)
+        sum_cfg["fix"].setdefault("len_scale", model.len_scale)
+    # use len_scale_<i> for fitting if len_scale not fixed
+    # use weights for fitting if len_scale fixed
+    # either way: len_scale not used for fitting in sum_model
+    para_select["len_scale"] = False
+    # check variance
+    if "var" in para_select:
+        sum_cfg["fix"].setdefault("var", model.var)
+    # use var_<i> for fitting if var not fixed
+    # use weights for fitting if var fixed
+    # either way: var not used for fitting in sum_model
+    para_select["var"] = False
+    # whether var and len_scale are fixed
+    var_fix = "var" in sum_cfg["fix"]
+    len_fix = "len_scale" in sum_cfg["fix"]
+    size = model.size
+    # check for fixed bounds
+    for i in range(size):
+        for par in [f"var_{i}", f"len_scale_{i}"]:
+            if not bnd[par][0] < bnd[par][1]:
+                para_select[par] = False
+    # check sub arguments (var_<i> and len_scale_<i>)
+    var_skip = []
+    len_skip = []
+    remove = []
+    for par in para_select:
+        if par.startswith("var_"):
+            var_skip.append(int(par[4:]))
+            # for fixed var, fit by weights
+            if var_fix:
+                remove.append(par)
+        if par.startswith("len_scale_"):
+            len_skip.append(int(par[10:]))
+    for par in remove:
+        para_select.pop(par)
+    var_skip.sort()
+    len_skip.sort()
+    var_fit = sorted(set(range(size)) - set(var_skip))
+    len_fit = sorted(set(range(size)) - set(len_skip))
+    # if all sub-vars fixed, total variance is fixed
+    if not var_fit:
+        para_select["var"] = False
+    # set values related to var and len_scale in sum-model
+    for i in var_skip:
+        sum_cfg["fix"].setdefault(f"var_{i}", model.vars[i])
+        setattr(model, f"var_{i}", sum_cfg["fix"][f"var_{i}"])
+    for i in len_skip:
+        sum_cfg["fix"].setdefault(f"len_scale_{i}", model.len_scales[i])
+        setattr(model, f"len_scale_{i}", sum_cfg["fix"][f"len_scale_{i}"])
+    if var_fix:
+        var_min = sum(sum_cfg["fix"][f"var_{i}"] for i in var_skip)
+        if var_min > sum_cfg["fix"]["var"]:
+            msg = "fit: fixed sub-variances greater than fixed total variance."
+            raise ValueError(msg)
+        model.set_var_weights(
+            np.ones_like(var_fit), var_skip, sum_cfg["fix"]["var"]
+        )
+        # if all sub_vars except one is fixed and total var as well, all are fixed
+        if len(var_fit) == 1:
+            var_fit = []
+            var_skip = list(range(size))
+    if len_fix:
+        model.len_scale = sum_cfg["fix"]["len_scale"]
+    # update config for the sum-model fitting
+    sum_cfg["var_skip"] = var_skip  # sub vars to skip
+    sum_cfg["len_skip"] = len_skip  # sub lens to skip
+    sum_cfg["var_fit"] = var_fit  # sub vars to fit
+    sum_cfg["len_fit"] = len_fit  # sub lens to fit
+    sum_cfg["var_fix"] = var_fix  # total variance fixed
+    sum_cfg["len_fix"] = len_fix  # total len-scale fixed
+    # number or sub-var parameters
+    sum_cfg["var_size"] = max(size - len(var_skip) - int(var_fix), 0)
+    # number or sub-len parameters
+    sum_cfg["len_size"] = max(size - len(len_skip) - int(len_fix), 0)
 
 
 def _pre_init_guess(model, init_guess, mean_x=1.0, mean_y=1.0):
@@ -276,15 +403,17 @@ def _pre_init_guess(model, init_guess, mean_x=1.0, mean_y=1.0):
     if default_guess not in ["default", "current"]:
         raise ValueError(f"fit_variogram: unknown def. guess: {default_guess}")
     default = default_guess == "default"
+    is_sum = hasattr(model, "sub_arg")
+    sub_args = getattr(model, "sub_arg", [])
+    valid_args = model.iso_arg + sub_args + ["anis"]
     # check invalid names for given init guesses
-    invalid_para = set(init_guess) - set(model.iso_arg + ["anis"])
+    invalid_para = set(init_guess) - set(valid_args)
     if invalid_para:
         raise ValueError(f"fit_variogram: unknown init guess: {invalid_para}")
     bnd = model.arg_bounds
     # default length scale is mean of given bin centers (respecting "rescale")
-    init_guess.setdefault(
-        "len_scale", mean_x * model.rescale if default else model.len_scale
-    )
+    len_def = mean_x * model.rescale
+    init_guess.setdefault("len_scale", len_def if default else model.len_scale)
     # init guess for variance and nugget is mean of given variogram
     for par in ["var", "nugget"]:
         init_guess.setdefault(par, mean_y if default else getattr(model, par))
@@ -304,8 +433,17 @@ def _pre_init_guess(model, init_guess, mean_x=1.0, mean_y=1.0):
                 else getattr(model, opt)
             ),
         )
+    # SumModel: check for var_<i> and len_scale_<i> and set defaults
+    if is_sum:
+        for i in range(model.size):
+            init_guess.setdefault(
+                f"len_scale_{i}", len_def if default else model.len_scales[i]
+            )
+            init_guess.setdefault(
+                f"var_{i}", mean_y if default else model.len_scales[i]
+            )
     # convert all init guesses to float (except "anis")
-    for arg in model.iso_arg:
+    for arg in model.iso_arg + sub_args:
         init_guess[arg] = float(init_guess[arg])
     return init_guess
 
@@ -351,24 +489,59 @@ def _set_weights(model, weights, x_data, curve_fit_kwargs, is_dir_vario):
         curve_fit_kwargs["absolute_sigma"] = True
 
 
-def _init_curve_fit_para(model, para, init_guess, constrain_sill, sill, anis):
+def _init_curve_fit_para(model, para, init_guess, sill, anis, sum_cfg):
     """Create initial guess and bounds for fitting."""
+    is_sum = hasattr(model, "sub_arg")
     low_bounds = []
     top_bounds = []
     init_guess_list = []
+    bnd = model.arg_bounds
     for par in model.iso_arg:
         if para[par]:
-            low_bounds.append(model.arg_bounds[par][0])
-            if par == "var" and constrain_sill:  # var <= sill in this case
+            low_bounds.append(bnd[par][0])
+            if par == "var" and sill is not None:  # var <= sill in this case
                 top_bounds.append(sill)
             else:
-                top_bounds.append(model.arg_bounds[par][1])
+                top_bounds.append(bnd[par][1])
             init_guess_list.append(
                 _init_guess(
                     bounds=[low_bounds[-1], top_bounds[-1]],
                     default=init_guess[par],
                 )
             )
+    if is_sum:
+        if sum_cfg["var_fix"]:
+            for _ in range(sum_cfg["var_size"]):
+                low_bounds.append(0.0)
+                top_bounds.append(1.0)
+                init_guess_list.append(0.5)
+        else:
+            for i in sum_cfg["var_fit"]:
+                par = f"var_{i}"
+                low_bounds.append(bnd[par][0])
+                top_bounds.append(bnd[par][1])
+                init_guess_list.append(
+                    _init_guess(
+                        bounds=[low_bounds[-1], top_bounds[-1]],
+                        default=init_guess[par],
+                    )
+                )
+        if sum_cfg["len_fix"]:
+            for _ in range(sum_cfg["len_size"]):
+                low_bounds.append(0.0)
+                top_bounds.append(1.0)
+                init_guess_list.append(0.5)
+        else:
+            for i in sum_cfg["len_fit"]:
+                par = f"len_scale_{i}"
+                low_bounds.append(bnd[par][0])
+                top_bounds.append(bnd[par][1])
+                init_guess_list.append(
+                    _init_guess(
+                        bounds=[low_bounds[-1], top_bounds[-1]],
+                        default=init_guess[par],
+                    )
+                )
     if anis:
         for i in range(model.dim - 1):
             low_bounds.append(model.anis_bounds[0])
@@ -389,8 +562,9 @@ def _init_guess(bounds, default):
     return default_arg_from_bounds(bounds)
 
 
-def _get_curve(model, para, constrain_sill, sill, anis, is_dir_vario):
+def _get_curve(model, para, sill, anis, is_dir_vario, sum_cfg):
     """Create the curve for scipys curve_fit."""
+    is_sum = hasattr(model, "sub_arg")
 
     # we need arg1, otherwise curve_fit throws an error (bug?!)
     def curve(x, arg1, *args):
@@ -401,7 +575,34 @@ def curve(x, arg1, *args):
             if para[par]:
                 setattr(model, par, args[para_skip])
                 para_skip += 1
-        if constrain_sill:
+        # set var and len-scale ratios in sum-models
+        if is_sum:
+            if sum_cfg["var_size"] > 0:
+                var_vals = args[para_skip : para_skip + sum_cfg["var_size"]]
+                para_skip += sum_cfg["var_size"]
+                if sum_cfg["var_fix"]:
+                    model.set_var_weights(
+                        stick_breaking_uniform(var_vals),
+                        sum_cfg["var_skip"],
+                        sum_cfg["fix"]["var"],
+                    )
+                else:
+                    for i, val in zip(sum_cfg["var_fit"], var_vals):
+                        setattr(model, f"var_{i}", val)
+            if sum_cfg["len_size"] > 0:
+                len_vals = args[para_skip : para_skip + sum_cfg["len_size"]]
+                para_skip += sum_cfg["len_size"]
+                if sum_cfg["len_fix"]:
+                    model.set_len_weights(
+                        stick_breaking_uniform(len_vals),
+                        sum_cfg["len_skip"],
+                        sum_cfg["fix"]["len_scale"],
+                    )
+                else:
+                    for i, val in zip(sum_cfg["len_fit"], len_vals):
+                        setattr(model, f"len_scale_{i}", val)
+        # handle sill
+        if sill is not None and para["var"]:
             nugget_tmp = sill - model.var
             # punishment, if resulting nugget out of range for fixed sill
             if check_arg_in_bounds(model, "nugget", nugget_tmp) > 0:
@@ -421,8 +622,9 @@ def curve(x, arg1, *args):
     return curve
 
 
-def _post_fitting(model, para, popt, anis, is_dir_vario):
+def _post_fitting(model, para, popt, sill, anis, is_dir_vario, sum_cfg):
     """Postprocess fitting results and application to model."""
+    is_sum = hasattr(model, "sub_arg")
     fit_para = {}
     para_skip = 0
     for par in model.iso_arg:
@@ -432,6 +634,40 @@ def _post_fitting(model, para, popt, anis, is_dir_vario):
             para_skip += 1
         else:
             fit_para[par] = getattr(model, par)
+    # set var and len-scale ratios in sum-models
+    if is_sum:
+        if sum_cfg["var_size"] > 0:
+            var_vals = popt[para_skip : para_skip + sum_cfg["var_size"]]
+            para_skip += sum_cfg["var_size"]
+            if sum_cfg["var_fix"]:
+                model.set_var_weights(
+                    stick_breaking_uniform(var_vals),
+                    sum_cfg["var_skip"],
+                    sum_cfg["fix"]["var"],
+                )
+            else:
+                for i, val in zip(sum_cfg["var_fit"], var_vals):
+                    setattr(model, f"var_{i}", val)
+        if sum_cfg["len_size"] > 0:
+            len_vals = popt[para_skip : para_skip + sum_cfg["len_size"]]
+            para_skip += sum_cfg["len_size"]
+            if sum_cfg["len_fix"]:
+                model.set_len_weights(
+                    stick_breaking_uniform(len_vals),
+                    sum_cfg["len_skip"],
+                    sum_cfg["fix"]["len_scale"],
+                )
+            else:
+                for i, val in zip(sum_cfg["len_fit"], len_vals):
+                    setattr(model, f"len_scale_{i}", val)
+        for i in range(model.size):
+            fit_para[f"var_{i}"] = model.vars[i]
+            fit_para[f"len_scale_{i}"] = model.len_scales[i]
+    # handle sill
+    if sill is not None and para["var"]:
+        nugget = sill - model.var
+        fit_para["nugget"] = nugget
+        model.nugget = nugget
     if is_dir_vario:
         if anis:
             model.anis = popt[1 - model.dim :]
@@ -484,3 +720,34 @@ def func(x_data):
         return 1.0 / (1.0 + np.exp(growth * (x_mean - x_data)))
 
     return func
+
+
+def stick_breaking_uniform(u):
+    """
+    Generate a single sample (x_1, ..., x_n) uniformly from the (n-1)-simplex.
+
+    This is using Beta transforms of uniform samples. The x_i will sum to 1.
+
+    Parameters
+    ----------
+    u : array-like of shape (n-1,)
+        Uniform(0,1) random values
+
+    Returns
+    -------
+    x : ndarray of shape (n,)
+        A random sample in the (n-1)-simplex.
+    """
+    n = len(u) + 1
+    x = np.zeros(n)
+    leftover = 1.0
+    for i in range(n - 1):
+        # 2) Compute the inverse CDF of Beta(1, b) = Beta(1, n-1-i)
+        fraction = 1.0 - (1.0 - u[i]) ** (1.0 / (n - 1 - i))
+        # 3) Break that fraction of the current leftover
+        x[i] = leftover * fraction
+        # 4) Subtract that from leftover
+        leftover -= x[i]
+    # Last coordinate
+    x[-1] = leftover
+    return x

From b120d935a9a2d55672bfca1db66938a81b79df58 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Sebastian=20M=C3=BCller?= <mueller.seb@posteo.de>
Date: Mon, 6 Jan 2025 15:32:47 +0100
Subject: [PATCH 34/58] tests: update tests to check for new var/nugget/sill
 handling

---
 tests/test_covmodel.py | 16 ++++++----------
 1 file changed, 6 insertions(+), 10 deletions(-)

diff --git a/tests/test_covmodel.py b/tests/test_covmodel.py
index c42e47746..62ea8a1c2 100644
--- a/tests/test_covmodel.py
+++ b/tests/test_covmodel.py
@@ -239,20 +239,16 @@ def test_fitting(self):
 
         # treatment of sill/var/nugget by fitting
         model = Stable()
-        model.fit_variogram(
-            self.gamma_x, self.gamma_y, nugget=False, var=False, sill=2
-        )
-        self.assertAlmostEqual(model.var, 1)
-        self.assertAlmostEqual(model.nugget, 1)
+        with self.assertRaises(ValueError):
+            model.fit_variogram(
+                self.gamma_x, self.gamma_y, nugget=False, var=False, sill=2
+            )
         model.fit_variogram(self.gamma_x, self.gamma_y, var=2, sill=3)
         self.assertAlmostEqual(model.var, 2)
         self.assertAlmostEqual(model.nugget, 1)
         model.var = 3
-        model.fit_variogram(
-            self.gamma_x, self.gamma_y, nugget=False, var=False, sill=2
-        )
-        self.assertAlmostEqual(model.var, 2)
-        self.assertAlmostEqual(model.nugget, 0)
+        with self.assertRaises(ValueError):
+            model.fit_variogram(self.gamma_x, self.gamma_y, var=False, sill=2)
         model.fit_variogram(self.gamma_x, self.gamma_y, weights="inv")
         len_save = model.len_scale
         model.fit_variogram(

From 302a056fb7d0c2f78d61579bffabcafa5ccefb8d Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Sebastian=20M=C3=BCller?= <mueller.seb@posteo.de>
Date: Mon, 6 Jan 2025 15:51:50 +0100
Subject: [PATCH 35/58] CI: update to include py3.13

---
 .github/workflows/main.yml | 11 ++++++-----
 1 file changed, 6 insertions(+), 5 deletions(-)

diff --git a/.github/workflows/main.yml b/.github/workflows/main.yml
index 9fc071d59..e22d3d52c 100644
--- a/.github/workflows/main.yml
+++ b/.github/workflows/main.yml
@@ -22,10 +22,10 @@ jobs:
         with:
           fetch-depth: '0'
 
-      - name: Set up Python 3.9
+      - name: Set up Python 3.11
         uses: actions/setup-python@v5
         with:
-          python-version: 3.9
+          python-version: 3.11
 
       - name: Install dependencies
         run: |
@@ -75,7 +75,7 @@ jobs:
           fetch-depth: '0'
 
       - name: Build wheels
-        uses: pypa/cibuildwheel@v2.18.0
+        uses: pypa/cibuildwheel@v2.22.0
         with:
           output-dir: dist
 
@@ -97,7 +97,8 @@ jobs:
           - {py: '3.10', np: '==1.21.6', sp: '==1.7.2'}
           - {py: '3.11', np: '==1.23.2', sp: '==1.9.2'}
           - {py: '3.12', np: '==1.26.2', sp: '==1.11.2'}
-          - {py: '3.12', np: '>=2.0.0rc1', sp: '>=1.13.0'}
+          - {py: '3.13', np: '==2.1.0', sp: '==1.14.1'}
+          - {py: '3.13', np: '>=2.1.0', sp: '>=1.14.1'}
         exclude:
           - os: macos-14
             ver: {py: '3.8', np: '==1.20.0', sp: '==1.5.4'}
@@ -137,7 +138,7 @@ jobs:
           python -m build --sdist --outdir dist .
 
       - uses: actions/upload-artifact@v3
-        if: matrix.os == 'ubuntu-latest' && matrix.ver.py == '3.9'
+        if: matrix.os == 'ubuntu-latest' && matrix.ver.py == '3.11'
         with:
           path: dist/*.tar.gz
 

From 9c2673eac3aa8a27a68b42d19959dd2223bb7b9b Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Sebastian=20M=C3=BCller?= <mueller.seb@posteo.de>
Date: Mon, 6 Jan 2025 15:52:47 +0100
Subject: [PATCH 36/58] Doc: use myst_parser instead of m2r2 and update docs

---
 AUTHORS.md                |  6 +++---
 README.md                 |  4 ++--
 docs/Makefile             |  2 +-
 docs/source/authors.rst   |  2 ++
 docs/source/changelog.rst |  3 ++-
 docs/source/conf.py       | 26 +++++++++++++++-----------
 docs/source/contents.rst  |  1 +
 docs/source/index.rst     |  2 +-
 pyproject.toml            |  4 ++--
 src/gstools/__init__.py   |  2 +-
 10 files changed, 30 insertions(+), 22 deletions(-)
 create mode 100644 docs/source/authors.rst

diff --git a/AUTHORS.md b/AUTHORS.md
index 63601d87c..2feeba611 100644
--- a/AUTHORS.md
+++ b/AUTHORS.md
@@ -1,10 +1,10 @@
-# GSTools
+# Authors
 
 GSTools is available on [GitHub](https://github.com/GeoStat-Framework/GSTools)
 and was created by following people.
 
 
-## Main Authors
+## Core developers
 
 - Sebastian Müller, GitHub: [@MuellerSeb](https://github.com/MuellerSeb), Email:  <sebastian@geostat-framework.org>
 - Lennart Schüler, GitHub: [@LSchueler](https://github.com/LSchueler), Email:  <lennart@geostat-framework.org>
@@ -14,4 +14,4 @@ and was created by following people.
 
 - Falk Heße, GitHub: [@fhesze](https://github.com/fhesze), Email: <falk.hesse@ufz.de>
 - Bane Sullivan, GitHub: [@banesullivan](https://github.com/banesullivan)
-- Tobias Glaubach, GitHub: [@TobiasGlaubach](https://github.com/TobiasGlaubach)
\ No newline at end of file
+- Tobias Glaubach, GitHub: [@TobiasGlaubach](https://github.com/TobiasGlaubach)
diff --git a/README.md b/README.md
index b90ad1294..58feb3f56 100644
--- a/README.md
+++ b/README.md
@@ -33,7 +33,7 @@
 
 <img align="right" width="450" src="https://raw.githubusercontent.com/GeoStat-Framework/GSTools/main/docs/source/pics/demonstrator.png" alt="">
 
-GeoStatTools provides geostatistical tools for various purposes:
+GSTools provides geostatistical tools for various purposes:
 - random field generation, including periodic boundaries
 - simple, ordinary, universal and external drift kriging
 - conditioned field generation
@@ -369,7 +369,7 @@ You can contact us via <info@geostat-framework.org>.
 
 ## License
 
-[LGPLv3][license_link] © 2018-2024
+[LGPLv3][license_link] © 2018-2025
 
 [pip_link]: https://pypi.org/project/gstools
 [conda_link]: https://docs.conda.io/en/latest/miniconda.html
diff --git a/docs/Makefile b/docs/Makefile
index 39bd7f46e..d1a288eba 100644
--- a/docs/Makefile
+++ b/docs/Makefile
@@ -4,7 +4,7 @@
 # You can set these variables from the command line.
 SPHINXOPTS    =
 SPHINXBUILD   = python3 -msphinx
-SPHINXPROJ    = GeoStatTools
+SPHINXPROJ    = GSTools
 SOURCEDIR     = source
 BUILDDIR      = build
 
diff --git a/docs/source/authors.rst b/docs/source/authors.rst
new file mode 100644
index 000000000..49b950323
--- /dev/null
+++ b/docs/source/authors.rst
@@ -0,0 +1,2 @@
+.. include:: ../../AUTHORS.md
+   :parser: myst_parser.docutils_
diff --git a/docs/source/changelog.rst b/docs/source/changelog.rst
index ab37940f4..a01bcd9f8 100644
--- a/docs/source/changelog.rst
+++ b/docs/source/changelog.rst
@@ -1 +1,2 @@
-.. mdinclude:: ../../CHANGELOG.md
+.. include:: ../../CHANGELOG.md
+   :parser: myst_parser.docutils_
diff --git a/docs/source/conf.py b/docs/source/conf.py
index e89928fc9..f0dc1a1e6 100644
--- a/docs/source/conf.py
+++ b/docs/source/conf.py
@@ -1,6 +1,6 @@
 #!/usr/bin/env python3
 #
-# GeoStatTools documentation build configuration file, created by
+# GSTools documentation build configuration file, created by
 # sphinx-quickstart on Fri Jan  5 14:20:43 2018.
 #
 # This file is execfile()d with the current directory set to its
@@ -67,7 +67,7 @@ def setup(app):
     "sphinx.ext.napoleon",  # parameters look better than with numpydoc only
     "numpydoc",
     "sphinx_gallery.gen_gallery",
-    "m2r2",
+    "myst_parser",
     "sphinxcontrib.youtube",
 ]
 
@@ -98,7 +98,11 @@ def setup(app):
 
 # The suffix(es) of source filenames.
 # You can specify multiple suffix as a list of string:
-source_suffix = [".rst", ".md"]
+source_suffix = {
+    ".rst": "restructuredtext",
+    ".md": "markdown",
+}
+# source_suffix = [".rst", ".md"]
 # source_suffix = ".rst"
 
 # The master toctree document.
@@ -152,7 +156,7 @@ def setup(app):
     #    'canonical_url': '',
     #    'analytics_id': '',
     "logo_only": False,
-    "display_version": True,
+    "version_selector": True,
     "prev_next_buttons_location": "top",
     #    'style_external_links': False,
     #    'vcs_pageview_mode': '',
@@ -188,7 +192,7 @@ def setup(app):
 # -- Options for HTMLHelp output ------------------------------------------
 
 # Output file base name for HTML help builder.
-htmlhelp_basename = "GeoStatToolsdoc"
+htmlhelp_basename = "GSToolsdoc"
 # logos for the page
 html_logo = "pics/gstools_150.png"
 html_favicon = "pics/gstools.ico"
@@ -217,8 +221,8 @@ def setup(app):
 latex_documents = [
     (
         master_doc,
-        "GeoStatTools.tex",
-        "GeoStatTools Documentation",
+        "GSTools.tex",
+        "GSTools Documentation",
         "Sebastian Müller, Lennart Schüler",
         "manual",
     )
@@ -230,7 +234,7 @@ def setup(app):
 # One entry per manual page. List of tuples
 # (source start file, name, description, authors, manual section).
 man_pages = [
-    (master_doc, "geostattools", "GeoStatTools Documentation", [author], 1)
+    (master_doc, "GSTools", "GSTools Documentation", [author], 1)
 ]
 
 
@@ -242,10 +246,10 @@ def setup(app):
 texinfo_documents = [
     (
         master_doc,
-        "GeoStatTools",
-        "GeoStatTools Documentation",
+        "GSTools",
+        "GSTools Documentation",
         author,
-        "GeoStatTools",
+        "GSTools",
         "Geo-statistical toolbox.",
         "Miscellaneous",
     )
diff --git a/docs/source/contents.rst b/docs/source/contents.rst
index 3224356ee..b6979a5ce 100644
--- a/docs/source/contents.rst
+++ b/docs/source/contents.rst
@@ -9,4 +9,5 @@ Contents
    index
    tutorials
    api
+   authors
    changelog
diff --git a/docs/source/index.rst b/docs/source/index.rst
index 7abb1e9e9..bba2309e8 100644
--- a/docs/source/index.rst
+++ b/docs/source/index.rst
@@ -22,7 +22,7 @@ GSTools Quickstart
 Purpose
 =======
 
-GeoStatTools provides geostatistical tools for various purposes:
+GSTools provides geostatistical tools for various purposes:
 
 - random field generation, including periodic boundaries
 - simple, ordinary, universal and external drift kriging
diff --git a/pyproject.toml b/pyproject.toml
index 88185bf30..f076a8dca 100644
--- a/pyproject.toml
+++ b/pyproject.toml
@@ -56,7 +56,7 @@ dependencies = [
 
 [project.optional-dependencies]
 doc = [
-    "m2r2>=0.2.8",
+    "myst_parser",
     "matplotlib>=3.7",
     "meshzoo>=0.7",
     "numpydoc>=1.1",
@@ -64,7 +64,7 @@ doc = [
     "pyvista>=0.40",
     "sphinx>=7",
     "sphinx-gallery>=0.8",
-    "sphinx-rtd-theme>=2",
+    "sphinx-rtd-theme>=3",
     "sphinxcontrib-youtube>=1.1",
 ]
 plotting = [
diff --git a/src/gstools/__init__.py b/src/gstools/__init__.py
index 7488dc9ce..bd92c8329 100644
--- a/src/gstools/__init__.py
+++ b/src/gstools/__init__.py
@@ -2,7 +2,7 @@
 Purpose
 =======
 
-GeoStatTools is a library providing geostatistical tools
+GSTools is a library providing geostatistical tools
 for random field generation, conditioned field generation,
 kriging and variogram estimation
 based on a list of provided or even user-defined covariance models.

From d03146faf7b8f96c88f3a7cb28b11cc9ed2a6365 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Sebastian=20M=C3=BCller?= <mueller.seb@posteo.de>
Date: Mon, 6 Jan 2025 17:26:50 +0100
Subject: [PATCH 37/58] lint

---
 docs/source/conf.py               | 4 +---
 src/gstools/covmodel/fit.py       | 2 +-
 src/gstools/covmodel/sum_tools.py | 7 +++----
 3 files changed, 5 insertions(+), 8 deletions(-)

diff --git a/docs/source/conf.py b/docs/source/conf.py
index f0dc1a1e6..af4adc6ba 100644
--- a/docs/source/conf.py
+++ b/docs/source/conf.py
@@ -233,9 +233,7 @@ def setup(app):
 
 # One entry per manual page. List of tuples
 # (source start file, name, description, authors, manual section).
-man_pages = [
-    (master_doc, "GSTools", "GSTools Documentation", [author], 1)
-]
+man_pages = [(master_doc, "GSTools", "GSTools Documentation", [author], 1)]
 
 
 # -- Options for Texinfo output -------------------------------------------
diff --git a/src/gstools/covmodel/fit.py b/src/gstools/covmodel/fit.py
index 455b93b7b..f6363f51d 100755
--- a/src/gstools/covmodel/fit.py
+++ b/src/gstools/covmodel/fit.py
@@ -315,7 +315,7 @@ def _check_sum(model, para_select, sum_cfg, bnd):
             if par.startswith("len_scale_"):
                 msg = (
                     "fit: for sum-models you can only fix "
-                    f"'len_scale' or the sub-arguments 'len_scale_<i>', not both."
+                    "'len_scale' or the sub-arguments 'len_scale_<i>', not both."
                 )
                 raise ValueError(msg)
         sum_cfg["fix"].setdefault("len_scale", model.len_scale)
diff --git a/src/gstools/covmodel/sum_tools.py b/src/gstools/covmodel/sum_tools.py
index 11c31efbd..7e6c6b4f3 100644
--- a/src/gstools/covmodel/sum_tools.py
+++ b/src/gstools/covmodel/sum_tools.py
@@ -20,7 +20,6 @@
 # pylint: disable=W0212
 import numpy as np
 
-from gstools.covmodel.tools import check_arg_in_bounds
 from gstools.tools import RADIAN_SCALE
 from gstools.tools.misc import list_format
 
@@ -146,14 +145,14 @@ def sum_set_var_weights(summod, weights, skip=None, var=None):
         msg = "SumModel.set_var_weights: skipped variances to big."
         raise RatioError(msg)
     weights_sum = sum(weights)
-    vars = summod.vars
+    var_list = summod.vars
     j = 0
     for i in ids:
         if i in skip:
             continue
-        vars[i] = var_diff * weights[j] / weights_sum
+        var_list[i] = var_diff * weights[j] / weights_sum
         j += 1
-    summod.vars = vars
+    summod.vars = var_list
 
 
 def sum_set_len_weights(summod, weights, skip=None, len_scale=None):

From 5dc6b16eebde2aa31806a4ebed9fb4a03829f1f6 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Sebastian=20M=C3=BCller?= <mueller.seb@posteo.de>
Date: Mon, 6 Jan 2025 18:03:23 +0100
Subject: [PATCH 38/58] pylint: set all config in pyproject

---
 .github/workflows/main.yml | 2 +-
 pyproject.toml             | 1 +
 2 files changed, 2 insertions(+), 1 deletion(-)

diff --git a/.github/workflows/main.yml b/.github/workflows/main.yml
index 754218338..c78431fbb 100644
--- a/.github/workflows/main.yml
+++ b/.github/workflows/main.yml
@@ -47,7 +47,7 @@ jobs:
 
       - name: pylint check
         run: |
-          python -m pylint --max-positional-arguments 20 src/gstools/
+          python -m pylint src/gstools/
 
       - name: cython-lint check
         run: |
diff --git a/pyproject.toml b/pyproject.toml
index f076a8dca..f7fbc3fa2 100644
--- a/pyproject.toml
+++ b/pyproject.toml
@@ -157,6 +157,7 @@ target-version = [
     max-statements = 85
     max-attributes = 30
     max-public-methods = 80
+    max-positional-arguments=20
 
 [tool.cibuildwheel]
 # Switch to using build

From 64b102eb82a3beea13a2e93885eb866b9ff3ef95 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Sebastian=20M=C3=BCller?= <mueller.seb@posteo.de>
Date: Mon, 6 Jan 2025 18:03:59 +0100
Subject: [PATCH 39/58] Normalizer: implement simple derivative to remove scipy
 dependency

---
 src/gstools/normalizer/base.py | 8 ++++++--
 1 file changed, 6 insertions(+), 2 deletions(-)

diff --git a/src/gstools/normalizer/base.py b/src/gstools/normalizer/base.py
index 4a8477c60..52c880539 100644
--- a/src/gstools/normalizer/base.py
+++ b/src/gstools/normalizer/base.py
@@ -13,10 +13,14 @@
 import warnings
 
 import numpy as np
-import scipy.misc as spm
 import scipy.optimize as spo
 
 
+def _derivative(f, x, dx=1e-6):
+    """Central difference formula."""
+    return (f(x + dx) - f(x - dx)) / (2 * dx)
+
+
 class Normalizer:
     """Normalizer class.
 
@@ -57,7 +61,7 @@ def _normalize(self, data):
         return data
 
     def _derivative(self, data):
-        return spm.derivative(self._normalize, data, dx=self._dx)
+        return _derivative(self._normalize, data, dx=self._dx)
 
     def _loglikelihood(self, data):
         add = -0.5 * np.size(data) * (np.log(2 * np.pi) + 1)

From f475c1e85aaeb3b4f4f0d1e5f7ddad2c7c68438d Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Sebastian=20M=C3=BCller?= <mueller.seb@posteo.de>
Date: Mon, 6 Jan 2025 18:19:01 +0100
Subject: [PATCH 40/58] fit: correctly check for parameters to fit

---
 src/gstools/covmodel/fit.py | 10 +++++-----
 1 file changed, 5 insertions(+), 5 deletions(-)

diff --git a/src/gstools/covmodel/fit.py b/src/gstools/covmodel/fit.py
index f6363f51d..7765ca54e 100755
--- a/src/gstools/covmodel/fit.py
+++ b/src/gstools/covmodel/fit.py
@@ -172,12 +172,15 @@ def fit_variogram(
     # prepare variogram data
     # => concatenate directional variograms to have a 1D array for x and y
     x_data, y_data, is_dir_vario = _check_vario(model, x_data, y_data)
+    # only fit anisotropy if a directional variogram was given
+    anis &= is_dir_vario
+    sub_fitting = sum_cfg.get("var_size", 0) + sum_cfg.get("len_size", 0) > 0
+    if not (any(para.values()) or anis or sub_fitting):
+        raise ValueError("fit: no parameters selected for fitting.")
     # prepare init guess dictionary
     init_guess = _pre_init_guess(
         model, init_guess, np.mean(x_data), np.mean(y_data)
     )
-    # only fit anisotropy if a directional variogram was given
-    anis &= is_dir_vario
     # set weights
     _set_weights(model, weights, x_data, curve_fit_kwargs, is_dir_vario)
     # set the lower/upper boundaries for the variogram-parameters
@@ -245,9 +248,6 @@ def _pre_para(model, para_select, sill, anis):
     if not isinstance(anis, bool):
         model.anis = anis
         anis = False
-    sub_fitting = sum_cfg.get("var_size", 0) + sum_cfg.get("len_size", 0) > 0
-    if not (any(para.values()) or anis or sub_fitting):
-        raise ValueError("fit: no parameters selected for fitting.")
     return para, sill, anis, sum_cfg
 
 

From 8ea78f8c65d83eb9a88ae8b2aa086379b2393768 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Sebastian=20M=C3=BCller?= <mueller.seb@posteo.de>
Date: Mon, 6 Jan 2025 23:26:39 +0100
Subject: [PATCH 41/58] fit: fix init guess for var_i

---
 src/gstools/covmodel/fit.py | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/src/gstools/covmodel/fit.py b/src/gstools/covmodel/fit.py
index 7765ca54e..54c0e06c6 100755
--- a/src/gstools/covmodel/fit.py
+++ b/src/gstools/covmodel/fit.py
@@ -440,7 +440,7 @@ def _pre_init_guess(model, init_guess, mean_x=1.0, mean_y=1.0):
                 f"len_scale_{i}", len_def if default else model.len_scales[i]
             )
             init_guess.setdefault(
-                f"var_{i}", mean_y if default else model.len_scales[i]
+                f"var_{i}", mean_y / model.size if default else model.len_vars[i]
             )
     # convert all init guesses to float (except "anis")
     for arg in model.iso_arg + sub_args:

From 3eee1ddc06fe72e591d5c39ae99854ee863a38a2 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Sebastian=20M=C3=BCller?= <mueller.seb@posteo.de>
Date: Mon, 6 Jan 2025 23:36:43 +0100
Subject: [PATCH 42/58] fix typo

---
 src/gstools/covmodel/fit.py | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/src/gstools/covmodel/fit.py b/src/gstools/covmodel/fit.py
index 54c0e06c6..0384ccbf5 100755
--- a/src/gstools/covmodel/fit.py
+++ b/src/gstools/covmodel/fit.py
@@ -440,7 +440,7 @@ def _pre_init_guess(model, init_guess, mean_x=1.0, mean_y=1.0):
                 f"len_scale_{i}", len_def if default else model.len_scales[i]
             )
             init_guess.setdefault(
-                f"var_{i}", mean_y / model.size if default else model.len_vars[i]
+                f"var_{i}", mean_y / model.size if default else model.vars[i]
             )
     # convert all init guesses to float (except "anis")
     for arg in model.iso_arg + sub_args:

From e3b944e0b585c996d84d7a37d82646a17d9d09e5 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Sebastian=20M=C3=BCller?= <mueller.seb@posteo.de>
Date: Tue, 7 Jan 2025 00:08:40 +0100
Subject: [PATCH 43/58] Doc: add mini-gallery to tutorials

---
 docs/source/conf.py       | 2 +-
 docs/source/tutorials.rst | 9 +++++++++
 2 files changed, 10 insertions(+), 1 deletion(-)

diff --git a/docs/source/conf.py b/docs/source/conf.py
index af4adc6ba..03ad07218 100644
--- a/docs/source/conf.py
+++ b/docs/source/conf.py
@@ -322,7 +322,7 @@ def setup(app):
     # Sort gallery example by file name instead of number of lines (default)
     "within_subsection_order": FileNameSortKey,
     # directory where function granular galleries are stored
-    "backreferences_dir": None,
+    "backreferences_dir": "examples/backreferences",
     # Modules for which function level galleries are created.  In
     "doc_module": "gstools",
     # "first_notebook_cell": (
diff --git a/docs/source/tutorials.rst b/docs/source/tutorials.rst
index 3c25f597a..519b69b0b 100644
--- a/docs/source/tutorials.rst
+++ b/docs/source/tutorials.rst
@@ -30,3 +30,12 @@ explore its whole beauty and power.
 
    .. youtube:: qZBJ-AZXq6Q
       :width: 100%
+
+   |
+
+   Gallery
+   =======
+
+   .. minigallery::
+
+      ../../examples/**/*.py

From 8a97d0d5ee5fbf3f74e32eebf95b6bdecd01892c Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Sebastian=20M=C3=BCller?= <mueller.seb@posteo.de>
Date: Tue, 7 Jan 2025 11:54:24 +0100
Subject: [PATCH 44/58] typo

---
 src/gstools/covmodel/base.py | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/src/gstools/covmodel/base.py b/src/gstools/covmodel/base.py
index f633ff10b..4893b215c 100644
--- a/src/gstools/covmodel/base.py
+++ b/src/gstools/covmodel/base.py
@@ -1207,7 +1207,7 @@ class SumModel(CovModel):
     Parameters
     ----------
     *models
-        tuple of :any:`CovModel` instances of subclasses to sum.
+        tuple of :any:`CovModel` instances or subclasses to sum.
         All models will get a nugget of zero and the nugget will
         be set in the SumModel directly.
         Models need to have matching temporal setting,

From f0b113e336d7a6213bb67d4ec32870b9e47806c5 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Sebastian=20M=C3=BCller?= <mueller.seb@posteo.de>
Date: Tue, 7 Jan 2025 11:54:48 +0100
Subject: [PATCH 45/58] Docs: add first examples for sum models

---
 docs/source/conf.py                           |  2 +
 docs/source/tutorials.rst                     |  1 +
 examples/11_sum_model/00_simple_sum_model.py  | 59 ++++++++++++++
 examples/11_sum_model/01_fitting_sum_model.py | 76 +++++++++++++++++++
 examples/11_sum_model/README.rst              | 33 ++++++++
 5 files changed, 171 insertions(+)
 create mode 100644 examples/11_sum_model/00_simple_sum_model.py
 create mode 100644 examples/11_sum_model/01_fitting_sum_model.py
 create mode 100644 examples/11_sum_model/README.rst

diff --git a/docs/source/conf.py b/docs/source/conf.py
index 03ad07218..3a0e5decd 100644
--- a/docs/source/conf.py
+++ b/docs/source/conf.py
@@ -300,6 +300,7 @@ def setup(app):
         "../../examples/08_geo_coordinates/",
         "../../examples/09_spatio_temporal/",
         "../../examples/10_normalizer/",
+        "../../examples/11_sum_model/",
     ],
     # path where to save gallery generated examples
     "gallery_dirs": [
@@ -314,6 +315,7 @@ def setup(app):
         "examples/08_geo_coordinates/",
         "examples/09_spatio_temporal/",
         "examples/10_normalizer/",
+        "examples/11_sum_model/",
     ],
     # Pattern to search for example files
     "filename_pattern": r"\.py",
diff --git a/docs/source/tutorials.rst b/docs/source/tutorials.rst
index 519b69b0b..1df3aa7fc 100644
--- a/docs/source/tutorials.rst
+++ b/docs/source/tutorials.rst
@@ -22,6 +22,7 @@ explore its whole beauty and power.
    examples/08_geo_coordinates/index
    examples/09_spatio_temporal/index
    examples/10_normalizer/index
+   examples/11_sum_model/index
    examples/00_misc/index
 
 .. only:: html
diff --git a/examples/11_sum_model/00_simple_sum_model.py b/examples/11_sum_model/00_simple_sum_model.py
new file mode 100644
index 000000000..035d68571
--- /dev/null
+++ b/examples/11_sum_model/00_simple_sum_model.py
@@ -0,0 +1,59 @@
+r"""
+Creating a Sum Model
+--------------------
+
+This tutorial demonstrates how to create and use sum models in GSTools.
+We'll combine a Spherical and a Gaussian covariance model to construct
+a sum model, visualize its variogram, and generate spatial random fields.
+
+Let's start with importing GSTools setting up the domain size.
+"""
+
+import gstools as gs
+
+x = y = range(100)
+
+###############################################################################
+# First, we create two individual covariance models: a :any:`Spherical` model and a
+# :any:`Gaussian` model. The Spherical model will emphasize small-scale variability,
+# while the Gaussian model captures larger-scale patterns.
+
+m0 = gs.Spherical(dim=2, var=2.0, len_scale=5.0)
+m1 = gs.Gaussian(dim=2, var=1.0, len_scale=10.0)
+
+###############################################################################
+# Next, we create a sum model by adding these two models together.
+# Let's visualize the resulting variogram alongside the individual models.
+
+model = m0 + m1
+ax = model.plot(x_max=20)
+m0.plot(x_max=20, ax=ax)
+m1.plot(x_max=20, ax=ax)
+
+###############################################################################
+# As shown, the Spherical model controls the behavior at shorter distances,
+# while the Gaussian model dominates at longer distances.
+#
+# Using the sum model, we can generate a spatial random field. Let's visualize
+# the field created by the sum model.
+
+srf = gs.SRF(model, seed=20250107)
+srf.structured((x, y))
+srf.plot()
+
+###############################################################################
+# For comparison, we generate random fields using the individual models
+# to observe their contributions more clearly.
+
+srf0 = gs.SRF(m0, seed=20250107)
+srf0.structured((x, y))
+srf0.plot()
+
+srf1 = gs.SRF(m1, seed=20250107)
+srf1.structured((x, y))
+srf1.plot()
+
+###############################################################################
+# As seen, the Gaussian model introduces large-scale structures, while the
+# Spherical model influences the field's roughness. The sum model combines
+# these effects, resulting in a field that reflects multi-scale variability.
diff --git a/examples/11_sum_model/01_fitting_sum_model.py b/examples/11_sum_model/01_fitting_sum_model.py
new file mode 100644
index 000000000..b717107ba
--- /dev/null
+++ b/examples/11_sum_model/01_fitting_sum_model.py
@@ -0,0 +1,76 @@
+r"""
+Fitting a Sum Model
+--------------------
+
+In this tutorial, we demonstrate how to fit a sum model consisting of two
+covariance models to an empirical variogram.
+
+We will generate synthetic data, compute an empirical variogram, and fit a
+sum model combining a Spherical and Gaussian model to it.
+"""
+
+import gstools as gs
+
+x = y = range(100)
+
+###############################################################################
+# First, we create a synthetic random field based on a known sum model.
+# This will serve as the ground truth for fitting.
+
+# Define the true sum model
+m0 = gs.Spherical(dim=2, var=2.0, len_scale=5.0)
+m1 = gs.Gaussian(dim=2, var=1.0, len_scale=10.0)
+true_model = m0 + m1
+
+# Generate synthetic field
+srf = gs.SRF(true_model, seed=20250107)
+field = srf.structured((x, y))
+
+###############################################################################
+# Next, we calculate the empirical variogram from the synthetic data.
+
+# Compute empirical variogram
+bin_center, gamma = gs.vario_estimate((x, y), field)
+
+###############################################################################
+# Now we define a sum model to fit to the empirical variogram.
+# Initially, the parameters of the models are arbitrary.
+#
+# A sum model can also be created by a list of model classes together with
+# the common arguments (like dim in this case).
+
+fit_model = gs.SumModel(gs.Spherical, gs.Gaussian, dim=2)
+
+###############################################################################
+# We fit the sum model to the empirical variogram using GSTools' built-in
+# fitting capabilities. As seen in the representation, the variances and length
+# scales of the individual models can be accessed by the attributes
+# :any:`SumModel.vars` and :any:`SumModel.len_scales`.
+
+fit_model.fit_variogram(bin_center, gamma)
+print(f"{true_model=}")
+print(f" {fit_model=}")
+
+###############################################################################
+# The variance of a sum model is the sum of the sub variances
+# from the contained models. The length scale is a weighted sum of the sub
+# length scales where the weights are the ratios of the sub variances
+# to the total variance of the sum model.
+
+print(f"{true_model.var=:.2}, {true_model.len_scale=:.2}")
+print(f" {fit_model.var=:.2},  {fit_model.len_scale=:.2}")
+
+###############################################################################
+# After fitting, we can visualize the empirical variogram alongside the
+# fitted sum model and its components. A Sum Model is subscriptable to access
+# the individual models its contains.
+
+ax = fit_model.plot(x_max=max(bin_center))
+ax.scatter(bin_center, gamma)
+# Extract individual components
+fit_model[0].plot(x_max=max(bin_center), ax=ax)
+fit_model[1].plot(x_max=max(bin_center), ax=ax)
+
+###############################################################################
+# As we can see, the fitted sum model closely matches the empirical variogram,
+# demonstrating its ability to capture multi-scale variability effectively.
diff --git a/examples/11_sum_model/README.rst b/examples/11_sum_model/README.rst
new file mode 100644
index 000000000..d50c15314
--- /dev/null
+++ b/examples/11_sum_model/README.rst
@@ -0,0 +1,33 @@
+Summing Covariance Models
+=========================
+
+In geostatistics, the spatial variability of natural phenomena is often represented using covariance models,
+which describe how values of a property correlate over distance.
+A single covariance model may capture specific features of spatial correlation, such as smoothness or the range of influence.
+However, many real-world spatial processes are complex, involving multiple overlapping structures
+that cannot be adequately described by a single covariance model.
+
+This is where **sum models** come into play.
+A sum model combines multiple covariance models into a single representation,
+allowing for a more flexible and comprehensive description of spatial variability.
+By summing covariance models, we can:
+
+1. **Capture Multi-Scale Variability:** Many spatial phenomena exhibit variability at different scales.
+   For example, soil properties may have small-scale variation due to local heterogeneities and large-scale variation due to regional trends.
+   A sum model can combine short-range and long-range covariance models to reflect this behavior.
+2. **Incorporate Multiple Physical Processes:** Different processes may contribute to the observed spatial pattern.
+   For instance, in hydrology, the spatial distribution of groundwater levels might be influenced by both geological structures and human activities.
+   A sum model can represent these independent contributions.
+3. **Improve Model Fit and Prediction Accuracy:** By combining models, sum models can better match empirical variograms or other observed data,
+   leading to more accurate predictions in kriging or simulation tasks.
+4. **Enhance Interpretability:** Each component of a sum model can be associated with a specific spatial process or scale,
+   providing insights into the underlying mechanisms driving spatial variability.
+
+The new :any:`SumModel` introduced in GSTools makes it straightforward to define and work with such composite covariance structures.
+It allows users to combine any number of base models, each with its own parameters, in a way that is both intuitive and computationally efficient.
+
+In the following tutorials, we'll explore how to use the :any:`SumModel` in GSTools,
+including practical examples that demonstrate its utility in real-world scenarios.
+
+Examples
+--------

From af06e4961d901e6e0ad38fb5b7d106d817f1f154 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Sebastian=20M=C3=BCller?= <mueller.seb@posteo.de>
Date: Fri, 7 Feb 2025 22:15:04 +0100
Subject: [PATCH 46/58] SumModel: copy CovModel instances to prevent strange
 errors

---
 src/gstools/covmodel/base.py | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/src/gstools/covmodel/base.py b/src/gstools/covmodel/base.py
index 4893b215c..544aacf2f 100644
--- a/src/gstools/covmodel/base.py
+++ b/src/gstools/covmodel/base.py
@@ -1310,7 +1310,7 @@ def __init__(self, *models, **kwargs):
                 if to_init is not None and to_init:
                     raise ValueError(imsg)
                 to_init = False
-                self._models.append(mod)
+                self._models.append(copy.deepcopy(mod))
             elif isinstance(mod, type) and issubclass(mod, CovModel):
                 if to_init is not None and not to_init:
                     raise ValueError(imsg)

From fb0a01e1128af138b6a296e2144172a3a2da33d3 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Sebastian=20M=C3=BCller?= <mueller.seb@posteo.de>
Date: Fri, 7 Feb 2025 23:11:51 +0100
Subject: [PATCH 47/58] SumModel: also copy models of other summodels

---
 src/gstools/covmodel/base.py | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/src/gstools/covmodel/base.py b/src/gstools/covmodel/base.py
index 544aacf2f..f3453d859 100644
--- a/src/gstools/covmodel/base.py
+++ b/src/gstools/covmodel/base.py
@@ -1304,7 +1304,7 @@ def __init__(self, *models, **kwargs):
                 if to_init is not None and to_init:
                     raise ValueError(imsg)
                 to_init = False
-                self._models += mod.models
+                self._models += copy.deepcopy(mod.models)
                 add_nug += mod.nugget
             elif isinstance(mod, CovModel):
                 if to_init is not None and to_init:

From e8ea93c050ba1fe4d234a5035bbcd1d4c96fa635 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Sebastian=20M=C3=BCller?= <mueller.seb@posteo.de>
Date: Fri, 7 Feb 2025 23:17:43 +0100
Subject: [PATCH 48/58] make sum-models comparable

---
 src/gstools/covmodel/base.py      |  9 +++++++++
 src/gstools/covmodel/sum_tools.py | 18 ++++++++++++++++++
 2 files changed, 27 insertions(+)

diff --git a/src/gstools/covmodel/base.py b/src/gstools/covmodel/base.py
index f3453d859..c0aa1ea32 100644
--- a/src/gstools/covmodel/base.py
+++ b/src/gstools/covmodel/base.py
@@ -22,6 +22,7 @@
 from gstools.covmodel.sum_tools import (
     default_mod_kwargs,
     sum_check,
+    sum_compare,
     sum_default_arg_bounds,
     sum_default_opt_arg_bounds,
     sum_model_repr,
@@ -1166,6 +1167,8 @@ def __eq__(self, other):
         """Compare CovModels."""
         if not isinstance(other, CovModel):
             return False
+        if isinstance(other, SumModel):
+            return False
         return compare(self, other)
 
     def __setattr__(self, name, value):
@@ -1738,6 +1741,12 @@ def __getattr__(self, name):
                         return getattr(mod, opt)
         raise AttributeError(f"'{self.name}' object has no attribute '{name}'")
 
+    def __eq__(self, other):
+        """Compare SumModels."""
+        if not isinstance(other, SumModel):
+            return False
+        return sum_compare(self, other)
+
     def __repr__(self):
         """Return String representation."""
         return sum_model_repr(self)
diff --git a/src/gstools/covmodel/sum_tools.py b/src/gstools/covmodel/sum_tools.py
index 7e6c6b4f3..54e0533b1 100644
--- a/src/gstools/covmodel/sum_tools.py
+++ b/src/gstools/covmodel/sum_tools.py
@@ -10,6 +10,7 @@
    ARG_DEF
    default_mod_kwargs
    sum_check
+   sum_compare
    sum_default_arg_bounds
    sum_default_opt_arg_bounds
    sum_set_norm_var_ratios
@@ -28,6 +29,7 @@
     "ARG_DEF",
     "default_mod_kwargs",
     "sum_check",
+    "sum_compare",
     "sum_default_arg_bounds",
     "sum_default_opt_arg_bounds",
     "sum_set_var_weights",
@@ -210,6 +212,22 @@ def sum_set_len_weights(summod, weights, skip=None, len_scale=None):
     summod.len_scales = len_scales
 
 
+def sum_compare(this, that):
+    """
+    Compare SumModels.
+
+    Parameters
+    ----------
+    this / that : :any:`SumModel`
+        The sum models to compare.
+    """
+    if len(this) != len(that):
+        return False
+    if not np.isclose(this.nugget, that.nugget):
+        return False
+    return all([mod1 == mod2 for (mod1, mod2) in zip(this, that)])
+
+
 def sum_model_repr(summod):  # pragma: no cover
     """
     Generate the sum-model string representation.

From c64af49121563a31a5ad6390313c3eff0abd40b6 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Sebastian=20M=C3=BCller?= <mueller.seb@posteo.de>
Date: Sat, 8 Feb 2025 13:36:54 +0100
Subject: [PATCH 49/58] SumModel: remove RatioError, use ValueError instead

---
 src/gstools/covmodel/sum_tools.py | 10 ++--------
 1 file changed, 2 insertions(+), 8 deletions(-)

diff --git a/src/gstools/covmodel/sum_tools.py b/src/gstools/covmodel/sum_tools.py
index 54e0533b1..1f288244a 100644
--- a/src/gstools/covmodel/sum_tools.py
+++ b/src/gstools/covmodel/sum_tools.py
@@ -6,7 +6,6 @@
 The following classes and functions are provided
 
 .. autosummary::
-   RatioError
    ARG_DEF
    default_mod_kwargs
    sum_check
@@ -25,7 +24,6 @@
 from gstools.tools.misc import list_format
 
 __all__ = [
-    "RatioError",
     "ARG_DEF",
     "default_mod_kwargs",
     "sum_check",
@@ -38,10 +36,6 @@
 ]
 
 
-class RatioError(Exception):
-    """Error for invalid ratios in SumModel."""
-
-
 ARG_DEF = {
     "dim": 3,
     "latlon": False,
@@ -145,7 +139,7 @@ def sum_set_var_weights(summod, weights, skip=None, var=None):
     var_diff = var - var_sum
     if var_diff < 0:
         msg = "SumModel.set_var_weights: skipped variances to big."
-        raise RatioError(msg)
+        raise ValueError(msg)
     weights_sum = sum(weights)
     var_list = summod.vars
     j = 0
@@ -200,7 +194,7 @@ def sum_set_len_weights(summod, weights, skip=None, len_scale=None):
     len_diff = len_scale - len_sum
     if len_diff < 0:
         msg = "SumModel.set_len_weights: skipped length scales to big."
-        raise RatioError(msg)
+        raise ValueError(msg)
     weights_sum = sum(weights)
     len_scales = summod.len_scales
     j = 0

From c885952ed9afbf472bf722bf8f87d552d2e34414 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Sebastian=20M=C3=BCller?= <mueller.seb@posteo.de>
Date: Sat, 8 Feb 2025 13:37:06 +0100
Subject: [PATCH 50/58] SumModel: add tests

---
 tests/data/variogram.txt |   2 +
 tests/test_sum_model.py  | 179 +++++++++++++++++++++++++++++++++++++++
 2 files changed, 181 insertions(+)
 create mode 100644 tests/data/variogram.txt
 create mode 100644 tests/test_sum_model.py

diff --git a/tests/data/variogram.txt b/tests/data/variogram.txt
new file mode 100644
index 000000000..cad2a2be3
--- /dev/null
+++ b/tests/data/variogram.txt
@@ -0,0 +1,2 @@
+2.500000000000000000e-01 7.500000000000000000e-01 1.250000000000000000e+00 1.750000000000000000e+00 2.250000000000000000e+00 2.750000000000000000e+00 3.250000000000000000e+00 3.750000000000000000e+00 4.250000000000000000e+00 4.750000000000000000e+00 5.250000000000000000e+00 5.750000000000000000e+00 6.250000000000000000e+00 6.750000000000000000e+00 7.250000000000000000e+00 7.750000000000000000e+00 8.250000000000000000e+00 8.750000000000000000e+00 9.250000000000000000e+00 9.750000000000000000e+00 1.025000000000000000e+01 1.075000000000000000e+01 1.125000000000000000e+01 1.175000000000000000e+01 1.225000000000000000e+01 1.275000000000000000e+01 1.325000000000000000e+01 1.375000000000000000e+01 1.425000000000000000e+01 1.475000000000000000e+01 1.525000000000000000e+01 1.575000000000000000e+01 1.625000000000000000e+01 1.675000000000000000e+01 1.725000000000000000e+01 1.775000000000000000e+01 1.825000000000000000e+01 1.875000000000000000e+01 1.925000000000000000e+01 1.975000000000000000e+01 2.025000000000000000e+01 2.075000000000000000e+01 2.125000000000000000e+01 2.175000000000000000e+01 2.225000000000000000e+01 2.275000000000000000e+01 2.325000000000000000e+01 2.375000000000000000e+01 2.425000000000000000e+01 2.475000000000000000e+01 2.525000000000000000e+01 2.575000000000000000e+01 2.625000000000000000e+01 2.675000000000000000e+01 2.725000000000000000e+01 2.775000000000000000e+01 2.825000000000000000e+01 2.875000000000000000e+01 2.925000000000000000e+01 2.975000000000000000e+01 3.025000000000000000e+01 3.075000000000000000e+01 3.125000000000000000e+01 3.175000000000000000e+01 3.225000000000000000e+01 3.275000000000000000e+01 3.325000000000000000e+01 3.375000000000000000e+01 3.425000000000000000e+01 3.475000000000000000e+01 3.525000000000000000e+01 3.575000000000000000e+01 3.625000000000000000e+01 3.675000000000000000e+01 3.725000000000000000e+01 3.775000000000000000e+01 3.825000000000000000e+01 3.875000000000000000e+01 3.925000000000000000e+01 3.975000000000000000e+01 4.025000000000000000e+01 4.075000000000000000e+01 4.125000000000000000e+01 4.175000000000000000e+01 4.225000000000000000e+01 4.275000000000000000e+01 4.325000000000000000e+01 4.375000000000000000e+01 4.425000000000000000e+01 4.475000000000000000e+01 4.525000000000000000e+01 4.575000000000000000e+01 4.625000000000000000e+01 4.675000000000000000e+01 4.725000000000000000e+01 4.775000000000000000e+01 4.825000000000000000e+01 4.875000000000000000e+01 4.925000000000000000e+01 4.975000000000000000e+01
+9.341287109808403821e-01 8.453134048372161757e-01 1.097876163640169045e+00 1.452675652987244215e+00 1.717011001508209400e+00 1.982144808494284316e+00 2.116570368080190168e+00 2.354956050093370834e+00 2.256700499821684858e+00 2.608287292075601727e+00 3.019563428869475263e+00 3.005582884425558987e+00 3.026403707124349474e+00 2.920331468705290145e+00 2.739513460596619598e+00 2.983946141180413569e+00 3.002935993020715788e+00 3.098828622420069401e+00 3.252632934536162423e+00 2.972322302660843629e+00 3.162730649439238206e+00 3.104794978590492249e+00 3.129840716346046658e+00 3.374897140494525605e+00 3.281178715947507207e+00 3.414929159332046993e+00 3.366128007332230609e+00 3.350479282803380254e+00 3.274833816819866961e+00 3.146498371924009607e+00 3.098730655077019946e+00 3.203013412974681895e+00 3.386168318931909393e+00 3.349129048195039093e+00 3.403803954611331228e+00 3.228624506901822677e+00 3.320464283280209816e+00 3.312375916258147424e+00 3.476587847562206068e+00 3.293083753810556846e+00 3.526290649886577366e+00 3.576052045142498859e+00 3.440057124681327405e+00 3.388464225293542853e+00 3.350244356710671667e+00 3.502962399475493704e+00 3.519679087894428626e+00 3.563334300898933993e+00 3.540069270604898843e+00 3.743566023712460389e+00 3.586432982353473964e+00 3.820556981535015773e+00 3.821630249625350473e+00 3.674450536968612013e+00 3.935179864707886388e+00 3.829663738092724312e+00 3.767186362291683466e+00 3.829516094024750128e+00 3.939988520639752245e+00 3.788358175483068191e+00 3.929699245388040385e+00 3.794618804853678196e+00 3.732527204313444535e+00 3.691846735135440483e+00 3.504551184978119682e+00 3.866643330145578261e+00 3.768275820398073073e+00 3.645363755293978603e+00 3.835573383488326105e+00 3.811532952179032208e+00 3.609181746819485781e+00 3.652959035004343491e+00 3.766941214399454285e+00 3.884337832265876589e+00 3.818059727339796705e+00 3.957999650541710324e+00 3.862889397483654896e+00 3.649107551259477056e+00 3.703662969618971346e+00 3.775887102684075014e+00 3.666801729040440438e+00 3.798614706003481167e+00 3.592569765506981483e+00 3.617445665128535826e+00 3.746142254086667567e+00 3.634920395024769935e+00 3.685908603892140256e+00 3.552532538928491412e+00 3.704709972576425869e+00 3.510318544036457933e+00 3.557027380970486874e+00 3.431417669026475270e+00 3.625033767218668324e+00 3.646369536691860258e+00 3.645946180745845311e+00 3.629191872454593071e+00 3.668750066947414457e+00 3.626446691977532222e+00 3.651811593215607665e+00 3.623437445821843017e+00
diff --git a/tests/test_sum_model.py b/tests/test_sum_model.py
new file mode 100644
index 000000000..c36f2ffb2
--- /dev/null
+++ b/tests/test_sum_model.py
@@ -0,0 +1,179 @@
+"""
+This is the unittest of the SumModel class.
+"""
+
+import unittest
+from pathlib import Path
+
+import numpy as np
+
+import gstools as gs
+
+
+class TestSumModel(unittest.TestCase):
+    def test_init(self):
+        s1 = gs.SumModel(dim=3)
+        s2 = gs.Nugget(dim=3)
+        self.assertTrue(s1 == s2)
+
+        m1 = gs.Nugget(nugget=1)
+        s1 = m1 + m1 + m1
+        s2 = gs.Nugget(nugget=3)
+        self.assertTrue(s1 == s2)
+        self.assertFalse(m1 == s2)
+
+        s1 = gs.Gaussian(var=1) + gs.Exponential(var=2)
+        s2 = gs.Exponential(var=1) + gs.Gaussian(var=2)
+        self.assertFalse(s1 == gs.Nugget(dim=3))
+        self.assertFalse(s1 == s2)
+
+        # check that models get copied
+        m1 = gs.Gaussian()
+        s1 = m1 + m1
+        var = [1.0, 2.0]
+        s1.vars = var
+        np.testing.assert_array_almost_equal(s1.vars, var)
+
+        m1 = gs.Gaussian(dim=2, var=1.0, len_scale=1.0)
+        m2 = gs.Matern(dim=2, var=2.0, len_scale=2.0, nu=2.0)
+        m3 = gs.Integral(dim=2, var=3.0, len_scale=3.0, nu=3.0)
+        s1 = gs.SumModel(m1, m2, m3)
+        s2 = m1 + m2 + m3
+        s3 = gs.SumModel(
+            *(gs.Gaussian, gs.Matern, gs.Integral),
+            dim=2,
+            vars=[1.0, 2.0, 3.0],
+            len_scales=[1.0, 2.0, 3.0],
+            nu_1=2.0,
+            nu_2=3.0,
+        )
+        self.assertTrue(s1 == s2)
+        self.assertTrue(s2 == s3)
+        self.assertTrue(s1 == s3)
+
+        with self.assertRaises(ValueError):
+            gs.SumModel(gs.Exponential, gs.Gaussian(dim=2))
+
+        with self.assertRaises(ValueError):
+            gs.SumModel(gs.Nugget, gs.Nugget(nugget=2))
+
+        with self.assertRaises(ValueError):
+            model = gs.Spherical() + gs.Spherical()
+            model[0].dim = 2
+            model.check()
+
+        with self.assertRaises(ValueError):
+            model = gs.Spherical() + gs.Spherical()
+            model[0].nugget = 2
+            model.check()
+
+        with self.assertRaises(ValueError):
+            gs.Spherical(latlon=True) + gs.Spherical()
+
+        with self.assertRaises(ValueError):
+            gs.Spherical(geo_scale=2) + gs.Spherical()
+
+        with self.assertRaises(ValueError):
+            gs.Spherical(temporal=True) + gs.Spherical()
+
+        with self.assertRaises(ValueError):
+            gs.Spherical(anis=0.5) + gs.Spherical()
+
+        with self.assertRaises(ValueError):
+            gs.Spherical(angles=0.5) + gs.Spherical()
+
+    def test_generate(self):
+        x = np.random.RandomState(19970221).rand(1000) * 100.0
+        y = np.random.RandomState(20011012).rand(1000) * 100.0
+        m1 = gs.Spherical(dim=2, var=2, len_scale=5)
+        m2 = gs.Spherical(dim=2, var=1, len_scale=10)
+        m3 = gs.Gaussian(dim=2, var=1, len_scale=20)
+        model = m1 + m2 + m3
+        srf = gs.SRF(model, mean=0, seed=199702212)
+        field = srf((x, y))
+        self.assertAlmostEqual(np.var(field), 3.7, places=1)
+        # used for test_fit (see below)
+        # bin_center, gamma = gs.vario_estimate((x, y), field, max_dist=50, bin_no=100)
+
+    def test_fit(self):
+        here = Path(__file__).parent
+        bin_center, gamma = np.loadtxt(here / "data" / "variogram.txt")
+        s2 = gs.SumModel(gs.Gaussian, gs.Spherical, gs.Spherical, dim=2)
+        res1, _ = s2.fit_variogram(bin_center, gamma, nugget=False)
+        res2, _ = s2.fit_variogram(
+            bin_center, gamma, nugget=False, len_scale_2=5, var_0=1
+        )
+        res3, _ = s2.fit_variogram(bin_center, gamma, nugget=False, var=3.7)
+        res4, _ = s2.fit_variogram(
+            bin_center, gamma, nugget=False, len_scale=15
+        )
+        res5, _ = s2.fit_variogram(
+            bin_center, gamma, len_scale=15, sill=3.7, var_1=1
+        )
+        res6, _ = s2.fit_variogram(
+            bin_center, gamma, len_scale=15, sill=3.7, var_1=1, nugget=0.1
+        )
+        res7, _ = s2.fit_variogram(
+            bin_center,
+            gamma,
+            len_scale=15,
+            sill=3.7,
+            var_1=1,
+            var_2=1,
+            nugget=0.1,
+        )
+
+        self.assertAlmostEqual(res1["var"], 3.7, places=1)
+        self.assertAlmostEqual(res2["var"], 3.7, places=1)
+        self.assertAlmostEqual(res3["var"], 3.7, places=5)
+        self.assertAlmostEqual(res4["len_scale"], 15.0, places=5)
+        self.assertAlmostEqual(res5["var"] + res5["nugget"], 3.7, places=2)
+        self.assertAlmostEqual(res5["len_scale"], 15.0, places=5)
+        self.assertAlmostEqual(res6["var"] + res6["nugget"], 3.7, places=2)
+        self.assertAlmostEqual(res7["var_0"], 1.6, places=5)
+
+        mod_n = gs.Nugget(dim=2)
+        res, _ = mod_n.fit_variogram(bin_center, gamma)
+        self.assertAlmostEqual(res["nugget"], 3.4, places=1)
+        # nothing to fit
+        with self.assertRaises(ValueError):
+            mod_n.fit_variogram(bin_center, gamma, nugget=False)
+        # fixed sub-vars greated than fixed total var
+        with self.assertRaises(ValueError):
+            s2.fit_variogram(bin_center, gamma, var=1, var_0=1, var_1=1)
+        # fixed len_scale and sub-len-scale not possible
+        with self.assertRaises(ValueError):
+            s2.fit_variogram(bin_center, gamma, len_scale=15, len_scale_0=5)
+
+    def test_sum_weights(self):
+        mod = gs.Gaussian() + gs.Gaussian() + gs.Gaussian()
+        # fixed sub-vars too big
+        with self.assertRaises(ValueError):
+            mod.set_var_weights([1], skip=[1, 2], var=1)
+        # too many ids
+        with self.assertRaises(ValueError):
+            mod.set_var_weights([1, 1], skip=[1, 2])
+        # wrong skip
+        with self.assertRaises(ValueError):
+            mod.set_var_weights([1, 1], skip=[10])
+
+        mod = gs.Gaussian() + gs.Gaussian() + gs.Gaussian()
+        # fixed sub-lens too big
+        with self.assertRaises(ValueError):
+            mod.set_len_weights([1], skip=[1, 2], len_scale=0.1)
+        # too many ids
+        with self.assertRaises(ValueError):
+            mod.set_len_weights([1, 1], skip=[1, 2])
+        # wrong skip
+        with self.assertRaises(ValueError):
+            mod.set_len_weights([1, 1], skip=[10])
+        # check setting with skipping
+        mod.set_len_weights([1, 1], skip=[0], len_scale=10)
+        np.testing.assert_array_almost_equal(mod.len_scales, [1, 14.5, 14.5])
+        # check setting
+        mod.set_len_weights([1, 1, 2], len_scale=10)
+        np.testing.assert_array_almost_equal(mod.len_scales, [7.5, 7.5, 15])
+
+
+if __name__ == "__main__":
+    unittest.main()

From c66d72283fcc136e2e63e8cb19def31a22bb9522 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Sebastian=20M=C3=BCller?= <mueller.seb@posteo.de>
Date: Sat, 8 Feb 2025 13:37:19 +0100
Subject: [PATCH 51/58] add test data to manifest

---
 MANIFEST.in | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/MANIFEST.in b/MANIFEST.in
index 24184482a..7d3d0bb9e 100644
--- a/MANIFEST.in
+++ b/MANIFEST.in
@@ -1,4 +1,4 @@
 prune **
-recursive-include tests *.py
+recursive-include tests *.py *.txt
 recursive-include src/gstools *.py *.pyx
 include AUTHORS.md LICENSE README.md pyproject.toml setup.py

From fedb294ba8018ea9cb67edb3313050f6538105f1 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Sebastian=20M=C3=BCller?= <mueller.seb@posteo.de>
Date: Sat, 8 Feb 2025 14:04:08 +0100
Subject: [PATCH 52/58] pylint fix

---
 src/gstools/covmodel/sum_tools.py | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/src/gstools/covmodel/sum_tools.py b/src/gstools/covmodel/sum_tools.py
index 1f288244a..74604d381 100644
--- a/src/gstools/covmodel/sum_tools.py
+++ b/src/gstools/covmodel/sum_tools.py
@@ -219,7 +219,7 @@ def sum_compare(this, that):
         return False
     if not np.isclose(this.nugget, that.nugget):
         return False
-    return all([mod1 == mod2 for (mod1, mod2) in zip(this, that)])
+    return all(mod1 == mod2 for (mod1, mod2) in zip(this, that))
 
 
 def sum_model_repr(summod):  # pragma: no cover

From 48ec88cd6eae2c14551b4d4cef28a4ed69da4d1b Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Sebastian=20M=C3=BCller?= <mueller.seb@posteo.de>
Date: Sat, 8 Feb 2025 14:04:16 +0100
Subject: [PATCH 53/58] SumModel: more tests

---
 tests/test_sum_model.py | 55 +++++++++++++++++++++++++++++++++++++++++
 1 file changed, 55 insertions(+)

diff --git a/tests/test_sum_model.py b/tests/test_sum_model.py
index c36f2ffb2..faba3be22 100644
--- a/tests/test_sum_model.py
+++ b/tests/test_sum_model.py
@@ -22,10 +22,40 @@ def test_init(self):
         self.assertTrue(s1 == s2)
         self.assertFalse(m1 == s2)
 
+        # Nugget cant set positive var or len-scale
+        with self.assertRaises(ValueError):
+            m1.var = 1
+        with self.assertRaises(ValueError):
+            m1.len_scale = 1
+
+        s1 = gs.Gaussian() + gs.Gaussian() + gs.Gaussian()
+        with self.assertRaises(ValueError):
+            s1.vars = [1, 2]
+        with self.assertRaises(ValueError):
+            s1.len_scales = [1, 2]
+
+        s1.integral_scale = 10
+        self.assertAlmostEqual(s1.integral_scale_0, 10)
+        self.assertAlmostEqual(s1.integral_scale_1, 10)
+        self.assertAlmostEqual(s1.integral_scale_2, 10)
+
+        s1.var = 2
+        s1.ratios = [0.2, 0.2, 0.6]
+        self.assertAlmostEqual(s1.vars[0], 0.4)
+        self.assertAlmostEqual(s1.vars[1], 0.4)
+        self.assertAlmostEqual(s1.vars[2], 1.2)
+
+        with self.assertRaises(ValueError):
+            s1.ratios = [0.3, 0.2, 0.6]
+        with self.assertRaises(ValueError):
+            s1.ratios = [0.3, 0.2]
+
         s1 = gs.Gaussian(var=1) + gs.Exponential(var=2)
         s2 = gs.Exponential(var=1) + gs.Gaussian(var=2)
         self.assertFalse(s1 == gs.Nugget(dim=3))
         self.assertFalse(s1 == s2)
+        self.assertFalse(gs.Exponential() == (gs.Exponential() + gs.Nugget()))
+        self.assertFalse((gs.Exponential() + gs.Nugget()) == gs.Exponential())
 
         # check that models get copied
         m1 = gs.Gaussian()
@@ -34,6 +64,19 @@ def test_init(self):
         s1.vars = var
         np.testing.assert_array_almost_equal(s1.vars, var)
 
+        s1 = gs.SumModel(gs.Exponential, gs.Exponential, var=3)
+        np.testing.assert_array_almost_equal(s1.vars, [1.5, 1.5])
+        self.assertFalse(gs.Gaussian() in s1)
+
+        s1 = gs.SumModel(gs.Exponential, gs.Exponential, len_scale=10)
+        np.testing.assert_array_almost_equal(s1.len_scales, [10, 10])
+
+        s1 = gs.SumModel(
+            gs.Exponential, gs.Exponential, temporal=True, spatial_dim=2
+        )
+        self.assertTrue(all(mod.temporal for mod in s1))
+        self.assertTrue(all(mod.dim == 3 for mod in s1))
+
         m1 = gs.Gaussian(dim=2, var=1.0, len_scale=1.0)
         m2 = gs.Matern(dim=2, var=2.0, len_scale=2.0, nu=2.0)
         m3 = gs.Integral(dim=2, var=3.0, len_scale=3.0, nu=3.0)
@@ -51,12 +94,24 @@ def test_init(self):
         self.assertTrue(s2 == s3)
         self.assertTrue(s1 == s3)
 
+        with self.assertRaises(ValueError):
+            gs.Exponential() + 1
+
+        with self.assertRaises(ValueError):
+            1 + gs.Exponential()
+
         with self.assertRaises(ValueError):
             gs.SumModel(gs.Exponential, gs.Gaussian(dim=2))
 
         with self.assertRaises(ValueError):
             gs.SumModel(gs.Nugget, gs.Nugget(nugget=2))
 
+        with self.assertRaises(ValueError):
+            gs.SumModel(gs.Gaussian(dim=2), gs.Exponential)
+
+        with self.assertRaises(ValueError):
+            gs.SumModel(gs.Nugget(nugget=2), gs.Nugget)
+
         with self.assertRaises(ValueError):
             model = gs.Spherical() + gs.Spherical()
             model[0].dim = 2

From 4010e866a54f473cd3f92937a0aed315fa910303 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Sebastian=20M=C3=BCller?= <mueller.seb@posteo.de>
Date: Tue, 22 Apr 2025 16:43:24 +0200
Subject: [PATCH 54/58] update sum-model examples after review

---
 examples/11_sum_model/00_simple_sum_model.py  | 10 ++++++----
 examples/11_sum_model/01_fitting_sum_model.py | 20 +++++++++++--------
 examples/11_sum_model/README.rst              | 13 +++++-------
 3 files changed, 23 insertions(+), 20 deletions(-)

diff --git a/examples/11_sum_model/00_simple_sum_model.py b/examples/11_sum_model/00_simple_sum_model.py
index 035d68571..a9ded98ed 100644
--- a/examples/11_sum_model/00_simple_sum_model.py
+++ b/examples/11_sum_model/00_simple_sum_model.py
@@ -6,7 +6,7 @@
 We'll combine a Spherical and a Gaussian covariance model to construct
 a sum model, visualize its variogram, and generate spatial random fields.
 
-Let's start with importing GSTools setting up the domain size.
+Let's start with importing GSTools and setting up the domain size.
 """
 
 import gstools as gs
@@ -15,8 +15,9 @@
 
 ###############################################################################
 # First, we create two individual covariance models: a :any:`Spherical` model and a
-# :any:`Gaussian` model. The Spherical model will emphasize small-scale variability,
-# while the Gaussian model captures larger-scale patterns.
+# :any:`Gaussian` model. The Spherical model with its short length scale
+# will emphasize small-scale variability, while the Gaussian model with a larger length scale
+# captures larger-scale patterns.
 
 m0 = gs.Spherical(dim=2, var=2.0, len_scale=5.0)
 m1 = gs.Gaussian(dim=2, var=1.0, len_scale=10.0)
@@ -32,7 +33,8 @@
 
 ###############################################################################
 # As shown, the Spherical model controls the behavior at shorter distances,
-# while the Gaussian model dominates at longer distances.
+# while the Gaussian model dominates at longer distances. The ratio of influence
+# is thereby controlled by the provided variances of the individual models.
 #
 # Using the sum model, we can generate a spatial random field. Let's visualize
 # the field created by the sum model.
diff --git a/examples/11_sum_model/01_fitting_sum_model.py b/examples/11_sum_model/01_fitting_sum_model.py
index b717107ba..7a0dfa294 100644
--- a/examples/11_sum_model/01_fitting_sum_model.py
+++ b/examples/11_sum_model/01_fitting_sum_model.py
@@ -23,7 +23,7 @@
 true_model = m0 + m1
 
 # Generate synthetic field
-srf = gs.SRF(true_model, seed=20250107)
+srf = gs.SRF(true_model, seed=20250405)
 field = srf.structured((x, y))
 
 ###############################################################################
@@ -43,19 +43,18 @@
 
 ###############################################################################
 # We fit the sum model to the empirical variogram using GSTools' built-in
-# fitting capabilities. As seen in the representation, the variances and length
-# scales of the individual models can be accessed by the attributes
-# :any:`SumModel.vars` and :any:`SumModel.len_scales`.
+# fitting capabilities. We deactivate the nugget fitting to not overparameterize
+# our model.
 
-fit_model.fit_variogram(bin_center, gamma)
+fit_model.fit_variogram(bin_center, gamma, nugget=False)
 print(f"{true_model=}")
 print(f" {fit_model=}")
 
 ###############################################################################
-# The variance of a sum model is the sum of the sub variances
+# The variance of a sum model is the sum of the sub variances (:any:`SumModel.vars`)
 # from the contained models. The length scale is a weighted sum of the sub
-# length scales where the weights are the ratios of the sub variances
-# to the total variance of the sum model.
+# length scales (:any:`SumModel.len_scales`) where the weights are the ratios
+# of the sub variances to the total variance of the sum model.
 
 print(f"{true_model.var=:.2}, {true_model.len_scale=:.2}")
 print(f" {fit_model.var=:.2},  {fit_model.len_scale=:.2}")
@@ -70,7 +69,12 @@
 # Extract individual components
 fit_model[0].plot(x_max=max(bin_center), ax=ax)
 fit_model[1].plot(x_max=max(bin_center), ax=ax)
+# True models
+true_model.plot(x_max=max(bin_center), ax=ax, ls="--", c="C0", label="")
+true_model[0].plot(x_max=max(bin_center), ax=ax, ls="--", c="C1", label="")
+true_model[1].plot(x_max=max(bin_center), ax=ax, ls="--", c="C2", label="")
 
 ###############################################################################
 # As we can see, the fitted sum model closely matches the empirical variogram,
 # demonstrating its ability to capture multi-scale variability effectively.
+# The "true" variograms are shown with dashed lines for comparison.
diff --git a/examples/11_sum_model/README.rst b/examples/11_sum_model/README.rst
index d50c15314..c4ccb2306 100644
--- a/examples/11_sum_model/README.rst
+++ b/examples/11_sum_model/README.rst
@@ -1,9 +1,9 @@
 Summing Covariance Models
 =========================
 
-In geostatistics, the spatial variability of natural phenomena is often represented using covariance models,
+In geostatistics, the spatial relations of natural phenomena is often represented using covariance models,
 which describe how values of a property correlate over distance.
-A single covariance model may capture specific features of spatial correlation, such as smoothness or the range of influence.
+A single covariance model may capture specific features of the spatial correlation, such as smoothness or the range of influence.
 However, many real-world spatial processes are complex, involving multiple overlapping structures
 that cannot be adequately described by a single covariance model.
 
@@ -15,19 +15,16 @@ By summing covariance models, we can:
 1. **Capture Multi-Scale Variability:** Many spatial phenomena exhibit variability at different scales.
    For example, soil properties may have small-scale variation due to local heterogeneities and large-scale variation due to regional trends.
    A sum model can combine short-range and long-range covariance models to reflect this behavior.
-2. **Incorporate Multiple Physical Processes:** Different processes may contribute to the observed spatial pattern.
-   For instance, in hydrology, the spatial distribution of groundwater levels might be influenced by both geological structures and human activities.
-   A sum model can represent these independent contributions.
-3. **Improve Model Fit and Prediction Accuracy:** By combining models, sum models can better match empirical variograms or other observed data,
+2. **Improve Model Fit and Prediction Accuracy:** By combining models, sum models can better match empirical variograms or other observed data,
    leading to more accurate predictions in kriging or simulation tasks.
-4. **Enhance Interpretability:** Each component of a sum model can be associated with a specific spatial process or scale,
+3. **Enhance Interpretability:** Each component of a sum model can be associated with a specific spatial process or scale,
    providing insights into the underlying mechanisms driving spatial variability.
 
 The new :any:`SumModel` introduced in GSTools makes it straightforward to define and work with such composite covariance structures.
 It allows users to combine any number of base models, each with its own parameters, in a way that is both intuitive and computationally efficient.
 
 In the following tutorials, we'll explore how to use the :any:`SumModel` in GSTools,
-including practical examples that demonstrate its utility in real-world scenarios.
+including practical examples that demonstrate its utility in different scenarios.
 
 Examples
 --------

From 342efbec8ed238f3ab1ca4176202266e57327ae4 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Sebastian=20M=C3=BCller?= <mueller.seb@posteo.de>
Date: Tue, 22 Apr 2025 20:58:12 +0200
Subject: [PATCH 55/58] tests: better structure sum-model tests

---
 tests/test_sum_model.py | 11 +++++++++++
 1 file changed, 11 insertions(+)

diff --git a/tests/test_sum_model.py b/tests/test_sum_model.py
index faba3be22..a1bd181f1 100644
--- a/tests/test_sum_model.py
+++ b/tests/test_sum_model.py
@@ -28,6 +28,7 @@ def test_init(self):
         with self.assertRaises(ValueError):
             m1.len_scale = 1
 
+    def test_attr_set(self):
         s1 = gs.Gaussian() + gs.Gaussian() + gs.Gaussian()
         with self.assertRaises(ValueError):
             s1.vars = [1, 2]
@@ -50,6 +51,7 @@ def test_init(self):
         with self.assertRaises(ValueError):
             s1.ratios = [0.3, 0.2]
 
+    def test_compare(self):
         s1 = gs.Gaussian(var=1) + gs.Exponential(var=2)
         s2 = gs.Exponential(var=1) + gs.Gaussian(var=2)
         self.assertFalse(s1 == gs.Nugget(dim=3))
@@ -57,6 +59,7 @@ def test_init(self):
         self.assertFalse(gs.Exponential() == (gs.Exponential() + gs.Nugget()))
         self.assertFalse((gs.Exponential() + gs.Nugget()) == gs.Exponential())
 
+    def test_copy(self):
         # check that models get copied
         m1 = gs.Gaussian()
         s1 = m1 + m1
@@ -64,19 +67,26 @@ def test_init(self):
         s1.vars = var
         np.testing.assert_array_almost_equal(s1.vars, var)
 
+    def test_var_dist(self):
         s1 = gs.SumModel(gs.Exponential, gs.Exponential, var=3)
         np.testing.assert_array_almost_equal(s1.vars, [1.5, 1.5])
+
+    def test_presence(self):
+        s1 = gs.SumModel(gs.Exponential, gs.Exponential)
         self.assertFalse(gs.Gaussian() in s1)
 
+    def test_len_dist(self):
         s1 = gs.SumModel(gs.Exponential, gs.Exponential, len_scale=10)
         np.testing.assert_array_almost_equal(s1.len_scales, [10, 10])
 
+    def test_temporal(self):
         s1 = gs.SumModel(
             gs.Exponential, gs.Exponential, temporal=True, spatial_dim=2
         )
         self.assertTrue(all(mod.temporal for mod in s1))
         self.assertTrue(all(mod.dim == 3 for mod in s1))
 
+    def test_magic(self):
         m1 = gs.Gaussian(dim=2, var=1.0, len_scale=1.0)
         m2 = gs.Matern(dim=2, var=2.0, len_scale=2.0, nu=2.0)
         m3 = gs.Integral(dim=2, var=3.0, len_scale=3.0, nu=3.0)
@@ -94,6 +104,7 @@ def test_init(self):
         self.assertTrue(s2 == s3)
         self.assertTrue(s1 == s3)
 
+    def test_exceptions(self):
         with self.assertRaises(ValueError):
             gs.Exponential() + 1
 

From e7f28c8bd1e0e772efb1464d16e9752fccce142c Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Sebastian=20M=C3=BCller?= <mueller.seb@posteo.de>
Date: Tue, 22 Apr 2025 21:23:33 +0200
Subject: [PATCH 56/58] sum-model: update doc-strings and error messages from
 review

---
 src/gstools/covmodel/base.py      | 18 +++++++++++++++++-
 src/gstools/covmodel/fit.py       |  4 ++++
 src/gstools/covmodel/models.py    |  3 +++
 src/gstools/covmodel/sum_tools.py | 18 ++++++++++++++----
 4 files changed, 38 insertions(+), 5 deletions(-)

diff --git a/src/gstools/covmodel/base.py b/src/gstools/covmodel/base.py
index c0aa1ea32..53a925f09 100644
--- a/src/gstools/covmodel/base.py
+++ b/src/gstools/covmodel/base.py
@@ -150,6 +150,7 @@ class CovModel:
     """
 
     _add_doc = True
+    """Flag to append the above doc-string about parameters to the model implementation."""
 
     def __init__(
         self,
@@ -499,7 +500,11 @@ def percentile_scale(self, per=0.9):
 
     @property
     def needs_fourier_transform(self):
-        """Whether the model doesn't implement an analytical spectral density."""
+        """
+        Flag whether the model needs a fourier transform to calculate
+        the spectral density from the covariance function or if
+        it implements an analytical spectral density.
+        """
         base_method = getattr(CovModel, "spectral_density")
         instance_method = getattr(self.__class__, "spectral_density")
         return base_method == instance_method
@@ -1430,6 +1435,11 @@ def set_var_weights(self, weights, skip=None, var=None):
         """
         Set variances of contained models by weights.
 
+        The variances of the sub-models act as ratios for the sum-model.
+        This function enables to keep the total variance of the sum-model
+        and reset the individual variances of the contained sub-models
+        by the given weights.
+
         Parameters
         ----------
         weights : iterable
@@ -1450,6 +1460,12 @@ def set_len_weights(self, weights, skip=None, len_scale=None):
         """
         Set length scales of contained models by weights.
 
+        The variances of the sub-models act as ratios for the sub-model length
+        scales to determine the total length scale of the sum-model.
+        This function enables to keep the total length scale of the sum-model as
+        well as the variances of the sub-models and reset the individual length scales
+        of the contained sub-models by the given weights.
+
         Parameters
         ----------
         weights : iterable
diff --git a/src/gstools/covmodel/fit.py b/src/gstools/covmodel/fit.py
index 0384ccbf5..8ad510de4 100755
--- a/src/gstools/covmodel/fit.py
+++ b/src/gstools/covmodel/fit.py
@@ -252,6 +252,10 @@ def _pre_para(model, para_select, sill, anis):
 
 
 def _check_sill(model, para_select, sill, bnd, sum_cfg):
+    """
+    This functions checks if the selected values for
+    variance, nugget and sill are valid.
+    """
     is_sum = hasattr(model, "sub_arg")
     sill_low = bnd["var"][0] + bnd["nugget"][0]
     sill_up = bnd["var"][1] + bnd["nugget"][1]
diff --git a/src/gstools/covmodel/models.py b/src/gstools/covmodel/models.py
index 7c1cd836d..cb495751f 100644
--- a/src/gstools/covmodel/models.py
+++ b/src/gstools/covmodel/models.py
@@ -53,6 +53,9 @@
 class Nugget(SumModel):
     r"""Pure nugget model.
 
+    This model has no correlated variability and represents pure noise.
+    The length scale of the model will be zero.
+
     Parameters
     ----------
     dim : :class:`int`, optional
diff --git a/src/gstools/covmodel/sum_tools.py b/src/gstools/covmodel/sum_tools.py
index 74604d381..899917919 100644
--- a/src/gstools/covmodel/sum_tools.py
+++ b/src/gstools/covmodel/sum_tools.py
@@ -132,13 +132,18 @@ def sum_set_var_weights(summod, weights, skip=None, var=None):
         raise ValueError(msg)
     ids = range(len(summod))
     if fail := set(skip) - set(ids):
-        msg = f"SumModel.set_var_weights: skip ids not valid: {fail}"
+        msg = (
+            f"SumModel.set_var_weights: ids given by 'skip' not valid: {fail}"
+        )
         raise ValueError(msg)
     var = summod.var if var is None else float(var)
     var_sum = sum(summod.models[i].var for i in skip)
     var_diff = var - var_sum
     if var_diff < 0:
-        msg = "SumModel.set_var_weights: skipped variances to big."
+        msg = (
+            "SumModel.set_var_weights: summed variances selected "
+            "with 'skip' already to big to keep total variance."
+        )
         raise ValueError(msg)
     weights_sum = sum(weights)
     var_list = summod.vars
@@ -175,7 +180,9 @@ def sum_set_len_weights(summod, weights, skip=None, len_scale=None):
         raise ValueError(msg)
     ids = range(len(summod))
     if fail := set(skip) - set(ids):
-        msg = f"SumModel.set_len_weights: skip ids not valid: {fail}"
+        msg = (
+            f"SumModel.set_len_weights: ids given by 'skip' not valid: {fail}"
+        )
         raise ValueError(msg)
     len_scale = summod.len_scale if len_scale is None else float(len_scale)
     # also skip models with no variance (not contributing to total len scale)
@@ -193,7 +200,10 @@ def sum_set_len_weights(summod, weights, skip=None, len_scale=None):
     len_sum = sum(summod[i].len_scale * summod.ratios[i] for i in skip)
     len_diff = len_scale - len_sum
     if len_diff < 0:
-        msg = "SumModel.set_len_weights: skipped length scales to big."
+        msg = (
+            "SumModel.set_len_weights: summed length scales "
+            "selected with 'skip' already to big to keep total length scale."
+        )
         raise ValueError(msg)
     weights_sum = sum(weights)
     len_scales = summod.len_scales

From b11af601ab72a6e2f70d820f9de8ef983f021553 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Sebastian=20M=C3=BCller?= <mueller.seb@posteo.de>
Date: Mon, 28 Apr 2025 07:54:12 +0200
Subject: [PATCH 57/58] fix spelling

---
 src/gstools/covmodel/sum_tools.py | 4 ++--
 1 file changed, 2 insertions(+), 2 deletions(-)

diff --git a/src/gstools/covmodel/sum_tools.py b/src/gstools/covmodel/sum_tools.py
index 899917919..253ee5a4b 100644
--- a/src/gstools/covmodel/sum_tools.py
+++ b/src/gstools/covmodel/sum_tools.py
@@ -142,7 +142,7 @@ def sum_set_var_weights(summod, weights, skip=None, var=None):
     if var_diff < 0:
         msg = (
             "SumModel.set_var_weights: summed variances selected "
-            "with 'skip' already to big to keep total variance."
+            "with 'skip' already too big to keep total variance."
         )
         raise ValueError(msg)
     weights_sum = sum(weights)
@@ -202,7 +202,7 @@ def sum_set_len_weights(summod, weights, skip=None, len_scale=None):
     if len_diff < 0:
         msg = (
             "SumModel.set_len_weights: summed length scales "
-            "selected with 'skip' already to big to keep total length scale."
+            "selected with 'skip' already too big to keep total length scale."
         )
         raise ValueError(msg)
     weights_sum = sum(weights)

From ee622a727b982dd5ea168934e0e5dbfe05aa42d5 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Sebastian=20M=C3=BCller?= <mueller.seb@posteo.de>
Date: Mon, 28 Apr 2025 07:54:30 +0200
Subject: [PATCH 58/58] move sum examples to number 12

---
 docs/source/conf.py                                         | 6 +++---
 docs/source/tutorials.rst                                   | 2 +-
 .../{11_sum_model => 12_sum_model}/00_simple_sum_model.py   | 0
 .../{11_sum_model => 12_sum_model}/01_fitting_sum_model.py  | 0
 examples/{11_sum_model => 12_sum_model}/README.rst          | 0
 5 files changed, 4 insertions(+), 4 deletions(-)
 rename examples/{11_sum_model => 12_sum_model}/00_simple_sum_model.py (100%)
 rename examples/{11_sum_model => 12_sum_model}/01_fitting_sum_model.py (100%)
 rename examples/{11_sum_model => 12_sum_model}/README.rst (100%)

diff --git a/docs/source/conf.py b/docs/source/conf.py
index 2b35fe59f..aaddc29f3 100644
--- a/docs/source/conf.py
+++ b/docs/source/conf.py
@@ -199,7 +199,7 @@ def setup(app):
 
 # -- Options for LaTeX output ---------------------------------------------
 # latex_engine = 'lualatex'
-# logo to big
+# logo too big
 latex_logo = "pics/gstools_150.png"
 
 # latex_show_urls = 'footnote'
@@ -300,8 +300,8 @@ def setup(app):
         "../../examples/08_geo_coordinates/",
         "../../examples/09_spatio_temporal/",
         "../../examples/10_normalizer/",
-        "../../examples/11_sum_model/",
         "../../examples/11_plurigaussian/",
+        "../../examples/12_sum_model/",
     ],
     # path where to save gallery generated examples
     "gallery_dirs": [
@@ -316,8 +316,8 @@ def setup(app):
         "examples/08_geo_coordinates/",
         "examples/09_spatio_temporal/",
         "examples/10_normalizer/",
-        "examples/11_sum_model/",
         "examples/11_plurigaussian/",
+        "examples/12_sum_model/",
     ],
     # Pattern to search for example files
     "filename_pattern": r"\.py",
diff --git a/docs/source/tutorials.rst b/docs/source/tutorials.rst
index 5e9c27c47..e3b20db1b 100644
--- a/docs/source/tutorials.rst
+++ b/docs/source/tutorials.rst
@@ -22,8 +22,8 @@ explore its whole beauty and power.
    examples/08_geo_coordinates/index
    examples/09_spatio_temporal/index
    examples/10_normalizer/index
-   examples/11_sum_model/index
    examples/11_plurigaussian/index
+   examples/12_sum_model/index
    examples/00_misc/index
 
 .. only:: html
diff --git a/examples/11_sum_model/00_simple_sum_model.py b/examples/12_sum_model/00_simple_sum_model.py
similarity index 100%
rename from examples/11_sum_model/00_simple_sum_model.py
rename to examples/12_sum_model/00_simple_sum_model.py
diff --git a/examples/11_sum_model/01_fitting_sum_model.py b/examples/12_sum_model/01_fitting_sum_model.py
similarity index 100%
rename from examples/11_sum_model/01_fitting_sum_model.py
rename to examples/12_sum_model/01_fitting_sum_model.py
diff --git a/examples/11_sum_model/README.rst b/examples/12_sum_model/README.rst
similarity index 100%
rename from examples/11_sum_model/README.rst
rename to examples/12_sum_model/README.rst