sfepy.discrete.fem.mesh module

class sfepy.discrete.fem.mesh.Mesh(name='mesh', filename=None, prefix_dir=None, **kwargs)[source]

Contains the FEM mesh together with all utilities related to it.

Input and output is handled by the MeshIO class and subclasses. The Mesh class only contains the real mesh - nodes, connectivity, regions, plus methods for doing operations on this mesh.

Example of creating and working with a mesh:

In [1]: from sfepy.discrete.fem import Mesh
In [2]: m = Mesh.from_file("meshes/3d/cylinder.vtk")
sfepy: reading mesh (meshes/3d/cylinder.vtk)...
sfepy: ...done in 0.04 s

In [3]: m.coors
Out[3]:
array([[  1.00000000e-01,   2.00000000e-02,  -1.22460635e-18],
       [  1.00000000e-01,   1.80193774e-02,   8.67767478e-03],
       [  1.00000000e-01,   1.24697960e-02,   1.56366296e-02],
       ...,
       [  8.00298527e-02,   5.21598617e-03,  -9.77772215e-05],
       [  7.02544004e-02,   3.61610291e-04,  -1.16903153e-04],
       [  3.19633596e-02,  -1.00335972e-02,   9.60460305e-03]])

In [4]: m.ngroups
Out[4]: array([0, 0, 0, ..., 0, 0, 0])

In [5]: m.conns
Out[5]:
[array([[ 28,  60,  45,  29],
       [ 28,  60,  57,  45],
       [ 28,  57,  27,  45],
       ...,
       [353, 343, 260, 296],
       [353, 139, 181, 140],
       [353, 295, 139, 140]])]

In [6]: m.mat_ids
Out[6]: [array([6, 6, 6, ..., 6, 6, 6])]

In [7]: m.descs
Out[7]: ['3_4']

In [8]: m
Out[8]: Mesh:meshes/3d/cylinder

In [9]: print m
Mesh:meshes/3d/cylinder
  conns:
    [array([[ 28,  60,  45,  29],
           [ 28,  60,  57,  45],
           [ 28,  57,  27,  45],
           ...,
           [353, 343, 260, 296],
           [353, 139, 181, 140],
           [353, 295, 139, 140]])]
  coors:
    [[  1.00000000e-01   2.00000000e-02  -1.22460635e-18]
     [  1.00000000e-01   1.80193774e-02   8.67767478e-03]
     [  1.00000000e-01   1.24697960e-02   1.56366296e-02]
     ...,
     [  8.00298527e-02   5.21598617e-03  -9.77772215e-05]
     [  7.02544004e-02   3.61610291e-04  -1.16903153e-04]
     [  3.19633596e-02  -1.00335972e-02   9.60460305e-03]]
  descs:
    ['3_4']
  dim:
    3
  el_offsets:
    [   0 1348]
  io:
    None
  mat_ids:
    [array([6, 6, 6, ..., 6, 6, 6])]
  n_e_ps:
    [4]
  n_el:
    1348
  n_els:
    [1348]
  n_nod:
    354
  name:
    meshes/3d/cylinder
  ngroups:
    [0 0 0 ..., 0 0 0]
  setup_done:
    0

The Mesh().coors is an array of node coordinates and Mesh().conns is the list of elements of each type (see Mesh().desc), so for example if you want to know the coordinates of the nodes of the fifth finite element of the type 3_4 do:

In [10]: m.descs
Out[10]: ['3_4']

So now you know that the finite elements of the type 3_4 are in a.conns[0]:

In [11]: m.coors[m.conns[0][4]]
Out[11]:
array([[  1.00000000e-01,   1.80193774e-02,  -8.67767478e-03],
       [  1.00000000e-01,   1.32888539e-02,  -4.35893200e-04],
       [  1.00000000e-01,   2.00000000e-02,  -1.22460635e-18],
       [  9.22857574e-02,   1.95180454e-02,  -4.36416134e-03]])

The element ids are of the form “<dimension>_<number of nodes>”, i.e.:

  • 2_2 ... line
  • 2_3 ... triangle
  • 2_4 ... quadrangle
  • 3_2 ... line
  • 3_4 ... tetrahedron
  • 3_8 ... hexahedron
copy(name=None)[source]

Make a deep copy of self.

Parameters:

name : str

Name of the copied mesh.

create_conn_graph(verbose=True)[source]

Create a graph of mesh connectivity.

Returns:

graph : csr_matrix

The mesh connectivity graph as a SciPy CSR matrix.

explode_groups(eps, return_emap=False)[source]

Explode the mesh element groups by eps, i.e. split group interface nodes and shrink each group towards its centre by eps.

Parameters:

eps : float in [0.0, 1.0]

The group shrinking factor.

return_emap : bool, optional

If True, also return the mapping against original mesh coordinates that result in the exploded mesh coordinates. The mapping can be used to map mesh vertex data to the exploded mesh vertices.

Returns:

mesh : Mesh

The new mesh with exploded groups.

emap : spmatrix, optional

The maping for exploding vertex values. Only provided if return_emap is True.

static from_data(name, coors, ngroups, conns, mat_ids, descs, igs=None, nodal_bcs=None)[source]

Create a mesh from mesh data.

static from_file(filename=None, io='auto', prefix_dir=None, omit_facets=False)[source]

Read a mesh from a file.

Parameters:

filename : string or function or MeshIO instance or Mesh instance

The name of file to read the mesh from. For convenience, a mesh creation function or a MeshIO instance or directly a Mesh instance can be passed in place of the file name.

io : *MeshIO instance

Passing *MeshIO instance has precedence over filename.

prefix_dir : str

If not None, the filename is relative to that directory.

omit_facets : bool

If True, do not read cells of lower dimension than the space dimension (faces and/or edges). Only some MeshIO subclasses support this!

static from_region(region, mesh_in, save_edges=False, save_faces=False, localize=False, is_surface=False)[source]

Create a mesh corresponding to a given region.

static from_surface(surf_faces, mesh_in)[source]

Create a mesh given a set of surface faces and the original mesh.

get_bounding_box()[source]
get_element_coors(ig=None)[source]

Get the coordinates of vertices elements in group ig.

Parameters:

ig : int, optional

The element group. If None, the coordinates for all groups are returned, filled with zeros at places of missing vertices, i.e. where elements having less then the full number of vertices (n_ep_max) are.

Returns:

coors : array

The coordinates in an array of shape (n_el, n_ep_max, dim).

localize(inod)[source]

Strips nodes not in inod and remaps connectivities. Omits elements where remap[conn] contains -1...

transform_coors(mtx_t, ref_coors=None)[source]

Transform coordinates of the mesh by the given transformation matrix.

Parameters:

mtx_t : array

The transformation matrix T (2D array). It is applied depending on its shape:

  • (dim, dim): x = T * x
  • (dim, dim + 1): x = T[:, :-1] * x + T[:, -1]

ref_coors : array, optional

Alternative coordinates to use for the transformation instead of the mesh coordinates, with the same shape as self.coors.

write(filename=None, io=None, coors=None, igs=None, out=None, float_format=None, **kwargs)[source]

Write mesh + optional results in out to a file.

Parameters:

filename : str, optional

The file name. If None, the mesh name is used instead.

io : MeshIO instance or ‘auto’, optional

Passing ‘auto’ respects the extension of filename.

coors : array, optional

The coordinates that can be used instead of the mesh coordinates.

igs : array_like, optional

Passing a list of group ids selects only those groups for writing.

out : dict, optional

The output data attached to the mesh vertices and/or cells.

float_format : str, optional

The format string used to print floats in case of a text file format.

**kwargs : dict, optional

Additional arguments that can be passed to the MeshIO instance.

sfepy.discrete.fem.mesh.find_map(x1, x2, eps=1e-08, allow_double=False, join=True)[source]

Find a mapping between common coordinates in x1 and x2, such that x1[cmap[:,0]] == x2[cmap[:,1]]

sfepy.discrete.fem.mesh.fix_double_nodes(coor, ngroups, conns, eps)[source]

Detect and attempt fixing double nodes in a mesh.

The double nodes are nodes having the same coordinates w.r.t. precision given by eps.

sfepy.discrete.fem.mesh.get_min_edge_size(coor, conns)[source]

Get the smallest edge length.

sfepy.discrete.fem.mesh.get_min_vertex_distance(coor, guess)[source]

Can miss the minimum, but is enough for our purposes.

sfepy.discrete.fem.mesh.get_min_vertex_distance_naive(coor)[source]
sfepy.discrete.fem.mesh.make_inverse_connectivity(conns, n_nod, ret_offsets=True)[source]

For each mesh node referenced in the connectivity conns, make a list of elements it belongs to.

sfepy.discrete.fem.mesh.make_mesh(coor, ngroups, conns, mesh_in)[source]

Create a mesh reusing mat_ids and descs of mesh_in.

sfepy.discrete.fem.mesh.make_point_cells(indx, dim)[source]
sfepy.discrete.fem.mesh.merge_mesh(x1, ngroups1, conns1, mat_ids1, x2, ngroups2, conns2, mat_ids2, cmap, eps=1e-08)[source]

Merge two meshes in common coordinates found in x1, x2.

Notes

Assumes the same number and kind of element groups in both meshes!