David Lambert/mesh

From J Wiki
Jump to navigation Jump to search

This work contributes to the ultimate goal to produce meshes suitable for wire frame graphics and finite element programs, such as CalculiX. The verbs mesh2, mesh3 and plot_mesh provide the main functionality to produce and display uniform rectilinear grids which can later be warped, repositioned, and merged. The meshes are scaled from _1 to 1 on all axes. Index position serves as node or element number, and otherwise, with so many good meshing programs available, this is an yet another effort---in j.

Ultimately---a nice mesher might take a bounded 3D figure, use Schwarz–Christoffel transform to mesh the 2D (possibly warped) surfaces, then mesh the interiors with Delaunay simplexes. Missing, a method to choose interior nodes. I've use such a 2d remesh algorithm in which uniformly spaced points on complex circle worked out nicely on the inverse transformation.

require'plot'

NB. The mesh data structure is elements ; nodes
NB. elements are the connectivity of nodes by index

NB. extraction, although most common is
NB. 'nodes elements'=. MESH
elements =: 0&{::  NB. elements MESH
nodes =: 1&{::     NB. nodes MESH

Note 'element node order'
 Node order is non-standard,
 The n-dimensional linear element node order given by
 coordinates
    ('x';'x y';'x y z'),: |."1@#:@:i.&.>2^1 2 3
 +-+---+-----+
 |x|x y|x y z|
 +-+---+-----+
 |0|0 0|0 0 0|
 |1|1 0|1 0 0|
 | |0 1|0 1 0|
 | |1 1|1 1 0|
 | |   |0 0 1|
 | |   |1 0 1|
 | |   |0 1 1|
 | |   |1 1 1|
 +-+---+-----+
 These are the corner positions of linear elements 0 1, bilinear 0 1 2 3, and trilinear i.8 .
)

Monad=: :[:
Dyad=: [: :

cross =: 1 |. ([ * 1 |. ]) - ] * 1 |. [ NB. j forum

odometer6=: (4 $. $.)@:($&1)  NB. contributed by David Ward Lambert
grid=: ([: <: odometer6 (*"1) 2&%@:<:)&.:(|."1) Monad NB.  grid Nx [Ny [Nz]]
NB.grid =: ([: <: odometer6 (*"1) 2&%@:<:) Monad
square_grid =: ([: grid 2 # ]) Monad

e =: ,: NB. itemize a vector
egen =: ([: ,/ (+/~ ({: * [: i. {.)))~ Dyad NB.  (repetition, nodal increment) egen elements
egrid2 =: (<:@:(1&{) , {.) egen ((1 ,~ <:@:{.) egen [: e 0 1 , 0 1 + {.) Monad
egrid3=: ((<:@:(2&{) , */@:(2&{.)) egen (,. (+ >:@:((<_1 _1)&{)))@:egrid2) Monad

NB.egrid3 =: (] ,. */@:(2&{.)@:[ + ]) egrid2

mesh2=: (egrid2 ; grid) Monad  NB. mesh Nx, Ny
mesh3=: (egrid3 ; grid) Monad  NB. mesh Nx, Ny, Nz

NB. segments ELEMENTS  NB. constructs a list of segments to draw
NB. (>0 1;1 3;3 2;2 0) are corners of mesh
segments =: (>0 1;1 3;3 2;2 0)&{"1  NB. per element arrays of faces
unique_segments =: [: ~. [: /:"1~ [: ,/ segments NB. nub of sorted faces 
plot_mesh =: 'color black' plot j./"1@:(({~ unique_segments)&>~/)

echo 'example: plot_mesh mesh2 4 9'