A simple mesh generator for non-convex regions in n-D space
An unstructured simplex requires a choice of mesh points (vertex nodes) and
a triangulation. This is a simple and short algorithm that improves the
quality of a mesh by relocating the mesh points according to a relaxation
scheme of forces in a truss structure. The topology of the truss is reset
using Delaunay triangulation. A (sufficiently smooth) user supplied signed
distance function (fd
) indicates if a given node is inside or
outside the region. Points outside the region are projected back to the
boundary.
distmeshnd(fdist, fh, h, box, pfix = array(dim = c(0, ncol(box))), ..., ptol = 0.001, ttol = 0.1, deltat = 0.1, geps = 0.1 * h, deps = sqrt(.Machine$double.eps) * h)
fdist |
Vectorized signed distance function, for example
|
fh |
Vectorized function, for example |
h |
Initial distance between mesh nodes. |
box |
|
pfix |
|
... |
parameters that are passed to |
ptol |
Algorithm stops when all node movements are smaller than
|
ttol |
Controls how far the points can move (relatively) before a
retriangulation with |
deltat |
Size of the time step in Euler's method. |
geps |
Tolerance in the geometry evaluations. |
deps |
Stepsize Δ x in numerical derivative computation for distance function. |
This is an implementation of original Matlab software of Per-Olof Persson.
Excerpt (modified) from the reference below:
‘The algorithm is based on a mechanical analogy between a triangular mesh and a n-D truss structure. In the physical model, the edges of the Delaunay triangles of a set of points correspond to bars of a truss. Each bar has a force-displacement relationship F(L,L0) depending on its current length L and its unextended length L0.’
‘External forces on the structure come at the boundaries, on which
external forces have normal orientations. These external forces are just
large enough to prevent nodes from moving outside the boundary. The
position of the nodes are the unknowns, and are found by solving for a
static force equilibrium. The hope is that (when fh = function(p)
return(rep(1,nrow(p)))
), the lengths of all the bars at equilibrium will
be nearly equal, giving a well-shaped triangular mesh.’
See the references below for all details. Also, see the comments in the
source file of distmesh2d
.
m
-by-n
matrix with node positions.
Implement in C/Fortran
Translate other functions of the Matlab package
Raoul Grasman; translated from original Matlab sources of Per-Olof Persson.
P.-O. Persson, G. Strang, A Simple Mesh Generator in MATLAB. SIAM Review, Volume 46 (2), pp. 329-345, June 2004
## Not run: # examples distmeshnd require(rgl) fd = function(p, ...) sqrt((p^2)%*%c(1,1,1)) - 1 # also predefined as `mesh.dsphere' fh = function(p,...) rep(1,nrow(p)) # also predefined as `mesh.hunif' bbox = matrix(c(-1,1),2,3) p = distmeshnd(fd,fh,0.2,bbox, maxiter=100) # this may take a while: # press Esc to get result of current iteration ## End(Not run)
Please choose more modern alternatives, such as Google Chrome or Mozilla Firefox.