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

using Markdown
using InteractiveUtils

# This Pluto notebook uses @bind for interactivity. When running this notebook outside of Pluto, the following 'mock version' of @bind gives bound variables a default value (instead of an error).
macro bind(def, element)
    quote
        local iv = try Base.loaded_modules[Base.PkgId(Base.UUID("6e696c72-6542-2067-7265-42206c756150"), "AbstractPlutoDingetjes")].Bonds.initial_value catch; b -> missing; end
        local el = $(esc(element))
        global $(esc(def)) = Core.applicable(Base.get, el) ? Base.get(el) : iv(el)
        el
    end
end

# ╔═╡ 8fa836a6-1133-4a54-b996-a02083fc6bba
begin
    using Pkg
	Pkg.activate(joinpath(Pkg.devdir(), "MLCourse"))
    using Random, Statistics, DataFrames, Plots, MLJ, MLJLinearModels, MLCourse
	import MLCourse: fitted_linear_func
end

# ╔═╡ 269a609c-74af-4e7e-86df-e2279096a7a6
using NearestNeighborModels

# ╔═╡ 12942f34-efb1-11eb-3eb4-c1a38396cfb8
begin
    using PlutoUI
    PlutoUI.TableOfContents()
end

# ╔═╡ 12942f50-efb1-11eb-01c0-055b6be166e0
md"# When Linear Methods Are Not Flexible Enough

We will use the function ``f(x) = 0.3 \sin(10x) + 0.7 x`` to design two artificial
data sets: one for regression with a single predictor and one for classification
with two predictors.

In this section we will see that standard linear regression or classification is not flexible enough to fit the artifical data sets.
"

# ╔═╡ 12942f5a-efb1-11eb-399a-a1300d636217
begin
    f(x) = .3 * sin(10x) + .7x
    # f(x) = .0015*(12x-6)^4 -.035(12x-6)^2 + .5x + .2  # an alternative
	σ(x) = 1 / (1 + exp(-x))
    function regression_data_generator(; n, rng = Random.GLOBAL_RNG)
        x = rand(rng, n)
        DataFrame(x = x, y = f.(x) .+ .1*randn(rng, n))
    end
    function classification_data_generator(; n, rng = Random.GLOBAL_RNG)
        X1 = rand(rng, n)
        X2 = rand(rng, n)
        df = DataFrame(X1 = X1, X2 = X2,
                       y = categorical(σ.(20(f.(X1) .- X2)) .> rand(rng, n),
					                   levels = [false, true], ordered = true))
    end
end;

# ╔═╡ 12942f62-efb1-11eb-3f59-f981bc32f308
regression_data = regression_data_generator(n = 50, rng = MersenneTwister(3))

# ╔═╡ 12942f6e-efb1-11eb-0a49-01a6a2d0196f
m1 = machine(LinearRegressor(), select(regression_data, :x), regression_data.y) |> fit!;

# ╔═╡ 12942f6e-efb1-11eb-35c0-19a7c4d6b44e
begin
    scatter(regression_data.x, regression_data.y, label = "training data")
    plot!(fitted_linear_func(m1), w = 3, xlabel = "X", ylabel = "Y",
          label = "linear regression", legend = :topleft)
end

# ╔═╡ 12942f94-efb1-11eb-234e-51492b51e583
classification_data = classification_data_generator(n = 400, rng = MersenneTwister(8))

# ╔═╡ 12942fa0-efb1-11eb-01c8-6dae80c55fb8
m2 = machine(LogisticClassifier(penalty = :none),
             select(classification_data, Not(:y)),
             classification_data.y) |> fit!;

# ╔═╡ 12942f94-efb1-11eb-2c48-a3418b53b886
begin
    xgrid = MLCourse.grid(0:.01:1, 0:.02:1, names = (:X1, :X2))
    scatter(xgrid.X1, xgrid.X2, color = coerce(predict_mode(m2, xgrid), Count),
            markersize = 2, label = nothing, markerstrokewidth = 0)
    scatter!(classification_data.X1, classification_data.X2,
		     xlabel = "X₁", ylabel = "X₂",
             color = coerce(classification_data.y, Count), label = nothing)
end

# ╔═╡ 2ae86454-1877-4972-9cf6-24ef9350a296
md"The small dots in the background indicate how logistic regression would
would classify test data (with decision threshold 0.5; blue = false, red = true)
and the large points are the training data."

# ╔═╡ d487fcd9-8b45-4237-ab2c-21f82ddf7f7c
md"## Polynomial Regression

One way to increase the flexibility is to fit a polynomial. This can be achieved by running linear regression on a transformed data set. The transformation consists of computing higher powers of the original input. In this course we do this transformation with the static transformation machine `Polynomial`.
"

# ╔═╡ c50ed135-68d5-43bc-9c26-9c265702a1f0
MLJ.transform(machine(Polynomial(degree = 5)), regression_data)

# ╔═╡ 75c5471d-926d-4b05-9002-28b14e1dd428
md"To create a machine that transforms the input to a polynomial representation and run linear regression on this transformed representation, we use an `MLJ.Pipeline`.
Pipelines in MLJ can be written explicitely with `Pipeline(Polynomial(degree = 5), LinearRegressor())` or with the pipe operator `|>`, i.e. `Polynomial(degree = 5) |> LinearRegressor()`.
"

# ╔═╡ 710d3104-e197-44c1-a10b-de1098d57dd6
md"degree = $(@bind degree Slider(1:17, default = 4, show_value = true))"

# ╔═╡ ad50d244-c644-4f61-bd8b-995d0110811d
m3 = machine(Polynomial(degree = degree) |> LinearRegressor(),
             select(regression_data, Not(:y)),
             regression_data.y);

# ╔═╡ e9b0ea86-a5f5-43fa-aa16-5a0240f298dd
fit!(m3, verbosity = 0);

# ╔═╡ 0272460a-5b9f-4728-a531-2b497b26c512
function compute_loss(mach, data, lossfunc; operation = predict)
    pred = operation(mach, select(data, Not(:y)))
    lossfunc(pred, data.y)
end;

# ╔═╡ cfcb8f61-af91-40dd-951a-09e8dbf17e30
begin
    mse(ŷ, y) = mean((ŷ .- y).^2)
    regression_test_data = regression_data_generator(n = 10^4)
    plosses = hcat([let m = fit!(machine(Polynomial(degree = d) |> LinearRegressor(),
                                         select(regression_data, Not(:y)),
                                         regression_data.y), verbosity = 0)
                        [compute_loss(m, regression_data, mse),
                         compute_loss(m, regression_test_data, mse)]
                   end
                   for d in 1:17]...)
end;

# ╔═╡ 6fa9b644-d4a6-4c53-9146-9d978207bfd0
begin
    scatter(regression_data.x, regression_data.y, label = "training data")
    p5 = plot!(0:.01:1, predict(m3, (x = 0:.01:1,)), w = 3,
		       xlabel = "X", ylabel = "Y",
               label = "$degree-polynomial regression", legend = :topleft)
	plot!(f, color = :green, label = "data generator", w = 2)
    p6 = plot(1:17, plosses[1, :], color = :blue, label = "training loss")
    plot!(1:17, plosses[2, :], color = :red, label = "test loss")
    hline!([minimum(plosses[2, :])], color = :red,
           label = "minimal test loss", style = :dash)
    hline!([.1^2], color = :green, style = :dash, label = "irreducible error")
    scatter!([degree], [plosses[1, degree]], color = :blue, label = nothing)
    scatter!([degree], [plosses[2, degree]], color = :red, label = nothing)
    vline!([argmin(plosses[2, :])], color = :red,
                label = nothing, style = :dash, xlabel = "flexibility (degree)", ylabel = "mse")
    plot(p5, p6, layout = (1, 2), size = (700, 400))
end

# ╔═╡ edfb269d-677e-4687-8bff-0aa9ae6e64c3
md"## Polynomial Classification"

# ╔═╡ 5ea1b31d-91e5-4c8f-93d6-5d31816fdbf5
MLJ.transform(machine(Polynomial(degree = 2, predictors = (:X1, :X2))),
	          classification_data)

# ╔═╡ 59acced5-16eb-49b8-8cf2-0c43a88d838e
md"degree = $(@bind degree2 Slider(1:17, default = 3, show_value = true))"

# ╔═╡ 2fa54070-e261-462d-bd63-c225b92fa876
m4 = machine(Polynomial(degree = degree2, predictors = (:X1, :X2)) |> LogisticClassifier(penalty = :none),
             select(classification_data, Not(:y)),
             classification_data.y);

# ╔═╡ e0acbf00-f6de-483b-902b-31db99298da7
fit!(m4, verbosity = 0);

# ╔═╡ 16f0d1b3-bd97-407d-9a79-25b0fb05bbeb
begin
    classification_test_data = classification_data_generator(n = 10^4)
    cplosses = hcat([let m = fit!(machine(Polynomial(degree = d,
                                                     predictors = (:X1, :X2)) |>
                                          LogisticClassifier(penalty = :none),
                                          select(classification_data, Not(:y)),
                                          classification_data.y), verbosity = 0)
                         [mean(compute_loss(m, classification_data, log_loss)),
                          mean(compute_loss(m, classification_test_data, log_loss))]
                   end
                   for d in 1:17]...)
    c_irred_error = let data = classification_test_data,
                        p = σ.(20(f.(data.X1) .- data.X2))
        mean(log_loss(UnivariateFinite([false, true], p,
			                            augment = true, pool = missing),
			 data.y))
    end
end;

# ╔═╡ ed62cb94-8187-4d50-a6f9-6967893dd021
begin
    scatter(xgrid.X1, xgrid.X2, color = coerce(predict_mode(m4, xgrid), Count),
            markersize = 2, label = nothing, markerstrokewidth = 0,
            xlabel = "X1", ylabel = "X2")
    p7 = scatter!(classification_data.X1, classification_data.X2,
                  color = coerce(classification_data.y, Count), label = nothing)
    plot(1:17, cplosses[1, :], color = :blue, label = "training loss")
    plot!(1:17, cplosses[2, :], color = :red, label = "test loss")
    hline!([c_irred_error], color = :green, style = :dash, label = "irreducible error")
    hline!([minimum(cplosses[2, :])], color = :red,
           label = "minimal test loss", style = :dash)
    scatter!([degree2], [cplosses[1, degree2]], color = :blue, label = nothing)
    scatter!([degree2], [cplosses[2, degree2]], color = :red, label = nothing)
    p8 = vline!([argmin(cplosses[2, :])], color = :red,
                label = nothing, style = :dash, xlabel = "flexibility (degree)", ylabel = "negative loglikelihood")
    plot(p7, p8, layout = (1, 2), size = (700, 400))
end

# ╔═╡ 5c984615-f123-47fc-8330-66694ab1cb9f
md"# K-Nearest-Neighbor Regression

K = $(@bind K Slider(1:50, show_value = true))
"

# ╔═╡ 12942f82-efb1-11eb-2827-df957759b02c
m12 = machine(KNNRegressor(K = K),
	          select(regression_data, :x),
	          regression_data.y);

# ╔═╡ e1477620-dc57-4e9a-b342-8798cf6aeffe
fit!(m12);

# ╔═╡ 8ba77b77-1016-4f5d-9f9e-76b2ad1f9eac
begin
    losses = hcat([let m = fit!(machine(KNNRegressor(K = k),
                                        select(regression_data, Not(:y)),
                                        regression_data.y), verbosity = 0)
                       [compute_loss(m, regression_data, mse),
                        compute_loss(m, regression_test_data, mse)]
                   end
                   for k in 1:50]...)
end;

# ╔═╡ 542327e9-1599-4a14-805f-5d2a3c4eae14
begin
    scatter(regression_data.x, regression_data.y, label = "training data")
    p1 = plot!(0:.001:1, predict(m12, (x = 0:.001:1,)), w = 3,
               label = "$K-NN fit", legend = :topleft, xlabel = "x", ylabel = "y")
    plot(1 ./ (1:50), losses[1, :], color = :blue, label = "training loss")
    plot!(1 ./ (1:50), losses[2, :], color = :red, label = "test loss")
    hline!([minimum(losses[2, :])], color = :red,
           label = "minimal test loss", style = :dash)
    hline!([.1^2], color = :green, style = :dash, label = "irreducible error")
    scatter!([1 / K], [losses[1, K]], color = :blue, label = nothing)
    scatter!([1 / K], [losses[2, K]], color = :red, label = nothing)
    p2 = vline!([1 / argmin(losses[2, :])], color = :red,
                label = nothing, style = :dash, xlabel = "flexibility 1/K", ylabel = "mse")
    plot(p1, p2, layout = (1, 2), size = (700, 400))
end

# ╔═╡ 12942f8c-efb1-11eb-284c-393f6a694818
md"# K-Nearest-Neighbor Classification

K = $(@bind Kc Slider(1:100, show_value = true))"

# ╔═╡ 12942fc8-efb1-11eb-3180-dff1921c5bf9
m14 = machine(KNNClassifier(K = Kc),
             select(classification_data, Not(:y)),
             classification_data.y);

# ╔═╡ 87468757-49c0-474f-8bea-bd5b37d10161
fit!(m14);

# ╔═╡ 0a57f15b-c292-4c64-986d-f046260da66e
begin
    closses = hcat([let m = fit!(machine(KNNClassifier(K = k),
                                         select(classification_data, Not(:y)),
                                         classification_data.y), verbosity = 0)
                        [mean(compute_loss(m, classification_data, log_loss)),
                         mean(compute_loss(m, classification_test_data, log_loss))]
                   end
                   for k in 1:100]...)
end;

# ╔═╡ 12942fc8-efb1-11eb-0c02-0150ef55ae98
begin
    scatter(xgrid.X1, xgrid.X2, color = coerce(predict_mode(m14, xgrid), Count),
            markersize = 2, label = nothing, markerstrokewidth = 0,
            xlabel = "X1", ylabel = "X2")
    p3 = scatter!(classification_data.X1, classification_data.X2,
                  color = coerce(classification_data.y, Count), label = nothing)
    plot(1 ./ (1:100), closses[1, :], color = :blue, label = "training loss")
    plot!(1 ./ (1:100), closses[2, :], color = :red, label = "test loss")
    hline!([c_irred_error], color = :green, style = :dash, label = "irreducible error")
    hline!([minimum(closses[2, :])], color = :red,
           label = "minimal test loss", style = :dash)
    scatter!([1 / Kc], [closses[1, Kc]], color = :blue, label = nothing)
    scatter!([1 / Kc], [closses[2, Kc]], color = :red, label = nothing)
    p4 = vline!([1 / argmin(closses[2, :])], color = :red,
                label = nothing, style = :dash, xlabel = "flexibility 1/K", ylabel = "negative loglikelihood")
    plot(p3, p4, layout = (1, 2), size = (700, 400))
end

# ╔═╡ cc8ed1de-beab-43e5-979e-e83df23f96ae
md"## Application to Handwritten Digit Recognition (MNIST)

In the following we are fitting a first nearest-neigbor classifier to the MNIST
data set.

WARNING: The following code takes more than 10 minutes to run.
Especially the prediction is slow, because for every test image the closest out
of 60'000 training images has to be found.

```julia
using OpenML
mnist_x, mnist_y = let df = OpenML.load(554) |> DataFrame
    coerce!(df, :class => Multiclass)
    coerce!(df, Count => Continuous)
    df[:, 1:end-1] ./ 255,
    df.class
end
m5 = fit!(machine(KNNClassifier(K = 1), mnist_x[1:60000, :], mnist_y[1:60000]))
mnist_errorrate = mean(predict_mode(m5, mnist_x[60001:70000, :]) .!= mnist_y[60001:70000])
```
We find a misclassification rate of approximately 3%. This is clearly better than the
approximately 7.8% obtained with Multinomial Logistic Regression.
"

# ╔═╡ f6093c98-7e89-48ba-95c9-4d1f60a25033
md"# Bias-Variance Decomposition

Below we define a function that fits a polynomial regression of a given `degree` to some `training_data` and evaluates it on some `test_data`. The result is returned in form of a `DataFrame` with a single row.
"

# ╔═╡ 376d62fc-2859-422d-9944-bd9c7929f942
function fit_and_evaluate(degree, training_data, test_data)
    m = fit!(machine(Polynomial(; degree) |> LinearRegressor(),
                     select(training_data, Not(:y)),
                     training_data.y), verbosity = 0)
    ŷ = predict(m, select(test_data, Not(:y)))
    DataFrame(degree = degree,
              training_loss = compute_loss(m, training_data, mse),
              test_loss = compute_loss(m, test_data, mse),
              prediction = Ref(ŷ)) # we use Ref to store the reference to the vector of predictions instead of inserting the predicted values as different rows into the DataFrame.
end;

# ╔═╡ 021f812a-b52a-46a3-b81a-fa1bfefff295
md"Now we will run this function on 9 different polynomial degrees and 100 different training sets. We use always the same (large) test set."

# ╔═╡ 1fd15c67-a196-4324-94f6-9e4ccdc92482
results = vcat([fit_and_evaluate(degree,
                         regression_data_generator(n = 50,
                                                   rng = MersenneTwister(seed)),
                         regression_test_data)
                for degree in 1:9, seed in 1:100]...)

# ╔═╡ bdaa48d6-30bb-4d04-b5dd-cc1825cb69e3
result = combine(groupby(results, :degree),
	             :training_loss => mean,
                 :test_loss => mean,
                 :prediction => (x -> mean(var(x))) => :variance,
                 :prediction => (x -> mean((mean(x) .- f.(regression_test_data.x)).^2)) => :bias)

# ╔═╡ 5ac9bcb4-c5a5-4b53-bb45-73dde74b7f60
begin
    plot(result.degree, result.training_loss_mean,
         label = "training loss", w = 2, ylabel = "loss",
         xlabel = "flexibility (degree)")
	plot!(result.degree, result.test_loss_mean, label = "test loss", w = 2)
    plot!(result.degree, result.bias, label = "bias^2", w = 2)
    plot!(result.degree, result.variance, label = "variance", w = 2)
    plot!(result.degree, result.variance .+ result.bias .+ .1^2, label = "bias^2 + variance + Var(ϵ)", color = :black, w = 2, linestyle = :dash)
    hline!([.1^2], label = "Var(ϵ)")
end

# ╔═╡ ae86ee9c-3645-43d2-9a42-79658521c3fb
md"# Exercises

## Conceptual
1. For each of example below, indicate whether we would generally expect the performance of a flexible statistical learning method to be better or worse than an inflexible method. Justify your answer.
   + The sample size n is extremely large, and the number of predictors ``p`` is small.
   + The number of predictors ``p`` is extremely large, and the number of observations ```n`` is small.
   + The relationship between the predictors and response is highly non-linear.
   + The variance of the error terms, i.e. ``\sigma^2 = Var(\epsilon)``, is extremely high.
1. The table below provides a training data set containing six observations, three predictors, and one qualitative response variable. Suppose we wish to use this data set to make a prediction for ``Y`` when ``X_1 = X_2 = X_3 = 0`` using K-nearest neighbors.
   - Compute the Euclidean distance between each observation and the test point, X1 = X2 = X3 = 0.
   - What is our prediction with K = 1? Why?
   - What is our prediction with K = 3? Why?
   - If the decision boundary at threshold 0.5 in this problem is highly non- linear, then would we expect the best value for K to be large or small? Why?
|Obs.|``X_1``|``X_2``|``X_3``|``Y``|
|:----|-------|-------|-------|-----:|
|1 | 0 | 3 | 0 | Red|
|2 | 2 | 0 | 0 | Red|
|3 | 0 | 1 | 3 | Red|
|4 | 0 | 1 | 2 | Green|
|5 | -1| 0 | 1 | Green|
|6 | 1 | 1 | 1 | Red|
3. Try without looking at the figure above:
   - Provide a sketch of typical (squared) bias, variance, training error, test error, and Bayes (or irreducible) error curves, on a single plot, as we go from less flexible statistical learning methods towards more flexible approaches. The x-axis should represent the amount of flexibility in the method, and the y-axis should represent the values for each curve. There should be five curves. Make sure to label each one.
   - Explain why each of the five curves has its particular shape.
2. Suppose that we take a data set with mutually distinct inputs ``x_i\neq x_j`` for ``i\neq j``, divide it into equally-sized training and test sets, and then try out two different classification procedures. First we use logistic regression and get an error rate of 20% on the training data and 30% on the test data. Next we use 1-nearest neighbors (i.e. ``K = 1``) and get an average error rate (averaged over both test and training data sets) of 18%. Based on these results, which method should we prefer to use for classification of new observations? Why?

## Applied
1. Apply K-nearest neighbors regression to the weather data set. Use as input all predictors except `:time` and `:LUZ_wind_peak`.
    * Compute the training and the test loss for ``K = 5, 10, 20, 50, 100``.
    * Which value of the hyper-parameter ``K`` should we prefer to make predictions on new data?
    * Should we prefer K-nearest neighbors with optimal ``K`` or multiple linear regression to make predictions on new data? *Hint*: Remember that we found a training error (RMSE) of approximately 8.1 and a test error of 8.9.
2. In this exercise we review the error-decomposition and the bias-variance decomposition.
* Write a data generator where the mean of the output depends through the non-linear function ``f(x) = x^2 * \sin(x) + 4 * \tanh(10x)`` on the input and normally distributed noise ``\epsilon`` with mean 0 and standard deviation 1.5.
    * Take the linear function ``\hat f(x) = 2x`` and estimate its reducible error at input point ``x = 0`` and at input point ``x = 2`` in two ways:
        * Using directly ``f``.
        * Using ``10^5`` samples from the data generator and your knowledge about the irreducible error.
    * Generate ``10^4`` training sets of 100 data points with input ``x`` normally distributed with standard deviation 2 and mean 0 and estimate the bias of linear regression at ``x = 4``  in two ways:
        * Using directly ``f``.
        * Using ``10^4`` samples from the data generator, your knowledge about the irreducible error and your estimate of the variance of linear regression.
"

# ╔═╡ efda845a-4390-40bc-bdf2-89e555d3b1b2
MLCourse.list_notebooks(@__FILE__)

# ╔═╡ 2320f424-7652-4e9f-83ef-fc011b722dcc
MLCourse.footer()

# ╔═╡ Cell order:
# ╟─12942f50-efb1-11eb-01c0-055b6be166e0
# ╠═8fa836a6-1133-4a54-b996-a02083fc6bba
# ╠═12942f5a-efb1-11eb-399a-a1300d636217
# ╠═12942f62-efb1-11eb-3f59-f981bc32f308
# ╠═12942f6e-efb1-11eb-0a49-01a6a2d0196f
# ╟─12942f6e-efb1-11eb-35c0-19a7c4d6b44e
# ╠═12942f94-efb1-11eb-234e-51492b51e583
# ╠═12942fa0-efb1-11eb-01c8-6dae80c55fb8
# ╟─12942f94-efb1-11eb-2c48-a3418b53b886
# ╟─2ae86454-1877-4972-9cf6-24ef9350a296
# ╟─d487fcd9-8b45-4237-ab2c-21f82ddf7f7c
# ╠═c50ed135-68d5-43bc-9c26-9c265702a1f0
# ╟─75c5471d-926d-4b05-9002-28b14e1dd428
# ╟─710d3104-e197-44c1-a10b-de1098d57dd6
# ╠═ad50d244-c644-4f61-bd8b-995d0110811d
# ╠═e9b0ea86-a5f5-43fa-aa16-5a0240f298dd
# ╟─6fa9b644-d4a6-4c53-9146-9d978207bfd0
# ╠═0272460a-5b9f-4728-a531-2b497b26c512
# ╠═cfcb8f61-af91-40dd-951a-09e8dbf17e30
# ╟─edfb269d-677e-4687-8bff-0aa9ae6e64c3
# ╠═5ea1b31d-91e5-4c8f-93d6-5d31816fdbf5
# ╟─59acced5-16eb-49b8-8cf2-0c43a88d838e
# ╠═2fa54070-e261-462d-bd63-c225b92fa876
# ╠═e0acbf00-f6de-483b-902b-31db99298da7
# ╟─ed62cb94-8187-4d50-a6f9-6967893dd021
# ╠═16f0d1b3-bd97-407d-9a79-25b0fb05bbeb
# ╟─5c984615-f123-47fc-8330-66694ab1cb9f
# ╠═269a609c-74af-4e7e-86df-e2279096a7a6
# ╠═12942f82-efb1-11eb-2827-df957759b02c
# ╠═e1477620-dc57-4e9a-b342-8798cf6aeffe
# ╟─542327e9-1599-4a14-805f-5d2a3c4eae14
# ╟─8ba77b77-1016-4f5d-9f9e-76b2ad1f9eac
# ╟─12942f8c-efb1-11eb-284c-393f6a694818
# ╠═12942fc8-efb1-11eb-3180-dff1921c5bf9
# ╠═87468757-49c0-474f-8bea-bd5b37d10161
# ╟─12942fc8-efb1-11eb-0c02-0150ef55ae98
# ╟─0a57f15b-c292-4c64-986d-f046260da66e
# ╟─cc8ed1de-beab-43e5-979e-e83df23f96ae
# ╟─f6093c98-7e89-48ba-95c9-4d1f60a25033
# ╠═376d62fc-2859-422d-9944-bd9c7929f942
# ╟─021f812a-b52a-46a3-b81a-fa1bfefff295
# ╠═1fd15c67-a196-4324-94f6-9e4ccdc92482
# ╠═bdaa48d6-30bb-4d04-b5dd-cc1825cb69e3
# ╟─5ac9bcb4-c5a5-4b53-bb45-73dde74b7f60
# ╟─ae86ee9c-3645-43d2-9a42-79658521c3fb
# ╟─efda845a-4390-40bc-bdf2-89e555d3b1b2
# ╟─12942f34-efb1-11eb-3eb4-c1a38396cfb8
# ╟─2320f424-7652-4e9f-83ef-fc011b722dcc

max_stars_repo_name: baltisberger/MLCourse

max_issues_repo_name: baltisberger/MLCourse

max_forks_repo_name: baltisberger/MLCourse


text:
using Dates

const filename = "NOTES.TXT"

if length(ARGS) == 0
    fp = open(filename, "r")
    println(read(fp, String))
else
    fp = open(filename, "a+")
    write(fp, string(DateTime(now()), "\n\t", join(ARGS, " "), "\n"))
end
close(fp)

ext: jl

max_stars_repo_name: ethansaxenian/RosettaDecode

max_issues_repo_name: ethansaxenian/RosettaDecode

max_forks_repo_name: ethansaxenian/RosettaDecode


text:
#DA algorythm

#recursive version
function call_match(m, n, m_prefs, f_prefs)
  males = 1:m
  females = 1:n

  pointer_males = ones(m)
  temp_matched_males = zeros(n)

  da_match(m, n, m_prefs, f_prefs, pointer_males, temp_matched_males)
end

function da_match(m, n, m_prefs, f_prefs, pointer_males, temp_matched_males)#m offers to f
  candidate_males = Array{Int32}(m, n)
  for i in 1:m
    push!(candidate_males[m_prefs[i, pointer_males[i]]], i)
  end

  choose_best_male!(n, f_prefs, candidate_males, temp_matched_males)
  if all(pointer_males .== n)
    return temp_matched_males
  else
    da_match(m, n, m_prefs, f_prefs, pointer_males, temp_matched_males)
  end
end

function choose_best_male(n, f_prefs, candidate_males, temp_matched_males)
  for j in 1:n
    for i in f_prefs[j]
      if i in candidate_males[:, j]
        temp_matched_males[j] = i
        break
      end
    end
    temp_matched_males[j] = 0
  end
end

call_match(2, 2, [1 2; 2 1], [1 2; 2 1])

max_stars_repo_path: Julia/Untitled Folder/2-DA.jl

max_issues_repo_path: Julia/Untitled Folder/2-DA.jl

max_forks_repo_path: Julia/Untitled Folder/2-DA.jl


text:
Base.show(io::IO, F::AbstractFiber) = print(io, typeof(F), "()")

function Base.show(io::IO, mime::MIME"text/html", NP::NewtonPolygon)
    if inIJulia() && plotsdefined()
        show(io, mime, Main.Plots.plot(NP))
    else
        show(io, NP)
    end
end
function Base.show(io::IO, NP::NewtonPolygon)
    #println(typeof(A), ":")
    println(io, "Newton polygon")
    println(io, "* lattices: ", lattices(NP))
    println(io, "* vertices: ", vertices(NP))
end

function Base.show(io::IO, B::Bitmap2D)
    n, m = size(B)
    print(io, "Bitmap2D of size ", n, "×", m)
end
function Base.show(io::IO, B::Bitmap3D)
    n1, n2, n3 = size(B)
    print(io, "Bitmap2D of size ", n1, "×", n2, "×", n3)
end
function Base.show(io::IO, mime::MIME"text/html", M::AbstractBitmap)
    if inIJulia() && plotsdefined()
        show(io, mime, Main.Plots.plot(M))
    else
        show(io, M)
    end
end


text:
using BinaryBuilder

name = "PETSc"
version = v"3.13.4"

# Collection of sources required to build PETSc. Avoid using the git repository, it will
# require building SOWING which fails in all non-linux platforms.
sources = [
    ArchiveSource("https://www.mcs.anl.gov/petsc/mirror/release-snapshots/petsc-3.13.4.tar.gz",
    "2cf8abd70ec332510b17bff7dc8a465bf1eb053e17e2a0451a4a0256385657d8"),
    DirectorySource("./bundled"),
]

# Bash recipe for building across all platforms
script = raw"""
cd $WORKSPACE/srcdir/petsc*
includedir="${prefix}/include"
atomic_patch -p1 $WORKSPACE/srcdir/patches/petsc_name_mangle.patch

BLAS_LAPACK_LIB="${libdir}/libopenblas.${dlext}"

if [[ "${target}" == *-mingw* ]]; then
    #atomic_patch -p1 $WORKSPACE/srcdir/patches/fix-header-cases.patch
    MPI_LIBS="${libdir}/msmpi.${dlext}"
else
    MPI_LIBS="[${libdir}/libmpi.${dlext},libmpifort.${dlext}]"
fi

./configure --prefix=${prefix} \
    CC=${CC} \
    FC=${FC} \
    CXX=${CXX} \
    COPTFLAGS='-O3' \
    CXXOPTFLAGS='-O3' \
    FOPTFLAGS='-O3' \
    --with-64-bit-indices=0 \
    --with-debugging=0 \
    --with-batch \
    --PETSC_ARCH=$target \
    --with-blaslapack-lib=$BLAS_LAPACK_LIB \
    --with-blaslapack-suffix="" \
    --known-64-bit-blas-indices=0 \
    --with-mpi-lib="${MPI_LIBS}" \
    --known-mpi-int64_t=0 \
    --with-mpi-include="${includedir}" \
    --with-sowing=0

if [[ "${target}" == *-mingw* ]]; then
    export CPPFLAGS="-Dpetsc_EXPORTS"
elif [[ "${target}" == powerpc64le-* ]]; then
    export CFLAGS="-fPIC"
    export FFLAGS="-fPIC"
fi

make -j${nproc} \
    PETSC_DIR="${PWD}" \
    PETSC_ARCH="${target}" \
    CPPFLAGS="${CPPFLAGS}" \
    CFLAGS="${CFLAGS}" \
    FFLAGS="${FFLAGS}" \
    DEST_DIR="${prefix}" \
    all

make PETSC_DIR=$PWD PETSC_ARCH=$target DEST_DIR=$prefix install

if [[ "${target}" == *-mingw* ]]; then
    # Move library to ${libdir} on Windows,
    # changing the extension from so to dll.
    mv ${prefix}/lib/libpetsc.so.*.*.* "${libdir}/libpetsc.${dlext}"
    # Remove useless links
    rm ${prefix}/lib/libpetsc.*
fi
"""

# We attempt to build for all defined platforms
platforms = expand_gfortran_versions(supported_platforms(exclude=[Windows(:i686)]))

products = [
    LibraryProduct("libpetsc", :libpetsc),
]

dependencies = [
    Dependency("OpenBLAS32_jll"),
    Dependency("MPICH_jll"),
    Dependency("MicrosoftMPI_jll"),
    Dependency("CompilerSupportLibraries_jll"),
]

# Build the tarballs.
build_tarballs(ARGS, name, version, sources, script, platforms, products, dependencies; preferred_gcc_version = v"6")


text:
using Delta
import Delta.paramsof
import Delta.forward
import Delta.predict
import Base.length
import Base.getindex

# 1. define model struct and its length and indexing method
mutable struct Model
    blocks::Vector
    function Model(featdims::Int)
        c1 = PlainConv1d(featdims,256,8,stride=3)
        c2 = PlainConv1d(256,256,3,stride=2)
        c3 = PlainConv1d(256,512,3,stride=1)
        c4 = PlainConv1d(512,512,3,stride=1)
        c5 = PlainConv1d(512,512,3,stride=1)
        r1 = IndRNN(512,512)
        r2 = IndRNN(512,512)
        f1 = Linear(512,140)
        new([c1,c2,c3,c4,c5, Chain(r1,r2,f1)])
    end
end

Base.length(m::Model) = length(m.blocks)
Base.getindex(m::Model, k...) = m.blocks[k...]

# 2. define how to extrac model's params
function paramsof(m::Model)
    params = Vector{Variable}(undef,0)
    for i = 1:length(m)
        append!(params, paramsof(m[i]))
    end
    return params
end

# 3. define model's forward calculation
function forward(m::Model, input::Variable)
    x1 = forward(m[1], input)
    x1 = relu(x1)
    x2 = forward(m[2], x1)
    x2 = relu(x2)
    x3 = forward(m[3], x2)
    x3 = relu(x3)
    x4 = forward(m[4], x3)
    x4 = relu(x4)
    x5 = forward(m[5], x4)
    x5 = relu(x5)
    return PackedSeqForward(m[6], x5)
end

# 4. experimental constants
featdims  = 64
timesteps = 124
batchsize = 16

# 5. pseudo inputs & labels
x = Variable(randn(featdims,timesteps,batchsize));
seqlabels = Vector(undef,batchsize)
for i=1:batchsize
    seqlabels[i] = [2, 3, 4, 5]
end

# 6. instantiate an model object
#    extrac its parameters
#    choose an optimizer
asr = Model(featdims);
asrparams = paramsof(asr);
optim = Momentum(asrparams);
epoch = 37

# 7. train this model
for e = 1:epoch
    y = forward(asr, x)
    loglikely = CRNN_Batch_CTC_With_Softmax(y, seqlabels);
    println("iter $e loss = $loglikely");
    backward()
    update!(optim, asrparams, clipvalue=20.0)
    zerograds!(asrparams)
end

# 8. define predict calculation
function predict(m::Model, input::AbstractArray)
    x1 = predict(m[1], input)
    x1 = relu(x1)
    x2 = predict(m[2], x1)
    x2 = relu(x2)
    x3 = predict(m[3], x2)
    x3 = relu(x3)
    x4 = predict(m[4], x3)
    x4 = relu(x4)
    x5 = predict(m[5], x4)
    x5 = relu(x5)
    return PackedSeqPredict(m[6], x5)
end

# 9. predict
y = predict(asr, randn(featdims,timesteps,batchsize));
r = CTCGreedySearch(softmax(y[:,:,10],dims=1));
println("decoding result: ", r);

# for more humanlity iteration style
Base.lastindex(m::Model)  = length(m)
Base.firstindex(m::Model) = 1
Base.iterate(m::Model, i=firstindex(m)) = i>length(m) ? nothing : (m[i], i+1)

# so we could iterate like `for x in X`
function Delta.nparamsof(model::Model)
    nparams = 0
    for m in model
        nparams += nparamsof(m)
    end
    return nparams
end

# of course we could also calculate how many bytes it uses
function bytesof(model::Model, unit::String="MB")
    n = nparamsof(model)
    u = lowercase(unit)
    if u == "kb" return n * sizeof(eltype(model[1].w)) / 1024 end
    if u == "mb" return n * sizeof(eltype(model[1].w)) / 1048576 end
    if u == "gb" return n * sizeof(eltype(model[1].w)) / 1073741824 end
    if u == "tb" return n * sizeof(eltype(model[1].w)) / 1099511627776 end
end

max_stars_repo_name: sonosole/Delta

max_issues_repo_name: sonosole/Delta

max_forks_repo_name: sonosole/Delta


text:
let
    PerfectHeuristic()
    SoftConflictHeuristic()
end
let
    G = initialize_regular_grid_graph(;n_obstacles_x=1,n_obstacles_y=1)
    start_times = [0,-2,2]
    starts = [1,2,3]
    goals = [4,5,6]
    let
        h = PerfectHeuristic(G,starts,goals)
        @test get_heuristic_cost(h,goals[1],starts[1]) == gdistances(G,starts[1])[goals[1]]
    end
    # let
        # h = SoftConflictHeuristic(length(starts),num(G),start_times,starts,goals)
        # @test get_heuristic_cost(h,starts[1],1) >= 1.0
    # end
    let
        h1 = PerfectHeuristic(G,starts,goals)
        # h2 = SoftConflictHeuristic(G,start_times,starts,goals)
    end
    let
        h = construct_composite_heuristic(
            NullHeuristic(),
            PerfectHeuristic(G,starts,goals)
        )
        cost = get_heuristic_cost(h,goals[1],starts[1])
        @test cost[1] == 0.0
        @test cost[2] ==  gdistances(G,starts[1])[goals[1]]
    end
end
# DeadlineHeuristic
# let
#     G = Graph(3)
#     for v in 1:nv(G)-1
#         add_edge!(G,v,v+1)
#     end
#     starts = [1]
#     goals = [nv(G)]
#     h1 = PerfectHeuristic(G,starts,goals)
#     h = DeadlineHeuristic(t_max=ne(G)+1.0, h=h1)
#     @test get_heuristic_cost(h, 0.0, goals[1], starts[1]) == 0.0
#     @test get_heuristic_cost(h, 1.0, goals[1], starts[1]) == 0.0
#     @test get_heuristic_cost(h, 2.0, goals[1], starts[1]) == 1.0
# end
# HardConflictHeuristic
let
    G = Graph(5)
    for v in 1:nv(G)-1
        add_edge!(G,v,v+1)
    end
    T = 10
    num_agents = 2
    h = HardConflictHeuristic(G,T,num_agents)
    set_path!(h,num_agents,[1,2],2)
    # Agent 1 has a cost of 1 whenever it conflicts with agent 2's path
    @test get_heuristic_cost(h,1,1,2) == 1.0
    @test get_heuristic_cost(h,1,2,3) == 1.0
    @test get_heuristic_cost(h,1,2,1) == 0.0
    @test get_heuristic_cost(h,1,3,4) == 0.0
    # Agent 2 should have zero cost even when it conflicts with its own path
    @test get_heuristic_cost(h,2,1,2) == 0.0
    # Now that agent 2's path has been updated, the cost should be zero for a
    # vtx no longer on the path
    set_path!(h,num_agents,[1,2],3)
    @test get_heuristic_cost(h,1,1,2) == 0.0
end
# Full CompositeCost and CompositeHeuristic for TaskGraphs route planner
let
    G = Graph(3)
    for v in 1:nv(G)-1
        add_edge!(G,v,v+1)
    end
    start_times = [0]
    starts = [1]
    goals = [nv(G)]
    deadline = 10
    num_agents = 1
    cost_model = construct_composite_cost_model(
        FullDeadlineCost(DeadlineCost(deadline)),
        FullCostModel(sum,NullCost()),
        SumOfTravelTime()
    )
    heuristic_model = construct_composite_heuristic(
        PerfectHeuristic(G,starts,goals),
        HardConflictHeuristic(G,ne(G),num_agents),
        PerfectHeuristic(G,starts,goals),
    )
end
let
    mapf = init_mapf_1()
    env = build_env(mapf,initialize_root_node(mapf),1)
    d = get_distance(env,get_start(mapf,1),get_goal(mapf,1))
    c = get_heuristic_cost(EnvDistanceHeuristic(),env,get_start(mapf,1))
    @test d == c
end
# multi stage env distance
let
    mapf = init_mapf_1()
    env = build_env(mapf,initialize_root_node(mapf),1)
    goals = Dict(1=>[1,2,3],2=>[3,4,5])

    h = CRCBS.construct_multi_stage_env_distance_heuristic(env.graph,goals)

    for (idx,goal_sequence) in goals
        for stage in 1:length(goal_sequence)
            d = 0
            for (i,v) in enumerate(goal_sequence[stage:end-1])
                d += get_distance(env.graph,v,goal_sequence[stage+i])
            end
            @test CRCBS.cost_from_stage(h,idx,stage) == d
        end
    end
end


text:
using CEnum

#! format: off

const CeedInt = Int32

const CeedScalar = Cdouble

mutable struct Ceed_private end

const Ceed = Ptr{Ceed_private}

mutable struct CeedRequest_private end

const CeedRequest = Ptr{CeedRequest_private}

mutable struct CeedVector_private end

const CeedVector = Ptr{CeedVector_private}

mutable struct CeedElemRestriction_private end

const CeedElemRestriction = Ptr{CeedElemRestriction_private}

mutable struct CeedBasis_private end

const CeedBasis = Ptr{CeedBasis_private}

mutable struct CeedQFunction_private end

const CeedQFunction = Ptr{CeedQFunction_private}

mutable struct CeedQFunctionContext_private end

const CeedQFunctionContext = Ptr{CeedQFunctionContext_private}

mutable struct CeedOperator_private end

const CeedOperator = Ptr{CeedOperator_private}

function CeedRegistryGetList(n, resources, array)
    ccall((:CeedRegistryGetList, libceed), Cint, (Ptr{Csize_t}, Ptr{Ptr{Ptr{Cchar}}}, Ptr{Ptr{CeedInt}}), n, resources, array)
end

function CeedInit(resource, ceed)
    ccall((:CeedInit, libceed), Cint, (Ptr{Cchar}, Ptr{Ceed}), resource, ceed)
end

function CeedReferenceCopy(ceed, ceed_copy)
    ccall((:CeedReferenceCopy, libceed), Cint, (Ceed, Ptr{Ceed}), ceed, ceed_copy)
end

function CeedGetResource(ceed, resource)
    ccall((:CeedGetResource, libceed), Cint, (Ceed, Ptr{Ptr{Cchar}}), ceed, resource)
end

function CeedIsDeterministic(ceed, is_deterministic)
    ccall((:CeedIsDeterministic, libceed), Cint, (Ceed, Ptr{Bool}), ceed, is_deterministic)
end

function CeedView(ceed, stream)
    ccall((:CeedView, libceed), Cint, (Ceed, Ptr{Libc.FILE}), ceed, stream)
end

function CeedDestroy(ceed)
    ccall((:CeedDestroy, libceed), Cint, (Ptr{Ceed},), ceed)
end

# automatic type deduction for variadic arguments may not be what you want, please use with caution
@generated function CeedErrorImpl(arg1, arg2, arg3, arg4, arg5, arg6, va_list...)
        :(@ccall(libceed.CeedErrorImpl(arg1::Ceed, arg2::Ptr{Cchar}, arg3::Cint, arg4::Ptr{Cchar}, arg5::Cint, arg6::Ptr{Cchar}; $(to_c_type_pairs(va_list)...))::Cint))
    end

# typedef int ( * CeedErrorHandler ) ( Ceed , const char * , int , const char * , int , const char * , va_list * )
const CeedErrorHandler = Ptr{Cvoid}

function CeedSetErrorHandler(ceed, eh)
    ccall((:CeedSetErrorHandler, libceed), Cint, (Ceed, CeedErrorHandler), ceed, eh)
end

function CeedGetErrorMessage(arg1, err_msg)
    ccall((:CeedGetErrorMessage, libceed), Cint, (Ceed, Ptr{Ptr{Cchar}}), arg1, err_msg)
end

function CeedResetErrorMessage(arg1, err_msg)
    ccall((:CeedResetErrorMessage, libceed), Cint, (Ceed, Ptr{Ptr{Cchar}}), arg1, err_msg)
end

function CeedGetVersion(major, minor, patch, release)
    ccall((:CeedGetVersion, libceed), Cint, (Ptr{Cint}, Ptr{Cint}, Ptr{Cint}, Ptr{Bool}), major, minor, patch, release)
end

@cenum CeedErrorType::Int32 begin
    CEED_ERROR_SUCCESS = 0
    CEED_ERROR_MINOR = 1
    CEED_ERROR_DIMENSION = 2
    CEED_ERROR_INCOMPLETE = 3
    CEED_ERROR_INCOMPATIBLE = 4
    CEED_ERROR_ACCESS = 5
    CEED_ERROR_MAJOR = -1
    CEED_ERROR_BACKEND = -2
    CEED_ERROR_UNSUPPORTED = -3
end

@cenum CeedMemType::UInt32 begin
    CEED_MEM_HOST = 0
    CEED_MEM_DEVICE = 1
end

function CeedGetPreferredMemType(ceed, type)
    ccall((:CeedGetPreferredMemType, libceed), Cint, (Ceed, Ptr{CeedMemType}), ceed, type)
end

@cenum CeedCopyMode::UInt32 begin
    CEED_COPY_VALUES = 0
    CEED_USE_POINTER = 1
    CEED_OWN_POINTER = 2
end

@cenum CeedNormType::UInt32 begin
    CEED_NORM_1 = 0
    CEED_NORM_2 = 1
    CEED_NORM_MAX = 2
end

function CeedVectorCreate(ceed, len, vec)
    ccall((:CeedVectorCreate, libceed), Cint, (Ceed, CeedInt, Ptr{CeedVector}), ceed, len, vec)
end

function CeedVectorReferenceCopy(vec, vec_copy)
    ccall((:CeedVectorReferenceCopy, libceed), Cint, (CeedVector, Ptr{CeedVector}), vec, vec_copy)
end

function CeedVectorSetArray(vec, mem_type, copy_mode, array)
    ccall((:CeedVectorSetArray, libceed), Cint, (CeedVector, CeedMemType, CeedCopyMode, Ptr{CeedScalar}), vec, mem_type, copy_mode, array)
end

function CeedVectorSetValue(vec, value)
    ccall((:CeedVectorSetValue, libceed), Cint, (CeedVector, CeedScalar), vec, value)
end

function CeedVectorSyncArray(vec, mem_type)
    ccall((:CeedVectorSyncArray, libceed), Cint, (CeedVector, CeedMemType), vec, mem_type)
end

function CeedVectorTakeArray(vec, mem_type, array)
    ccall((:CeedVectorTakeArray, libceed), Cint, (CeedVector, CeedMemType, Ptr{Ptr{CeedScalar}}), vec, mem_type, array)
end

function CeedVectorGetArray(vec, mem_type, array)
    ccall((:CeedVectorGetArray, libceed), Cint, (CeedVector, CeedMemType, Ptr{Ptr{CeedScalar}}), vec, mem_type, array)
end

function CeedVectorGetArrayRead(vec, mem_type, array)
    ccall((:CeedVectorGetArrayRead, libceed), Cint, (CeedVector, CeedMemType, Ptr{Ptr{CeedScalar}}), vec, mem_type, array)
end

function CeedVectorRestoreArray(vec, array)
    ccall((:CeedVectorRestoreArray, libceed), Cint, (CeedVector, Ptr{Ptr{CeedScalar}}), vec, array)
end

function CeedVectorRestoreArrayRead(vec, array)
    ccall((:CeedVectorRestoreArrayRead, libceed), Cint, (CeedVector, Ptr{Ptr{CeedScalar}}), vec, array)
end

function CeedVectorNorm(vec, type, norm)
    ccall((:CeedVectorNorm, libceed), Cint, (CeedVector, CeedNormType, Ptr{CeedScalar}), vec, type, norm)
end

function CeedVectorScale(x, alpha)
    ccall((:CeedVectorScale, libceed), Cint, (CeedVector, CeedScalar), x, alpha)
end

function CeedVectorAXPY(y, alpha, x)
    ccall((:CeedVectorAXPY, libceed), Cint, (CeedVector, CeedScalar, CeedVector), y, alpha, x)
end

function CeedVectorPointwiseMult(w, x, y)
    ccall((:CeedVectorPointwiseMult, libceed), Cint, (CeedVector, CeedVector, CeedVector), w, x, y)
end

function CeedVectorReciprocal(vec)
    ccall((:CeedVectorReciprocal, libceed), Cint, (CeedVector,), vec)
end

function CeedVectorView(vec, fp_fmt, stream)
    ccall((:CeedVectorView, libceed), Cint, (CeedVector, Ptr{Cchar}, Ptr{Libc.FILE}), vec, fp_fmt, stream)
end

function CeedVectorGetLength(vec, length)
    ccall((:CeedVectorGetLength, libceed), Cint, (CeedVector, Ptr{CeedInt}), vec, length)
end

function CeedVectorDestroy(vec)
    ccall((:CeedVectorDestroy, libceed), Cint, (Ptr{CeedVector},), vec)
end

function CeedRequestWait(req)
    ccall((:CeedRequestWait, libceed), Cint, (Ptr{CeedRequest},), req)
end

@cenum CeedTransposeMode::UInt32 begin
    CEED_NOTRANSPOSE = 0
    CEED_TRANSPOSE = 1
end

function CeedElemRestrictionCreate(ceed, num_elem, elem_size, num_comp, comp_stride, l_size, mem_type, copy_mode, offsets, rstr)
    ccall((:CeedElemRestrictionCreate, libceed), Cint, (Ceed, CeedInt, CeedInt, CeedInt, CeedInt, CeedInt, CeedMemType, CeedCopyMode, Ptr{CeedInt}, Ptr{CeedElemRestriction}), ceed, num_elem, elem_size, num_comp, comp_stride, l_size, mem_type, copy_mode, offsets, rstr)
end

function CeedElemRestrictionCreateStrided(ceed, num_elem, elem_size, num_comp, l_size, strides, rstr)
    ccall((:CeedElemRestrictionCreateStrided, libceed), Cint, (Ceed, CeedInt, CeedInt, CeedInt, CeedInt, Ptr{CeedInt}, Ptr{CeedElemRestriction}), ceed, num_elem, elem_size, num_comp, l_size, strides, rstr)
end

function CeedElemRestrictionCreateBlocked(ceed, num_elem, elem_size, blk_size, num_comp, comp_stride, l_size, mem_type, copy_mode, offsets, rstr)
    ccall((:CeedElemRestrictionCreateBlocked, libceed), Cint, (Ceed, CeedInt, CeedInt, CeedInt, CeedInt, CeedInt, CeedInt, CeedMemType, CeedCopyMode, Ptr{CeedInt}, Ptr{CeedElemRestriction}), ceed, num_elem, elem_size, blk_size, num_comp, comp_stride, l_size, mem_type, copy_mode, offsets, rstr)
end

function CeedElemRestrictionCreateBlockedStrided(ceed, num_elem, elem_size, blk_size, num_comp, l_size, strides, rstr)
    ccall((:CeedElemRestrictionCreateBlockedStrided, libceed), Cint, (Ceed, CeedInt, CeedInt, CeedInt, CeedInt, CeedInt, Ptr{CeedInt}, Ptr{CeedElemRestriction}), ceed, num_elem, elem_size, blk_size, num_comp, l_size, strides, rstr)
end

function CeedElemRestrictionReferenceCopy(rstr, rstr_copy)
    ccall((:CeedElemRestrictionReferenceCopy, libceed), Cint, (CeedElemRestriction, Ptr{CeedElemRestriction}), rstr, rstr_copy)
end

function CeedElemRestrictionCreateVector(rstr, lvec, evec)
    ccall((:CeedElemRestrictionCreateVector, libceed), Cint, (CeedElemRestriction, Ptr{CeedVector}, Ptr{CeedVector}), rstr, lvec, evec)
end

function CeedElemRestrictionApply(rstr, t_mode, u, ru, request)
    ccall((:CeedElemRestrictionApply, libceed), Cint, (CeedElemRestriction, CeedTransposeMode, CeedVector, CeedVector, Ptr{CeedRequest}), rstr, t_mode, u, ru, request)
end

function CeedElemRestrictionApplyBlock(rstr, block, t_mode, u, ru, request)
    ccall((:CeedElemRestrictionApplyBlock, libceed), Cint, (CeedElemRestriction, CeedInt, CeedTransposeMode, CeedVector, CeedVector, Ptr{CeedRequest}), rstr, block, t_mode, u, ru, request)
end

function CeedElemRestrictionGetCompStride(rstr, comp_stride)
    ccall((:CeedElemRestrictionGetCompStride, libceed), Cint, (CeedElemRestriction, Ptr{CeedInt}), rstr, comp_stride)
end

function CeedElemRestrictionGetNumElements(rstr, num_elem)
    ccall((:CeedElemRestrictionGetNumElements, libceed), Cint, (CeedElemRestriction, Ptr{CeedInt}), rstr, num_elem)
end

function CeedElemRestrictionGetElementSize(rstr, elem_size)
    ccall((:CeedElemRestrictionGetElementSize, libceed), Cint, (CeedElemRestriction, Ptr{CeedInt}), rstr, elem_size)
end

function CeedElemRestrictionGetLVectorSize(rstr, l_size)
    ccall((:CeedElemRestrictionGetLVectorSize, libceed), Cint, (CeedElemRestriction, Ptr{CeedInt}), rstr, l_size)
end

function CeedElemRestrictionGetNumComponents(rstr, num_comp)
    ccall((:CeedElemRestrictionGetNumComponents, libceed), Cint, (CeedElemRestriction, Ptr{CeedInt}), rstr, num_comp)
end

function CeedElemRestrictionGetNumBlocks(rstr, num_blk)
    ccall((:CeedElemRestrictionGetNumBlocks, libceed), Cint, (CeedElemRestriction, Ptr{CeedInt}), rstr, num_blk)
end

function CeedElemRestrictionGetBlockSize(rstr, blk_size)
    ccall((:CeedElemRestrictionGetBlockSize, libceed), Cint, (CeedElemRestriction, Ptr{CeedInt}), rstr, blk_size)
end

function CeedElemRestrictionGetMultiplicity(rstr, mult)
    ccall((:CeedElemRestrictionGetMultiplicity, libceed), Cint, (CeedElemRestriction, CeedVector), rstr, mult)
end

function CeedElemRestrictionView(rstr, stream)
    ccall((:CeedElemRestrictionView, libceed), Cint, (CeedElemRestriction, Ptr{Libc.FILE}), rstr, stream)
end

function CeedElemRestrictionDestroy(rstr)
    ccall((:CeedElemRestrictionDestroy, libceed), Cint, (Ptr{CeedElemRestriction},), rstr)
end

@cenum CeedEvalMode::UInt32 begin
    CEED_EVAL_NONE = 0
    CEED_EVAL_INTERP = 1
    CEED_EVAL_GRAD = 2
    CEED_EVAL_DIV = 4
    CEED_EVAL_CURL = 8
    CEED_EVAL_WEIGHT = 16
end

@cenum CeedQuadMode::UInt32 begin
    CEED_GAUSS = 0
    CEED_GAUSS_LOBATTO = 1
end

@cenum CeedElemTopology::UInt32 begin
    CEED_LINE = 65536
    CEED_TRIANGLE = 131073
    CEED_QUAD = 131074
    CEED_TET = 196611
    CEED_PYRAMID = 196612
    CEED_PRISM = 196613
    CEED_HEX = 196614
end

function CeedBasisCreateTensorH1Lagrange(ceed, dim, num_comp, P, Q, quad_mode, basis)
    ccall((:CeedBasisCreateTensorH1Lagrange, libceed), Cint, (Ceed, CeedInt, CeedInt, CeedInt, CeedInt, CeedQuadMode, Ptr{CeedBasis}), ceed, dim, num_comp, P, Q, quad_mode, basis)
end

function CeedBasisCreateTensorH1(ceed, dim, num_comp, P_1d, Q_1d, interp_1d, grad_1d, q_ref_1d, q_weight_1d, basis)
    ccall((:CeedBasisCreateTensorH1, libceed), Cint, (Ceed, CeedInt, CeedInt, CeedInt, CeedInt, Ptr{CeedScalar}, Ptr{CeedScalar}, Ptr{CeedScalar}, Ptr{CeedScalar}, Ptr{CeedBasis}), ceed, dim, num_comp, P_1d, Q_1d, interp_1d, grad_1d, q_ref_1d, q_weight_1d, basis)
end

function CeedBasisCreateH1(ceed, topo, num_comp, num_nodes, nqpts, interp, grad, q_ref, q_weights, basis)
    ccall((:CeedBasisCreateH1, libceed), Cint, (Ceed, CeedElemTopology, CeedInt, CeedInt, CeedInt, Ptr{CeedScalar}, Ptr{CeedScalar}, Ptr{CeedScalar}, Ptr{CeedScalar}, Ptr{CeedBasis}), ceed, topo, num_comp, num_nodes, nqpts, interp, grad, q_ref, q_weights, basis)
end

function CeedBasisReferenceCopy(basis, basis_copy)
    ccall((:CeedBasisReferenceCopy, libceed), Cint, (CeedBasis, Ptr{CeedBasis}), basis, basis_copy)
end

function CeedBasisView(basis, stream)
    ccall((:CeedBasisView, libceed), Cint, (CeedBasis, Ptr{Libc.FILE}), basis, stream)
end

function CeedBasisApply(basis, num_elem, t_mode, eval_mode, u, v)
    ccall((:CeedBasisApply, libceed), Cint, (CeedBasis, CeedInt, CeedTransposeMode, CeedEvalMode, CeedVector, CeedVector), basis, num_elem, t_mode, eval_mode, u, v)
end

function CeedBasisGetDimension(basis, dim)
    ccall((:CeedBasisGetDimension, libceed), Cint, (CeedBasis, Ptr{CeedInt}), basis, dim)
end

function CeedBasisGetTopology(basis, topo)
    ccall((:CeedBasisGetTopology, libceed), Cint, (CeedBasis, Ptr{CeedElemTopology}), basis, topo)
end

function CeedBasisGetNumComponents(basis, num_comp)
    ccall((:CeedBasisGetNumComponents, libceed), Cint, (CeedBasis, Ptr{CeedInt}), basis, num_comp)
end

function CeedBasisGetNumNodes(basis, P)
    ccall((:CeedBasisGetNumNodes, libceed), Cint, (CeedBasis, Ptr{CeedInt}), basis, P)
end

function CeedBasisGetNumNodes1D(basis, P_1d)
    ccall((:CeedBasisGetNumNodes1D, libceed), Cint, (CeedBasis, Ptr{CeedInt}), basis, P_1d)
end

function CeedBasisGetNumQuadraturePoints(basis, Q)
    ccall((:CeedBasisGetNumQuadraturePoints, libceed), Cint, (CeedBasis, Ptr{CeedInt}), basis, Q)
end

function CeedBasisGetNumQuadraturePoints1D(basis, Q_1d)
    ccall((:CeedBasisGetNumQuadraturePoints1D, libceed), Cint, (CeedBasis, Ptr{CeedInt}), basis, Q_1d)
end

function CeedBasisGetQRef(basis, q_ref)
    ccall((:CeedBasisGetQRef, libceed), Cint, (CeedBasis, Ptr{Ptr{CeedScalar}}), basis, q_ref)
end

function CeedBasisGetQWeights(basis, q_weights)
    ccall((:CeedBasisGetQWeights, libceed), Cint, (CeedBasis, Ptr{Ptr{CeedScalar}}), basis, q_weights)
end

function CeedBasisGetInterp(basis, interp)
    ccall((:CeedBasisGetInterp, libceed), Cint, (CeedBasis, Ptr{Ptr{CeedScalar}}), basis, interp)
end

function CeedBasisGetInterp1D(basis, interp_1d)
    ccall((:CeedBasisGetInterp1D, libceed), Cint, (CeedBasis, Ptr{Ptr{CeedScalar}}), basis, interp_1d)
end

function CeedBasisGetGrad(basis, grad)
    ccall((:CeedBasisGetGrad, libceed), Cint, (CeedBasis, Ptr{Ptr{CeedScalar}}), basis, grad)
end

function CeedBasisGetGrad1D(basis, grad_1d)
    ccall((:CeedBasisGetGrad1D, libceed), Cint, (CeedBasis, Ptr{Ptr{CeedScalar}}), basis, grad_1d)
end

function CeedBasisDestroy(basis)
    ccall((:CeedBasisDestroy, libceed), Cint, (Ptr{CeedBasis},), basis)
end

function CeedGaussQuadrature(Q, q_ref_1d, q_weight_1d)
    ccall((:CeedGaussQuadrature, libceed), Cint, (CeedInt, Ptr{CeedScalar}, Ptr{CeedScalar}), Q, q_ref_1d, q_weight_1d)
end

function CeedLobattoQuadrature(Q, q_ref_1d, q_weight_1d)
    ccall((:CeedLobattoQuadrature, libceed), Cint, (CeedInt, Ptr{CeedScalar}, Ptr{CeedScalar}), Q, q_ref_1d, q_weight_1d)
end

function CeedQRFactorization(ceed, mat, tau, m, n)
    ccall((:CeedQRFactorization, libceed), Cint, (Ceed, Ptr{CeedScalar}, Ptr{CeedScalar}, CeedInt, CeedInt), ceed, mat, tau, m, n)
end

function CeedSymmetricSchurDecomposition(ceed, mat, lambda, n)
    ccall((:CeedSymmetricSchurDecomposition, libceed), Cint, (Ceed, Ptr{CeedScalar}, Ptr{CeedScalar}, CeedInt), ceed, mat, lambda, n)
end

function CeedSimultaneousDiagonalization(ceed, mat_A, mat_B, x, lambda, n)
    ccall((:CeedSimultaneousDiagonalization, libceed), Cint, (Ceed, Ptr{CeedScalar}, Ptr{CeedScalar}, Ptr{CeedScalar}, Ptr{CeedScalar}, CeedInt), ceed, mat_A, mat_B, x, lambda, n)
end

# typedef int ( * CeedQFunctionUser ) ( void * ctx , const CeedInt Q , const CeedScalar * const * in , CeedScalar * const * out )
const CeedQFunctionUser = Ptr{Cvoid}

function CeedQFunctionCreateInterior(ceed, vec_length, f, source, qf)
    ccall((:CeedQFunctionCreateInterior, libceed), Cint, (Ceed, CeedInt, CeedQFunctionUser, Ptr{Cchar}, Ptr{CeedQFunction}), ceed, vec_length, f, source, qf)
end

function CeedQFunctionCreateInteriorByName(ceed, name, qf)
    ccall((:CeedQFunctionCreateInteriorByName, libceed), Cint, (Ceed, Ptr{Cchar}, Ptr{CeedQFunction}), ceed, name, qf)
end

function CeedQFunctionCreateIdentity(ceed, size, in_mode, out_mode, qf)
    ccall((:CeedQFunctionCreateIdentity, libceed), Cint, (Ceed, CeedInt, CeedEvalMode, CeedEvalMode, Ptr{CeedQFunction}), ceed, size, in_mode, out_mode, qf)
end

function CeedQFunctionReferenceCopy(qf, qf_copy)
    ccall((:CeedQFunctionReferenceCopy, libceed), Cint, (CeedQFunction, Ptr{CeedQFunction}), qf, qf_copy)
end

function CeedQFunctionAddInput(qf, field_name, size, eval_mode)
    ccall((:CeedQFunctionAddInput, libceed), Cint, (CeedQFunction, Ptr{Cchar}, CeedInt, CeedEvalMode), qf, field_name, size, eval_mode)
end

function CeedQFunctionAddOutput(qf, field_name, size, eval_mode)
    ccall((:CeedQFunctionAddOutput, libceed), Cint, (CeedQFunction, Ptr{Cchar}, CeedInt, CeedEvalMode), qf, field_name, size, eval_mode)
end

function CeedQFunctionSetContext(qf, ctx)
    ccall((:CeedQFunctionSetContext, libceed), Cint, (CeedQFunction, CeedQFunctionContext), qf, ctx)
end

function CeedQFunctionView(qf, stream)
    ccall((:CeedQFunctionView, libceed), Cint, (CeedQFunction, Ptr{Libc.FILE}), qf, stream)
end

function CeedQFunctionApply(qf, Q, u, v)
    ccall((:CeedQFunctionApply, libceed), Cint, (CeedQFunction, CeedInt, Ptr{CeedVector}, Ptr{CeedVector}), qf, Q, u, v)
end

function CeedQFunctionDestroy(qf)
    ccall((:CeedQFunctionDestroy, libceed), Cint, (Ptr{CeedQFunction},), qf)
end

function CeedQFunctionContextCreate(ceed, ctx)
    ccall((:CeedQFunctionContextCreate, libceed), Cint, (Ceed, Ptr{CeedQFunctionContext}), ceed, ctx)
end

function CeedQFunctionContextReferenceCopy(ctx, ctx_copy)
    ccall((:CeedQFunctionContextReferenceCopy, libceed), Cint, (CeedQFunctionContext, Ptr{CeedQFunctionContext}), ctx, ctx_copy)
end

function CeedQFunctionContextSetData(ctx, mem_type, copy_mode, size, data)
    ccall((:CeedQFunctionContextSetData, libceed), Cint, (CeedQFunctionContext, CeedMemType, CeedCopyMode, Csize_t, Ptr{Cvoid}), ctx, mem_type, copy_mode, size, data)
end

function CeedQFunctionContextTakeData(ctx, mem_type, data)
    ccall((:CeedQFunctionContextTakeData, libceed), Cint, (CeedQFunctionContext, CeedMemType, Ptr{Cvoid}), ctx, mem_type, data)
end

function CeedQFunctionContextGetData(ctx, mem_type, data)
    ccall((:CeedQFunctionContextGetData, libceed), Cint, (CeedQFunctionContext, CeedMemType, Ptr{Cvoid}), ctx, mem_type, data)
end

function CeedQFunctionContextRestoreData(ctx, data)
    ccall((:CeedQFunctionContextRestoreData, libceed), Cint, (CeedQFunctionContext, Ptr{Cvoid}), ctx, data)
end

function CeedQFunctionContextView(ctx, stream)
    ccall((:CeedQFunctionContextView, libceed), Cint, (CeedQFunctionContext, Ptr{Libc.FILE}), ctx, stream)
end

function CeedQFunctionContextDestroy(ctx)
    ccall((:CeedQFunctionContextDestroy, libceed), Cint, (Ptr{CeedQFunctionContext},), ctx)
end

function CeedOperatorCreate(ceed, qf, dqf, dqfT, op)
    ccall((:CeedOperatorCreate, libceed), Cint, (Ceed, CeedQFunction, CeedQFunction, CeedQFunction, Ptr{CeedOperator}), ceed, qf, dqf, dqfT, op)
end

function CeedCompositeOperatorCreate(ceed, op)
    ccall((:CeedCompositeOperatorCreate, libceed), Cint, (Ceed, Ptr{CeedOperator}), ceed, op)
end

function CeedOperatorReferenceCopy(op, op_copy)
    ccall((:CeedOperatorReferenceCopy, libceed), Cint, (CeedOperator, Ptr{CeedOperator}), op, op_copy)
end

function CeedOperatorSetField(op, field_name, r, b, v)
    ccall((:CeedOperatorSetField, libceed), Cint, (CeedOperator, Ptr{Cchar}, CeedElemRestriction, CeedBasis, CeedVector), op, field_name, r, b, v)
end

function CeedCompositeOperatorAddSub(composite_op, sub_op)
    ccall((:CeedCompositeOperatorAddSub, libceed), Cint, (CeedOperator, CeedOperator), composite_op, sub_op)
end

function CeedOperatorLinearAssembleQFunction(op, assembled, rstr, request)
    ccall((:CeedOperatorLinearAssembleQFunction, libceed), Cint, (CeedOperator, Ptr{CeedVector}, Ptr{CeedElemRestriction}, Ptr{CeedRequest}), op, assembled, rstr, request)
end

function CeedOperatorLinearAssembleDiagonal(op, assembled, request)
    ccall((:CeedOperatorLinearAssembleDiagonal, libceed), Cint, (CeedOperator, CeedVector, Ptr{CeedRequest}), op, assembled, request)
end

function CeedOperatorLinearAssembleAddDiagonal(op, assembled, request)
    ccall((:CeedOperatorLinearAssembleAddDiagonal, libceed), Cint, (CeedOperator, CeedVector, Ptr{CeedRequest}), op, assembled, request)
end

function CeedOperatorLinearAssemblePointBlockDiagonal(op, assembled, request)
    ccall((:CeedOperatorLinearAssemblePointBlockDiagonal, libceed), Cint, (CeedOperator, CeedVector, Ptr{CeedRequest}), op, assembled, request)
end

function CeedOperatorLinearAssembleAddPointBlockDiagonal(op, assembled, request)
    ccall((:CeedOperatorLinearAssembleAddPointBlockDiagonal, libceed), Cint, (CeedOperator, CeedVector, Ptr{CeedRequest}), op, assembled, request)
end

function CeedOperatorLinearAssembleSymbolic(op, num_entries, rows, cols)
    ccall((:CeedOperatorLinearAssembleSymbolic, libceed), Cint, (CeedOperator, Ptr{CeedInt}, Ptr{Ptr{CeedInt}}, Ptr{Ptr{CeedInt}}), op, num_entries, rows, cols)
end

function CeedOperatorLinearAssemble(op, values)
    ccall((:CeedOperatorLinearAssemble, libceed), Cint, (CeedOperator, CeedVector), op, values)
end

function CeedOperatorMultigridLevelCreate(op_fine, p_mult_fine, rstr_coarse, basis_coarse, op_coarse, op_prolong, op_restrict)
    ccall((:CeedOperatorMultigridLevelCreate, libceed), Cint, (CeedOperator, CeedVector, CeedElemRestriction, CeedBasis, Ptr{CeedOperator}, Ptr{CeedOperator}, Ptr{CeedOperator}), op_fine, p_mult_fine, rstr_coarse, basis_coarse, op_coarse, op_prolong, op_restrict)
end

function CeedOperatorMultigridLevelCreateTensorH1(op_fine, p_mult_fine, rstr_coarse, basis_coarse, interp_c_to_f, op_coarse, op_prolong, op_restrict)
    ccall((:CeedOperatorMultigridLevelCreateTensorH1, libceed), Cint, (CeedOperator, CeedVector, CeedElemRestriction, CeedBasis, Ptr{CeedScalar}, Ptr{CeedOperator}, Ptr{CeedOperator}, Ptr{CeedOperator}), op_fine, p_mult_fine, rstr_coarse, basis_coarse, interp_c_to_f, op_coarse, op_prolong, op_restrict)
end

function CeedOperatorMultigridLevelCreateH1(op_fine, p_mult_fine, rstr_coarse, basis_coarse, interp_c_to_f, op_coarse, op_prolong, op_restrict)
    ccall((:CeedOperatorMultigridLevelCreateH1, libceed), Cint, (CeedOperator, CeedVector, CeedElemRestriction, CeedBasis, Ptr{CeedScalar}, Ptr{CeedOperator}, Ptr{CeedOperator}, Ptr{CeedOperator}), op_fine, p_mult_fine, rstr_coarse, basis_coarse, interp_c_to_f, op_coarse, op_prolong, op_restrict)
end

function CeedOperatorCreateFDMElementInverse(op, fdm_inv, request)
    ccall((:CeedOperatorCreateFDMElementInverse, libceed), Cint, (CeedOperator, Ptr{CeedOperator}, Ptr{CeedRequest}), op, fdm_inv, request)
end

function CeedOperatorSetNumQuadraturePoints(op, num_qpts)
    ccall((:CeedOperatorSetNumQuadraturePoints, libceed), Cint, (CeedOperator, CeedInt), op, num_qpts)
end

function CeedOperatorView(op, stream)
    ccall((:CeedOperatorView, libceed), Cint, (CeedOperator, Ptr{Libc.FILE}), op, stream)
end

function CeedOperatorApply(op, in, out, request)
    ccall((:CeedOperatorApply, libceed), Cint, (CeedOperator, CeedVector, CeedVector, Ptr{CeedRequest}), op, in, out, request)
end

function CeedOperatorApplyAdd(op, in, out, request)
    ccall((:CeedOperatorApplyAdd, libceed), Cint, (CeedOperator, CeedVector, CeedVector, Ptr{CeedRequest}), op, in, out, request)
end

function CeedOperatorDestroy(op)
    ccall((:CeedOperatorDestroy, libceed), Cint, (Ptr{CeedOperator},), op)
end

function CeedIntPow(base, power)
    ccall((:CeedIntPow, libceed), CeedInt, (CeedInt, CeedInt), base, power)
end

function CeedIntMin(a, b)
    ccall((:CeedIntMin, libceed), CeedInt, (CeedInt, CeedInt), a, b)
end

function CeedIntMax(a, b)
    ccall((:CeedIntMax, libceed), CeedInt, (CeedInt, CeedInt), a, b)
end

function CeedRegisterAll()
    ccall((:CeedRegisterAll, libceed), Cint, ())
end

function CeedQFunctionRegisterAll()
    ccall((:CeedQFunctionRegisterAll, libceed), Cint, ())
end

function CeedQFunctionSetCUDAUserFunction(qf, f)
    ccall((:CeedQFunctionSetCUDAUserFunction, libceed), Cint, (CeedQFunction, Cint), qf, f)
end

# automatic type deduction for variadic arguments may not be what you want, please use with caution
@generated function CeedDebugImpl(arg1, arg2, va_list...)
        :(@ccall(libceed.CeedDebugImpl(arg1::Ceed, arg2::Ptr{Cchar}; $(to_c_type_pairs(va_list)...))::Cvoid))
    end

# automatic type deduction for variadic arguments may not be what you want, please use with caution
@generated function CeedDebugImpl256(arg1, arg2, arg3, va_list...)
        :(@ccall(libceed.CeedDebugImpl256(arg1::Ceed, arg2::Cuchar, arg3::Ptr{Cchar}; $(to_c_type_pairs(va_list)...))::Cvoid))
    end

mutable struct CeedTensorContract_private end

const CeedTensorContract = Ptr{CeedTensorContract_private}

function CeedMallocArray(n, unit, p)
    ccall((:CeedMallocArray, libceed), Cint, (Csize_t, Csize_t, Ptr{Cvoid}), n, unit, p)
end

function CeedCallocArray(n, unit, p)
    ccall((:CeedCallocArray, libceed), Cint, (Csize_t, Csize_t, Ptr{Cvoid}), n, unit, p)
end

function CeedReallocArray(n, unit, p)
    ccall((:CeedReallocArray, libceed), Cint, (Csize_t, Csize_t, Ptr{Cvoid}), n, unit, p)
end

function CeedFree(p)
    ccall((:CeedFree, libceed), Cint, (Ptr{Cvoid},), p)
end

mutable struct CeedQFunctionField_private end

const CeedQFunctionField = Ptr{CeedQFunctionField_private}

mutable struct CeedOperatorField_private end

const CeedOperatorField = Ptr{CeedOperatorField_private}

function CeedRegister(prefix, init, priority)
    ccall((:CeedRegister, libceed), Cint, (Ptr{Cchar}, Ptr{Cvoid}, Cuint), prefix, init, priority)
end

function CeedIsDebug(ceed, is_debug)
    ccall((:CeedIsDebug, libceed), Cint, (Ceed, Ptr{Bool}), ceed, is_debug)
end

function CeedGetParent(ceed, parent)
    ccall((:CeedGetParent, libceed), Cint, (Ceed, Ptr{Ceed}), ceed, parent)
end

function CeedGetDelegate(ceed, delegate)
    ccall((:CeedGetDelegate, libceed), Cint, (Ceed, Ptr{Ceed}), ceed, delegate)
end

function CeedSetDelegate(ceed, delegate)
    ccall((:CeedSetDelegate, libceed), Cint, (Ceed, Ceed), ceed, delegate)
end

function CeedGetObjectDelegate(ceed, delegate, obj_name)
    ccall((:CeedGetObjectDelegate, libceed), Cint, (Ceed, Ptr{Ceed}, Ptr{Cchar}), ceed, delegate, obj_name)
end

function CeedSetObjectDelegate(ceed, delegate, obj_name)
    ccall((:CeedSetObjectDelegate, libceed), Cint, (Ceed, Ceed, Ptr{Cchar}), ceed, delegate, obj_name)
end

function CeedGetOperatorFallbackResource(ceed, resource)
    ccall((:CeedGetOperatorFallbackResource, libceed), Cint, (Ceed, Ptr{Ptr{Cchar}}), ceed, resource)
end

function CeedSetOperatorFallbackResource(ceed, resource)
    ccall((:CeedSetOperatorFallbackResource, libceed), Cint, (Ceed, Ptr{Cchar}), ceed, resource)
end

function CeedGetOperatorFallbackParentCeed(ceed, parent)
    ccall((:CeedGetOperatorFallbackParentCeed, libceed), Cint, (Ceed, Ptr{Ceed}), ceed, parent)
end

function CeedSetDeterministic(ceed, is_deterministic)
    ccall((:CeedSetDeterministic, libceed), Cint, (Ceed, Bool), ceed, is_deterministic)
end

function CeedSetBackendFunction(ceed, type, object, func_name, f)
    ccall((:CeedSetBackendFunction, libceed), Cint, (Ceed, Ptr{Cchar}, Ptr{Cvoid}, Ptr{Cchar}, Ptr{Cvoid}), ceed, type, object, func_name, f)
end

function CeedGetData(ceed, data)
    ccall((:CeedGetData, libceed), Cint, (Ceed, Ptr{Cvoid}), ceed, data)
end

function CeedSetData(ceed, data)
    ccall((:CeedSetData, libceed), Cint, (Ceed, Ptr{Cvoid}), ceed, data)
end

function CeedReference(ceed)
    ccall((:CeedReference, libceed), Cint, (Ceed,), ceed)
end

function CeedVectorGetCeed(vec, ceed)
    ccall((:CeedVectorGetCeed, libceed), Cint, (CeedVector, Ptr{Ceed}), vec, ceed)
end

function CeedVectorGetState(vec, state)
    ccall((:CeedVectorGetState, libceed), Cint, (CeedVector, Ptr{UInt64}), vec, state)
end

function CeedVectorAddReference(vec)
    ccall((:CeedVectorAddReference, libceed), Cint, (CeedVector,), vec)
end

function CeedVectorGetData(vec, data)
    ccall((:CeedVectorGetData, libceed), Cint, (CeedVector, Ptr{Cvoid}), vec, data)
end

function CeedVectorSetData(vec, data)
    ccall((:CeedVectorSetData, libceed), Cint, (CeedVector, Ptr{Cvoid}), vec, data)
end

function CeedVectorReference(vec)
    ccall((:CeedVectorReference, libceed), Cint, (CeedVector,), vec)
end

function CeedElemRestrictionGetCeed(rstr, ceed)
    ccall((:CeedElemRestrictionGetCeed, libceed), Cint, (CeedElemRestriction, Ptr{Ceed}), rstr, ceed)
end

function CeedElemRestrictionGetStrides(rstr, strides)
    ccall((:CeedElemRestrictionGetStrides, libceed), Cint, (CeedElemRestriction, Ptr{NTuple{3, CeedInt}}), rstr, strides)
end

function CeedElemRestrictionGetOffsets(rstr, mem_type, offsets)
    ccall((:CeedElemRestrictionGetOffsets, libceed), Cint, (CeedElemRestriction, CeedMemType, Ptr{Ptr{CeedInt}}), rstr, mem_type, offsets)
end

function CeedElemRestrictionRestoreOffsets(rstr, offsets)
    ccall((:CeedElemRestrictionRestoreOffsets, libceed), Cint, (CeedElemRestriction, Ptr{Ptr{CeedInt}}), rstr, offsets)
end

function CeedElemRestrictionIsStrided(rstr, is_strided)
    ccall((:CeedElemRestrictionIsStrided, libceed), Cint, (CeedElemRestriction, Ptr{Bool}), rstr, is_strided)
end

function CeedElemRestrictionHasBackendStrides(rstr, has_backend_strides)
    ccall((:CeedElemRestrictionHasBackendStrides, libceed), Cint, (CeedElemRestriction, Ptr{Bool}), rstr, has_backend_strides)
end

function CeedElemRestrictionGetELayout(rstr, layout)
    ccall((:CeedElemRestrictionGetELayout, libceed), Cint, (CeedElemRestriction, Ptr{NTuple{3, CeedInt}}), rstr, layout)
end

function CeedElemRestrictionSetELayout(rstr, layout)
    ccall((:CeedElemRestrictionSetELayout, libceed), Cint, (CeedElemRestriction, Ptr{CeedInt}), rstr, layout)
end

function CeedElemRestrictionGetData(rstr, data)
    ccall((:CeedElemRestrictionGetData, libceed), Cint, (CeedElemRestriction, Ptr{Cvoid}), rstr, data)
end

function CeedElemRestrictionSetData(rstr, data)
    ccall((:CeedElemRestrictionSetData, libceed), Cint, (CeedElemRestriction, Ptr{Cvoid}), rstr, data)
end

function CeedElemRestrictionReference(rstr)
    ccall((:CeedElemRestrictionReference, libceed), Cint, (CeedElemRestriction,), rstr)
end

function CeedBasisGetCollocatedGrad(basis, colo_grad_1d)
    ccall((:CeedBasisGetCollocatedGrad, libceed), Cint, (CeedBasis, Ptr{CeedScalar}), basis, colo_grad_1d)
end

function CeedHouseholderApplyQ(A, Q, tau, t_mode, m, n, k, row, col)
    ccall((:CeedHouseholderApplyQ, libceed), Cint, (Ptr{CeedScalar}, Ptr{CeedScalar}, Ptr{CeedScalar}, CeedTransposeMode, CeedInt, CeedInt, CeedInt, CeedInt, CeedInt), A, Q, tau, t_mode, m, n, k, row, col)
end

function CeedBasisGetCeed(basis, ceed)
    ccall((:CeedBasisGetCeed, libceed), Cint, (CeedBasis, Ptr{Ceed}), basis, ceed)
end

function CeedBasisIsTensor(basis, is_tensor)
    ccall((:CeedBasisIsTensor, libceed), Cint, (CeedBasis, Ptr{Bool}), basis, is_tensor)
end

function CeedBasisGetData(basis, data)
    ccall((:CeedBasisGetData, libceed), Cint, (CeedBasis, Ptr{Cvoid}), basis, data)
end

function CeedBasisSetData(basis, data)
    ccall((:CeedBasisSetData, libceed), Cint, (CeedBasis, Ptr{Cvoid}), basis, data)
end

function CeedBasisReference(basis)
    ccall((:CeedBasisReference, libceed), Cint, (CeedBasis,), basis)
end

function CeedBasisGetTopologyDimension(topo, dim)
    ccall((:CeedBasisGetTopologyDimension, libceed), Cint, (CeedElemTopology, Ptr{CeedInt}), topo, dim)
end

function CeedBasisGetTensorContract(basis, contract)
    ccall((:CeedBasisGetTensorContract, libceed), Cint, (CeedBasis, Ptr{CeedTensorContract}), basis, contract)
end

function CeedBasisSetTensorContract(basis, contract)
    ccall((:CeedBasisSetTensorContract, libceed), Cint, (CeedBasis, CeedTensorContract), basis, contract)
end

function CeedTensorContractCreate(ceed, basis, contract)
    ccall((:CeedTensorContractCreate, libceed), Cint, (Ceed, CeedBasis, Ptr{CeedTensorContract}), ceed, basis, contract)
end

function CeedTensorContractApply(contract, A, B, C, J, t, t_mode, Add, u, v)
    ccall((:CeedTensorContractApply, libceed), Cint, (CeedTensorContract, CeedInt, CeedInt, CeedInt, CeedInt, Ptr{CeedScalar}, CeedTransposeMode, CeedInt, Ptr{CeedScalar}, Ptr{CeedScalar}), contract, A, B, C, J, t, t_mode, Add, u, v)
end

function CeedTensorContractGetCeed(contract, ceed)
    ccall((:CeedTensorContractGetCeed, libceed), Cint, (CeedTensorContract, Ptr{Ceed}), contract, ceed)
end

function CeedTensorContractGetData(contract, data)
    ccall((:CeedTensorContractGetData, libceed), Cint, (CeedTensorContract, Ptr{Cvoid}), contract, data)
end

function CeedTensorContractSetData(contract, data)
    ccall((:CeedTensorContractSetData, libceed), Cint, (CeedTensorContract, Ptr{Cvoid}), contract, data)
end

function CeedTensorContractReference(contract)
    ccall((:CeedTensorContractReference, libceed), Cint, (CeedTensorContract,), contract)
end

function CeedTensorContractDestroy(contract)
    ccall((:CeedTensorContractDestroy, libceed), Cint, (Ptr{CeedTensorContract},), contract)
end

function CeedQFunctionRegister(arg1, arg2, arg3, arg4, init)
    ccall((:CeedQFunctionRegister, libceed), Cint, (Ptr{Cchar}, Ptr{Cchar}, CeedInt, CeedQFunctionUser, Ptr{Cvoid}), arg1, arg2, arg3, arg4, init)
end

function CeedQFunctionSetFortranStatus(qf, status)
    ccall((:CeedQFunctionSetFortranStatus, libceed), Cint, (CeedQFunction, Bool), qf, status)
end

function CeedQFunctionGetCeed(qf, ceed)
    ccall((:CeedQFunctionGetCeed, libceed), Cint, (CeedQFunction, Ptr{Ceed}), qf, ceed)
end

function CeedQFunctionGetVectorLength(qf, vec_length)
    ccall((:CeedQFunctionGetVectorLength, libceed), Cint, (CeedQFunction, Ptr{CeedInt}), qf, vec_length)
end

function CeedQFunctionGetNumArgs(qf, num_input_fields, num_output_fields)
    ccall((:CeedQFunctionGetNumArgs, libceed), Cint, (CeedQFunction, Ptr{CeedInt}, Ptr{CeedInt}), qf, num_input_fields, num_output_fields)
end

function CeedQFunctionGetSourcePath(qf, source)
    ccall((:CeedQFunctionGetSourcePath, libceed), Cint, (CeedQFunction, Ptr{Ptr{Cchar}}), qf, source)
end

function CeedQFunctionGetUserFunction(qf, f)
    ccall((:CeedQFunctionGetUserFunction, libceed), Cint, (CeedQFunction, Ptr{CeedQFunctionUser}), qf, f)
end

function CeedQFunctionGetContext(qf, ctx)
    ccall((:CeedQFunctionGetContext, libceed), Cint, (CeedQFunction, Ptr{CeedQFunctionContext}), qf, ctx)
end

function CeedQFunctionGetInnerContext(qf, ctx)
    ccall((:CeedQFunctionGetInnerContext, libceed), Cint, (CeedQFunction, Ptr{CeedQFunctionContext}), qf, ctx)
end

function CeedQFunctionIsIdentity(qf, is_identity)
    ccall((:CeedQFunctionIsIdentity, libceed), Cint, (CeedQFunction, Ptr{Bool}), qf, is_identity)
end

function CeedQFunctionGetData(qf, data)
    ccall((:CeedQFunctionGetData, libceed), Cint, (CeedQFunction, Ptr{Cvoid}), qf, data)
end

function CeedQFunctionSetData(qf, data)
    ccall((:CeedQFunctionSetData, libceed), Cint, (CeedQFunction, Ptr{Cvoid}), qf, data)
end

function CeedQFunctionReference(qf)
    ccall((:CeedQFunctionReference, libceed), Cint, (CeedQFunction,), qf)
end

function CeedQFunctionGetFields(qf, input_fields, output_fields)
    ccall((:CeedQFunctionGetFields, libceed), Cint, (CeedQFunction, Ptr{Ptr{CeedQFunctionField}}, Ptr{Ptr{CeedQFunctionField}}), qf, input_fields, output_fields)
end

function CeedQFunctionFieldGetName(qf_field, field_name)
    ccall((:CeedQFunctionFieldGetName, libceed), Cint, (CeedQFunctionField, Ptr{Ptr{Cchar}}), qf_field, field_name)
end

function CeedQFunctionFieldGetSize(qf_field, size)
    ccall((:CeedQFunctionFieldGetSize, libceed), Cint, (CeedQFunctionField, Ptr{CeedInt}), qf_field, size)
end

function CeedQFunctionFieldGetEvalMode(qf_field, eval_mode)
    ccall((:CeedQFunctionFieldGetEvalMode, libceed), Cint, (CeedQFunctionField, Ptr{CeedEvalMode}), qf_field, eval_mode)
end

function CeedQFunctionContextGetCeed(cxt, ceed)
    ccall((:CeedQFunctionContextGetCeed, libceed), Cint, (CeedQFunctionContext, Ptr{Ceed}), cxt, ceed)
end

function CeedQFunctionContextGetState(ctx, state)
    ccall((:CeedQFunctionContextGetState, libceed), Cint, (CeedQFunctionContext, Ptr{UInt64}), ctx, state)
end

function CeedQFunctionContextGetContextSize(ctx, ctx_size)
    ccall((:CeedQFunctionContextGetContextSize, libceed), Cint, (CeedQFunctionContext, Ptr{Csize_t}), ctx, ctx_size)
end

function CeedQFunctionContextGetBackendData(ctx, data)
    ccall((:CeedQFunctionContextGetBackendData, libceed), Cint, (CeedQFunctionContext, Ptr{Cvoid}), ctx, data)
end

function CeedQFunctionContextSetBackendData(ctx, data)
    ccall((:CeedQFunctionContextSetBackendData, libceed), Cint, (CeedQFunctionContext, Ptr{Cvoid}), ctx, data)
end

function CeedQFunctionContextReference(ctx)
    ccall((:CeedQFunctionContextReference, libceed), Cint, (CeedQFunctionContext,), ctx)
end

function CeedOperatorGetCeed(op, ceed)
    ccall((:CeedOperatorGetCeed, libceed), Cint, (CeedOperator, Ptr{Ceed}), op, ceed)
end

function CeedOperatorGetNumElements(op, num_elem)
    ccall((:CeedOperatorGetNumElements, libceed), Cint, (CeedOperator, Ptr{CeedInt}), op, num_elem)
end

function CeedOperatorGetNumQuadraturePoints(op, num_qpts)
    ccall((:CeedOperatorGetNumQuadraturePoints, libceed), Cint, (CeedOperator, Ptr{CeedInt}), op, num_qpts)
end

function CeedOperatorGetNumArgs(op, num_args)
    ccall((:CeedOperatorGetNumArgs, libceed), Cint, (CeedOperator, Ptr{CeedInt}), op, num_args)
end

function CeedOperatorIsSetupDone(op, is_setup_done)
    ccall((:CeedOperatorIsSetupDone, libceed), Cint, (CeedOperator, Ptr{Bool}), op, is_setup_done)
end

function CeedOperatorGetQFunction(op, qf)
    ccall((:CeedOperatorGetQFunction, libceed), Cint, (CeedOperator, Ptr{CeedQFunction}), op, qf)
end

function CeedOperatorIsComposite(op, is_composite)
    ccall((:CeedOperatorIsComposite, libceed), Cint, (CeedOperator, Ptr{Bool}), op, is_composite)
end

function CeedOperatorGetNumSub(op, num_suboperators)
    ccall((:CeedOperatorGetNumSub, libceed), Cint, (CeedOperator, Ptr{CeedInt}), op, num_suboperators)
end

function CeedOperatorGetSubList(op, sub_operators)
    ccall((:CeedOperatorGetSubList, libceed), Cint, (CeedOperator, Ptr{Ptr{CeedOperator}}), op, sub_operators)
end

function CeedOperatorGetData(op, data)
    ccall((:CeedOperatorGetData, libceed), Cint, (CeedOperator, Ptr{Cvoid}), op, data)
end

function CeedOperatorSetData(op, data)
    ccall((:CeedOperatorSetData, libceed), Cint, (CeedOperator, Ptr{Cvoid}), op, data)
end

function CeedOperatorReference(op)
    ccall((:CeedOperatorReference, libceed), Cint, (CeedOperator,), op)
end

function CeedOperatorSetSetupDone(op)
    ccall((:CeedOperatorSetSetupDone, libceed), Cint, (CeedOperator,), op)
end

function CeedOperatorGetFields(op, input_fields, output_fields)
    ccall((:CeedOperatorGetFields, libceed), Cint, (CeedOperator, Ptr{Ptr{CeedOperatorField}}, Ptr{Ptr{CeedOperatorField}}), op, input_fields, output_fields)
end

function CeedOperatorFieldGetElemRestriction(op_field, rstr)
    ccall((:CeedOperatorFieldGetElemRestriction, libceed), Cint, (CeedOperatorField, Ptr{CeedElemRestriction}), op_field, rstr)
end

function CeedOperatorFieldGetBasis(op_field, basis)
    ccall((:CeedOperatorFieldGetBasis, libceed), Cint, (CeedOperatorField, Ptr{CeedBasis}), op_field, basis)
end

function CeedOperatorFieldGetVector(op_field, vec)
    ccall((:CeedOperatorFieldGetVector, libceed), Cint, (CeedOperatorField, Ptr{CeedVector}), op_field, vec)
end

function CeedMatrixMultiply(ceed, mat_A, mat_B, mat_C, m, n, kk)
    ccall((:CeedMatrixMultiply, libceed), Cint, (Ceed, Ptr{CeedScalar}, Ptr{CeedScalar}, Ptr{CeedScalar}, CeedInt, CeedInt, CeedInt), ceed, mat_A, mat_B, mat_C, m, n, kk)
end

# Skipping MacroDefinition: CEED_EXTERN extern

# Skipping MacroDefinition: CEED_QFUNCTION_HELPER static inline

const CEED_VERSION_MAJOR = 0

const CEED_VERSION_MINOR = 9

const CEED_VERSION_PATCH = 0

const CEED_VERSION_RELEASE = false

# Skipping MacroDefinition: CEED_INTERN CEED_EXTERN __attribute__ ( ( visibility ( "hidden" ) ) )

# Skipping MacroDefinition: CEED_UNUSED __attribute__ ( ( unused ) )

const CEED_MAX_RESOURCE_LEN = 1024

const CEED_MAX_BACKEND_PRIORITY = UINT_MAX

const CEED_ALIGN = 64

const CEED_COMPOSITE_MAX = 16

const CEED_EPSILON = 1.0e-16

# Skipping MacroDefinition: CeedPragmaOptimizeOff _Pragma ( "clang optimize off" )

# Skipping MacroDefinition: CeedPragmaOptimizeOn _Pragma ( "clang optimize on" )

const CEED_DEBUG_COLOR = 0

max_stars_repo_name: jedbrown/libCEED

max_issues_repo_name: jedbrown/libCEED

max_forks_repo_name: jedbrown/libCEED


text:
using MAT

useproc("GPU") # GPU is much faster than CPU
#useproc("CPU")

Ntrain=60000
Ntest=10000
BatchSize=100
TrainingIts=5000 # number of training iterations
include("loadmnist.jl")
images,label=loadmnist()
label=convert(Array{Int,2},label)
r=randperm(size(images,2))
data=images[:,r]; label=label[r]
class=zeros(10,70000)+0.00001
for i=1:70000
    class[label[i]+1,i]=0.9999
end
traindata=data[:,1:Ntrain]
testdata=data[:,Ntrain+1:Ntrain+Ntest]
trainclass=class[:,1:Ntrain]
testclass=class[:,Ntrain+1:Ntrain+Ntest]

# Construct the DAG function:
H=[784 3000 200 100 50 10] # number of units in each layer

L=length(H) # number of hidden layers
# node indices:
bias=Array(ADnode,L) # weight index
w=Array(ADnode,L) # weight index
h=Array(ADnode,L) # hidden layer index (note that I call the input layer h[1])
#hh=Array(ADnode,L) # hidden layer index (note that I call the input layer h[1])
dropout=Array(ADnode,L) # hidden layer index (note that I call the input layer h[1])
StartCode()
class=ADnode() # node has no parents
x=ADnode()
dropout[1]=ADnode()
h[1]=diagm(dropout[1])*x
for i=2:L-1
    dropout[i]=ADnode()
    w[i]=ADvariable()
    bias[i]=ADvariable()
    h[i]=diagm(dropout[i])*rectlinAXplusBias(w[i],h[i-1],bias[i])
end
w[L]=ADvariable()
h[L]= w[L]*h[L-1]
loss=KLsoftmax(class,h[L]) # the loss we are minimising must be the final node in the graph.
net=EndCode()

# instantiate parameter nodes and inputs:
net.value[x]=traindata[:,1:BatchSize]
net.value[class]=trainclass[:,1:BatchSize]
for i=2:L
    net.value[w[i]]=sign(randn(H[i],H[i-1]))/sqrt(H[i-1])
end
net.value[dropout[1]]=ones(H[1])
for i=2:L-1
    net.value[bias[i]]=zeros(H[i])
    net.value[dropout[i]]=ones(H[i],1)
end

net=compile(net) # compile the DAG and preallocate memory
#gradcheck(net)

@gpu CUDArt.init([0])
@gpu net=convert(net,"GPU")

tmp=net.value[1]

ADforward!(net)

dropoutprob=zeros(L)
dropoutprob[1]=0
#for i=2:L-1
#    dropoutprob[i]=0.2
#end
dropoutprob[2]=0.5
println("Training: using $(net.gpu==true? "GPU" : "CPU") ")
error=Array(Float64,0)
ParsToUpdate=Parameters(net)
meangrad=zeros(TrainingIts)
meanv=zeros(TrainingIts)
velo=NesterovInit(net)
minibatchstart=1 # starting datapoint for the minibatch
for iter=1:TrainingIts
    LearningRate=0.1/(1+iter/1500)
    minibatchstart,minibatch=GetBatch(minibatchstart,BatchSize,Ntrain)
    net.value[x]=cArray(PROC,traindata[:,minibatch]) # select batch of data
    net.value[class]=cArray(PROC,trainclass[:,minibatch]) # select batch of data
    for i=2:L-1
        net.value[dropout[i]]=cArray(PROC,1.0*(rand(H[i],1).>dropoutprob[i]))
    end
    ADforward!(net)
    ADbackward!(net)
    push!(error,extract(net.value[net.FunctionNode])[1])
    printover("iteration $iter: training loss = $(error[iter])")
    for par in ParsToUpdate
        #GradientDescentUpdate!(net.value[par],net.gradient[par],LearningRate)
        NesterovGradientDescentUpdate!(net.value[par],net.gradient[par],velo[par],LearningRate,iter)
    end
end
#For prediction scaling the weights is equivalent to using a fixed dropout value:
for i=1:L-1
    net.value[dropout[i]]=cArray((1-dropoutprob[i])*ones(H[i],1))
end
ADforward!(net)

net=convert(net,"CPU") # convert to CPU for easy results analysis

net.value[h[L]]
classpredtrain=argcolmax(softmax(net.value[h[L]]))
classtrain=argcolmax(net.value[class])
println("\nlast batch train accuracy = $(mean(classpredtrain.==classtrain))")

println("Testing: using $(net.gpu==true? "GPU" : "CPU") ")
# evaluate the test predictions:
net.value[x]=testdata
ForwardPassList!(net,ExcludeNodes=[class])
# ADforward! by default computes values inplace and does no additional memory allocation on the graph. If we change the input value dimensions (which happens here since there are a different number of test points to train points), we need to reallocatte memory based on the shapes of the new values on the graph. This only needs to be done once.
ADforward!(net,AllocateMemory=true) # reallocate since shape of net.value[x] changed

classpredtest=argcolmax(softmax(net.value[h[L]]))
classtest=argcolmax(testclass)
println("test accuracy = $(mean(classpredtest.==classtest))")

#end

max_stars_repo_name: davidbarber/AutoDiff

max_issues_repo_name: davidbarber/AutoDiff

max_forks_repo_name: davidbarber/AutoDiff


text:
import JuliaParser.Parser

using AbstractTrees
using Base.Test

if length(ARGS) != 1
    error("need to specify filepath: ./diffast.jl [filepath]")
elseif !ispath(ARGS[1])
    error("need to specify valid filepath, got: $(ARGS[1])")
end

include("ast.jl")

const src = IOBuffer()

# wrap source in toplevel block
write(src, "module X\n")
write(src, open(readstring, ARGS[1]))
write(src, "\nend")

const jlsrc = String(src)

tmp1 = tempname()
tmp2 = tempname()

global ast1, ast2

#ast = let
    #local ast1::Expr
    open("$tmp1", "w+") do io
        global ast1
        ast1 = Parser.parse(jlsrc) |> norm_ast
        AbstractTrees.print_tree(io, ast1)
    end

    #local ast2::Expr
    open("$tmp2", "w+") do io
        global ast2
        ast2 = parse(jlsrc) |> norm_ast
        AbstractTrees.print_tree(io, ast2)
    end

    if ast1 == ast2
        print_with_color(:green, "OK\n")
    else
        print_with_color(:red, "FAILED\n")
    end

    ast1
#end

#topast = Expr(:toplevel)
#append!(topast.args, ast.args)

#eval(Main, topast)

run(`vimdiff $tmp1 $tmp2`)


text:
#######################################################################
# States can be a simple Vector (for <:Lattice) or a SupercellState
#######################################################################
struct SupercellState{V,S<:Supercell,A<:OffsetArray{V}} <: AbstractVector{V}
    vector::A
    supercell::S
end

SupercellState(lat::Superlattice{E,L,T};
      type::Type{Tv} = Complex{T},
      vector = OffsetArray{orbitaltype(lat, Tv)}(undef, maskranges(lat))) where {E,L,T,Tv} =
    SupercellState(vector, lat.supercell)

displayeltype(s::SupercellState{V}) where {V<:Number} = V
displayeltype(s::SupercellState{V}) where {T,N,V<:SVector{N,T}} = T

displayorbsize(s::SupercellState{V}) where {V<:Number} = 1
displayorbsize(s::SupercellState{V}) where {T,N,V<:SVector{N,T}} = N

Base.show(io::IO, s::SupercellState) = show(io, MIME("text/plain"), s)
Base.show(io::IO, ::MIME"text/plain", s::SupercellState{V,S}) where {V,L,L´,S<:Supercell{L,L´}} =
    print(io,
"SupercellState{$L} : state of an $(L)D Superlattice
  Element type     : $(displayeltype(s))
  Max orbital size : $(displayorbsize(s))
  Sites            : $(nsites(s.supercell))")

function floatsorbs(V::Type{<:SVector}, T)
    n, r = divrem(sizeof(eltype(V)), sizeof(T))
    N = length(V) # orbitals
    r == 0 || throw(
        error("Unexpected error: cannot reinterpret orbital type $V as a number of floats"))
    return n, N
end

function floatsorbs(V::Type{<:Number}, T)
    n, r = divrem(sizeof(V), sizeof(T))
    N = 1 # orbitals
    r == 0 || throw(
        error("Unexpected error: cannot reinterpret orbital type $V as a number of floats"))
    return n, N
end

function maybe_wrapstate(sv, lat::Superlattice)
    o = OffsetArray(sv, maskranges(lat))
    return SupercellState(o, lat.supercell)
end
maybe_wrapstate(sv, lat::Lattice) = sv

# External API #

# build a random state, with zeroed-out padding orbitals, throughout the supercell/unitcell
function randomstate(h::Hamiltonian{LA}; type::Type{Tv} = Complex{T}) where {E,L,T,Tv,LA<:AbstractLattice{E,L,T}}
    lat = h.lattice
    masked = ismasked(lat)
    V = orbitaltype(h.orbitals, type)
    n, N = floatsorbs(V, T)
    masksize = length.(maskranges(lat))
    norbs = length.(h.orbitals)
    v = rand(T, n * N, masksize...) # for performance, use n×N Floats to build an S
    if !all(x -> x == norbs[1], norbs)  # zero out missing orbitals
        @inbounds for c in CartesianIndices(masksize)
            site = first(Tuple(c))
            insupercell = !masked || parent(lat.supercell.mask)[c]
            norb = norbs[sublat(lat, site)] * insupercell
            for j in 1:N, i in 1:n
                v[i + (j-1)*n, Tuple(c)...] =
                    (v[i + (j-1)*n, Tuple(c)...] - T(0.5)) * (j <= norb)
            end
        end
    end
    rmul!(v, inv(norm(v))) # normalize! without needing to cast v as vector
    rv = reinterpret(V, v)
    sv = V[rv[1, c] for c in CartesianIndices(masksize)]
    return maybe_wrapstate(sv, lat)
end

Base.copy!(t::S, s::S) where {S<:SupercellState} = SupercellState(copy!(t.vector, s.vector), s.supercell)
Base.copy(s::SupercellState) = SupercellState(copy(s.vector), s.supercell)
Base.similar(s::SupercellState) = SupercellState(similar(s.vector),  s.supercell)
Base.size(s::SupercellState, i...) = size(s.vector, i...)
Base.length(s::SupercellState) = length(s.vector)
Base.getindex(s::SupercellState, i...) = getindex(s.vector, i...)

######################################################################
# mul!
######################################################################
function SparseArrays.mul!(t::S, hb::SupercellBloch, s::S, α::Number = true, β::Number = false) where {V,S<:SupercellState{V}}
    C = t.vector
    B = s.vector
    ham = hb.hamiltonian
    phases = hb.phases
    cells = s.supercell.cells
    cols = 1:size(first(ham.harmonics).h, 2)
    pinvint = pinvmultiple(s.supercell.matrix)
    zeroV = zero(V)
    # Scale target by β
    if β != 1
        β != 0 ? rmul!(C, β) : fill!(C, zeroV)
    end
    # Add α * blochphase * h * source to target
    Threads.@threads for ic in cells
        i = Tuple(ic)
        # isemptycell(s, i) && continue # good for performance? No much
        for h in ham.harmonics
            olddn = h.dn + SVector(i)
            newdn = new_dn(olddn, pinvint)
            j = Tuple(wrap_dn(olddn, newdn, s.supercell.matrix))
            CartesianIndex(j) in cells || continue # boundaries in unwrapped directions
            α´ = α * cis(phases' * newdn)
            nzv = nonzeros(h.h)
            rv = rowvals(h.h)
            for col in cols
                @inbounds αxj = B[col, i...] * α´
                for p in nzrange(h.h, col)
                    @inbounds C[rv[p], j...] += applyfield(ham.field, nzv[p], rv[p], col, h.dn) * αxj
                end
            end
        end
    end
    # Filter out sites not in supercell
    @simd for j in eachindex(t.vector)
        @inbounds isinmask(s.supercell, j) || (t.vector[j] = zeroV)
    end
    return t
end


text:
#-----------------------------------------------------------------------
#   RTC21
#-----------------------------------------------------------------------

#-----------------------------------------------------------------------
#   FIVE POINT LOCATOR FUNCTION (RTC21)
#-----------------------------------------------------------------------
"""
    locator5RTC21( X::RealMtx )

The function presented here implements the five point relativistic
location formula given by Ruggiero, Tartaglia, and Casalino in Ruggiero
et al., arxiv:2111.13423 (2021).

"""     #---------------------------------------------------------------
function locator5RTC21( X::RealMtx )
    tpfl=typeof(X[1,1])
    l = size(X)

    if l[2] >= 5 || l[1] != 4
        (t1,t2,t3,t4,t5) = (X[1,1],X[1,2],X[1,3],X[1,4],X[1,5])
        (x1,x2,x3,x4,x5) = (X[2,1],X[2,2],X[2,3],X[2,4],X[2,5])
        (y1,y2,y3,y4,y5) = (X[3,1],X[3,2],X[3,3],X[3,4],X[3,5])
        (z1,z2,z3,z4,z5) = (X[4,1],X[4,2],X[4,3],X[4,4],X[4,5])

        C  = [  t1-t2   x2-x1   y2-y1   z2-z1    ;
                t2-t3   x3-x2   y3-y2   z3-z2    ;
                t3-t4   x4-x3   y4-y3   z4-z3    ;
                t4-t5   x5-x4   y5-y4   z5-z4    ]

        B  = [  x2^2 - x1^2 + y2^2 - y1^2 + z2^2 - z1^2 + t1^2 - t2^2  ;
                x3^2 - x2^2 + y3^2 - y2^2 + z3^2 - z2^2 + t2^2 - t3^2  ;
                x4^2 - x3^2 + y4^2 - y3^2 + z4^2 - z3^2 + t3^2 - t4^2  ;
                x5^2 - x4^2 + y5^2 - y4^2 + z5^2 - z4^2 + t4^2 - t5^2  ]

        return inv(2 .* C) * B
    else
        return zeros(tpfl,4)
    end
end      #---------------------------------------------------------------


text:
module ParallelBenchmarks

include(joinpath(dirname(@__FILE__), "..", "utils", "RandUtils.jl"))

using .RandUtils
using BenchmarkTools
using Compat

if VERSION >= v"0.7.0-DEV.2954"
    using Distributed
end

const SUITE = BenchmarkGroup()

#################################################
# Echoing data between processes (Issue #14467) #
#################################################

if nprocs() > 1
    function otherid(id)
        workers = procs()
        return workers[findfirst(w -> w != id, workers)]
    end
    g = addgroup!(SUITE, "remotecall", ["io", "remotecall_fetch"])
    for s in (2, 64, 512, 1024, 4096)
        z = zeros(UInt8, s)
        i = otherid(myid())
        g["identity", s] = @benchmarkable remotecall_fetch(identity, $i, $z)
    end
end

##################
# Multithreading #
##################

# Threading is broken for now, so we don't have tuned parameters for this yet
# if VERSION >= v"0.5.0-dev+923" && Base.Threads.nthreads() > 1
#     include("Laplace3D.jl")
#     include("ThreadedStockCorr.jl")
#     include("LatticeBoltzmann.jl")
#     g = addgroup!(SUITE, "multithread", ["thread", "laplace", "laplacian"])
#     g["laplace3d"] = @benchmarkable Laplace3D.perf_laplace3d()
#     g["pstockcorr"] = @benchmarkable ThreadedStockCorr.perf_pstockcorr(10^4)
#     g["lattice_boltzmann"] = @benchmarkable LatticeBoltzmann.perf_lattice_boltzmann(36)
# end

end # module


text:
export get_processed_data, get_test_data

using Images, Flux, MLDatasets
using Flux: onehotbatch
using Zygote: @nograd
using Base.Iterators: partition
using Random
using Statistics

function get_test_data(get_proportion = 1f0, batchsize = 0, random_seed = -1)
	if random_seed > -1
	    Random.seed!(random_seed)
	end
    test_x, test_y = MLDatasets.CIFAR10.testdata(Float32)
	total_img = Int(floor(length(test_y)*get_proportion))
	order = 1:total_img
	X = test_x[:,:,:,order]
	y = test_y[order]
	imgs = [X[:,:,:,i] for i in 1:total_img]
	labels = Matrix(onehotbatch([y[i] for i in 1:total_img],0:9))
	if batchsize == 0
	    test = (imgs,[labels[:,i] for i in 1:total_img])
	else
		test = [(cat(imgs[i]..., dims = 4), labels[:,i]) for i in partition(1:total_img, batchsize)]
	end
    return test
end

function get_processed_data(splitr = 0.5f0, batchsize = 64, mini = 1f0, val_batchsize = 0, random_seed = -1)
	if random_seed > -1
	    Random.seed!(random_seed)
	end
	train_x, train_y = MLDatasets.CIFAR10.traindata(Float32)
    total_img = Int(floor(length(train_y)*mini))
	order = shuffle(1:total_img)
    X = train_x[:,:,:,order]
	y = train_y[order]
	imgs = [X[:,:,:,i] for i in 1:total_img]
    labels = Matrix(onehotbatch([y[i] for i in 1:total_img],0:9))
    train_pop = Int(floor((1-splitr)* total_img))
    train = [(cat(imgs[i]..., dims = 4), labels[:,i]) for i in partition(1:train_pop, batchsize)]
    if train_pop < total_img
        if val_batchsize == 0
	    	val_batchsize = batchsize
        end
        val = [(cat(imgs[i]..., dims = 4), labels[:,i]) for i in partition(train_pop+1:total_img, val_batchsize)]
    else
        val = []
    end
    return train, val
end

max_stars_repo_name: huangyxi/DifferentiableNAS

max_issues_repo_name: huangyxi/DifferentiableNAS

max_forks_repo_name: huangyxi/DifferentiableNAS


text:
# """
#     MOI.add_constraint(model::Optimizer, sgvar::MOI.VectorOfVariables, set::VariablesEquality)

# Interface function which add a constraint to the model (which is himself an Optimizer)
# This constraints a single variable to be different or equal to another variable.
# """
function MOI.add_constraint(model::Optimizer, vectOfVar::MOI.VectorOfVariables, ::NotEqualSet)
    id1, id2 = vectOfVar.variables[1], vectOfVar.variables[2]

    ci = MOI.ConstraintIndex{MOI.VectorOfVariables, NotEqualSet}(length(model.moimodel.constraints) + 1)
    moiConstraint = MOIConstraint(NotEqualSet, (id1, id2), ci)
    push!(model.moimodel.constraints, moiConstraint)

    return ci
end

function create_CPConstraint(moiconstraint::MOIConstraint{NotEqualSet}, optimizer::Optimizer)
    id1, id2 = moiconstraint.args
    x, y = get_cp_variable(optimizer, id1), get_cp_variable(optimizer, id2)
    NotEqual(x, y, optimizer.cpmodel.trailer)
end

"""
    MOI.add_constraint(model::Optimizer, sgvar::MOI.SingleVariable, set::MOI.LessThan)

Interface function which add a constraint to the model (which is himself an Optimizer).
This constraints a single variable to be less than a given Integer.
"""
function MOI.add_constraint(model::Optimizer, sgvar::MOI.SingleVariable, set::MOI.LessThan{Float64})
    # get the VariableIndex and convert it to string
    id = sgvar.variable.value
    if isnothing(model.moimodel.variables[id].max) || model.moimodel.variables[id].max >= set.upper
        model.moimodel.variables[id].max = set.upper
    end
    return MOI.ConstraintIndex{MOI.SingleVariable, MOI.LessThan{Float64}}(1)
end

function MOI.add_constraint(model::Optimizer, saf::MOI.ScalarAffineFunction{Float64}, set::MOI.LessThan{Float64})
    if length(saf.terms) != 2 || saf.terms[1].coefficient != 1. || saf.terms[2].coefficient != -1. || saf.constant != set.upper
        moiaff = MOIAffineFunction(nothing, saf)
        push!(model.moimodel.affines, moiaff)

        newId = length(model.moimodel.constraints) + 1
        ci = MOI.ConstraintIndex{MOI.ScalarAffineFunction{Float64}, MOI.LessThan{Float64}}(newId)
        constraint = MOIConstraint(SeaPearl.LessOrEqualConstant, (AffineIndex(length(model.moimodel.affines)), set.upper), ci)
        push!(model.moimodel.constraints, constraint)

        return ci

    end

    saf.constant -= set.upper

    newId = length(model.moimodel.constraints) + 1
    ci = MOI.ConstraintIndex{MOI.ScalarAffineFunction{Float64}, MOI.LessThan{Float64}}(newId)
    constraint = MOIConstraint(SeaPearl.LessOrEqual, (saf.terms[1].variable_index, saf.terms[2].variable_index), ci)
    push!(model.moimodel.constraints, constraint)

    return ci
end

function create_CPConstraint(moiconstraint::MOIConstraint{SeaPearl.LessOrEqual}, optimizer::Optimizer)
    if !isa(moiconstraint.ci, MOI.ConstraintIndex{MOI.ScalarAffineFunction{Float64}, MOI.LessThan{Float64}})
        return nothing
    end
    x = get_cp_variable(optimizer, moiconstraint.args[1])
    y = get_cp_variable(optimizer, moiconstraint.args[2])

    LessOrEqual(x, y, optimizer.cpmodel.trailer)
end

function create_CPConstraint(moiconstraint::MOIConstraint{SeaPearl.LessOrEqualConstant}, optimizer::Optimizer)
    x = get_cp_variable(optimizer, moiconstraint.args[1])
    c = Int(floor(moiconstraint.args[2]))

    removeAbove!(x.domain, c)
    nothing
end

"""
    MOI.add_constraint(model::Optimizer, sgvar::MOI.SingleVariable, set::MOI.GreaterThan)

Give an upper bound to `sgvar`
"""
function MOI.add_constraint(model::Optimizer, sgvar::MOI.SingleVariable, set::MOI.GreaterThan{Float64})
    # get the VariableIndex and convert it to string
    id = sgvar.variable.value
    if isnothing(model.moimodel.variables[id].min) || model.moimodel.variables[id].min <= set.lower
        model.moimodel.variables[id].min = set.lower
    end
    return MOI.ConstraintIndex{MOI.SingleVariable, MOI.GreaterThan{Float64}}(1)
end

function MOI.add_constraint(model::Optimizer, saf::MOI.ScalarAffineFunction{Float64}, set::MOI.EqualTo{Float64})
    saf.constant -= set.value
    moiaff = MOIAffineFunction(nothing, saf)
    push!(model.moimodel.affines, moiaff)

    newId = length(model.moimodel.constraints) + 1
    ci = MOI.ConstraintIndex{MOI.ScalarAffineFunction{Float64}, MOI.EqualTo{Float64}}(newId)
    constraint = MOIConstraint(MOI.EqualTo, (AffineIndex(length(model.moimodel.affines)),), ci)
    push!(model.moimodel.constraints, constraint)

    return ci
end

function create_CPConstraint(moiconstraint::MOIConstraint{MOI.EqualTo}, optimizer::Optimizer)
    x = get_cp_variable(optimizer, moiconstraint.args[1])

    @assert 0 in x.domain "Infeasible problem"
    assign!(x, 0)
    nothing
end


text:
# include("Viability.jl")
@everywhere include("Viability.jl")

@everywhere using Viability:
    COORDINATE_TYPE,
    STATE_TYPE,
    create_normalized_2d_grid,
    map_over_points,
    rhs2stepfct,
    sp_mparallel,
    sp_parallel,
    sp_viability_kernel,
    spatial_normalization,
    temporal_homogenization,
    points2bounds,
    write_result_file

OUTFILENAME = "consum.txt"

dim = 2
num_per_dim = unsigned(201)
stepsize=1/(num_per_dim-1)
println("stepsize $stepsize")
delta_t = 7*stepsize
println("delta_t $delta_t")
delta_ball = 1.5*stepsize
println("delta_ball $delta_ball")
num_total = num_per_dim ^ 2
bounds = COORDINATE_TYPE[0 2; 0 3]

estimation_args = Dict{String,Any}(
    "num-per-dim" => num_per_dim,
    "bounds" => bounds,
    "delta_t" => delta_t,
    "delta_ball" => delta_ball
)

model_params = Dict{String,Any}(
    "u+" => 0.5,
    "u-" => -0.5
)

model_info = Dict{String,Any}(
    "model" => "consum",
    "dimension" => 2,
    "estimation-args" => estimation_args,
    "model-params" => model_params
)

@everywhere RET_ARRAY=COORDINATE_TYPE[0, 0]
@everywhere function consum_rhs(xy, u)
    RET_ARRAY[1] = xy[1] - xy[2]
    RET_ARRAY[2] = u
    return RET_ARRAY
end

rhss = [x->consum_rhs(x, model_params["u-"]), x->consum_rhs(x,  model_params["u+"])]
rhss_transformed =
    map(temporal_homogenization,
        map(
            f->spatial_normalization(f,bounds),
            rhss
        )
    )
step_functions = map(f->rhs2stepfct(f, delta_t), rhss_transformed)

function normalized_is_sunny(x)
    # set only the boundary dark
    lower = stepsize / 2
    upper = 1-lower
    x[1] > lower && x[2] > lower && x[1] < upper && x[2] < upper
end

mode = sp_mparallel
# mode = sp_parallel

points = create_normalized_2d_grid(num_per_dim)
states = map_over_points(STATE_TYPE, normalized_is_sunny, points)
@time ret = sp_viability_kernel(
    points,
    states,
    step_functions,
    delta=delta_ball,
    computation_mode=mode,
    r=delta_t
)

points = create_normalized_2d_grid(num_per_dim)
states = map_over_points(STATE_TYPE, normalized_is_sunny, points)
@time ret = sp_viability_kernel(
    points,
    states,
    step_functions,
    delta=delta_ball,
    computation_mode=mode,
    r=delta_t
)

actual_points = points2bounds(points, bounds)

print("writing output to '$OUTFILENAME' ... ")
write_result_file(OUTFILENAME, model_info, actual_points, states)
println("done")


text:
module Variable

import MathOptInterface
const MOI = MathOptInterface
const MOIU = MOI.Utilities
const MOIB = MOI.Bridges

import ComplexOptInterface
const COI = ComplexOptInterface

include("psd.jl")

end


text:
using GLAbstraction, GeometryTypes, ColorTypes, GLFW, Reactive, ModernGL
const GLA = GLAbstraction

const window = GLFW.Window(name="Example Geometry Shader")
GLFW.MakeContextCurrent(window)
GLA.set_context!(window)

const vert = GLA.vert"""
#version 330
in vec2 pos;
out vec2 g_position;

void main() {
    g_position = pos;
}
"""

const frag = GLA.frag"""
#version 330

out vec4 outColor;

void main() {
    outColor = vec4(1,0,0,1);
}
"""

const geom = GLA.geom"""
#version 330

layout(points) in;
layout(triangle_strip, max_vertices = 4) out;
in vec2 g_position[];

void main(void)
{
  vec2 p0 = g_position[0];

  gl_Position = vec4(p0, 0, 1);
  EmitVertex();

  gl_Position = vec4(p0+vec2(0,0.1), 0, 1);
  EmitVertex();

  gl_Position = vec4(p0+vec2(0.1,0), 0, 1);
  EmitVertex();

  gl_Position = vec4(p0+vec2(0.1), 0, 1);
  EmitVertex();

  EndPrimitive();
}
"""

const b = Point2f0[(-0.5,0),(0.0, 0.0),(0.4, 0.3)]

program = GLA.Program(vert, geom, frag)
robj = GLA.VertexArray(GLA.generate_buffers(program, GLA.GEOMETRY_DIVISOR, pos=b))

glClearColor(0,0,0,1)

GLA.bind(program)
while !GLFW.WindowShouldClose(window)
    glClear(GL_COLOR_BUFFER_BIT | GL_DEPTH_BUFFER_BIT)
    GLA.bind(robj)
    GLA.draw(robj)
    GLFW.SwapBuffers(window)
    GLFW.PollEvents()
    if GLFW.GetKey(window, GLFW.KEY_ESCAPE) == GLFW.PRESS
        GLFW.SetWindowShouldClose(window, true)
    end
end
GLFW.DestroyWindow(window)