porepy.geometry.half_space module

This module contains functions for computations relating to half spaces.

half_space_interior_point(n, x0, pts, recompute=True)

Find an interior point for the halfspaces.

Note

We use linear programming to find one interior point for the half spaces. Assume, num_planes halfspaces defined by

aj*x1+bj*x2+cj*x3+dj<=0, j=1..num_planes.

Perform the following linear program:

max(x5) aj*x1+bj*x2+cj*x3+dj*x4+x5<=0, j=1..num_planes

Then, if [x1,x2,x3,x4,x5] is an optimal solution with x4>0 and x5>0 we get:

aj*(x1/x4)+bj*(x2/x4)+cj*(x3/x4)+dj<=(-x5/x4) j=1..num_planes

and

(-x5/x4)<0,

and conclude that the point [x1/x4,x2/x4,x3/x4] is in the interior of all the halfspaces. Since x5 is optimal, this point is “way in” the interior (good for precision errors). For more information, see http://www.qhull.org/html/qhalf.htm#notes

Parameters

• n (ndarray) –

  shape=(3, num_planes)

  This is the normal vectors of the half planes. The normal vectors are assumed to be coherently oriented for all the half spaces (inward or outward).

• x0 (ndarray) –

  shape=(3, num_planes)

  Point on the boundary of the half-spaces. Half space i is given by all points satisfying (x - x0[:, i]) * n[:, i] <= 0.

• pts (ndarray) –

  shape=(3, np)

  Points used to bound the search space for interior point. The optimum solution will be sought within (pts.min(axis=1), pts.max(axis=1)).

• recompute (bool) –

  default=True

  If the algorithm fails, try again with flipped normals.

Raises

ValueError – If the inequalities define an empty half space.

Returns

Interior point of the halfspaces with shape=(np, ).

Return type

ndarray

point_inside_half_space_intersection(n, x0, pts)

Find the points that lie in the intersection of half spaces (in 3D).

Examples

>>> import numpy as np
>>> n = np.array([[0, 1], [1, 0], [0, 0]])
>>> x0 = np.array([[0, -1], [0, 0], [0, 0]])
>>> pts = np.array([[-1 ,-1 ,4], [2, -2, -2], [0, 0, 0]])
>>> half_space_int(n,x0,pts)
array([False,  True, False], dtype=bool)

Parameters

• n (ndarray) –

  shape=(3, num_planes)

  The normal vectors of the half planes. The normal vectors are assumed to point out of the half spaces.

• x0 (ndarray) –

  shape=(3, num_planes)

  Point on the boundary of the half-spaces. Half space i is given by all points satisfying (x - x0[:,i])*n[:,i]<=0.

• pts (ndarray) –

  shape=(3, np)

  The points to be tested if they are in the intersection of all half-spaces or not.

Raises

• ValueError – If either of the parameters n, x0 or pts are not three dimensional.

• ValueError – If the number of columns in n and x0 are not equal.

Returns

A logical array with shape=(np, ). out[i] is True if pts[:, i] is in all half-spaces.

Return type

ndarray

vertexes_of_convex_domain(A, b)

Find the vertexes of a convex domain specified as an intersection of half spaces.

Note

The function assumes the domain is defined by inequalities on the form

A * x + b <= 0

For more information, see scipy.spatial functions HalfspaceIntersection.

The function has been tested for 2d and 3d domains.

Parameters

• A (ndarray) –

  shape=(num_planes, nd)

  Matrix of normal vectors (in rows) for the half planes. Should be oriented so that A * x + b < 0

• b (ndarray) –

  shape=(num_planes,)

  Constants used to define inequalities of the half spaces. Should be scaled so that A * x + b < 0.

Raises

QhullError – QhullError: If A and b are not set up right (e.g. sign errors that imply that the inequalities do not form a closed domain).

Returns

Vertexes of a convex domain.

Return type

ndarray