Update contrib/python/contourpy to 1.2.0

13 files changed, 88 insertions, 17 deletions

diff --git a/contrib/python/contourpy/src/base.h b/contrib/python/contourpy/src/base.h
index 7faf9335ec1..e626e994318 100644
--- a/contrib/python/contourpy/src/base.h
+++ b/contrib/python/contourpy/src/base.h
@@ -202,6 +202,7 @@ private:
     bool _direct_line_offsets;        // Whether line offsets array is written direct to Python.
     bool _direct_outer_offsets;       // Whether outer offsets array is written direct to Python.
     bool _outer_offsets_into_points;  // Otherwise into line offsets.  Only used if _identify_holes.
+    bool _nan_separated;              // Whether adjacent lines' points are separated by nans.
     unsigned int _return_list_count;
 };
 
diff --git a/contrib/python/contourpy/src/base_impl.h b/contrib/python/contourpy/src/base_impl.h
index 61b9c461322..6175ad5c567 100644
--- a/contrib/python/contourpy/src/base_impl.h
+++ b/contrib/python/contourpy/src/base_impl.h
@@ -3,6 +3,7 @@
 
 #include "base.h"
 #include "converter.h"
+#include "util.h"
 #include <iostream>
 
 namespace contourpy {
@@ -121,6 +122,7 @@ BaseContourGenerator<Derived>::BaseContourGenerator(
       _direct_line_offsets(false),
       _direct_outer_offsets(false),
       _outer_offsets_into_points(false),
+      _nan_separated(false),
       _return_list_count(0)
 {
     if (_x.ndim() != 2 || _y.ndim() != 2 || _z.ndim() != 2)
@@ -351,8 +353,8 @@ LineType BaseContourGenerator<Derived>::default_line_type()
 template <typename Derived>
 py::sequence BaseContourGenerator<Derived>::filled(double lower_level, double upper_level)
 {
-    if (lower_level > upper_level)
-        throw std::invalid_argument("upper and lower levels are the wrong way round");
+    if (lower_level >= upper_level)
+        throw std::invalid_argument("upper_level must be larger than lower_level");
 
     _filled = true;
     _lower_level = lower_level;
@@ -367,6 +369,7 @@ py::sequence BaseContourGenerator<Derived>::filled(double lower_level, double up
     _direct_outer_offsets = (_fill_type == FillType::ChunkCombinedCodeOffset ||
                              _fill_type == FillType::ChunkCombinedOffsetOffset);
     _outer_offsets_into_points = (_fill_type == FillType::ChunkCombinedCodeOffset);
+    _nan_separated = false;
     _return_list_count = (_fill_type == FillType::ChunkCombinedCodeOffset ||
                           _fill_type == FillType::ChunkCombinedOffsetOffset) ? 3 : 2;
 
@@ -1911,6 +1914,12 @@ void BaseContourGenerator<Derived>::line(const Location& start_location, ChunkLo
     Location location = start_location;
     count_t point_count = 0;
 
+    // Insert nan if required before start of new line.
+    if (_nan_separated and local.pass > 0 && local.line_count > 0) {
+        *local.points.current++ = Util::nan;
+        *local.points.current++ = Util::nan;
+    }
+
     // finished == true indicates closed line loop.
     bool finished = follow_interior(location, start_location, local, point_count);
 
@@ -1942,7 +1951,12 @@ py::sequence BaseContourGenerator<Derived>::lines(double level)
     _direct_line_offsets = (_line_type == LineType::ChunkCombinedOffset);
     _direct_outer_offsets = false;
     _outer_offsets_into_points = false;
-    _return_list_count = (_line_type == LineType::Separate) ? 1 : 2;
+    _return_list_count = (_line_type == LineType::Separate ||
+                          _line_type == LineType::ChunkCombinedNan) ? 1 : 2;
+    _nan_separated = (_line_type == LineType::ChunkCombinedNan);
+
+    if (_nan_separated)
+        Util::ensure_nan_loaded();
 
     return march_wrapper();
 }
@@ -2073,6 +2087,14 @@ void BaseContourGenerator<Derived>::march_chunk(
         if (j_final_start < local.jend)
             _cache[local.istart + (j_final_start+1)*_nx] |= MASK_NO_MORE_STARTS;
 
+        if (_nan_separated && local.line_count > 1) {
+            // If _nan_separated, each line after the first has an extra nan to separate it from the
+            // previous line's points. If we were returning line offsets to the caller then this
+            // would need to occur in line() where the line_count is incremented. But as we are not
+            // returning line offsets it is faster just to add the extra here all at once.
+            local.total_point_count += local.line_count - 1;
+        }
+
         if (local.pass == 0) {
             if (local.total_point_count == 0) {
                 local.points.clear();
@@ -2169,8 +2191,13 @@ py::sequence BaseContourGenerator<Derived>::march_wrapper()
 
     // Return to python objects.
     if (_return_list_count == 1) {
-        assert(!_filled && _line_type == LineType::Separate);
-        return return_lists[0];
+        assert(!_filled);
+        if (_line_type == LineType::Separate)
+            return return_lists[0];
+        else {
+            assert(_line_type == LineType::ChunkCombinedNan);
+            return py::make_tuple(return_lists[0]);
+        }
     }
     else if (_return_list_count == 2)
         return py::make_tuple(return_lists[0], return_lists[1]);
@@ -2384,6 +2411,7 @@ bool BaseContourGenerator<Derived>::supports_line_type(LineType line_type)
         case LineType::SeparateCode:
         case LineType::ChunkCombinedCode:
         case LineType::ChunkCombinedOffset:
+        case LineType::ChunkCombinedNan:
             return true;
         default:
             return false;
diff --git a/contrib/python/contourpy/src/chunk_local.h b/contrib/python/contourpy/src/chunk_local.h
index f9b68ea6bed..0298ea4814b 100644
--- a/contrib/python/contourpy/src/chunk_local.h
+++ b/contrib/python/contourpy/src/chunk_local.h
@@ -21,7 +21,7 @@ struct ChunkLocal
     int pass;
 
     // Data for whole pass.
-    count_t total_point_count;
+    count_t total_point_count;           // Includes nan separators if used.
     count_t line_count;                  // Count of all lines
     count_t hole_count;                  // Count of holes only.
 
diff --git a/contrib/python/contourpy/src/fill_type.h b/contrib/python/contourpy/src/fill_type.h
index f77c765d28c..20d41aeb3b2 100644
--- a/contrib/python/contourpy/src/fill_type.h
+++ b/contrib/python/contourpy/src/fill_type.h
@@ -9,10 +9,10 @@ namespace contourpy {
 // C++11 scoped enum, must be fully qualified to use.
 enum class FillType
 {
-    OuterCode= 201,
+    OuterCode = 201,
     OuterOffset = 202,
     ChunkCombinedCode = 203,
-    ChunkCombinedOffset= 204,
+    ChunkCombinedOffset = 204,
     ChunkCombinedCodeOffset = 205,
     ChunkCombinedOffsetOffset = 206,
 };
diff --git a/contrib/python/contourpy/src/line_type.cpp b/contrib/python/contourpy/src/line_type.cpp
index 5ef2d68dac6..fbd97ea78ce 100644
--- a/contrib/python/contourpy/src/line_type.cpp
+++ b/contrib/python/contourpy/src/line_type.cpp
@@ -18,6 +18,9 @@ std::ostream &operator<<(std::ostream &os, const LineType& line_type)
         case LineType::ChunkCombinedOffset:
             os << "ChunkCombinedOffset";
             break;
+        case LineType::ChunkCombinedNan:
+            os << "ChunkCombinedNan";
+            break;
     }
     return os;
 }
diff --git a/contrib/python/contourpy/src/line_type.h b/contrib/python/contourpy/src/line_type.h
index afc6edf0c36..c1e5dda702e 100644
--- a/contrib/python/contourpy/src/line_type.h
+++ b/contrib/python/contourpy/src/line_type.h
@@ -12,7 +12,8 @@ enum class LineType
     Separate = 101,
     SeparateCode = 102,
     ChunkCombinedCode = 103,
-    ChunkCombinedOffset= 104,
+    ChunkCombinedOffset = 104,
+    ChunkCombinedNan = 105,
 };
 
 std::ostream &operator<<(std::ostream &os, const LineType& line_type);
diff --git a/contrib/python/contourpy/src/mpl2005.cpp b/contrib/python/contourpy/src/mpl2005.cpp
index dfa5c546192..a62d6e86ce4 100644
--- a/contrib/python/contourpy/src/mpl2005.cpp
+++ b/contrib/python/contourpy/src/mpl2005.cpp
@@ -48,8 +48,8 @@ Mpl2005ContourGenerator::~Mpl2005ContourGenerator()
 
 py::tuple Mpl2005ContourGenerator::filled(const double& lower_level, const double& upper_level)
 {
-    if (lower_level > upper_level)
-        throw std::invalid_argument("upper and lower levels are the wrong way round");
+    if (lower_level >= upper_level)
+        throw std::invalid_argument("upper_level must be larger than lower_level");
 
     double levels[2] = {lower_level, upper_level};
     return cntr_trace(_site, levels, 2);
diff --git a/contrib/python/contourpy/src/mpl2014.cpp b/contrib/python/contourpy/src/mpl2014.cpp
index 99107ad9ea5..d5012de6fc4 100644
--- a/contrib/python/contourpy/src/mpl2014.cpp
+++ b/contrib/python/contourpy/src/mpl2014.cpp
@@ -485,8 +485,8 @@ void Mpl2014ContourGenerator::edge_interp(
 py::tuple Mpl2014ContourGenerator::filled(
     const double& lower_level, const double& upper_level)
 {
-    if (lower_level > upper_level)
-        throw std::invalid_argument("upper and lower levels are the wrong way round");
+    if (lower_level >= upper_level)
+        throw std::invalid_argument("upper_level must be larger than lower_level");
 
     init_cache_levels(lower_level, upper_level);
 
diff --git a/contrib/python/contourpy/src/serial.cpp b/contrib/python/contourpy/src/serial.cpp
index 08d851ca349..dd69471ef64 100644
--- a/contrib/python/contourpy/src/serial.cpp
+++ b/contrib/python/contourpy/src/serial.cpp
@@ -115,6 +115,10 @@ void SerialContourGenerator::export_lines(
             // return_lists[0][local.chunk] already contains points.
             // return_lists[1][local.chunk] already contains line offsets.
             break;
+        case LineType::ChunkCombinedNan:
+            assert(has_direct_points());
+            // return_lists[0][local.chunk] already contains points.
+            break;
     }
 }
 
diff --git a/contrib/python/contourpy/src/threaded.cpp b/contrib/python/contourpy/src/threaded.cpp
index e27cd55d519..9f6eec1beda 100644
--- a/contrib/python/contourpy/src/threaded.cpp
+++ b/contrib/python/contourpy/src/threaded.cpp
@@ -208,6 +208,10 @@ void ThreadedContourGenerator::export_lines(
             // return_lists[0][local.chunk] already contains points.
             // return_lists[1][local.chunk] already contains line offsets.
             break;
+        case LineType::ChunkCombinedNan:
+            assert(has_direct_points());
+            // return_lists[0][local.chunk] already contains points.
+            break;
     }
 }
 
diff --git a/contrib/python/contourpy/src/util.cpp b/contrib/python/contourpy/src/util.cpp
index d41f0edef94..59b6e4da558 100644
--- a/contrib/python/contourpy/src/util.cpp
+++ b/contrib/python/contourpy/src/util.cpp
@@ -3,6 +3,19 @@
 
 namespace contourpy {
 
+bool Util::_nan_loaded = false;
+
+double Util::nan = 0.0;
+
+void Util::ensure_nan_loaded()
+{
+    if (!_nan_loaded) {
+        auto numpy = py::module_::import("numpy");
+        nan = numpy.attr("nan").cast<double>();
+        _nan_loaded = true;
+    }
+}
+
 index_t Util::get_max_threads()
 {
     return static_cast<index_t>(std::thread::hardware_concurrency());
diff --git a/contrib/python/contourpy/src/util.h b/contrib/python/contourpy/src/util.h
index 22dcaca7fd4..f4a1abda79a 100644
--- a/contrib/python/contourpy/src/util.h
+++ b/contrib/python/contourpy/src/util.h
@@ -8,7 +8,20 @@ namespace contourpy {
 class Util
 {
 public:
+    static void ensure_nan_loaded();
+
     static index_t get_max_threads();
+
+    // This is the NaN used internally and returned to calling functions. The value is taken from
+    // numpy rather than the standard C++ approach so that it is guaranteed to work with
+    // numpy.isnan(). The value is actually the same for many platforms, but this approach
+    // guarantees it works for all platforms that numpy supports.
+    //
+    // ensure_nan_loaded() must be called before this value is read.
+    static double nan;
+
+private:
+    static bool _nan_loaded;
 };
 
 } // namespace contourpy
diff --git a/contrib/python/contourpy/src/wrap.cpp b/contrib/python/contourpy/src/wrap.cpp
index 2c796419908..b784f814a5f 100644
--- a/contrib/python/contourpy/src/wrap.cpp
+++ b/contrib/python/contourpy/src/wrap.cpp
@@ -20,8 +20,7 @@ static contourpy::FillType mpl20xx_fill_type = contourpy::FillType::OuterCode;
 PYBIND11_MODULE(_contourpy, m) {
     m.doc() =
         "C++11 extension module wrapped using `pybind11`_.\n\n"
-        ".. note::\n"
-        "   It should not be necessary to access classes and functions in this extension module "
+        "It should not be necessary to access classes and functions in this extension module "
         "directly. Instead, :func:`contourpy.contour_generator` should be used to create "
         ":class:`~contourpy.ContourGenerator` objects, and the enums "
         "(:class:`~contourpy.FillType`, :class:`~contourpy.LineType` and "
@@ -53,11 +52,13 @@ PYBIND11_MODULE(_contourpy, m) {
     py::enum_<contourpy::LineType>(m, "LineType",
         "Enum used for ``line_type`` keyword argument in :func:`~contourpy.contour_generator`.\n\n"
         "This controls the format of contour line data returned from "
-        ":meth:`~contourpy.ContourGenerator.lines`.")
+        ":meth:`~contourpy.ContourGenerator.lines`.\n\n"
+        "``LineType.ChunkCombinedNan`` added in version 1.2.0")
         .value("Separate", contourpy::LineType::Separate)
         .value("SeparateCode", contourpy::LineType::SeparateCode)
         .value("ChunkCombinedCode", contourpy::LineType::ChunkCombinedCode)
         .value("ChunkCombinedOffset", contourpy::LineType::ChunkCombinedOffset)
+        .value("ChunkCombinedNan", contourpy::LineType::ChunkCombinedNan)
         .export_values();
 
     py::enum_<contourpy::ZInterp>(m, "ZInterp",
@@ -93,7 +94,10 @@ PYBIND11_MODULE(_contourpy, m) {
         "    upper_level (float): Upper z-level of the filled contours.\n"
         "Return:\n"
         "    Filled contour polygons as one or more sequences of numpy arrays. The exact format is "
-        "determined by the ``fill_type`` used by the ``ContourGenerator``.";
+        "determined by the ``fill_type`` used by the ``ContourGenerator``.\n\n"
+        "Raises a ``ValueError`` if ``lower_level >= upper_level``.\n\n"
+        "To return filled contours below a ``level`` use ``filled(-np.inf, level)``.\n"
+        "To return filled contours above a ``level`` use ``filled(level, np.inf)``";
     const char* line_type_doc = "Return the ``LineType``.";
     const char* lines_doc =
         "Calculate and return contour lines at a particular level.\n\n"