Skip to content
Open
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
90 changes: 62 additions & 28 deletions src/primitives/surfacemeshes/mesh_cuboid.jl
Original file line number Diff line number Diff line change
Expand Up @@ -8,8 +8,8 @@ returns Mesh(vertices, faces)

Function returns a simplicial areal mesh of a cuboid. It takes kwarg - generator
:compsciencemeshes - is default, and gives a structured mesh of
a cuboid, the dimensions of the cuboid are approximated by multiples
of edge length -
a cuboid, the edge lengths are approximated to fit the dimensions
of the geometry -

Number of faces = 2*(2*m*n) + 2*(2*n*p) + 2*(2*m*p) -> front, back; top, bottom;
left, right
Expand Down Expand Up @@ -65,10 +65,10 @@ function meshcuboid(len::F, breadth::F, width::F, edge_len::F;
if generator == :gmsh
msh = gmshcuboid(len, breadth, width, edge_len)
elseif generator == :compsciencemeshes
# @info "Generating a structured mesh: The dimensions of the cuboid are
# approximated by multiples of edge length.
# @info "Generating a structured mesh: The edge lengths are
# approximated to fit the dimensions of the geometry.
# For exact dimensions/ unstructured grid, use kwarg - generator = :gmsh"
msh = mesh_cuboid(len, width, width, edge_len)
msh = mesh_cuboid(len, breadth, width, edge_len)
else
@error "generators are gmsh and compsciencemeshes only"
end
Expand Down Expand Up @@ -107,17 +107,51 @@ returns an areal structured mesh of a cuboid.

"""
@generated function mesh_cuboid(a, b, c, h)
Core.println("Generating a structured mesh: The dimensions of the cuboid are
approximated by multiples of edge length. For exact dimensions and
Core.println("Generating a structured mesh: The edge lengths are
approximated to fit the dimensions of the geometry. For exact edge length and
unstructured grids, use kwarg - generator = :gmsh")
return :(mesh_cuboid_impl(a, b, c, h))
end

function mesh_cuboid_impl(a::F, b::F, c::F, h::F) where F
# if isapprox(a%h, F(0)) && isapprox(b%h, F(0) && isapprox(c%h, F(0)))
n = Int(round(a/h)) # number of elements along a
m = Int(round(b/h)) # number of elements along b
p = Int(round(c/h)) #number of elements along c
n_u = Int(ceil(a/h)) # number of elements along a
n_d = Int(floor(a/h))
m_u = Int(ceil(b/h)) # number of elements along b
m_d = Int(floor(b/h))
p_u = Int(ceil(c/h)) #number of elements along c
p_d = Int(floor(c/h))

#determining edge lengths closest to given edge length
h_1u = F((a/n_u))
h_1d = F((a/n_d))
if abs(h_1u - h) < abs(h_1d - h)
h_1 = h_1u
n = n_u
else
h_1 = h_1d
n = n_d
end

h_2u = F(b/m_u)
h_2d = F(b/m_d)
if abs(h_2u - h) < abs(h_2d - h)
h_2 = h_2u
m = m_u
else
h_2 = h_2d
m = m_d
end

h_3u = F(c/p_u)
h_3d = F(c/p_d)
if abs(h_3u - h) < abs(h_3d - h)
h_3 = h_3u
p = p_u
else
h_3 = h_3d
p = p_d
end

nodes = Vector{SVector{3, F}}(undef, 2*(m*n + m*p + n*p + 1))
faces = Vector{SVector{3, Int64}}(undef, 4*m*n + 4*n*p + 4*m*p)
Expand All @@ -130,12 +164,12 @@ function mesh_cuboid_impl(a::F, b::F, c::F, h::F) where F
for iy in 1:m
#nodes
#front
nodes[(ix - 1)*(m + 1) + iy] = SVector((ix - 1)*h, (iy - 1)*h, F(0))
nodes[(ix - 1)*(m + 1) + iy] = SVector((ix - 1)*h_1, (iy - 1)*h_2, F(0))
#back
nodes[back_node + (ix - 1)*(m + 1) + iy] = SVector(
(ix - 1)*h,
(iy - 1)*h,
p*h
(ix - 1)*h_1,
(iy - 1)*h_2,
p*h_3
)
#faces
#front
Expand Down Expand Up @@ -163,18 +197,18 @@ function mesh_cuboid_impl(a::F, b::F, c::F, h::F) where F
end
# for the mth element in y-direction
#front
nodes[ix*(m + 1)] = SVector((ix - 1)*h, m*h, F(0))
nodes[ix*(m + 1)] = SVector((ix - 1)*h_1, m*h_2, F(0))
#back
nodes[back_node + ix*(m + 1)] = SVector(
(ix - 1)*h, m*h, p*h)
(ix - 1)*h_1, m*h_2, p*h_3)
end
# for ix = n
for iy in 0:m
#front
nodes[n*(m + 1) + iy + 1] = SVector(n*h, (iy*h), F(0))
nodes[n*(m + 1) + iy + 1] = SVector(n*h_1, (iy*h_2), F(0))
#back
nodes[back_node + n*(m + 1) + iy + 1] = SVector(
n*h, (iy*h), p*h)
n*h_1, (iy*h_2), p*h_3)
end

# along y - z
Expand All @@ -187,7 +221,7 @@ function mesh_cuboid_impl(a::F, b::F, c::F, h::F) where F
left_face = 4*m*n + 4*n*p + 2*m*p + (iz - 2)*2*m
#ix in 1 -> left nodes
nodes[left_node_front + iy] = SVector(
F(0), (iy - 1)*h, (iz - 1)*h)
F(0), (iy - 1)*h_2, (iz - 1)*h_3)
if iz == 2 && iy != (m + 1)
#left faces
faces[left_face + (2*iy - 1)] = SVector(
Expand Down Expand Up @@ -237,9 +271,9 @@ function mesh_cuboid_impl(a::F, b::F, c::F, h::F) where F
end
#ix in n -> right nodes
nodes[left_node_front + (m + 2*n - 1) + iy] = SVector(
a,
(iy - 1)*h,
(iz - 1)*h
n*h_1,
(iy - 1)*h_2,
(iz - 1)*h_3
)
end
#iz = p + 1
Expand Down Expand Up @@ -372,10 +406,10 @@ function mesh_cuboid_impl(a::F, b::F, c::F, h::F) where F
#nodes for iz = p
#iy = 1 -> bottom nodes
nodes[(n + 1)*(m + 1) + (p - 2)*(2*(m + n)) + (m + 1) + (ix - 2)*2 + 1] = SVector(
(ix - 1)*h, F(0), (p - 1)*h)
(ix - 1)*h_1, F(0), (p - 1)*h_3)
#iy = m -> top nodes
nodes[(n + 1)*(m + 1) + (p - 2)*(2*(m + n)) + (m + 1) + (ix - 1)*2] = SVector(
(ix - 1)*h, b, (p - 1)*h)
(ix - 1)*h_1, m*h_2, (p - 1)*h_3)
if (ix != n)
# iz = 1
#bottom faces
Expand Down Expand Up @@ -475,10 +509,10 @@ function mesh_cuboid_impl(a::F, b::F, c::F, h::F) where F
#nodes
#iy = 1 -> bottom nodes
nodes[(n + 1)*(m + 1) + (iz - 2)*(2*(m + n)) + (m + 1) + (ix - 2)*2 + 1] = SVector(
(ix - 1)*h, F(0), (iz - 1)*h)
(ix - 1)*h_1, F(0), (iz - 1)*h_3)
#iy = m -> top nodes
nodes[(n + 1)*(m + 1) + (iz - 2)*(2*(m + n)) + (m + 1) + (ix - 1)*2] = SVector(
(ix - 1)*h, b, (iz - 1)*h)
(ix - 1)*h_1, m*h_2, (iz - 1)*h_3)
else
#for iy = 1 -> bottom faces
faces[2*m*n + (ix - 1)*2*p + (2*iz - 1)] = SVector(
Expand All @@ -505,10 +539,10 @@ function mesh_cuboid_impl(a::F, b::F, c::F, h::F) where F
#nodes
#iy = 1 -> bottom nodes
nodes[(n + 1)*(m + 1) + (iz - 2)*(2*(m + n)) + (m + 1) + (ix - 2)*2 + 1] = SVector(
(ix - 1)*h, F(0), (iz - 1)*h)
(ix - 1)*h_1, F(0), (iz - 1)*h_3)
#iy = m -> top nodes
nodes[(n + 1)*(m + 1) + (iz - 2)*(2*(m + n)) + (m + 1) + (ix - 1)*2] = SVector(
(ix - 1)*h, b, (iz - 1)*h)
(ix - 1)*h_1, m*h_2, (iz - 1)*h_3)
end
end
end
Expand Down
47 changes: 35 additions & 12 deletions src/primitives/surfacemeshes/mesh_rectangle.jl
Original file line number Diff line number Diff line change
Expand Up @@ -53,8 +53,8 @@ function meshrectangle(len::F, breadth::F, edge_length::F, udim = 3;
@assert udim==3
msh = meshrectangle_quadrilateral(len, breadth, edge_length)
elseif generator == :compsciencemeshes && element == :triangle
# @info "Generating a structured mesh: The dimensions of the rectangle are
# approximated by multiples of edge length.
# @info "Generating a structured mesh: The edge lengths are
# approximated to fit the dimensions of the geometry.
# For exact dimensions/ unstructured grid, use kwarg - generator = :gmsh"
msh = mesh_rectangle(len, breadth, edge_length, udim)
else
Expand All @@ -74,8 +74,8 @@ end
structured mesh.
"""
@generated function mesh_rectangle(a, b, h, udim)
Core.println("Generating a structured mesh: The dimensions of the cuboid are
approximated by multiples of edge length. For exact dimensions and
Core.println("Generating a structured mesh: The edge lengths are
approximated to fit the dimensions of the geometry. For exact dimensions and
unstructured grids, use kwarg - generator = :gmsh")
return :(mesh_rectangle_impl(a, b, h, udim))
end
Expand All @@ -84,18 +84,41 @@ end
function mesh_rectangle_impl(a::F, b::F, h::F, udim) where F
@assert udim == 2 || udim == 3 "Universal dimension can only be 2 or 3"
#structured mesh: isapprox(a%h, F(0)) && isapprox(b%h, F(0))
n = Int(round(a/h)) # number of elements along a
m = Int(round(b/h)) # number of elements along b
n_u = Int(ceil(a/h)) # number of elements along a
n_d = Int(floor(a/h))
m_u = Int(ceil(b/h)) # number of elements along b
m_d = Int(floor(b/h))

#determining edge lengths closest to given edge length
h_1u = F((a/n_u))
h_1d = F((a/n_d))
if abs(h_1u - h) < abs(h_1d - h)
h_1 = h_1u
n = n_u
else
h_1 = h_1d
n = n_d
end

h_2u = F(b/m_u)
h_2d = F(b/m_d)
if abs(h_2u - h) < abs(h_2d - h)
h_2 = h_2u
m = m_u
else
h_2 = h_2d
m = m_d
end

nodes = zeros(SVector{udim, F}, (m + 1)*(n + 1))
faces = Vector{SVector{3, Int64}}(undef, 2*m*n)

for ix in 0 : n - 1
for iy in 1 : m
if udim == 3
node = SVector((ix)*h, (iy - 1)*h, F(0))
node = SVector((ix)*h_1, (iy - 1)*h_2, F(0))
else
node = SVector((ix)*h, (iy - 1)*h)
node = SVector((ix)*h_1, (iy - 1)*h_2)
end

nodes[(ix)*(m + 1) + iy] = node
Expand All @@ -115,18 +138,18 @@ function mesh_rectangle_impl(a::F, b::F, h::F, udim) where F
end
# for the mth element in y-direction
if udim == 3
nodes[(ix + 1)*(m + 1)] = SVector((ix)*h, m*h, F(0))
nodes[(ix + 1)*(m + 1)] = SVector((ix)*h_1, m*h_2, F(0))
else
nodes[(ix + 1)*(m + 1)] = SVector((ix)*h, m*h)
nodes[(ix + 1)*(m + 1)] = SVector((ix)*h_1, m*h_2)
end

end
# for ix = n
for iy in 0 : m
if udim == 3
nodes[n*(m + 1) + iy + 1] = SVector(n*h, (iy*h), F(0))
nodes[n*(m + 1) + iy + 1] = SVector(n*h_1, (iy*h_2), F(0))
else
nodes[n*(m + 1) + iy + 1] = SVector(n*h, (iy*h))
nodes[n*(m + 1) + iy + 1] = SVector(n*h_1, (iy*h_2))
end
end

Expand Down
Loading
Loading