text:
### A Pluto.jl notebook ###
# v0.17.0

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 el = $(esc(element))
        global $(esc(def)) = Core.applicable(Base.get, el) ? Base.get(el) : missing
        el
    end
end

# ╔═╡ 104b6690-212b-11ec-3a2b-25c420f62534
begin
	using Plots
	using LinearAlgebra
	using Kronecker
	using PlutoUI
	using SparseArrays
	#https://juliapackages.com/p/kronecker
	md"""## Packages"""
end

# ╔═╡ 21712488-f636-4861-913c-012693cb524d
md"""# Order of tensor product basis (TPB)"""

# ╔═╡ 89b58dd0-c70c-4566-b991-2d739990f2ff
md"""
Note, that throughout the thesis a colexicographical ordering of the tensor product is preferred:
$\ket{a_1 b_1}, \ket{a_2 b_1}, \ket{a_1 b_2}, \ket{a_2 b_2}$. In this case the tensor product is realized via the conventional Kronecker product but using flipped arguments:
$A \otimes B \rightarrow \boldsymbol{B} \otimes \boldsymbol{A}$.

Here you can choose the ordering of the tensor product basis of the tensor product ``\mathcal{H}_{\mathcal{A}} \otimes \mathcal{H}_{\mathcal{B}}``:

👉 $(@bind order Select(["colex" => "Colexicographical order", "lex" => "Lexicographical order"]))
"""

# ╔═╡ 23c019dd-53d2-461e-a6ba-f1b5da2f6710
if order == "colex"
md"""
``\begin{pmatrix}{\color{magenta}a}&{\color{red}C}\\ \boldsymbol{b} & D   \end{pmatrix} \overset{\text{c}} {\otimes}\begin{pmatrix} 1 & \large 3 \\ {\color{blue}2} & \scriptsize {\color{green}4}\end{pmatrix} = \begin{pmatrix}{\color{magenta}a}1 & {\color{red}C}1 & {\color{magenta}a} \large 3 & {\color{red}C} \large 3 \\\boldsymbol{b}1 & D1 & \boldsymbol{b} \large 3 & D \large 3 \\ {\color{magenta}a}{\color{blue}2} & {\color{red}C}{\color{blue}2} & {\color{magenta}a}\scriptsize {\color{green}4} & {\color{red}C}\scriptsize {\color{green}4} \\ \boldsymbol{b}{\color{blue}2} & D{\color{blue}2} & \boldsymbol{b}\scriptsize {\color{green}4} & D\scriptsize {\color{green}4} \end{pmatrix}, \quad A \overset{\text{c}}{\otimes} B:= \textbf{B} \otimes \textbf{A}``
"""
else
md"""
``\begin{pmatrix}{\color{magenta}a}&{\color{red}C}\\ \boldsymbol{b} & D   \end{pmatrix} {\otimes}\begin{pmatrix} 1 & \large 3 \\ {\color{blue}2} & \scriptsize {\color{green}4}\end{pmatrix}  = \begin{pmatrix}{\color{magenta}a}1 & {\color{magenta}a} \large 3 & {\color{red}C}1 & {\color{red}C} \large 3 \\{\color{magenta}a}{\color{blue}2} & {\color{magenta}a}\scriptsize {\color{green}4} & {\color{red}C}{\color{blue}2} & {\color{red}C}\scriptsize {\color{green}4} \\  \boldsymbol{b}1 & \boldsymbol{b} \large 3 & D1 & D \large 3 \\  \boldsymbol{b}{\color{blue}2} & \boldsymbol{b}\scriptsize {\color{green}4} & D{\color{blue}2} & D\scriptsize {\color{green}4} \end{pmatrix}``
"""
end

# ╔═╡ f5992b42-ff99-479a-b093-70a9bfd85655
@show order

# ╔═╡ 0a95d032-8576-492f-bee1-ef67303bcb3a
begin
	if order == "colex"
		⊗(A,B) = kron(B,A)
	else
		⊗(A,B) = kron(A,B)
	end
end

# ╔═╡ 7458474c-3c4d-4a99-9640-4bfd2e6465ea
@bind 🍎 Scrubbable(1:10)

# ╔═╡ edf8bfba-ca6b-4adf-9a7d-e618ca511a80
@bind time_r Clock()

# ╔═╡ d7f3dc01-1201-41de-8a63-c3e2bc020cf6
σ₄ = 5

# ╔═╡ 2e2854f6-8cab-4252-a580-7abc71efd2a0
rand(2,2)  ⊗ I(2)

# ╔═╡ a646342f-3177-4579-8625-cb7d7ea27a49
🍎 + σ₄ + time_r

# ╔═╡ 792748ec-199c-41d8-a801-6cd97881ee72
🍎*4

# ╔═╡ 85bca31d-c3be-46b6-b330-b80648472378
md"""
# System Hamiltonian of bipartite system
"""

# ╔═╡ ebd5efea-94b6-419e-8ff1-9dcf624ddaba
md"""
The system Hamiltonian for a bipartite system $\mathcal{H}_{\mathcal{A}} \otimes \mathcal{H}_{\mathcal{B}}$ is given by

$H = \sum_{i = \{\mathcal{A}, \mathcal{B}\}} \varepsilon_i \hat{n}_i + \sum_{i=1}^3 \alpha_i (\mathbb{1}_{\mathcal{A}} \otimes {\sigma}_i) + \sum_{i=1}^3 \beta_i ({\sigma}_i \otimes \mathbb{1}_{\mathcal{B}}) + \sum_{ij=1}^3 \gamma_{ij} ({\sigma}_i \otimes {\sigma}_j)$

with Pauli matrices
$\sigma_1 = \sigma_x = \begin{pmatrix} 0 & 1 \\ 1 & 0  \end{pmatrix}$,
$\sigma_2 = \sigma_y = \begin{pmatrix} 0 & -i \\ i & 0  \end{pmatrix}$,
$\sigma_3 = \sigma_z = \begin{pmatrix} 1 & 0 \\ 0 & -1\end{pmatrix}$

"""

# ╔═╡ 913d405a-e525-49d6-9021-227543b37953
md"""
# Initial state ``ρ(t_0)``
"""

# ╔═╡ a36e5515-e7f9-4c34-a48e-52a004d80f23
md"""
Switch from different initial states 👉 $(@bind initial_state Select(["thesis" => "Example in 📖", "example" => "Another example 📍","slider"=> "Choose your own initial state 📐"]))
"""

# ╔═╡ 7479e663-e533-442c-8324-c60f4849c470
md"""
# Dynamical map for different initial states
"""

# ╔═╡ 0872a46c-4f4b-49fa-8152-9193b71f332b
md"""
## Initial product state
"""

# ╔═╡ 5e1377b7-c784-4e42-8190-ffbf48f0ecc1
md"""
##  Initial classical-quantum correlated state
"""

# ╔═╡ 822ea5ab-4a3a-45ee-903c-ff9bad9d5f16
md"""
## Example of non--positive dynamical map
"""

# ╔═╡ 4cdab1a4-1784-4430-ae4f-50a4459d83a4
md"""
# Example of a non-accessible map
The transposition superoperator

"""

# ╔═╡ a558b847-3f09-473c-9ec5-70badcda7589
md"""
# Appendix: Functions and Initializations
"""

# ╔═╡ 8407fc8c-95b4-4536-bbda-b3f4945518b0
TableOfContents()

# ╔═╡ 7f43573f-010f-455a-82c7-f2f509cce3d3
begin
	time_slider = @bind time Slider(0:0.1:20, default = 1, show_value = true)
	md""" ## Interactive widgets
	Definition of Sliders"""
end

# ╔═╡ e8c989d6-9ba6-40d1-ad36-d768172d5f5e
time_slider

# ╔═╡ 002504f4-e74a-4e8c-80fa-f4c5f2e90df0
time_slider

# ╔═╡ d28b8274-7435-484a-bde7-61443d7e4ae2
begin
	scrub_1 = @bind s1 Scrubbable(-0.5:0.01:0.5, format=".03")
	scrub_2 = @bind s2 Scrubbable(-0.5:0.01:0.5, format=".03")
	scrub_3 = @bind s3 Scrubbable(-0.5:0.01:0.5, format=".03")
	scrub_4 = @bind s4 Scrubbable(-0.5:0.01:0.5, format=".03")
	scrub_5 = @bind s5 Scrubbable(-0.5:0.01:0.5, format=".03")
	scrub_6 = @bind s6 Scrubbable(-0.5:0.01:0.5, format=".03")
	scrub_7 = @bind s7 Scrubbable(-0.5:0.01:0.5, format=".03", default = 0.25)
	scrub_8 = @bind s8 Scrubbable(-0.5:0.01:0.5, format=".03")
	scrub_9 = @bind s9 Scrubbable(-0.50:0.01:0.50, format=".03")

	scrub_a1 = @bind a1 Scrubbable(-0.50:0.01:0.50, format=".03", default = -0.5)
	scrub_a2 = @bind a2 Scrubbable(-0.50:0.01:0.50, format=".03")
	scrub_a3 = @bind a3 Scrubbable(-0.50:0.01:0.50, format=".03")

	scrub_b1 = @bind b1 Scrubbable(-0.50:0.01:0.50, format=".03", default = -0.5)
	scrub_b2 = @bind b2 Scrubbable(-0.50:0.01:0.50, format=".03")
	scrub_b3 = @bind b3 Scrubbable(-0.50:0.01:0.50, format=".03")

	scrub_i1 = @bind i1 Scrubbable(-0.50:0.01:0.50, format=".03", default = -0.5)
	scrub_i2 = @bind i2 Scrubbable(-0.50:0.01:0.50, format=".03")
	scrub_i3 = @bind i3 Scrubbable(-0.50:0.01:0.50, format=".03")

	md" Definition of Scrubbables"
end

# ╔═╡ fc12024b-2785-4f89-a457-95a2deb5045f
begin
	if initial_state == "slider"
		md"""
		You can choose the parameters of the initial state:

		$\rho(t_0) =\frac{1}{4}\mathbb{1}_{\mathcal{A}} \otimes \mathbb{1}_{\mathcal{B}}+ \sum_{i=1}^3 \alpha_i (\mathbb{1}_{\mathcal{A}} \otimes {\sigma}_i) + \sum_{i=1}^3 \beta_i ({\sigma}_i \otimes \mathbb{1}_{\mathcal{B}}) + \sum_{ij=1}^3 \gamma_{ij} ({\sigma}_i \otimes {\sigma}_j)$

		``\alpha_1`` 👉 $(scrub_a1),  ``\alpha_2`` 👉 $(scrub_a2), ``\alpha_3`` 👉 $(scrub_a3)

		``\beta_1`` 👉 $(scrub_b1),  ``\beta_2`` 👉 $(scrub_b2), ``\beta_3`` 👉 $(scrub_b3)

		``\gamma_{11}`` 👉 $(scrub_1),  ``\gamma_{12}`` 👉 $(scrub_2), ``\gamma_{13}`` 👉 $(scrub_3)

		``\gamma_{21}`` 👉 $(scrub_4),  ``\gamma_{22}`` 👉 $(scrub_5), ``\gamma_{23}`` 👉 $(scrub_6)

		``\gamma_{31}`` 👉 $(scrub_7),  ``\gamma_{32}`` 👉 $(scrub_8), ``\gamma_{33}`` 👉 $(scrub_9)

		"""

	end
end

# ╔═╡ 1e4adc4d-52cf-41df-ad62-cab016f55397
[scrub_i1, scrub_i2, scrub_i3]

# ╔═╡ 0c09a3b6-920d-4732-af67-86b36bd88292
scrub_a1, scrub_a2, scrub_a3

# ╔═╡ dc95c17d-de0f-4b72-b57e-0dab473a69cd
[scrub_1,scrub_2,scrub_3,scrub_4, scrub_5, scrub_6, scrub_7, scrub_8, scrub_9]

# ╔═╡ 75ebae9f-1c69-4d92-9cdf-a0f4b0971576
begin
σ = [[0 1; 1 0],  [0 -im; im 0],  [1 0; 0 -1]]

	function  qu_bit_system(epsilon  = 1,
			alpha = [0,0,0],
			beta= [0,0,0],
			gamma = zeros(3,3))

	return H = epsilon * kron(I(2),I(2)) +
		sum([(I(2) ⊗ σ[i]) * alpha[i] for i = 1:3]) +
		sum([(σ[i] ⊗ I(2)) * beta[i] for i=1:3]) +
		sum([gamma[i,j] * (σ[i] ⊗ σ[j]) for i=1:3, j=1:3])
	end

	function density_pauli(α)
		return α[1]/2 * I(2) + sum([α[i+1] * σ[i] for i = 1:3])
	end

	function density_matrix(n, neg = 0, real_valued = false)
		#trace 1 remains preserved!
		dia = rand(n)
		if neg > 0
			dia[rand(1:n, neg)] .*= -1
		end
		dia = dia/sum(dia)
		if real_valued
			temp = rand(n,n)
		else
			temp = rand(n,n) + rand(n,n) * im
		end
		ew, ev = eigen(Hermitian(temp'*temp))
		d_mat = Hermitian(ev* Diagonal(dia) * ev')
		return d_mat
	end

	function unitary(n)
		A = (rand(n,n) .-0.5) + im*(rand(n,n) .-0.5)
		return eigvecs(A'+A)
	end

	function unity(i,n)
		E = spzeros(Bool,n,n)
		E[i] = true
    	return E
	end
	function partial_trace_matrix(dim_1, dim_2, traced_system, order="lex")
# traces out the system indicated by traced_system
# for C = kron(A,B), tracing over the first system (traced_system = 1);
# yields: reshape(p1*C[:], 4,4)/tr(A)  == B
#
# when colex ordered:
# flip traced_system

		if order == "colex"
			traced_system = 3 - traced_system
		end

		if traced_system == 1
			T = [kron(sparse(I, dim_1,dim_1), unity(i,dim_2))[:] for i = 1:dim_2^2]
			col = reduce(vcat,[findnz(T[i])[1] for i = 1:dim_2^2])
			row = kron(1:dim_2^2, ones(Int64,dim_1))
			return sparse(row, col, trues(dim_2^2*dim_1))
		else
			T = [kron(unity(i, dim_1), sparse(I, dim_2, dim_2))[:] for i = 1:dim_1^2]
			col = reduce(vcat, [findnz(T[i])[1] for i = 1:dim_1^2])
			row = kron(1:dim_1^2, ones(Int64,dim_2))
			return sparse(row, col, trues(dim_1^2*dim_2))
			#return    Tr_B = [reshape(kron(reshape(unity(i,dim_1), dim_1, dim_1), I(dim_2)), 1, dim_2^2 * dim_1^2)[1,j] for i = 1:dim_1^2, j =1:dim_1^2*dim_2^2]
		end
	end

	function choi(A, dim_1)
		mm, nn = size(A)
		n = mm ÷ dim_1
		return reshape(permutedims(reshape(A,n,dim_1,n,dim_1),(1,3,2,4)), mm,mm)
	end

	function choi_rect(A,dim_1, dim_2, dim_3, dim_4)
		return reshape(permutedims(reshape(A, dim_1, dim_2, dim_3, dim_4), (1,3,2,4)), dim_1*dim_3, dim_2*dim_4)
	end

	function choi_index(i, n)
		return  (i.-1) .% n[1] .+ 1 .+
		n[1] .* ((i.-1) .% prod(n[1:3]) .÷ prod(n[1:2]) .+
		n[3] .* ((i.-1) .% prod(n[1:2]) .÷ n[1] .+
		n[2] .* ((i.-1) .÷ prod(n[1:3]))))
	end

	function ptr(ρ, dim_1, traced_system = 1, order="lex")
		# get partial trace superoperator
		N,M = size(ρ)
		dim_2 = N÷dim_1
		PT = partial_trace_matrix(dim_1, dim_2, traced_system, order)
		new_dim = traced_system == 1 ? dim_2 : dim_1
		# apply superoperator to vectorized density matrix and reshape to new dimensions
		σ = reshape(PT * ρ[:], new_dim, new_dim)
		return σ
	end

	function factorize(A, dim_1, dim_2, order="lex")
		σ = ptr(A, dim_1, 2, order)
		ρ = ptr(A, dim_1, 1, order)
		A_tensor_prod = (σ ⊗ ρ)
		rho_c = A - A_tensor_prod
		return A_tensor_prod, rho_c
	end

	function normed_vector(k)
		v = (rand(k) .- 0.5) + im*(rand(k) .- 0.5)
		return v/norm(v)
	end

	function prob_dist(k)
		p = rand(k)
		return p/sum(p)
	end

	function PPT(rho,dim_1)
		dim_2 = size(rho,1)÷dim_1
		if dim_1 *dim_2 > 6
			error("positive partial trace does not work for too large dimensions")
		end
		return eigvals(reshape(permutedims(reshape(rho,dim_1,dim_2,dim_1,dim_2), [1,4,3,2]),dim_1*dim_2, dim_1*dim_2))

	end

	function sep_state(k=1,n=2)
		VVV = sum([(normed_vector(n) ⊗ normed_vector(n)) for _=1:k] .* (prob_dist(k)))
		return VVV / norm(VVV)
	end

	#test for k-block positivity
	# strategy: generate random matrices in 2x2 and test whether Choi matrix maps to positive operators, if true for all, then positive.
	# alternatively: check k-positivity: k as optional parameter. Create states with Schmidt rank k
	#Schmidt rank:
	function pos_test(choi, k = 1, n = 2)

		runs = 10000
		test = zeros(runs)
		vector = zeros(ComplexF64, n^2)
		for i = 1:runs

			#sum([density_matrix(2) for _=1:k] * p_rand/sum(p_rand))

			VVN = sep_state(k,n)

			test[i] = real(VVN'* choi *VVN)

			if real(test[i]) < 0
				vector = VVN
			end
		end

		return sum(test .> 0) / runs

	end

md"""
# Functions
"""
end

# ╔═╡ c8e02ee8-2c21-4a37-9875-c4a463e0b462
begin
	ϵ = [1,2]
	α = [0,0,0]
	β = [0,0,0]
	γ = [1 0 0; 0 0 0; 0 0 0]

	H = qu_bit_system(0, α, β, γ) + I(2) ⊗ Diagonal([0,ϵ[2]]) +
		Diagonal([0,ϵ[1]]) ⊗ I(2)

end

# ╔═╡ a5e38ca3-8184-4160-a6b1-4c1783c6e185
begin

	if initial_state == "thesis"
		rho_0 = qu_bit_system(1/4, [0.0, 0.00, +0.0],
							  [0.00, -0.0, -0.0], [0 0 0; 0 0 0; 1/4 0 0])
	elseif initial_state == "example"
		rho_0 = qu_bit_system(1/4, [0.05, 0.00, +0.0],
							  [0.00, 0.05, -0.0], [0 0 0.05; 0.05 0 0; 0.05 0 0])
	else
		rho_0 = qu_bit_system(1/4, [a1, a2, a3],
								[b1, b2, b3], [s1 s2 s3; s4 s5 s6;s7 s8 s9])
	end

	# alternative
	Dim = 2
	RHO_b = [density_matrix(Dim) for _=1:Dim]
	RHO_s = density_matrix(Dim)
	evals, Pro = eigen(RHO_s)

	prob = prob_dist(Dim)

	#rho_0 = sum([kron(Pro[:,i] * Pro[:,i]', RHO_b[i]) for i = 1:Dim] .*prob)

	rho_t_0, rho_c = factorize(rho_0,2,2, order)
end

# ╔═╡ 4f87569b-777e-4e97-ac80-8b7b76115d2a
rho_0

# ╔═╡ b595349b-b572-4cc8-92b1-07631cc8fad7
eigvals(rho_0) # eigenvalues of the initial state

# ╔═╡ 8cd0438c-c480-4b73-93af-79e3f3a5dfa9
eigvals(rho_0)

# ╔═╡ 596a9128-2fff-47e8-8b6f-e9fd8686a279
[ptr(rho_0,2,2,order), ptr(rho_0,2,1,order)] #σ(t_0) and ρ_B(t_0)

# ╔═╡ c578581e-467e-482c-9cce-2ad55d32267f
PPT(rho_0,2) #initial state is separable, but is it classical quantum correlated?

# ╔═╡ 8c779399-7db0-421e-b1b2-659d28fc8c41
begin
	#using the product state rho_t_0
	rho_B = ptr(rho_t_0,2,1, order) # the bath denstiy matrix
	sigma_0 = ptr(rho_t_0,2,2, order) # the system density matrix
	compare = sigma_0 ⊗ rho_B - rho_t_0 #check if this is really a product state

	U = exp(-im*H*time)
	#NOTE: the choi matrix notation needs a transpose the match the lexicographical order
	if order == "colex"
		Λ = choi(choi(U,2) * kron(rho_B, I(2)) * choi(U,2)',2)#
	else
		Λ = transpose(choi(transpose(choi(U,2)) * kron(I(2), rho_B) * transpose(choi(U,2))',2)) # this is a real Kronecker product
	end

	Λ_c = ptr(U * rho_c * U', 2, 2, order)[:] * collect(I(2)[:])'

	test = density_matrix(2)
	dyn_1 = reshape(Λ * test[:], 2,2)

	dyn_2 = ptr(U * (test ⊗ rho_B) * U',2,2, order)
	dyn_1 - dyn_2
end

# ╔═╡ 9330d1cf-e0e1-41c1-a876-5fae11633ce8
Λ_c + Λ

# ╔═╡ 04583ab3-2364-4c32-be8e-d3e65ca5321e
round.(Λ + Λ_c, digits=2)

# ╔═╡ 5d450f8d-c8c4-4a58-b2a6-97bfe109f91f
begin
	#check

	norm(ptr(U*rho_0*U',2,2, order)- reshape((Λ+Λ_c) * ptr(rho_0,2,2, order)[:],2,2))

end

# ╔═╡ 4fb4c614-ed74-411e-91db-641c331c75fa
eigvals(choi(Λ + Λ_c,2))
# map is not completely positive

# ╔═╡ 36df576e-5111-4d66-9170-2eb9573cea32
pos_test(choi(Λ + Λ_c,2),1,2)
# dynamical map is not positive

# ╔═╡ 2cab05f7-bc59-4d26-a535-495b5ab104c0
begin
	sig_n = 1/2 * I(2) + i1 * σ[1] + i2 * σ[2] + i3*σ[3]
	#pure state: a1 = -0.5
	eigvals(sig_n)
	#sig_n is not compatible:
	eigvals((sig_n ⊗ rho_B) + rho_c)
end

# ╔═╡ faf456bc-9145-4d87-ad4d-6ecc8c78bd41
eigvals(reshape((Λ + Λ_c)* sig_n[:],2,2))
# yields negative result

# ╔═╡ e63d17ba-767c-4136-8dc3-3da609f8dda5
reshape((Λ + Λ_c)* sig_n[:],2,2)

# ╔═╡ 5695cdc2-0da6-4ef1-9aa9-d553f9b194e4
begin

K_p = [1/2*(I(2) + σ[3]), 1/2*(I(2) - σ[3]), 1/sqrt(2)*σ[1]]
K_m = 1/sqrt(2)*[im * σ[2]]
Λ_t = sum([kron(conj.(K_p[i]), K_p[i]) for i = 1:3]) - kron(conj.(K_m[1]), K_m[1])

choi(Λ_t,2)

C_t = sum([K_p[i][:] * K_p[i][:]' for i = 1:3]) - K_m[1][:] * K_m[1][:]'

C_t - kron(I(2), density_pauli([0,a1,a2,a3]))
end

# ╔═╡ 3140729e-d516-434d-ba32-3ff9705ee83b
begin
keep_working(text=md"The answer is not quite right.", title="Keep working on it!") = Markdown.MD(Markdown.Admonition("danger", title, [text]));

almost(text, title="Almost there!") = Markdown.MD(Markdown.Admonition("warning", title, [text]));

hint(text, title ="Hint") = Markdown.MD(Markdown.Admonition("hint", title, [text]));

correct(text=md"Great! You got the right answer! Let's move on to the next section.", title="Got it!") = Markdown.MD(Markdown.Admonition("correct", title, [text]));
md" Definition of Boxes"
end

# ╔═╡ 5afc42f7-3253-438d-96a8-d665a1b9a1cd
if any(eigvals(rho_0) .< 0)
	keep_working(md"Tune the parameters such, that the initial state becomes positive: ``\rho(t_0)> 0``.

The current eigenvalues are:
		$(string(round.(eigvals(rho_0), digits = 2)))" ,"Initial state  not valid!")
end

# ╔═╡ 00000000-0000-0000-0000-000000000001
PLUTO_PROJECT_TOML_CONTENTS = """
[deps]
Kronecker = "2c470bb0-bcc8-11e8-3dad-c9649493f05e"
LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
Plots = "91a5bcdd-55d7-5caf-9e0b-520d859cae80"
PlutoUI = "7f904dfe-b85e-4ff6-b463-dae2292396a8"
SparseArrays = "2f01184e-e22b-5df5-ae63-d93ebab69eaf"

[compat]
Kronecker = "~0.5.0"
Plots = "~1.22.3"
PlutoUI = "~0.7.12"
"""

# ╔═╡ 00000000-0000-0000-0000-000000000002
PLUTO_MANIFEST_TOML_CONTENTS = """
# This file is machine-generated - editing it directly is not advised

[[AbstractFFTs]]
deps = ["LinearAlgebra"]
git-tree-sha1 = "485ee0867925449198280d4af84bdb46a2a404d0"
uuid = "621f4979-c628-5d54-868e-fcf4e3e8185c"
version = "1.0.1"

[[Adapt]]
deps = ["LinearAlgebra"]
git-tree-sha1 = "84918055d15b3114ede17ac6a7182f68870c16f7"
uuid = "79e6a3ab-5dfb-504d-930d-738a2a938a0e"
version = "3.3.1"

[[ArgTools]]
uuid = "0dad84c5-d112-42e6-8d28-ef12dabb789f"

[[Artifacts]]
uuid = "56f22d72-fd6d-98f1-02f0-08ddc0907c33"

[[Base64]]
uuid = "2a0f44e3-6c83-55bd-87e4-b1978d98bd5f"

[[Bzip2_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "19a35467a82e236ff51bc17a3a44b69ef35185a2"
uuid = "6e34b625-4abd-537c-b88f-471c36dfa7a0"
version = "1.0.8+0"

[[Cairo_jll]]
deps = ["Artifacts", "Bzip2_jll", "Fontconfig_jll", "FreeType2_jll", "Glib_jll", "JLLWrappers", "LZO_jll", "Libdl", "Pixman_jll", "Pkg", "Xorg_libXext_jll", "Xorg_libXrender_jll", "Zlib_jll", "libpng_jll"]
git-tree-sha1 = "f2202b55d816427cd385a9a4f3ffb226bee80f99"
uuid = "83423d85-b0ee-5818-9007-b63ccbeb887a"
version = "1.16.1+0"

[[ChainRulesCore]]
deps = ["Compat", "LinearAlgebra", "SparseArrays"]
git-tree-sha1 = "e8a30e8019a512e4b6c56ccebc065026624660e8"
uuid = "d360d2e6-b24c-11e9-a2a3-2a2ae2dbcce4"
version = "1.7.0"

[[ColorSchemes]]
deps = ["ColorTypes", "Colors", "FixedPointNumbers", "Random"]
git-tree-sha1 = "a851fec56cb73cfdf43762999ec72eff5b86882a"
uuid = "35d6a980-a343-548e-a6ea-1d62b119f2f4"
version = "3.15.0"

[[ColorTypes]]
deps = ["FixedPointNumbers", "Random"]
git-tree-sha1 = "024fe24d83e4a5bf5fc80501a314ce0d1aa35597"
uuid = "3da002f7-5984-5a60-b8a6-cbb66c0b333f"
version = "0.11.0"

[[Colors]]
deps = ["ColorTypes", "FixedPointNumbers", "Reexport"]
git-tree-sha1 = "417b0ed7b8b838aa6ca0a87aadf1bb9eb111ce40"
uuid = "5ae59095-9a9b-59fe-a467-6f913c188581"
version = "0.12.8"

[[Compat]]
deps = ["Base64", "Dates", "DelimitedFiles", "Distributed", "InteractiveUtils", "LibGit2", "Libdl", "LinearAlgebra", "Markdown", "Mmap", "Pkg", "Printf", "REPL", "Random", "SHA", "Serialization", "SharedArrays", "Sockets", "SparseArrays", "Statistics", "Test", "UUIDs", "Unicode"]
git-tree-sha1 = "31d0151f5716b655421d9d75b7fa74cc4e744df2"
uuid = "34da2185-b29b-5c13-b0c7-acf172513d20"
version = "3.39.0"

[[CompilerSupportLibraries_jll]]
deps = ["Artifacts", "Libdl"]
uuid = "e66e0078-7015-5450-92f7-15fbd957f2ae"

[[Contour]]
deps = ["StaticArrays"]
git-tree-sha1 = "9f02045d934dc030edad45944ea80dbd1f0ebea7"
uuid = "d38c429a-6771-53c6-b99e-75d170b6e991"
version = "0.5.7"

[[CovarianceEstimation]]
deps = ["LinearAlgebra", "Statistics", "StatsBase"]
git-tree-sha1 = "bc3930158d2be029e90b7c40d1371c4f54fa04db"
uuid = "587fd27a-f159-11e8-2dae-1979310e6154"
version = "0.2.6"

[[DataAPI]]
git-tree-sha1 = "cc70b17275652eb47bc9e5f81635981f13cea5c8"
uuid = "9a962f9c-6df0-11e9-0e5d-c546b8b5ee8a"
version = "1.9.0"

[[DataStructures]]
deps = ["Compat", "InteractiveUtils", "OrderedCollections"]
git-tree-sha1 = "7d9d316f04214f7efdbb6398d545446e246eff02"
uuid = "864edb3b-99cc-5e75-8d2d-829cb0a9cfe8"
version = "0.18.10"

[[DataValueInterfaces]]
git-tree-sha1 = "bfc1187b79289637fa0ef6d4436ebdfe6905cbd6"
uuid = "e2d170a0-9d28-54be-80f0-106bbe20a464"
version = "1.0.0"

[[Dates]]
deps = ["Printf"]
uuid = "ade2ca70-3891-5945-98fb-dc099432e06a"

[[DelimitedFiles]]
deps = ["Mmap"]
uuid = "8bb1440f-4735-579b-a4ab-409b98df4dab"

[[Distributed]]
deps = ["Random", "Serialization", "Sockets"]
uuid = "8ba89e20-285c-5b6f-9357-94700520ee1b"

[[Downloads]]
deps = ["ArgTools", "LibCURL", "NetworkOptions"]
uuid = "f43a241f-c20a-4ad4-852c-f6b1247861c6"

[[EarCut_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "3f3a2501fa7236e9b911e0f7a588c657e822bb6d"
uuid = "5ae413db-bbd1-5e63-b57d-d24a61df00f5"
version = "2.2.3+0"

[[Expat_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "b3bfd02e98aedfa5cf885665493c5598c350cd2f"
uuid = "2e619515-83b5-522b-bb60-26c02a35a201"
version = "2.2.10+0"

[[FFMPEG]]
deps = ["FFMPEG_jll"]
git-tree-sha1 = "b57e3acbe22f8484b4b5ff66a7499717fe1a9cc8"
uuid = "c87230d0-a227-11e9-1b43-d7ebe4e7570a"
version = "0.4.1"

[[FFMPEG_jll]]
deps = ["Artifacts", "Bzip2_jll", "FreeType2_jll", "FriBidi_jll", "JLLWrappers", "LAME_jll", "Libdl", "Ogg_jll", "OpenSSL_jll", "Opus_jll", "Pkg", "Zlib_jll", "libass_jll", "libfdk_aac_jll", "libvorbis_jll", "x264_jll", "x265_jll"]
git-tree-sha1 = "d8a578692e3077ac998b50c0217dfd67f21d1e5f"
uuid = "b22a6f82-2f65-5046-a5b2-351ab43fb4e5"
version = "4.4.0+0"

[[FixedPointNumbers]]
deps = ["Statistics"]
git-tree-sha1 = "335bfdceacc84c5cdf16aadc768aa5ddfc5383cc"
uuid = "53c48c17-4a7d-5ca2-90c5-79b7896eea93"
version = "0.8.4"

[[Fontconfig_jll]]
deps = ["Artifacts", "Bzip2_jll", "Expat_jll", "FreeType2_jll", "JLLWrappers", "Libdl", "Libuuid_jll", "Pkg", "Zlib_jll"]
git-tree-sha1 = "21efd19106a55620a188615da6d3d06cd7f6ee03"
uuid = "a3f928ae-7b40-5064-980b-68af3947d34b"
version = "2.13.93+0"

[[Formatting]]
deps = ["Printf"]
git-tree-sha1 = "8339d61043228fdd3eb658d86c926cb282ae72a8"
uuid = "59287772-0a20-5a39-b81b-1366585eb4c0"
version = "0.4.2"

[[FreeType2_jll]]
deps = ["Artifacts", "Bzip2_jll", "JLLWrappers", "Libdl", "Pkg", "Zlib_jll"]
git-tree-sha1 = "87eb71354d8ec1a96d4a7636bd57a7347dde3ef9"
uuid = "d7e528f0-a631-5988-bf34-fe36492bcfd7"
version = "2.10.4+0"

[[FriBidi_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "aa31987c2ba8704e23c6c8ba8a4f769d5d7e4f91"
uuid = "559328eb-81f9-559d-9380-de523a88c83c"
version = "1.0.10+0"

[[GLFW_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Libglvnd_jll", "Pkg", "Xorg_libXcursor_jll", "Xorg_libXi_jll", "Xorg_libXinerama_jll", "Xorg_libXrandr_jll"]
git-tree-sha1 = "0c603255764a1fa0b61752d2bec14cfbd18f7fe8"
uuid = "0656b61e-2033-5cc2-a64a-77c0f6c09b89"
version = "3.3.5+1"

[[GR]]
deps = ["Base64", "DelimitedFiles", "GR_jll", "HTTP", "JSON", "Libdl", "LinearAlgebra", "Pkg", "Printf", "Random", "Serialization", "Sockets", "Test", "UUIDs"]
git-tree-sha1 = "c2178cfbc0a5a552e16d097fae508f2024de61a3"
uuid = "28b8d3ca-fb5f-59d9-8090-bfdbd6d07a71"
version = "0.59.0"

[[GR_jll]]
deps = ["Artifacts", "Bzip2_jll", "Cairo_jll", "FFMPEG_jll", "Fontconfig_jll", "GLFW_jll", "JLLWrappers", "JpegTurbo_jll", "Libdl", "Libtiff_jll", "Pixman_jll", "Pkg", "Qt5Base_jll", "Zlib_jll", "libpng_jll"]
git-tree-sha1 = "ef49a187604f865f4708c90e3f431890724e9012"
uuid = "d2c73de3-f751-5644-a686-071e5b155ba9"
version = "0.59.0+0"

[[GeometryBasics]]
deps = ["EarCut_jll", "IterTools", "LinearAlgebra", "StaticArrays", "StructArrays", "Tables"]
git-tree-sha1 = "58bcdf5ebc057b085e58d95c138725628dd7453c"
uuid = "5c1252a2-5f33-56bf-86c9-59e7332b4326"
version = "0.4.1"

[[Gettext_jll]]
deps = ["Artifacts", "CompilerSupportLibraries_jll", "JLLWrappers", "Libdl", "Libiconv_jll", "Pkg", "XML2_jll"]
git-tree-sha1 = "9b02998aba7bf074d14de89f9d37ca24a1a0b046"
uuid = "78b55507-aeef-58d4-861c-77aaff3498b1"
version = "0.21.0+0"

[[Glib_jll]]
deps = ["Artifacts", "Gettext_jll", "JLLWrappers", "Libdl", "Libffi_jll", "Libiconv_jll", "Libmount_jll", "PCRE_jll", "Pkg", "Zlib_jll"]
git-tree-sha1 = "7bf67e9a481712b3dbe9cb3dac852dc4b1162e02"
uuid = "7746bdde-850d-59dc-9ae8-88ece973131d"
version = "2.68.3+0"

[[Graphite2_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "344bf40dcab1073aca04aa0df4fb092f920e4011"
uuid = "3b182d85-2403-5c21-9c21-1e1f0cc25472"
version = "1.3.14+0"

[[Grisu]]
git-tree-sha1 = "53bb909d1151e57e2484c3d1b53e19552b887fb2"
uuid = "42e2da0e-8278-4e71-bc24-59509adca0fe"
version = "1.0.2"

[[HTTP]]
deps = ["Base64", "Dates", "IniFile", "Logging", "MbedTLS", "NetworkOptions", "Sockets", "URIs"]
git-tree-sha1 = "14eece7a3308b4d8be910e265c724a6ba51a9798"
uuid = "cd3eb016-35fb-5094-929b-558a96fad6f3"
version = "0.9.16"

[[HarfBuzz_jll]]
deps = ["Artifacts", "Cairo_jll", "Fontconfig_jll", "FreeType2_jll", "Glib_jll", "Graphite2_jll", "JLLWrappers", "Libdl", "Libffi_jll", "Pkg"]
git-tree-sha1 = "8a954fed8ac097d5be04921d595f741115c1b2ad"
uuid = "2e76f6c2-a576-52d4-95c1-20adfe4de566"
version = "2.8.1+0"

[[HypertextLiteral]]
git-tree-sha1 = "72053798e1be56026b81d4e2682dbe58922e5ec9"
uuid = "ac1192a8-f4b3-4bfe-ba22-af5b92cd3ab2"
version = "0.9.0"

[[IOCapture]]
deps = ["Logging", "Random"]
git-tree-sha1 = "f7be53659ab06ddc986428d3a9dcc95f6fa6705a"
uuid = "b5f81e59-6552-4d32-b1f0-c071b021bf89"
version = "0.2.2"

[[IniFile]]
deps = ["Test"]
git-tree-sha1 = "098e4d2c533924c921f9f9847274f2ad89e018b8"
uuid = "83e8ac13-25f8-5344-8a64-a9f2b223428f"
version = "0.5.0"

[[InteractiveUtils]]
deps = ["Markdown"]
uuid = "b77e0a4c-d291-57a0-90e8-8db25a27a240"

[[IterTools]]
git-tree-sha1 = "05110a2ab1fc5f932622ffea2a003221f4782c18"
uuid = "c8e1da08-722c-5040-9ed9-7db0dc04731e"
version = "1.3.0"

[[IteratorInterfaceExtensions]]
git-tree-sha1 = "a3f24677c21f5bbe9d2a714f95dcd58337fb2856"
uuid = "82899510-4779-5014-852e-03e436cf321d"
version = "1.0.0"

[[JLLWrappers]]
deps = ["Preferences"]
git-tree-sha1 = "642a199af8b68253517b80bd3bfd17eb4e84df6e"
uuid = "692b3bcd-3c85-4b1f-b108-f13ce0eb3210"
version = "1.3.0"

[[JSON]]
deps = ["Dates", "Mmap", "Parsers", "Unicode"]
git-tree-sha1 = "8076680b162ada2a031f707ac7b4953e30667a37"
uuid = "682c06a0-de6a-54ab-a142-c8b1cf79cde6"
version = "0.21.2"

[[JpegTurbo_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "d735490ac75c5cb9f1b00d8b5509c11984dc6943"
uuid = "aacddb02-875f-59d6-b918-886e6ef4fbf8"
version = "2.1.0+0"

[[Kronecker]]
deps = ["LinearAlgebra", "NamedDims", "SparseArrays", "StatsBase"]
git-tree-sha1 = "f6e3cc35572a6be64308ffb0a56d70be36dc6a85"
uuid = "2c470bb0-bcc8-11e8-3dad-c9649493f05e"
version = "0.5.0"

[[LAME_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "f6250b16881adf048549549fba48b1161acdac8c"
uuid = "c1c5ebd0-6772-5130-a774-d5fcae4a789d"
version = "3.100.1+0"

[[LZO_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "e5b909bcf985c5e2605737d2ce278ed791b89be6"
uuid = "dd4b983a-f0e5-5f8d-a1b7-129d4a5fb1ac"
version = "2.10.1+0"

[[LaTeXStrings]]
git-tree-sha1 = "c7f1c695e06c01b95a67f0cd1d34994f3e7db104"
uuid = "b964fa9f-0449-5b57-a5c2-d3ea65f4040f"
version = "1.2.1"

[[Latexify]]
deps = ["Formatting", "InteractiveUtils", "LaTeXStrings", "MacroTools", "Markdown", "Printf", "Requires"]
git-tree-sha1 = "a4b12a1bd2ebade87891ab7e36fdbce582301a92"
uuid = "23fbe1c1-3f47-55db-b15f-69d7ec21a316"
version = "0.15.6"

[[LibCURL]]
deps = ["LibCURL_jll", "MozillaCACerts_jll"]
uuid = "b27032c2-a3e7-50c8-80cd-2d36dbcbfd21"

[[LibCURL_jll]]
deps = ["Artifacts", "LibSSH2_jll", "Libdl", "MbedTLS_jll", "Zlib_jll", "nghttp2_jll"]
uuid = "deac9b47-8bc7-5906-a0fe-35ac56dc84c0"

[[LibGit2]]
deps = ["Base64", "NetworkOptions", "Printf", "SHA"]
uuid = "76f85450-5226-5b5a-8eaa-529ad045b433"

[[LibSSH2_jll]]
deps = ["Artifacts", "Libdl", "MbedTLS_jll"]
uuid = "29816b5a-b9ab-546f-933c-edad1886dfa8"

[[Libdl]]
uuid = "8f399da3-3557-5675-b5ff-fb832c97cbdb"

[[Libffi_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "761a393aeccd6aa92ec3515e428c26bf99575b3b"
uuid = "e9f186c6-92d2-5b65-8a66-fee21dc1b490"
version = "3.2.2+0"

[[Libgcrypt_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Libgpg_error_jll", "Pkg"]
git-tree-sha1 = "64613c82a59c120435c067c2b809fc61cf5166ae"
uuid = "d4300ac3-e22c-5743-9152-c294e39db1e4"
version = "1.8.7+0"

[[Libglvnd_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg", "Xorg_libX11_jll", "Xorg_libXext_jll"]
git-tree-sha1 = "7739f837d6447403596a75d19ed01fd08d6f56bf"
uuid = "7e76a0d4-f3c7-5321-8279-8d96eeed0f29"
version = "1.3.0+3"

[[Libgpg_error_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "c333716e46366857753e273ce6a69ee0945a6db9"
uuid = "7add5ba3-2f88-524e-9cd5-f83b8a55f7b8"
version = "1.42.0+0"

[[Libiconv_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "42b62845d70a619f063a7da093d995ec8e15e778"
uuid = "94ce4f54-9a6c-5748-9c1c-f9c7231a4531"
version = "1.16.1+1"

[[Libmount_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "9c30530bf0effd46e15e0fdcf2b8636e78cbbd73"
uuid = "4b2f31a3-9ecc-558c-b454-b3730dcb73e9"
version = "2.35.0+0"

[[Libtiff_jll]]
deps = ["Artifacts", "JLLWrappers", "JpegTurbo_jll", "Libdl", "Pkg", "Zlib_jll", "Zstd_jll"]
git-tree-sha1 = "340e257aada13f95f98ee352d316c3bed37c8ab9"
uuid = "89763e89-9b03-5906-acba-b20f662cd828"
version = "4.3.0+0"

[[Libuuid_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "7f3efec06033682db852f8b3bc3c1d2b0a0ab066"
uuid = "38a345b3-de98-5d2b-a5d3-14cd9215e700"
version = "2.36.0+0"

[[LinearAlgebra]]
deps = ["Libdl"]
uuid = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"

[[Logging]]
uuid = "56ddb016-857b-54e1-b83d-db4d58db5568"

[[MacroTools]]
deps = ["Markdown", "Random"]
git-tree-sha1 = "5a5bc6bf062f0f95e62d0fe0a2d99699fed82dd9"
uuid = "1914dd2f-81c6-5fcd-8719-6d5c9610ff09"
version = "0.5.8"

[[Markdown]]
deps = ["Base64"]
uuid = "d6f4376e-aef5-505a-96c1-9c027394607a"

[[MbedTLS]]
deps = ["Dates", "MbedTLS_jll", "Random", "Sockets"]
git-tree-sha1 = "1c38e51c3d08ef2278062ebceade0e46cefc96fe"
uuid = "739be429-bea8-5141-9913-cc70e7f3736d"
version = "1.0.3"

[[MbedTLS_jll]]
deps = ["Artifacts", "Libdl"]
uuid = "c8ffd9c3-330d-5841-b78e-0817d7145fa1"

[[Measures]]
git-tree-sha1 = "e498ddeee6f9fdb4551ce855a46f54dbd900245f"
uuid = "442fdcdd-2543-5da2-b0f3-8c86c306513e"
version = "0.3.1"

[[Missings]]
deps = ["DataAPI"]
git-tree-sha1 = "bf210ce90b6c9eed32d25dbcae1ebc565df2687f"
uuid = "e1d29d7a-bbdc-5cf2-9ac0-f12de2c33e28"
version = "1.0.2"

[[Mmap]]
uuid = "a63ad114-7e13-5084-954f-fe012c677804"

[[MozillaCACerts_jll]]
uuid = "14a3606d-f60d-562e-9121-12d972cd8159"

[[NaNMath]]
git-tree-sha1 = "bfe47e760d60b82b66b61d2d44128b62e3a369fb"
uuid = "77ba4419-2d1f-58cd-9bb1-8ffee604a2e3"
version = "0.3.5"

[[NamedDims]]
deps = ["AbstractFFTs", "ChainRulesCore", "CovarianceEstimation", "LinearAlgebra", "Pkg", "Requires", "Statistics"]
git-tree-sha1 = "fb4530603a1e62aa5ed7569f283d4b47c2a92f61"
uuid = "356022a1-0364-5f58-8944-0da4b18d706f"
version = "0.2.38"

[[NetworkOptions]]
uuid = "ca575930-c2e3-43a9-ace4-1e988b2c1908"

[[Ogg_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "7937eda4681660b4d6aeeecc2f7e1c81c8ee4e2f"
uuid = "e7412a2a-1a6e-54c0-be00-318e2571c051"
version = "1.3.5+0"

[[OpenSSL_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "15003dcb7d8db3c6c857fda14891a539a8f2705a"
uuid = "458c3c95-2e84-50aa-8efc-19380b2a3a95"
version = "1.1.10+0"

[[Opus_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "51a08fb14ec28da2ec7a927c4337e4332c2a4720"
uuid = "91d4177d-7536-5919-b921-800302f37372"
version = "1.3.2+0"

[[OrderedCollections]]
git-tree-sha1 = "85f8e6578bf1f9ee0d11e7bb1b1456435479d47c"
uuid = "bac558e1-5e72-5ebc-8fee-abe8a469f55d"
version = "1.4.1"

[[PCRE_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "b2a7af664e098055a7529ad1a900ded962bca488"
uuid = "2f80f16e-611a-54ab-bc61-aa92de5b98fc"
version = "8.44.0+0"

[[Parsers]]
deps = ["Dates"]
git-tree-sha1 = "a8709b968a1ea6abc2dc1967cb1db6ac9a00dfb6"
uuid = "69de0a69-1ddd-5017-9359-2bf0b02dc9f0"
version = "2.0.5"

[[Pixman_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "b4f5d02549a10e20780a24fce72bea96b6329e29"
uuid = "30392449-352a-5448-841d-b1acce4e97dc"
version = "0.40.1+0"

[[Pkg]]
deps = ["Artifacts", "Dates", "Downloads", "LibGit2", "Libdl", "Logging", "Markdown", "Printf", "REPL", "Random", "SHA", "Serialization", "TOML", "Tar", "UUIDs", "p7zip_jll"]
uuid = "44cfe95a-1eb2-52ea-b672-e2afdf69b78f"

[[PlotThemes]]
deps = ["PlotUtils", "Requires", "Statistics"]
git-tree-sha1 = "a3a964ce9dc7898193536002a6dd892b1b5a6f1d"
uuid = "ccf2f8ad-2431-5c83-bf29-c5338b663b6a"
version = "2.0.1"

[[PlotUtils]]
deps = ["ColorSchemes", "Colors", "Dates", "Printf", "Random", "Reexport", "Statistics"]
git-tree-sha1 = "2537ed3c0ed5e03896927187f5f2ee6a4ab342db"
uuid = "995b91a9-d308-5afd-9ec6-746e21dbc043"
version = "1.0.14"

[[Plots]]
deps = ["Base64", "Contour", "Dates", "Downloads", "FFMPEG", "FixedPointNumbers", "GR", "GeometryBasics", "JSON", "Latexify", "LinearAlgebra", "Measures", "NaNMath", "PlotThemes", "PlotUtils", "Printf", "REPL", "Random", "RecipesBase", "RecipesPipeline", "Reexport", "Requires", "Scratch", "Showoff", "SparseArrays", "Statistics", "StatsBase", "UUIDs"]
git-tree-sha1 = "cfbd033def161db9494f86c5d18fbf874e09e514"
uuid = "91a5bcdd-55d7-5caf-9e0b-520d859cae80"
version = "1.22.3"

[[PlutoUI]]
deps = ["Base64", "Dates", "HypertextLiteral", "IOCapture", "InteractiveUtils", "JSON", "Logging", "Markdown", "Random", "Reexport", "UUIDs"]
git-tree-sha1 = "d1fb76655a95bf6ea4348d7197b22e889a4375f4"
uuid = "7f904dfe-b85e-4ff6-b463-dae2292396a8"
version = "0.7.14"

[[Preferences]]
deps = ["TOML"]
git-tree-sha1 = "00cfd92944ca9c760982747e9a1d0d5d86ab1e5a"
uuid = "21216c6a-2e73-6563-6e65-726566657250"
version = "1.2.2"

[[Printf]]
deps = ["Unicode"]
uuid = "de0858da-6303-5e67-8744-51eddeeeb8d7"

[[Qt5Base_jll]]
deps = ["Artifacts", "CompilerSupportLibraries_jll", "Fontconfig_jll", "Glib_jll", "JLLWrappers", "Libdl", "Libglvnd_jll", "OpenSSL_jll", "Pkg", "Xorg_libXext_jll", "Xorg_libxcb_jll", "Xorg_xcb_util_image_jll", "Xorg_xcb_util_keysyms_jll", "Xorg_xcb_util_renderutil_jll", "Xorg_xcb_util_wm_jll", "Zlib_jll", "xkbcommon_jll"]
git-tree-sha1 = "ad368663a5e20dbb8d6dc2fddeefe4dae0781ae8"
uuid = "ea2cea3b-5b76-57ae-a6ef-0a8af62496e1"
version = "5.15.3+0"

[[REPL]]
deps = ["InteractiveUtils", "Markdown", "Sockets", "Unicode"]
uuid = "3fa0cd96-eef1-5676-8a61-b3b8758bbffb"

[[Random]]
deps = ["Serialization"]
uuid = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c"

[[RecipesBase]]
git-tree-sha1 = "44a75aa7a527910ee3d1751d1f0e4148698add9e"
uuid = "3cdcf5f2-1ef4-517c-9805-6587b60abb01"
version = "1.1.2"

[[RecipesPipeline]]
deps = ["Dates", "NaNMath", "PlotUtils", "RecipesBase"]
git-tree-sha1 = "7ad0dfa8d03b7bcf8c597f59f5292801730c55b8"
uuid = "01d81517-befc-4cb6-b9ec-a95719d0359c"
version = "0.4.1"

[[Reexport]]
git-tree-sha1 = "45e428421666073eab6f2da5c9d310d99bb12f9b"
uuid = "189a3867-3050-52da-a836-e630ba90ab69"
version = "1.2.2"

[[Requires]]
deps = ["UUIDs"]
git-tree-sha1 = "4036a3bd08ac7e968e27c203d45f5fff15020621"
uuid = "ae029012-a4dd-5104-9daa-d747884805df"
version = "1.1.3"

[[SHA]]
uuid = "ea8e919c-243c-51af-8825-aaa63cd721ce"

[[Scratch]]
deps = ["Dates"]
git-tree-sha1 = "0b4b7f1393cff97c33891da2a0bf69c6ed241fda"
uuid = "6c6a2e73-6563-6170-7368-637461726353"
version = "1.1.0"

[[Serialization]]
uuid = "9e88b42a-f829-5b0c-bbe9-9e923198166b"

[[SharedArrays]]
deps = ["Distributed", "Mmap", "Random", "Serialization"]
uuid = "1a1011a3-84de-559e-8e89-a11a2f7dc383"

[[Showoff]]
deps = ["Dates", "Grisu"]
git-tree-sha1 = "91eddf657aca81df9ae6ceb20b959ae5653ad1de"
uuid = "992d4aef-0814-514b-bc4d-f2e9a6c4116f"
version = "1.0.3"

[[Sockets]]
uuid = "6462fe0b-24de-5631-8697-dd941f90decc"

[[SortingAlgorithms]]
deps = ["DataStructures"]
git-tree-sha1 = "b3363d7460f7d098ca0912c69b082f75625d7508"
uuid = "a2af1166-a08f-5f64-846c-94a0d3cef48c"
version = "1.0.1"

[[SparseArrays]]
deps = ["LinearAlgebra", "Random"]
uuid = "2f01184e-e22b-5df5-ae63-d93ebab69eaf"

[[StaticArrays]]
deps = ["LinearAlgebra", "Random", "Statistics"]
git-tree-sha1 = "3240808c6d463ac46f1c1cd7638375cd22abbccb"
uuid = "90137ffa-7385-5640-81b9-e52037218182"
version = "1.2.12"

[[Statistics]]
deps = ["LinearAlgebra", "SparseArrays"]
uuid = "10745b16-79ce-11e8-11f9-7d13ad32a3b2"

[[StatsAPI]]
git-tree-sha1 = "1958272568dc176a1d881acb797beb909c785510"
uuid = "82ae8749-77ed-4fe6-ae5f-f523153014b0"
version = "1.0.0"

[[StatsBase]]
deps = ["DataAPI", "DataStructures", "LinearAlgebra", "Missings", "Printf", "Random", "SortingAlgorithms", "SparseArrays", "Statistics", "StatsAPI"]
git-tree-sha1 = "8cbbc098554648c84f79a463c9ff0fd277144b6c"
uuid = "2913bbd2-ae8a-5f71-8c99-4fb6c76f3a91"
version = "0.33.10"

[[StructArrays]]
deps = ["Adapt", "DataAPI", "StaticArrays", "Tables"]
git-tree-sha1 = "2ce41e0d042c60ecd131e9fb7154a3bfadbf50d3"
uuid = "09ab397b-f2b6-538f-b94a-2f83cf4a842a"
version = "0.6.3"

[[TOML]]
deps = ["Dates"]
uuid = "fa267f1f-6049-4f14-aa54-33bafae1ed76"

[[TableTraits]]
deps = ["IteratorInterfaceExtensions"]
git-tree-sha1 = "c06b2f539df1c6efa794486abfb6ed2022561a39"
uuid = "3783bdb8-4a98-5b6b-af9a-565f29a5fe9c"
version = "1.0.1"

[[Tables]]
deps = ["DataAPI", "DataValueInterfaces", "IteratorInterfaceExtensions", "LinearAlgebra", "TableTraits", "Test"]
git-tree-sha1 = "1162ce4a6c4b7e31e0e6b14486a6986951c73be9"
uuid = "bd369af6-aec1-5ad0-b16a-f7cc5008161c"
version = "1.5.2"

[[Tar]]
deps = ["ArgTools", "SHA"]
uuid = "a4e569a6-e804-4fa4-b0f3-eef7a1d5b13e"

[[Test]]
deps = ["InteractiveUtils", "Logging", "Random", "Serialization"]
uuid = "8dfed614-e22c-5e08-85e1-65c5234f0b40"

[[URIs]]
git-tree-sha1 = "97bbe755a53fe859669cd907f2d96aee8d2c1355"
uuid = "5c2747f8-b7ea-4ff2-ba2e-563bfd36b1d4"
version = "1.3.0"

[[UUIDs]]
deps = ["Random", "SHA"]
uuid = "cf7118a7-6976-5b1a-9a39-7adc72f591a4"

[[Unicode]]
uuid = "4ec0a83e-493e-50e2-b9ac-8f72acf5a8f5"

[[Wayland_jll]]
deps = ["Artifacts", "Expat_jll", "JLLWrappers", "Libdl", "Libffi_jll", "Pkg", "XML2_jll"]
git-tree-sha1 = "3e61f0b86f90dacb0bc0e73a0c5a83f6a8636e23"
uuid = "a2964d1f-97da-50d4-b82a-358c7fce9d89"
version = "1.19.0+0"

[[Wayland_protocols_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg", "Wayland_jll"]
git-tree-sha1 = "2839f1c1296940218e35df0bbb220f2a79686670"
uuid = "2381bf8a-dfd0-557d-9999-79630e7b1b91"
version = "1.18.0+4"

[[XML2_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Libiconv_jll", "Pkg", "Zlib_jll"]
git-tree-sha1 = "1acf5bdf07aa0907e0a37d3718bb88d4b687b74a"
uuid = "02c8fc9c-b97f-50b9-bbe4-9be30ff0a78a"
version = "2.9.12+0"

[[XSLT_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Libgcrypt_jll", "Libgpg_error_jll", "Libiconv_jll", "Pkg", "XML2_jll", "Zlib_jll"]
git-tree-sha1 = "91844873c4085240b95e795f692c4cec4d805f8a"
uuid = "aed1982a-8fda-507f-9586-7b0439959a61"
version = "1.1.34+0"

[[Xorg_libX11_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg", "Xorg_libxcb_jll", "Xorg_xtrans_jll"]
git-tree-sha1 = "5be649d550f3f4b95308bf0183b82e2582876527"
uuid = "4f6342f7-b3d2-589e-9d20-edeb45f2b2bc"
version = "1.6.9+4"

[[Xorg_libXau_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "4e490d5c960c314f33885790ed410ff3a94ce67e"
uuid = "0c0b7dd1-d40b-584c-a123-a41640f87eec"
version = "1.0.9+4"

[[Xorg_libXcursor_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg", "Xorg_libXfixes_jll", "Xorg_libXrender_jll"]
git-tree-sha1 = "12e0eb3bc634fa2080c1c37fccf56f7c22989afd"
uuid = "935fb764-8cf2-53bf-bb30-45bb1f8bf724"
version = "1.2.0+4"

[[Xorg_libXdmcp_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "4fe47bd2247248125c428978740e18a681372dd4"
uuid = "a3789734-cfe1-5b06-b2d0-1dd0d9d62d05"
version = "1.1.3+4"

[[Xorg_libXext_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg", "Xorg_libX11_jll"]
git-tree-sha1 = "b7c0aa8c376b31e4852b360222848637f481f8c3"
uuid = "1082639a-0dae-5f34-9b06-72781eeb8cb3"
version = "1.3.4+4"

[[Xorg_libXfixes_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg", "Xorg_libX11_jll"]
git-tree-sha1 = "0e0dc7431e7a0587559f9294aeec269471c991a4"
uuid = "d091e8ba-531a-589c-9de9-94069b037ed8"
version = "5.0.3+4"

[[Xorg_libXi_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg", "Xorg_libXext_jll", "Xorg_libXfixes_jll"]
git-tree-sha1 = "89b52bc2160aadc84d707093930ef0bffa641246"
uuid = "a51aa0fd-4e3c-5386-b890-e753decda492"
version = "1.7.10+4"

[[Xorg_libXinerama_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg", "Xorg_libXext_jll"]
git-tree-sha1 = "26be8b1c342929259317d8b9f7b53bf2bb73b123"
uuid = "d1454406-59df-5ea1-beac-c340f2130bc3"
version = "1.1.4+4"

[[Xorg_libXrandr_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg", "Xorg_libXext_jll", "Xorg_libXrender_jll"]
git-tree-sha1 = "34cea83cb726fb58f325887bf0612c6b3fb17631"
uuid = "ec84b674-ba8e-5d96-8ba1-2a689ba10484"
version = "1.5.2+4"

[[Xorg_libXrender_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg", "Xorg_libX11_jll"]
git-tree-sha1 = "19560f30fd49f4d4efbe7002a1037f8c43d43b96"
uuid = "ea2f1a96-1ddc-540d-b46f-429655e07cfa"
version = "0.9.10+4"

[[Xorg_libpthread_stubs_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "6783737e45d3c59a4a4c4091f5f88cdcf0908cbb"
uuid = "14d82f49-176c-5ed1-bb49-ad3f5cbd8c74"
version = "0.1.0+3"

[[Xorg_libxcb_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg", "XSLT_jll", "Xorg_libXau_jll", "Xorg_libXdmcp_jll", "Xorg_libpthread_stubs_jll"]
git-tree-sha1 = "daf17f441228e7a3833846cd048892861cff16d6"
uuid = "c7cfdc94-dc32-55de-ac96-5a1b8d977c5b"
version = "1.13.0+3"

[[Xorg_libxkbfile_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg", "Xorg_libX11_jll"]
git-tree-sha1 = "926af861744212db0eb001d9e40b5d16292080b2"
uuid = "cc61e674-0454-545c-8b26-ed2c68acab7a"
version = "1.1.0+4"

[[Xorg_xcb_util_image_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg", "Xorg_xcb_util_jll"]
git-tree-sha1 = "0fab0a40349ba1cba2c1da699243396ff8e94b97"
uuid = "12413925-8142-5f55-bb0e-6d7ca50bb09b"
version = "0.4.0+1"

[[Xorg_xcb_util_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg", "Xorg_libxcb_jll"]
git-tree-sha1 = "e7fd7b2881fa2eaa72717420894d3938177862d1"
uuid = "2def613f-5ad1-5310-b15b-b15d46f528f5"
version = "0.4.0+1"

[[Xorg_xcb_util_keysyms_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg", "Xorg_xcb_util_jll"]
git-tree-sha1 = "d1151e2c45a544f32441a567d1690e701ec89b00"
uuid = "975044d2-76e6-5fbe-bf08-97ce7c6574c7"
version = "0.4.0+1"

[[Xorg_xcb_util_renderutil_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg", "Xorg_xcb_util_jll"]
git-tree-sha1 = "dfd7a8f38d4613b6a575253b3174dd991ca6183e"
uuid = "0d47668e-0667-5a69-a72c-f761630bfb7e"
version = "0.3.9+1"

[[Xorg_xcb_util_wm_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg", "Xorg_xcb_util_jll"]
git-tree-sha1 = "e78d10aab01a4a154142c5006ed44fd9e8e31b67"
uuid = "c22f9ab0-d5fe-5066-847c-f4bb1cd4e361"
version = "0.4.1+1"

[[Xorg_xkbcomp_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg", "Xorg_libxkbfile_jll"]
git-tree-sha1 = "4bcbf660f6c2e714f87e960a171b119d06ee163b"
uuid = "35661453-b289-5fab-8a00-3d9160c6a3a4"
version = "1.4.2+4"

[[Xorg_xkeyboard_config_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg", "Xorg_xkbcomp_jll"]
git-tree-sha1 = "5c8424f8a67c3f2209646d4425f3d415fee5931d"
uuid = "33bec58e-1273-512f-9401-5d533626f822"
version = "2.27.0+4"

[[Xorg_xtrans_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "79c31e7844f6ecf779705fbc12146eb190b7d845"
uuid = "c5fb5394-a638-5e4d-96e5-b29de1b5cf10"
version = "1.4.0+3"

[[Zlib_jll]]
deps = ["Libdl"]
uuid = "83775a58-1f1d-513f-b197-d71354ab007a"

[[Zstd_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "cc4bf3fdde8b7e3e9fa0351bdeedba1cf3b7f6e6"
uuid = "3161d3a3-bdf6-5164-811a-617609db77b4"
version = "1.5.0+0"

[[libass_jll]]
deps = ["Artifacts", "Bzip2_jll", "FreeType2_jll", "FriBidi_jll", "HarfBuzz_jll", "JLLWrappers", "Libdl", "Pkg", "Zlib_jll"]
git-tree-sha1 = "5982a94fcba20f02f42ace44b9894ee2b140fe47"
uuid = "0ac62f75-1d6f-5e53-bd7c-93b484bb37c0"
version = "0.15.1+0"

[[libfdk_aac_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "daacc84a041563f965be61859a36e17c4e4fcd55"
uuid = "f638f0a6-7fb0-5443-88ba-1cc74229b280"
version = "2.0.2+0"

[[libpng_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg", "Zlib_jll"]
git-tree-sha1 = "94d180a6d2b5e55e447e2d27a29ed04fe79eb30c"
uuid = "b53b4c65-9356-5827-b1ea-8c7a1a84506f"
version = "1.6.38+0"

[[libvorbis_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Ogg_jll", "Pkg"]
git-tree-sha1 = "c45f4e40e7aafe9d086379e5578947ec8b95a8fb"
uuid = "f27f6e37-5d2b-51aa-960f-b287f2bc3b7a"
version = "1.3.7+0"

[[nghttp2_jll]]
deps = ["Artifacts", "Libdl"]
uuid = "8e850ede-7688-5339-a07c-302acd2aaf8d"

[[p7zip_jll]]
deps = ["Artifacts", "Libdl"]
uuid = "3f19e933-33d8-53b3-aaab-bd5110c3b7a0"

[[x264_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "4fea590b89e6ec504593146bf8b988b2c00922b2"
uuid = "1270edf5-f2f9-52d2-97e9-ab00b5d0237a"
version = "2021.5.5+0"

[[x265_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "ee567a171cce03570d77ad3a43e90218e38937a9"
uuid = "dfaa095f-4041-5dcd-9319-2fabd8486b76"
version = "3.5.0+0"

[[xkbcommon_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg", "Wayland_jll", "Wayland_protocols_jll", "Xorg_libxcb_jll", "Xorg_xkeyboard_config_jll"]
git-tree-sha1 = "ece2350174195bb31de1a63bea3a41ae1aa593b6"
uuid = "d8fb68d0-12a3-5cfd-a85a-d49703b185fd"
version = "0.9.1+5"
"""

# ╔═╡ Cell order:
# ╟─21712488-f636-4861-913c-012693cb524d
# ╟─89b58dd0-c70c-4566-b991-2d739990f2ff
# ╟─23c019dd-53d2-461e-a6ba-f1b5da2f6710
# ╠═f5992b42-ff99-479a-b093-70a9bfd85655
# ╠═0a95d032-8576-492f-bee1-ef67303bcb3a
# ╠═7458474c-3c4d-4a99-9640-4bfd2e6465ea
# ╠═edf8bfba-ca6b-4adf-9a7d-e618ca511a80
# ╠═d7f3dc01-1201-41de-8a63-c3e2bc020cf6
# ╠═2e2854f6-8cab-4252-a580-7abc71efd2a0
# ╠═a646342f-3177-4579-8625-cb7d7ea27a49
# ╠═792748ec-199c-41d8-a801-6cd97881ee72
# ╟─85bca31d-c3be-46b6-b330-b80648472378
# ╟─ebd5efea-94b6-419e-8ff1-9dcf624ddaba
# ╟─c8e02ee8-2c21-4a37-9875-c4a463e0b462
# ╟─913d405a-e525-49d6-9021-227543b37953
# ╟─4f87569b-777e-4e97-ac80-8b7b76115d2a
# ╟─a36e5515-e7f9-4c34-a48e-52a004d80f23
# ╟─fc12024b-2785-4f89-a457-95a2deb5045f
# ╟─5afc42f7-3253-438d-96a8-d665a1b9a1cd
# ╠═a5e38ca3-8184-4160-a6b1-4c1783c6e185
# ╠═596a9128-2fff-47e8-8b6f-e9fd8686a279
# ╠═b595349b-b572-4cc8-92b1-07631cc8fad7
# ╠═c578581e-467e-482c-9cce-2ad55d32267f
# ╟─7479e663-e533-442c-8324-c60f4849c470
# ╟─0872a46c-4f4b-49fa-8152-9193b71f332b
# ╠═8c779399-7db0-421e-b1b2-659d28fc8c41
# ╠═e8c989d6-9ba6-40d1-ad36-d768172d5f5e
# ╠═9330d1cf-e0e1-41c1-a876-5fae11633ce8
# ╠═5d450f8d-c8c4-4a58-b2a6-97bfe109f91f
# ╟─5e1377b7-c784-4e42-8190-ffbf48f0ecc1
# ╟─822ea5ab-4a3a-45ee-903c-ff9bad9d5f16
# ╠═04583ab3-2364-4c32-be8e-d3e65ca5321e
# ╠═4fb4c614-ed74-411e-91db-641c331c75fa
# ╠═8cd0438c-c480-4b73-93af-79e3f3a5dfa9
# ╠═36df576e-5111-4d66-9170-2eb9573cea32
# ╠═1e4adc4d-52cf-41df-ad62-cab016f55397
# ╠═2cab05f7-bc59-4d26-a535-495b5ab104c0
# ╠═faf456bc-9145-4d87-ad4d-6ecc8c78bd41
# ╠═e63d17ba-767c-4136-8dc3-3da609f8dda5
# ╠═002504f4-e74a-4e8c-80fa-f4c5f2e90df0
# ╟─4cdab1a4-1784-4430-ae4f-50a4459d83a4
# ╠═5695cdc2-0da6-4ef1-9aa9-d553f9b194e4
# ╠═0c09a3b6-920d-4732-af67-86b36bd88292
# ╟─a558b847-3f09-473c-9ec5-70badcda7589
# ╟─8407fc8c-95b4-4536-bbda-b3f4945518b0
# ╠═104b6690-212b-11ec-3a2b-25c420f62534
# ╠═7f43573f-010f-455a-82c7-f2f509cce3d3
# ╠═d28b8274-7435-484a-bde7-61443d7e4ae2
# ╠═dc95c17d-de0f-4b72-b57e-0dab473a69cd
# ╠═75ebae9f-1c69-4d92-9cdf-a0f4b0971576
# ╟─3140729e-d516-434d-ba32-3ff9705ee83b
# ╟─00000000-0000-0000-0000-000000000001
# ╟─00000000-0000-0000-0000-000000000002


text:
function slide!(cur_player::AbstractPlayer, cur_piece::AbstractPiece, cur_grid::AbstractGrid, cur_direction::Integer)
  can_move(cur_player) || return false

  cur_blocks = map(
    cur_block -> calc_block_coords(cur_block),
    cur_piece.blocks
  )

  did_move = false

  new_direction = 0

  for tmp_direction in sign(cur_direction):sign(cur_direction):cur_direction

    for (cur_index, cur_block) in enumerate(cur_blocks)
      cur_blocks[cur_index] = (
        first(cur_block) + 0,
        last(cur_block) + sign(cur_direction)
      )
    end

    check(cur_blocks, cur_grid) || break

    did_move = true

    new_direction = tmp_direction

  end

  did_move || return false

  cur_piece.j += new_direction

  reset_shadow!(cur_player)

  move!(cur_player)

  true
end

function slide!(cur_player::AbstractPlayer, cur_piece::Nullable{AbstractPiece}, cur_grid::AbstractGrid, cur_direction::Integer)
  false
end

ext: jl


text:
module Wynn
    using SymPy

    # data structure to store the epsilon table
    # ϵ_ij indexed by ints i & j, stored within dict.
    struct EpsilonTable{T}
        series::T
        terms::Vector{T}
        etable::Dict{Tuple{Int,Int},T}
    end

    function EpsilonTable(terms::Vector{T}; simplified::Bool = true) where T<:Union{Float64,Sym}
        # maximum occuring i & j index
        max_ind = length(terms)

        # setting base case j = -1
        etable = Dict((i, -1) => eltype(terms)(0) for i in 0:max_ind)

        # setting base case j = 0, i = 0:max_ind-1
        merge!(etable, Dict((i-1, 0) => sum(terms[1:i]) for i in 1:max_ind))

        # setting base case i = -j-1, j even
        merge!(etable, Dict((-j-1, 2j) => eltype(terms)(0) for j in 0:max_ind))

        # setting base case for odd j
        merge!(etable, Dict((-j-1, 2j-1) => eltype(terms)(0) for j in 0:max_ind))

        # recursive calculations, j>=1, i = floor(Int,-j/2):(max_ind-j-1)
        for j in 1:2*(max_ind-1), i in floor(Int,-j/2):(max_ind-j-1)#j-1
            ϵ_ij = etable[i+1, j-2] + 1 / (etable[i+1, j-1] - etable[i, j-1])
            eltype(terms) <: Sym && simplified ? push!(etable, (i, j) =>
                simplify(ϵ_ij)) : push!(etable, (i, j) => ϵ_ij)
        end
        series = sum(terms)
        EpsilonTable(series, terms, etable)
    end

    EpsilonTable(terms::Vector{Int}; simplified::Bool = true) = EpsilonTable(convert(Vector{Float64}, terms); simplified = simplified)

    # recursive function to generate just one term of the epsilon table, much more efficient.
    function epsilon(terms::Vector{T}, order::Tuple{Int,Int}; simplified::Bool = true) where T<:Union{Real,Sym}
        max_ind = length(terms)
        if order[2] == -1
            return 0
        elseif (0 <= order[1] < max_ind) && (order[2] == 0)
            return sum(terms[1:order[1]+1])
        elseif iseven(order[2]) && (order[1] == -order[2]/2-1)
            return 0
        elseif isodd(order[2]) && (order[1] == -(order[2]+1)/2-1)
            return 0
        end
        max_ind = length(terms)
        ϵ_order = epsilon(terms, (order[1]+1, order[2]-2)) + 1 / (epsilon(terms, (order[1]+1, order[2]-1)) - epsilon(terms, (order[1], order[2]-1)))
        if (eltype(terms) <: Sym) && simplified
            return simplify(ϵ_order)
        else
            return ϵ_order
        end
    end

    export EpsilonTable, epsilon
end

max_stars_repo_name: J-Revell/Wynn


text:
@testset "basic building blocks" begin
    mt = UltimateMortality([0.5,0.5])

    # whole life insurance
    ins = Insurance(
            SingleLife(mort = mt,issue_age = 0),
            Yields.Constant(0.05)
    )

    @test all(survival(SingleLife(mort = mt,issue_age = 0)) .== [1.0, 0.5,0.25])

    @test timepoints(ins) == [1.0,2.0]
    @test survival(ins) == [1.0,0.5]
    @test discount(ins) == [1.0 / 1.05, 1 / 1.05^2]
    @test benefit(ins) == [1.0,1.0]
    @test probability(ins) == [0.5,0.25]
    @test cashflows(ins) == [0.5,0.25]
    @test cashflows(ins) == benefit(ins) .* probability(ins)
    @test present_value(ins)  ≈ 0.5 / 1.05 + 0.5 * 0.5 / 1.05 ^ 2
    @test present_value(ins) ≈ LifeContingencies.A(ins.life,ins.int)

    # basic life contingency tests
    @test survival(LifeContingency(ins),0,1) ≈ 0.5
    @test survival(LifeContingency(ins),0,2) ≈ 0.25
    @test discount(LifeContingency(ins),0,1) ≈ 1/1.05
    @test discount(LifeContingency(ins),0,2) ≈ 1/1.05^2

    # term life insurance
    ins = Insurance(
            SingleLife(mort = mt,issue_age = 0),
            Yields.Constant(0.05),
            n=1
    )

    @test survival(ins) == [1.0]
    @test discount(ins) == [1.0 / 1.05]
    @test benefit(ins) == [1.0]
    @test probability(ins) == [0.5]
    @test cashflows(ins) == [0.5]
    @test cashflows(ins) == benefit(ins) .* probability(ins)
    @test timepoints(ins) == [1.0]
    @test present_value(ins)  ≈ 0.5 / 1.05

    # annuity due
    ins = AnnuityDue(
            SingleLife(mort = mt,issue_age = 0),
            Yields.Constant(0.05)
    )

    @test survival(ins) == [1.0,0.5,.25]
    @test discount(ins) == [1.0, 1 / 1.05^1, 1 / 1.05^2]
    @test benefit(ins) == [1.0,1.0,1.0]
    @test timepoints(ins) == [0.0,1.0,2.0]
    @test probability(ins) == [1.,0.5,0.25]
    @test cashflows(ins) == [1.0,0.5,0.25]
    @test cashflows(ins) == benefit(ins) .* probability(ins)
    @test present_value(ins)  ≈ 1 + 1 * .5 / 1.05 +  1 * .25 / 1.05 ^2

    ins = AnnuityDue(
        SingleLife(mort = mt,issue_age = 0),
        Yields.Constant(0.05),
        n = 2
)
    @test present_value(ins)  ≈ 1 + 1 * .5 / 1.05

    ins = AnnuityImmediate(
        SingleLife(mort = mt,issue_age = 0),
        Yields.Constant(0.05),
        n = 1
    )

    @test present_value(ins)  ≈ 1 * .5 * 1 / 1.05

    ins = AnnuityDue(
        SingleLife(mort = mt,issue_age = 0),
        Yields.Constant(0.05),
        n = 1,
        start_time=1
    )

    @test timepoints(ins)  == [1.0]
    @test probability(ins)  == [0.5]
    @test present_value(ins)  ≈ 1 * .5 * 1 / 1.05

end

@testset "one  year" begin
    ins = LifeContingency(
        SingleLife(
            mort = UltimateMortality([0.5]),
            issue_age = 0
        ),
        Yields.Constant(0.05)
    )

    @test survival(ins,1) ≈ 0.5
    @test survival(ins,0) ≈ 1.0
    @test survival(ins.life,0,0.5) ≈ 1 - 0.5 * 0.5
    @test survival(ins,0.5) ≈ 1 - 0.5 * 0.5

    @test omega(ins) ≈ 1
    @test present_value(AnnuityDue(ins)) ≈ 1 + 1 * .5 / 1.05
    @test present_value(AnnuityDue(ins,n=1)) ≈ 1
    @test present_value(AnnuityDue(ins,n=0)) == 0
    @test present_value(AnnuityImmediate(ins,n=1)) ≈ 1 * .5 / 1.05
    @test present_value(AnnuityImmediate(ins,n=0)) == 0

    @test present_value(Insurance(ins)) ≈ 0.5 / 1.05
    @test present_value(Insurance(ins,n=1)) ≈ 0.5 / 1.05
    @test present_value(Insurance(ins,n=0)) ≈ 0

    ins_jl = LifeContingency(
        JointLife((ins.life,ins.life),LastSurvivor(),Frasier()),
        Yields.Constant(0.05)
    )

    @test omega(ins_jl) ≈ 1
    @test present_value(AnnuityDue(ins_jl)) ≈ 1 + 1 * .75 / 1.05
    @test present_value(AnnuityDue(ins_jl,n=1)) ≈ 1
    @test present_value(AnnuityImmediate(ins_jl,n=1)) ≈ 1 * .75 / 1.05
    @test present_value(AnnuityDue(ins_jl,n=2)) ≈ 1 + 1 * .75 / 1.05
    @test present_value(AnnuityDue(ins_jl,n=2;certain=2)) ≈ 1 + 1 / 1.05
    @test present_value(AnnuityDue(ins_jl,n=0)) == 0

    @test survival(ins_jl,1) ≈ .5 + .5 - .5 * .5
    @test all(survival(JointLife((ins.life,ins.life),LastSurvivor(),Frasier())) .== [1.,.75])
    @test present_value(Insurance(ins_jl)) ≈ .25 / 1.05
    @test present_value(Insurance(ins_jl,n=1)) ≈ 0.25 / 1.05
    @test present_value(Insurance(ins_jl,n=0)) ≈ 0
end

@testset "two year no discount" begin
    ins = LifeContingency(
        SingleLife(
            mort = UltimateMortality([0.0,0.0]),
            issue_age = 0
        ),
        Yields.Constant(0.00)
    )

    @test omega(ins) ≈ 2
    @test present_value(AnnuityDue(ins)) ≈ 3
    @test present_value(AnnuityDue(ins;start_time=1)) ≈ 2
    @test present_value(AnnuityDue(ins,n=1)) ≈ 1
    @test present_value(AnnuityDue(ins,n=2)) ≈ 2
    @test_throws BoundsError present_value(AnnuityDue(ins,n=2;start_time=2))
    @test present_value(AnnuityDue(ins,n=3)) ≈ 3
    @test present_value(AnnuityDue(ins,n=0)) == 0
    @test present_value(AnnuityImmediate(ins,n=0)) == 0
    @test present_value(AnnuityImmediate(ins,n=1)) ≈ 1
    @test present_value(AnnuityImmediate(ins,n=1;start_time=1)) ≈ 1
    @test_throws BoundsError  present_value(AnnuityImmediate(ins;start_time=2))
    @test present_value(AnnuityImmediate(ins,n=2)) ≈ 2
    @test present_value(AnnuityImmediate(ins)) ≈ 2
    @test present_value(AnnuityImmediate(ins;certain=0)) ≈ 2
    @test present_value(AnnuityImmediate(ins;certain=2)) ≈ 2

    @test present_value(Insurance(ins)) ≈ 0
    @test present_value(Insurance(ins,n=1)) ≈ 0
    @test present_value(Insurance(ins,n=0)) ≈ 0

    ins_jl = LifeContingency(
        JointLife(
            lives = (ins.life,ins.life),
            contingency = LastSurvivor(),
            joint_assumption = Frasier()),
        Yields.Constant(0.00)
    )

    @test survival(ins_jl,0) ≈ 1.0
    @test survival(ins_jl,1) ≈ 1
    @test survival(ins_jl,2) ≈ 1

end

@testset "two year with interest" begin
    ins = LifeContingency(
        SingleLife(
            mort = UltimateMortality([0.5,0.5]),
            issue_age = 0
        ),
        Yields.Constant(0.05)
    )

    @test omega(ins) ≈ 2
    @test present_value(AnnuityDue(ins)) ≈ 1 + 1 * .5 * 1 / 1.05 +  1 * .25 / 1.05 ^2
    @test present_value(AnnuityDue(ins,n=1)) ≈ 1
    @test present_value(AnnuityDue(ins,n=2)) ≈ 1 + 1 * .5 * 1 / 1.05
    @test present_value(AnnuityDue(ins,n=3)) ≈ 1 + 1 * .5 * 1 / 1.05 +  1 * .25 / 1.05 ^2
    @test present_value(AnnuityDue(ins,n=0)) ≈ 0
    @test present_value(AnnuityImmediate(ins,n=0)) ≈ 0
    @test present_value(AnnuityImmediate(ins,n=1)) ≈ 1 * .5 * 1 / 1.05

    @test present_value(Insurance(ins)) ≈ 0.5 / 1.05 + 0.5 * 0.5 / 1.05 ^ 2
    @test present_value(Insurance(ins,n=1)) ≈ 0.5 / 1.05
    @test present_value(Insurance(ins,n=0)) ≈ 0

    ins_jl = LifeContingency(
        JointLife(
            lives = (ins.life,ins.life),
            contingency = LastSurvivor(),
            joint_assumption = Frasier()),
        Yields.Constant(0.05)
    )

    @test survival(ins_jl,0) ≈ 1.0
    @test survival(ins_jl,1) ≈ .5 + .5 - .5 * .5
    @test survival(ins_jl,2) ≈ (.25 + .25 - .25 * .25)

end

@testset "two years" begin
    ins = LifeContingency(
        SingleLife(
            mort = UltimateMortality([0.5,0.5]),
            issue_age = 0
        ),
        Yields.Constant(0.05)
    )

    @test omega(ins) ≈ 2
    @test present_value(AnnuityDue(ins,n=1)) ≈ 1
    @test present_value(AnnuityDue(ins,n=2)) ≈ 1 + 1 * .5 * 1 / 1.05
    @test present_value(AnnuityDue(ins,n=3)) ≈ 1 + 1 * .5 * 1 / 1.05 + 1 * .25 * 1 / 1.05 ^ 2
    @test present_value(AnnuityDue(ins,n=0)) == 0
    @test present_value(AnnuityImmediate(ins,n=0)) == 0
    @test present_value(AnnuityImmediate(ins,n=1)) ≈ 1 * .5 * 1 / 1.05

    @test present_value(Insurance(ins)) ≈ 0.5 / 1.05 + 0.5 * 0.5 / 1.05 ^ 2
    @test present_value(Insurance(ins,n=1)) ≈ 0.5 / 1.05
    @test present_value(Insurance(ins,n=0)) ≈ 0

    ins_jl = LifeContingency(
        JointLife(
            lives = (ins.life,ins.life),
            contingency = LastSurvivor(),
            joint_assumption = Frasier()),
        Yields.Constant(0.05)
    )

    @test omega(ins_jl) ≈ 2
    @test present_value(AnnuityDue(ins_jl,n=1)) ≈ 1
    @test present_value(AnnuityDue(ins_jl,n=2)) ≈ 1 + 1 * .75 /1.05
    @test present_value(AnnuityDue(ins_jl,n=3)) ≈ 1 + 1 * .75 /1.05 + 1 * survival(ins_jl,2) / 1.05 ^ 2
    @test present_value(AnnuityDue(ins_jl,n=0)) == 0
    @test present_value(AnnuityImmediate(ins_jl,n=0)) == 0
    @test present_value(AnnuityImmediate(ins_jl,n=1)) ≈ 1 * .75 /1.05

end

# assumes embedded 'testMort' table in combination with MortalityTables
t = UltimateMortality(maleMort)

@testset "demo mortality" begin
    i = Yields.Constant(0.05)
    ins = LifeContingency(SingleLife(mort = t, issue_age = 0), i)

    @test l(ins, 0) ≈ 1.0
    @test l(ins, 1) ≈ 0.993010000000000
    @test l(ins, 2) ≈ 0.992566124530000

    @test D(ins, 0) ≈ 1.0
    @test D(ins, 1) ≈ 0.945723809523809
    @test D(ins, 2) ≈ 0.900286734267574

    @test N(ins, 0) ≈ 20.113017073119200
    @test N(ins, 1) ≈ 19.113017073119200
    @test N(ins, 2) ≈ 18.167293263595400

    @test C(ins, 0) ≈ 0.006657142857142910
    @test C(ins, 1) ≈ 0.000402608136054441
    @test C(ins, 2) ≈ 0.000258082197156658

    @test M(ins, 0) ≈ 0.0422372822324182
    @test M(ins, 1) ≈ 0.0355801393752753
    @test M(ins, 2) ≈ 0.0351775312392208

    @test present_value(Insurance(ins   )) ≈ 0.04223728223
    @test present_value(Insurance(ins, n=0)) ≈ 0.0
    @test present_value(Insurance(ins, n=1)) ≈ 0.0066571428571429
    @test present_value(AnnuityDue(ins, n=0)) ≈ 0.0
    @test present_value(AnnuityDue(ins, n=1)) ≈ 1.0
    @test present_value(AnnuityDue(ins, n=2)) ≈ 1.0 + survival(ins,1) / 1.05
    @test present_value(AnnuityImmediate(ins, n=0)) ≈ 0.0
    @test present_value(AnnuityImmediate(ins, n=1)) ≈ survival(ins,1) / 1.05

    @test present_value(Insurance(ins, n=30)) ≈ 0.0137761089686975

    @test N(ins, 26) ≈ 5.156762988852310
    @test D(ins, 26) ≈ 0.275358702015970
    @test present_value(AnnuityDue(ins, n=26)) ≈ 14.9562540842669

end


text:
# Copyright (c) 2018 Maarten Pronk, MIT license
#
# @references:
# 1. `https://github.com/evetion/GeoArrays.jl/blob/master/src/geoutils.jl`

function get_affine_map(ds::ArchGDAL.Dataset)
    # ArchGDAL fails hard on datasets without
    # an affinemap. GDAL documents that on fail
    # a default affinemap should be returned.
    try
        global gt = ArchGDAL.getgeotransform(ds)
    catch y
        @warn y.msg
        global gt = [0.0, 1.0, 0.0, 0.0, 0.0, 1.0]
    end
    geotransform_to_affine(SVector{6,Float64}(gt))
end

function geotransform_to_affine(gt::SVector{6,Float64})
    # See https://lists.osgeo.org/pipermail/gdal-dev/2011-July/029449.html
    # for an explanation of the geotransform format
    AffineMap(SMatrix{2,2}([gt[2] gt[3]; gt[5] gt[6]]), SVector{2}([gt[1], gt[4]]))
end
geotransform_to_affine(A::Vector{Float64}) = geotransform_to_affine(SVector{6}(A))

function affine_to_geotransform(am::AffineMap)
    l = am.linear
    t = am.translation
    (length(l) == 4 && length(t) == 2) || error("AffineMap has wrong dimensions.")
    [t[1], l[1], l[3], t[2], l[2], l[4]]
end

function unitrange_to_affine(x::AbstractRange, y::AbstractRange)
    δx, δy = step(x), step(y)
    AffineMap(
        SMatrix{2,2}(δx, 0, 0, δy),
        SVector(x[1] - δx / 2, y[1] - δy / 2)
    )
end

"""
    bbox(xmin, ymin, xmax, ymax)
    bbox(;xmin, ymin, xmax, ymax)
    bbox2tuple(b::bbox)
    bbox2vec(b::bbox)
    bbox2affine(size::Tuple{Integer, Integer}, b::bbox)

Spatial bounding box
"""
Base.@kwdef struct bbox
    xmin::Float64
    ymin::Float64
    xmax::Float64
    ymax::Float64
end

function bbox2dims(b::bbox, cellsize)
    lon = b.xmin+cellsize/2:cellsize:b.xmax
    lat = b.ymin+cellsize/2:cellsize:b.ymax
    lon, lat
end

bbox2range(b::bbox) = [b.xmin, b.xmax, b.ymin, b.ymax]

bbox2tuple(b::bbox) = (xmin=b.xmin, ymin=b.ymin, xmax=b.xmax, ymax=b.ymax)
bbox2vec(b::bbox) = [b.xmin, b.ymin, b.xmax, b.ymax]

function bbox2affine(size::Tuple{Integer,Integer}, b::bbox)
    AffineMap(
        SMatrix{2,2}((b.xmax - b.xmin) / size[1], 0, 0, -(b.ymax - b.ymin) / size[2]),
        SVector(b.xmin, b.ymax))
end

const default_affinemap = geotransform_to_affine(SVector(0.0, 1.0, 0.0, 0.0, 0.0, 1.0))
const TYPE_CRS = Union{WellKnownText,AbstractString}

crs2wkt(crs::AbstractString="") = GFT.WellKnownText(GFT.CRS(), crs)
crs2wkt(crs::WellKnownText) = crs

WGS84_str = "GEOGCS[\"WGS 84\",DATUM[\"WGS_1984\",SPHEROID[\"WGS 84\",6378137,298.257223563,AUTHORITY[\"EPSG\",\"7030\"]],AUTHORITY[\"EPSG\",\"6326\"]],PRIMEM[\"Greenwich\",0,AUTHORITY[\"EPSG\",\"8901\"]],UNIT[\"degree\",0.0174532925199433,AUTHORITY[\"EPSG\",\"9122\"]],AXIS[\"Latitude\",NORTH],AXIS[\"Longitude\",EAST],AUTHORITY[\"EPSG\",\"4326\"]]"

WGS84_wkt = crs2wkt(WGS84_str)

export affine_to_geotransform, geotransform_to_affine,
    crs2wkt, bbox,
    bbox2range, bbox2vec, bbox2dims;


text:
#! /usr/bin/julia

# Rosetta Code, Sorting algorithms/Counting sort

function countsort{T<:Integer}(a::Array{T,1})
    (lo, hi) = extrema(a)
    b = zeros(T, length(a))
    cnt = zeros(Int, hi - lo + 1)
    for i in a
        cnt[i - lo + 1] += 1
    end
    z = one(Int)
    for i in lo:hi
        while cnt[i - lo + 1] > 0
            b[z] = i
            z += 1
            cnt[i - lo + 1] -= 1
        end
    end
    return b
end

a = Uint8[rand(1:typemax(Uint8)) for i in 1:20]
println("Sort of Unsigned Bytes:")
println("    Before Sort:")
println("    ", a)
a = countsort(a)
println("\n    After Sort:")
println("    ", a, "\n")

a = [rand(1:2^10) for i in 1:20]
println("Sort of Integers:")
println("    Before Sort:")
println("    ", a)
a = countsort(a)
println("\n    After Sort:")
println("    ", a)

max_stars_repo_name: MichaeLeroy/rosetta-code

max_issues_repo_name: MichaeLeroy/rosetta-code

max_forks_repo_name: MichaeLeroy/rosetta-code


text:
__precompile__(true)
"""
InternedStrings

Copyright 2017-2018 Lyndon White (and other contributors)
Licensed under MIT License, see LICENSE.md
"""
module InternedStrings

export @i_str, intern

Base.@deprecate_binding(InternedString, String, true)

using Base: ht_keyindex2!

########################
# The pool/interning lookup core code

# This forces the type to be inferred (I don't know that the @noinline is reqired or even good)
@noinline getvalue(::Type{K}, wk) where {K} = wk.value::K

# NOTE: This code is carefully optimised.
# Do not tweak it (for readability or otherwise) without benchmarking
@inline function intern!(wkd::WeakKeyDict{K}, key)::K where {K}
    lock(wkd.lock)
    # hand positioning the locks and unlocks
    # (rather than do block or try finally, seems to be faster)
    index = ht_keyindex2!(wkd.ht, key) # returns index if present, or -index if not
    # note hash of weakref is equal to the hash of value, so avoid constructing it if not required
    if index > 0
        # found it
        @inbounds found_key = wkd.ht.keys[index]
        unlock(wkd.lock)
        return getvalue(K, found_key) # return the strong ref
    else
        # Not found, so add it,
        # and mark it as a reference we track to delete!
        kk::K = convert(K, key)
        finalizer(wkd.finalizer, kk) # finalizer is set on the strong ref
        @inbounds Base._setindex!(wkd.ht, nothing, WeakRef(kk), -index)
        unlock(wkd.lock)
        return kk # Return the strong ref
    end
end
#####################################################
# Setup for types

const pool = Dict{DataType, WeakKeyDict}()

@inline function get_pool(::Type{T})::WeakKeyDict{T, Nothing} where {T}
    get!(pool, T) do
        WeakKeyDict{T, Nothing}()
    end
end

###################################

"""
    intern(s::T)

Return a reference to a interned instance of `s`,
adding it to the interning pool if it did not already exist.
"""
function intern(s::T)::T where {T}
    intern(T, s)
end

"""
    intern(::Type{T}, s)

Intern `s` as if it were type `T`, converting it if required.
Note that this will lead to unexpected behavour if the type of `s`, and `T`,
do not have equivalent equality and hash functions
(i.e. this is not safe if `hash(s) != hash(convert(T, s))`).
"""
function intern(::Type{T}, s)::T where {T}
    intern!(get_pool(T), s)
end

"""
Substrings are interned as their parent string type
"""
function intern(substr::SubString{T})::T where {T}
    intern(T, substr)
end

#############################

macro i_str(s)
    # This is done to get interpolation to work correctly
    true_string_expr = esc(Meta.parse(string('"', unescape_string(s), '"')))
    Expr(:call, intern, true_string_expr)
end

end # module InternedStrings


text:
# Functions for Bayesian logit model with a Normal prior N(0, priorVar*I)

function logPrior(pars::Vector{Float64}, nPars::Int, data::Dict{Any, Any})
  return (-dot(pars,pars)/data["priorVar"]
    -nPars*log(2*pi*data["priorVar"]))/2
end

function logLikelihood(pars::Vector{Float64}, nPars::Int, data::Dict{Any, Any})
  XPars = data["X"]*pars
  return (XPars'*data["y"]-sum(log(1+exp(XPars))))[1]
end

function gradLogPosterior(pars::Vector{Float64}, nPars::Int,
  data::Dict{Any, Any})
  return (data["X"]'*(data["y"]-1./(1+exp(-data["X"]*pars)))
    -pars/data["priorVar"])
end

function tensor(pars::Vector{Float64}, nPars::Int, data::Dict{Any, Any})
  p = 1./(1+exp(-data["X"]*pars))
  return ((data["X"]'.*repmat((p.*(1-p))', nPars, 1))*data["X"]
    +(eye(nPars)/data["priorVar"]))
end

function derivTensor(pars::Vector{Float64}, nPars::Int, data::Dict{Any, Any})
  matrix = Array(Float64, data["nData"], nPars)
  output = Array(Float64, nPars, nPars, nPars)

  p = 1./(1+exp(-data["X"]*pars))

  for i = 1:nPars
    for j =1:nPars
      matrix[:, j] = data["X"][:, j].*((p.*(1-p)).*((1-2*p).*data["X"][:, i]))
    end

    output[:, :, i] = matrix'*data["X"]
  end

  return output
end

function randPrior(nPars::Int, data::Dict{Any, Any})
  return rand(Normal(0.0, sqrt(data["priorVar"])), nPars)
end


text:
#FileIO2: Method of categorizing/dispatching on data formats & files.
__precompile__(true)
module FileIO2

import Base: open, read, write, close
import Base: eof
import InteractiveUtils: subtypes

include("codegen.jl")
include("base.jl")
include("textfiles.jl")
include("mime.jl")

const rootpath = realpath(joinpath(dirname(realpath(@__FILE__)),"../."))

export File
export Extensions
export IOOptions, IOOptionsRead, IOOptionsWrite

#Useful file types/aliases:
#TODO: Stop exporting these aliases to avoid namespace pollution.
#TODO: Find a new way to export types (maybe collect them in a different submodule).
export TextFmt, CSVFmt
export HTMLFmt, XMLFmt, MarkdownFmt, AsciiDocFmt
export BMPFmt, PNGFmt, GIFFmt, JPEGFmt, TIFFFmt
export SVGFmt, CGMFmt, EPSFmt, EMFFmt, STLImgFmt

export AbstractDataIO, AbstractReader, AbstractWriter, AbstractDataIORW
#==Create your own reader/writer state machine using the following pattern:
struct MyDataType <: FileIO2.DataFormat; end
mutable struct MyReader <: AbstractReader{MyDataType}
...
end
==#

end #FileIO2


text:
##################################################################
#  html to pdf conversion example
##################################################################

using Wkhtmltox

# Initialize converter
img_init(1)

exfile = joinpath(dirname(@__FILE__),"../examples/example.html")

# Create and populate the settings object
is = ImgSettings("in" => exfile,
                 "out" => outfile,
                 "fmt" => "png")  # png format output

# create the converter object
conv = ImgConverter(is)

# do conversion
run(conv) == 1

# free objects
conv = nothing
img_deinit()


text:
# This file was generated, do not modify it. # hide
Normal() ≪ Lebesgue(ℝ)

max_stars_repo_name: cscherrer/InformativePrior

max_issues_repo_name: cscherrer/InformativePrior

max_forks_repo_name: cscherrer/InformativePrior


text:
"""
    get_cmpt_volumes(mesh)

Get compartment volumes.
"""
function get_cmpt_volumes(mesh)
    map((p, e) -> sum(get_mesh_volumes(p, e)[1]), mesh.points, mesh.elements)
end

"""
    get_mesh_volumes(points, elements)

Get total volume and volumes of each tetrahedron of mesh.
"""
function get_mesh_volumes(points, elements)

    # Sizes
    nelement = size(elements, 2)

    # Elements
    x = reshape(points[:, elements], 3, 4, nelement)

    # Element centers
    centers = mean(x, dims = 2)[:, 1, :]

    # Element volumes
    areavectors =
        [cross(x[:, 2, i] - x[:, 4, i], x[:, 3, i] - x[:, 4, i]) for i = 1:nelement]

    volumes = [1 / 6 * abs(dot(x[:, 1, i] - x[:, 4, i], areavectors[i])) for i = 1:nelement]

    volumes, centers
end


text:
struct Exp <: _Function end
(::Exp)(input::Variable)::Variable = _function(x -> exp.(x), input)


text: println("CHEM 324 Tutorials")


text:
function st_prim_budgeted_lagrangian_megiddo(i::SpanningTreeInstance{T}, budget::Int) where T
  # Initial set of values for λ.
  λlow = 0
  λhigh = Inf

  # TODO: add a fake edge in the graph for the remaining term.

  # TODO: maybe use https://pdfs.semanticscholar.org/a11c/7d19cc2f06fc4d1f8d5b8dcbf8e87e372147.pdf for implementation?
  # Or other references from https://www2.cs.duke.edu/courses/spring07/cps296.2/scribe_notes/p412-agarwal.pdf.

  # Start Prim's algorithm.
  remaining_edges = PriorityQueue{Edge{T}, Float64}() # Easy retrieval of minimum-cost edge.
  node_done = falses(nv(graph(i)))
  node_done[first(vertices(graph(i)))] = true

  # Helper methods.
  edges_around(g, v) = (edgetype(g)(v, x) for x in neighbors(g, v))

  # Initialise with the source node (arbitrarily, the first one).
  first_node = first(vertices(graph(i)))
  for e in edges_around(graph(i), first_node)
    enqueue!(remaining_edges, e => cost(i, e))
  end

  # Build the spanning tree.
  solution = Edge{T}[]
  current_node = first_node
  while length(solution) < nv(graph(i)) - 1 # While there is still an edge to add...
    # Find an admissible edge, starting with those of minimum cost.
    e = nothing
    while true
      e = dequeue!(remaining_edges)
      if (node_done[src(e)] && ! node_done[dst(e)]) || (! node_done[src(e)] && node_done[dst(e)])
        break
      end
    end

    # Update the current node.
    if src(e) == current_node
      current_node = dst(e)
    else
      current_node = src(e)
    end

    # Use this edge in the solution.
    push!(solution, e)
    node_done[current_node] = true

    # Prepare to use the neighbours of this new node.
    for e in edges_around(graph(i), current_node)
      enqueue!(remaining_edges, e => cost(i, e))
    end
  end

  return SpanningTreeSolution(i, solution)
end


text:
export MSet

mutable struct MSet{T} <: AbstractSet{T}
  dict::Dict{T, Int}
  MSet{T}() where {T} = new(Dict{T,Int}())
  MSet{T}(itr) where {T} = union!(new(Dict{T,Int}()), itr)
end

MSet() = MSet{Any}()
MSet(itr) = MSet{eltype(itr)}(itr)

Base.eltype(::Type{MSet{T}}) where {T} = T
Base.similar(s::MSet{T}) where {T} = MSet{T}()
Base.similar(s::MSet, T::Type) = MSet{T}()

#TODO: compact printing, remove trailing , ... the works...
function Base.show(io::IO, s::MSet)
    print(io,"MSet")
    if isempty(s)
        print(io,"{",eltype(s),"}()")
        return
    end
    print(io,"(")
    for (k,v) = s.dict
      if v > 1
        print(io, "$k : $v, ") # ugly
      else
        print(io, "$k, ") # ugly
      end
    end
    print(io,")")
end

Base.isempty(s::MSet) = isempty(s.dict)
Base.length(s::MSet)  = BigInt(sum(values(s.dict)))
Base.in(x, s::MSet) = haskey(s.dict, x)
function Base.push!(s::MSet, x)
  if haskey(s.dict, x)
    s.dict[x] += 1
  else
    s.dict[x] = 1
  end
end

function Base.pop!(s::MSet, x)
  s.dict[x] -= 1
  if s.dict[x] == 0
    delete!(s.dict, x)
  end
  return x
end
Base.pop!(s::MSet, x, deflt) = x in s ? pop!(s, x) : deflt
Base.pop!(s::MSet) = (idx = start(s.dict); val = s.dict.keys[idx]; pop!(s, val))
Base.delete!(s::MSet, x) = (delete!(s.dict, x); s)

Base.copy(s::MSet) = union!(similar(s), s)

Base.start(s::MSet)       = (start(s.dict), 1)
Base.done(s::MSet, state) = done(s.dict, state[1])
function Base.next(s::MSet, state)
  if state[2] < s.dict.vals[state[1]]
    return s.dict.keys[state[1]], (state[1], state[2]+1)
  else
    val, st = next(s.dict, state[1])
    return (val[1], (st, 1))
  end
end

union(s::MSet) = copy(s)
function union(s::MSet, sets::Set...)
    u = MSet{join_eltype(s, sets...)}()
    union!(u,s)
    for t in sets
        union!(u,t)
    end
    return u
end
const ∪ = union
union!(s::MSet, xs) = (for x=xs; push!(s,x); end; s)
union!(s::MSet, xs::AbstractArray) = (for x=xs; push!(s,x); end; s)

function multiplicities(s::MSet)
  return values(s.dict)
end

function Base.unique(s::MSet)
  return collect(keys(s.dict))
end

#... to be completed from base/Set.jl ...


text:
"""
    MultivariateGamma(α,θ)

Construct a distribution that correspond to several comonotonous gammas with common shapes α and respectives scales θ[i].
The distribution can then be used through several methods,
following the Distributions.jl standard, mainly to obtain samples.

The code is type stable and handles any <:Real types, given by the parameters.

Note that the supprot of this distribution is a half-line in the d-dimensional space, [0,+∞], with angle given by θ.

# Examples
```julia-repl
julia> dist = MultivariateGamma(2,[3,4,5]);
julia> sample = zeros(Float64,(3,10));
julia> Random.rand!(dist,sample);
```
"""
struct MultivariateGamma{T<:Real, V<:AbstractVector{T}, M<:Distributions.Gamma{T} } <: Distributions.ContinuousMultivariateDistribution where T
    α::T
    θ::V
    Γα::M
end

#### Constructor :
function MultivariateGamma(α::T1,θ::AbstractVector{T2}) where {T1 <: Real, T2<:Real}
     T = Base.promote_eltype([α],θ,[1.0])
     α = T(α)
     θ = T.(θ)
     MultivariateGamma(α,θ,Distributions.Gamma(α,T(1)))
end
MultivariateGamma(α::T1,θ::T2) where {T1<:Real, T2<:Real} = Distributions.Gamma(promote(α,θ)...)

#### eltype, length, support
Base.eltype(d::MultivariateGamma) = typeof(d.α)
Base.length(d::MultivariateGamma) = length(d.θ)
function Distributions.insupport(d::MultivariateGamma,x::AbstractVector{T}) where {T <: Real}
    # for theta that are 0, x must be 0
    # for other theta, x/theta must be a constant.
    j = 0
    for i in 1:length(x)
        if d.θ[i] == T(0)
            if not x[i] == T(0)
                return false
            end
        else
            j = i
        end
    end
    if j == 0
        # then all theta are 0 and all x are 0
        return true
    end
    constant = x[j]/d.θ[j]
    for i in 1:length(x)
        if d.θ[i] != 0
            # for stability issue, we doubled the tolerence of the type. Otherwise even simulations are not inbounds, which causes pdf == 0 on simulated data !
            if !isapprox(x[i]/d.θ[i],constant,atol=2*eps(T),rtol=T(0))
                return false
            end
        end
    end
    return true
end

#### Sampling
struct MGSPL <: Distributions.Sampleable{Distributions.Multivariate,Distributions.Continuous}
    α
    θ
    Γα
end
Distributions.sampler(d::MultivariateGamma) = MGSPL(d.α,d.θ,d.Γα)
Base.length(s::MGSPL) = length(s.θ)
function Distributions._rand!(rng::Distributions.AbstractRNG, s::MGSPL, x::AbstractVector{T}) where T<:Real
    x .= s.θ .* rand(rng,s.Γα) # simulate a gamma(alpha,1) and multiply by theta
end

#### Logpdf
function Distributions._logpdf(d::MultivariateGamma{T}, x::AbstractArray) where T<:Real
    if !insupport(d,x)
        return T(-Inf)
    end
    for i in 1:length(x)
        if d.θ[i] != 0
            return log(Distributions.pdf(d.Γα,x[i]/d.θ[i])/d.θ[i])
        end
    end
    print("we go there")
    return all(x .== T(0)) ? T(0) : T(-Inf) # if all theta are zero, then the distribution is a diract in (0,...,0)
end

mean(d::MultivariateGamma) = d.α .* d.θ
var(d::MultivariateGamma) = d.α .* d.θ .* d.θ
Statistics.cov(d::MultivariateGamma) = d.α .* d.θ[na,:] .* d.θ[:,na]

# # Okay this should be enough to obtain a pdf
#
# dist = MultivariateGamma(2,[3,4,5])
# sample = zeros(Float64,(3,10))
# import Random
# Random.rand!(dist,sample)
# display(sample)
#
#
# display(Distributions.pdf(dist,sample))
# display(Distributions.pdf(Distributions.Gamma(dist.α,dist.θ[1]),sample[1,:]))


text:
@testset "benchmark_integration" begin
    integrators = [Path(), SplitPath(), Horn(), Frankot()]
    benchmark_integration(integrators, [0, 0.1], SynthSphere(50), debug_mode = true)
    @test true
end


text:
#!/usr/bin/env julia

module PE22

using DelimitedFiles

const pnames = readdlm(joinpath(@__DIR__,"p022_names.txt"),',',String)[1,:];

score(str::AbstractString,index)=index*sum(Int.(collect(str)).-Int('A').+1)

end

export pe22

pe22(data::Vector{<:AbstractString}=PE22.pnames)=sum(PE22.score(val,i) for (i,val) in enumerate(sort(data)))

if !haskey(ENV,"PROJECT_EULER_WITHOUT_RESULT")
    pe22()|>print
end

max_stars_repo_path: 22. Names scores/22.jl

max_stars_repo_name: miRoox/ProjectEuler

max_issues_repo_path: 22. Names scores/22.jl

max_issues_repo_name: miRoox/ProjectEuler

max_forks_repo_path: 22. Names scores/22.jl

max_forks_repo_name: miRoox/ProjectEuler


text:
### A Pluto.jl notebook ###
# v0.12.10

using Markdown
using InteractiveUtils

# ╔═╡ e077ab16-22e5-11eb-15ed-5f0654bf7328
using Pkg, DrWatson

# ╔═╡ e75a9a06-22e5-11eb-2e36-8d4f62b830ed
begin
	@quickactivate "StatisticsWithJuliaPlutoNotebooks"
	using DataFrames, CSV, Dates, StatsBase, Plots, TimeSeries
end

# ╔═╡ ca5fbdc8-22e5-11eb-0a60-cb0a127a3ca5
md"## Listing4.20"

# ╔═╡ fa68607e-22e5-11eb-0558-c9a4d9f77426
begin
	data = CSV.read(datadir("temperatures.csv"),copycols = true)
	brisbane = data.Brisbane
	dates = [Date(
				Year(data.Year[i]),
				Month(data.Month[i]),
				Day(data.Day[i])
			) for i in 1:nrow(data)]

	window1, window2 = 7, 14
	d1 = values(moving(mean,TimeArray(dates,brisbane),window1))
	d2 = values(moving(mean,TimeArray(dates,brisbane),window2))
end;

# ╔═╡ 56bc1334-250c-11eb-3c15-f98b6097c3a6
begin
	grid = (2pi:-2pi/365:0) .+ pi/2
	monthsNames = Dates.monthname.(dates[1:31:365])

	plot(grid, d1[1:366],
		c=:blue, proj=:polar, label="Brisbane weekly average temp.")
	plot!(grid, d2[1:366],
		xticks=([mod.((11pi/6:-pi/6:0) .+ pi/2,2pi) ;], monthsNames),
		c=:red, proj=:polar,
		label="Brisbane fortnightly average temp.", legend=:outerbottom)
end

# ╔═╡ 475ff888-22e6-11eb-2354-09f8f40a8e12
md"## End of listing4.20"

# ╔═╡ Cell order:
# ╟─ca5fbdc8-22e5-11eb-0a60-cb0a127a3ca5
# ╠═e077ab16-22e5-11eb-15ed-5f0654bf7328
# ╠═e75a9a06-22e5-11eb-2e36-8d4f62b830ed
# ╠═fa68607e-22e5-11eb-0558-c9a4d9f77426
# ╠═56bc1334-250c-11eb-3c15-f98b6097c3a6
# ╟─475ff888-22e6-11eb-2354-09f8f40a8e12


text:
# Part 1: programs

import Base: convert

abstract Exp

immutable Sym <: Exp
    name::Symbol
end

immutable Num <: Exp
    val::Int
end

type Call <: Exp
    head::Sym
    args::Vector{Exp}
end

abstract Stmt

type Assign <: Stmt
    lhs::Sym
    rhs::Exp
end

type Goto <: Stmt
    label::Int
end

type GotoIf <: Stmt
    label::Int
    cond::Exp
end

type Ret <: Stmt
end

convert(::Type{Exp}, s::Symbol) = Sym(s)
convert(::Type{Sym}, s::Symbol) = Sym(s)
convert(::Type{Exp}, n::Int)    = Num(n)
convert(::Type{Num}, n::Int)    = Num(n)

# Part 2: lattices

# partial order, meet, join, top, bottom, and their identities

import Base: <=, ==, <, show

abstract LatticeElement

# Note: == and < are defined such that future LatticeElements only
# need to implement <=

<=(x::LatticeElement, y::LatticeElement) = x===y

==(x::LatticeElement, y::LatticeElement) = x<=y && y<=x

<(x::LatticeElement, y::LatticeElement) = x<=y && !(y<=x)

immutable TopElement <: LatticeElement; end
immutable BotElement <: LatticeElement; end

const ⊤ = TopElement()
const ⊥ = BotElement()

show(io::IO, ::TopElement) = print(io, "⊤")
show(io::IO, ::BotElement) = print(io, "⊥")

<=(::BotElement, ::TopElement) = true
<=(::BotElement, ::LatticeElement) = true
<=(::LatticeElement, ::TopElement) = true

# join
⊔(x::LatticeElement, y::LatticeElement) = (x <= y ? y : y <= x ? x : ⊤)

# meet
⊓(x::LatticeElement, y::LatticeElement) = (x <= y ? x : y <= x ? y : ⊥)

# Note: the above definitions are such that we get flat lattices
# "for free" by wrapping any simple julia value in an immutable
# LatticeElement.

# Part 3: dataflow analysis

# Note: the paper uses U+1D56E MATHEMATICAL BOLD FRAKTUR CAPITAL C for this
typealias AbstractValue Dict{Symbol,LatticeElement}

# Here we extend lattices of values to lattices of mappings of variables
# to values. meet and join operate elementwise, and from there we only
# need equality on dictionaries to get <= and <.

⊔(X::AbstractValue, Y::AbstractValue) = [ v => X[v] ⊔ Y[v] for v in keys(X) ]
⊓(X::AbstractValue, Y::AbstractValue) = [ v => X[v] ⊓ Y[v] for v in keys(X) ]

<=(X::AbstractValue, Y::AbstractValue) = X⊓Y == X
< (X::AbstractValue, Y::AbstractValue) = X!=Y && X<=Y

function max_fixed_point(P::Vector, a₁::AbstractValue, eval)
    n = length(P)
    bot = AbstractValue([ v => ⊥ for v in keys(a₁) ])
    s = [ a₁; [ bot for i = 2:n ] ]
    W = IntSet(1)

    while !isempty(W)
        pc = first(W)
        while pc != n+1
            delete!(W, pc)
            I = P[pc]
            new = s[pc]
            if isa(I, Assign)
                # for an assignment, outgoing value is different from incoming
                new = copy(new)
                new[I.lhs.name] = eval(I.rhs, new)
            end
            if isa(I, Goto)
                pc´ = I.label
            else
                pc´ = pc+1
                if isa(I, GotoIf)
                    l = I.label
                    if !(new <= s[l])
                        push!(W, l)
                        s[l] = s[l] ⊔ new
                    end
                end
            end
            if pc´<=n && !(new <= s[pc´])
                s[pc´] = s[pc´] ⊔ new
                pc = pc´
            else
                pc = n+1
            end
        end
    end
    s
end

# example problem - find uses of undefined variables

# flat lattice of variable definedness

immutable IsDefined <: LatticeElement
    is::Bool
end

const undef = IsDefined(false)
const def   = IsDefined(true)

abstract_eval(x::Sym, s::AbstractValue) = get(s, x.name, ⊥)

abstract_eval(x::Num, s::AbstractValue) = def

function abstract_eval(x::Call, s::AbstractValue)
    if any(a->(abstract_eval(a,s) == ⊥), x.args)
        return ⊥
    end
    return def
end

prog1 = [Assign(:x, 0),                       # 1
         GotoIf(5, Call(:randbool, Exp[])),   # 2
         Assign(:y, 1),                       # 3
         Goto(5),                             # 4
         Assign(:z, Call(:pair, Exp[:x,:y])), # 5
         Ret()]

# variables initially undefined
l = AbstractValue(:x => undef, :y => undef, :z => undef)

max_fixed_point(prog1, l, abstract_eval)


text:
using .DataFrames

function anova(df::DataFrame, observationscolumn::Symbol, factorcolumns::Vector{Symbol}, factortypes::Vector{FactorType} = FactorType[]; factornames::Vector{<:AbstractString} = String[])
    observations = df[!, observationscolumn]
    length(observations) > 0 || return
    eltype(observations) <: Number || error("Obervations must be numeric")
    isempty(factornames) && (factornames = [String(col) for col ∈ factorcolumns])
    anova(observations, [df[!, x] for x ∈ factorcolumns], factortypes, factornames = factornames)
end


text:
# Script to replicate figure 4 and 5

using StatsBase, Statistics, LinearAlgebra, StatsPlots, XLSX, PrettyTables, GLM

# -------------- Figure 4 --------------------------------------------------------------------------------------------------------------------------------

t_fig_4 = DataFrame(XLSX.readtable("clean/t_fig_2_3.xlsx", "Sheet1")...) # load directly from file
transform!(t_fig_4, :country .=> ByRow(string) .=> :country, renamecols=false) # change type to String
transform!(t_fig_4, :year .=> ByRow(Int) .=> :year, renamecols=false) # change type to Int64
transform!(t_fig_4, [:FD_GO, :VA_GO, :U, :D] .=> ByRow(Float64) .=> [:FD_GO, :VA_GO, :U, :D], renamecols=false) # change type to Float64

function plot_fig4(df::DataFrame, a, b, year::Int64)

    lims = ifelse(a == :D, (1.4,3.1), (0.28,0.72)) # different x, y axis limits
    df = subset(df, :year => ByRow(x -> x == year))
    # unfortunately unable to change size of annotations
    p = @df df scatter(cols(a), cols(b), xlims=lims, ylims=lims, legend=:none, series_annotations=:country,
    xlabel="$(String(a)) ($(year))", ylabel="$(String(b)) ($(year))") # need to use cols(a) to pass from function
    plot!(0.1:0.1:3.1,0.1:0.1:3.1, col=:red, lw=2) # 45 degree line

    return p
end

l = @layout [a b ; c d] # 2×2 layout
p1 = plot_fig4(t_fig_4, :VA_GO, :FD_GO, 1995)
p2 = plot_fig4(t_fig_4, :VA_GO, :FD_GO, 2011)
p3 = plot_fig4(t_fig_4, :D, :U, 1995)
p4 = plot_fig4(t_fig_4, :D, :U, 2011)
title = "Figure 4: GVC Measures and their Correlation over time"
p_fig4 = plot(p1, p2, p3, p4, layout=l, plot_title=title, plot_titlefontsize=10)
savefig(p_fig4, "images/figure4.png") # export image to folder

# -------------- Figure 5 and 6 --------------------------------------------------------------------------------------------------------------------------------

# data for upper panel
gdf = groupby(t_fig_4, :year)
corr_FD_VA = combine(gdf, [:FD_GO, :VA_GO] => ((x, y) -> cor(x, y)) => :correlation)
corr_U_D = combine(gdf, [:U, :D] => ((x, y) -> cor(x, y)) => :correlation)

# data for lower panel: simple OLS grouped per year
corr_confint_FD_VA = DataFrame(year=Int[], beta=Float64[], lower=Float64[], upper=Float64[]) # initiatlizing DataFrame
corr_confint_U_D = DataFrame(year=Int[], beta=Float64[], lower=Float64[], upper=Float64[])

for i in years

    df = subset(t_fig_4, :year => ByRow(x -> x == i))

    ols = lm(@formula(FD_GO ~ VA_GO), df)
    β = coef(ols)[2]
    conf = confint(ols)[2,:]
    push!(corr_confint_FD_VA, [i β conf[1] conf[2]]) # push data to the dataframe

    ols = lm(@formula(U ~ D), df)
    β = coef(ols)[2]
    conf = confint(ols)[2,:]
    push!(corr_confint_U_D, [i β conf[1] conf[2]]) # push data to the dataframe

end

# figure 5
p1 = @df corr_FD_VA plot(:year, :correlation, xlims=(1994,2012), ylims=(0.79,1.01),
 lw=2, color=:blue, label=:none, ylabel="Correlation: FD on VA")

p2 = @df corr_confint_FD_VA scatter(:year, [:beta, :lower, :upper], xlims=(1994,2012), ylims=(0.68,1.32), label=:none,
 ylabel="Slope coeff and conf int")
plot!(corr_confint_FD_VA.year, corr_confint_FD_VA.beta, color=:blue, lw=2, label=:none)

l = @layout [a ; b] # 2×2 layout
title = "Figure 5: FU/GO and VA/GO over time"
p_fig5 = plot(p1, p2, layout=l, plot_title=title, plot_titlefontsize=10)
savefig(p_fig5, "images/figure5.png") # export image to folder

# figure 6
p1 = @df corr_U_D plot(:year, :correlation, xlims=(1994,2012), ylims=(0.79,1.01),
 lw=2, color=:blue, label=:none, ylabel="Correlation: U on D")

p2 = @df corr_confint_U_D scatter(:year, [:beta, :lower, :upper], xlims=(1994,2012), ylims=(0.68,1.32), label=:none,
 ylabel="Slope coeff and conf int")
plot!(corr_confint_U_D.year, corr_confint_U_D.beta, color=:blue, lw=2, label=:none)

l = @layout [a ; b] # 2×2 layout
title = "Figure 6: U and D over time"
p_fig6 = plot(p1, p2, layout=l, plot_title=title, plot_titlefontsize=10)
savefig(p_fig6, "images/figure6.png") # export image to folder


text:
" Boundary condition expressing `u = value` at a set of nodes "
type DirichletBC <: NodalBC
    " The variable the BC will be applied to "
    u::Variable

    " The boundary ID to apply this BC to "
    bids::Array{Int64}

    " The value on the boundary "
    value::Float64
end

function computeResidual(bc::DirichletBC)
    u = bc.u

    return u.nodal_value - bc.value
end

function computeJacobian(bc::DirichletBC, v::Variable{Float64})
    u = bc.u

    if u.id == v.id
        return 1.
    end

    return 0.
end


text:
using Pkg
Pkg.activate(".")
Pkg.add(PackageSpec(url="https://github.com/AlgebraicJulia/Decapods.jl.git"))
Pkg.instantiate()

max_stars_repo_name: AlgebraicJulia/DECAPODES-Benchmarks

max_issues_repo_name: AlgebraicJulia/DECAPODES-Benchmarks

max_forks_repo_name: AlgebraicJulia/DECAPODES-Benchmarks


text:
function rollout!(p_data::PolicyData, m_data::ModelData; α = 1.0)
    x = m_data.x
    u = m_data.u
    w = m_data.w
    h = m_data.h
    T = m_data.T
    x̄ = m_data.x̄
    ū = m_data.ū
    model = m_data.model

    K = p_data.K
    k = p_data.k

    # initial state
    x[1] = copy(x̄[1])

    # rollout
    for t = 1:T-1
        u[t] .= ū[t] + K[t] * (x[t] - x̄[t]) + α * k[t]
        x[t+1] .= fd(model, x[t], u[t], w[t], h, t)
    end
end

function rollout(model, x1, u, w, h, T)
    x_hist = [x1]

    for t = 1:T-1
        push!(x_hist, fd(model, x_hist[end], u[t], w[t], h, t))
    end

    return x_hist
end


text:
# ------------------------------------------------------------------------------
#
#  Gmsh Julia tutorial 3
#
#  Extruded meshes, options
#
# ------------------------------------------------------------------------------
import gmsh

gmsh.initialize()

gmsh.model.add("t3")

# Copied from t1.jl...
lc = 1e-2
gmsh.model.geo.addPoint(0, 0, 0, lc, 1)
gmsh.model.geo.addPoint(.1, 0,  0, lc, 2)
gmsh.model.geo.addPoint(.1, .3, 0, lc, 3)
gmsh.model.geo.addPoint(0, .3, 0, lc, 4)
gmsh.model.geo.addLine(1, 2, 1)
gmsh.model.geo.addLine(3, 2, 2)
gmsh.model.geo.addLine(3, 4, 3)
gmsh.model.geo.addLine(4, 1, 4)
gmsh.model.geo.addCurveLoop([4, 1, -2, 3], 1)
gmsh.model.geo.addPlaneSurface([1], 1)
gmsh.model.geo.synchronize()
gmsh.model.addPhysicalGroup(0, [1, 2], 1)
gmsh.model.addPhysicalGroup(1, [1, 2], 2)
gmsh.model.addPhysicalGroup(2, [1], 6)
gmsh.model.setPhysicalName(2, 6, "My surface")

# As in `t2.jl', we plan to perform an extrusion along the z axis.  But here,
# instead of only extruding the geometry, we also want to extrude the 2D
# mesh. This is done with the same `extrude()' function, but by specifying
# element 'Layers' (2 layers in this case, the first one with 8 subdivisions and
# the second one with 2 subdivisions, both with a height of h/2). The number of
# elements for each layer and the (end) height of each layer are specified in
# two vectors:
h = 0.1
angle = 90.
ov = gmsh.model.geo.extrude([(2,1)], 0, 0, h, [8,2], [0.5,1])

# The extrusion can also be performed with a rotation instead of a translation,
# and the resulting mesh can be recombined into prisms (we use only one layer
# here, with 7 subdivisions). All rotations are specified by an an axis point
# (-0.1, 0, 0.1), an axis direction (0, 1, 0), and a rotation angle (-Pi/2):
ov = gmsh.model.geo.revolve([(2,28)], -0.1,0,0.1, 0,1,0, -pi/2, [7])

# Using the built-in geometry kernel, only rotations with angles < Pi are
# supported. To do a full turn, you will thus need to apply at least 3
# rotations. The OpenCASCADE geometry kernel does not have this limitation.

# A translation (-2*h, 0, 0) and a rotation ((0,0.15,0.25), (1,0,0), Pi/2) can
# also be combined to form a "twist".  The last (optional) argument for the
# extrude() and twist() functions specifies whether the extruded mesh should be
# recombined or not.
ov = gmsh.model.geo.twist([(2,50)], 0,0.15,0.25, -2*h,0,0, 1,0,0,
                          angle*pi/180., [10], [], true)

gmsh.model.geo.synchronize()

# All the extrusion functions return a vector of extruded entities: the "top" of
# the extruded surface (in `ov[0]'), the newly created volume (in `ov[1]') and
# the tags of the lateral surfaces (in `ov[2]', `ov[3]', ...).

# We can then define a new physical volume (with tag 101) to group all the
# elementary volumes:
gmsh.model.addPhysicalGroup(3, [1, 2, ov[2][2]], 101)

gmsh.model.mesh.generate(3)
gmsh.write("t3.msh")

gmsh.finalize()


text:
using CommonCrawl

@assert isdefined(:CommonCrawl)
@assert typeof(CommonCrawl) === Module