Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Enhance Mesh: Mesh tools as methods #425

Merged
merged 9 commits into from
Apr 5, 2023
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
2 changes: 2 additions & 0 deletions CHANGELOG.md
Original file line number Diff line number Diff line change
Expand Up @@ -10,6 +10,8 @@ All notable changes to this project will be documented in this file. The format
### Changed
- Change `einsumt` from an optional to a required dependency.
- Vectorize implementations of `MultiPointConstraint` and `MultiPointContact` and re-implement both as `scipy.sparse.lil_matrix()`.
- Rename old `Mesh` to `DiscreteGeometry` and rebase new `Mesh` on `DiscreteGeometry`.
- Simplify the usage of explicit mesh-related tools by adding them as methods to `Mesh`, i.e. `mesh.tools.expand(Rectangle())` is equivalent to `Rectangle().expand()`.

### Fixed
- Fix `tools.project()` for higher-order quad- and hexahedron elements.
Expand Down
2 changes: 1 addition & 1 deletion src/felupe/mesh/_convert.py
Original file line number Diff line number Diff line change
Expand Up @@ -27,8 +27,8 @@

import numpy as np

from ._discrete_geometry import DiscreteGeometry
from ._helpers import mesh_or_data
from ._mesh import Mesh


@mesh_or_data
Expand Down
102 changes: 102 additions & 0 deletions src/felupe/mesh/_discrete_geometry.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,102 @@
# -*- coding: utf-8 -*-
"""
_______ _______ ___ __ __ _______ _______
| || || | | | | || || |
| ___|| ___|| | | | | || _ || ___|
| |___ | |___ | | | |_| || |_| || |___
| ___|| ___|| |___ | || ___|| ___|
| | | |___ | || || | | |___
|___| |_______||_______||_______||___| |_______|

This file is part of felupe.

Felupe is free software: you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation, either version 3 of the License, or
(at your option) any later version.

Felupe 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.

You should have received a copy of the GNU General Public License
along with Felupe. If not, see <http://www.gnu.org/licenses/>.

"""

import numpy as np


class DiscreteGeometry:
"""A discrete geometry with points, cells and optional a specified cell type.

Parameters
----------
points : ndarray
Point coordinates.
cells : ndarray
Point-connectivity of cells.
cell_type : str or None, optional
An optional string in VTK-convention that specifies the cell type (default is None). Necessary when a discrete geometry is saved to a file.

Attributes
----------
points : ndarray
Point coordinates.
cells : ndarray
Point-connectivity of cells.
cell_type : str or None
A string in VTK-convention that specifies the cell type.
npoints : int
Amount of points.
dim : int
Dimension of mesh point coordinates.
ndof : int
Amount of degrees of freedom.
ncells : int
Amount of cells.
points_with_cells : array
Array with points connected to cells.
points_without_cells : array
Array with points not connected to cells.
cells_per_point : array
Array which counts connected cells per point. Used for averaging results.

"""

def __init__(self, points, cells, cell_type=None):
self.points = np.array(points)
self.cells = np.array(cells)
self.cell_type = cell_type

self.update(self.cells)

def update(self, cells, cell_type=None):
"Update the cell and dimension attributes with a given cell array."
self.cells = cells

if cell_type is not None:
self.cell_type = cell_type

# obtain dimensions
self.npoints, self.dim = self.points.shape
self.ndof = self.points.size
self.ncells = self.cells.shape[0]

# get number of cells per point
points_in_cell, self.cells_per_point = np.unique(cells, return_counts=True)

# check if there are points without cells
if self.npoints != len(self.cells_per_point):
self.point_has_cell = np.isin(np.arange(self.npoints), points_in_cell)
# update "cells_per_point" ... cells per point
cells_per_point = -np.ones(self.npoints, dtype=int)
cells_per_point[points_in_cell] = self.cells_per_point
self.cells_per_point = cells_per_point

self.points_without_cells = np.arange(self.npoints)[~self.point_has_cell]
self.points_with_cells = np.arange(self.npoints)[self.point_has_cell]
else:
self.points_without_cells = np.array([], dtype=int)
self.points_with_cells = np.arange(self.npoints)
15 changes: 8 additions & 7 deletions src/felupe/mesh/_helpers.py
Original file line number Diff line number Diff line change
Expand Up @@ -27,12 +27,10 @@

from functools import wraps

from ._mesh import Mesh


def mesh_or_data(meshfun):
"""If a ``Mesh`` is passed to a mesh function, extract ``points`` and
``cells`` arrays along with the ``cell_type`` and return a ``Mesh``
"""If a ``DiscreteGeometry`` is passed to a mesh function, extract ``points`` and
``cells`` arrays along with the ``cell_type`` and return a ``DiscreteGeometry``
as a result."""

@wraps(meshfun)
Expand All @@ -44,8 +42,8 @@ def mesh_or_points_cells_type(*args, **kwargs):
# check if unnamed args are passed
if len(args) > 0:

# meshfun(Mesh)
if isinstance(args[0], Mesh):
# meshfun(DiscreteGeometry)
if hasattr(args[0], "__mesh__"):

# set mesh flag
is_mesh = True
Expand All @@ -55,6 +53,9 @@ def mesh_or_points_cells_type(*args, **kwargs):
cells = args[0].cells
cell_type = args[0].cell_type

# get mesh class
Mesh = args[0].__mesh__

# remove Mesh from args
args = args[1:]

Expand Down Expand Up @@ -104,7 +105,7 @@ def mesh_or_points_cells_type(*args, **kwargs):
# call mesh manipulation function
points, cells, cell_type = meshfun(points, cells, cell_type, *args, **kwargs)

# return a Mesh if a Mesh was passed
# return a DiscreteGeometry if a DiscreteGeometry was passed
if is_mesh:
return Mesh(points=points, cells=cells, cell_type=cell_type)

Expand Down
161 changes: 112 additions & 49 deletions src/felupe/mesh/_mesh.py
Original file line number Diff line number Diff line change
Expand Up @@ -26,11 +26,29 @@
"""

from copy import deepcopy
from functools import wraps

import numpy as np

from ._convert import (
add_midpoints_edges,
add_midpoints_faces,
add_midpoints_volumes,
collect_edges,
collect_faces,
collect_volumes,
convert,
)
from ._discrete_geometry import DiscreteGeometry
from ._tools import expand, mirror, revolve, rotate, runouts, sweep, triangulate

class Mesh:

def as_mesh(obj):
"Convert a ``DiscreteGeometry`` object to a ``Mesh`` object."
return Mesh(points=obj.points, cells=obj.cells, cell_type=obj.cell_type)


class Mesh(DiscreteGeometry):
"""A mesh with points, cells and optional a specified cell type.

Parameters
Expand All @@ -50,20 +68,6 @@ class Mesh:
Point-connectivity of cells.
cell_type : str or None
A string in VTK-convention that specifies the cell type.
npoints : int
Amount of points.
dim : int
Dimension of mesh point coordinates.
ndof : int
Amount of degrees of freedom.
ncells : int
Amount of cells.
points_with_cells : array
Array with points connected to cells.
points_without_cells : array
Array with points not connected to cells.
cells_per_point : array
Array which counts connected cells per point. Used for averging results.

"""

Expand All @@ -72,36 +76,20 @@ def __init__(self, points, cells, cell_type=None):
self.cells = np.array(cells)
self.cell_type = cell_type

self.update(self.cells)

def update(self, cells, cell_type=None):
"Update the cell and dimension attributes with a given cell array."
self.cells = cells
super().__init__(points=points, cells=cells, cell_type=cell_type)

if cell_type is not None:
self.cell_type = cell_type
self.__mesh__ = Mesh

# obtain dimensions
self.npoints, self.dim = self.points.shape
self.ndof = self.points.size
self.ncells = self.cells.shape[0]

# get number of cells per point
points_in_cell, self.cells_per_point = np.unique(cells, return_counts=True)
def __repr__(self):
header = "<felupe mesh object>"
points = f" Number of points: {len(self.points)}"
cells_header = " Number of cells:"
cells = [f" {self.cell_type}: {self.ncells}"]

# check if there are points without cells
if self.npoints != len(self.cells_per_point):
self.point_has_cell = np.isin(np.arange(self.npoints), points_in_cell)
# update "cells_per_point" ... cells per point
cells_per_point = -np.ones(self.npoints, dtype=int)
cells_per_point[points_in_cell] = self.cells_per_point
self.cells_per_point = cells_per_point
return "\n".join([header, points, cells_header, *cells])

self.points_without_cells = np.arange(self.npoints)[~self.point_has_cell]
self.points_with_cells = np.arange(self.npoints)[self.point_has_cell]
else:
self.points_without_cells = np.array([], dtype=int)
self.points_with_cells = np.arange(self.npoints)
def __str__(self):
return self.__repr__()

def disconnect(self, points_per_cell=None, calc_points=True):
"""Return a new instance of a Mesh with disconnected cells. Optionally, the
Expand Down Expand Up @@ -160,13 +148,88 @@ def copy(self):

return deepcopy(self)

def __repr__(self):
header = "<felupe mesh object>"
points = f" Number of points: {len(self.points)}"
cells_header = " Number of cells:"
cells = [f" {self.cell_type}: {self.ncells}"]
@wraps(expand)
def expand(self, n=11, z=1):
return as_mesh(expand(self, n=n, z=z))

@wraps(rotate)
def rotate(self, angle_deg, axis, center=None):
return as_mesh(rotate(self, angle_deg=angle_deg, axis=axis, center=center))

@wraps(revolve)
def revolve(self, n=11, phi=180, axis=0):
return as_mesh(revolve(self, n=n, phi=phi, axis=axis))

@wraps(sweep)
def sweep(self, decimals=None):
return as_mesh(sweep(self, decimals=decimals))

@wraps(mirror)
def mirror(self, normal=[1, 0, 0], centerpoint=[0, 0, 0], axis=None):
return as_mesh(mirror(self, normal=normal, centerpoint=centerpoint, axis=axis))

@wraps(triangulate)
def triangulate(self, mode=3):
return as_mesh(triangulate(self, mode=mode))

@wraps(runouts)
def add_runouts(
self,
values=[0.1, 0.1],
centerpoint=[0, 0, 0],
axis=0,
exponent=5,
mask=slice(None),
):
return as_mesh(
runouts(
self,
values=values,
centerpoint=centerpoint,
axis=axis,
exponent=exponent,
mask=mask,
)
)

@wraps(convert)
def convert(
self,
order=0,
calc_points=False,
calc_midfaces=False,
calc_midvolumes=False,
):
return as_mesh(
convert(
self,
order=order,
calc_points=calc_points,
calc_midfaces=calc_midfaces,
calc_midvolumes=calc_midvolumes,
)
)

return "\n".join([header, points, cells_header, *cells])
@wraps(collect_edges)
def collect_edges(self):
return collect_edges

def __str__(self):
return self.__repr__()
@wraps(collect_faces)
def collect_faces(self):
return collect_faces

@wraps(collect_volumes)
def collect_volumes(self):
return collect_volumes

@wraps(add_midpoints_edges)
def add_midpoints_edges(self):
return add_midpoints_edges

@wraps(add_midpoints_faces)
def add_midpoints_faces(self):
return add_midpoints_faces

@wraps(add_midpoints_volumes)
def add_midpoints_volumes(self):
return add_midpoints_volumes
3 changes: 2 additions & 1 deletion src/felupe/mesh/_tools.py
Original file line number Diff line number Diff line change
Expand Up @@ -29,7 +29,6 @@

from ..math import rotation_matrix
from ._helpers import mesh_or_data
from ._mesh import Mesh


@mesh_or_data
Expand Down Expand Up @@ -332,6 +331,8 @@ def mirror(
def concatenate(meshes):
"Join a sequence of meshes with identical cell types."

Mesh = meshes[0].__mesh__

points = np.vstack([mesh.points for mesh in meshes])
offsets = np.cumsum(np.insert([mesh.npoints for mesh in meshes][:-1], 0, 0))
cells = np.vstack([offset + mesh.cells for offset, mesh in zip(offsets, meshes)])
Expand Down
Loading