Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

move GeometryBasics to Makie requires block #1322

Merged
merged 1 commit into from
Jan 6, 2023
Merged
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
2 changes: 0 additions & 2 deletions Project.toml
Original file line number Diff line number Diff line change
Expand Up @@ -9,7 +9,6 @@ ConstructionBase = "187b0558-2788-49d3-abe0-74a17ed4e7c9"
EllipsisNotation = "da5c29d0-fa7d-589e-88eb-ea29b0a81949"
FillArrays = "1a297f60-69ca-5386-bcde-b61e274b549b"
ForwardDiff = "f6369f11-7733-5829-9624-2563aa707210"
GeometryBasics = "5c1252a2-5f33-56bf-86c9-59e7332b4326"
HDF5 = "f67ccb44-e63f-5c2f-98bd-6dc0ccc4ba2f"
IfElse = "615f187c-cbe4-4ef1-ba3b-2fcf58d6d173"
LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
Expand Down Expand Up @@ -46,7 +45,6 @@ ConstructionBase = "1.3"
EllipsisNotation = "1.0"
FillArrays = "0.13.2"
ForwardDiff = "0.10.18"
GeometryBasics = "0.3, 0.4"
HDF5 = "0.14, 0.15, 0.16"
IfElse = "0.1"
LinearMaps = "2.7, 3.0"
Expand Down
3 changes: 1 addition & 2 deletions src/Trixi.jl
Original file line number Diff line number Diff line change
Expand Up @@ -44,7 +44,6 @@ using LoopVectorization: LoopVectorization, @turbo, indices
using LoopVectorization.ArrayInterface: static_length
using MPI: MPI
using MuladdMacro: @muladd
using GeometryBasics: GeometryBasics
using Octavian: Octavian, matmul!
using Polyester: @batch # You know, the cheapest threads you can find...
using OffsetArrays: OffsetArray, OffsetVector
Expand Down Expand Up @@ -253,7 +252,7 @@ function __init__()

@require Makie="ee78f7c6-11fb-53f2-987a-cfe4a2b5a57a" begin
include("visualization/recipes_makie.jl")
using .Makie: Makie
using .Makie: Makie, GeometryBasics
export iplot, iplot! # interactive plot
end

Expand Down
107 changes: 107 additions & 0 deletions src/visualization/recipes_makie.jl
Original file line number Diff line number Diff line change
Expand Up @@ -4,6 +4,113 @@
# See https://ranocha.de/blog/Optimizing_EC_Trixi for further details.
@muladd begin

# First some utilities
# Given a reference plotting triangulation, this function generates a plotting triangulation for
# the entire global mesh. The output can be plotted using `Makie.mesh`.
function global_plotting_triangulation_makie(pds::PlotDataSeries{<:PlotData2DTriangulated};
set_z_coordinate_zero = false)
@unpack variable_id = pds
pd = pds.plot_data
@unpack x, y, data, t = pd

makie_triangles = Makie.to_triangles(t)

# trimesh[i] holds GeometryBasics.Mesh containing plotting information on the ith element.
# Note: Float32 is required by GeometryBasics
num_plotting_nodes, num_elements = size(x)
trimesh = Vector{GeometryBasics.Mesh{3, Float32}}(undef, num_elements)
coordinates = zeros(Float32, num_plotting_nodes, 3)
for element in Base.OneTo(num_elements)
for i in Base.OneTo(num_plotting_nodes)
coordinates[i, 1] = x[i, element]
coordinates[i, 2] = y[i, element]
if set_z_coordinate_zero == false
coordinates[i, 3] = data[i, element][variable_id]
end
end
trimesh[element] = GeometryBasics.normal_mesh(Makie.to_vertices(coordinates), makie_triangles)
end
plotting_mesh = merge([trimesh...]) # merge meshes on each element into one large mesh
return plotting_mesh
end

# Returns a list of `Makie.Point`s which can be used to plot the mesh, or a solution "wireframe"
# (e.g., a plot of the mesh lines but with the z-coordinate equal to the value of the solution).
function convert_PlotData2D_to_mesh_Points(pds::PlotDataSeries{<:PlotData2DTriangulated};
set_z_coordinate_zero = false)
@unpack variable_id = pds
pd = pds.plot_data
@unpack x_face, y_face, face_data = pd

if set_z_coordinate_zero
# plot 2d surface by setting z coordinate to zero.
# Uses `x_face` since `face_data` may be `::Nothing`, as it's not used for 2D plots.
sol_f = zeros(eltype(first(x_face)), size(x_face))
else
sol_f = StructArrays.component(face_data, variable_id)
end

# This line separates solution lines on each edge by NaNs to ensure that they are rendered
# separately. The coordinates `xf`, `yf` and the solution `sol_f`` are assumed to be a matrix
# whose columns correspond to different elements. We add NaN separators by appending a row of
# NaNs to this matrix. We also flatten (e.g., apply `vec` to) the result, as this speeds up
# plotting.
xyz_wireframe = GeometryBasics.Point.(map(x->vec(vcat(x, fill(NaN, 1, size(x, 2)))), (x_face, y_face, sol_f))...)

return xyz_wireframe
end

# Creates a GeometryBasics triangulation for the visualization of a ScalarData2D plot object.
function global_plotting_triangulation_makie(pd::PlotData2DTriangulated{<:ScalarData};
set_z_coordinate_zero = false)
@unpack x, y, data, t = pd

makie_triangles = Makie.to_triangles(t)

# trimesh[i] holds GeometryBasics.Mesh containing plotting information on the ith element.
# Note: Float32 is required by GeometryBasics
num_plotting_nodes, num_elements = size(x)
trimesh = Vector{GeometryBasics.Mesh{3, Float32}}(undef, num_elements)
coordinates = zeros(Float32, num_plotting_nodes, 3)
for element in Base.OneTo(num_elements)
for i in Base.OneTo(num_plotting_nodes)
coordinates[i, 1] = x[i, element]
coordinates[i, 2] = y[i, element]
if set_z_coordinate_zero == false
coordinates[i, 3] = data.data[i, element]
end
end
trimesh[element] = GeometryBasics.normal_mesh(Makie.to_vertices(coordinates), makie_triangles)
end
plotting_mesh = merge([trimesh...]) # merge meshes on each element into one large mesh
return plotting_mesh
end

# Returns a list of `GeometryBasics.Point`s which can be used to plot the mesh, or a solution "wireframe"
# (e.g., a plot of the mesh lines but with the z-coordinate equal to the value of the solution).
function convert_PlotData2D_to_mesh_Points(pd::PlotData2DTriangulated{<:ScalarData};
set_z_coordinate_zero = false)
@unpack x_face, y_face, face_data = pd

if set_z_coordinate_zero
# plot 2d surface by setting z coordinate to zero.
# Uses `x_face` since `face_data` may be `::Nothing`, as it's not used for 2D plots.
sol_f = zeros(eltype(first(x_face)), size(x_face))
else
sol_f = face_data
end

# This line separates solution lines on each edge by NaNs to ensure that they are rendered
# separately. The coordinates `xf`, `yf` and the solution `sol_f`` are assumed to be a matrix
# whose columns correspond to different elements. We add NaN separators by appending a row of
# NaNs to this matrix. We also flatten (e.g., apply `vec` to) the result, as this speeds up
# plotting.
xyz_wireframe = GeometryBasics.Point.(map(x->vec(vcat(x, fill(NaN, 1, size(x, 2)))), (x_face, y_face, sol_f))...)

return xyz_wireframe
end


# We set the Makie default colormap to match Plots.jl, which uses `:inferno` by default.
default_Makie_colormap() = :inferno

Expand Down
105 changes: 0 additions & 105 deletions src/visualization/utilities.jl
Original file line number Diff line number Diff line change
Expand Up @@ -1377,111 +1377,6 @@ function reference_node_coordinates_2d(dg::DGSEM)
end


# Given a reference plotting triangulation, this function generates a plotting triangulation for
# the entire global mesh. The output can be plotted using `Makie.mesh`.
function global_plotting_triangulation_makie(pds::PlotDataSeries{<:PlotData2DTriangulated};
set_z_coordinate_zero = false)
@unpack variable_id = pds
pd = pds.plot_data
@unpack x, y, data, t = pd

makie_triangles = Makie.to_triangles(t)

# trimesh[i] holds GeometryBasics.Mesh containing plotting information on the ith element.
# Note: Float32 is required by GeometryBasics
num_plotting_nodes, num_elements = size(x)
trimesh = Vector{GeometryBasics.Mesh{3, Float32}}(undef, num_elements)
coordinates = zeros(Float32, num_plotting_nodes, 3)
for element in Base.OneTo(num_elements)
for i in Base.OneTo(num_plotting_nodes)
coordinates[i, 1] = x[i, element]
coordinates[i, 2] = y[i, element]
if set_z_coordinate_zero == false
coordinates[i, 3] = data[i, element][variable_id]
end
end
trimesh[element] = GeometryBasics.normal_mesh(Makie.to_vertices(coordinates), makie_triangles)
end
plotting_mesh = merge([trimesh...]) # merge meshes on each element into one large mesh
return plotting_mesh
end

# Returns a list of `Makie.Point`s which can be used to plot the mesh, or a solution "wireframe"
# (e.g., a plot of the mesh lines but with the z-coordinate equal to the value of the solution).
function convert_PlotData2D_to_mesh_Points(pds::PlotDataSeries{<:PlotData2DTriangulated};
set_z_coordinate_zero = false)
@unpack variable_id = pds
pd = pds.plot_data
@unpack x_face, y_face, face_data = pd

if set_z_coordinate_zero
# plot 2d surface by setting z coordinate to zero.
# Uses `x_face` since `face_data` may be `::Nothing`, as it's not used for 2D plots.
sol_f = zeros(eltype(first(x_face)), size(x_face))
else
sol_f = StructArrays.component(face_data, variable_id)
end

# This line separates solution lines on each edge by NaNs to ensure that they are rendered
# separately. The coordinates `xf`, `yf` and the solution `sol_f`` are assumed to be a matrix
# whose columns correspond to different elements. We add NaN separators by appending a row of
# NaNs to this matrix. We also flatten (e.g., apply `vec` to) the result, as this speeds up
# plotting.
xyz_wireframe = GeometryBasics.Point.(map(x->vec(vcat(x, fill(NaN, 1, size(x, 2)))), (x_face, y_face, sol_f))...)

return xyz_wireframe
end

# Creates a GeometryBasics triangulation for the visualization of a ScalarData2D plot object.
function global_plotting_triangulation_makie(pd::PlotData2DTriangulated{<:ScalarData};
set_z_coordinate_zero = false)
@unpack x, y, data, t = pd

makie_triangles = Makie.to_triangles(t)

# trimesh[i] holds GeometryBasics.Mesh containing plotting information on the ith element.
# Note: Float32 is required by GeometryBasics
num_plotting_nodes, num_elements = size(x)
trimesh = Vector{GeometryBasics.Mesh{3, Float32}}(undef, num_elements)
coordinates = zeros(Float32, num_plotting_nodes, 3)
for element in Base.OneTo(num_elements)
for i in Base.OneTo(num_plotting_nodes)
coordinates[i, 1] = x[i, element]
coordinates[i, 2] = y[i, element]
if set_z_coordinate_zero == false
coordinates[i, 3] = data.data[i, element]
end
end
trimesh[element] = GeometryBasics.normal_mesh(Makie.to_vertices(coordinates), makie_triangles)
end
plotting_mesh = merge([trimesh...]) # merge meshes on each element into one large mesh
return plotting_mesh
end

# Returns a list of `GeometryBasics.Point`s which can be used to plot the mesh, or a solution "wireframe"
# (e.g., a plot of the mesh lines but with the z-coordinate equal to the value of the solution).
function convert_PlotData2D_to_mesh_Points(pd::PlotData2DTriangulated{<:ScalarData};
set_z_coordinate_zero = false)
@unpack x_face, y_face, face_data = pd

if set_z_coordinate_zero
# plot 2d surface by setting z coordinate to zero.
# Uses `x_face` since `face_data` may be `::Nothing`, as it's not used for 2D plots.
sol_f = zeros(eltype(first(x_face)), size(x_face))
else
sol_f = face_data
end

# This line separates solution lines on each edge by NaNs to ensure that they are rendered
# separately. The coordinates `xf`, `yf` and the solution `sol_f`` are assumed to be a matrix
# whose columns correspond to different elements. We add NaN separators by appending a row of
# NaNs to this matrix. We also flatten (e.g., apply `vec` to) the result, as this speeds up
# plotting.
xyz_wireframe = GeometryBasics.Point.(map(x->vec(vcat(x, fill(NaN, 1, size(x, 2)))), (x_face, y_face, sol_f))...)

return xyz_wireframe
end


# Find element and triangle ids containing coordinates given as a matrix [ndims, npoints]
function get_ids_by_coordinates!(ids, coordinates, pd)
Expand Down