Gist of one of my FHN MakieLayout figures

fitzhugh_nagumo_makielayout.jl

using Makie, MakieTeX, CairoMakie
# Makie v0.17 and associated packages
using DifferentialEquations, Polynomials

# define the default theme

Makie.set_theme!(Theme(
    font = "CMU Serif",
    fontsize = 12, # 12,
    Axis = (
        xticksize = 6,
        yticksize = 6,
        xtickalign = 1,
        ytickalign = 1,
        spinewidth = 0.7,
        gridwidth = 0.5,
    ),
    Heatmap = (rasterize = 5,),
    Surface = (rasterize = 5,),
))

# Returns a Vector{Point2f0} of points which define the isocline
function isoclines(xs::AbstractVector{<: Real}, ys::AbstractVector{<: Real}, zs::AbstractMatrix{<: Real}, isoval::Real = 0.0)
		return Makie.Contours.contour(xs, ys, zs, eltype(xs)(isoval)) |> Makie.Contours.lines |> x -> map(y -> Point2f.(y.vertices), x) .|>  Makie.GeometryBasics.LineString
end

# extract isoclines from a multivariate function
function isoclines(xs::Makie.IntervalSets.ClosedInterval{<: Real}, ys::Makie.IntervalSets.ClosedInterval{<: Real}, f::Function, isoval::Real = 0.0; samples = 500, ind = 1)
  	xr = LinRange(xs, samples)
    yr = LinRange(ys, samples)
    return isoclines(xr, yr, getindex.(f.(Point2f.(xr, yr')), ind), isoval)
end


function fitzhugh_nagumo(x::Point2{T}; e = 0.08, b = 0.4, g = 0.8, i_ext = 0.2) where T
    v, w = x
    return Point2{T}(
        v - v^3/3 - w + i_ext,    # dv₀
        e * (v - g * w + b)       # dw₀
    )
end


ode_fhn(u, p, t) = fitzhugh_nagumo(u; e = p[1], b = p[2], g = p[3], i_ext = p[4])

u0 = Point2{Float64}(1.0, -0.51)

p_exc = [#=e=# 0.08, #=b=# 0.4, #=g=# 0.8, #=i_ext=# -0.2]
p_osc = [#=e=# 0.08, #=b=# 0.4, #=g=# 0.8, #=i_ext=#  0.5]
tspan = (0.0, 150.0)

prob_exc = ODEProblem(ode_fhn, u0, tspan, p_exc)
sol_exc = solve(prob_exc, Tsit5(), reltol=1e-8, abstol = 1e-8)
plot_sol_exc = sol_exc(LinRange(tspan..., 1000))
# see solution
lines(plot_sol_exc.u; color = plot_sol_exc.t, colormap = :linear_bgy_10_95_c74_n256, linewidth = 4)

prob_osc = ODEProblem(ode_fhn, u0, tspan, p_osc)
sol_osc = solve(prob_osc, Tsit5(), reltol=1e-8, abstol = 1e-8)
plot_sol_osc = sol_osc(LinRange(tspan..., 1000))
# see solution
lines(plot_sol_osc.u; color = plot_sol_osc.t, colormap = :linear_bgy_10_95_c74_n256, linewidth = 4)


fh_exc = (x...) -> ode_fhn(x, p_exc, 0)
fh_osc = (x...) -> ode_fhn(x, p_osc, 0)

plot_interval = -2..2

iso_vosc = isoclines(plot_interval, plot_interval, fh_osc; ind=1)
iso_wosc = isoclines(plot_interval, plot_interval, fh_osc; ind=2)

iso_vexc = isoclines(plot_interval, plot_interval, fh_exc; ind=1)
iso_wexc = isoclines(plot_interval, plot_interval, fh_exc; ind=2)


fig = Figure(resolution = (6, 6) .* 72, figure_padding = 5) # resolution in points)

axs = [
            Axis(fig[1, i];
                xlabel = "𝑣",
                ylabel = "𝑤",
                xlabelpadding = 1,
                ylabelpadding = 1,
                tellwidth = true,
                tellheight = false,
            )
            for i in 1:2
        ]
rowsize!(fig.layout, 1, Aspect(1, 1))
hidexdecorations!.(axs); hideydecorations!.(axs)
setproperty!.(axs, :xlabelvisible, true)
setproperty!.(axs, :ylabelvisible, true)
axs[2].yaxisposition = :right

# setproperty!.(axs, :titlesize, 12)
# `axs[1]` is the excitatory axis.
axs[1].title = "Excitable mode"

vexc_lines = lines!(axs[1], iso_vexc; linewidth = 1, color = Makie.wong_colors()[1])
wexc_lines = lines!(axs[1], iso_wexc; linewidth = 1, color = Makie.wong_colors()[2])

# `axs[2]` is the oscillatory axis.
axs[2].title = "Oscillatory mode"

vosc_lines = lines!(axs[2], iso_vosc; linewidth = 1, color = Makie.wong_colors()[1])
wosc_lines = lines!(axs[2], iso_wosc; linewidth = 1, color = Makie.wong_colors()[2])

labels = [Label(fig[1, i, TopLeft()], string('A'-1+i); font = "CMU Serif Bold") for i in 1:2]

sp1 = streamplot!(axs[1], fh_exc,  plot_interval, plot_interval; linewidth = 0.5, arrow_size = 4, colormap = :linear_bmy_10_95_c71_n256, figure = (resolution = (4, 4) .* 72, figure_padding = 0), rasterize = 15)
sp2 = streamplot!(axs[2], fh_osc,  plot_interval, plot_interval; linewidth = 0.5, arrow_size = 4, colormap = :linear_bmy_10_95_c71_n256, figure = (resolution = (4, 4) .* 72, figure_padding = 0), rasterize = 15)

sp1.density = 0.9
sp2.density = 0.9


sc1 = scatter!(axs[1], [u0]; markersize = 5, markerspace = :pixel, markercolor = :blue)
sc2 = scatter!(axs[2], [u0]; markersize = 5, markerspace = :pixel, markercolor = :blue)

ts_axs = [Axis(fig[2, i]; xlabel = L"t", width = @lift(Fixed($(pixelarea(axs[i].scene)).widths[1]))) for i in 1:2]

setproperty!.(ts_axs, :xlabelpadding, 0)

hidexdecorations!.(ts_axs; ticks = true, ticklabels = false, label = false)
hideydecorations!.(ts_axs; ticks = true, ticklabels = true, label = true)


line_traj_exc = lines!(axs[1], plot_sol_exc.u; color = plot_sol_exc.t, colormap = :linear_bgy_10_95_c74_n256)
line_traj_osc = lines!(axs[2], plot_sol_osc.u; color = plot_sol_osc.t, colormap = :linear_bgy_10_95_c74_n256)

ts_exc_1 = lines!(ts_axs[1], plot_sol_exc.t, first.(plot_sol_exc.u); color = Makie.wong_colors()[3], linewidth = 1.2)
ts_exc_2 = lines!(ts_axs[1], plot_sol_exc.t, last.(plot_sol_exc.u); color = Makie.wong_colors()[4], linewidth = 1.2)

ts_osc_1 = lines!(ts_axs[2], plot_sol_osc.t, first.(plot_sol_osc.u); color = Makie.wong_colors()[3], linewidth = 1.2)
ts_osc_2 = lines!(ts_axs[2], plot_sol_osc.t, last.(plot_sol_osc.u); color = Makie.wong_colors()[4], linewidth = 1.2)
# Now, we create the legend.

leg_elements = [
    vexc_lines, ts_osc_1, [LineElement(points = Point2f.(LinRange(0, 1, 100), 0.5), color = 2:101, colormap = :linear_bmy_10_95_c71_n256), MarkerElement(markersize = 4, points = Point2f[(0.5, 0.5)], marker = '▶')], LineElement(points = Point2f.(LinRange(0, 1, 100), 0.5), color = 1:100, colormap = Reverse(:viridis)),
    wexc_lines, ts_osc_2, sc1,
]

leg = Legend(
    fig[3, 1:2],
    leg_elements,
    fill("", length(leg_elements));
    tellwidth = false,
    tellheight = true,
    rowgap = 0,
    framewidth = 0.25,
    padding = (0, 0, 0, 0),
    framevisible = false,
    nbanks = 4,
    orientation = :vertical,
)

leg.nbanks = 4
leg.orientation = :vertical

leg.blockscene.children[1].plots[end-3].colormap = :linear_bgy_10_95_c74_n256
leg.blockscene.children[1].plots[end-5].colormap = :linear_bmy_10_95_c71_n256
# The legend doesn't need much space...80 points should be enough.
rowgap!(fig.layout, 1, 4)
rowgap!(fig.layout, 2, 4)

# Now comes the shady hacking section.  We add layouted TeX to the legend's internal gridlayout.

ltexs = LTeX.(
    Ref(fig), [
        L"dv = 0", L"v", "Trajectories", "Trajectory from \$(v_0, w_0)\$",
        L"dw = 0", L"w", "Initial point"
    ];
    halign = :left,
)
translate!.(getproperty.(ltexs, :blockscene), 0, 0, 100)

grid_indices = vcat(CartesianIndex.(1, 2:2:8), CartesianIndex.(2, 2:2:6))

for (ind, ltex) in zip(grid_indices, ltexs)
    leg.grid[1, 1][ind[1], ind[2]] = ltex
end


rowsize!(fig.layout, 2, Aspect(2, 0.2))
resize_to_layout!(fig)

save("../model/dynamics/FHN_nullclines_exc_and_osc.pdf", fig; pt_per_unit = 1)
save("fhn_null_show.png", fig; px_per_unit = 3)