text:
# final scheme. adaptive kernel.

import Printf
import PyPlot
import Random
import Optim

using LinearAlgebra

using FFTW

import Statistics

import Distributions
#import HCubature
import Interpolations
import SpecialFunctions

# import SignalTools
#
# import Utilities # https://gitlab.com/RoyCCWang/utilities
#
# import VisualizationTools # https://gitlab.com/RoyCCWang/visualizationtools

#import Calculus
#import ForwardDiff

import Utilities # https://gitlab.com/RoyCCWang/utilities
import ProbabilityUtilities # https://gitlab.com/RoyCCWang/probabilityutilities
import RiemannianOptim # https://gitlab.com/RoyCCWang/riemannianoptim

include("../src/misc/declarations.jl")

include("../src/RKHS/RKHS.jl")
include("../src/RKHS/kernel.jl")
include("../src/RKHS/interpolators.jl")

include("../src/warp_map/Rieszanalysis.jl")
include("../src/misc/utilities.jl")

include("../src/optimization/sdp.jl")

include("../src/misc/front_end.jl")

PyPlot.close("all")
fig_num = 1

Random.seed!(25)

function visualizemeshgridpcolor(   x_ranges::Vector{LinRange{T}},
                                    Y::Matrix{T},
                                    marker_locations::Vector,
                                    marker_symbol::String,
                                    fig_num::Int,
                                    title_string::String;
                                    x1_title_string::String = "Dimension 1",
                                    x2_title_string::String = "Dimension 2",
                                    cmap = "Greens_r") where T <: Real
    #
    @assert length(x_ranges) == 2
    x_coords = collect( collect(x_ranges[d]) for d = 1:2 )

    PyPlot.figure(fig_num)
    fig_num += 1
    PyPlot.pcolormesh(x_coords[2], x_coords[1], Y, cmap = cmap)
    PyPlot.xlabel(x2_title_string)
    PyPlot.ylabel(x1_title_string)
    PyPlot.title(title_string)

    for i = 1:length(marker_locations)
        pt = reverse(marker_locations[i])
        PyPlot.annotate(marker_symbol, xy=pt, xycoords="data")
    end

    PyPlot.plt.colorbar()

    return fig_num
end

### I am here. write test script for fitnDdensityRiemannian. Then, move onto finite KR, helmt.

D = 2

N_components = 3
N_realizations = 100
N_array = [200; 200]

### set up coordinates.
τ = 15.0
limit_a = [-τ; -τ]
limit_b = [τ; τ]

x_ranges = collect( LinRange(limit_a[d], limit_b[d], N_array[d]) for d = 1:D )
X_nD = Utilities.ranges2collection(x_ranges, Val(D))
X = vec(X_nD)

### set up oracle density.
mixture_weights = rand(N_components)
mixture_weights = mixture_weights ./ sum(mixture_weights)

f, oracle_dist = ProbabilityUtilities.generaterandomGMM(mixture_weights, D)

# visualize full joint density.
f_X = f.(X)
f_X_nD = f.(X_nD)

fig_num = visualizemeshgridpcolor(x_ranges, f_X_nD, [], ".", fig_num, "oracle density, full joint, D = 2")

##### RKHS for density

zero_tol = 1e-13
max_iters = 500
σ² = 1e-3

optim_max_iters = 50000
show_trace_flag = false

########## second dimension.
println("Start working on d_select = 2")
d_select = 2

#######

𝑋 = collect( rand(oracle_dist) for n = 1:N_realizations)
𝑌 = f.(𝑋)

### fit density.
𝑌_nD = f.(X_nD)
σ² = 1e-5

# # warp map.
# amplification_factor = 1.0
# attenuation_factor_at_cut_off = 2
# N_bands = 5
# reciprocal_cut_off_percentages = ones(N_bands) ./collect(LinRange(1.0,0.2,N_bands))
# ω_set = collect( π/(reciprocal_cut_off_percentages[i]*sqrt(2*log(attenuation_factor_at_cut_off))) for i = 1:length(reciprocal_cut_off_percentages) )
# pass_band_factor = abs(ω_set[1]-ω_set[2])*0.2
#
# ϕ = getRieszwarpmapsamples(𝑌_nD, Val(:simple), Val(:uniform), ω_set, pass_band_factor)
#
# ϕ_map_func_2D, d_ϕ_map_func_2D, d2_ϕ_map_func_2D = Utilities.setupcubicitp(ϕ, x_ranges[1:d_select], amplification_factor)

θ_canonical = RationalQuadraticKernelType(0.7)
# θ_a_2D = AdaptiveKernelType(θ_canonical, ϕ_map_func_2D)
#
# println("Starting density fit: canonical kernel")
# @time η, α_SDP = fitnDdensity(vec(𝑌), vec(𝑋), σ², θ_canonical, zero_tol, max_iters)
# println("finished density fit")

# ## save data.
# import BSON
#
# y = vec(𝑌)
# X = vec(𝑋)
# println("construct matrix")
# @time K = constructkernelmatrix(X, θ_canonical)
# μ = σ²
# BSON.bson("density_fit_data1.bson", K = K,
#                                     y = y,
#                                     X = X,
#                                     μ = μ,
#                                     α_SDP = α_SDP)
# @assert 1==2

# ζ, α_optim, final_results,
#   initial_results_minimizer, ζ_costfunc = fitdensityoptim( vec(𝑌), vec(𝑋), σ², θ_canonical;
#                                         zero_tol = zero_tol,
#                                         max_iters = optim_max_iters,
#                                         show_trace = show_trace_flag )

#optim_max_iters = 5000
println("Starting density fit: canonical kernel, Riemannian")

# lower, faster, less accurate.
avg_Δf_tol = 1e-5 # 1e=10 gives better score than α_SDP.

@time ζ, α_optim = fitnDdensityRiemannian( vec(𝑌), vec(𝑋), σ²,
                                          θ_canonical, zero_tol, optim_max_iters;
                                          avg_Δf_tol = avg_Δf_tol )
println("finished density fit")
#
# discrepancy = norm(α_SDP - α_optim)
# println("SDP vs. Optim l-2 discrepancy: ", discrepancy)
#
K_costfunc = constructkernelmatrix(vec(𝑋), θ_canonical)
costfunc = aa->RiemannianOptim.RKHSfitdensitycostfunc(aa, K_costfunc, vec(𝑌), σ²)
#
# println("θ_canonical: costfunc(α_SDP) ", costfunc(α_SDP) )
println("θ_canonical: costfunc(α_optim) ", costfunc(α_optim) )
println()

# println("[α_SDP α_optim] = ")
# display([α_SDP α_optim])
# println()

# I am here. save another dataset, and tune the zero_tol of the retraction.
# finish up with swarm optim.

# TODO future works, or if insuccificent material for publication, fit conditional densities.

@assert 1==2

println("Starting density fit: adaptive kernel")
σ² = 1e-5
@time η_a, α_a_SDP = fitnDdensity(vec(𝑌), vec(𝑋), σ², θ_a_2D, zero_tol, max_iters)
println("finished density fit")

ζ_a, α_a_optim, final_results_a,
  initial_results_minimizer_a, ζ_costfunc_a = fitdensityoptim( vec(𝑌), vec(𝑋), σ², θ_a_2D;
                                            zero_tol = zero_tol,
                                            max_iters = optim_max_iters,
                                            show_trace = show_trace_flag )

discrepancy = norm(α_SDP - α_optim)
println("SDP vs. Optim l-2 discrepancy: ", discrepancy)
println("θ_canonical: costfunc(α_a_SDP) ", ζ_costfunc_a(α_a_SDP) )
println("θ_canonical: costfunc(α_a_optim) ", ζ_costfunc_a(α_a_optim) )
println()

# next, determined the kernel centers.

fq = xx->sum( η.c[n]*evalkernel(xx,η.X[n],η.θ) for n = 1:length(η.c) )

#
fq_a = xx->sum( η_a.c[n]*evalkernel(xx, η_a.X[n], η_a.θ) for n = 1:length(η_a.c) )

### visualize. Show adaptive kernel is useful.
# replace [] with 𝑋 to see kernel centers.
𝑌q = fq.(X_nD)
𝑌q_a = fq_a.(X_nD)

fig_num = visualizemeshgridpcolor(x_ranges[1:d_select], 𝑌q, [], "x",
                fig_num, "canonical query", "x1", "x2")

fig_num = visualizemeshgridpcolor(x_ranges[1:d_select], 𝑌q_a, [], "x",
                fig_num, "adaptive query", "x1", "x2")

###

gq = xx->sum( ζ.c[n]*evalkernel(xx,ζ.X[n],ζ.θ) for n = 1:length(ζ.c) )

#
gq_a = xx->sum( ζ_a.c[n]*evalkernel(xx, ζ_a.X[n], ζ_a.θ) for n = 1:length(ζ_a.c) )

### visualize. Show adaptive kernel is useful.
# replace [] with 𝑋 to see kernel centers.
gq_X = gq.(X_nD)
gq_a_X = gq_a.(X_nD)

fig_num = visualizemeshgridpcolor(x_ranges[1:d_select], gq_X,
                [], "x", fig_num, "canonical query: gq", "x1", "x2")

fig_num = visualizemeshgridpcolor(x_ranges[1:d_select], gq_a_X,
                [], "x", fig_num, "adaptive query: gq", "x1", "x2")

max_stars_repo_name: RoyCCWang/RKHSRegularization

max_issues_repo_name: RoyCCWang/RKHSRegularization

max_forks_repo_name: RoyCCWang/RKHSRegularization


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

using Markdown
using InteractiveUtils

# ╔═╡ de9d25b0-9ecc-4752-8c4d-e452dcf82c57
begin
	# Packages
	using PlutoUI, LightGraphs, GraphPlot, Clustering, SparseArrays
	using LinearAlgebra, Arpack, Plots, Random
end

# ╔═╡ 0f6cd296-6cb8-4a9d-a21f-02b689dfe4e9
TableOfContents(title="📚 Table of Contents", aside=true)

# ╔═╡ 89299cb7-a4eb-48f3-b3f1-ab2b5c6cfa2e
md"""
# Spectral Partitioning of Bipartite Graphs

Typical example of bipartite graph is a graph obtained from a collection of documents presented as a _term $\times$ document_  matrix.

__Prerequisites__

The reader should be familiar with k-means algorithm and spectral graph partitioning theory and algorithms.

__Competences__

The reader should be able to apply spectral partitioning of bipartite graphs to data clustering problems.

__Credits.__ The notebook was initially derived from M.Sc. Thesis of Ivančica Mirošević.
"""

# ╔═╡ a92d0e70-479d-444e-b06d-e47d2e385b10
md"""
## Definitions

__Undirected bipartite graph__ $G$ is a triplet $G=(T,D,E)$, where $T=\{t_{1},\cdots ,t_{m}\}$ and $D=\{d_{1},...,d_{n}\}$ are two sets of vertices and $E=\{(t_{i},d_{j}):t_{i}\in R,d_{j}\in D\}$, is a set of edges.

 $G$ is __weighted__ if there is weight $\omega(e)$ associated with each edge $e\in E$.

For example, $D$ is a set of documents, $T$ is a set of terms (words) and  edge $e=(t_{i},d_{j})$ exists if document $d_{j}$ contains term $t_{i}$. Weight $\omega(e)$ can be number of appearances of the term $t_i$ in the document $d_j$.

A __term-by-document-matrix__ is a matrix $A\in\mathbb{R}^{m\times n}$ with $A_{ij}=\omega((t_i,d_j))$.
"""

# ╔═╡ fddcc2a7-f718-4900-a135-764eb36496a7
md"""
## Facts

1. The weight matrix of $G$ is $W=\begin{bmatrix}0 & A \\ A^{T} & 0 \end{bmatrix}$.

2. The Laplacian matrix of $G$ is

$$L=\begin{bmatrix} \Delta_{1} & -A \\ -A^{T} & \Delta_{2}\end{bmatrix},$$

where $\Delta_1$ and $\Delta_2$ are diagonal matrices with elements

$\Delta_{1,ii}=\sum\limits_{j=1}^n A_{ij},\quad i=1,\ldots,m,$

$\Delta_{2,jj}=\sum\limits_{i=1}^m A_{ij},\quad j=1,\ldots,n.$

3. The normalized Laplacian matrix of $G$ is

$$L_n=\begin{bmatrix}
I & -\Delta_{1}^{-\frac{1}{2}}A\Delta_{2}^{-\frac{1}{2}} \\
-\Delta_{2}^{-\frac{1}{2}}A^T\Delta_{1}^{-\frac{1}{2}} & I
\end{bmatrix} \equiv
\begin{bmatrix} I & -A_n \\ -A_n^T & I \end{bmatrix}.$$
"""

# ╔═╡ b3b44560-203a-11eb-2782-170e4c12627a
md"""
4. Let $\lambda$ be an eigenvalue of $L_n$ with an eigenvector $w=\begin{bmatrix} u \\ v\end{bmatrix}$, where $u\in \mathbb{R}^{m}$ $v\in\mathbb{R}^{n}$. Then $L_n w=\lambda w$ implies $A_n v =(1-\lambda)u$ and $A_n^T u=(1-\lambda)v$. Vice versa, if $(u,\sigma,v)$ is a singular triplet of $A_n$, then $1-\sigma$ is an eigenvalue of $L_n$ with (non-unit) eigenvector $w=\begin{bmatrix} u \\ v\end{bmatrix}$.

5. The second largest singular value of $A_n$ corresponds to the second smallest eigenvalue of $L_n$, and computing the former is numerically more stable.
"""

# ╔═╡ 1be90f30-203b-11eb-34da-c3116523ec3a
md"

# Algorithms

## Bipartitioning algorithm

1. For given $A$ compute $A_{n}$.
2. Compute singular vectors of $A_{n}$, $u^{[2]}$ and $v^{[2]}$, which correspond to the second largest singular value, $\sigma_2(A_n)$.
3. Assign the partitions $T=\{T_1,T_2\}$ and $D=\{D_1,D_2\}$ according to the signs of $u^{[2]}$ and $v^{[2]}$. The pair $(T,D)$ is now partitioned as $\{(T_1,D_1),(T_2,D_2)\}$.
"

# ╔═╡ 5211ce2e-203b-11eb-0b6f-55e809177fb9
md"
## Recursive bipartitioning algorithm

1. Compute the bipartition $\pi=\{(T_1,D_1),(T_2,D_2)\}$ of $(T,D)$. Set the counter $c=2$.

2. While $c<k$ repeat

   - compute bipartitions of each of the subpartitions of $(T,D)$,
   - among all $(c+1)$-subpartitions, choose the one with the smallest $\mathop{\mathrm{pcut}}(\pi_{c+1})$ or $\mathop{\mathrm{ncut}}(\pi_{c+1})$, respectively.
   - Set $c=c+1$

3. Stop
"

# ╔═╡ 38333a20-203c-11eb-2b28-59e07b3e50ee
md"
## Multipartitioning algorithm
1. For given $A$ compute $A_{n}$.
2. Compute $k$ left and right singular vectors, $u^{[1]},\ldots,u^{[k]}$ and $v^{[1]},\ldots,v^{[k]}$, which correspond to $k$ largest singular values $\sigma_1\geq \cdots \geq \sigma_k$ of $A_n$.
3. Partition the rows of matrices $\Delta_{1}^{-\frac{1}{2}}\begin{bmatrix} u^{[1]} & \ldots & u^{[k]}\end{bmatrix}$ and $\Delta_{2}^{-\frac{1}{2}}\begin{bmatrix} v^{[1]} & \ldots & v^{[k]}\end{bmatrix}$ with the k-means algorithm.
"

# ╔═╡ 1bc20329-3a6e-435d-a25c-8cce79fb7d9d
md"
# Examples

## Small term-by- document matrix
"

# ╔═╡ b9d0830a-0a33-40c1-a08a-49538cf7a879
begin
	# Make a nicer spy function
	import Plots.spy
	spy(A)=heatmap(A,yflip=true,color=:RdBu,aspectratio=1,clim=(-1,1.0))
end

# ╔═╡ f6492dfe-faac-407f-9a0e-5dcf321a4f2d
begin
	# Sources, targets, and weights
	dn=[6,6,7,6,7,7]
	tn=[1,2,2,3,4,5]
	wn=[3,1,3,2,2,3]
	[dn tn wn]
end

# ╔═╡ 1a8b6f7b-4098-4cd2-adad-8ee788b47333
mynames=["Term 1";"Term 2";"Term 3";"Term 4";"Term 5";"Doc 1";"Doc 2"]

# ╔═╡ 4079aa8e-2bb5-4d51-bca8-04d0997db2bf
begin
	G=Graph(7)
	for i=1:length(dn)
	    add_edge!(G,tn[i],dn[i])
	end
	gplot(G, nodelabel=mynames, edgelabel=wn)
end

# ╔═╡ b6d41ab1-24f2-417d-b604-3e12384e65f8
begin
	# We define some functions
	function WeightMatrix(src::Array,dst::Array,weights::Array)
	    n=nv(G)
	    sparse([src;dst],[dst;src],[weights;weights],n,n)
	end

	Laplacian(W::AbstractMatrix)=spdiagm(0=>vec(sum(W,dims=2)))-W

	function NormalizedLaplacian(L::AbstractMatrix)
	    D=1.0./sqrt.(diag(L))
	    n=length(D)
	    [L[i,j]*(D[i]*D[j]) for i=1:n, j=1:n]
	end
end

# ╔═╡ d5d8ee5d-9cce-45fc-a31b-b2ed8237cc86
W=WeightMatrix(tn,dn,wn)

# ╔═╡ 86fd4f5a-af57-490d-a17f-2339a9414057
Matrix(W)

# ╔═╡ 1f53b5cf-ab6a-4107-b53b-9ce96afcdb7f
begin
	L=Laplacian(W)
	Matrix(L)
end

# ╔═╡ 63403411-cca4-4d71-8a76-8dea8638f69f
Lₙ=NormalizedLaplacian(L)

# ╔═╡ 3cf471bc-d316-4de9-af63-71db9b035c40
begin
	A=W[1:5,6:7]
	Δ₁=sqrt.(sum(A,dims=2))
	Δ₂=sqrt.(sum(A,dims=1))
	Aₙ=[A[i,j]/(Δ₁[i]*Δ₂[j]) for i=1:size(A,1), j=1:size(A,2)]
end

# ╔═╡ ac0b5eda-299c-4ae2-ae78-bda92b378f19
# The partitioning - explain the results!
U,σ,V=svd(Aₙ)

# ╔═╡ 9bf6bbca-f96c-4b38-a4fd-5188ec29a2a2
U[:,2]

# ╔═╡ aa51688c-21a0-40f4-a6d6-01d441bfd4b7
V[:,2]

# ╔═╡ 5f27f1a2-80de-41d1-b1e2-da3a99493217
md"
## Sets of points
"

# ╔═╡ 7cd3af86-f459-4e4b-8e35-c29bc6966eed
begin
	# Define sizes
	m=[200,100,100]
	n=[100,200,100]
	density=[0.5,0.7,0.4]
	Pts=Array{Any}(undef,3)
	Random.seed!(421)
	for i=1:3
	    # Generate sparse random, not with sprandn
	    Pts[i]=sprand(m[i],n[i],density[i])
	end
	B=blockdiag(Pts[1],Pts[2],Pts[3])
end

# ╔═╡ e1149b05-d594-4d75-ae56-1ce18bd0f60d
spy(Matrix(B))

# ╔═╡ b95d4070-c127-4db0-987b-5fc17421b667
# The structure of singular vectors reflects the block
S,rest₀=svds(B,nsv=3);

# ╔═╡ 57d472c1-548d-4678-b884-85f863771aa7
# S is a structure
S.S

# ╔═╡ 4806378b-63cb-445a-8c71-583458fce846
begin
	# Plot the first three left singular vectors
	k=size(B,1)
	x=collect(1:k)
	scatter(x,S.U[:,1],title="Left Singular Vectors",label="U[:,1]")
	scatter!(x,S.U[:,2],label="U[:,2]",legend=:topleft)
	scatter!(x,S.U[:,3],label="U[:,3]")
end

# ╔═╡ a08c99a3-9c10-47be-9dc6-987d92a85087
begin
	# Plot the first three right singular vectors
	scatter(x,S.Vt[1,:],title="Right Singular Vectors",label="V[:,1]")
	scatter!(x,S.Vt[2,:],label="V[:,2]")
	scatter!(x,S.Vt[3,:],label="V[:,3]")
end

# ╔═╡ 6888e8a9-f306-4490-9147-2460fbda4889
begin
	# Add random noise
	noise=sprand(k,k,0.3)
	C=B+noise
	spy(Matrix(C))
end

# ╔═╡ 3bd5ad15-31cd-4efd-b4e3-7cf7ede3edae
begin
	# Apply random permutation to rows and columns of C
	D=C[randperm(k),randperm(k)]
	spy(Matrix(D))
end

# ╔═╡ e4487f35-5de0-492e-96f1-69318d7d8e8b
md"""
__Question.__ Given D, can we recover C?

__Answer.__ Yes (with spectral partitioning)!
"""

# ╔═╡ 81d99f24-82b9-4216-8d17-6a18317beb77
Sₙ,rest=svds(D,nsv=3);

# ╔═╡ 8f925ff7-90b4-4fe4-a869-86dec693c0f4
# K-means on rows of U
outU=kmeans(Matrix(transpose(Sₙ.U)),3)

# ╔═╡ d446724a-97d1-403f-baaf-d0619203f351
# K-means on Vt
outV=kmeans(Sₙ.Vt,3)

# ╔═╡ 7572c241-dbcc-4094-97e3-ee3172255294
sortperm(outU.assignments)

# ╔═╡ 22e90c11-ba4e-4073-9f78-f3b7beb24e69
begin
	# RECOVERY of B
	E=D[sortperm(outU.assignments),sortperm(outV.assignments)]
	spy(Matrix(E))
end

# ╔═╡ 00000000-0000-0000-0000-000000000001
PLUTO_PROJECT_TOML_CONTENTS = """
[deps]
Arpack = "7d9fca2a-8960-54d3-9f78-7d1dccf2cb97"
Clustering = "aaaa29a8-35af-508c-8bc3-b662a17a0fe5"
GraphPlot = "a2cc645c-3eea-5389-862e-a155d0052231"
LightGraphs = "093fc24a-ae57-5d10-9952-331d41423f4d"
LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
Plots = "91a5bcdd-55d7-5caf-9e0b-520d859cae80"
PlutoUI = "7f904dfe-b85e-4ff6-b463-dae2292396a8"
Random = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c"
SparseArrays = "2f01184e-e22b-5df5-ae63-d93ebab69eaf"

[compat]
Arpack = "~0.5.3"
Clustering = "~0.14.2"
GraphPlot = "~0.4.4"
LightGraphs = "~1.3.5"
Plots = "~1.21.3"
PlutoUI = "~0.7.9"
"""

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

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

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

[[ArnoldiMethod]]
deps = ["LinearAlgebra", "Random", "StaticArrays"]
git-tree-sha1 = "f87e559f87a45bece9c9ed97458d3afe98b1ebb9"
uuid = "ec485272-7323-5ecc-a04f-4719b315124d"
version = "0.1.0"

[[Arpack]]
deps = ["Arpack_jll", "Libdl", "LinearAlgebra", "Logging"]
git-tree-sha1 = "91ca22c4b8437da89b030f08d71db55a379ce958"
uuid = "7d9fca2a-8960-54d3-9f78-7d1dccf2cb97"
version = "0.5.3"

[[Arpack_jll]]
deps = ["Libdl", "OpenBLAS_jll", "Pkg"]
git-tree-sha1 = "e214a9b9bd1b4e1b4f15b22c0994862b66af7ff7"
uuid = "68821587-b530-5797-8361-c406ea357684"
version = "3.5.0+3"

[[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"

[[Clustering]]
deps = ["Distances", "LinearAlgebra", "NearestNeighbors", "Printf", "SparseArrays", "Statistics", "StatsBase"]
git-tree-sha1 = "75479b7df4167267d75294d14b58244695beb2ac"
uuid = "aaaa29a8-35af-508c-8bc3-b662a17a0fe5"
version = "0.14.2"

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

[[ColorTypes]]
deps = ["FixedPointNumbers", "Random"]
git-tree-sha1 = "32a2b8af383f11cbb65803883837a149d10dfe8a"
uuid = "3da002f7-5984-5a60-b8a6-cbb66c0b333f"
version = "0.10.12"

[[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 = "727e463cfebd0c7b999bbf3e9e7e16f254b94193"
uuid = "34da2185-b29b-5c13-b0c7-acf172513d20"
version = "3.34.0"

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

[[Compose]]
deps = ["Base64", "Colors", "DataStructures", "Dates", "IterTools", "JSON", "LinearAlgebra", "Measures", "Printf", "Random", "Requires", "Statistics", "UUIDs"]
git-tree-sha1 = "c6461fc7c35a4bb8d00905df7adafcff1fe3a6bc"
uuid = "a81c6b42-2e10-5240-aca2-a61377ecd94b"
version = "0.9.2"

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

[[DataAPI]]
git-tree-sha1 = "bec2532f8adb82005476c141ec23e921fc20971b"
uuid = "9a962f9c-6df0-11e9-0e5d-c546b8b5ee8a"
version = "1.8.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"

[[Distances]]
deps = ["LinearAlgebra", "Statistics", "StatsAPI"]
git-tree-sha1 = "9f46deb4d4ee4494ffb5a40a27a2aced67bdd838"
uuid = "b4f34e82-e78d-54a5-968a-f98e89d6e8f7"
version = "0.10.4"

[[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 = "dba1e8614e98949abfa60480b13653813d8f0157"
uuid = "0656b61e-2033-5cc2-a64a-77c0f6c09b89"
version = "3.3.5+0"

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

[[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"

[[GraphPlot]]
deps = ["ArnoldiMethod", "ColorTypes", "Colors", "Compose", "DelimitedFiles", "LightGraphs", "LinearAlgebra", "Random", "SparseArrays"]
git-tree-sha1 = "dd8f15128a91b0079dfe3f4a4a1e190e54ac7164"
uuid = "a2cc645c-3eea-5389-862e-a155d0052231"
version = "0.4.4"

[[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 = "60ed5f1643927479f845b0135bb369b031b541fa"
uuid = "cd3eb016-35fb-5094-929b-558a96fad6f3"
version = "0.9.14"

[[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"

[[Inflate]]
git-tree-sha1 = "f5fc07d4e706b84f72d54eedcc1c13d92fb0871c"
uuid = "d25df0c9-e2be-5dd7-82c8-3ad0b3e990b9"
version = "0.1.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"

[[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"

[[LightGraphs]]
deps = ["ArnoldiMethod", "DataStructures", "Distributed", "Inflate", "LinearAlgebra", "Random", "SharedArrays", "SimpleTraits", "SparseArrays", "Statistics"]
git-tree-sha1 = "432428df5f360964040ed60418dd5601ecd240b6"
uuid = "093fc24a-ae57-5d10-9952-331d41423f4d"
version = "1.3.5"

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

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

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

[[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 = "2ca267b08821e86c5ef4376cffed98a46c2cb205"
uuid = "e1d29d7a-bbdc-5cf2-9ac0-f12de2c33e28"
version = "1.0.1"

[[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"

[[NearestNeighbors]]
deps = ["Distances", "StaticArrays"]
git-tree-sha1 = "16baacfdc8758bc374882566c9187e785e85c2f0"
uuid = "b8a86587-4115-5ab1-83bc-aa920d37bbce"
version = "0.4.9"

[[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"

[[OpenBLAS_jll]]
deps = ["Artifacts", "CompilerSupportLibraries_jll", "Libdl"]
uuid = "4536629a-c528-5b80-bd46-f80d51c5b363"

[[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 = "438d35d2d95ae2c5e8780b330592b6de8494e779"
uuid = "69de0a69-1ddd-5017-9359-2bf0b02dc9f0"
version = "2.0.3"

[[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 = "9ff1c70190c1c30aebca35dc489f7411b256cd23"
uuid = "995b91a9-d308-5afd-9ec6-746e21dbc043"
version = "1.0.13"

[[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 = "2dbafeadadcf7dadff20cd60046bba416b4912be"
uuid = "91a5bcdd-55d7-5caf-9e0b-520d859cae80"
version = "1.21.3"

[[PlutoUI]]
deps = ["Base64", "Dates", "InteractiveUtils", "JSON", "Logging", "Markdown", "Random", "Reexport", "Suppressor"]
git-tree-sha1 = "44e225d5837e2a2345e69a1d1e01ac2443ff9fcb"
uuid = "7f904dfe-b85e-4ff6-b463-dae2292396a8"
version = "0.7.9"

[[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 = "d4491becdc53580c6dadb0f6249f90caae888554"
uuid = "01d81517-befc-4cb6-b9ec-a95719d0359c"
version = "0.4.0"

[[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"

[[SimpleTraits]]
deps = ["InteractiveUtils", "MacroTools"]
git-tree-sha1 = "5d7e3f4e11935503d3ecaf7186eac40602e7d231"
uuid = "699a6c99-e7fa-54fc-8d76-47d257e15c1d"
version = "0.9.4"

[[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 = "1700b86ad59348c0f9f68ddc95117071f947072d"
uuid = "09ab397b-f2b6-538f-b94a-2f83cf4a842a"
version = "0.6.1"

[[Suppressor]]
git-tree-sha1 = "a819d77f31f83e5792a76081eee1ea6342ab8787"
uuid = "fd094767-a336-5f1f-9728-57cf17d0bbfb"
version = "0.2.0"

[[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 = "d0c690d37c73aeb5ca063056283fde5585a41710"
uuid = "bd369af6-aec1-5ad0-b16a-f7cc5008161c"
version = "1.5.0"

[[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:
# ╠═de9d25b0-9ecc-4752-8c4d-e452dcf82c57
# ╠═0f6cd296-6cb8-4a9d-a21f-02b689dfe4e9
# ╟─89299cb7-a4eb-48f3-b3f1-ab2b5c6cfa2e
# ╟─a92d0e70-479d-444e-b06d-e47d2e385b10
# ╟─fddcc2a7-f718-4900-a135-764eb36496a7
# ╟─b3b44560-203a-11eb-2782-170e4c12627a
# ╟─1be90f30-203b-11eb-34da-c3116523ec3a
# ╟─5211ce2e-203b-11eb-0b6f-55e809177fb9
# ╟─38333a20-203c-11eb-2b28-59e07b3e50ee
# ╟─1bc20329-3a6e-435d-a25c-8cce79fb7d9d
# ╠═b9d0830a-0a33-40c1-a08a-49538cf7a879
# ╠═b6d41ab1-24f2-417d-b604-3e12384e65f8
# ╠═f6492dfe-faac-407f-9a0e-5dcf321a4f2d
# ╠═1a8b6f7b-4098-4cd2-adad-8ee788b47333
# ╠═4079aa8e-2bb5-4d51-bca8-04d0997db2bf
# ╠═d5d8ee5d-9cce-45fc-a31b-b2ed8237cc86
# ╠═86fd4f5a-af57-490d-a17f-2339a9414057
# ╠═1f53b5cf-ab6a-4107-b53b-9ce96afcdb7f
# ╠═63403411-cca4-4d71-8a76-8dea8638f69f
# ╠═3cf471bc-d316-4de9-af63-71db9b035c40
# ╠═ac0b5eda-299c-4ae2-ae78-bda92b378f19
# ╠═9bf6bbca-f96c-4b38-a4fd-5188ec29a2a2
# ╠═aa51688c-21a0-40f4-a6d6-01d441bfd4b7
# ╠═5f27f1a2-80de-41d1-b1e2-da3a99493217
# ╠═7cd3af86-f459-4e4b-8e35-c29bc6966eed
# ╠═e1149b05-d594-4d75-ae56-1ce18bd0f60d
# ╠═b95d4070-c127-4db0-987b-5fc17421b667
# ╠═57d472c1-548d-4678-b884-85f863771aa7
# ╠═4806378b-63cb-445a-8c71-583458fce846
# ╠═a08c99a3-9c10-47be-9dc6-987d92a85087
# ╠═6888e8a9-f306-4490-9147-2460fbda4889
# ╠═3bd5ad15-31cd-4efd-b4e3-7cf7ede3edae
# ╟─e4487f35-5de0-492e-96f1-69318d7d8e8b
# ╠═81d99f24-82b9-4216-8d17-6a18317beb77
# ╠═8f925ff7-90b4-4fe4-a869-86dec693c0f4
# ╠═d446724a-97d1-403f-baaf-d0619203f351
# ╠═7572c241-dbcc-4094-97e3-ee3172255294
# ╠═22e90c11-ba4e-4073-9f78-f3b7beb24e69
# ╟─00000000-0000-0000-0000-000000000001
# ╟─00000000-0000-0000-0000-000000000002

ext: jl

max_stars_repo_path: Lectures/C4 Spectral Partitioning of Bipartite Graphs.jl


text:
@testset "DirichletProblem" begin

    @testset "constant solution" begin
    prob = DirichletProblem(identity, i=1, b=1)

    θ = solve(prob)

    ϕ = range(0, 20, length=100)

    @test all(@inferred θ.(ϕ) .≈ θ.i)
    @test all(@inferred θ.(ϕ) .≈ θ.b)
    @test all(@inferred d_dϕ.(θ, ϕ) .≈ 0)
    @test θ.iterations == 0
    @test isnan(@inferred θ(-1))
    end

    @testset "exact" begin
    # Reference: Philip (1960) Table 1, No. 13
    # https://doi.org/10.1071/PH600001
    prob = DirichletProblem(θ -> 0.5*(1 - log(θ)), i=0, b=1)

    θ = solve(prob)

    ϕ = range(0, 20, length=100)

    @test all(@. isapprox(θ(ϕ), exp(-ϕ), atol=1e-3))
    @test all(@. isapprox(d_dϕ(θ,ϕ), -exp(-ϕ), atol=1e-3))
    @test θ.iterations > 0
    end

    @testset "HF135" begin

    r = [0.    , 0.0025, 0.005 , 0.0075, 0.01  , 0.0125, 0.015 , 0.0175,
         0.02  , 0.0225, 0.025 , 0.0275, 0.03  , 0.0325, 0.035 , 0.0375,
         0.04  , 0.0425, 0.045 , 0.0475, 0.05  ]

    t = 60

    # Data obtained with the porousMultiphaseFoam toolbox, version 1906
    # https://github.com/phorgue/porousMultiphaseFoam
    θ_pmf = [0.945   , 0.944845, 0.944188, 0.942814, 0.940517, 0.937055,
             0.93214 , 0.925406, 0.916379, 0.904433, 0.888715, 0.868016,
             0.840562, 0.803597, 0.752494, 0.678493, 0.560999, 0.314848,
             0.102755, 0.102755, 0.102755]
    U_pmf = [2.66135e-04,  2.66133e-04,  2.66111e-04,  2.66038e-04,
             2.65869e-04,  2.65542e-04,  2.64975e-04,  2.64060e-04,
             2.62644e-04,  2.60523e-04,  2.57404e-04,  2.52863e-04,
             2.46269e-04,  2.36619e-04,  2.22209e-04,  1.99790e-04,
             1.61709e-04,  7.64565e-05, -2.45199e-21, -7.35598e-21,
             0.00000e+00]

    ϵ = 1e-7

    # Wetting of an HF135 membrane, Van Genuchten model
    # Data from Buser (PhD thesis, 2016)
    # http://hdl.handle.net/1773/38064
    θr = 0.0473
    θs = 0.945
    k = 5.50e-13  # m**2
    α = 0.2555  # 1/m
    n = 2.3521
    θi = 0.102755  # Computed from P0

    θb = θs - ϵ

    model = VanGenuchten(n=n, α=α, k=k, θr=θr, θs=θs)

    prob = DirichletProblem(model, i=θi, b=θb)

    θ = solve(prob, itol=1e-7)

    @test all(@. isapprox(θ(r,t), θ_pmf, atol=1e-3))
    @test all(@. isapprox(flux(θ,r,t), U_pmf, atol=1e-6))
    @test θ.iterations > 0
    @test isnan(@inferred θ(-1,t))
    end

    @testset "bad argument" begin
    prob = DirichletProblem(identity, i=0, b=1)
    @test_throws ArgumentError solve(prob, maxiter=-1)
    end

    @testset "unsolved" begin
    prob = DirichletProblem(identity, i=0, b=1)
    @test_throws SolvingError solve(prob, maxiter=0)
    end

    @testset "unsolvable" begin
    prob = DirichletProblem(identity, i=-1e-3, b=1)
    @test_throws DomainError solve(prob)

    prob = DirichletProblem(identity, i=0, b=-1)
    @test_throws DomainError solve(prob)
    end

end


text:
module GModel

using DocStringExtensions

using Random
using Distributions
using LinearAlgebra
using FFTW

export run_G
export run_G_ensemble
export lorenz_forward

"""
    GSettings{FT<:AbstractFloat, KT, D}

Structure to hold all information to run the forward model G

# Fields
$(DocStringExtensions.FIELDS)
"""
struct GSettings{FT<:AbstractFloat, KT, D}
    "model output"
    out::Array{FT, 1}
    "time period over which to run the model, e.g., `(0, 1)`"
    tspan::Tuple{FT, FT}
    "x domain"
    x::Array{FT, 1}
    "y domain"
    y::Array{FT, 1}
end

struct LSettings
	# Stationary or transient dynamics
	dynamics::Int32
	# G model statistics type
	stats_type::Int32
	# Integration time start
	t_start::Float64
	# Data collection length
	T::Float64
	# Integration length
	Ts::Float64
	# Duration of polynomial fit, number of samples of Ts to fit over
	Tfit::Float64
	# Initial perturbation
	Fp::Array{Float64}
	# Number of longitude steps
	N::Int32
	# Timestep
	dt::Float64
	# Simulation end time
	tend::Float64
	# For stats_type=2, number of frequencies to consider
	kmax::Int32
	# Should ω be learned using CES?
	ω_fixed::Bool
	# Truth ω
	ω_true::Float64
end

struct LParams
	# Mean forcing
	F::Float64
	# Forcing frequency for transient terms
	ω::Float64
	# Forcing amplitude for transient terms
	A::Float64
end

"""
    run_G_ensemble(params::Array{FT, 2},
                   settings::GSettings{FT},
                   update_params,
                   moment,
                   get_src;
                   rng_seed=42) where {FT<:AbstractFloat}

Run the forward model G for an array of parameters by iteratively
calling run_G for each of the N_ensemble parameter values.
Return g_ens, an array of size N_data x N_ensemble, where
g_ens[:,j] = G(params[:,j])

 - `params` - array of size N_parameters x N_ensemble containing the
              parameters for which G will be run
 - `settings` - a GSetttings struct

"""
function run_G_ensemble(params::Array{FT, 2},
                        settings::LSettings,
                        rng_seed=42) where {FT<:AbstractFloat}
    # Initilize ensemble
    N_ens = size(params, 2) # params is N_params x N_ens
    if settings.stats_type == 1
	    nd = 2
    elseif settings.stats_type == 2
	    nd = 1 + (2*settings.kmax)
    elseif settings.stats_type == 3
	    nd = 3
    elseif settings.stats_type == 4
	    nd = Int64(2*(settings.T/settings.Ts/settings.Tfit))
    elseif settings.stats_type == 5
	    nd = Int64(settings.T/settings.Ts/settings.Tfit)
    end
    g_ens = zeros(nd, N_ens)
    # Lorenz parameters
    #Fp = rand(Normal(0.0, 0.01), N); # Initial perturbation
    F = params[1, :] # Forcing
    if settings.dynamics==2
        A = params[2, :] # Transience amplitude
    else
        A = F .* 0. # Set to zero, it is not used
    end
    if settings.ω_fixed==false;
    	ω = params[3, :] # Transience frequency
    end

    Random.seed!(rng_seed)
    for i in 1:N_ens
        # run the model with the current parameters, i.e., map θ to G(θ)
	if settings.ω_fixed==false
	    lorenz_params = GModel.LParams(F[i], ω[i], A[i])
	elseif settings.ω_fixed==true
	    lorenz_params = GModel.LParams(F[i], settings.ω_true, A[i])
	end
	g_ens[:, i] = lorenz_forward(settings, lorenz_params)
    end

    return g_ens
end

# Forward pass of the Lorenz 96 model
# Inputs: settings: structure with stats_type, t_start, T
# F: scalar forcing Float64(1), Fp: initial perturbation Float64(N)
# N: number of longitude steps Int64(1), dt: time step Float64(1)
# tend: End of simulation Float64(1), nstep:
function lorenz_forward(settings::LSettings, params::LParams)
	# run the Lorenz simulation
	xn, t = lorenz_solve(settings, params)
	# Get statistics
	gt = stats(settings, xn, t)
	return gt
end

# Solve the Lorenz 96 system
# Inputs: F: scalar forcing Float64(1), Fp: initial perturbation Float64(N)
# N: number of longitude steps Int64(1), dt: time step Float64(1)
# tend: End of simulation Float64(1), nstep:
function lorenz_solve(settings::LSettings, params::LParams)
	# Initialize
	nstep = Int32(ceil(settings.tend/settings.dt));
	xn = zeros(Float64, settings.N, nstep)
	t = zeros(Float64, nstep)
	# Initial perturbation
	X = fill(Float64(0.), settings.N)
	X = X + settings.Fp
	# March forward in time
	for j in 1:nstep
		t[j] = settings.dt*j
		if settings.dynamics==1
		    	X = RK4(X,settings.dt,settings.N,params.F)
	        elseif settings.dynamics==2
			F_local = params.F + params.A*sin(params.ω*t[j])
		    	X = RK4(X,settings.dt,settings.N,F_local)
	        end
		xn[:,j] = X
	end
	# Output
	return xn, t

end

# Lorenz 96 system
# f = dx/dt
# Inputs: x: state, N: longitude steps, F: forcing
function f(x,N,F)
	f = zeros(Float64, N)
	# Loop over N positions
	for i in 3:N-1
		f[i] = -x[i-2]*x[i-1] + x[i-1]*x[i+1] - x[i] + F
	end
	# Periodic boundary conditions
	f[1] = -x[N-1]*x[N] + x[N]*x[2] - x[1] + F
	f[2] = -x[N]*x[1] + x[1]*x[3] - x[2] + F
	f[N] = -x[N-2]*x[N-1] + x[N-1]*x[1] - x[N] + F
	# Output
	return f
end

# RK4 solve
function RK4(xold, dt, N, F)
	# Predictor steps
	k1 = f(xold, N, F)
	k2 = f(xold + k1*dt/2., N, F)
	k3 = f(xold + k2*dt/2., N, F)
	k4 = f(xold + k3*dt, N, F)
	# Step
	xnew = xold + (dt/6.0)*(k1 + 2.0*k2 + 2.0*k3 + k4)
	# Output
	return xnew
end

function regression(X,y)
	beta_tilde = [ ones(length(X),1) X ] \ y;
	v = beta_tilde[1];  beta = beta_tilde[2];
	return beta, v
end

function spectra(signal, Ts, t)
	# Init
	Nt = length(t)
	if mod(Nt,2)!=0
		t=t[1:Nt-1]; signal = signal[1:Nt-1]; Nt = Nt-1;
	end
	# Remove mean
	u = signal .- mean(signal,dims=1)
	# FFT
	F = fft(u) |> fftshift
	freqs = fftfreq(length(t), 1.0/Ts) |> fftshift
	#k = gen_k(Nt,Ts)  |> fftshift
	T = (Nt-1)*Ts;
	A = Ts^2 / (2*pi*T);
	uhat = A*abs.(F.^2);
	# Output
	mid = Int32(Nt/2)+1;
	f = freqs[mid:Nt]
	Fp = uhat[mid:Nt]
	return f, Fp
end

###############################
## Extract statistics
###############################
function stats(settings, xn, t)
# Define averaging indices range
	indices = findall(x -> (x>settings.t_start)
			  && (x<settings.t_start+settings.T), t)
	# Define statistics of interest
	if settings.stats_type == 1 # Mean
		# Average in time and over longitude
		gtm = mean(mean(xn[:,indices], dims=2), dims=1)
		# Variance
		gtc = mean(mean((xn[:,indices].-gtm).^2, dims=2), dims=1)
		# Combine statistics
		gt = vcat(gtm, gtc)
	elseif settings.stats_type == 2 #What about an integral under parts of the spectrum
		# Average in time and over longitude
		gtm = mean(mean(xn[:,indices], dims=2), dims=1)
		# Power spectra
		# Get size
		f, Fp = spectra(xn[1,indices]', settings.dt, t[indices]);
		# Compute temporal spectrum for each longitude
		Fp = zeros(size(Fp,1), settings.N)
		for i in 1:settings.N
			f, Fp[:,i] = spectra(xn[i,indices]', settings.dt, t[indices]);
		end
		# Average spectra over periodic directions
		Fp = dropdims(mean(Fp, dims=2), dims=2)
		ys = partialsortperm(Fp, 1:settings.kmax; rev=true)
		gt = vcat(gtm, Fp[ys]..., f[ys]...)
	elseif settings.stats_type == 3 # Structure function
		# Average in time and over longitude
		gtm = mean(mean(xn[:,indices], dims=2), dims=1)
		# Maximum
		mxval, mxind = findmax(xn[:,indices]; dims=2);
		mxval_out = mean(mxval, dims=1);
		# Power spectra
		# Get size
		f, Fp = spectra(xn[1,indices]', settings.dt, t[indices]);
		# Compute temporal spectrum for each longitude
		Fp = zeros(size(Fp,1), settings.N)
		for i in 1:settings.N
			f, Fp[:,i] = spectra(xn[i,indices]', settings.dt, t[indices]);
		end
		# Average spectra over periodic directions
		Fp = dropdims(mean(Fp, dims=2), dims=2)
		ys = partialsortperm(Fp, 1; rev=true)
		# Period
		T = 1. / f[ys]
		mxval, r = findmin(abs.(t.-T); dims=1); r = r[1]
                # Structure function
		xp = xn .- gtm
		st = zeros(settings.N)
		for i in 1:settings.N
			st[i] = mean( ( xp[1,indices[1]:indices[length(indices)]-r] .-
				       xp[1,indices[1]+r:indices[length(indices)]] ).^2 );
		end
		# Combine
		gt = vcat(gtm, mean(st)..., mxval_out...)
	elseif settings.stats_type == 4 # Variance sub-samples, linear fit over window
		T = settings.T
		# Calculate variance on subsamples
		# Variances over sliding windows
		Ts = settings.Ts
		N = Int64(T/Ts); # number of sliding windows
		nt = Int64(Ts/settings.dt)-1 # length of sliding windows
		var_slide = zeros(N); t_slide = zeros(nt,N);
		for i in 1:N;
			var_slide[i] = var(vcat(xn[:,indices[(i-1)*nt+1:i*nt]]...))
			t_slide[:,i] = t[(i-1)*nt+1:i*nt];
		end
		# Polynomial fit over a batch of sliding windows
		Tfit = Int64(settings.Tfit);
		NT = Int64(N/Tfit); a1 = zeros(NT); a2 = zeros(NT); t_start = zeros(NT)
		y = zeros(Tfit,NT); ty = zeros(Tfit,NT); ym = zeros(Tfit,NT)
		for i in 1:NT
			ind_low = (i-1)*Tfit+1; ind_high = i*Tfit;
			t_start[i] = t_slide[1,ind_low]
			ty[:,i] = t_slide[1, ind_low:ind_high] .- t_start[i]
			a1[i], a2[i] = regression(ty[:,i], var_slide[ind_low:ind_high])
			y[:,i] = a1[i].*ty[:,i] .+ a2[i]
			ym[i] = mean(var_slide[ind_low:ind_high])
		end
		# Combine
		gt = vcat(a1..., a2...)
	elseif settings.stats_type == 5 # Variance sub-samples, mean over window
		T = settings.T
		# Calculate variance on subsamples
		# Variances over sliding windows
		Ts = settings.Ts
		N = Int64(T/Ts); # number of sliding windows
		nt = Int64(Ts/settings.dt)-1 # length of sliding windows
		var_slide = zeros(N); t_slide = zeros(nt,N);
		for i in 1:N;
			var_slide[i] = var(vcat(xn[:,indices[(i-1)*nt+1:i*nt]]...))
			t_slide[:,i] = t[(i-1)*nt+1:i*nt];
		end
		# Polynomial fit over a batch of sliding windows
		Tfit = Int64(settings.Tfit);
		NT = Int64(N/Tfit); a1 = zeros(NT); a2 = zeros(NT); t_start = zeros(NT)
		y = zeros(Tfit,NT); ty = zeros(Tfit,NT); ym = zeros(NT)
		for i in 1:NT
			ind_low = (i-1)*Tfit+1; ind_high = i*Tfit;
			t_start[i] = t_slide[1,ind_low]
			ty[:,i] = t_slide[1, ind_low:ind_high] .- t_start[i]
			a1[i], a2[i] = regression(ty[:,i], var_slide[ind_low:ind_high])
			y[:,i] = a1[i].*ty[:,i] .+ a2[i]
			ym[i] = mean(var_slide[ind_low:ind_high])
		end
		# Combine
		gt = vcat(ym...)
	else
		ArgumentError("Setting "*string(settings.stats_type)*" not implemented.")
	end
	return gt
end

## Forward pass of the Lorenz 96 model
## Inputs: settings: structure with stats_type, t_start, T
## F: scalar forcing Float64(1), Fp: initial perturbation Float64(N)
## N: number of longitude steps Int64(1), dt: time step Float64(1)
## tend: End of simulation Float64(1), nstep:
#function lorenz_forward(settings::LSettings, F)
#	# run the Lorenz simulation
#	xn, t = lorenz_solve(settings, F)
#	# Define averaging indices range
#	indices = findall(x -> (x>settings.t_start)
#			  && (x<settings.t_start+settings.T), t)
#	# Define statistics of interest
#	if settings.stats_type == 1 # Mean
#		# Average in time and over longitude
#		gtm = mean(mean(xn[:,indices], dims=2), dims=1)
#		# Variance
#		gtc = mean(mean((xn[:,indices].-gtm).^2, dims=2), dims=1)
#		# Combine statistics
#		gt = vcat(gtm, gtc)
#	elseif settings.stats_type == 2 #What about an integral under parts of the spectrum
#		# Average in time and over longitude
#		gtm = mean(mean(xn[:,indices], dims=2), dims=1)
#		# Power spectra
#		# Get size
#		f, Fp = spectra(xn[1,indices]', settings.dt, t[indices]);
#		# Compute temporal spectrum for each longitude
#		Fp = zeros(size(Fp,1), settings.N)
#		for i in 1:settings.N
#			f, Fp[:,i] = spectra(xn[i,indices]', settings.dt, t[indices]);
#		end
#		# Average spectra over periodic directions
#		Fp = dropdims(mean(Fp, dims=2), dims=2)
#		ys = partialsortperm(Fp, 1:settings.kmax; rev=true)
#		gt = vcat(gtm, Fp[ys]..., f[ys]...)
#	elseif settings.stats_type == 3 # Structure function
#		# Average in time and over longitude
#		gtm = mean(mean(xn[:,indices], dims=2), dims=1)
#		# Maximum
#		mxval, mxind = findmax(xn[:,indices]; dims=2);
#		mxval_out = mean(mxval, dims=1);
#		# Power spectra
#		# Get size
#		f, Fp = spectra(xn[1,indices]', settings.dt, t[indices]);
#		# Compute temporal spectrum for each longitude
#		Fp = zeros(size(Fp,1), settings.N)
#		for i in 1:settings.N
#			f, Fp[:,i] = spectra(xn[i,indices]', settings.dt, t[indices]);
#		end
#		# Average spectra over periodic directions
#		Fp = dropdims(mean(Fp, dims=2), dims=2)
#		ys = partialsortperm(Fp, 1; rev=true)
#		# Period
#		T = 1. / f[ys]
#		mxval, r = findmin(abs.(t.-T); dims=1); r = r[1]
#                # Structure function
#		xp = xn .- gtm
#		st = zeros(settings.N)
#		for i in 1:settings.N
#			st[i] = mean( ( xp[1,indices[1]:indices[length(indices)]-r] .-
#				       xp[1,indices[1]+r:indices[length(indices)]] ).^2 );
#		end
#		# Combine
#		gt = vcat(gtm, mean(st)..., mxval_out...)
#	else
#		ArgumentError("Setting "*string(settings.stats_type)*" not implemented.")
#	end
#	return gt
#end
#

end # module


text:
abstract type RobertsRosenthalMCTune <: MCTunerState end

mutable struct UnvRobertsRosenthalMCTune <: RobertsRosenthalMCTune
  logσ::Real # Standard deviation of underlying normal
  δ::Real # Increase or decrease of logσ at each batch
  batch::Integer
  accepted::Integer # Number of accepted MCMC samples during current tuning period
  proposed::Integer # Number of proposed MCMC samples during current tuning period
  totproposed::Integer # Total number of proposed MCMC samples during burnin
  rate::Real # Observed acceptance rate over current tuning period

  function UnvRobertsRosenthalMCTune(
    logσ::Real, δ::Real, batch::Integer, accepted::Integer, proposed::Integer, totproposed::Integer, rate::Real
  )
    @assert batch >= 0 "Batch should be non-negative"
    @assert accepted >= 0 "Number of accepted MCMC samples should be non-negative"
    @assert proposed >= 0 "Number of proposed MCMC samples should be non-negative"
    @assert totproposed >= 0 "Total number of proposed MCMC samples should be non-negative"
    if !isnan(rate)
      @assert 0 <= rate <= 1 "Observed acceptance rate should be in [0, 1]"
    end
    new(logσ, NaN, batch, accepted, proposed, totproposed, rate)
  end
end

UnvRobertsRosenthalMCTune(; logσ::Real=0., accepted::Integer=0, proposed::Integer=0, totproposed::Integer=0) =
  UnvRobertsRosenthalMCTune(logσ, NaN, 0, accepted, proposed, totproposed, NaN)

mutable struct MuvRobertsRosenthalMCTune <: RobertsRosenthalMCTune
  logσ::RealVector # Standard deviation of underlying normal
  δ::Real # Increase or decrease of logσ at each batch
  batch::Integer
  accepted::IntegerVector # Number of accepted MCMC samples during current tuning period
  proposed::Integer # Number of proposed MCMC samples during current tuning period
  totproposed::Integer # Total number of proposed MCMC samples during burnin
  rate::RealVector # Observed acceptance rate over current tuning period

  function MuvRobertsRosenthalMCTune(
    logσ::RealVector,
    δ::Real,
    batch::Integer,
    accepted::IntegerVector,
    proposed::Integer,
    totproposed::Integer,
    rate::RealVector
  )
    @assert batch >= 0 "Batch should be non-negative"
    @assert all(i -> i >= 0, accepted) "Number of accepted MCMC samples should be non-negative"
    @assert proposed >= 0 "Number of proposed MCMC samples should be non-negative"
    @assert totproposed >= 0 "Total number of proposed MCMC samples should be non-negative"
    if all(i -> !isnan(i), rate)
      @assert all(i -> 0 <= i <= 1, rate) "Observed acceptance rates should be in [0, 1]"
    end
    new(logσ, NaN, batch, accepted, proposed, totproposed, rate)
  end
end

MuvRobertsRosenthalMCTune(
  logσ::RealVector; accepted::IntegerVector=fill(0, length(logσ)), proposed::Integer=0, totproposed::Integer=0
) =
  MuvRobertsRosenthalMCTune(logσ, NaN, 0, accepted, proposed, totproposed, fill(NaN, length(logσ)))

MuvRobertsRosenthalMCTune(
  d::Integer; logσ::Real=0., accepted::IntegerVector=fill(0, d), proposed::Integer=0, totproposed::Integer=0
) =
  MuvRobertsRosenthalMCTune(fill(logσ, d), NaN, 0, accepted, proposed, totproposed, fill(NaN, d))

reset!(tune::UnvRobertsRosenthalMCTune) = reset_burnin!(tune::MCTunerState)

function reset!(tune::MuvRobertsRosenthalMCTune, i::Integer)
  tune.accepted[i] = 0
  tune.rate[i] = NaN
end

function reset!(tune::MuvRobertsRosenthalMCTune)
  tune.totproposed += tune.proposed
  tune.proposed = 0
end

rate!(tune::MuvRobertsRosenthalMCTune, i::Integer) = (tune.rate[i] = tune.accepted[i]/tune.proposed)

### RobertsRosenthalMCTuner

struct RobertsRosenthalMCTuner <: MCTuner
  targetrate::Real # Target acceptance rate
  period::Integer # Tuning period over which acceptance rate is computed
  verbose::Bool # If verbose=false then the tuner is silent, else it is verbose

  function RobertsRosenthalMCTuner(targetrate::Real, period::Integer, verbose::Bool)
    @assert 0 < targetrate < 1 "Target acceptance rate should be between 0 and 1"
    @assert period > 0 "Tuning period should be positive"
    new(targetrate, period, verbose)
  end
end

RobertsRosenthalMCTuner(; targetrate::Real=0.44, period::Integer=50, verbose::Bool=false) =
  RobertsRosenthalMCTuner(targetrate, period, verbose)

set_batch!(tune::RobertsRosenthalMCTune, tuner::RobertsRosenthalMCTuner) = (tune.batch = tune.totproposed/tuner.period)

set_delta!(tune::RobertsRosenthalMCTune, tuner::RobertsRosenthalMCTuner) = (tune.δ = min(0.01, tune.batch^-0.5))

tune!(tune::UnvRobertsRosenthalMCTune, tuner::RobertsRosenthalMCTuner) =
  tune.logσ += (tune.rate < tuner.targetrate ? -tune.δ : tune.δ)

tune!(tune::MuvRobertsRosenthalMCTune, tuner::RobertsRosenthalMCTuner, i::Integer) =
  tune.logσ[i] += (tune.rate[i] < tuner.targetrate ? -tune.δ : tune.δ)

show(io::IO, tuner::RobertsRosenthalMCTuner) =
  print(
    io,
    string(
      "RobertsRosenthalMCTuner: target rate = ",
      tuner.targetrate,
      ", period = ",
      tuner.period,
      ", verbose = ",
      tuner.verbose
    )
  )


text:
# Ridge Regression (Tikhonov regularization)

#### auxiliary

function lreg_chkdims(X::AbstractMatrix, Y::AbstractVecOrMat, trans::Bool)
    mX, nX = size(X)
    dX = ifelse(trans, mX, nX)
    dY = ifelse(trans, nX, mX)
    size(Y, 1) == dY || throw(DimensionMismatch("Dimensions of X and Y mismatch."))
    return dX
end

lrsoltype(::AbstractVector{T}) where T = Vector{T}
lrsoltype(::AbstractMatrix{T}) where T = Matrix{T}

_vaug(X::AbstractMatrix{T}) where T = vcat(X, ones(T, 1, size(X,2)))::Matrix{T}
_haug(X::AbstractMatrix{T}) where T = hcat(X, ones(T, size(X,1), 1))::Matrix{T}

## linear least square

function llsq(X::AbstractMatrix{T}, Y::AbstractVecOrMat{T};
              trans::Bool=false, bias::Bool=true) where T<:AbstractFloat
    if trans
        mX, nX = size(X)
        size(Y, 1) == nX || throw(DimensionMismatch("Dimensions of X and Y mismatch."))
        mX <= nX || error("mX <= nX is required when trans is false.")
    else
        mX, nX = size(X)
        size(Y, 1) == mX || throw(DimensionMismatch("Dimensions of X and Y mismatch."))
        mX >= nX || error("mX >= nX is required when trans is false.")
    end
    _ridge(X, Y, zero(T), trans, bias)
end

## ridge regression

function ridge(X::AbstractMatrix{T}, Y::AbstractVecOrMat{T}, r::Real;
               trans::Bool=false, bias::Bool=true) where T<:AbstractFloat
    lreg_chkdims(X, Y, trans)
    r >= zero(r) || error("r must be non-negative.")
    _ridge(X, Y, convert(T, r), trans, bias)
end

function ridge(X::AbstractMatrix{T}, Y::AbstractVecOrMat{T}, r::AbstractVector{T};
               trans::Bool=false, bias::Bool=true) where T<:AbstractFloat
    d = lreg_chkdims(X, Y, trans)
    length(r) == d || throw(DimensionMismatch("Incorrect length of r."))
    _ridge(X, Y, r, trans, bias)
end

function ridge(X::AbstractMatrix{T}, Y::AbstractVecOrMat{T}, r::AbstractMatrix{T};
               trans::Bool=false, bias::Bool=true) where T<:AbstractFloat
    d = lreg_chkdims(X, Y, trans)
    size(r) == (d, d) || throw(DimensionMismatch("Incorrect size of r."))
    _ridge(X, Y, r, trans, bias)
end

## implementation

function _ridge(X::AbstractMatrix{T}, Y::AbstractVecOrMat{T},
                r::Union{Real, AbstractVecOrMat}, trans::Bool, bias::Bool) where T<:AbstractFloat
    if bias
        if trans
            X_ = _vaug(X)
            A = cholesky!(Hermitian(_ridge_reg!(X_ * transpose(X_), r, bias))) \ (X_ * Y)
        else
            X_ = _haug(X)
            A = cholesky!(Hermitian(_ridge_reg!(X_'X_, r, bias))) \ (X_'Y)
        end
    else
        if trans
            A = cholesky!(Hermitian(_ridge_reg!(X * X', r, bias))) \ (X * Y)
        else
            A = cholesky!(Hermitian(_ridge_reg!(X'X, r, bias))) \ (X'Y)
        end
    end
    return A::lrsoltype(Y)
end

function _ridge_reg!(Q::Matrix, r::Real, bias::Bool)
    if r > zero(r)
        n = size(Q, 1) - Int(bias)
        for i = 1:n
            @inbounds Q[i,i] += r
        end
    end
    return Q
end

function _ridge_reg!(Q::AbstractMatrix, r::AbstractVector, bias::Bool)
    n = size(Q, 1) - Int(bias)
    @assert length(r) == n
    for i = 1:n
        @inbounds Q[i,i] += r[i]
    end
    return Q
end

function _ridge_reg!(Q::AbstractMatrix, r::AbstractMatrix, bias::Bool)
    n = size(Q, 1) - Int(bias)
    @assert size(r) == (n, n)
    for j = 1:n, i = 1:n
        @inbounds Q[i,j] += r[i,j]
    end
    return Q
end