PyCheeger

pycheeger.compute_cheeger.compute_cheeger(eta, max_tri_area_fm=0.002, max_iter_fm=10000, plot_results_fm=False, num_boundary_vertices_ld=50, max_tri_area_ld=0.005, step_size_ld=0.01, max_iter_ld=500, convergence_tol_ld=0.0001, num_iter_resampling_ld=None, plot_results_ld=False)

Compute the Cheeger set associated to the weight function eta

Parameters

• eta (function) – Function to be integrated. f must handle array inputs with shape (N, 2)

• max_tri_area_fm (float) – Fixed mesh step parameter. Maximum triangle area allowed for the domain mesh

• max_iter_fm (int) – Fixed mesh step parameter. Maximum number of iterations for the primal dual algorithm

• plot_results_fm (bool) – Fixed mesh step parameter. Whether to plot the results of the fixed mesh step or not

• num_boundary_vertices_ld (int) – Local descent step parameter. Number of boundary vertices used to represent the simple set

• max_tri_area_ld (float) – Local descent step parameter. Maximum triangle area allowed for the inner mesh of the simple set

• step_size_ld (float) – Local descent step parameter. Step size used in the local descent

• max_iter_ld (int) – Local descent step parameter. Maximum number of iterations allowed for the local descent

• convergence_tol_ld (float) – Local descent step parameter. Convergence tol for the local descent

• num_iter_resampling_ld (None or int) – Local descent step parameter. Number of iterations between two resampling of the boundary curve (None for no resampling)

• plot_results_ld (bool) – Local descent step parameter. Whether to plot the results of the local descent step or not

Returns

• simplet_set (SimpleSet) – Cheeger set

• obj_tab (array, shape (n_iter_ld,)) – Values of the objective over the course of the local descent

• grad_norm_tab (array, shape (n_iter_ld,)) – Values of the objective gradient norm over the course of the local descent

class pycheeger.simple_set.SimpleSet(boundary_vertices, max_tri_area=None)

Simple set

boundary_vertices

Each row contains the two coordinates of a boundary vertex

Type

array, shape (N, 2)

num_boundary_vertices

Number of boundary vertices

Type

int

is_clockwise

Whether the list of boundary vertices is clockwise ordered or not

Type

bool

mesh_vertices

Each row contains the two coordinates of a vertex of the simplet set inner mesh

Type

array, shape (M, 2)

mesh_faces

Each row contains the three indices describing a face of the simple set inner mesh

Type

array, shape (K, 3)

mesh_boundary_faces_indices

One dimensional array, contains the indices of the boundary faces of the simple set inner mesh

Type

array

compute_perimeter()

Compute the perimeter of the set

Returns

The perimeter

Return type

float

compute_perimeter_gradient()

Compute the “gradient” of the perimeter

Returns

Each row contains the two coordinates of the translation to apply at each boundary vertex

Return type

array, shape (N, 2)

Notes

See [1]_ (first variation of the perimeter)

References

1

Maggi, F. (2012). Sets of finite perimeter and geometric variational problems: an introduction to Geometric Measure Theory (No. 135). Cambridge University Press.

compute_weighted_area(f)

Compute the integral of f over the set

Parameters

f (function) – Function to be integrated. f must handle array inputs with shape (N, 2). It can be vector valued

Returns

Value computed for the integral of f over the set (if f takes values in dimension D, the result will be an array of shape (D,))

Return type

float or array of shape (D,)

compute_weighted_area_gradient(f)

Compute the “gradient” of the weighted area, for a given weight function

Parameters

f (function) – Function to be integrated. f must handle array inputs with shape (N, 2). It can be vector valued

Returns

Return type

array, shape

Notes

Vectorized computations are really nasty here, mainly because f can be vector valued.

compute_weighted_area_tab(f, boundary_faces_only=False)

Compute the integral of f on each face of the inner mesh

Parameters

• f (function) – Function to be integrated. f must handle array inputs with shape (N, 2). It can be vector valued

• boundary_faces_only (bool) – Whether to compute weighted areas only on boundary faces, defaut False

Returns

Value computed for the integral of f on each of the N triangles (if f takes values in dimension D, the shape of the resulting array is (N, D))

Return type

array, shape (N,) or (N,D)

contains(x)

Whether a given point x is inside the set or not

Parameters

x (array, shape (2,)) – The input point

Returns

Whether x is in the set or not

Return type

bool

create_mesh(boundary_vertices, max_tri_area)

Create the inner mesh of the set

Parameters

• boundary_vertices (array, shape (N, 2)) – Each row contains the two coordinates of a boundary vertex

• max_tri_area (float) – Maximum triangle area for the inner mesh

pycheeger.simple_set.disk(center, radius, num_vertices=20, max_tri_area=None)

Create a SimpleSet with boundary vertices regularly spaced on a circle

Parameters

• center (array, shape (2,)) – Coordinates of the center

• radius (float) – Radius of the circles

• num_vertices (int) – Number of boundary vertices (on the disk)

• max_tri_area (None or float) – Maximum area allowed for triangles, see Shewchuk’s Triangle mesh generator, defaut None (no constraint)

Returns

Return type

SimpleSet

class pycheeger.mesh.Mesh(raw_mesh)

Custom mesh class

vertices

Each row contains the two coordinates of a vertex

Type

array, shape (N, 2)

faces

Each row contains the three indices of the face’s vertices. For convenience reasons, the values in each row are sorted, i.e. the face whose vertices have indices 0, 1, 2 will always be stored as [0, 1, 2] (and not [1, 0, 2])

Type

array, shape (M, 3)

edges

Each row contains the two indices of the edge’s vertices. For convenience reasons, the values in each row are sorted, i.e. the edge whose vertices have indices 0, 1 will always be stored as [0, 1] (and not [1, 0])

Type

array, shape (K, 2)

num_edges

Equals K

Type

int

num_faces

Equals M

Type

int

num_vertices

Equals N

Type

int

grad_mat

Linear operator mapping a vector of length M (the values taken by a piecewise constant function on the mesh) to a vector of length K (the jumps / values taken by the gradient on each edge of the mesh)

Type

scipy.sparse.csr_matrix

build_grad_matrix()

Compute the gradient matrix associated to the mesh

The gradient is the linear operator mapping a vector of length M (the values taken by a piecewise constant function on the mesh) to a vector of length K (the jumps / values taken by the gradient on each edge of the mesh), M being the number of faces in the mesh, and K being the number of edges.

Returns

The gradient matrix

Return type

scipy.sparse.csr_matrix

get_edge_adjacent_faces(edge)

Find the index of all faces to which a given edge belongs

Parameters

edge (array, shape (2,)) – Array containing two integers, which are the indices of the edge’s vertices. Does not need to be sorted.

Returns

One dimensional integer valued array containing the indices of the relevant faces

Return type

array

get_edge_index(edge)

Find the index of an edge given by its two vertices in the array of edges

Parameters

edge (array, shape (2,)) – Array containing two integers, which are the indices of the edge’s vertices. Does not need to be sorted.

Returns

The index of the input edge (raises an error if the edge is not found or found multiple times)

Return type

int

get_edge_length(edge)

Compute the euclidean distance between the two vertices of an edge

Parameters

edge (array, shape(2,)) – Array containing two integers, which are the indices of the edge’s vertices. Does not need to be sorted.

Returns

The length of the edge

Return type

float

get_face_area(face_index)

Compute the area of a face

Parameters

face_index (int) – Index of the face in self.faces

Returns

Area of the input face

Return type

float

get_face_edges(face_index)

Get the list of edges which belong to a given face

Parameters

face_index (int) – The index of the face in self.faces

Returns

Each row contains the two indices of an edge’s vertices (sorted).

Return type

array, shape (3, 2)

get_orientation(face_index, edge)

Compute the orientation of an edge with respect to a given face.

The orientation of an edge with respect to a face is defined as +1 if the normal is the outward normal and -1 otherwise. The normal is computed by applying a counter clockwise rotation to the edge vector, which is himself given by v2 - v1 where v1 and v2 are the two vertices of the edge, and the index of v1 in self.vertices is smaller than the one of v2 (edges are always stored sorted). This definition of the orientation is of course arbitrary, but the only thing we need is to keep one that is consistent all along.

Parameters

• face_index (int) – Index of the face in self.faces

• edge (array, shape (2,)) – Array containing two integers, which are the indices of the edge’s vertices. Does not need to be sorted.

Returns

+1 or -1 if the input face contains the input edge (see the explanations above), and 0 otherwise

Return type

int

integrate(f)

Compute the integral of a given function on each face of the mesh

Parameters

f (function) – Function to be integrated. f must handle array inputs with shape (N, 2), N being the number of faces of the mesh. It can be vector valued

Returns

Value computed for the integral of f on each of the N triangles (if f takes values in dimension D, the shape of the resulting array is (N, D))

Return type

array, shape (N,) or (N, D)

pycheeger.tools.find_threshold(x)

Compute the value of the Lagrange multiplier involved in the projection of x into the unit l1 ball

Parameters

x (array, shape (N,)) – Vector to be projected

Returns

res – Value of the Lagrange multiplier

Return type

float

Notes

Sorting based algorithm. See [1]_ for a detailed explanation of the computations and alternative algorithms.

References

1

L. Condat, Fast Projection onto the Simplex and the l1 Ball, Mathematical Programming, Series A, Springer, 2016, 158 (1), pp.575-585.

pycheeger.tools.integrate_on_triangles(f, triangles)

Numerical integration of f on a list of triangles

Parameters

• f (function) – Function to be integrated. f must handle array inputs with shape (N, 2). It can be vector valued

• triangles (array, shape (N, 3, 2)) – triangles[i, j] contains the coordinates of the j-th vertex of the i-th triangle

Returns

Value computed for the integral of f on each of the N triangles (if f takes values in dimension D, the shape of the resulting array is (N, D))

Return type

array, shape (N,) or (N, D)

Notes

Here, quadpy is only used to extract the scheme’s characteristics, in order to speed up computations (in quadpy, the code handles arbitrary dimensions and is too generic)

pycheeger.tools.postprocess_indicator(x, grad_mat)

Post process a piecewise constant function on a mesh to get an indicator function of a union of cells

Parameters

• x (array, shape (N,)) – Values describing the piecewise constant function to be processed

• grad_mat (array, shape (M, N)) – Matrix representing the linear operator which maps the values describing a piecewise constant function to its jumps on each edge of the mesh

Returns

Values of the indicator function on each cell of the mesh

Return type

array, shape (N,)

pycheeger.tools.proj_unit_l1_ball(x)

Projection onto the l1 unit ball

Parameters

x (array, shape (N,)) – Vector to be projected

Returns

res – Projection

Return type

array, shape (N,)

Notes

See [1]_ for a detailed explanation of the computations and alternative algorithms.

References

1

L. Condat, Fast Projection onto the Simplex and the l1 Ball, Mathematical Programming, Series A, Springer, 2016, 158 (1), pp.575-585.

pycheeger.tools.prox_dot_prod(x, tau, a)

Proximal map of the inner product between x and a

Parameters

• x (array, shape (N,)) –

• tau (float) –

• a (array, shape (N,)) –

Returns

Return type

array, shape (N,)

pycheeger.tools.prox_inf_norm(x, tau)

Proximal map of the l-infinity norm

Parameters

• x (array, shape (N,)) –

• tau (float) –

Returns

Return type

array, shape (N,)

Notes

\[prox_{\tau \, ||.||_{\infty}}(x) = x - \tau ~ \text{proj}_{\{||.||_{\infty}\leq 1\}}(x / \tau)\]

pycheeger.tools.resample(curve, num_points)

Resample a closed polygonal chain

Parameters

• curve (array, shape (N, 2)) – Curve to be resampled, described by the list of its vertices. The last vertex should not be equal to the first one (the input array should be a minimal description of the closed polygonal chain)

• num_points (int) – Number of vertices of the output curve

Returns

Resampled curve

Return type

array, shape (num_points, 2)

pycheeger.tools.run_primal_dual(mesh, eta_bar, max_iter, grad_mat_norm, verbose=False)

Solves the “fixed mesh weighted Cheeger problem” by running a primal dual algorithm

Parameters

• mesh (Mesh) – Triangle mesh made of N triangles and M edges

• eta_bar (array, shape (N, 2)) – Integral of the weight function on each triangle

• max_iter (integer) – Maximum number of iterations (for now, exact number of iterations, since no convergence criterion is implemented yet)

• grad_mat_norm (float) – Norm of the gradient operator for piecewise constant functions on the mesh

• verbose (bool, defaut False) – Whether to print some information at the end of the algorithm or not

Returns

Values describing a piecewise constant function on the mesh, which solves the fixed mesh weighted Cheeger problem

Return type

array, shape (N, 2)

pycheeger.tools.triangulate(vertices, max_triangle_area=None, split_boundary=False, plot_result=False)

Triangulate the interior of a closed polygonal curve using the triangle library (Python wrapper around Shewchuk’s Triangle mesh generator)

Parameters

• vertices (array, shape (N, 2)) – Coordinates of the curve’s vertices

• max_triangle_area (None or float) – Maximum area allowed for triangles, see Shewchuk’s Triangle mesh generator, defaut None (no constraint)

• split_boundary (bool) – Whether to allow boundary segments to be splitted or not, defaut False

• plot_result (bool) – If True, the resulting triangulation is shown along with the input, defaut False

Returns

raw_mesh – Output mesh, see the documentation of the triangle library

Return type

dict

pycheeger.tools.winding(x, vertices)

Compute the winding number of a closed polygonal curve described by its vertices around the point x. This number is zero if and only if x is outside the polygon.

Parameters

• x (array, shape (2,)) – Point around which the winding number is computed

• vertices (array, shape (N, 2)) – Coordinates of the curve’s vertices

Returns

wn – The winding number

Return type

float