Add support for arbitrary dimension in Surface.is_valid

[dhermes]

Currently only 2D is supported.

Right now, we just take the 2x2 Jacobian J = J(s, t) and convert it into the determinant polynomial j(s, t) = det(J) and then make sure that polynomial always takes the same sign on the reference triangle. In the 2D case, if the surface is degree n, the components of J are degree n - 1 hence j(s, t) is degree 2(n - 1).

However, in order for a surface to be "valid", (I think) we only care that the Jacobian J has full-rank throughout (the inverse fn. thm. doesn't really cover the non-square case, I'm drawing a blank right now). So if the dimension is d, then J is d x 2 and full rank means it'll have rank 2. Instead of checking the rank of J, we can apply our determinant trick to J^T J (which is 2 x 2). In this case, each component of J^T J is degree 2(n - 1) hence the degree of jj(s, t) = det(J^T J) is 4(n - 1).

[dhermes]

When the Bezier object is the same "type" as the ambient space, this works best (due to topology / differential geometry reasons). So for example a 2D Bezier curve would be embedded as an edge of a Bezier surface/triangle in 2D:

>>> import matplotlib.pyplot as plt
>>> import numpy as np
>>> import seaborn
>>> 
>>> import bezier
>>> 
>>> nodes = np.array([
...     [-1.5, 1.625],
...     [ 2.0, 1.0  ],
...     [-2.0, 0.0  ],
...     [ 0.0, 2.0  ],
...     [-1.5, 2.125],
...     [-1.0, 2.0  ],
...     [-0.5, 2.375],
...     [-1.5, 2.625],
...     [-1.0, 2.75 ],
...     [-1.5, 3.125],
... ])
>>> 
>>> curve = bezier.Curve(nodes[:4, :], 3)
>>> surface = bezier.Surface(nodes, 3)
>>> 
>>> figure, (ax1, ax2) = plt.subplots(1, 2)
>>> 
>>> curve.plot(256, ax=ax1)
<matplotlib.axes._subplots.AxesSubplot object at 0x7fe7233dc050>
>>> _ = ax1.axis('scaled')
>>> ax1.set_xlim(-1.575, 0.075)
(-1.575, 0.075)
>>> ax1.set_ylim(0.625, 3.25)
(0.625, 3.25)
>>> 
>>> surface.plot(256, ax=ax2)
<matplotlib.axes._subplots.AxesSubplot object at 0x7fe7222e0fd0>
>>> _ = ax2.axis('scaled')
>>> ax2.set_xlim(-1.575, 0.075)
(-1.575, 0.075)
>>> ax2.set_ylim(0.625, 3.25)
(0.625, 3.25)
>>> 
>>> figure.show()

>>> surface.is_valid
False

H/T to Persson for telling me this / reminding me of an idea I had previously had.

---

$ git log -1 --pretty=%H
c74557dd593e31a16846beea19b0c600476d512e