adtzlr / felupe

:mag: finite element analysis for continuum mechanics of solid bodies
https://felupe.readthedocs.io/
GNU General Public License v3.0
66 stars 10 forks source link

Import `pyvista.UnstructuredGrid` as `MeshContainer`? #830

Open adtzlr opened 3 weeks ago

adtzlr commented 3 weeks ago

Take the points and cells_dict of an unstructured grid and create a MeshContainer. This should simplify things like

import felupe as fem
import pyvista as pv
import skgmsh as sg

source = pv.Polygon(n_sides=6, radius=1, fill=False)
delaunay_2d = sg.frontal_delaunay_2d(edge_source=source, target_sizes=0.1)

mesh = (
    fem.Mesh(
        points=delaunay_2d.points[:, :2],  # in-plane points
        cells=delaunay_2d.cells_dict[pv.CellType.TRIANGLE],  # triangle cells
        cell_type="triangle",
    )
    .flip()
    .rotate(90, axis=2)
)
region = fem.RegionTriangle(mesh)
field = fem.FieldContainer([fem.FieldPlaneStrain(region, dim=2)])
boundaries, loadcase = fem.dof.shear(
    field,
    sym=False,
    moves=(0.75, 0.0, 0.0),
)
solid = fem.SolidBody(umat=fem.LinearElastic(E=1, nu=0.4), field=field)
step = fem.Step(items=[solid], boundaries=boundaries)
job = fem.Job(steps=[step]).evaluate()
ax = solid.imshow("Principal Values of Logarithmic Strain", project=fem.project)

image

to


source = pv.Polygon(n_sides=6, radius=1, fill=False)
delaunay_2d = sg.frontal_delaunay_2d(edge_source=source, target_sizes=0.1)

mesh = fem.mesh.from_unstructured_grid(delaunay_2d)[0].rotate(90, axis=2)
mesh.update(points=mesh.points[:, :2])
# ...

Tasks

tkoyama010 commented 3 weeks ago

Thank you for considering this. Currently, I am considering the following points in the new version of scikit-gmsh.

  1. Defines a new class, FrontalDelaunay2D, and recommends using that.
  2. The PyVista's PolyData cannot guarantee the integrity of the geometry in . Therefore, the FrontalDelaunay2D class uses shapely to define the geometry.
  3. Defines a new class, FrontalDelaunay3D, and recommends using that.
  4. Delaunay3D uses PyVista's PolyData as geometry (same as delaunay_3d function).
adtzlr commented 3 weeks ago

Thanks for the information! I'll start with Task 1 (importing an unstructured grid) anyway and when a new version of scikit-gmsh is ready I'll add an example where the mesh is created by scikit-gmsh.

adtzlr commented 1 week ago

Updated Examples

These code-blocks use felupe.MeshContainer.from_unstructured_grid(). The geometries / meshes are based on the examples in https://scikit-gmsh.readthedocs.io/.

2D Plane Strain

import skgmsh as sg
import felupe as fem

shell = [(0, 0, 0), (0, 10, 0), (10, 10, 0), (10, 0, 0), (0, 0, 0)]
holes = [[(2, 2, 0), (2, 4, 0), (4, 4, 0), (4, 2, 0), (2, 2, 0)]]
alg = sg.Delaunay2D(shell=shell, holes=holes)
alg.cell_size = 0.3
grid = alg.mesh.split_bodies()[0]

mesh = fem.MeshContainer.from_unstructured_grid(grid, dim=2)[0]
mesh.update(cells=mesh.cells[:, ::-1])
region = fem.RegionTriangle(mesh)
field = fem.FieldContainer([fem.FieldPlaneStrain(region, dim=2)])
boundaries, loadcase = fem.dof.shear(
    field,
    sym=False,
    moves=(3, 0, 0),
)
solid = fem.SolidBody(umat=fem.LinearElasticLargeStrain(E=1, nu=0.3), field=field)
step = fem.Step(items=[solid], boundaries=boundaries)
job = fem.Job(steps=[step]).evaluate()
ax = solid.imshow("Principal Values of Logarithmic Strain", project=fem.project)

grafik

3D

import pyvista as pv
import skgmsh as sg

edge_source = pv.Cylinder(resolution=16)
edge_source.merge(pv.PolyData(edge_source.points), merge_points=True, inplace=True)
alg = sg.Delaunay3D(edge_source)
alg.cell_size = 0.2

mesh = fem.MeshContainer.from_unstructured_grid(alg.mesh, dim=3)[0]
mesh.update(cells=mesh.cells[:, ::-1])
region = fem.RegionTetra(mesh)
field = fem.FieldContainer([fem.Field(region, dim=3)])
boundaries, loadcase = fem.dof.shear(
    field,
    axes=(1, 0),
    sym=False,
    moves=(1.5, -0.5, 0.5),
)
solid = fem.SolidBody(umat=fem.LinearElasticLargeStrain(E=1, nu=0.3), field=field)
step = fem.Step(items=[solid], boundaries=boundaries)
job = fem.Job(steps=[step]).evaluate()
solid.plot("Principal Values of Logarithmic Strain", project=fem.project).show()

grafik