docs: def plot(
          self,
          show_coedges=True,
          show_centroids=True,
          mesh_color="k",
          nondelaunay_edge_color="#d62728",  # mpl 2.0 default red
          boundary_edge_color=None,
          comesh_color=(0.8, 0.8, 0.8),
          show_axes=True,
      ):
          """Show the mesh using matplotlib.
          """
          # Importing matplotlib takes a while, so don't do that at the header.
          from matplotlib import pyplot as plt
          from matplotlib.collections import LineCollection
  
          # from mpl_toolkits.mplot3d import Axes3D
          fig = plt.figure()
          ax = fig.gca()
          # plt.axis("equal")
          if not show_axes:
              ax.set_axis_off()
  
          xmin = numpy.amin(self.node_coords[:, 0])
          xmax = numpy.amax(self.node_coords[:, 0])
          ymin = numpy.amin(self.node_coords[:, 1])
          ymax = numpy.amax(self.node_coords[:, 1])
  
          width = xmax - xmin
          xmin -= 0.1 * width
          xmax += 0.1 * width
  
          height = ymax - ymin
          ymin -= 0.1 * height
          ymax += 0.1 * height
  
          # ax.set_xlim(xmin, xmax)
          # ax.set_ylim(ymin, ymax)
  
          if self.edges is None:
              self.create_edges()
  
          # Get edges, cut off z-component.
          e = self.node_coords[self.edges["nodes"]][:, :, :2]
          # Plot regular edges, mark those with negative ce-ratio red.
          ce_ratios = self.ce_ratios_per_interior_edge
          pos = ce_ratios >= 0
  
          is_pos = numpy.zeros(len(self.edges["nodes"]), dtype=bool)
          is_pos[self._edge_to_edge_gid[2][pos]] = True
  
          # Mark Delaunay-conforming boundary edges
          is_pos_boundary = self.ce_ratios[self.is_boundary_edge] >= 0
          is_pos[self._edge_to_edge_gid[1][is_pos_boundary]] = True
  
          line_segments0 = LineCollection(e[is_pos], color=mesh_color)
          ax.add_collection(line_segments0)
          #
          line_segments1 = LineCollection(e[~is_pos], color=nondelaunay_edge_color)
          ax.add_collection(line_segments1)
  
          if show_coedges:
              # Connect all cell circumcenters with the edge midpoints
              cc = self.cell_circumcenters
  
              edge_midpoints = 0.5 * (
                  self.node_coords[self.edges["nodes"][:, 0]]
                  + self.node_coords[self.edges["nodes"][:, 1]]
              )
  
              # Plot connection of the circumcenter to the midpoint of all three
              # axes.
              a = numpy.stack(
                  [cc[:, :2], edge_midpoints[self.cells["edges"][:, 0], :2]], axis=1
              )
              b = numpy.stack(
                  [cc[:, :2], edge_midpoints[self.cells["edges"][:, 1], :2]], axis=1
              )
              c = numpy.stack(
                  [cc[:, :2], edge_midpoints[self.cells["edges"][:, 2], :2]], axis=1
              )
  
              line_segments = LineCollection(
                  numpy.concatenate([a, b, c]), color=comesh_color
              )
              ax.add_collection(line_segments)
  
          if boundary_edge_color:
              e = self.node_coords[self.edges["nodes"][self.is_boundary_edge_individual]][
                  :, :, :2
              ]
              line_segments1 = LineCollection(e, color=boundary_edge_color)
              ax.add_collection(line_segments1)
  
          if show_centroids:
              centroids = self.control_volume_centroids
              ax.plot(
                  centroids[:, 0],
                  centroids[:, 1],
                  linestyle="",
                  marker=".",
                  color="#d62728",
              )
  
          return fig