Generating a nice mesh on a sphere S^2

Specific Domains Visualization
1

I would like to generate a nice mesh/grid (i.e., as homogeneous/isotropic as possible) on the unit sphere S^2 in R^3. In particular, I would like to generate two possible grid systems: 1) based on triangles; or 2) based on the Goldberg polyhedra (i.e., the so-called Geodesic Grid). I am wondering whether there are Julia packages to generate those meshes/grids.
Thanks a lot in advance!

2

Gmsh.jl can generate triangular meshes for spheres (via the addSphere(x,y,z,radius) function).

2 Likes
3

Thanks a lot, @stevengj ! I’ll try that package.

4

Do you want quadrilaterals or triangles? Either can be generated with FinEtools. For example https://github.com/PetrKryslUCSD/FinEtoolsMeshing.jl/blob/8f2460a191e46fa3779fd2ac6a6662654aae9a8e/test/test_meshing.jl#L2963

1 Like
5

Here is the GMSH geo file I use to generate spherical (and ellipsoidal meshes) (borrowed from http://bempp.com/)

// from bempp shapes.py
//
// Return an ellipsoid grid.
//
// Parameters
// ----------
// r1 : float
// Radius of first major axis
// r2 : float
// Radius of second major axis
// r3 : float
// Radius of third major axis
// origin : tuple
// Tuple specifying the origin of the ellipsoid
// h : float
// Element size.
// “”"
// stub = “”"

max_freq = 2000;
C = 343.2;

wav = C / max_freq;

lc = wav / 6; // Characteristic Length – i.e. basic mesh size parameter

r1 = 0.25;
r2 = 0.25;
r3 = 0.25;

orig0 = 0.0;
orig1 = 0.0;
orig2 = 0.0;

cl = lc;
cl2 = lc;

Point(1) = {orig0,orig1,orig2,cl};
Point(2) = {orig0+r1,orig1,orig2,cl};
Point(3) = {orig0,orig1+r2,orig2,cl2};
Point(4) = {orig0-r1,orig1,orig2,cl};
Point(5) = {orig0,orig1-r2,orig2,cl2};
Point(6) = {orig0,orig1,orig2-r3,cl2};
Point(7) = {orig0,orig1,orig2+r3,cl2};

Ellipse(1) = {2,1,2,3};
Ellipse(2) = {3,1,4,4};
Ellipse(3) = {4,1,4,5};
Ellipse(4) = {5,1,2,2};
Ellipse(5) = {3,1,3,6};
Ellipse(6) = {6,1,5,5};
Ellipse(7) = {5,1,5,7};
Ellipse(8) = {7,1,3,3};
Ellipse(9) = {2,1,2,7};
Ellipse(10) = {7,1,4,4};
Ellipse(11) = {4,1,4,6};
Ellipse(12) = {6,1,2,2};

Line Loop(13) = {2,8,-10};
Ruled Surface(14) = {13};

Line Loop(15) = {10,3,7};
Ruled Surface(16) = {15};

Line Loop(17) = {-8,-9,1};
Ruled Surface(18) = {17};

Line Loop(19) = {-11,-2,5};
Ruled Surface(20) = {19};

Line Loop(21) = {-5,-12,-1};
Ruled Surface(22) = {21};

Line Loop(23) = {-3,11,6};
Ruled Surface(24) = {23};

Line Loop(25) = {-7,4,9};
Ruled Surface(26) = {25};

Line Loop(27) = {-4,12,-6};
Ruled Surface(28) = {27};

Surface Loop(29) = {28,26,16,14,20,24,22,18};

// GMSH did NOT like this compound constraint…
// Compound Surface{28,26,16,14,20,24,22,18};

Volume(30) = {29};
Physical Surface(1) = {14};
Physical Surface(2) = {16};
Physical Surface(3) = {18};
Physical Surface(4) = {20};
Physical Surface(5) = {22};
Physical Surface(6) = {24};
Physical Surface(7) = {26};
Physical Surface(8) = {28};

// To generate a curvilinear mesh and optimize it to produce provably valid
// curved elements (see A. Johnen, J.-F. Remacle and C. Geuzaine. Geometric
// validity of curvilinear finite elements. Journal of Computational Physics
// 233, pp. 359-372, 2013; and T. Toulorge, C. Geuzaine, J.-F. Remacle,
// J. Lambrechts. Robust untangling of curvilinear meshes. Journal of
// Computational Physics 254, pp. 8-26, 2013), you can uncomment the following
// lines:
//
// Mesh.ElementOrder = 2;
// Mesh.HighOrderOptimize = 2;

Mesh.Algorithm = 2;
Mesh.ElementOrder = 2;
Mesh.HighOrderOptimize = 2;

2 Likes
6

gmshsphere
gmshsphere770Ă—951 173 KB

3 Likes
7

image

2 Likes
8

Here is a solution in pure Julia:

using Meshes

# sphere geometry
sphere = Sphere((0.,0.,0.), 1.)

# polygonal mesh with quadrangles + triangles
mesh = discretize(sphere, RegularDiscretization(30, 30))

# triangulate the mesh so that it contains only triangles
tmesh = triangulate(mesh)

# you can refine a mesh further with triangles
ttmesh = refine(tmesh, TriRefinement())

# or refine it further with quadrangles
qmesh = refine(tmesh, QuadRefinement())

Here is a visualization of all these meshes:

using MeshViz

# choose a Makie backend
import GLMakie as Mke

fig = Mke.Figure()
viz(fig[1,1], mesh, showfacets = true)
viz(fig[1,2], tmesh, showfacets = true)
viz(fig[1,3], ttmesh, showfacets = true)
viz(fig[1,4], qmesh, showfacets = true)
fig

image
image1893Ă—767 147 KB

More discretization methods are available, check the docs: Discretization · Meshes.jl

We are always happy to review PRs with new methods.

5 Likes
9

Thank you everyone for suggesting all these packages! I’ll check them all, and decide which one I would like to use.