Hot-keys on this page
r m x p toggle line displays
j k next/prev highlighted chunk
0 (zero) top of page
1 (one) first highlighted chunk
r""" Graphics 3D Object for Representing and Triangulating Isosurfaces.
AUTHORS:
- Robert Hanson (2007): initial Java version, in Jmol. - Carl Witty (2009-01): first Cython version. - Bill Cauchois (2009): improvements for inclusion into Sage. """
#***************************************************************************** # Copyright (C) 2009 Carl Witty <Carl.Witty@gmail.com> # # Distributed under the terms of the GNU General Public License (GPL) # # This code is distributed in the hope that it will be useful, # but WITHOUT ANY WARRANTY; without even the implied warranty of # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU # General Public License for more details. # # The full text of the GPL is available at: # # http://www.gnu.org/licenses/ #*****************************************************************************
# Pieces of this file are strongly based on the marching cubes # implementation in Jmol located at src/org/jmol/jvxl/calc/MarchingCubes.java. # The data tables are in fact directly copied from there.
#***************************************************************************** # This copyright is inherited from MarchingCubes.java in the Jmol source code. # # * Copyright (C) 2007 Miguel, Bob, Jmol Development # * # * Contact: hansonr@stolaf.edu # * # * This library is free software; you can redistribute it and/or # * modify it under the terms of the GNU Lesser General Public # * License as published by the Free Software Foundation; either # * version 2.1 of the License, or (at your option) any later version. # * # * This library is distributed in the hope that it will be useful, # * but WITHOUT ANY WARRANTY; without even the implied warranty of # * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU # * Lesser General License for more details. # * # * You should have received a copy of the GNU Lesser General Public # * License along with this library; if not, write to the Free Software # * Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA. #*****************************************************************************
# There's a framework in here for computing multiple isosurfaces of a # single function. Currently, it's only used for a single # implicit_plot3d where contour=... is given a list, but it's ready to # be extended. I think the best UI would be if animate(...) and # show(...) had a prepass where they went through their rendering # trees, found ImplicitSurface objects with the same function, # bounding box, and plot_points (other parameters, such as contours, # hole, color, would be allowed to be different), # and arranged to have them all plotted simultaneously. These # prepasses have not been written yet. Combining multiple # ImplicitSurface plots would be particularly advantageous for animate(...), # but for a big animation, watch out for memory usage.
# If you have a reasonably simple surface (not a space-filling fractal), # then if n is your resolution, we have n^3 evaluations of the main # function, about n^2 evaluations of auxiliary functions (hole, gradient, # color), and output of size about n^2.
# With this in mind, we pay particular attention to optimizing the n^3 # function evaluations. (But keep in mind that n may be as small as 20 # or so, so we shouldn't ignore the efficiency of the n^2 parts.)
# Also, we make sure not to ever allocate O(n^3) memory; we do the # computation a slice at a time to avoid this. (Jmol always allocates # n^3 memory when it reads the JVXL file, but that might be on a different # computer; Tachyon would only allocate memory proportional to the # output size.) from __future__ import absolute_import
cimport numpy as np
from sage.plot.plot3d.transform cimport point_c, face_c, color_c, Transformation from sage.plot.plot3d.base cimport PrimitiveObject from sage.plot.plot3d.index_face_set cimport IndexFaceSet
from sage.libs.gsl.math cimport gsl_isnan
include "point_c.pxi"
# The default value for plot_points (a keyword argument to implicit_plot3d()), # assumed when plot_points is set to "automatic".
cdef inline bint marching_has_edge(double a, double b, double contour, double *f, bint *has_nan): # XXX Would be nicer to use isnan(), because it's inlined. # Is it portable enough? has_nan[0] = True return False
# Returns 0 or 1 cdef inline int marching_is_inside(double v, double contour):
cdef void interpolate_point_c(point_c *result, double frac, point_c *inputs):
cdef class VertexInfo: # The point in "integer space" cdef point_c pt
# The gradient at this point in "evaluation space" cdef point_c gradient
# (R,G,B) color ; a color triple cdef color_c color
# This point in "evaluation space" cdef point_c eval_pt
cdef void update_eval_pt(self, point_c *eval_min, point_c *eval_scale): """ Use eval_min and eval_scale to transform self.pt into evaluation space and store the result in self.eval_pt. """
def __repr__(self): """ TESTS::
sage: from sage.plot.plot3d.implicit_surface import VertexInfo sage: VertexInfo() <0.0, 0.0, 0.0> """
cdef mk_VertexInfo(double x, double y, double z, point_c *eval_min, point_c *eval_scale): cdef VertexInfo v
cdef class MarchingCubes: r""" Handles marching cube rendering.
Protocol:
1. Create the class. 2. Call ``process_slice`` once for each X slice, from self.nx > x >= 0. 3. Call ``finish()``, which returns a list of strings.
.. NOTE::
Actually, only 4 slices ever exist; the caller will re-use old storage. """
cdef readonly object xrange cdef readonly object yrange cdef readonly object zrange cdef readonly double contour cdef readonly int nx cdef readonly int ny cdef readonly int nz cdef readonly Transformation transform cdef readonly object region cdef readonly object gradient cdef readonly bint smooth cdef readonly object color_function cdef readonly object colormap cdef readonly object results
# We deal with three coordinate systems. We do most of our work # in an integral coordinate system where x ranges from # 0 <= x < self.nx, etc.; we do function evaluations where # self.xrange[0] <= x <= self.xrange[1], etc.; and output # is in a third coordinate system. # (Note that in "integer space", function evaluations of the main # function happen on integer coordinates; but function evaluations # of the auxiliary functions will have one non-integer coordinate.)
# These parameters convert from integer space to function-evaluation # space: eval_coord = eval_min + int_coord * eval_scale cdef point_c eval_min cdef point_c eval_scale
# The componentwise reciprocal of eval_scale; just used to change # some divisions into multiplications cdef point_c eval_scale_inv
cdef point_c out_origin, out_plus_x, out_plus_y, out_plus_z
def __init__(self, xrange, yrange, zrange, contour, plot_points, transform=None, region=None, gradient=None, smooth=True, color_function=None, colormap=None): """ TESTS:
Marching cubes is an abstract base class, you can't do anything with it::
sage: from sage.plot.plot3d.implicit_surface import MarchingCubes sage: cube_marcher = MarchingCubes((0, 1), (0, 1), (0, 1), 1, (10, 10, 10)) """
cdef point_c zero_pt, origin, plus_x, plus_y, plus_z self.transform.transform_point_c(&self.out_origin, origin) self.transform.transform_point_c(&self.out_plus_x, plus_x) self.transform.transform_point_c(&self.out_plus_y, plus_y) self.transform.transform_point_c(&self.out_plus_z, plus_z) else:
def finish(self): """ Return the results of the marching cubes algorithm as a list.
The format is specific to the subclass implementing this method.
TESTS:
By default, it returns an empty list::
sage: from sage.plot.plot3d.implicit_surface import MarchingCubes sage: cube_marcher = MarchingCubes((0, 1), (0, 1), (0, 1), 1, (10, 10, 10), None) sage: cube_marcher.finish() [] """
cdef class MarchingCubesTriangles(MarchingCubes): """ A subclass of MarchingCubes that returns its results as a list of triangles, including their vertices and normals (if ``smooth=True``).
And also their vertex colors if a vertex coloring function is given. """
cdef readonly np.ndarray x_vertices cdef readonly np.ndarray y_vertices cdef readonly np.ndarray z_vertices
cdef readonly np.ndarray y_vertices_swapped cdef readonly np.ndarray z_vertices_swapped
cdef readonly slices
def __init__(self, *args, **kwargs): """ TESTS::
sage: from sage.plot.plot3d.implicit_surface import MarchingCubesTriangles sage: cube_marcher = MarchingCubesTriangles((0, 1), (0, 1), (0, 1), 1, (10, 10, 10)) """
# Create new arrays that share the same underlying data, but # have 0 and 1 reversed for the first coordinate.
def process_slice(self, unsigned int x, np.ndarray slice): """ Process a single slice of function evaluations at the specified `x` coordinate.
EXAMPLES::
sage: from sage.plot.plot3d.implicit_surface import MarchingCubesTriangles sage: import numpy as np sage: cube_marcher = MarchingCubesTriangles((-2, 2), (-2, 2), (-2, 2), 4, (10,)*3, smooth=False) sage: f = lambda x, y, z: x^2 + y^2 + z^2 sage: slices = np.zeros((10, 10, 10), dtype=np.double) sage: for x in reversed(range(0, 10)): ....: for y in range(0, 10): ....: for z in range(0, 10): ....: slices[x, y, z] = f(*[a * (4 / 9) -2 for a in (x, y, z)]) ....: cube_marcher.process_slice(x, slices[x, :, :]) sage: faces = cube_marcher.finish() sage: faces[0][0] {'x': 1.555555555555..., 'y': -1.111111111111..., 'z': -0.555555555555...}
We render the isosurface using IndexFaceSet::
sage: from sage.plot.plot3d.index_face_set import IndexFaceSet sage: IndexFaceSet([tuple((p['x'], p['y'], p['z']) for p in face) for face in faces]) Graphics3d Object """ # We go to a lot of effort to avoid repeating computations. # (I hope that the related bookkeeping is actually faster # than repeating the computations, but I haven't tested.)
# We get evaluation points, one slice at a time. But we # don't want to process slices of evaluation points; # we want to process slices of cubes. (There is one fewer slice # of cubes than there is of evaluation points.)
# Without any caching, we would repeat this for each cube: # Find which vertices are inside/outside the surface. # For edges which cross the surface, interpolate to # find the exact crossing location; these will be the vertices # of triangles. # Compute the color and surface normal for each of these vertices. # Assemble the vertices into triangles, using the "marching cubes" # tables.
# The problem with this is that each vertex is actually shared # among four neighbor cubes (except on the edge of the array), # so we would be repeating work (potentially a lot of work, if # the user-specified gradient or color routines are expensive) # four times.
# So we cache the information associated with each vertex. # Let's call an edge an X, Y, or Z edge depending on which # axis it is parallel to; and a vertex is an X, Y, or Z # vertex depending on which kind of edge it lies on.
# X vertices are shared between cubes that are all within a # single cube slice. However, Y and Z vertices are shared # between adjacent slices, so we need to keep those caches # around. But we still need only two slices of vertex cache # at a time (the two that are adjacent to the current cube slice).
# To reduce memory allocations, we allocate these caches at the # beginning of the run (they are ndarrays of objects). We # set a VertexInfo object in the cache when we first compute it; # then we look in the cache when we need it for the other three # cubes.
# The X vertex cache is a 2-D ndarray. The Y and Z vertex caches # are 3-D, where the first (x) dimension is indexed by 0 or 1.
# The Cython buffer interface (that we use for fast access to # ndarray's) has to reinitialize itself separately in each # function. For this reason, we use fewer, bigger functions; # in particular, we don't call any functions per-cube that # require buffer access.
# To compute interpolated gradients using central differencing, # we need up to 4 slices. Consider a gradient computed at # an X vertex, between slices 1 and 2. This is interpolated # between the gradient at slice 1 and the gradient at slice 2. # Computing the gradient at slice 1 requires slices 0 and 2; # computing the gradient at slice 2 requires slices 1 and 3.
# This means we need to queue up slices and process them # in a somewhat strange order. This function only does # the queuing, and then passes all the real work off to # other functions.
None,
None,
cpdef _update_yz_vertices(self, int x, np.ndarray _prev, np.ndarray _cur, np.ndarray _next): """ TESTS::
sage: from sage.plot.plot3d.implicit_surface import MarchingCubesTriangles sage: import numpy as np sage: cube_marcher = MarchingCubesTriangles((0, 1), (0, 1), (0, 1), 0, (2,)*3, smooth=False) sage: prev_slice = next_slice = np.ones((2, 2), dtype=np.double) sage: cur_slice = prev_slice.copy() sage: cur_slice[0, 0] = -1 # Seed the slice data with an "interesting" value. sage: cube_marcher._update_yz_vertices(1, prev_slice, cur_slice, next_slice) sage: cube_marcher.z_vertices.tolist() [[[<1.0, 0.0, 0.5>], [None]], [[None], [None]]] sage: cube_marcher.y_vertices.tolist() [[[<1.0, 0.5, 0.0>, None]], [[None, None]]] sage: cube_marcher.x_vertices.any() # This shouldn't affect the X vertices. """
cdef int y cdef int z cdef VertexInfo v cdef double frac cdef bint has_nan cdef point_c gradients[2] cdef int i # We must compute a gradient. self.apply_point_func(&v.gradient, self.gradient, v) else: # Use central differencing. x, y+i, z, else:
# OK, that updated the Y vertices. Now we do almost exactly # the same thing to update Z vertices. # We must compute a gradient. self.apply_point_func(&v.gradient, self.gradient, v) else: # Use central differencing. x, y, z+i, else:
cpdef _update_x_vertices(self, int x, np.ndarray _prev, np.ndarray _left, np.ndarray _right, np.ndarray _next): """ TESTS::
sage: from sage.plot.plot3d.implicit_surface import MarchingCubesTriangles sage: import numpy as np sage: cube_marcher = MarchingCubesTriangles((0, 1), (0, 1), (0, 1), 0, (4, 2, 2), smooth=False) sage: prev_slice = right_slice = next_slice = np.ones((2, 2), dtype=np.double) sage: left_slice = prev_slice.copy() sage: left_slice[1, 1] = -1 sage: cube_marcher._update_x_vertices(1, prev_slice, left_slice, right_slice, next_slice) sage: cube_marcher.x_vertices.tolist() [[None, None], [None, <1.5, 1.0, 1.0>]] sage: cube_marcher.y_vertices.any() or cube_marcher.z_vertices.any() # This shouldn't affect the Y or Z vertices. """
cdef int y cdef int z
cdef VertexInfo v cdef double frac cdef bint has_nan cdef point_c gradients[2] # We must compute a gradient. self.apply_point_func(&v.gradient, self.gradient, v) else: # Use central differencing. x, y, z, x, y, z, else:
cdef bint in_region(self, VertexInfo v):
cdef apply_point_func(self, point_c *pt, fn, VertexInfo v): if isinstance(fn, tuple): pt[0].x = fn[0](v.eval_pt.x, v.eval_pt.y, v.eval_pt.z) pt[0].y = fn[1](v.eval_pt.x, v.eval_pt.y, v.eval_pt.z) pt[0].z = fn[2](v.eval_pt.x, v.eval_pt.y, v.eval_pt.z) else: t = fn(v.eval_pt.x, v.eval_pt.y, v.eval_pt.z) pt[0].x = t[0] pt[0].y = t[1] pt[0].z = t[2]
cdef apply_color_func(self, color_c *pt, fn, cm, VertexInfo v):
cdef get_gradient(self, point_c *g, int x, int y, int z, double center, double lx, double ux, double ly, double uy, double lz, double uz): # What a mess! It would be much nicer-looking code to pass slices # in here and do the subscripting in here. Unfortunately, # that would also be much slower, because we'd have to re-initialize # the Cython buffer interface on each call.
cpdef process_cubes(self, np.ndarray _left, np.ndarray _right): """ TESTS::
sage: from sage.plot.plot3d.implicit_surface import MarchingCubesTriangles sage: import numpy as np sage: cube_marcher = MarchingCubesTriangles((0, 1), (0, 1), (0, 1), 0, (3, 2, 2), smooth=False) sage: slices = [np.ones((2, 2), dtype=np.double) for i in range(3)] sage: slices[0][1, 1] = -1 sage: cube_marcher._update_yz_vertices(0, None, slices[0], slices[1]) sage: cube_marcher._update_x_vertices(0, None, slices[0], slices[1], slices[2]) sage: cube_marcher.process_cubes(slices[0], slices[1]) sage: cube_marcher.finish() [({'x': 0.0, 'y': 1.0, 'z': 0.5}, {'x': 0.25, 'y': 1.0, 'z': 1.0}, {'x': 0.0, 'y': 0.5, 'z': 1.0})] """
cdef int y cdef int z
# based on generateSurfaceData in MarchingCubes.java cdef int insideMask
# Cool ASCII art from MarchingCubes.java: # * Y # * 4 --------4--------- 5 # * /| /| # * / | / | # * / | / | # * 7 8 5 | # * / | / 9 # * / | / | # * 7 --------6--------- 6 | # * | | | | # * | 0 ---------0--|----- 1 X # * | / | / # * 11 / 10 / # * | 3 | 1 # * | / | / # * | / | / # * 3 ---------2-------- 2 # * Z
# We see from the above that vertices are labeled 0 to 7, and # edges are labeled 0 to 11.
cdef tuple my_triangles
cdef int i
# For each vertex (0 to 7), set the corresponding bit # of insideMask iff the vertex is inside the surface.
# OK, we have a cube on the surface. Copy all of the vertex # info into an array for easier reference.
# In wireframe mode, my_triangles[i+3] specifies # whether or not to draw the corresponding edge. # See MarchingCubes.java for details.
cpdef add_triangle(self, VertexInfo v1, VertexInfo v2, VertexInfo v3): """ Called when a new triangle is generated by the marching cubes algorithm to update the results array.
TESTS::
sage: from sage.plot.plot3d.implicit_surface import MarchingCubesTriangles, VertexInfo sage: cube_marcher = MarchingCubesTriangles((0, 1), (0, 1), (0, 1), 0, (10,)*3, smooth=False) sage: cube_marcher.add_triangle(VertexInfo(), VertexInfo(), VertexInfo()) sage: cube_marcher.finish() [({'x': 0.0, 'y': 0.0, 'z': 0.0}, {'x': 0.0, 'y': 0.0, 'z': 0.0}, {'x': 0.0, 'y': 0.0, 'z': 0.0})] """ # This happens if there is a NaN nearby, or if a hole was # specified here.
cdef: point_c v1_ev_pt, v2_ev_pt, v3_ev_pt point_c n1_ev_vec, n2_ev_vec, n3_ev_vec color_c v1_col, v2_col, v3_col
self.transform.transform_point_c(&v1_ev_pt, v1.eval_pt) self.transform.transform_point_c(&v2_ev_pt, v2.eval_pt) self.transform.transform_point_c(&v3_ev_pt, v3.eval_pt) else:
# XXX I believe this is wrong for non-uniform transforms self.transform.transform_vector_c(&n1_ev_vec, v1.gradient) self.transform.transform_vector_c(&n2_ev_vec, v2.gradient) self.transform.transform_vector_c(&n3_ev_vec, v3.gradient) else:
cpdef render_implicit(f, xrange, yrange, zrange, plot_points, cube_marchers): """ INPUT:
- ``f`` - a (fast!) callable function
- ``xrange`` - a 2-tuple (x_min, x_max)
- ``yrange`` - a 2-tuple (y_min, y_may)
- ``zrange`` - a 2-tuple (z_min, z_maz)
- ``plot_points`` - a triple of integers indicating the number of function evaluations in each direction.
- ``cube_marchers`` - a list of cube marchers, one for each contour.
OUTPUT:
A representation of the isosurface, in the format specified by the individual cube marchers.
TESTS::
sage: from sage.plot.plot3d.implicit_surface import render_implicit, MarchingCubesTriangles sage: plot_points, f = (40,)*3, lambda x, y, z: x + y + z sage: cube_marcher = MarchingCubesTriangles((0, 1), (0, 1), (0, 1), 1, (10,)*3) sage: results = render_implicit(lambda x, y, z: x + y + z, ....: (0, 1), (0, 1), (0, 1), (10,)*3, [cube_marcher]) sage: results[0][0] {'x': 1.0, 'y': 0.0, 'z': 0.0} """
cdef point_c eval_min, eval_scale
# A possible future extension would be to allow passing in "f" # as a numpy ndarray. If that were done, we could slightly modify # the following code to just pass slices of f to the renderers # (no copying of the data would be required).
# The current marching cube renderers need only at most four slices at # a time.
cdef np.ndarray[double, ndim=2] slice
cdef unsigned int x, y, z cdef int n
cdef double eval_x, eval_y, eval_z
cdef MarchingCubes marcher
cdef class ImplicitSurface(IndexFaceSet): cdef readonly object f cdef readonly object vars cdef readonly tuple xrange cdef readonly tuple yrange cdef readonly tuple zrange cdef readonly list contours cdef readonly object region cdef readonly bint smooth cdef readonly object gradient cdef readonly object color_function cdef readonly object colormap cdef readonly tuple plot_points
def __init__(self, f, xrange, yrange, zrange, contour=0, plot_points="automatic", region=None, smooth=True, gradient=None, **kwds): """ TESTS::
sage: from sage.plot.plot3d.implicit_surface import ImplicitSurface sage: var('x,y,z') (x, y, z) sage: G = ImplicitSurface(x^2 + y^2 + z^2, (x,-2, 2), (y,-2, 2), (z,-2, 2), contour=4) sage: show(G)
A colored case::
sage: t = (1-sin(2*x*y+3*z)**2).function(x,y,z) sage: cm = colormaps.autumn sage: G = ImplicitSurface(x^2 + y^2 + z^2, (x,-2, 2), (y,-2, 2), (z,-2, 2), contour=4, color=(t,cm)) sage: G.show(viewer='tachyon') """ except (TypeError, AttributeError): pass # case of a global color else: # case of a color depending on parameters
else:
# Comments from Carl Witty, who first wrote this some of this code # See Trac 9483 # When I first wrote the code, I had the idea to create a # direct-to-tachyon backend that would use vertex normals # to create much nicer-looking plots with smaller numbers # of plot_points, and laid the groundwork for this backend # with the gradient and smooth arguments. But I abandoned the # code without writing this backend (and leaving many other parts # of the code unfinished). # When William Cauchois took over and finished the code (thank you, # William!), he only wrote an IndexFaceSet backend, that can't # (currently) make use of vertex normals. So the gradient code is # currently useless. # But it's still open for somebody to write a direct-to-tachyon backend, # or to extend IndexFaceSet to support vertex normals.
# Since IndexFaceSet doesn't even support smooth shading, we overwrite # the passed-in smooth parameter.
try: gradient = (orig_f.diff(self.vars[0]), orig_f.diff(self.vars[1]), orig_f.diff(self.vars[2])) except Exception: # Would be nice to have more nuanced error handling here.
# If anything goes wrong, we'll just use central differencing. pass else: raise ValueError
def bounding_box(self): """ Return a bounding box for the ``ImplicitSurface``, as a tuple of two 3-dimensional points.
EXAMPLES:
Note that the bounding box corresponds exactly to the x-, y-, and z- range::
sage: from sage.plot.plot3d.implicit_surface import ImplicitSurface sage: G = ImplicitSurface(0, (0, 1), (0, 1), (0, 1)) sage: G.bounding_box() ((0.0, 0.0, 0.0), (1.0, 1.0, 1.0)) """
def obj_repr(self, render_params): """ Return a representation of this object in the .obj format.
TESTS:
We graph a simple plane::
sage: from sage.plot.plot3d.implicit_surface import ImplicitSurface sage: var('x,y,z') (x, y, z) sage: G = ImplicitSurface(x + y + z, (x,-1, 1), (y,-1, 1), (z,-1, 1)) sage: obj = G.obj_repr(G.default_render_params()) sage: vertices = obj[2]
The number of vertices in the OBJ representation should equal the number of vertices in the face set::
sage: len(vertices) == len(G.vertex_list()) True
The vertices in the OBJ representation should also be approximately equal to the vertices in the face set -- the small error is due to rounding which occurs during output (we test only the first 20 points for the sake of speed)::
sage: def points_equal(a, b, epsilon=(1e-5)): ....: return all(abs(x0-x1) < epsilon for x0, x1 in zip(a, b)) sage: checklist = [] sage: assert len(vertices) >= 20 # I should hope so, we're rendering at the default resolution! sage: for vertex, surf_vertex in list(zip(vertices, G.vertex_list()))[0:20]: ....: checklist.append(points_equal(map(float, vertex.split(' ')[1:]), surf_vertex)) sage: all(checklist) True """
def tachyon_repr(self, render_params): """ Return a representation of this object suitable for use with the Tachyon renderer.
TESTS::
sage: from sage.plot.plot3d.implicit_surface import ImplicitSurface sage: var('x,y,z') (x, y, z) sage: G = ImplicitSurface(x + y + z, (x,-1, 1), (y,-1, 1), (z,-1, 1)) sage: G.tachyon_repr(G.default_render_params())[0].startswith('TRI') True """
def jmol_repr(self, render_params): """ Return a representation of this object suitable for use with the Jmol renderer.
TESTS::
sage: from sage.plot.plot3d.implicit_surface import ImplicitSurface sage: var('x,y,z') (x, y, z) sage: G = ImplicitSurface(x + y + z, (x,-1, 1), (y,-1, 1), (z,-1, 1)) sage: show(G, viewer='jmol') # indirect doctest """
def json_repr(self, render_params): """ Return a representation of this object in JavaScript Object Notation (JSON).
TESTS::
sage: from sage.plot.plot3d.implicit_surface import ImplicitSurface sage: var('x,y,z') (x, y, z) sage: G = ImplicitSurface(x + y + z, (x,-1, 1), (y,-1, 1), (z,-1, 1)) sage: G.json_repr(G.default_render_params())[0].startswith('{"vertices":') True """
def triangulate(self, force=False): """ The IndexFaceSet will be empty until you call this method, which generates the faces and vertices according to the parameters specified in the constructor for ImplicitSurface.
Note that if you call this method more than once, subsequent invocations will have no effect (this is an optimization to avoid repeated work) unless you specify ``force=True`` in the keywords.
EXAMPLES::
sage: from sage.plot.plot3d.implicit_surface import ImplicitSurface sage: var('x,y,z') (x, y, z) sage: G = ImplicitSurface(x + y + z, (x,-1, 1), (y,-1, 1), (z,-1, 1)) sage: len(G.vertex_list()), len(G.face_list()) (0, 0) sage: G.triangulate() sage: len(G.vertex_list()) > 0, len(G.face_list()) > 0 (True, True) sage: G.show() # This should be fast, since the mesh is already triangulated.
""" # The mesh is already triangulated
cdef: face_c* dest_face point_c* dest_vertex
# ct is the mean of the colors of vertices
# Data table (courtesy of MarchingCubes.java) ( 0, 8, 3, 7 ), ( 0, 1, 9, 7 ), ( 1, 8, 3, 6, 9, 8, 1, 5 ), ( 1, 2, 10, 7 ), ( 0, 8, 3, 7, 1, 2, 10, 7 ), ( 9, 2, 10, 6, 0, 2, 9, 5 ), ( 2, 8, 3, 6, 2, 10, 8, 1, 10, 9, 8, 3 ), ( 3, 11, 2, 7 ), ( 0, 11, 2, 6, 8, 11, 0, 5 ), ( 1, 9, 0, 7, 2, 3, 11, 7 ), ( 1, 11, 2, 6, 1, 9, 11, 1, 9, 8, 11, 3 ), ( 3, 10, 1, 6, 11, 10, 3, 5 ), ( 0, 10, 1, 6, 0, 8, 10, 1, 8, 11, 10, 3 ), ( 3, 9, 0, 6, 3, 11, 9, 1, 11, 10, 9, 3 ), ( 9, 8, 10, 5, 10, 8, 11, 6 ), ( 4, 7, 8, 7 ), ( 4, 3, 0, 6, 7, 3, 4, 5 ), ( 0, 1, 9, 7, 8, 4, 7, 7 ), ( 4, 1, 9, 6, 4, 7, 1, 1, 7, 3, 1, 3 ), ( 1, 2, 10, 7, 8, 4, 7, 7 ), ( 3, 4, 7, 6, 3, 0, 4, 3, 1, 2, 10, 7 ), ( 9, 2, 10, 6, 9, 0, 2, 3, 8, 4, 7, 7 ), ( 2, 10, 9, 3, 2, 9, 7, 0, 2, 7, 3, 6, 7, 9, 4, 6 ), ( 8, 4, 7, 7, 3, 11, 2, 7 ), ( 11, 4, 7, 6, 11, 2, 4, 1, 2, 0, 4, 3 ), ( 9, 0, 1, 7, 8, 4, 7, 7, 2, 3, 11, 7 ), ( 4, 7, 11, 3, 9, 4, 11, 1, 9, 11, 2, 2, 9, 2, 1, 6 ), ( 3, 10, 1, 6, 3, 11, 10, 3, 7, 8, 4, 7 ), ( 1, 11, 10, 6, 1, 4, 11, 0, 1, 0, 4, 3, 7, 11, 4, 5 ), ( 4, 7, 8, 7, 9, 0, 11, 1, 9, 11, 10, 6, 11, 0, 3, 6 ), ( 4, 7, 11, 3, 4, 11, 9, 4, 9, 11, 10, 6 ), ( 9, 5, 4, 7 ), ( 9, 5, 4, 7, 0, 8, 3, 7 ), ( 0, 5, 4, 6, 1, 5, 0, 5 ), ( 8, 5, 4, 6, 8, 3, 5, 1, 3, 1, 5, 3 ), ( 1, 2, 10, 7, 9, 5, 4, 7 ), ( 3, 0, 8, 7, 1, 2, 10, 7, 4, 9, 5, 7 ), ( 5, 2, 10, 6, 5, 4, 2, 1, 4, 0, 2, 3 ), ( 2, 10, 5, 3, 3, 2, 5, 1, 3, 5, 4, 2, 3, 4, 8, 6 ), ( 9, 5, 4, 7, 2, 3, 11, 7 ), ( 0, 11, 2, 6, 0, 8, 11, 3, 4, 9, 5, 7 ), ( 0, 5, 4, 6, 0, 1, 5, 3, 2, 3, 11, 7 ), ( 2, 1, 5, 3, 2, 5, 8, 0, 2, 8, 11, 6, 4, 8, 5, 5 ), ( 10, 3, 11, 6, 10, 1, 3, 3, 9, 5, 4, 7 ), ( 4, 9, 5, 7, 0, 8, 1, 5, 8, 10, 1, 2, 8, 11, 10, 3 ), ( 5, 4, 0, 3, 5, 0, 11, 0, 5, 11, 10, 6, 11, 0, 3, 6 ), ( 5, 4, 8, 3, 5, 8, 10, 4, 10, 8, 11, 6 ), ( 9, 7, 8, 6, 5, 7, 9, 5 ), ( 9, 3, 0, 6, 9, 5, 3, 1, 5, 7, 3, 3 ), ( 0, 7, 8, 6, 0, 1, 7, 1, 1, 5, 7, 3 ), ( 1, 5, 3, 5, 3, 5, 7, 6 ), ( 9, 7, 8, 6, 9, 5, 7, 3, 10, 1, 2, 7 ), ( 10, 1, 2, 7, 9, 5, 0, 5, 5, 3, 0, 2, 5, 7, 3, 3 ), ( 8, 0, 2, 3, 8, 2, 5, 0, 8, 5, 7, 6, 10, 5, 2, 5 ), ( 2, 10, 5, 3, 2, 5, 3, 4, 3, 5, 7, 6 ), ( 7, 9, 5, 6, 7, 8, 9, 3, 3, 11, 2, 7 ), ( 9, 5, 7, 3, 9, 7, 2, 0, 9, 2, 0, 6, 2, 7, 11, 6 ), ( 2, 3, 11, 7, 0, 1, 8, 5, 1, 7, 8, 2, 1, 5, 7, 3 ), ( 11, 2, 1, 3, 11, 1, 7, 4, 7, 1, 5, 6 ), ( 9, 5, 8, 5, 8, 5, 7, 6, 10, 1, 3, 3, 10, 3, 11, 6 ), ( 5, 7, 0, 1, 5, 0, 9, 6, 7, 11, 0, 1, 1, 0, 10, 5, 11, 10, 0, 1 ), ( 11, 10, 0, 1, 11, 0, 3, 6, 10, 5, 0, 1, 8, 0, 7, 5, 5, 7, 0, 1 ), ( 11, 10, 5, 3, 7, 11, 5, 5 ), ( 10, 6, 5, 7 ), ( 0, 8, 3, 7, 5, 10, 6, 7 ), ( 9, 0, 1, 7, 5, 10, 6, 7 ), ( 1, 8, 3, 6, 1, 9, 8, 3, 5, 10, 6, 7 ), ( 1, 6, 5, 6, 2, 6, 1, 5 ), ( 1, 6, 5, 6, 1, 2, 6, 3, 3, 0, 8, 7 ), ( 9, 6, 5, 6, 9, 0, 6, 1, 0, 2, 6, 3 ), ( 5, 9, 8, 3, 5, 8, 2, 0, 5, 2, 6, 6, 3, 2, 8, 5 ), ( 2, 3, 11, 7, 10, 6, 5, 7 ), ( 11, 0, 8, 6, 11, 2, 0, 3, 10, 6, 5, 7 ), ( 0, 1, 9, 7, 2, 3, 11, 7, 5, 10, 6, 7 ), ( 5, 10, 6, 7, 1, 9, 2, 5, 9, 11, 2, 2, 9, 8, 11, 3 ), ( 6, 3, 11, 6, 6, 5, 3, 1, 5, 1, 3, 3 ), ( 0, 8, 11, 3, 0, 11, 5, 0, 0, 5, 1, 6, 5, 11, 6, 6 ), ( 3, 11, 6, 3, 0, 3, 6, 1, 0, 6, 5, 2, 0, 5, 9, 6 ), ( 6, 5, 9, 3, 6, 9, 11, 4, 11, 9, 8, 6 ), ( 5, 10, 6, 7, 4, 7, 8, 7 ), ( 4, 3, 0, 6, 4, 7, 3, 3, 6, 5, 10, 7 ), ( 1, 9, 0, 7, 5, 10, 6, 7, 8, 4, 7, 7 ), ( 10, 6, 5, 7, 1, 9, 7, 1, 1, 7, 3, 6, 7, 9, 4, 6 ), ( 6, 1, 2, 6, 6, 5, 1, 3, 4, 7, 8, 7 ), ( 1, 2, 5, 5, 5, 2, 6, 6, 3, 0, 4, 3, 3, 4, 7, 6 ), ( 8, 4, 7, 7, 9, 0, 5, 5, 0, 6, 5, 2, 0, 2, 6, 3 ), ( 7, 3, 9, 1, 7, 9, 4, 6, 3, 2, 9, 1, 5, 9, 6, 5, 2, 6, 9, 1 ), ( 3, 11, 2, 7, 7, 8, 4, 7, 10, 6, 5, 7 ), ( 5, 10, 6, 7, 4, 7, 2, 1, 4, 2, 0, 6, 2, 7, 11, 6 ), ( 0, 1, 9, 7, 4, 7, 8, 7, 2, 3, 11, 7, 5, 10, 6, 7 ), ( 9, 2, 1, 6, 9, 11, 2, 2, 9, 4, 11, 1, 7, 11, 4, 5, 5, 10, 6, 7 ), ( 8, 4, 7, 7, 3, 11, 5, 1, 3, 5, 1, 6, 5, 11, 6, 6 ), ( 5, 1, 11, 1, 5, 11, 6, 6, 1, 0, 11, 1, 7, 11, 4, 5, 0, 4, 11, 1 ), ( 0, 5, 9, 6, 0, 6, 5, 2, 0, 3, 6, 1, 11, 6, 3, 5, 8, 4, 7, 7 ), ( 6, 5, 9, 3, 6, 9, 11, 4, 4, 7, 9, 5, 7, 11, 9, 1 ), ( 10, 4, 9, 6, 6, 4, 10, 5 ), ( 4, 10, 6, 6, 4, 9, 10, 3, 0, 8, 3, 7 ), ( 10, 0, 1, 6, 10, 6, 0, 1, 6, 4, 0, 3 ), ( 8, 3, 1, 3, 8, 1, 6, 0, 8, 6, 4, 6, 6, 1, 10, 6 ), ( 1, 4, 9, 6, 1, 2, 4, 1, 2, 6, 4, 3 ), ( 3, 0, 8, 7, 1, 2, 9, 5, 2, 4, 9, 2, 2, 6, 4, 3 ), ( 0, 2, 4, 5, 4, 2, 6, 6 ), ( 8, 3, 2, 3, 8, 2, 4, 4, 4, 2, 6, 6 ), ( 10, 4, 9, 6, 10, 6, 4, 3, 11, 2, 3, 7 ), ( 0, 8, 2, 5, 2, 8, 11, 6, 4, 9, 10, 3, 4, 10, 6, 6 ), ( 3, 11, 2, 7, 0, 1, 6, 1, 0, 6, 4, 6, 6, 1, 10, 6 ), ( 6, 4, 1, 1, 6, 1, 10, 6, 4, 8, 1, 1, 2, 1, 11, 5, 8, 11, 1, 1 ), ( 9, 6, 4, 6, 9, 3, 6, 0, 9, 1, 3, 3, 11, 6, 3, 5 ), ( 8, 11, 1, 1, 8, 1, 0, 6, 11, 6, 1, 1, 9, 1, 4, 5, 6, 4, 1, 1 ), ( 3, 11, 6, 3, 3, 6, 0, 4, 0, 6, 4, 6 ), ( 6, 4, 8, 3, 11, 6, 8, 5 ), ( 7, 10, 6, 6, 7, 8, 10, 1, 8, 9, 10, 3 ), ( 0, 7, 3, 6, 0, 10, 7, 0, 0, 9, 10, 3, 6, 7, 10, 5 ), ( 10, 6, 7, 3, 1, 10, 7, 1, 1, 7, 8, 2, 1, 8, 0, 6 ), ( 10, 6, 7, 3, 10, 7, 1, 4, 1, 7, 3, 6 ), ( 1, 2, 6, 3, 1, 6, 8, 0, 1, 8, 9, 6, 8, 6, 7, 6 ), ( 2, 6, 9, 1, 2, 9, 1, 6, 6, 7, 9, 1, 0, 9, 3, 5, 7, 3, 9, 1 ), ( 7, 8, 0, 3, 7, 0, 6, 4, 6, 0, 2, 6 ), ( 7, 3, 2, 3, 6, 7, 2, 5 ), ( 2, 3, 11, 7, 10, 6, 8, 1, 10, 8, 9, 6, 8, 6, 7, 6 ), ( 2, 0, 7, 1, 2, 7, 11, 6, 0, 9, 7, 1, 6, 7, 10, 5, 9, 10, 7, 1 ), ( 1, 8, 0, 6, 1, 7, 8, 2, 1, 10, 7, 1, 6, 7, 10, 5, 2, 3, 11, 7 ), ( 11, 2, 1, 3, 11, 1, 7, 4, 10, 6, 1, 5, 6, 7, 1, 1 ), ( 8, 9, 6, 1, 8, 6, 7, 6, 9, 1, 6, 1, 11, 6, 3, 5, 1, 3, 6, 1 ), ( 0, 9, 1, 7, 11, 6, 7, 7 ), ( 7, 8, 0, 3, 7, 0, 6, 4, 3, 11, 0, 5, 11, 6, 0, 1 ), ( 7, 11, 6, 7 ), ( 7, 6, 11, 7 ), ( 3, 0, 8, 7, 11, 7, 6, 7 ), ( 0, 1, 9, 7, 11, 7, 6, 7 ), ( 8, 1, 9, 6, 8, 3, 1, 3, 11, 7, 6, 7 ), ( 10, 1, 2, 7, 6, 11, 7, 7 ), ( 1, 2, 10, 7, 3, 0, 8, 7, 6, 11, 7, 7 ), ( 2, 9, 0, 6, 2, 10, 9, 3, 6, 11, 7, 7 ), ( 6, 11, 7, 7, 2, 10, 3, 5, 10, 8, 3, 2, 10, 9, 8, 3 ), ( 7, 2, 3, 6, 6, 2, 7, 5 ), ( 7, 0, 8, 6, 7, 6, 0, 1, 6, 2, 0, 3 ), ( 2, 7, 6, 6, 2, 3, 7, 3, 0, 1, 9, 7 ), ( 1, 6, 2, 6, 1, 8, 6, 0, 1, 9, 8, 3, 8, 7, 6, 3 ), ( 10, 7, 6, 6, 10, 1, 7, 1, 1, 3, 7, 3 ), ( 10, 7, 6, 6, 1, 7, 10, 4, 1, 8, 7, 2, 1, 0, 8, 3 ), ( 0, 3, 7, 3, 0, 7, 10, 0, 0, 10, 9, 6, 6, 10, 7, 5 ), ( 7, 6, 10, 3, 7, 10, 8, 4, 8, 10, 9, 6 ), ( 6, 8, 4, 6, 11, 8, 6, 5 ), ( 3, 6, 11, 6, 3, 0, 6, 1, 0, 4, 6, 3 ), ( 8, 6, 11, 6, 8, 4, 6, 3, 9, 0, 1, 7 ), ( 9, 4, 6, 3, 9, 6, 3, 0, 9, 3, 1, 6, 11, 3, 6, 5 ), ( 6, 8, 4, 6, 6, 11, 8, 3, 2, 10, 1, 7 ), ( 1, 2, 10, 7, 3, 0, 11, 5, 0, 6, 11, 2, 0, 4, 6, 3 ), ( 4, 11, 8, 6, 4, 6, 11, 3, 0, 2, 9, 5, 2, 10, 9, 3 ), ( 10, 9, 3, 1, 10, 3, 2, 6, 9, 4, 3, 1, 11, 3, 6, 5, 4, 6, 3, 1 ), ( 8, 2, 3, 6, 8, 4, 2, 1, 4, 6, 2, 3 ), ( 0, 4, 2, 5, 4, 6, 2, 3 ), ( 1, 9, 0, 7, 2, 3, 4, 1, 2, 4, 6, 6, 4, 3, 8, 6 ), ( 1, 9, 4, 3, 1, 4, 2, 4, 2, 4, 6, 6 ), ( 8, 1, 3, 6, 8, 6, 1, 0, 8, 4, 6, 3, 6, 10, 1, 3 ), ( 10, 1, 0, 3, 10, 0, 6, 4, 6, 0, 4, 6 ), ( 4, 6, 3, 1, 4, 3, 8, 6, 6, 10, 3, 1, 0, 3, 9, 5, 10, 9, 3, 1 ), ( 10, 9, 4, 3, 6, 10, 4, 5 ), ( 4, 9, 5, 7, 7, 6, 11, 7 ), ( 0, 8, 3, 7, 4, 9, 5, 7, 11, 7, 6, 7 ), ( 5, 0, 1, 6, 5, 4, 0, 3, 7, 6, 11, 7 ), ( 11, 7, 6, 7, 8, 3, 4, 5, 3, 5, 4, 2, 3, 1, 5, 3 ), ( 9, 5, 4, 7, 10, 1, 2, 7, 7, 6, 11, 7 ), ( 6, 11, 7, 7, 1, 2, 10, 7, 0, 8, 3, 7, 4, 9, 5, 7 ), ( 7, 6, 11, 7, 5, 4, 10, 5, 4, 2, 10, 2, 4, 0, 2, 3 ), ( 3, 4, 8, 6, 3, 5, 4, 2, 3, 2, 5, 1, 10, 5, 2, 5, 11, 7, 6, 7 ), ( 7, 2, 3, 6, 7, 6, 2, 3, 5, 4, 9, 7 ), ( 9, 5, 4, 7, 0, 8, 6, 1, 0, 6, 2, 6, 6, 8, 7, 6 ), ( 3, 6, 2, 6, 3, 7, 6, 3, 1, 5, 0, 5, 5, 4, 0, 3 ), ( 6, 2, 8, 1, 6, 8, 7, 6, 2, 1, 8, 1, 4, 8, 5, 5, 1, 5, 8, 1 ), ( 9, 5, 4, 7, 10, 1, 6, 5, 1, 7, 6, 2, 1, 3, 7, 3 ), ( 1, 6, 10, 6, 1, 7, 6, 2, 1, 0, 7, 1, 8, 7, 0, 5, 9, 5, 4, 7 ), ( 4, 0, 10, 1, 4, 10, 5, 6, 0, 3, 10, 1, 6, 10, 7, 5, 3, 7, 10, 1 ), ( 7, 6, 10, 3, 7, 10, 8, 4, 5, 4, 10, 5, 4, 8, 10, 1 ), ( 6, 9, 5, 6, 6, 11, 9, 1, 11, 8, 9, 3 ), ( 3, 6, 11, 6, 0, 6, 3, 4, 0, 5, 6, 2, 0, 9, 5, 3 ), ( 0, 11, 8, 6, 0, 5, 11, 0, 0, 1, 5, 3, 5, 6, 11, 3 ), ( 6, 11, 3, 3, 6, 3, 5, 4, 5, 3, 1, 6 ), ( 1, 2, 10, 7, 9, 5, 11, 1, 9, 11, 8, 6, 11, 5, 6, 6 ), ( 0, 11, 3, 6, 0, 6, 11, 2, 0, 9, 6, 1, 5, 6, 9, 5, 1, 2, 10, 7 ), ( 11, 8, 5, 1, 11, 5, 6, 6, 8, 0, 5, 1, 10, 5, 2, 5, 0, 2, 5, 1 ), ( 6, 11, 3, 3, 6, 3, 5, 4, 2, 10, 3, 5, 10, 5, 3, 1 ), ( 5, 8, 9, 6, 5, 2, 8, 0, 5, 6, 2, 3, 3, 8, 2, 5 ), ( 9, 5, 6, 3, 9, 6, 0, 4, 0, 6, 2, 6 ), ( 1, 5, 8, 1, 1, 8, 0, 6, 5, 6, 8, 1, 3, 8, 2, 5, 6, 2, 8, 1 ), ( 1, 5, 6, 3, 2, 1, 6, 5 ), ( 1, 3, 6, 1, 1, 6, 10, 6, 3, 8, 6, 1, 5, 6, 9, 5, 8, 9, 6, 1 ), ( 10, 1, 0, 3, 10, 0, 6, 4, 9, 5, 0, 5, 5, 6, 0, 1 ), ( 0, 3, 8, 7, 5, 6, 10, 7 ), ( 10, 5, 6, 7 ), ( 11, 5, 10, 6, 7, 5, 11, 5 ), ( 11, 5, 10, 6, 11, 7, 5, 3, 8, 3, 0, 7 ), ( 5, 11, 7, 6, 5, 10, 11, 3, 1, 9, 0, 7 ), ( 10, 7, 5, 6, 10, 11, 7, 3, 9, 8, 1, 5, 8, 3, 1, 3 ), ( 11, 1, 2, 6, 11, 7, 1, 1, 7, 5, 1, 3 ), ( 0, 8, 3, 7, 1, 2, 7, 1, 1, 7, 5, 6, 7, 2, 11, 6 ), ( 9, 7, 5, 6, 9, 2, 7, 0, 9, 0, 2, 3, 2, 11, 7, 3 ), ( 7, 5, 2, 1, 7, 2, 11, 6, 5, 9, 2, 1, 3, 2, 8, 5, 9, 8, 2, 1 ), ( 2, 5, 10, 6, 2, 3, 5, 1, 3, 7, 5, 3 ), ( 8, 2, 0, 6, 8, 5, 2, 0, 8, 7, 5, 3, 10, 2, 5, 5 ), ( 9, 0, 1, 7, 5, 10, 3, 1, 5, 3, 7, 6, 3, 10, 2, 6 ), ( 9, 8, 2, 1, 9, 2, 1, 6, 8, 7, 2, 1, 10, 2, 5, 5, 7, 5, 2, 1 ), ( 1, 3, 5, 5, 3, 7, 5, 3 ), ( 0, 8, 7, 3, 0, 7, 1, 4, 1, 7, 5, 6 ), ( 9, 0, 3, 3, 9, 3, 5, 4, 5, 3, 7, 6 ), ( 9, 8, 7, 3, 5, 9, 7, 5 ), ( 5, 8, 4, 6, 5, 10, 8, 1, 10, 11, 8, 3 ), ( 5, 0, 4, 6, 5, 11, 0, 0, 5, 10, 11, 3, 11, 3, 0, 3 ), ( 0, 1, 9, 7, 8, 4, 10, 1, 8, 10, 11, 6, 10, 4, 5, 6 ), ( 10, 11, 4, 1, 10, 4, 5, 6, 11, 3, 4, 1, 9, 4, 1, 5, 3, 1, 4, 1 ), ( 2, 5, 1, 6, 2, 8, 5, 0, 2, 11, 8, 3, 4, 5, 8, 5 ), ( 0, 4, 11, 1, 0, 11, 3, 6, 4, 5, 11, 1, 2, 11, 1, 5, 5, 1, 11, 1 ), ( 0, 2, 5, 1, 0, 5, 9, 6, 2, 11, 5, 1, 4, 5, 8, 5, 11, 8, 5, 1 ), ( 9, 4, 5, 7, 2, 11, 3, 7 ), ( 2, 5, 10, 6, 3, 5, 2, 4, 3, 4, 5, 2, 3, 8, 4, 3 ), ( 5, 10, 2, 3, 5, 2, 4, 4, 4, 2, 0, 6 ), ( 3, 10, 2, 6, 3, 5, 10, 2, 3, 8, 5, 1, 4, 5, 8, 5, 0, 1, 9, 7 ), ( 5, 10, 2, 3, 5, 2, 4, 4, 1, 9, 2, 5, 9, 4, 2, 1 ), ( 8, 4, 5, 3, 8, 5, 3, 4, 3, 5, 1, 6 ), ( 0, 4, 5, 3, 1, 0, 5, 5 ), ( 8, 4, 5, 3, 8, 5, 3, 4, 9, 0, 5, 5, 0, 3, 5, 1 ), ( 9, 4, 5, 7 ), ( 4, 11, 7, 6, 4, 9, 11, 1, 9, 10, 11, 3 ), ( 0, 8, 3, 7, 4, 9, 7, 5, 9, 11, 7, 2, 9, 10, 11, 3 ), ( 1, 10, 11, 3, 1, 11, 4, 0, 1, 4, 0, 6, 7, 4, 11, 5 ), ( 3, 1, 4, 1, 3, 4, 8, 6, 1, 10, 4, 1, 7, 4, 11, 5, 10, 11, 4, 1 ), ( 4, 11, 7, 6, 9, 11, 4, 4, 9, 2, 11, 2, 9, 1, 2, 3 ), ( 9, 7, 4, 6, 9, 11, 7, 2, 9, 1, 11, 1, 2, 11, 1, 5, 0, 8, 3, 7 ), ( 11, 7, 4, 3, 11, 4, 2, 4, 2, 4, 0, 6 ), ( 11, 7, 4, 3, 11, 4, 2, 4, 8, 3, 4, 5, 3, 2, 4, 1 ), ( 2, 9, 10, 6, 2, 7, 9, 0, 2, 3, 7, 3, 7, 4, 9, 3 ), ( 9, 10, 7, 1, 9, 7, 4, 6, 10, 2, 7, 1, 8, 7, 0, 5, 2, 0, 7, 1 ), ( 3, 7, 10, 1, 3, 10, 2, 6, 7, 4, 10, 1, 1, 10, 0, 5, 4, 0, 10, 1 ), ( 1, 10, 2, 7, 8, 7, 4, 7 ), ( 4, 9, 1, 3, 4, 1, 7, 4, 7, 1, 3, 6 ), ( 4, 9, 1, 3, 4, 1, 7, 4, 0, 8, 1, 5, 8, 7, 1, 1 ), ( 4, 0, 3, 3, 7, 4, 3, 5 ), ( 4, 8, 7, 7 ), ( 9, 10, 8, 5, 10, 11, 8, 3 ), ( 3, 0, 9, 3, 3, 9, 11, 4, 11, 9, 10, 6 ), ( 0, 1, 10, 3, 0, 10, 8, 4, 8, 10, 11, 6 ), ( 3, 1, 10, 3, 11, 3, 10, 5 ), ( 1, 2, 11, 3, 1, 11, 9, 4, 9, 11, 8, 6 ), ( 3, 0, 9, 3, 3, 9, 11, 4, 1, 2, 9, 5, 2, 11, 9, 1 ), ( 0, 2, 11, 3, 8, 0, 11, 5 ), ( 3, 2, 11, 7 ), ( 2, 3, 8, 3, 2, 8, 10, 4, 10, 8, 9, 6 ), ( 9, 10, 2, 3, 0, 9, 2, 5 ), ( 2, 3, 8, 3, 2, 8, 10, 4, 0, 1, 8, 5, 1, 10, 8, 1 ), ( 1, 10, 2, 7 ), ( 1, 3, 8, 3, 9, 1, 8, 5 ), ( 0, 9, 1, 7 ), ( 0, 3, 8, 7 ), None ) |