From 12075a1496b3d2eece036760f991c35e922a60c8 Mon Sep 17 00:00:00 2001
From: Bas Nijholt <bas@nijho.lt>
Date: Fri, 19 Mar 2021 10:13:59 +0000
Subject: [PATCH 1/7] evaluate coeffs lazily

---
 adaptive/learner/integrator_coeffs.py  | 138 ++++++++++++++++++-------
 adaptive/learner/integrator_learner.py |  54 ++++------
 adaptive/tests/test_cquad.py           |   4 +-
 3 files changed, 126 insertions(+), 70 deletions(-)

diff --git a/adaptive/learner/integrator_coeffs.py b/adaptive/learner/integrator_coeffs.py
index 7719f601f..09f211f34 100644
--- a/adaptive/learner/integrator_coeffs.py
+++ b/adaptive/learner/integrator_coeffs.py
@@ -141,39 +141,105 @@ def calc_V(x, n):
     return np.array(V).T
 
 
-eps = np.spacing(1)
-
-# the nodes and Newton polynomials
-ns = (5, 9, 17, 33)
-xi = [-np.cos(np.linspace(0, np.pi, n)) for n in ns]
-
-# Make `xi` perfectly anti-symmetric, important for splitting the intervals
-xi = [(row - row[::-1]) / 2 for row in xi]
-
-# Compute the Vandermonde-like matrix and its inverse.
-V = [calc_V(x, n) for x, n in zip(xi, ns)]
-V_inv = list(map(scipy.linalg.inv, V))
-Vcond = [scipy.linalg.norm(a, 2) * scipy.linalg.norm(b, 2) for a, b in zip(V, V_inv)]
-
-# Compute the shift matrices.
-T_left, T_right = [V_inv[3] @ calc_V((xi[3] + a) / 2, ns[3]) for a in [-1, 1]]
-
-# If the relative difference between two consecutive approximations is
-# lower than this value, the error estimate is considered reliable.
-# See section 6.2 of Pedro Gonnet's thesis.
-hint = 0.1
-
-# Smallest acceptable relative difference of points in a rule.  This was chosen
-# such that no artifacts are apparent in plots of (i, log(a_i)), where a_i is
-# the sequence of estimates of the integral value of an interval and all its
-# ancestors..
-min_sep = 16 * eps
-
-ndiv_max = 20
-
-# set-up the downdate matrix
-k = np.arange(ns[3])
-alpha = np.sqrt((k + 1) ** 2 / (2 * k + 1) / (2 * k + 3))
-gamma = np.concatenate([[0, 0], np.sqrt(k[2:] ** 2 / (4 * k[2:] ** 2 - 1))])
-
-b_def = calc_bdef(ns)
+class Coefficients:
+    def __init__(self) -> None:
+        self.is_set = False
+        # The nodes
+        self.ns = (5, 9, 17, 33)
+        self.eps = np.spacing(1)
+
+        # If the relative difference between two consecutive approximations is
+        # lower than this value, the error estimate is considered reliable.
+        # See section 6.2 of Pedro Gonnet's thesis.
+        self.hint = 0.1
+
+        # Smallest acceptable relative difference of points in a rule.  This was chosen
+        # such that no artifacts are apparent in plots of (i, log(a_i)), where a_i is
+        # the sequence of estimates of the integral value of an interval and all its
+        # ancestors..
+        self.min_sep = 16 * self.eps
+
+        # Maximum amount of subdivisions
+        self.ndiv_max = 20
+
+    def set(self):
+        self.is_set = True
+        # The Newton polynomials
+        xi = [-np.cos(np.linspace(0, np.pi, n)) for n in self.ns]
+        # Make `xi` perfectly anti-symmetric, important for splitting the intervals
+        self._xi = [(row - row[::-1]) / 2 for row in xi]
+
+        # Compute the Vandermonde-like matrix and its inverse.
+        self._V = [calc_V(x, n) for x, n in zip(xi, self.ns)]
+        self._V_inv = list(map(scipy.linalg.inv, self._V))
+        self._Vcond = [
+            scipy.linalg.norm(a, 2) * scipy.linalg.norm(b, 2)
+            for a, b in zip(self._V, self._V_inv)
+        ]
+
+        # Compute the shift matrices.
+        self._T_left, self._T_right = [
+            self._V_inv[3] @ calc_V((xi[3] + a) / 2, self.ns[3]) for a in [-1, 1]
+        ]
+
+        # set-up the downdate matrix
+        k = np.arange(self.ns[3])
+        self._alpha = np.sqrt((k + 1) ** 2 / (2 * k + 1) / (2 * k + 3))
+        self._gamma = np.concatenate(
+            [[0, 0], np.sqrt(k[2:] ** 2 / (4 * k[2:] ** 2 - 1))]
+        )
+        self._b_def = calc_bdef(self.ns)
+
+    @property
+    def xi(self):
+        if not self.is_set:
+            self.set()
+        return self._xi
+
+    @property
+    def V(self):
+        if not self.is_set:
+            self.set()
+        return self._V
+
+    @property
+    def V_inv(self):
+        if not self.is_set:
+            self.set()
+        return self._V_inv
+
+    @property
+    def Vcond(self):
+        if not self.is_set:
+            self.set()
+        return self._Vcond
+
+    @property
+    def T_right(self):
+        if not self.is_set:
+            self.set()
+        return self._T_right
+
+    @property
+    def T_left(self):
+        if not self.is_set:
+            self.set()
+        return self._T_left
+
+    @property
+    def alpha(self):
+        if not self.is_set:
+            self.set()
+        return self._alpha
+
+    @property
+    def gamma(self):
+        if not self.is_set:
+            self.set()
+        return self._gamma
+
+    @property
+    def b_def(self):
+        if not self.is_set:
+            self.set()
+        return self._b_def
diff --git a/adaptive/learner/integrator_learner.py b/adaptive/learner/integrator_learner.py
index 52a930de2..46e5cc433 100644
--- a/adaptive/learner/integrator_learner.py
+++ b/adaptive/learner/integrator_learner.py
@@ -14,33 +14,21 @@
 from adaptive.notebook_integration import ensure_holoviews
 from adaptive.utils import cache_latest, restore
 
-from .integrator_coeffs import (
-    T_left,
-    T_right,
-    V_inv,
-    Vcond,
-    alpha,
-    b_def,
-    eps,
-    gamma,
-    hint,
-    min_sep,
-    ndiv_max,
-    ns,
-    xi,
-)
+from .integrator_coeffs import Coefficients
+
+c = Coefficients()
 
 
 def _downdate(c, nans, depth):
     # This is algorithm 5 from the thesis of Pedro Gonnet.
-    b = b_def[depth].copy()
-    m = ns[depth] - 1
+    b = c.b_def[depth].copy()
+    m = c.ns[depth] - 1
     for i in nans:
-        b[m + 1] /= alpha[m]
-        xii = xi[depth][i]
-        b[m] = (b[m] + xii * b[m + 1]) / alpha[m - 1]
+        b[m + 1] /= c.alpha[m]
+        xii = c.xi[depth][i]
+        b[m] = (b[m] + xii * b[m + 1]) / c.alpha[m - 1]
         for j in range(m - 1, 0, -1):
-            b[j] = (b[j] + xii * b[j + 1] - gamma[j + 1] * b[j + 2]) / alpha[j - 1]
+            b[j] = (b[j] + xii * b[j + 1] - c.gamma[j + 1] * b[j + 2]) / c.alpha[j - 1]
         b = b[1:]
 
         c[:m] -= c[m] / b[m] * b[:m]
@@ -62,7 +50,7 @@ def _zero_nans(fx):
 def _calc_coeffs(fx, depth):
     """Caution: this function modifies fx."""
     nans = _zero_nans(fx)
-    c_new = V_inv[depth] @ fx
+    c_new = c.V_inv[depth] @ fx
     if nans:
         fx[nans] = np.nan
         c_new = _downdate(c_new, nans, depth)
@@ -168,11 +156,11 @@ def T(self):
         left = self.a == self.parent.a
         right = self.b == self.parent.b
         assert left != right
-        return T_left if left else T_right
+        return c.T_left if left else c.T_right
 
     def refinement_complete(self, depth):
         """The interval has all the y-values to calculate the intergral."""
-        if len(self.data) < ns[depth]:
+        if len(self.data) < c.ns[depth]:
             return False
         return all(p in self.data for p in self.points(depth))
 
@@ -181,7 +169,7 @@ def points(self, depth=None):
             depth = self.depth
         a = self.a
         b = self.b
-        return (a + b) / 2 + (b - a) * xi[depth] / 2
+        return (a + b) / 2 + (b - a) * c.xi[depth] / 2
 
     def refine(self):
         self.depth += 1
@@ -234,7 +222,7 @@ def calc_ndiv(self):
         div = self.parent.c00 and self.c00 / self.parent.c00 > 2
         self.ndiv += div
 
-        if self.ndiv > ndiv_max and 2 * self.ndiv > self.rdepth:
+        if self.ndiv > c.ndiv_max and 2 * self.ndiv > self.rdepth:
             raise DivergentIntegralError
 
         if div:
@@ -243,7 +231,7 @@ def calc_ndiv(self):
 
     def update_ndiv_recursively(self):
         self.ndiv += 1
-        if self.ndiv > ndiv_max and 2 * self.ndiv > self.rdepth:
+        if self.ndiv > c.ndiv_max and 2 * self.ndiv > self.rdepth:
             raise DivergentIntegralError
 
         for child in self.children:
@@ -276,13 +264,13 @@ def complete_process(self, depth):
         if depth:
             # Refine
             c_diff = self.calc_err(c_old)
-            force_split = c_diff > hint * norm(self.c)
+            force_split = c_diff > c.hint * norm(self.c)
         else:
             # Split
             self.c00 = self.c[0]
 
             if self.parent.depth_complete is not None:
-                c_old = self.T[:, : ns[self.parent.depth_complete]] @ self.parent.c
+                c_old = self.T[:, : c.ns[self.parent.depth_complete]] @ self.parent.c
                 self.calc_err(c_old)
                 self.calc_ndiv()
 
@@ -290,7 +278,7 @@ def complete_process(self, depth):
                 if child.depth_complete is not None:
                     child.calc_ndiv()
                 if child.depth_complete == 0:
-                    c_old = child.T[:, : ns[self.depth_complete]] @ self.c
+                    c_old = child.T[:, : c.ns[self.depth_complete]] @ self.c
                     child.calc_err(c_old)
 
         if self.done_leaves is not None and not len(self.done_leaves):
@@ -320,7 +308,7 @@ def complete_process(self, depth):
                 ival.done_leaves -= old_leaves
                 ival = ival.parent
 
-        remove = self.err < (abs(self.igral) * eps * Vcond[depth])
+        remove = self.err < (abs(self.igral) * c.eps * c.Vcond[depth])
 
         return force_split, remove
 
@@ -496,8 +484,8 @@ def _fill_stack(self):
         points = ival.points()
 
         if (
-            points[1] - points[0] < points[0] * min_sep
-            or points[-1] - points[-2] < points[-2] * min_sep
+            points[1] - points[0] < points[0] * c.min_sep
+            or points[-1] - points[-2] < points[-2] * c.min_sep
         ):
             self.ivals.remove(ival)
         elif ival.depth == 3 or force_split:
diff --git a/adaptive/tests/test_cquad.py b/adaptive/tests/test_cquad.py
index b1b0b3be1..bd89813da 100644
--- a/adaptive/tests/test_cquad.py
+++ b/adaptive/tests/test_cquad.py
@@ -5,7 +5,7 @@
 import pytest
 
 from adaptive.learner import IntegratorLearner
-from adaptive.learner.integrator_coeffs import ns
+from adaptive.learner.integrator_coeffs import Coefficients
 from adaptive.learner.integrator_learner import DivergentIntegralError
 
 from .algorithm_4 import DivergentIntegralError as A4DivergentIntegralError
@@ -13,6 +13,8 @@
 
 eps = np.spacing(1)
 
+ns = Coefficients().ns
+
 
 def run_integrator_learner(f, a, b, tol, n):
     learner = IntegratorLearner(f, bounds=(a, b), tol=tol)

From 5be5a66a8a8aac33933d2e4160c257c9d8c63ffa Mon Sep 17 00:00:00 2001
From: Bas Nijholt <bas@nijho.lt>
Date: Fri, 19 Mar 2021 10:34:00 +0000
Subject: [PATCH 2/7] use cached_property

---
 adaptive/learner/integrator_coeffs.py | 97 ++++++++++-----------------
 1 file changed, 35 insertions(+), 62 deletions(-)

diff --git a/adaptive/learner/integrator_coeffs.py b/adaptive/learner/integrator_coeffs.py
index 09f211f34..b97b23923 100644
--- a/adaptive/learner/integrator_coeffs.py
+++ b/adaptive/learner/integrator_coeffs.py
@@ -1,5 +1,6 @@
 # Based on an adaptive quadrature algorithm by Pedro Gonnet
 
+import functools
 from collections import defaultdict
 from fractions import Fraction
 
@@ -141,6 +142,10 @@ def calc_V(x, n):
     return np.array(V).T
 
 
+def cached_property(f):
+    return property(functools.lru_cache(None)(f))
+
+
 class Coefficients:
     def __init__(self) -> None:
         self.is_set = False
@@ -162,84 +167,52 @@ def __init__(self) -> None:
         # Maximum amount of subdivisions
         self.ndiv_max = 20
 
-    def set(self):
-        self.is_set = True
+    @cached_property
+    def xi(self):
         # The Newton polynomials
         xi = [-np.cos(np.linspace(0, np.pi, n)) for n in self.ns]
         # Make `xi` perfectly anti-symmetric, important for splitting the intervals
-        self._xi = [(row - row[::-1]) / 2 for row in xi]
-
-        # Compute the Vandermonde-like matrix and its inverse.
-        self._V = [calc_V(x, n) for x, n in zip(xi, self.ns)]
-        self._V_inv = list(map(scipy.linalg.inv, self._V))
-        self._Vcond = [
-            scipy.linalg.norm(a, 2) * scipy.linalg.norm(b, 2)
-            for a, b in zip(self._V, self._V_inv)
-        ]
-
-        # Compute the shift matrices.
-        self._T_left, self._T_right = [
-            self._V_inv[3] @ calc_V((xi[3] + a) / 2, self.ns[3]) for a in [-1, 1]
-        ]
+        return [(row - row[::-1]) / 2 for row in xi]
 
-        # set-up the downdate matrix
-        k = np.arange(self.ns[3])
-        self._alpha = np.sqrt((k + 1) ** 2 / (2 * k + 1) / (2 * k + 3))
-        self._gamma = np.concatenate(
-            [[0, 0], np.sqrt(k[2:] ** 2 / (4 * k[2:] ** 2 - 1))]
-        )
-        self._b_def = calc_bdef(self.ns)
-
-    @property
-    def xi(self):
-        if not self.is_set:
-            self.set()
-        return self._xi
-
-    @property
+    @cached_property
     def V(self):
-        if not self.is_set:
-            self.set()
-        return self._V
+        # Compute the Vandermonde-like matrix and its inverse.
+        return [calc_V(x, n) for x, n in zip(self.xi, self.ns)]
 
-    @property
+    @cached_property
     def V_inv(self):
-        if not self.is_set:
-            self.set()
-        return self._V_inv
+        # Compute the inverse Vandermonde-like matrix
+        return list(map(scipy.linalg.inv, self.V))
 
-    @property
+    @cached_property
     def Vcond(self):
-        if not self.is_set:
-            self.set()
-        return self._Vcond
+        return [
+            scipy.linalg.norm(a, 2) * scipy.linalg.norm(b, 2)
+            for a, b in zip(self.V, self.V_inv)
+        ]
 
-    @property
+    @cached_property
+    def k(self):
+        return np.arange(self.ns[3])
+
+    @cached_property
     def T_right(self):
-        if not self.is_set:
-            self.set()
-        return self._T_right
+        return self.V_inv[3] @ calc_V((self.xi[3] + 1) / 2, self.ns[3])
 
-    @property
+    @cached_property
     def T_left(self):
-        if not self.is_set:
-            self.set()
-        return self._T_left
+        return self.V_inv[3] @ calc_V((self.xi[3] - 1) / 2, self.ns[3])
 
-    @property
+    @cached_property
     def alpha(self):
-        if not self.is_set:
-            self.set()
-        return self._alpha
+        return np.sqrt((self.k + 1) ** 2 / (2 * self.k + 1) / (2 * self.k + 3))
 
-    @property
+    @cached_property
     def gamma(self):
-        if not self.is_set:
-            self.set()
-        return self._gamma
+        return np.concatenate(
+            [[0, 0], np.sqrt(self.k[2:] ** 2 / (4 * self.k[2:] ** 2 - 1))]
+        )
 
-    @property
+    @cached_property
     def b_def(self):
-        if not self.is_set:
-            self.set()
-        return self._b_def
+        return calc_bdef(self.ns)

From d748cc70b1f35ec4151477da807432ba7d51719e Mon Sep 17 00:00:00 2001
From: Bas Nijholt <bas@nijho.lt>
Date: Fri, 19 Mar 2021 10:36:40 +0000
Subject: [PATCH 3/7] use cached_property is available

---
 adaptive/learner/integrator_coeffs.py | 10 +++++-----
 adaptive/utils.py                     |  4 ++++
 2 files changed, 9 insertions(+), 5 deletions(-)

diff --git a/adaptive/learner/integrator_coeffs.py b/adaptive/learner/integrator_coeffs.py
index b97b23923..b8a2126d8 100644
--- a/adaptive/learner/integrator_coeffs.py
+++ b/adaptive/learner/integrator_coeffs.py
@@ -1,12 +1,16 @@
 # Based on an adaptive quadrature algorithm by Pedro Gonnet
 
-import functools
 from collections import defaultdict
 from fractions import Fraction
 
 import numpy as np
 import scipy.linalg
 
+try:
+    from functools import cached_property
+except ImportError:
+    from adaptive.utils import cached_property
+
 
 def legendre(n):
     """Return the first n Legendre polynomials.
@@ -142,10 +146,6 @@ def calc_V(x, n):
     return np.array(V).T
 
 
-def cached_property(f):
-    return property(functools.lru_cache(None)(f))
-
-
 class Coefficients:
     def __init__(self) -> None:
         self.is_set = False
diff --git a/adaptive/utils.py b/adaptive/utils.py
index 035086205..97ed20f79 100644
--- a/adaptive/utils.py
+++ b/adaptive/utils.py
@@ -83,3 +83,7 @@ def __call__(self, *args, **kwargs):
                     msg = f"The attribute '{name}' should be of type {type_}, not {type(x)}."
                     raise TypeError(msg)
         return obj
+
+
+def cached_property(f):
+    return property(functools.lru_cache(None)(f))

From 4cf4a54227aa5c647f442b8457b3b0c7d4826df1 Mon Sep 17 00:00:00 2001
From: Bas Nijholt <bas@nijho.lt>
Date: Fri, 19 Mar 2021 10:41:31 +0000
Subject: [PATCH 4/7] delay setting of coefficients

---
 adaptive/learner/integrator_coeffs.py | 3 +--
 1 file changed, 1 insertion(+), 2 deletions(-)

diff --git a/adaptive/learner/integrator_coeffs.py b/adaptive/learner/integrator_coeffs.py
index b8a2126d8..2af5940b7 100644
--- a/adaptive/learner/integrator_coeffs.py
+++ b/adaptive/learner/integrator_coeffs.py
@@ -148,10 +148,8 @@ def calc_V(x, n):
 
 class Coefficients:
     def __init__(self) -> None:
-        self.is_set = False
         # The nodes
         self.ns = (5, 9, 17, 33)
-        self.eps = np.spacing(1)
 
         # If the relative difference between two consecutive approximations is
         # lower than this value, the error estimate is considered reliable.
@@ -162,6 +160,7 @@ def __init__(self) -> None:
         # such that no artifacts are apparent in plots of (i, log(a_i)), where a_i is
         # the sequence of estimates of the integral value of an interval and all its
         # ancestors..
+        self.eps = np.spacing(1)
         self.min_sep = 16 * self.eps
 
         # Maximum amount of subdivisions

From 878068ef16cdd877927b8ef309c9022b014b535b Mon Sep 17 00:00:00 2001
From: Bas Nijholt <bas@nijho.lt>
Date: Fri, 19 Mar 2021 10:50:04 +0000
Subject: [PATCH 5/7] rename c to coeffs

---
 adaptive/learner/integrator_learner.py | 42 ++++++++++++++------------
 1 file changed, 23 insertions(+), 19 deletions(-)

diff --git a/adaptive/learner/integrator_learner.py b/adaptive/learner/integrator_learner.py
index 46e5cc433..1b244f37e 100644
--- a/adaptive/learner/integrator_learner.py
+++ b/adaptive/learner/integrator_learner.py
@@ -16,19 +16,21 @@
 
 from .integrator_coeffs import Coefficients
 
-c = Coefficients()
+coeff = Coefficients()
 
 
 def _downdate(c, nans, depth):
     # This is algorithm 5 from the thesis of Pedro Gonnet.
-    b = c.b_def[depth].copy()
-    m = c.ns[depth] - 1
+    b = coeff.b_def[depth].copy()
+    m = coeff.ns[depth] - 1
     for i in nans:
-        b[m + 1] /= c.alpha[m]
-        xii = c.xi[depth][i]
-        b[m] = (b[m] + xii * b[m + 1]) / c.alpha[m - 1]
+        b[m + 1] /= coeff.alpha[m]
+        xii = coeff.xi[depth][i]
+        b[m] = (b[m] + xii * b[m + 1]) / coeff.alpha[m - 1]
         for j in range(m - 1, 0, -1):
-            b[j] = (b[j] + xii * b[j + 1] - c.gamma[j + 1] * b[j + 2]) / c.alpha[j - 1]
+            b[j] = (
+                b[j] + xii * b[j + 1] - coeff.gamma[j + 1] * b[j + 2]
+            ) / coeff.alpha[j - 1]
         b = b[1:]
 
         c[:m] -= c[m] / b[m] * b[:m]
@@ -50,7 +52,7 @@ def _zero_nans(fx):
 def _calc_coeffs(fx, depth):
     """Caution: this function modifies fx."""
     nans = _zero_nans(fx)
-    c_new = c.V_inv[depth] @ fx
+    c_new = coeff.V_inv[depth] @ fx
     if nans:
         fx[nans] = np.nan
         c_new = _downdate(c_new, nans, depth)
@@ -156,11 +158,11 @@ def T(self):
         left = self.a == self.parent.a
         right = self.b == self.parent.b
         assert left != right
-        return c.T_left if left else c.T_right
+        return coeff.T_left if left else coeff.T_right
 
     def refinement_complete(self, depth):
         """The interval has all the y-values to calculate the intergral."""
-        if len(self.data) < c.ns[depth]:
+        if len(self.data) < coeff.ns[depth]:
             return False
         return all(p in self.data for p in self.points(depth))
 
@@ -169,7 +171,7 @@ def points(self, depth=None):
             depth = self.depth
         a = self.a
         b = self.b
-        return (a + b) / 2 + (b - a) * c.xi[depth] / 2
+        return (a + b) / 2 + (b - a) * coeff.xi[depth] / 2
 
     def refine(self):
         self.depth += 1
@@ -222,7 +224,7 @@ def calc_ndiv(self):
         div = self.parent.c00 and self.c00 / self.parent.c00 > 2
         self.ndiv += div
 
-        if self.ndiv > c.ndiv_max and 2 * self.ndiv > self.rdepth:
+        if self.ndiv > coeff.ndiv_max and 2 * self.ndiv > self.rdepth:
             raise DivergentIntegralError
 
         if div:
@@ -231,7 +233,7 @@ def calc_ndiv(self):
 
     def update_ndiv_recursively(self):
         self.ndiv += 1
-        if self.ndiv > c.ndiv_max and 2 * self.ndiv > self.rdepth:
+        if self.ndiv > coeff.ndiv_max and 2 * self.ndiv > self.rdepth:
             raise DivergentIntegralError
 
         for child in self.children:
@@ -264,13 +266,15 @@ def complete_process(self, depth):
         if depth:
             # Refine
             c_diff = self.calc_err(c_old)
-            force_split = c_diff > c.hint * norm(self.c)
+            force_split = c_diff > coeff.hint * norm(self.c)
         else:
             # Split
             self.c00 = self.c[0]
 
             if self.parent.depth_complete is not None:
-                c_old = self.T[:, : c.ns[self.parent.depth_complete]] @ self.parent.c
+                c_old = (
+                    self.T[:, : coeff.ns[self.parent.depth_complete]] @ self.parent.c
+                )
                 self.calc_err(c_old)
                 self.calc_ndiv()
 
@@ -278,7 +282,7 @@ def complete_process(self, depth):
                 if child.depth_complete is not None:
                     child.calc_ndiv()
                 if child.depth_complete == 0:
-                    c_old = child.T[:, : c.ns[self.depth_complete]] @ self.c
+                    c_old = child.T[:, : coeff.ns[self.depth_complete]] @ self.c
                     child.calc_err(c_old)
 
         if self.done_leaves is not None and not len(self.done_leaves):
@@ -308,7 +312,7 @@ def complete_process(self, depth):
                 ival.done_leaves -= old_leaves
                 ival = ival.parent
 
-        remove = self.err < (abs(self.igral) * c.eps * c.Vcond[depth])
+        remove = self.err < (abs(self.igral) * coeff.eps * coeff.Vcond[depth])
 
         return force_split, remove
 
@@ -484,8 +488,8 @@ def _fill_stack(self):
         points = ival.points()
 
         if (
-            points[1] - points[0] < points[0] * c.min_sep
-            or points[-1] - points[-2] < points[-2] * c.min_sep
+            points[1] - points[0] < points[0] * coeff.min_sep
+            or points[-1] - points[-2] < points[-2] * coeff.min_sep
         ):
             self.ivals.remove(ival)
         elif ival.depth == 3 or force_split:

From 839c232be734b0b82867504a434e8b472763fb66 Mon Sep 17 00:00:00 2001
From: Bas Nijholt <bas@nijho.lt>
Date: Tue, 23 Mar 2021 10:14:55 +0100
Subject: [PATCH 6/7] implement integrator_coeffs using module.__getattr__

---
 adaptive/learner/integrator_coeffs.py  | 120 ++++++++++---------------
 adaptive/learner/integrator_learner.py |   5 +-
 adaptive/tests/test_cquad.py           |   8 +-
 adaptive/utils.py                      |   4 -
 4 files changed, 50 insertions(+), 87 deletions(-)

diff --git a/adaptive/learner/integrator_coeffs.py b/adaptive/learner/integrator_coeffs.py
index 2af5940b7..d590c8a29 100644
--- a/adaptive/learner/integrator_coeffs.py
+++ b/adaptive/learner/integrator_coeffs.py
@@ -2,15 +2,11 @@
 
 from collections import defaultdict
 from fractions import Fraction
+from functools import lru_cache
 
 import numpy as np
 import scipy.linalg
 
-try:
-    from functools import cached_property
-except ImportError:
-    from adaptive.utils import cached_property
-
 
 def legendre(n):
     """Return the first n Legendre polynomials.
@@ -146,72 +142,48 @@ def calc_V(x, n):
     return np.array(V).T
 
 
-class Coefficients:
-    def __init__(self) -> None:
-        # The nodes
-        self.ns = (5, 9, 17, 33)
-
-        # If the relative difference between two consecutive approximations is
-        # lower than this value, the error estimate is considered reliable.
-        # See section 6.2 of Pedro Gonnet's thesis.
-        self.hint = 0.1
-
-        # Smallest acceptable relative difference of points in a rule.  This was chosen
-        # such that no artifacts are apparent in plots of (i, log(a_i)), where a_i is
-        # the sequence of estimates of the integral value of an interval and all its
-        # ancestors..
-        self.eps = np.spacing(1)
-        self.min_sep = 16 * self.eps
-
-        # Maximum amount of subdivisions
-        self.ndiv_max = 20
-
-    @cached_property
-    def xi(self):
-        # The Newton polynomials
-        xi = [-np.cos(np.linspace(0, np.pi, n)) for n in self.ns]
-        # Make `xi` perfectly anti-symmetric, important for splitting the intervals
-        return [(row - row[::-1]) / 2 for row in xi]
-
-    @cached_property
-    def V(self):
-        # Compute the Vandermonde-like matrix and its inverse.
-        return [calc_V(x, n) for x, n in zip(self.xi, self.ns)]
-
-    @cached_property
-    def V_inv(self):
-        # Compute the inverse Vandermonde-like matrix
-        return list(map(scipy.linalg.inv, self.V))
-
-    @cached_property
-    def Vcond(self):
-        return [
-            scipy.linalg.norm(a, 2) * scipy.linalg.norm(b, 2)
-            for a, b in zip(self.V, self.V_inv)
-        ]
-
-    @cached_property
-    def k(self):
-        return np.arange(self.ns[3])
-
-    @cached_property
-    def T_right(self):
-        return self.V_inv[3] @ calc_V((self.xi[3] + 1) / 2, self.ns[3])
-
-    @cached_property
-    def T_left(self):
-        return self.V_inv[3] @ calc_V((self.xi[3] - 1) / 2, self.ns[3])
-
-    @cached_property
-    def alpha(self):
-        return np.sqrt((self.k + 1) ** 2 / (2 * self.k + 1) / (2 * self.k + 3))
-
-    @cached_property
-    def gamma(self):
-        return np.concatenate(
-            [[0, 0], np.sqrt(self.k[2:] ** 2 / (4 * self.k[2:] ** 2 - 1))]
-        )
-
-    @cached_property
-    def b_def(self):
-        return calc_bdef(self.ns)
+@lru_cache(maxsize=None)
+def _coefficients():
+    eps = np.spacing(1)
+
+    # the nodes and Newton polynomials
+    ns = (5, 9, 17, 33)
+    xi = [-np.cos(np.linspace(0, np.pi, n)) for n in ns]
+
+    # Make `xi` perfectly anti-symmetric, important for splitting the intervals
+    xi = [(row - row[::-1]) / 2 for row in xi]
+
+    # Compute the Vandermonde-like matrix and its inverse.
+    V = [calc_V(x, n) for x, n in zip(xi, ns)]
+    V_inv = list(map(scipy.linalg.inv, V))
+    Vcond = [
+        scipy.linalg.norm(a, 2) * scipy.linalg.norm(b, 2) for a, b in zip(V, V_inv)
+    ]
+
+    # Compute the shift matrices.
+    T_left, T_right = [V_inv[3] @ calc_V((xi[3] + a) / 2, ns[3]) for a in [-1, 1]]
+
+    # If the relative difference between two consecutive approximations is
+    # lower than this value, the error estimate is considered reliable.
+    # See section 6.2 of Pedro Gonnet's thesis.
+    hint = 0.1
+
+    # Smallest acceptable relative difference of points in a rule.  This was chosen
+    # such that no artifacts are apparent in plots of (i, log(a_i)), where a_i is
+    # the sequence of estimates of the integral value of an interval and all its
+    # ancestors..
+    min_sep = 16 * eps
+
+    ndiv_max = 20
+
+    # set-up the downdate matrix
+    k = np.arange(ns[3])
+    alpha = np.sqrt((k + 1) ** 2 / (2 * k + 1) / (2 * k + 3))
+    gamma = np.concatenate([[0, 0], np.sqrt(k[2:] ** 2 / (4 * k[2:] ** 2 - 1))])
+
+    b_def = calc_bdef(ns)
+    return locals()
+
+
+def __getattr__(attr):
+    return _coefficients()[attr]
diff --git a/adaptive/learner/integrator_learner.py b/adaptive/learner/integrator_learner.py
index 1b244f37e..15e3fcf2f 100644
--- a/adaptive/learner/integrator_learner.py
+++ b/adaptive/learner/integrator_learner.py
@@ -10,14 +10,11 @@
 from scipy.linalg import norm
 from sortedcontainers import SortedSet
 
+import adaptive.learner.integrator_coeffs as coeff
 from adaptive.learner.base_learner import BaseLearner
 from adaptive.notebook_integration import ensure_holoviews
 from adaptive.utils import cache_latest, restore
 
-from .integrator_coeffs import Coefficients
-
-coeff = Coefficients()
-
 
 def _downdate(c, nans, depth):
     # This is algorithm 5 from the thesis of Pedro Gonnet.
diff --git a/adaptive/tests/test_cquad.py b/adaptive/tests/test_cquad.py
index bd89813da..13e6552e3 100644
--- a/adaptive/tests/test_cquad.py
+++ b/adaptive/tests/test_cquad.py
@@ -4,8 +4,8 @@
 import numpy as np
 import pytest
 
+import adaptive.learner.integrator_coeffs as coeff
 from adaptive.learner import IntegratorLearner
-from adaptive.learner.integrator_coeffs import Coefficients
 from adaptive.learner.integrator_learner import DivergentIntegralError
 
 from .algorithm_4 import DivergentIntegralError as A4DivergentIntegralError
@@ -13,8 +13,6 @@
 
 eps = np.spacing(1)
 
-ns = Coefficients().ns
-
 
 def run_integrator_learner(f, a, b, tol, n):
     learner = IntegratorLearner(f, bounds=(a, b), tol=tol)
@@ -247,7 +245,7 @@ def test_removed_choose_mutiple_points_at_once():
     we should use the high-precision interval.
     """
     learner = IntegratorLearner(np.exp, bounds=(0, 1), tol=1e-15)
-    n = ns[-1] + 2 * (ns[0] - 2)  # first + two children (33+6=39)
+    n = coeff.ns[-1] + 2 * (coeff.ns[0] - 2)  # first + two children (33+6=39)
     xs, _ = learner.ask(n)
     for x in xs:
         learner.tell(x, learner.function(x))
@@ -259,7 +257,7 @@ def test_removed_ask_one_by_one():
         # This test should raise because integrating np.exp should be done
         # after the 33th point
         learner = IntegratorLearner(np.exp, bounds=(0, 1), tol=1e-15)
-        n = ns[-1] + 2 * (ns[0] - 2)  # first + two children (33+6=39)
+        n = coeff.ns[-1] + 2 * (coeff.ns[0] - 2)  # first + two children (33+6=39)
         for _ in range(n):
             xs, _ = learner.ask(1)
             for x in xs:
diff --git a/adaptive/utils.py b/adaptive/utils.py
index 97ed20f79..035086205 100644
--- a/adaptive/utils.py
+++ b/adaptive/utils.py
@@ -83,7 +83,3 @@ def __call__(self, *args, **kwargs):
                     msg = f"The attribute '{name}' should be of type {type_}, not {type(x)}."
                     raise TypeError(msg)
         return obj
-
-
-def cached_property(f):
-    return property(functools.lru_cache(None)(f))

From ebae237f8a97368970eb9a82cee82a01fc6defa6 Mon Sep 17 00:00:00 2001
From: Bas Nijholt <bas@nijho.lt>
Date: Tue, 23 Mar 2021 12:31:01 +0100
Subject: [PATCH 7/7] add doc-string to _coefficients

---
 adaptive/learner/integrator_coeffs.py | 1 +
 1 file changed, 1 insertion(+)

diff --git a/adaptive/learner/integrator_coeffs.py b/adaptive/learner/integrator_coeffs.py
index d590c8a29..33419c9bd 100644
--- a/adaptive/learner/integrator_coeffs.py
+++ b/adaptive/learner/integrator_coeffs.py
@@ -144,6 +144,7 @@ def calc_V(x, n):
 
 @lru_cache(maxsize=None)
 def _coefficients():
+    """Compute the coefficients on demand, in order to avoid doing linear algebra on import."""
     eps = np.spacing(1)
 
     # the nodes and Newton polynomials