2

Suppose I have a set of points defining a perimeter of a non rectangular shape in the 2d plane.

I need a function to create a triangular meshing where I can modify the number of triangle cells and return (x,y) coordinates of each cell.

Thank you.

RMPR
  • 3,368
  • 4
  • 19
  • 31
  • 1
    This is a big topic - automatic mesh generation. Generation of 2D triangular meshes has been an active area of research since the 1980s. I'd recommend you search to see what Python libraries are available today. Google is your friend: https://pypi.org/project/MeshPy/ – duffymo Mar 30 '20 at 13:34

2 Answers2

4

You should probably check out

(I'm the author of dmsh and pygmsh.)

Nico Schlömer
  • 53,797
  • 27
  • 201
  • 249
2

I poste what i did finally. Some problems are found when working with Python3 . I advice to run on GoogleColab.

Source: https://github.com/inducer/meshpy/blob/master/examples/test_triangle.py

Run in google colab if you have python3 and you do not want to downgrad to python2

!pip install Meshpy
!pip install pybind11
!pip install pyvtk

from __future__ import division
from __future__ import absolute_import

import meshpy.triangle as triangle
import numpy as np
import numpy.linalg as la
from six.moves import range
import matplotlib.pyplot as pt

def round_trip_connect(start, end):
    return [(i, i+1) for i in range(start, end)] + [(end, start)]



points = [(100,400),
          (400,400),
          (400,280),
          (320,270), 
          (320,160), 
          (120,160), 
          (120,250), 
          (100,250)]

facets = round_trip_connect(0, len(points)-1)

circ_start = len(points)
points.extend((3 * np.cos(angle), 3 * np.sin(angle)) for angle in np.linspace(0, 2*np.pi, 30, endpoint=False))

facets.extend(round_trip_connect(circ_start, len(points)-1))


def needs_refinement(vertices, area):
    bary = np.sum(np.array(vertices), axis=0)/3
    max_area = 5 + (la.norm(bary, np.inf)-1)*5
    return bool(area > max_area)

info = triangle.MeshInfo()
info.set_points(points)
info.set_holes([(0, 0)])
info.set_facets(facets)

mesh = triangle.build(info, refinement_func=needs_refinement)

mesh_points = np.array(mesh.points)
mesh_tris = np.array(mesh.elements)


pt.triplot(mesh_points[:, 0], mesh_points[:, 1], mesh_tris)
pt.show()

print (mesh_points)
Community
  • 1
  • 1