aboutsummaryrefslogtreecommitdiff
path: root/vquad
diff options
context:
space:
mode:
authorChristoph Groth <christoph.groth@cea.fr>2018-02-21 16:51:41 +0100
committerChristoph Groth <christoph.groth@cea.fr>2018-02-23 12:55:00 +0100
commit581bad7324c4c1164ba4adec69c7b11b0e25fb7e (patch)
treeeece21a8d282601ad08b2d3880d90be996e79afc /vquad
parentee6289434905c0b3674e0f071ab0e27f77dc1417 (diff)
move interval quality check into refine()
Diffstat (limited to 'vquad')
-rw-r--r--vquad/core.py31
1 files changed, 18 insertions, 13 deletions
diff --git a/vquad/core.py b/vquad/core.py
index 7435334..f6baedb 100644
--- a/vquad/core.py
+++ b/vquad/core.py
@@ -159,14 +159,26 @@ class Vquad:
return children
def refine(self, ival):
- """Increase degree of interval."""
+ """Increase degree of interval.
+
+ Returns True if the refined interval is OK, and False if it is
+ borderline and should not be refined/split further. This
+ happens when neigboring points can be only barely resolved by
+ floating point numbers, or when the estimated relative error is
+ already at the limit of numerical accuracy and cannot be reduced
+ further.
+ """
ival.level += 1
points = ival.points()
vals = np.empty(points.shape)
vals[0::2] = ival.vals
vals[1::2] = self.f(points[1::2])
ival.interpolate(vals, ival.coeffs)
- return points
+
+ return (points[1] - points[0] > points[0] * min_sep
+ and points[-1] - points[-2] > points[-2] * min_sep
+ and ival.err > (abs(ival.igral)
+ * eps * tbls.V_cond_nums[ival.level]))
def improve(self):
ival = self.ivals[-1]
@@ -174,21 +186,14 @@ class Vquad:
if ival.level == max_level:
split = True
else:
- points = self.refine(ival)
- split = ival.unreliable_err
-
- if (points[1] - points[0] < points[0] * min_sep
- or points[-1] - points[-2] < points[-2] * min_sep
- or ival.err < (abs(ival.igral)
- * eps * tbls.V_cond_nums[ival.level])):
- # Remove the interval (while remembering the excess integral
- # and error), since it is either too narrow, or the estimated
- # relative error is already at the limit of numerical accuracy
- # and cannot be reduced further.
+ if not self.refine(ival):
+ # Remove the interval but remember the excess integral and
+ # error.
self.err_excess += ival.err
self.igral_excess += ival.igral
self.attic.append(self.ivals.pop())
return
+ split = ival.unreliable_err
if split:
# Replace current interval by its children.