sfepy.discrete.common.global_interp module

Global interpolation functions.

sfepy.discrete.common.global_interp.get_potential_cells(coors, cmesh, centroids=None, extrapolate=True)[source]

Get cells that potentially contain points with the given physical coordinates.

Parameters:
coors : array

The physical coordinates.

cmesh : CMesh instance

The cmesh defining the cells.

centroids : array, optional

The centroids of the cells.

extrapolate : bool

If True, even the points that are surely outside of the cmesh are considered and assigned potential cells.

Returns:
potential_cells : array

The indices of the cells that potentially contain the points.

offsets : array

The offsets into potential_cells for each point: a point ip is potentially in cells potential_cells[offsets[ip]:offsets[ip+1]].

sfepy.discrete.common.global_interp.get_ref_coors(field, coors, strategy=’general’, close_limit=0.1, get_cells_fun=None, cache=None, verbose=False)[source]

Get reference element coordinates and elements corresponding to given physical coordinates.

Parameters:
field : Field instance

The field defining the approximation.

coors : array

The physical coordinates.

strategy : {‘general’, ‘convex’}, optional

The strategy for finding the elements that contain the coordinates. For convex meshes, the ‘convex’ strategy might be faster than the ‘general’ one.

close_limit : float, optional

The maximum limit distance of a point from the closest element allowed for extrapolation.

get_cells_fun : callable, optional

If given, a function with signature get_cells_fun(coors, cmesh, **kwargs) returning cells and offsets that potentially contain points with the coordinates coors. Applicable only when strategy is ‘general’. When not given, get_potential_cells() is used.

cache : Struct, optional

To speed up a sequence of evaluations, the field mesh and other data can be cached. Optionally, the cache can also contain the reference element coordinates as cache.ref_coors, cache.cells and cache.status, if the evaluation occurs in the same coordinates repeatedly. In that case the mesh related data are ignored.

verbose : bool

If False, reduce verbosity.

Returns:
ref_coors : array

The reference coordinates.

cells : array

The cell indices corresponding to the reference coordinates.

status : array

The status: 0 is success, 1 is extrapolation within close_limit, 2 is extrapolation outside close_limit, 3 is failure, 4 is failure due to non-convergence of the Newton iteration in tensor product cells. If close_limit is 0, then for the ‘general’ strategy the status 5 indicates points outside of the field domain that had no potential cells.

sfepy.discrete.common.global_interp.get_ref_coors_convex(field, coors, close_limit=0.1, cache=None, verbose=False)[source]

Get reference element coordinates and elements corresponding to given physical coordinates.

Parameters:
field : Field instance

The field defining the approximation.

coors : array

The physical coordinates.

close_limit : float, optional

The maximum limit distance of a point from the closest element allowed for extrapolation.

cache : Struct, optional

To speed up a sequence of evaluations, the field mesh and other data can be cached. Optionally, the cache can also contain the reference element coordinates as cache.ref_coors, cache.cells and cache.status, if the evaluation occurs in the same coordinates repeatedly. In that case the mesh related data are ignored.

verbose : bool

If False, reduce verbosity.

Returns:
ref_coors : array

The reference coordinates.

cells : array

The cell indices corresponding to the reference coordinates.

status : array

The status: 0 is success, 1 is extrapolation within close_limit, 2 is extrapolation outside close_limit, 3 is failure, 4 is failure due to non-convergence of the Newton iteration in tensor product cells.

Notes

Outline of the algorithm for finding xi such that X(xi) = P:

  1. make inverse connectivity - for each vertex have cells it is in.
  2. find the closest vertex V.
  3. choose initial cell: i0 = first from cells incident to V.
  4. while not P in C_i, change C_i towards P, check if P in new C_i.
sfepy.discrete.common.global_interp.get_ref_coors_general(field, coors, close_limit=0.1, get_cells_fun=None, cache=None, verbose=False)[source]

Get reference element coordinates and elements corresponding to given physical coordinates.

Parameters:
field : Field instance

The field defining the approximation.

coors : array

The physical coordinates.

close_limit : float, optional

The maximum limit distance of a point from the closest element allowed for extrapolation.

get_cells_fun : callable, optional

If given, a function with signature get_cells_fun(coors, cmesh, **kwargs) returning cells and offsets that potentially contain points with the coordinates coors. When not given, get_potential_cells() is used.

cache : Struct, optional

To speed up a sequence of evaluations, the field mesh and other data can be cached. Optionally, the cache can also contain the reference element coordinates as cache.ref_coors, cache.cells and cache.status, if the evaluation occurs in the same coordinates repeatedly. In that case the mesh related data are ignored.

verbose : bool

If False, reduce verbosity.

Returns:
ref_coors : array

The reference coordinates.

cells : array

The cell indices corresponding to the reference coordinates.

status : array

The status: 0 is success, 1 is extrapolation within close_limit, 2 is extrapolation outside close_limit, 3 is failure, 4 is failure due to non-convergence of the Newton iteration in tensor product cells. If close_limit is 0, then status 5 indicates points outside of the field domain that had no potential cells.