Hello, I’m trying to calculate the number of closed contours in a potential as well the area enclosed by closed,
i.e. In this image below at certain energy level there are 4 closed contours. I’m trying to figure a way to calculate the number of closed contours (4 here) and the area enclosed by them, given the energy levels.
sample code for above (adopted from Pluto.jl notebook example),
using GLMakie
len = 50
x = range(-1, 1, length=len)
y = range(-1, 1, length=len)
four_wells(x,y) = x^4 - x^2 + y^4 - y^2
z = [four_wells(x[i], y[j]) for i in 1:len, j in 1:len]
# energy
height = -0.4
f = Figure(;resolution=(680,400))
ax1 = Axis3(f[1, 1])
ax2 = Axis(f[1, 2])
# Draw the contour in the 3D graph
plane = [height for i in 1:length(x), j in 1:length(y)];
# Drawing the function and the contour on the 3G graph
surface!(ax1, x, y, z, colormap=:viridis)
surface!(ax1, x, y, plane; transparency=true, colormap=[:white,:white])
# Drawing the function and the countour on the heatmap
hm = contourf!(ax2, x, y, z)
Colorbar(f[:, 3], hm; vertical=true)
contour!(ax2, x, y, z; linewidth=3, levels=[height], color=:white)
f
Pluto notebook file for the same
### A Pluto.jl notebook ###
# v0.19.35
using Markdown
using InteractiveUtils
# This Pluto notebook uses @bind for interactivity. When running this notebook outside of Pluto, the following 'mock version' of @bind gives bound variables a default value (instead of an error).
macro bind(def, element)
quote
local iv = try Base.loaded_modules[Base.PkgId(Base.UUID("6e696c72-6542-2067-7265-42206c756150"), "AbstractPlutoDingetjes")].Bonds.initial_value catch; b -> missing; end
local el = $(esc(element))
global $(esc(def)) = Core.applicable(Base.get, el) ? Base.get(el) : iv(el)
el
end
end
# ╔═╡ ca18be6d-8a9b-4c53-bd2e-553ea01ebb13
using PlutoUI, GLMakie
# ╔═╡ 503aa0a0-ece8-4acc-b437-3119055d22bc
k = @bind height PlutoUI.Slider(-1:0.1:1; default=0, show_value=true)
# ╔═╡ 2dab2fcd-77bf-4651-8223-6854863f425a
md"""
# Appendix
"""
# ╔═╡ 58d43864-9a74-42da-96d5-afb59f22cd0f
begin
len = 50
x = range(-1, 1, length=len)
y = range(-1, 1, length=len)
four_wells(x,y) = x^4 - x^2 + y^4 - y^2
z = [four_wells(x[i], y[j]) for i in 1:len, j in 1:len]
end
# ╔═╡ e39373bd-8e12-469f-bd66-2d005ff6cff2
begin
# Setup of the figure
f = Figure(;resolution=(680,400))
ax1 = Axis3(f[1, 1])
ax2 = Axis(f[1, 2])
# Draw the contour in the 3D graph
plane = [height for i in 1:length(x), j in 1:length(y)];
# Drawing the function and the contour on the 3G graph
surface!(ax1, x, y, z, colormap=:viridis)
surface!(ax1, x, y, plane; transparency=true, colormap=[:white,:white])
# Drawing the function and the countour on the heatmap
hm = contourf!(ax2, x, y, z)
Colorbar(f[:, 3], hm; vertical=true)
contour!(ax2, x, y, z; linewidth=3, levels=[height], color=:white)
end;
# ╔═╡ 9c17c0d8-e4cb-41c0-b180-b18a3c87e265
# Displaying the function
f
# ╔═╡ Cell order:
# ╠═503aa0a0-ece8-4acc-b437-3119055d22bc
# ╟─9c17c0d8-e4cb-41c0-b180-b18a3c87e265
# ╠═e39373bd-8e12-469f-bd66-2d005ff6cff2
# ╟─2dab2fcd-77bf-4651-8223-6854863f425a
# ╠═ca18be6d-8a9b-4c53-bd2e-553ea01ebb13
# ╠═58d43864-9a74-42da-96d5-afb59f22cd0f
There’s a python implementation on stackoverflow, any pointers on which packages to look on Julia will be more helpful.
