text:
#=
Here I show the exact correspondence
between my spiking implementation and code by
Litwin-Kumar / Auguste Schulz
For one large random network
and inputs in currents (high freq input still missing)
=#

push!(LOAD_PATH, abspath(@__DIR__,".."))

using LinearAlgebra,Statistics,StatsBase,Distributions
using Plots,NamedColors ; theme(:dark)
using SparseArrays
using SpikingRNNs; const global S = SpikingRNNs
using BenchmarkTools
using ProgressMeter

function onesparsemat(w::Real)
  return sparse(fill(w,(1,1)))
end

##

using Random; Random.seed!(0)

function initweights(Ne, Ni, p, Jee0, Jei0, Jie, Jii)
	""" initialise the connectivity weight matrix based on initial
	E-to-E Jee0 initial weight e to e plastic
	I-to-E Jei0 initial weight i to e plastic
	E-to-I Jie0 constant weight e to i not plastic
	I-to-I Jii constant weight i to i not plastic"""
	Ncells = Ne+Ni

	# initialise weight matrix
	# w[i,j] is the weight from pre- i to postsynaptic j neuron
	weights = zeros(Float64,Ncells,Ncells)
	weights[1:Ne,1:Ne] .= Jee0
	weights[1:Ne,(1+Ne):Ncells] .= Jie
	weights[(1+Ne):Ncells,1:Ne] .= Jei0
	weights[(1+Ne):Ncells,(1+Ne):Ncells] .= Jii
	# set diagonal elements to 0
	# for cc = 1:Ncells
	# 	weights[cc,cc] = 0
	# end
	weights[diagind(weights)] .= 0.0
	# ensure that the connection probability is only p
	weights = weights.*(rand(Ncells,Ncells) .< p)
	return weights
end
function weightpars(;Ne = 4000, Ni = 1000, p = 0.2 )
	"""Ne, Ni number of excitatory, inhibitory neurons
	p initial connection probability"""
	Jee0 = 2.86 #initial weight e to e plastic
	Jei0 = 48.7 #initial weight i to e plastic
	Jie = 1.27 #constant weight e to i not plastic
	Jii = 16.2 #constant weight i to i not plastic
	return Ne,Ni,p,Jee0,Jei0,Jie,Jii
end
initweights() = initweights(weightpars()...)

# this is a bare-bone version of the Schulz code
function runsimulation_static_simplified(Ne::Integer,Ni::Integer,
  weights::Array{Float64,2},spiketimes::Array{<:Integer,2},
  v_start::Vector{Float64},input_constant::Vector{Float64}; dt = 0.1, T = 2000)

		# Naming convention
		# e corresponds to excitatory
		# i corresponds to inhibitory
		# x corresponds to external

		#membrane dynamics
		taue = 30 #e membrane time constant ms
		taui = 20 #i membrane time constant ms
		vreste = -60 #e resting potential mV
		vresti = -70 #i resting potential mV
		vpeak = 20 #cutoff for voltage.  when crossed, record a spike and reset mV
		eifslope = 2 #eif slope parameter mV
		C = 300 #capacitance pF
		erev = 0 #e synapse reversal potential mV
		irev = -75 #i synapse reversal potntial mV
		vth0 = -52 #initial spike voltage threshold mV
		thrchange = false # can be switched off to have vth constant at vth0
		ath = 10 #increase in threshold post spike mV
		tauth = 30 #threshold decay timescale ms
		vresete = vreste #reset potential mV
    vreseti = vresti
		taurefrac = 15 #absolute refractory period ms

		# total number of neurons
		Ncells = Ne+Ni

    # constant input
    in_const = input_constant

		# synaptic kernel
		tauerise = 1 #e synapse rise time
		tauedecay = 6 #e synapse decay time
		tauirise = .5 #i synapse rise time
		tauidecay = 2 #i synapse decay time

		# Initialisation of spike arrays ------------------------------
		# spike count and spike times
		# Scope problems - changing arrays are therefore defined here directly

		totalspikes = zeros(Int,Ncells)
		totsp::Int64 = 0; # total numberof spikes
		spmax::Int64 = size(spiketimes,1); # maximum recorded spikes

		# further arrays for storing inputs and synaptic integration
		forwardInputsE = zeros(Float64,Ncells) #sum of all incoming weights from excitatpry inputs both external and EE and IE
		forwardInputsI = zeros(Float64,Ncells) #sum of all incoming weights from inhibitory inputs both II and EI
		forwardInputsEPrev = zeros(Float64,Ncells) #as above, for previous timestep
		forwardInputsIPrev = zeros(Float64,Ncells)

		xerise = zeros(Float64,Ncells) #auxiliary variables for E/I currents (difference of exponentials)
		xedecay = zeros(Float64,Ncells)
		xirise = zeros(Float64,Ncells)
		xidecay = zeros(Float64,Ncells)

		v = zeros(Float64,Ncells) #membrane voltage

		# initialisation of membrane potentials and poisson inputs
    copy!(v,v_start)
		# for cc = 1:Ncells
		# 	v[cc] = vreset # + (vth0-vreset)*rand()
		# end

		vth = vth0*ones(Float64,Ncells) #adaptive threshold
		lastSpike = -100*ones(Ncells) #last time the neuron spiked

		# ---------------------------------- set up storing avg weights -----------------

		# total number of simulation steps, steps when to normalise, and save
		Nsteps = round(Int,T/dt)

		# true time
		t::Float64 = 0.0

		# bool counter if a neuron has just had a spike
		spiked = zeros(Bool,Ncells)
  v_out = Matrix{Float64}(undef,Ncells,Nsteps)
	@time	for tt = 1:Nsteps

				if mod(tt,Nsteps/100) == 1  #print percent complete
					print("\r",round(Int,100*tt/Nsteps))
				end

				forwardInputsE[:] .= 0.
				forwardInputsI[:] .= 0.
				t = dt*tt

				fill!(spiked,zero(Bool)) # reset spike bool without new memory allocation

				for cc = 1:Ncells
					xerise[cc] += -dt*xerise[cc]/tauerise + forwardInputsEPrev[cc]
					xedecay[cc] += -dt*xedecay[cc]/tauedecay + forwardInputsEPrev[cc]
					xirise[cc] += -dt*xirise[cc]/tauirise + forwardInputsIPrev[cc]
					xidecay[cc] += -dt*xidecay[cc]/tauidecay + forwardInputsIPrev[cc]

					if cc <= Ne # excitatory
						if thrchange
						vth[cc] += dt*(vth0 - vth[cc])/tauth;
						end
					end

					if t > (lastSpike[cc] + taurefrac) #not in refractory period
						# update membrane voltage

						ge = (xedecay[cc] - xerise[cc])/(tauedecay - tauerise);
						gi = (xidecay[cc] - xirise[cc])/(tauidecay - tauirise);

						if cc <= Ne #excitatory neuron (eif), has adaptation
							dv = (vreste - v[cc] + eifslope*exp((v[cc]-vth[cc])/eifslope))/taue +
                  ge*(erev-v[cc])/C + gi*(irev-v[cc])/C +
                  in_const[cc]/taue
							v[cc] += dt*dv
							if v[cc] > vpeak
								spiked[cc] = true
							end
						else
							dv = (vresti - v[cc])/taui +
                 ge*(erev-v[cc])/C + gi*(irev-v[cc])/C +
                  in_const[cc]/taui
							v[cc] += dt*dv
							if v[cc] > vth0
								spiked[cc] = true
							end
						end

						if spiked[cc] #spike occurred
							spiked[cc] = true;
							v[cc] = cc<=Ne ? vresete : vreseti;
							lastSpike[cc] = t;
							totalspikes[cc] += 1;
							totsp += 1;
							if totsp < spmax
								spiketimes[totsp,1] = tt; # time index as a sparse way to save spiketimes
								spiketimes[totsp,2] = cc; # cell id
							elseif totsp == spmax
								spiketimes[totsp,1] = tt; # time index
								spiketimes[totsp,2] = cc; # cell id

								totsp = 0 # reset counter total number of spikes
							end

							if cc <= Ne && thrchange # only change for excitatory cells and when thrchange == true
								vth[cc] = vth0 + ath;
							end

							#loop over synaptic projections
							for dd = 1:Ncells # postsynaptic cells dd  - cc presynaptic cells
								if cc <= Ne #excitatory synapse
									forwardInputsE[dd] += weights[cc,dd];
								else #inhibitory synapse
									forwardInputsI[dd] += weights[cc,dd];
								end
							end

						end #end if(spiked)
					end #end if(not refractory)
        v_out[cc,tt] = v[cc]
			end #end loop over cells

			forwardInputsEPrev = copy(forwardInputsE)
			forwardInputsIPrev = copy(forwardInputsI)

		end # tt loop over time
	print("\r")

	println("simulation finished")

	return totalspikes,v_out

end

##

const Ne = 4000
const Ni = 1000
const Ttot = 60.0
const all_in_e = (60.0-52.5)
const all_in_i = (70.0-53.0)

weights = permutedims(initweights(weightpars(;Ne=Ne,Ni=Ni)...))

v_start = let n=Ne+Ni,
  (vmin,vmax) = (-60.0,-51.0)
  rand(Uniform(vmin,vmax),n)
end
const n_spikes_count  = 100_000
test_spiketimes = fill(-1,(n_spikes_count,2))

totspikes,vtest = let dt=0.1,
  _weights = permutedims(weights),
  input = vcat(fill(all_in_e,Ne),fill(all_in_i,Ni)),
  T=round(Integer,Ttot*1E3)
  runsimulation_static_simplified(Ne,Ni,_weights,test_spiketimes,v_start,
    input;dt=dt, T=T)
end
rates_test = totspikes ./ Ttot

_ = let plt = plot(),nneu=22,
  ts = (1:size(vtest,2)).*0.1E-3
  plot!(plt,ts,vtest[nneu,:];linewidth=2,xlims=(0,1))
end

histogram(rates_test[1:Ne])

##
#=

		#membrane dynamics
		taue = 30 #e membrane time constant ms
		taui = 20 #i membrane time constant ms
		vreste = -60 #e resting potential mV
		vresti = -70 #i resting potential mV
		vpeak = 20 #cutoff for voltage.  when crossed, record a spike and reset mV
		eifslope = 2 #eif slope parameter mV
		C = 300 #capacitance pF
		erev = 0 #e synapse reversal potential mV
		irev = -75 #i synapse reversal potntial mV
		vth0 = -52 #initial spike voltage threshold mV
		thrchange = false # can be switched off to have vth constant at vth0
		ath = 10 #increase in threshold post spike mV
		tauth = 30 #threshold decay timescale ms
		vreset = -70 #reset potential mV
		taurefrac = 50 #absolute refractory period ms
		aw_adapt = 4 #adaptation parameter a nS conductance
		bwfactor = 100
		bw_adapt = bwfactor*0.805 #adaptation parameter b pA current

    # constant input
    in_const = fill(0.0,Ncells)
    in_const[2] = -51.0 - vreste

		# synaptic kernel
		tauerise = 1 #e synapse rise time
		tauedecay = 6 #e synapse decay time
		tauirise = .5 #i synapse rise time
		tauidecay = 2 #i synapse decay time

=#

##
# time in seconds, voltage in mV
const dt = 0.1E-3
myτe = 30E-3 # seconds
myτi = 20E-3 # seconds
τrefr= 15E-3 # refractoriness
vth_e = 20.   # mV
vthexp = -52.0 # actual threshold for spike-generation
vth_i = vthexp
eifslope = 2.0
v_rest_e = -60.0
v_rest_i = -70.0
v_rev_e = 0.0
v_rev_i = -75.0
v_leak_e = v_rest_e
v_leak_i = v_rest_i
v_reset_e = v_rest_e
v_reset_i = v_rest_i
Cap = 300.0 #capacitance mF
in_const_e = all_in_e*Cap/myτe
in_const_i = all_in_i*Cap/myτi

# synaptic kernel
tauerise = 1E-3 #e synapse rise time
tauedecay = 6E-3 #e synapse decay time
taueplus,taueminus  = tauedecay, tauerise
tauirise = 0.5E-3 #i synapse rise time
tauidecay = 2E-3 #i synapse decay time
tauiplus,tauiminus = tauidecay,tauirise

##

nt_e = let sker = S.SKExpDiff(taueplus,taueminus)
  sgen = S.SpikeGenEIF(vthexp,eifslope)
  S.NTLIFConductance(sker,sgen,myτe,Cap,
   vth_e,v_reset_e,v_rest_e,τrefr,v_rev_e)
end
ps_e = S.PSLIFConductance(nt_e,Ne)

nt_i = let sker = S.SKExpDiff(tauiplus,tauiminus)
  sgen = S.SpikeGenNone()
  S.NTLIFConductance(sker,sgen,myτi,Cap,
   vth_i,v_reset_i,v_rest_i,τrefr,v_rev_i)
end
ps_i = S.PSLIFConductance(nt_i,Ni)

# static inputs
ps_input_e = S.PSSimpleInput(S.InputSimpleOffset(in_const_e))
ps_input_i = S.PSSimpleInput(S.InputSimpleOffset(in_const_i))

## connections
conn_ee = let w_ee = weights[1:Ne,1:Ne]
  S.ConnGeneralIF2(sparse(w_ee))
end
conn_ii = let w_ii = weights[Ne+1:end,Ne+1:end]
  S.ConnGeneralIF2(sparse(w_ii))
end
conn_ei = let w_ei = weights[1:Ne,Ne+1:end]
  S.ConnGeneralIF2(sparse(w_ei))
end
conn_ie = let w_ie = weights[Ne+1:end,1:Ne]
  S.ConnGeneralIF2(sparse(w_ie))
end

# input populations
pop_e = S.Population(ps_e,(S.FakeConnection(),ps_input_e),
   (conn_ee,ps_e),(conn_ei,ps_i) )
pop_i = S.Population(ps_i,(S.FakeConnection(),ps_input_i),
   (conn_ii,ps_i),(conn_ie,ps_e) )

##
# that's it, let's make the network
myntw = S.RecurrentNetwork(dt,pop_e,pop_i)

# record spiketimes and internal potential
krec = 1
n_e_rec = 1000
n_i_rec = 1000
t_wup = 0.0
rec_state_e = S.RecStateNow(ps_e,krec,dt,Ttot;idx_save=collect(1:n_e_rec),t_warmup=t_wup)
rec_state_i = S.RecStateNow(ps_i,krec,dt,Ttot;idx_save=collect(1:n_i_rec),t_warmup=t_wup)
rec_spikes_e = S.RecSpikes(ps_e,5.0,Ttot;idx_save=collect(1:n_e_rec),t_warmup=t_wup)
rec_spikes_i = S.RecSpikes(ps_i,5.0,Ttot;idx_save=collect(1:n_i_rec),t_warmup=t_wup)

## Run

times = (0:myntw.dt:Ttot)
nt = length(times)
# clean up
S.reset!.([rec_state_e,rec_spikes_e])
S.reset!.([rec_state_i,rec_spikes_i])
S.reset!.([ps_e,ps_i])
S.reset!(conn_ei)
# initial conditions
ps_e.state_now .= v_start[1:Ne]
ps_i.state_now .= v_start[Ne+1:end]

@time begin
  @showprogress 5.0 "network simulation " for (k,t) in enumerate(times)
    rec_state_e(t,k,myntw)
    rec_state_i(t,k,myntw)
    rec_spikes_e(t,k,myntw)
    rec_spikes_i(t,k,myntw)
    S.dynamics_step!(t,myntw)
  end
end

S.add_fake_spikes!(1.0vth_e,rec_state_e,rec_spikes_e)
S.add_fake_spikes!(0.0,rec_state_i,rec_spikes_i)
##

rates_e = let rdic=S.get_mean_rates(rec_spikes_e,dt,Ttot)
  ret = fill(0.0,n_i_rec)
  for (k,v) in pairs(rdic)
    ret[k] = v
  end
  ret
end
rates_i = let rdic=S.get_mean_rates(rec_spikes_i,dt,Ttot)
  ret = fill(0.0,n_i_rec)
  for (k,v) in pairs(rdic)
    ret[k] = v
  end
  ret
end

_ = let plt=plot(;leg=false),
  netest = n_e_rec
  scatter!(rates_e,rates_test[1:netest];ratio=1,
    xlabel="Dylan's model",ylabel="Suchulz et al",title="E rates")
  plot!(plt,identity; linewidth=2)
end

_ = let plt=plot(;leg=false),
  ntest = n_i_rec
  scatter!(rates_i,rates_test[Ne+1:Ne+ntest];ratio=1,
    xlabel="Dylan's model",ylabel="Suchulz et al",title="I rates")
  plot!(plt,identity;linewidth=2)
end

##

_ = let neu = 202,
  plt=plot()
  plot!(plt,rec_state_e.times,rec_state_e.state_now[neu,:];linewidth=2,leg=false,
    xlims=(0,1),ylims=(-70,-40))
  ts = (1:size(vtest,2)).*0.1E-3
  plot!(plt,ts,vtest[neu,:];linewidth=2,linestyle=:dash)
end
_ = let plt=plot(),neu=16
  plot!(plt,rec_state_i.times,rec_state_i.state_now[neu,:];linewidth=2,leg=false,
    xlims=(0,1))
  ts = (1:size(vtest,2)).*0.1E-3
  plot!(plt,ts,vtest[Ne+neu,:];linewidth=2,linestyle=:dash)
end


text:
using Printf

function makesofn(n::Integer)
    buf = Vector{typeof(n)}(0)
    i = 0
    return function sofn(item)
        i += 1
        if i ≤ n
            push!(buf, item)
        else
            j = rand(1:i)
            if j ≤ n buf[j] = item end
        end
        return buf
    end
end

nhist = zeros(Int, 10)
for _ in 1:10^5
    kas = makesofn(3)
    for j in 0:8 kas(j) end
    for k in kas(9) nhist[k+1] += 1 end
end

println("Simulating sof3(0:9) 100000 times:")
for (i, c) in enumerate(nhist)
    @printf("%5d → %5d\n", i-1, c)
end

ext: jl

max_stars_repo_name: ethansaxenian/RosettaDecode

max_issues_repo_name: ethansaxenian/RosettaDecode

max_forks_repo_name: ethansaxenian/RosettaDecode


text:
module OSF

using HTTP
using JSON

export osfget, osfput, osfcontents, osfprint, osfheaders, osftree, osfupload

const osfurl = "https://api.osf.io/v2/nodes"
const default_osfheaders = ["Content-Type" => "application/vnd.api+json"]

###############################################################################
# Returns a dictionary mapping each remote file with urls
function osftree(id, headers = default_osfheaders; osfurl = osfurl)
    url = joinpath(osfurl, id)
    storage = osfget(joinpath(url, "files", "osfstorage"), headers)
    storagetree = crawler(storage)

    rootstorage = osfget(joinpath(url, "files"), headers)
    storagetree["/"] = getlinks(first(rootstorage["data"]))
    return storagetree
end

# Used by osftree() to crawl over remote files
function crawler(storage)
    contents = Dict{String, Dict{String, String}}()
    files = filter(isfile, storage["data"])
    folders = filter(isfolder, storage["data"])
    for file in files
      contents[getname(file)] = getlinks(file)
    end
    for folder in folders
        contents[getname(folder)] = getlinks(folder)
        merge!(contents, crawler(osfget(getsubfolder(folder))))
    end
    nextpage = storage["links"]["next"]
    if !isnothing(nextpage)
        merge!(contents, crawler(osfget(nextpage)))
    end
    return contents
end

###############################################################################
# Uploads `files` to osf repository
function osfupload(files, remote_id; headers = default_osfheaders, rootdir = "")
    remotetree = osftree(remote_id)
    return osfupload(files, remote_id, remotetree; headers, rootdir)
end

function osfupload(files, remote_id, remotetree; headers = default_osfheaders, rootdir = "")
    for file in files
        remotefile = replace(file, rootdir => "")
        if remotefile ∈ getfilenames(remotetree)
            println("$remotefile already exists. Skipping.")
            continue
        end
        createfolder(dirname(remotefile), remotetree, headers)
        createfile(file, remotefile, remotetree, headers)
    end
    return remotetree
end

function createfolder(folder, remotetree, headers)
    folder = confirmpath(folder)
    subfolders = splitpath(folder)
    subpaths = cumulativepath(subfolders)

    lastpath = "/"
    for (subpath, subfolder) in zip(subpaths, subfolders)
        remotepaths = getfoldernames(remotetree)
        if subpath ∉ remotepaths
            println("Not found. Will create $subfolder")
            url = remotetree[lastpath]["new_folder"] * "&name=$subfolder"
            response = osfput(url, headers)
            remotetree[subpath] = response["data"]["links"]
        end
        lastpath = subpath
    end
end

function createfile(path, remotepath, remotetree, headers)
    folder = confirmpath(dirname(remotepath))
    filename = basename(remotepath)
    uploadurl = remotetree[folder]["upload"] * "?kind=file&name=$filename"
    response = open(path, "r") do io
        osfput(uploadurl, headers, io)
    end
    remotetree[remotepath] = response["data"]["links"]
end

function cumulativepath(pathparts)
    # If input is: ["/", "a", "b", "c"]
    # Then output is: ["/", "/a/, "/a/b/", "/a/b/c/"]
    out = String[first(pathparts)]
    for part in pathparts[2:end]
        newentry = joinpath(last(out), part)
        newentry = confirmpath(newentry)
        push!(out, newentry)
    end
    return out
end

###############################################################################
# Helper functions
osfget(url, headers = default_osfheaders) = JSON.parse(String(HTTP.request("GET", url, headers).body))
osfput(url, headers = default_osfheaders) = JSON.parse(String(HTTP.request("PUT", url, headers).body))
osfput(url, headers, io) = JSON.parse(String(HTTP.request("PUT", url, headers, io).body))
osfcontents(id, headers = default_osfheaders) = osf_get("$(api_url)/nodes/$id/files", headers)
osfprint(json) = JSON.print(json, 2)

isfile(response::Dict) = response["attributes"]["kind"] == "file"
isfile(path::AbstractString) = !isdirpath(path)
isfolder(response::Dict) = response["attributes"]["kind"] == "folder"
isfolder(path::AbstractString) = isdirpath(path)

confirmpath(path) = isdirpath(path) ? path : path * "/"
getname(response::Dict) = response["attributes"]["materialized_path"]
getlinks(response::Dict) = response["links"]
getsubfolder(response::Dict) = response["relationships"]["files"]["links"]["related"]["href"]
getfoldernames(storagetree::Dict) = (filter(isfolder, keys(storagetree)))
getfilenames(storagetree::Dict) = (filter(isfile, keys(storagetree)))

function get_names(response)
    names = String[]
    if isempty(response)
        return names
    end
    for item in response
        if haskey(item["attributes"], "materialized_path")
            push!(names, item["attributes"]["materialized_path"])
        end
    end
    return names
end

# Creates a header for HTTP requests, conains the authorization token
function osfheaders(tokenpath)
    token = open(tokenpath) do file
        chomp(read(file, String))
    end
    return ["Authorization" => "Bearer $token", "Content-Type" => "application/vnd.api+json"]
end

end


text:
# This is a special case and can't use @ch_struct because we don't
# know the server revision before reading this packet
struct ServerInfo
    server_name::String
    server_major_ver::VarUInt
    server_minor_ver::VarUInt
    server_rev::VarUInt

    # DBMS_MIN_REVISION_WITH_SERVER_TIMEZONE
    server_timezone::String
    # DBMS_MIN_REVISION_WITH_SERVER_DISPLAY_NAME
    server_display_name::String
    # DBMS_MIN_REVISION_WITH_VERSION_PATCH
    server_version_patch::VarUInt
end

function chread(sock::ClickHouseSock, ::Type{ServerInfo})
    server_name = chread(sock, String)
    server_major_ver = chread(sock, VarUInt)
    server_minor_ver = chread(sock, VarUInt)
    server_rev = chread(sock, VarUInt)
    rev = UInt64(server_rev)
    server_timezone = has_server_timezone(rev) ?
                        chread(sock, String) : "UTC"
    server_display_name = has_server_display_name(rev) ?
                        chread(sock, String) : ""
    server_version_patch = has_version_patch(rev) ?
                        chread(sock, VarUInt) : VarUInt(0)
    return ServerInfo(
            server_name,
            server_major_ver,
            server_minor_ver,
            server_rev,
            server_timezone,
            server_display_name,
            server_version_patch
            )
end

function chwrite(sock::ClickHouseSock, info::ServerInfo)
    chwrite(sock, info.server_name)
    chwrite(sock, info.server_major_ver)
    chwrite(sock, info.server_minor_ver)
    chread(sock, info.server_rev)
    rev = UInt64(info.server_rev)
    has_server_timezone(rev) && chwrite(sock, info.server_timezone)
    has_server_display_name(rev) && chwrite(sock, info.server_display_name)
    has_version_patch(rev) && chwrite(sock, info.server_version_patch)
end

@ch_struct struct ServerPong
end

@ch_struct struct ServerProgress
    rows::VarUInt
    bytes::VarUInt
    total_rows::VarUInt

    @has_client_write_info written_rows::VarUInt = VarUInt(0)
    @has_client_write_info written_bytes::VarUInt = VarUInt(0)
end

@ch_struct struct ServerProfileInfo
    rows::VarUInt
    blocks::VarUInt
    bytes::VarUInt
    applied_limit::Bool
    rows_before_limit::VarUInt
    calc_rows_before_limit::Bool
end

@ch_struct struct ServerEndOfStream
end

@ch_struct struct ServerTableColumns
    external_table_name::String
    columns::String
    sample_block::Block
end

@ch_struct struct ServerData
    data::Block
end

@ch_struct struct ServerTotals
    data::Block
end

@ch_struct struct ServerExtremes
    data::Block
end

struct ServerException
    code::UInt32
    name::String
    message::String
    stack_trace::String
    nested::Union{Nothing, ServerException}
end

@ch_struct struct ServerExceptionBase
    code::UInt32
    name::String
    message::String
    stack_trace::String
    has_nested::Bool
end

function chread(sock::ClickHouseSock, ::Type{ServerException})::ServerException
    base = chread(sock, ServerExceptionBase)
    nested = base.has_nested ? chread(sock, ServerException) : nothing
    ServerException(base.code, base.name, base.message, base.stack_trace, nested)
end

function cwrite(sock::ClickHouseSock, x::ServerException)::ServerException
    has_nested = !isnothing(x.nested)
    base = ServerExceptionBase(x.code, x.name, x.message, x.stack_trace, has_nested)
    chwrite(sock, base)
    has_nested && chwrite(sock, x.nested)
end


text:
using Pkg
using SafeTestsets

const GROUP = get(ENV, "GROUP", "All")

const is_APPVEYOR = Sys.iswindows() && haskey(ENV,"APPVEYOR")

const is_TRAVIS = haskey(ENV,"TRAVIS")

const is_CI = haskey(ENV,"CI")

@time begin
  if GROUP == "All" || GROUP == "NNODE"
      @time @safetestset "NNODE" begin include("NNODE_tests.jl") end
  end
  if !is_APPVEYOR && (GROUP == "All" || GROUP == "NNPDEHan")
      @time @safetestset "NNPDEHan" begin include("NNPDEHan_tests.jl") end
  end
  if GROUP == "All" || GROUP == "NNPDENS"
      @time @safetestset "NNPDENS" begin include("NNPDENS_tests.jl") end
  end
  if GROUP == "All" || GROUP == "NNPDE"
      @time @safetestset "NNPDE" begin include("NNPDE_tests.jl") end
  end
  if GROUP == "All" || GROUP == "NNKOLMOGOROV"
      @time @safetestset "NNKolmogorov" begin include("NNKolmogorov_tests.jl") end
  end
  if GROUP == "All" || GROUP == "NNSTOPPINGTIME"
      @time @safetestset "NNStopping" begin include("Stopping_tests.jl") end
  end
  if GROUP == "All" || GROUP == "NNRODE"
        @time @safetestset "NNRODE" begin include("NNRODE_tests.jl") end
  end

  if !is_APPVEYOR && GROUP == "GPU"
     @safetestset "NNPDE_gpu" begin include("NNPDE_tests_gpu.jl") end
 end
end


text:
################################################################################
#
#       DedekindCriterion.jl : Dedekinds Criterion for maximality
#
# This file is part of Hecke.
#
# Copyright (c) 2015, 2016: Claus Fieker, Tommy Hofmann
# All rights reserved.
#
# Redistribution and use in source and binary forms, with or without
# modification, are permitted provided that the following conditions are met:
# * Redistributions of source code must retain the above copyright notice, this
#   list of conditions and the following disclaimer.
#
# * Redistributions in binary form must reproduce the above copyright notice,
#   this list of conditions and the following disclaimer in the documentation
#   and/or other materials provided with the distribution.
#
# THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
# AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
# IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
# DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE LIABLE
# FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
# DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR
# SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER
# CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY,
# OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
# OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
#
#
#  Copyright (C) 2015, 2016 Tommy Hofmann
#
################################################################################

function dedekind_ispmaximal(O::NfOrd, p::fmpz)
  !isequationorder(O) && error("Order must be an equation order")
  if rem(discriminant(O), p) != 0
    return true, O
  end

  pmaximal = true

  Zy, y = PolynomialRing(ZZ, "y")
  Kx, x = PolynomialRing(ResidueRing(ZZ, p, cached=false), "x", cached=false)

  f = nf(O).pol

  Zyf = Zy(f)

  fmodp = Kx(Zyf)

  fac = factor(fmodp)

  g = Zy(1)

  # first build 1/p ( f - prod_{g in fac} g^?)
  for (fi, ei) in fac
    g *= lift(Zy, fi)^ei
  end

  g = Zyf - g

  for i in 0:degree(g)
    q,r = divrem(coeff(g,i),p)
    @assert r == 0
    setcoeff!(g,i,q)
  end

  gmodp = Kx(g)

  U = one(Kx)

  for (fi,ei) in fac
    if ei != 1 && rem(gmodp, fi) == zero(Kx)
      U *= fi
      return false
    end
  end

  return true
end

function dedekind_ispmaximal(O::NfOrd, p::Integer)
  return dedekind_ispmaximal(O, ZZ(p))
end

function dedekind_poverorder(O::NfOrd, p::fmpz)
  _, O = dedekind_test(O, p)
  return O
end

function dedekind_poverorder(O::NfOrd, p::Integer)
  return dedekind_poverorder(O, ZZ(p))
end

function dedekind_test(O::NfOrd, p::fmpz)
  !isequationorder(O) && error("Order must be an equation order")

  if rem(discriminant(O), p) != 0
    return true, O
  end

  pmaximal = true

  Zy, y = PolynomialRing(ZZ, "y")
  Kx, x = PolynomialRing(ResidueRing(ZZ, p, cached=false), "x", cached=false)

  f = nf(O).pol

  Zyf = Zy(f)

  fmodp = Kx(Zyf)

  fac = factor(fmodp)

  g = Zy(1)

  # first build 1/p ( f - prod_{g in fac} g^?)
  for (fi,ei) in fac
    g *= lift(Zy, fi)^ei
  end

  g = Zyf - g

  for i in 0:degree(g)
    q,r = divrem(coeff(g,i),p)
    @assert r == 0
    setcoeff!(g,i,q)
  end

  gmodp = Kx(g)

  U = one(Kx)

  for (fi, ei) in fac
    if ei != 1 && rem(gmodp, fi) == zero(Kx)
      U *= fi
      pmaximal=false
    end
  end

  @assert rem(fmodp, U) == zero(Kx)
  U = divexact(fmodp, U)

  if pmaximal
    return true, O
  end

  @assert rem(discriminant(O),p^2) == 0

  alpha = nf(O)(parent(f)(lift(Zy,U)))

  # build the new basis matrix
  # we have to take the representation matrix of alpha!
  # concatenating the coefficient vector won't help
  n = vcat(num(basis_mat(O))*p,representation_mat(alpha))

  b = FakeFmpqMat(n,p)

  OO = Order(nf(O), sub(hnf(b),degree(O) + 1:2*degree(O), 1:degree(O)))

  OO.isequationorder = false

  OO.disc = divexact(discriminant(O),p^(2*(degree(O)-degree(U))))

  return false, OO
end

function dedekind_test(O::NfOrd, p::Integer)
  return dedekind_test(O, ZZ(p))
end


text:
# ------------ Choice sites ------------ #

@inline function (ctx::RegenerateContext)(call::typeof(trace),
                                          addr::T,
                                          d::Distribution{K}) where {T <: Address, K}
    visit!(ctx, addr)
    in_prev_chm = has_value(get_trace(ctx.prev), addr)
    in_sel = haskey(ctx.target, addr)

    if in_prev_chm
        prev = get_sub(get_trace(ctx.prev), addr)
    end

    if in_sel && in_prev_chm
        ret = rand(d)
        set_sub!(ctx.discard, addr, prev)
    elseif in_prev_chm
        ret = prev.val
    else
        ret = rand(d)
    end

    score = logpdf(d, ret)
    in_prev_chm && increment!(ctx, score - get_score(prev))
    add_choice!(ctx, addr, score, ret)
    return ret
end

# ------------ Learnable ------------ #

@inline function (ctx::RegenerateContext)(fn::typeof(learnable), addr::Address)
    visit!(ctx, addr)
    haskey(ctx.params, addr) && return getindex(ctx.params, addr)
    error("(learnable): parameter not provided at address $addr.")
end

# ------------ Call sites ------------ #

@inline function (ctx::RegenerateContext)(c::typeof(trace),
                                          addr::T,
                                          call::Function,
                                          args...) where T <: Address
    visit!(ctx, addr)
    ps = get_sub(ctx.params, addr)
    ss = get_sub(ctx.target, addr)
    if has_sub(get_trace(ctx.prev), addr)
        prev_call = get_prev(ctx, addr)
        ret, cl, w, retdiff, d = regenerate(ss, ps, prev_call, UnknownChange(), args...)
    else
        ret, cl, w = generate(ss, ps, call, args...)
    end
    add_call!(ctx, addr, cl)
    increment!(ctx, w)
    return ret
end

@inline function (ctx::RegenerateContext)(c::typeof(trace),
                                          addr::T,
                                          call::G,
                                          args...) where {G <: GenerativeFunction, T <: Address}
    visit!(ctx, addr)
    ps = get_sub(ctx.params, addr)
    ss = get_sub(ctx.target, addr)
    if has_sub(get_trace(ctx.prev), addr)
        prev_call = get_prev(ctx, addr)
        ret, cl, w, retdiff, d = regenerate(ss, ps, prev_call, UnknownChange(), args...)
    else
        ret, cl, w = generate(ss, ps, call.fn, args...)
    end
    add_call!(ctx, addr, cl)
    increment!(ctx, w)
    return ret
end

# ------------ Rewalk and update weights/prune choices ------------ #

function regenerate_projection_walk(tr::DynamicTrace,
                                    visited::Visitor)
    weight = 0.0
    for (k, v) in shallow_iterator(tr)
        if !(k in visited)
            weight += projection(v, SelectAll())[1]
        end
    end
    weight
end

function regenerate_discard_walk!(d::DynamicDiscard,
                                  visited::Visitor,
                                  prev::DynamicTrace)
    for (k, v) in shallow_iterator(prev)
        if !(k in visited)
            ss = get_sub(visited, k)
            if isempty(ss)
                set_sub!(d, k, v)
            else
                sd = get_sub(d, k)
                sd = isempty(sd) ? DynamicMap{Value}() : sd
                discard_walk!(sd, ss, v)
                set_sub!(d, k, sd)
            end
        end
    end
end

# ------------ Convenience ------------ #

function regenerate(ctx::RegenerateContext, cs::DynamicCallSite, args::NTuple{N, Diffed}) where N
    ret = ctx(cs.fn, tupletype(args...), args...)
    adj_w = regenerate_projection_walk(ctx.tr, ctx.visited)
    regenerate_discard_walk!(ctx.discard, ctx.visited, get_trace(cs))
    return ret, DynamicCallSite(ctx.tr, ctx.score - adj_w, cs.fn, map(a -> unwrap(a), args), ret), ctx.weight, UnknownChange(), ctx.discard
end

# Explicitly parametrize the compilation pipeline.
function regenerate(opt::J, sel::L, ps::P, cs::DynamicCallSite, args::Diffed...) where {J <: CompilationOptions, L <: AddressMap, P <: AddressMap}
    ctx = Regenerate(opt, sel, ps, cs, DynamicTrace(), DynamicDiscard())
    return regenerate(ctx, cs, args)
end

# Explicitly parametrize the compilation pipeline.
function regenerate(opt::J, sel::L, ps::P, cs::DynamicCallSite, args...) where {J <: CompilationOptions, L <: AddressMap, P <: AddressMap}
    ctx = Regenerate(opt, sel, ps, cs, DynamicTrace(), DynamicDiscard())
    args = map(args) do a
        a isa Diffed ? a : Δ(a, NoChange())
    end
    return regenerate(ctx, cs, args)
end

function regenerate(sel::L, ps::P, cs::DynamicCallSite, args...) where {L <: AddressMap, P <: AddressMap}
    ctx = Regenerate(DefaultPipeline(), sel, ps, cs, DynamicTrace(), DynamicDiscard())
    args = map(args) do a
        a isa Diffed ? a : Δ(a, NoChange())
    end
    return regenerate(ctx, cs, args)
end

function regenerate(sel::L, cs::DynamicCallSite, args...) where L <: AddressMap
    ctx = Regenerate(DefaultPipeline(), sel, Empty(), cs, DynamicTrace(), DynamicDiscard())
    args = map(args) do a
        a isa Diffed ? a : Δ(a, NoChange())
    end
    return regenerate(ctx, cs, args)
end

# ------------ Documentation ------------ #

@doc(
"""
```julia
mutable struct RegenerateContext{T <: Trace,
                                 L <: Target,
                                 P <: AddressMap} <: ExecutionContext
    prev::T
    tr::T
    target::L
    weight::Float64
    score::Float64
    discard::T
    visited::Visitor
    params::P
end
```

Inner constructors:

```julia
function RegenerateContext(tr::T, sel::Vector{Address}) where T <: Trace
    un_sel = targetion(sel)
    new{T, typeof(un_sel), EmptyAddressMap}(tr, Trace(), un_sel, 0.0, Trace(), Visitor(), AddressMap())
end
function RegenerateContext(tr::T, sel::L) where {T <: Trace, L <: Target}
    new{T, L, EmptyAddressMap}(tr, Trace(), sel, 0.0, Trace(), Visitor(), AddressMap())
end
```

Outer constructors:

```julia
Regenerate(tr::Trace, sel::Vector{Address}) = RegenerateContext(tr, sel)
Regenerate(tr::Trace, sel::Target) = RegenerateContext(tr, sel)
```

The `RegenerateContext` is used for MCMC algorithms, to propose new choices for addresses indicated by an `Target` in the `target` field.
""", RegenerateContext)

@doc(
"""
```julia
ret, cl = regenerate(sel::L, cs::DynamicCallSite, new_args...) where L <: Target
ret, cl = regenerate(sel::L, cs::DynamicCallSite) where L <: Target
```
`regenerate` is an API to the `RegenerateContext` execution context. `regenerate` requires that users provide an `Target`, an original call site, and possibly a set of new arguments to be used in the regeneration step. This context internally keeps track of the bookkeeping required to increment likelihood weights, as well as prune off parts of the trace which are invalid if a regenerated choice changes the shape of the trace (e.g. control flow), and returns a new return value `ret` as well as the modified call site `cl`.
""", regenerate)


text:
function tile(context::NamedTuple, embedding::NamedTuple, selections::NamedTuple, colors::NamedTuple, names::NamedTuple)
    request = context.request
	try
		response = HTTP.Response(200)
		HTTP.Messages.setheader(response,"Content-type" => "image/png")
		channel = match(r"(?<=channel=)[^&]*(?=&)",request.target)
		scale = match(r"(?<=scale=)[^&]*(?=&)",request.target)

		body = IOBuffer(response.body, write=true)
		t = @elapsed save( Stream(format"PNG",body), tile(request.target,embedding,selections,colors,names; channel=isnothing(channel) ? "" : channel.match, scale=isnothing(scale) ? 1.0 : parse(Float32,scale.match)) )

		@info """$(request.method) $(replace(request.target, r"\?.+" => "")) | $t seconds"""
		return response

	catch exception
		printstyled("ERROR $(replace(request.target, r"\?.+" => "")) | $exception\n",color=:red)
		return HTTP.Response(500)
	end
end

function tile(url::AbstractString, embedding::NamedTuple, selections::NamedTuple, colors::NamedTuple, names::NamedTuple; channel::AbstractString="", scale::Number=1.0)
    zxy = broadcast(parse, Int, match(r"(\d+)/(\d+)/(\d+)", url ).captures )
    return colorview(RGBA,tile(zxy...,embedding,selections,colors,names; channel=channel, scale=scale))
end

function tile( z::Int, x::Int, y::Int, embedding::NamedTuple, selections::NamedTuple, colors::NamedTuple, names::NamedTuple; channel::AbstractString="", padIndex::Int=16, scale::Number=1.0)
	(xmin,xmax), (ymin,ymax) = embedding.extrema
	imageSize = max(xmax-xmin,ymax-ymin)

    tileSize = imageSize / 2^z
    padding = tileSize/padIndex

    x = xmin + x*tileSize
    y = ymin + y*tileSize

	X,Y = embedding.array[1,:], embedding.array[2,:]
	selections = @. selections.rows & (x-padding < X < x+tileSize+padding) & (y-padding < Y < y+tileSize+padding)

    return rasterKernelCircle( scale*sqrt(z),
        rasterize( (256+2padIndex,256+2padIndex), embedding.array[:,selections],
			channel ∈ names.channels ? colors.channels.colors[:,colors.channels.rows[selections,channel]] : colors.labels.rows[:,selections],

            xlim=( x-padding, x+tileSize+padding),
            ylim=( y-padding, y+tileSize+padding)
		)
	)[:,padIndex+1:end-padIndex,padIndex+1:end-padIndex]
end

max_stars_repo_name: gszep/flow-explorer

max_issues_repo_name: gszep/flow-explorer

max_forks_repo_name: gszep/flow-explorer


text:
# UnitCommitment.jl: Optimization Package for Security-Constrained Unit Commitment
# Copyright (C) 2020, UChicago Argonne, LLC. All rights reserved.
# Released under the modified BSD license. See COPYING.md for more details.

"""
    function optimize!(model::JuMP.Model)::Nothing

Solve the given unit commitment model. Unlike JuMP.optimize!, this uses more
advanced methods to accelerate the solution process and to enforce transmission
and N-1 security constraints.
"""
function optimize!(model::JuMP.Model)::Nothing
    return UnitCommitment.optimize!(model, XavQiuWanThi2019.Method())
end


text:
# Autogenerated wrapper script for Birch_jll for x86_64-unknown-freebsd
export birch

using boost_jll
using LibYAML_jll
JLLWrappers.@generate_wrapper_header("Birch")
JLLWrappers.@declare_executable_product(birch)
function __init__()
    JLLWrappers.@generate_init_header(boost_jll, LibYAML_jll)
    JLLWrappers.@init_executable_product(
        birch,
        "bin/birch",
    )

    JLLWrappers.@generate_init_footer()
end  # __init__()


text:
using DelimitedFiles

function checkregionorder(model::Model, regions, file)
    regionaliases = Dict{AbstractString, Vector{AbstractString}}("EU" => [],
                                                                 "USA" => ["US"],
                                                                 "OECD" => ["OT"],
                                                                 "Africa" => ["AF"],
                                                                 "China" => ["CA"],
                                                                 "SEAsia" => ["IA"],
                                                                 "LatAmerica" => ["LA"],
                                                                 "USSR" => ["EE"])

    for ii in 1:length(regions)
        region_keys = Mimi.dim_keys(model.md, :region)
        if region_keys[ii] != regions[ii] && !in(regions[ii], regionaliases[region_keys[ii]])
            error("Region indices in $file do not match expectations: $(region_keys[ii]) <> $(regions[ii]).")
        end
    end
end

function checktimeorder(model::Model, times, file)
    for ii in 1:length(times)
        if Mimi.time_labels(model)[ii] != times[ii]
            error("Time indices in $file do not match expectations: $(Mimi.time_labels(model)[ii]) <> $(times[ii]).")
        end
    end
end

function readpagedata(model::Union{Model, Nothing}, filepath::AbstractString)
    # Handle relative paths
    if filepath[1] ∉ ['.', '/'] && !isfile(filepath)
        filepath = joinpath(@__DIR__, "..", "..", filepath)
    end

    content = readlines(filepath)

    firstline = chomp(content[1])
    if firstline == "# Index: region"
        data = readdlm(filepath, ',', header=true, comments=true)

        if model != nothing
            # Check that regions are in the right order
            checkregionorder(model, data[1][:, 1], basename(filepath))
        end

        return convert(Vector{Float64},vec(data[1][:, 2]))
    elseif firstline == "# Index: time"
        data = readdlm(filepath, ',', header=true, comments=true)

        if model != nothing
            # Check that the times are in the right order
            checktimeorder(model, data[1][:, 1], basename(filepath))
        end

        return convert(Vector{Float64}, vec(data[1][:, 2]))
    elseif firstline == "# Index: time, region"
        data = readdlm(filepath, ',', header=true, comments = true)

        if model != nothing
            # Check that both dimension match
            checktimeorder(model, data[1][:, 1], basename(filepath))
            checkregionorder(model, data[2][2:end], basename(filepath))
        end

        return convert(Array{Float64}, data[1][:, 2:end])
    else
        error("Unknown header in parameter file $filepath.")
    end
end

function load_parameters(model::Model)
    parameters = Dict{Any, Any}()

    parameter_directory = joinpath(dirname(@__FILE__), "..", "..", "data")
    for file in filter(q->splitext(q)[2]==".csv", readdir(parameter_directory))
        parametername = splitext(file)[1]
        filepath = joinpath(parameter_directory, file)

        parameters[parametername] = readpagedata(model, filepath)
    end
    return parameters
end

max_stars_repo_name: lrennels/PAGE-ICE

max_issues_repo_name: lrennels/PAGE-ICE

max_forks_repo_name: lrennels/PAGE-ICE


text:
module TestFunctionCompilation

using Test
using DataFramesMeta

@testset "function_compilation" begin
    # Lazy way of making sure all functions are pre-compiled.
    # @eval prevents julia from caching the intermediate anonymous functions.
    for _ in 1:2
        @eval begin
            df = DataFrame(a = [1], b = [2])

            testfun(x, y) = x .* y
            testdotfun(x, y) = x * y
            testnt(x) = (c = x,)

            @test @select(df, :c =  :a + :b) == DataFrame(c = [3])

            fasttime = @timed @select(df, :c =  :a + :b)
            slowtime = @timed select(df, [:a, :b] => ((a, b) -> a + b) => :c)
            (slowtime[2] > fasttime[2]) || @warn("Slow compilation")

            @test @select(df, :c =  begin :a + :b end) == DataFrame(c = [3])

            fasttime = @timed @select(df, :c =  begin :a + :b end)
            slowtime = @timed select(df, [:a, :b] => ((a, b) -> a + b) => :c)
            (slowtime[2] > fasttime[2]) || @warn("Slow compilation")

            @test @select(df, cols(:c) = :a + :b) == DataFrame(c = [3])

            fasttime = @timed @select(df, cols(:c) = :a + :b)
            slowtime = @timed select(df, [:a, :b] => ((a, b) -> a + b) => :c)
            (slowtime[2] > fasttime[2]) || @warn("Slow compilation")

            @test @select(df, :c =  :a .+ :b) == DataFrame(c = [3])

            fasttime = @timed @select(df, :c =  :a .+ :b)
            slowtime = @timed select(df, [:a, :b] => ((a, b) -> a .+ b) => :c)
            (slowtime[2] > fasttime[2]) || @warn("Slow compilation")

            @test @select(df, cols(:c) = :a .+ :b) == DataFrame(c = [3])

            fasttime = @timed @select(df, cols(:c) = :a .+ :b)
            slowtime = @timed select(df, [:a, :b] => ((a, b) -> a .+ b) => :c)
            (slowtime[2] > fasttime[2]) || @warn("Slow compilation")

            @test @select(df, :c =  cols(:a) + cols(:b)) == DataFrame(c = [3])

            fasttime = @timed @select(df, :c =  cols(:a) + cols(:b))
            slowtime = @timed select(df, [:a, :b] => ((a, b) -> a + b) => :c)
            (slowtime[2] > fasttime[2]) || @warn("Slow compilation")

            @test @select(df, cols(:c) = cols(:a) + cols(:b)) == DataFrame(c = [3])

            fasttime = @timed @select(df, cols(:c) = cols(:a) + cols(:b))
            slowtime = @timed select(df, [:a, :b] => ((a, b) -> a + b) => :c)
            (slowtime[2] > fasttime[2]) || @warn("Slow compilation")

            @test @select(df, :c =  cols(:a) .+ cols(:b)) == DataFrame(c = [3])

            fasttime = @timed @select(df, :c =  cols(:a) .+ cols(:b))
            slowtime = @timed select(df, [:a, :b] => ((a, b) -> a .+ b) => :c)
            (slowtime[2] > fasttime[2]) || @warn("Slow compilation")

            @test @select(df, cols(:c) = cols(:a) .+ cols(:b)) == DataFrame(c = [3])

            fasttime = @timed @select(df, cols(:c) = cols(:a) .+ cols(:b))
            slowtime = @timed select(df, [:a, :b] => ((a, b) -> a .+ b) => :c)
            (slowtime[2] > fasttime[2]) || @warn("Slow compilation")

            @test @select(df, :c =  :a) == DataFrame(c = [1])

            fasttime = @timed @select(df, :c =  :a)
            slowtime = @timed select(df, [:a] => (a -> identity(a)) => :c)
            (slowtime[2] > fasttime[2]) || @warn("Slow compilation")

            @test @select(df, cols(:c) = :a) == DataFrame(c = [1])

            fasttime = @timed @select(df, cols(:c) = :a)
            slowtime = @timed select(df, [:a] => (a -> identity(a)) => :c)
            (slowtime[2] > fasttime[2]) || @warn("Slow compilation")

            @test @select(df, :c =  cols(:a)) == DataFrame(c = [1])

            fasttime = @timed @select(df, :c =  cols(:a))
            slowtime = @timed select(df, [:a] => (a -> identity(a)) => :c)
            (slowtime[2] > fasttime[2]) || @warn("Slow compilation")

            @test @select(df, cols(:c) = cols(:a)) == DataFrame(c = [1])

            fasttime = @timed @select(df, cols(:c) = cols(:a))
            slowtime = @timed select(df, [:a] => (a -> identity(a)) => :c)
            (slowtime[2] > fasttime[2]) || @warn("Slow compilation")

            @test @select(df, :a) == df[:, [:a]]

            fasttime = @timed @select(df, :a)
            slowtime = @timed select(df, [:a] => (a -> identity(a)) => :a)
            (slowtime[2] > fasttime[2]) || @warn("Slow compilation")

            @test @select(df, cols(:a)) == df[:, [:a]]

            fasttime = @timed @select(df, cols(:a))
            slowtime = @timed select(df, [:a] => (a -> identity(a)) => :a)
            (slowtime[2] > fasttime[2]) || @warn("Slow compilation")

            @test @select(df, :c =  testfun(:a, :b)) == DataFrame(c = [2])

            fasttime = @timed @select(df, :c =  testfun(:a, :b))
            slowtime = @timed select(df, [:a, :b] => ((a, b) -> testfun(a, b)) => :c)
            (slowtime[2] > fasttime[2]) || @warn("Slow compilation")

            @test @select(df, cols(:c) = testfun(:a, :b)) == DataFrame(c = [2])

            fasttime = @timed @select(df, cols(:c) = testfun(:a, :b))
            slowtime = @timed select(df, [:a, :b] => ((a, b) -> testfun(a, b)) => :c)
            (slowtime[2] > fasttime[2]) || @warn("Slow compilation")

            @test @select(df, :c =  testfun(cols("a"), cols("b"))) == DataFrame(c = [2])

            fasttime = @timed @select(df, :c =  testfun(cols("a"), cols("b")))
            slowtime = @timed select(df, [:a, :b] => ((a, b) -> testfun(a, b)) => :c)
            (slowtime[2] > fasttime[2]) || @warn("Slow compilation")

            @test @select(df, cols(:c) = testfun(cols("a"), cols("b"))) == DataFrame(c = [2])

            fasttime = @timed @select(df, cols(:c) = testfun(cols("a"), cols("b")))
            slowtime = @timed select(df, [:a, :b] => ((a, b) -> testfun(a, b)) => :c)
            (slowtime[2] > fasttime[2]) || @warn("Slow compilation")

            @test @select(df, :c =  testdotfun.(:a, :b)) == DataFrame(c = [2])

            fasttime = @timed @select(df, :c =  testdotfun.(:a, :b))
            slowtime = @timed select(df, [:a, :b] => ((a, b) -> testdotfun.(a, b)) => :c)
            (slowtime[2] > fasttime[2]) || @warn("Slow compilation")

            @test @select(df, cols(:c) = testdotfun.(:a, :b)) == DataFrame(c = [2])

            fasttime = @timed @select(df, cols(:c) = testdotfun.(:a, :b))
            slowtime = @timed select(df, [:a, :b] => ((a, b) -> testdotfun.(a, b)) => :c)
            (slowtime[2] > fasttime[2]) || @warn("Slow compilation")

            @test @select(df, :c =  testdotfun.(cols("a"), cols("b"))) == DataFrame(c = [2])

            fasttime = @timed @select(df, :c =  testdotfun.(cols("a"), cols("b")))
            slowtime = @timed select(df, [:a, :b] => ((a, b) -> testdotfun.(a, b)) => :c)
            (slowtime[2] > fasttime[2]) || @warn("Slow compilation")

            @test @select(df, cols(:c) = testdotfun.(cols("a"), cols("b"))) == DataFrame(c = [2])

            fasttime = @timed @select(df, cols(:c) = testdotfun.(cols("a"), cols("b")))
            slowtime = @timed select(df, [:a, :b] => ((a, b) -> testdotfun.(a, b)) => :c)
            (slowtime[2] > fasttime[2]) || @warn("Slow compilation")

            gd = groupby(df, :a)

            @test @combine(gd, testnt(:b)) == DataFrame(a = [1], c =  [2])

            fasttime = @timed @combine(gd, testnt(:b))
            slowtime = @timed combine(gd, :b => (b -> testnt(b)) => AsTable)
            (slowtime[2] > fasttime[2]) || @warn("Slow compilation")

            @test @with df (:a + :b) == [3]

            fasttime = @timed @with df (:a + :b)
            slowtime = @timed @with df ((a, b) -> a + b)(df.a, df.b)
            (slowtime[2] > fasttime[2]) || @warn("Slow compilation")

            @test @with df (:a .* :b) == [2]

            @with df (:a .* :b)
            fasttime = @timed @with df (:a .* :b)
            slowtime = @timed @with df ((a, b) -> a .* b)(df.a, df.b)
            (slowtime[2] > fasttime[2]) || @warn("Slow compilation")

            @test @with df testfun(:a, :b) == [2]

            fasttime = @timed @with df testfun(:a, :b)
            slowtime = @timed @with df ((a, b) -> testfun(a, b))(df.a, df.b)
            (slowtime[2] > fasttime[2]) || @warn("Slow compilation")

            @test @with df testdotfun.(:a, :b) == [2]

            fasttime = @timed @with df testdotfun.(:a, :b)
            slowtime = @timed @with df ((a, b) -> testdotfun.(a, b))(df.a, df.b)
            (slowtime[2] > fasttime[2]) || @warn("Slow compilation")
        end
    end
end
end # module


text:
# ---
# title: 509. Fibonacci Number
# id: problem509
# author: Tian Jun
# date: 2020-10-31
# difficulty: Easy
# categories: Array
# link: <https://leetcode.com/problems/fibonacci-number/description/>
# hidden: true
# ---
#
# The  **Fibonacci numbers** , commonly denoted `F(n)` form a sequence, called
# the  **Fibonacci sequence** , such that each number is the sum of the two
# preceding ones, starting from `0` and `1`. That is,
#
#
#
#     F(0) = 0,   F(1) = 1
#     F(N) = F(N - 1) + F(N - 2), for N > 1.
#
#
# Given `N`, calculate `F(N)`.
#
#
#
# **Example 1:**
#
#
#
#     Input: 2
#     Output: 1
#     Explanation: F(2) = F(1) + F(0) = 1 + 0 = 1.
#
#
# **Example 2:**
#
#
#
#     Input: 3
#     Output: 2
#     Explanation: F(3) = F(2) + F(1) = 1 + 1 = 2.
#
#
# **Example 3:**
#
#
#
#     Input: 4
#     Output: 3
#     Explanation: F(4) = F(3) + F(2) = 2 + 1 = 3.
#
#
#
#
# **Note:**
#
# 0 ≤ `N` ≤ 30.
#
#
## @lc code=start
using LeetCode

## add your code here:
## @lc code=end


text:
# # Wiring Diagrams as Attributed C-Sets
# Catlab supports many different flavors of diagrammatic syntax. These support the different combinatorial data structures that we use for representing categorical constructions. We will discuss DirectedWiringDiagrams, UndirectedWiringDiagrams, and CPortGraphs in this document.
using Catlab, Catlab.Theories
using Catlab.CategoricalAlgebra
using Catlab.WiringDiagrams
using Catlab.Graphics
using Catlab.Graphics.Graphviz
using Catlab.Programs
using Catlab.WiringDiagrams

draw(d::WiringDiagram) = to_graphviz(d,
  orientation=LeftToRight,
  labels=true, label_attr=:xlabel,
  node_attrs=Graphviz.Attributes(
    :fontname => "Courier",
  ),
  edge_attrs=Graphviz.Attributes(
    :fontname => "Courier",
  )
)

draw(uwd::AbstractUWD) = to_graphviz(uwd, junction_labels=:variable, box_labels=:name)

# ## Directed Wiring Diagrams

# DWDs are used to represent the morphisms in a symmetric monoidal category. You can get started by presenting a `FreeSymmetricMonoidalCategory` with the `@present` macro.
@present P(FreeSymmetricMonoidalCategory) begin
  (A,B,C,D)::Ob
  f::Hom(A,B)
  g::Hom(B,A)
  h::Hom(A⊗B, C)
  k::Hom(C,D⊗A)
end
generators(P)

# These presentations are very syntactic objects and expose an API for manipulating expressions.
for g in generators(P)
  "$g is a $(gat_typeof(g)) with arguments $(gat_type_args(g))"
end

# The `gat_typeof` function computes the algebraic type of a term by analogy to `Base.typeof` which computes the Julia type of a value.

homs_P = filter(generators(P)) do g
  gat_typeof(g) == :Hom
end

# When the term is a Hom, you can get the domain and codomain of the morphism with the `dom` and `codom` functions.
map(homs_P) do f
  "$f: $(dom(f)) → $(codom(f))"
end

# With the presentation you can build up morphism expressions using a formula syntax.
P[:h]⋅P[:k]

# This syntactic API is useful for manipulating terms in an arbitrary GAT and is the formal language of Catlab for representing and manipulating algebraic structures.
# However, when we want to work with big expressions in an SMC, the tree structure inherent to formulas is too verbose, and we want to move to a port-graph structure called `DirectedWiringDiagrams`. This gives us the benefits of combinatorial data structures like graphs with the right expressional power for representing the morphisms in an SMC.
wd = @program P (a::A, b::B) begin
  c = h(a,b)
  return k(c)
end

draw(wd)

# Catlab gives you the tools for drawing wiring diagrams. Visualization of wiring diagrams is the oldest part of Catlab and the original motivation for its development. The `@program` macro allows you to define wiring diagrams using a syntax that feels like Julia code.
wd = @program P (a::A, b::B) begin
  c = h(a,b)
  d,a₁ = k(c)
  return d, f(a₁)
end
draw(wd)

# The input wires are declared as *arguments* to the program, and the output wires are declared as *returns* from the function. Variables that are not consumed or by another function or returned by the program are automatically dropped.
wd = @program P (a::A, b::B) begin
  c = h(a,b)
  d,a₁ = k(c)
  return f(a₁)
end
draw(wd)

# You can copy a value by using it more than once. This is visualized as a wire being split into two wires.
wd = @program P (b::B) begin
  h(g(b),b)
end
draw(wd)

# You can visualize the copy and delete morphisms explicitly with the `add_junctions` function. The dots with one wire input and multiple outputs are copying values and dots with no wires out are deletions (discarding values). Not all instances of a `SymmetricMonoidalCategory` support copy and delete, for example, in manufacturing you can't duplicate a resource, and in chemistry you can't discard species. Catlab would enforce that when you tried to interpret the wiring diagram in a specific SMC.
wd = @program P (b::B) begin
  c = h(g(b),b)
  d, a₁ = k(c)
  return f(a₁)
end
draw(add_junctions(wd))

# For more details about working with wiring diagrams in Catlab, you should look at the other vignettes in this section which explain how wiring diagrams interact with SMC expressions and the basics of constructing and manipulation wiring diagrams.

# ## Diagrams as C-Sets

# The underlying data of a wiring diagram is combinatorial. That means we can represent it as a C-Set
wd.diagram

# Ok, there is a lot in there. The columns with integer entries are the combinatorial data encoding the connectivity of the wiring diagram. The columns with Symbols in them are encoding the labels for the diagram, the `value` of a box is the content of the box. Ports have types and wires have values. When we define the wiring diagram with the `@program` macro, we get a diagram that has labels and types, but no values. These values are initialized to nothing, but could be filled with values to be carried down the wires in an application.

# The schema of for wiring diagrams is called TheoryAttributedWiringDiagrams and is a little overwhelming, so we can explore how to build it up with C-Set schema inheritance.
to_graphviz(WiringDiagrams.DirectedWiringDiagrams.TheoryAttributedWiringDiagram)

# From the file Catlab/src/WiringDiagrams/Directed.jl
# ```julia
# @present TheoryWiringDiagram(FreeSchema) begin
#   Box::Ob
#   (InPort, OutPort, OuterInPort, OuterOutPort)::Ob
#   (Wire, InWire, OutWire, PassWire)::Ob

#   src::Hom(Wire, OutPort)
#   tgt::Hom(Wire, InPort)
#   in_src::Hom(InWire, OuterInPort)
#   in_tgt::Hom(InWire, InPort)
#   out_src::Hom(OutWire, OutPort)
#   out_tgt::Hom(OutWire, OuterOutPort)
#   pass_src::Hom(PassWire, OuterInPort)
#   pass_tgt::Hom(PassWire, OuterOutPort)

#   in_port_box::Hom(InPort, Box)
#   out_port_box::Hom(OutPort, Box)
# end

# @abstract_acset_type AbstractWiringDiagram <: AbstractGraph

# @present TheoryTypedWiringDiagram <: TheoryWiringDiagram begin
#   PortValue::AttrType
#   in_port_type::Attr(InPort, PortValue)
#   out_port_type::Attr(OutPort, PortValue)
#   outer_in_port_type::Attr(OuterInPort, PortValue)
#   outer_out_port_type::Attr(OuterOutPort, PortValue)
# end

# @present TheoryAttributedWiringDiagram <: TheoryTypedWiringDiagram begin
#   WireValue::AttrType
#   BoxValue::AttrType
#   BoxType::AttrType

#   value::Attr(Box, BoxValue)
#   box_type::Attr(Box, BoxType)
#   wire_value::Attr(Wire, WireValue)
#   in_wire_value::Attr(InWire, WireValue)
#   out_wire_value::Attr(OutWire, WireValue)
#   pass_wire_value::Attr(PassWire, WireValue)
# end
# ```

# The bare minimum diagram language is:
to_graphviz(WiringDiagrams.DirectedWiringDiagrams.TheoryWiringDiagram)

# And then you can add back in the types.
to_graphviz(WiringDiagrams.DirectedWiringDiagrams.TheoryTypedWiringDiagram)

#Layout is hard, so if you want to understand the `TheoryAttributedWiringDiagrams`, you should do the layout by hand as an exercise.

# We can create our own version of the theory of DWDs to see how it works:
@present MyTheoryWiringDiagram(FreeSchema) begin
  Box::Ob
  (InPort, OutPort, OuterInPort, OuterOutPort)::Ob
  (Wire, InWire, OutWire, PassWire)::Ob

  src::Hom(Wire, OutPort)
  tgt::Hom(Wire, InPort)
  in_src::Hom(InWire, OuterInPort)
  in_tgt::Hom(InWire, InPort)
  out_src::Hom(OutWire, OutPort)
  out_tgt::Hom(OutWire, OuterOutPort)
  pass_src::Hom(PassWire, OuterInPort)
  pass_tgt::Hom(PassWire, OuterOutPort)

  in_port_box::Hom(InPort, Box)
  out_port_box::Hom(OutPort, Box)
end

# If your application of wiring diagrams needs to attach numeric or textual information to the boxes of a wiring diagram, you would extend the `TheoryWiringDiagram` with the attributes that you need. That will give you a custom data structure that has those fields. One of the goals of Catlab is to make it so much easier to generate custom data structures that interoperate, that you don't need to create generic structures that can be used for many purposes. Just snap your fingers and create a new structure perfectly tailored to your needs, when your needs change, snap again to get a new version of that structure.

# The `@acset_type` macro does the hard work of generating the data structure and accessors and mutators for you. The form of this call is `@acset_type NewStructName(Schema, index=[morphisms in Schema]) <: Supertype`. You should index any morphism where you need to use `incident` frequently. For wiring diagrams you will often want to know what all the wires that are incident to a port.
@acset_type MyWiringDiagramACSet(MyTheoryWiringDiagram,
  index=[:src, :tgt, :in_src, :in_tgt, :out_src, :out_tgt, :pass_src, :pass_tgt]) <: WiringDiagrams.DirectedWiringDiagrams.AbstractWiringDiagram

# We get the `@acset` macro from Catlab and can create DWDs by hand. It is very tedious, which is why the `@program` macro exists!
md = @acset MyWiringDiagramACSet begin
  Box = 3
  InPort = 6
  OutPort = 3
  Wire = 3
  src = [1,2,3]
  tgt = [3,4,5]
  in_port_box = [1,1,2,2,3,3]
  out_port_box = [1,2,3]
end

# ## Undirected Wiring Diagrams
# A much simpler structure than DWDs are known as undirected wiring diagrams. They are called undirected because ports boxes have one set of ports that aren't divided into inputs and outputs, and the wires are undirected. Wires connect junctions to ports (which live on boxes).
to_graphviz(WiringDiagrams.UndirectedWiringDiagrams.TheoryUWD)

# These UWDs are combinatorial syntax for relations. The junctions are variables and the boxes are the relations. A relation R ⊆ X × Y has two ports one for the value of X and one for the value of Y. The expression R(x:X, y:Y) says to connect the X port of R to the junction for the variable x, and the Y port of R to the y variable junction. If two ports are attached to the same junction, then you have a constraint that those values must be equal. The outer ports are the components of the final relation. For example the following UWD encodes the relation {(x,y,z) | R(x,y) and S(y,z) for all x∈X, y∈Y, z∈Z}.
uwd = @relation (x, y, z) begin
  R(x,y)
  S(y,z)
end

# These UWDs are drawn with circular boxes and undirected wires. Note that since all wires go from junction to port, they are not *symmetric* wires.
draw(uwd)

# By adding more relations we can get bigger relations.
uwd₂ = @relation (x, z) begin
  R(x,y)
  S(y,z)
  T(x,z,y)
end

# Not all of the junctions have to be exposed to the outside world. Note that there is no distinction between arguments and return values in the relation macro. This is because relations are inherently undirected, unlike functions.
draw(uwd₂)

# ## Circular Port Graphs
# CPGs are the natural data structure for representing process of interconnected systems that share information along wires, but send different information to their different neighbors.
to_graphviz(ThCPortGraph)

# They are also a kind of CSet, so we can use the `@acset` macro to construct them.
cpg = @acset CPortGraph begin
  Box = 3
  Port = 6
  Wire = 5

  box = [1,1,2,2,3,3]
  src = [1,2,5,2,6]
  tgt = [3,4,2,4,1]
end

# the layout for CPGs is not great with graphviz.
to_graphviz(cpg, port_labels=true, graph_attrs=Dict("nodesep"=>"2", "ranksep"=>"2"))

# Almost every application of graphs in computer science could be better served by using one of these extensions to the basic graph data structure.


text:
# This file is auto-generated by AWSMetadata.jl
using AWS
using AWS.AWSServices: shield

using AWS.Compat
using AWS.UUIDs
"""
    AssociateDRTLogBucket()

Authorizes the DDoS Response Team (DRT) to access the specified Amazon S3 bucket containing your AWS WAF logs. You can associate up to 10 Amazon S3 buckets with your subscription. To use the services of the DRT and make an AssociateDRTLogBucket request, you must be subscribed to the Business Support plan or the Enterprise Support plan.

# Required Parameters
- `LogBucket`: The Amazon S3 bucket that contains your AWS WAF logs.

"""

associate_drtlog_bucket(LogBucket; aws_config::AWSConfig=global_aws_config()) = shield("AssociateDRTLogBucket", Dict{String, Any}("LogBucket"=>LogBucket); aws_config=aws_config)
associate_drtlog_bucket(LogBucket, args::AbstractDict{String, <:Any}; aws_config::AWSConfig=global_aws_config()) = shield("AssociateDRTLogBucket", Dict{String, Any}(mergewith(_merge, Dict{String, Any}("LogBucket"=>LogBucket), args)); aws_config=aws_config)

"""
    AssociateDRTRole()

Authorizes the DDoS Response Team (DRT), using the specified role, to access your AWS account to assist with DDoS attack mitigation during potential attacks. This enables the DRT to inspect your AWS WAF configuration and create or update AWS WAF rules and web ACLs. You can associate only one RoleArn with your subscription. If you submit an AssociateDRTRole request for an account that already has an associated role, the new RoleArn will replace the existing RoleArn.  Prior to making the AssociateDRTRole request, you must attach the AWSShieldDRTAccessPolicy managed policy to the role you will specify in the request. For more information see Attaching and Detaching IAM Policies. The role must also trust the service principal  drt.shield.amazonaws.com. For more information, see IAM JSON Policy Elements: Principal. The DRT will have access only to your AWS WAF and Shield resources. By submitting this request, you authorize the DRT to inspect your AWS WAF and Shield configuration and create and update AWS WAF rules and web ACLs on your behalf. The DRT takes these actions only if explicitly authorized by you. You must have the iam:PassRole permission to make an AssociateDRTRole request. For more information, see Granting a User Permissions to Pass a Role to an AWS Service.  To use the services of the DRT and make an AssociateDRTRole request, you must be subscribed to the Business Support plan or the Enterprise Support plan.

# Required Parameters
- `RoleArn`: The Amazon Resource Name (ARN) of the role the DRT will use to access your AWS account. Prior to making the AssociateDRTRole request, you must attach the AWSShieldDRTAccessPolicy managed policy to this role. For more information see Attaching and Detaching IAM Policies.

"""

associate_drtrole(RoleArn; aws_config::AWSConfig=global_aws_config()) = shield("AssociateDRTRole", Dict{String, Any}("RoleArn"=>RoleArn); aws_config=aws_config)
associate_drtrole(RoleArn, args::AbstractDict{String, <:Any}; aws_config::AWSConfig=global_aws_config()) = shield("AssociateDRTRole", Dict{String, Any}(mergewith(_merge, Dict{String, Any}("RoleArn"=>RoleArn), args)); aws_config=aws_config)

"""
    AssociateHealthCheck()

Adds health-based detection to the Shield Advanced protection for a resource. Shield Advanced health-based detection uses the health of your AWS resource to improve responsiveness and accuracy in attack detection and mitigation.  You define the health check in Route 53 and then associate it with your Shield Advanced protection. For more information, see Shield Advanced Health-Based Detection in the AWS WAF and AWS Shield Developer Guide.

# Required Parameters
- `HealthCheckArn`: The Amazon Resource Name (ARN) of the health check to associate with the protection.
- `ProtectionId`: The unique identifier (ID) for the Protection object to add the health check association to.

"""

associate_health_check(HealthCheckArn, ProtectionId; aws_config::AWSConfig=global_aws_config()) = shield("AssociateHealthCheck", Dict{String, Any}("HealthCheckArn"=>HealthCheckArn, "ProtectionId"=>ProtectionId); aws_config=aws_config)
associate_health_check(HealthCheckArn, ProtectionId, args::AbstractDict{String, <:Any}; aws_config::AWSConfig=global_aws_config()) = shield("AssociateHealthCheck", Dict{String, Any}(mergewith(_merge, Dict{String, Any}("HealthCheckArn"=>HealthCheckArn, "ProtectionId"=>ProtectionId), args)); aws_config=aws_config)

"""
    AssociateProactiveEngagementDetails()

Initializes proactive engagement and sets the list of contacts for the DDoS Response Team (DRT) to use. You must provide at least one phone number in the emergency contact list.  After you have initialized proactive engagement using this call, to disable or enable proactive engagement, use the calls DisableProactiveEngagement and EnableProactiveEngagement.   This call defines the list of email addresses and phone numbers that the DDoS Response Team (DRT) can use to contact you for escalations to the DRT and to initiate proactive customer support. The contacts that you provide in the request replace any contacts that were already defined. If you already have contacts defined and want to use them, retrieve the list using DescribeEmergencyContactSettings and then provide it to this call.

# Required Parameters
- `EmergencyContactList`: A list of email addresses and phone numbers that the DDoS Response Team (DRT) can use to contact you for escalations to the DRT and to initiate proactive customer support.  To enable proactive engagement, the contact list must include at least one phone number.  The contacts that you provide here replace any contacts that were already defined. If you already have contacts defined and want to use them, retrieve the list using DescribeEmergencyContactSettings and then provide it here.

"""

associate_proactive_engagement_details(EmergencyContactList; aws_config::AWSConfig=global_aws_config()) = shield("AssociateProactiveEngagementDetails", Dict{String, Any}("EmergencyContactList"=>EmergencyContactList); aws_config=aws_config)
associate_proactive_engagement_details(EmergencyContactList, args::AbstractDict{String, <:Any}; aws_config::AWSConfig=global_aws_config()) = shield("AssociateProactiveEngagementDetails", Dict{String, Any}(mergewith(_merge, Dict{String, Any}("EmergencyContactList"=>EmergencyContactList), args)); aws_config=aws_config)

"""
    CreateProtection()

Enables AWS Shield Advanced for a specific AWS resource. The resource can be an Amazon CloudFront distribution, Elastic Load Balancing load balancer, AWS Global Accelerator accelerator, Elastic IP Address, or an Amazon Route 53 hosted zone. You can add protection to only a single resource with each CreateProtection request. If you want to add protection to multiple resources at once, use the AWS WAF console. For more information see Getting Started with AWS Shield Advanced and Add AWS Shield Advanced Protection to more AWS Resources.

# Required Parameters
- `Name`: Friendly name for the Protection you are creating.
- `ResourceArn`: The ARN (Amazon Resource Name) of the resource to be protected. The ARN should be in one of the following formats:   For an Application Load Balancer: arn:aws:elasticloadbalancing:region:account-id:loadbalancer/app/load-balancer-name/load-balancer-id     For an Elastic Load Balancer (Classic Load Balancer): arn:aws:elasticloadbalancing:region:account-id:loadbalancer/load-balancer-name     For an AWS CloudFront distribution: arn:aws:cloudfront::account-id:distribution/distribution-id     For an AWS Global Accelerator accelerator: arn:aws:globalaccelerator::account-id:accelerator/accelerator-id     For Amazon Route 53: arn:aws:route53:::hostedzone/hosted-zone-id     For an Elastic IP address: arn:aws:ec2:region:account-id:eip-allocation/allocation-id

"""

create_protection(Name, ResourceArn; aws_config::AWSConfig=global_aws_config()) = shield("CreateProtection", Dict{String, Any}("Name"=>Name, "ResourceArn"=>ResourceArn); aws_config=aws_config)
create_protection(Name, ResourceArn, args::AbstractDict{String, <:Any}; aws_config::AWSConfig=global_aws_config()) = shield("CreateProtection", Dict{String, Any}(mergewith(_merge, Dict{String, Any}("Name"=>Name, "ResourceArn"=>ResourceArn), args)); aws_config=aws_config)

"""
    CreateSubscription()

Activates AWS Shield Advanced for an account. When you initally create a subscription, your subscription is set to be automatically renewed at the end of the existing subscription period. You can change this by submitting an UpdateSubscription request.

"""

create_subscription(; aws_config::AWSConfig=global_aws_config()) = shield("CreateSubscription"; aws_config=aws_config)
create_subscription(args::AbstractDict{String, <:Any}; aws_config::AWSConfig=global_aws_config()) = shield("CreateSubscription", args; aws_config=aws_config)

"""
    DeleteProtection()

Deletes an AWS Shield Advanced Protection.

# Required Parameters
- `ProtectionId`: The unique identifier (ID) for the Protection object to be deleted.

"""

delete_protection(ProtectionId; aws_config::AWSConfig=global_aws_config()) = shield("DeleteProtection", Dict{String, Any}("ProtectionId"=>ProtectionId); aws_config=aws_config)
delete_protection(ProtectionId, args::AbstractDict{String, <:Any}; aws_config::AWSConfig=global_aws_config()) = shield("DeleteProtection", Dict{String, Any}(mergewith(_merge, Dict{String, Any}("ProtectionId"=>ProtectionId), args)); aws_config=aws_config)

"""
    DeleteSubscription()

Removes AWS Shield Advanced from an account. AWS Shield Advanced requires a 1-year subscription commitment. You cannot delete a subscription prior to the completion of that commitment.

"""

delete_subscription(; aws_config::AWSConfig=global_aws_config()) = shield("DeleteSubscription"; aws_config=aws_config)
delete_subscription(args::AbstractDict{String, <:Any}; aws_config::AWSConfig=global_aws_config()) = shield("DeleteSubscription", args; aws_config=aws_config)

"""
    DescribeAttack()

Describes the details of a DDoS attack.

# Required Parameters
- `AttackId`: The unique identifier (ID) for the attack that to be described.

"""

describe_attack(AttackId; aws_config::AWSConfig=global_aws_config()) = shield("DescribeAttack", Dict{String, Any}("AttackId"=>AttackId); aws_config=aws_config)
describe_attack(AttackId, args::AbstractDict{String, <:Any}; aws_config::AWSConfig=global_aws_config()) = shield("DescribeAttack", Dict{String, Any}(mergewith(_merge, Dict{String, Any}("AttackId"=>AttackId), args)); aws_config=aws_config)

"""
    DescribeDRTAccess()

Returns the current role and list of Amazon S3 log buckets used by the DDoS Response Team (DRT) to access your AWS account while assisting with attack mitigation.

"""

describe_drtaccess(; aws_config::AWSConfig=global_aws_config()) = shield("DescribeDRTAccess"; aws_config=aws_config)
describe_drtaccess(args::AbstractDict{String, <:Any}; aws_config::AWSConfig=global_aws_config()) = shield("DescribeDRTAccess", args; aws_config=aws_config)

"""
    DescribeEmergencyContactSettings()

A list of email addresses and phone numbers that the DDoS Response Team (DRT) can use to contact you if you have proactive engagement enabled, for escalations to the DRT and to initiate proactive customer support.

"""

describe_emergency_contact_settings(; aws_config::AWSConfig=global_aws_config()) = shield("DescribeEmergencyContactSettings"; aws_config=aws_config)
describe_emergency_contact_settings(args::AbstractDict{String, <:Any}; aws_config::AWSConfig=global_aws_config()) = shield("DescribeEmergencyContactSettings", args; aws_config=aws_config)

"""
    DescribeProtection()

Lists the details of a Protection object.

# Optional Parameters
- `ProtectionId`: The unique identifier (ID) for the Protection object that is described. When submitting the DescribeProtection request you must provide either the ResourceArn or the ProtectionID, but not both.
- `ResourceArn`: The ARN (Amazon Resource Name) of the AWS resource for the Protection object that is described. When submitting the DescribeProtection request you must provide either the ResourceArn or the ProtectionID, but not both.
"""

describe_protection(; aws_config::AWSConfig=global_aws_config()) = shield("DescribeProtection"; aws_config=aws_config)
describe_protection(args::AbstractDict{String, <:Any}; aws_config::AWSConfig=global_aws_config()) = shield("DescribeProtection", args; aws_config=aws_config)

"""
    DescribeSubscription()

Provides details about the AWS Shield Advanced subscription for an account.

"""

describe_subscription(; aws_config::AWSConfig=global_aws_config()) = shield("DescribeSubscription"; aws_config=aws_config)
describe_subscription(args::AbstractDict{String, <:Any}; aws_config::AWSConfig=global_aws_config()) = shield("DescribeSubscription", args; aws_config=aws_config)

"""
    DisableProactiveEngagement()

Removes authorization from the DDoS Response Team (DRT) to notify contacts about escalations to the DRT and to initiate proactive customer support.

"""

disable_proactive_engagement(; aws_config::AWSConfig=global_aws_config()) = shield("DisableProactiveEngagement"; aws_config=aws_config)
disable_proactive_engagement(args::AbstractDict{String, <:Any}; aws_config::AWSConfig=global_aws_config()) = shield("DisableProactiveEngagement", args; aws_config=aws_config)

"""
    DisassociateDRTLogBucket()

Removes the DDoS Response Team's (DRT) access to the specified Amazon S3 bucket containing your AWS WAF logs. To make a DisassociateDRTLogBucket request, you must be subscribed to the Business Support plan or the Enterprise Support plan. However, if you are not subscribed to one of these support plans, but had been previously and had granted the DRT access to your account, you can submit a DisassociateDRTLogBucket request to remove this access.

# Required Parameters
- `LogBucket`: The Amazon S3 bucket that contains your AWS WAF logs.

"""

disassociate_drtlog_bucket(LogBucket; aws_config::AWSConfig=global_aws_config()) = shield("DisassociateDRTLogBucket", Dict{String, Any}("LogBucket"=>LogBucket); aws_config=aws_config)
disassociate_drtlog_bucket(LogBucket, args::AbstractDict{String, <:Any}; aws_config::AWSConfig=global_aws_config()) = shield("DisassociateDRTLogBucket", Dict{String, Any}(mergewith(_merge, Dict{String, Any}("LogBucket"=>LogBucket), args)); aws_config=aws_config)

"""
    DisassociateDRTRole()

Removes the DDoS Response Team's (DRT) access to your AWS account. To make a DisassociateDRTRole request, you must be subscribed to the Business Support plan or the Enterprise Support plan. However, if you are not subscribed to one of these support plans, but had been previously and had granted the DRT access to your account, you can submit a DisassociateDRTRole request to remove this access.

"""

disassociate_drtrole(; aws_config::AWSConfig=global_aws_config()) = shield("DisassociateDRTRole"; aws_config=aws_config)
disassociate_drtrole(args::AbstractDict{String, <:Any}; aws_config::AWSConfig=global_aws_config()) = shield("DisassociateDRTRole", args; aws_config=aws_config)

"""
    DisassociateHealthCheck()

Removes health-based detection from the Shield Advanced protection for a resource. Shield Advanced health-based detection uses the health of your AWS resource to improve responsiveness and accuracy in attack detection and mitigation.  You define the health check in Route 53 and then associate or disassociate it with your Shield Advanced protection. For more information, see Shield Advanced Health-Based Detection in the AWS WAF and AWS Shield Developer Guide.

# Required Parameters
- `HealthCheckArn`: The Amazon Resource Name (ARN) of the health check that is associated with the protection.
- `ProtectionId`: The unique identifier (ID) for the Protection object to remove the health check association from.

"""

disassociate_health_check(HealthCheckArn, ProtectionId; aws_config::AWSConfig=global_aws_config()) = shield("DisassociateHealthCheck", Dict{String, Any}("HealthCheckArn"=>HealthCheckArn, "ProtectionId"=>ProtectionId); aws_config=aws_config)
disassociate_health_check(HealthCheckArn, ProtectionId, args::AbstractDict{String, <:Any}; aws_config::AWSConfig=global_aws_config()) = shield("DisassociateHealthCheck", Dict{String, Any}(mergewith(_merge, Dict{String, Any}("HealthCheckArn"=>HealthCheckArn, "ProtectionId"=>ProtectionId), args)); aws_config=aws_config)

"""
    EnableProactiveEngagement()

Authorizes the DDoS Response Team (DRT) to use email and phone to notify contacts about escalations to the DRT and to initiate proactive customer support.

"""

enable_proactive_engagement(; aws_config::AWSConfig=global_aws_config()) = shield("EnableProactiveEngagement"; aws_config=aws_config)
enable_proactive_engagement(args::AbstractDict{String, <:Any}; aws_config::AWSConfig=global_aws_config()) = shield("EnableProactiveEngagement", args; aws_config=aws_config)

"""
    GetSubscriptionState()

Returns the SubscriptionState, either Active or Inactive.

"""

get_subscription_state(; aws_config::AWSConfig=global_aws_config()) = shield("GetSubscriptionState"; aws_config=aws_config)
get_subscription_state(args::AbstractDict{String, <:Any}; aws_config::AWSConfig=global_aws_config()) = shield("GetSubscriptionState", args; aws_config=aws_config)

"""
    ListAttacks()

Returns all ongoing DDoS attacks or all DDoS attacks during a specified time period.

# Optional Parameters
- `EndTime`: The end of the time period for the attacks. This is a timestamp type. The sample request above indicates a number type because the default used by WAF is Unix time in seconds. However any valid timestamp format is allowed.
- `MaxResults`: The maximum number of AttackSummary objects to be returned. If this is left blank, the first 20 results will be returned. This is a maximum value; it is possible that AWS WAF will return the results in smaller batches. That is, the number of AttackSummary objects returned could be less than MaxResults, even if there are still more AttackSummary objects yet to return. If there are more AttackSummary objects to return, AWS WAF will always also return a NextToken.
- `NextToken`: The ListAttacksRequest.NextMarker value from a previous call to ListAttacksRequest. Pass null if this is the first call.
- `ResourceArns`: The ARN (Amazon Resource Name) of the resource that was attacked. If this is left blank, all applicable resources for this account will be included.
- `StartTime`: The start of the time period for the attacks. This is a timestamp type. The sample request above indicates a number type because the default used by WAF is Unix time in seconds. However any valid timestamp format is allowed.
"""

list_attacks(; aws_config::AWSConfig=global_aws_config()) = shield("ListAttacks"; aws_config=aws_config)
list_attacks(args::AbstractDict{String, <:Any}; aws_config::AWSConfig=global_aws_config()) = shield("ListAttacks", args; aws_config=aws_config)

"""
    ListProtections()

Lists all Protection objects for the account.

# Optional Parameters
- `MaxResults`: The maximum number of Protection objects to be returned. If this is left blank the first 20 results will be returned. This is a maximum value; it is possible that AWS WAF will return the results in smaller batches. That is, the number of Protection objects returned could be less than MaxResults, even if there are still more Protection objects yet to return. If there are more Protection objects to return, AWS WAF will always also return a NextToken.
- `NextToken`: The ListProtectionsRequest.NextToken value from a previous call to ListProtections. Pass null if this is the first call.
"""

list_protections(; aws_config::AWSConfig=global_aws_config()) = shield("ListProtections"; aws_config=aws_config)
list_protections(args::AbstractDict{String, <:Any}; aws_config::AWSConfig=global_aws_config()) = shield("ListProtections", args; aws_config=aws_config)

"""
    UpdateEmergencyContactSettings()

Updates the details of the list of email addresses and phone numbers that the DDoS Response Team (DRT) can use to contact you if you have proactive engagement enabled, for escalations to the DRT and to initiate proactive customer support.

# Optional Parameters
- `EmergencyContactList`: A list of email addresses and phone numbers that the DDoS Response Team (DRT) can use to contact you if you have proactive engagement enabled, for escalations to the DRT and to initiate proactive customer support. If you have proactive engagement enabled, the contact list must include at least one phone number.
"""

update_emergency_contact_settings(; aws_config::AWSConfig=global_aws_config()) = shield("UpdateEmergencyContactSettings"; aws_config=aws_config)
update_emergency_contact_settings(args::AbstractDict{String, <:Any}; aws_config::AWSConfig=global_aws_config()) = shield("UpdateEmergencyContactSettings", args; aws_config=aws_config)

"""
    UpdateSubscription()

Updates the details of an existing subscription. Only enter values for parameters you want to change. Empty parameters are not updated.

# Optional Parameters
- `AutoRenew`: When you initally create a subscription, AutoRenew is set to ENABLED. If ENABLED, the subscription will be automatically renewed at the end of the existing subscription period. You can change this by submitting an UpdateSubscription request. If the UpdateSubscription request does not included a value for AutoRenew, the existing value for AutoRenew remains unchanged.
"""

update_subscription(; aws_config::AWSConfig=global_aws_config()) = shield("UpdateSubscription"; aws_config=aws_config)
update_subscription(args::AbstractDict{String, <:Any}; aws_config::AWSConfig=global_aws_config()) = shield("UpdateSubscription", args; aws_config=aws_config)


text:
module ofMeshReader
#=
Functions to:
- Read polymesh in a OpenFOAM Case folder
- Finds and checks the polyMesh folder
- Parse files in polyMesh folder
- Generate OfMesh (ofmeshTypes.jl) from parsed entities

TODO:
- maybe do some file encoding checking, to avoid encoding errors?
- performance testing and optimization
- implement unpacking of  gzip files in RAM and procesing with ascii reading
- implement reading of binary dictionary files -> see also ofMeshDictTokenizer.jl
- tests against more differnt meshes
- is some additional masking of neighbour faces necessarry? or are the faces with no neighbour face_owner_cells allways the last faces in faces file, as asumbed now?
=#

include("ofMeshTypes.jl")
include("ofMeshDictParser.jl")

const OfMeshFileType = Dict(
    :ascii => 0,
    :binary => 1,
    :gz => 2
    )

export read_openfoam_mesh

"""
    read_openfoam_mesh()

read, check, parse openfoam mesh
"""
function read_openfoam_mesh(case_path::String)

    # toggle error if necessarry file in case_path/constant/polyMesh non
    # existent ...
    checkmesh(case_path)

    polymesh_path = joinpath(joinpath(case_path, "constant"), "polyMesh")

    # bf: boundary file
    bf = joinpath(polymesh_path,"boundary")
    # boundary file should be ascii every time ...
    bf_tokens = tokenize_ascii_ofmeshfile(bf)
    bf_preparse = parser(bf_tokens)
    bf_parsed = dfs_parse(bf_preparse, bf_tokens)
    meshinfo = getMeshInfo(bf_parsed, case_path)

    if meshinfo.format == "ascii"
        return parse_mesh(meshinfo, parse_ascii_file)
    else
        error("Binary and Gzip file formats not supported jet!")
    end

end

function checkmesh(case_path::String)

    polymesh_path = joinpath(joinpath(case_path, "constant"), "polyMesh")
    checkdir(polymesh_path)
    checkifmeshfilesexits(polymesh_path)
end

function checkdir(polymesh_path::String)

    if !ispath(polymesh_path) || !isdir(polymesh_path)
        error("Folder ", polymesh_path, " not found!")
    end
    return true
end

function checkifmeshfilesexits(polymesh_path::String)

    if !isfile(joinpath(polymesh_path, "boundary"))
        error("File boundary in constant/polyMesh is missing or is not a " *
              "valid file")
    end
    if !isfile(joinpath(polymesh_path, "faces"))
        error("File faces in constant/polyMesh is missing or is not a " *
              "valid file")
    end
    if !isfile(joinpath(polymesh_path, "neighbour"))
        error("File neighbour in constant/polyMesh is missing or is not a" *
              " valid file")
    end
    if !isfile(joinpath(polymesh_path, "owner"))
        error("File owner in constant/polyMesh is missing or is not a " *
              "valid file")
    end
    if !isfile(joinpath(polymesh_path, "points"))
        error("File points in constant/polyMesh is missing or is not a " *
              "valid file")
    end
    return true
end

"""
    getMeshInfo(bf_parsed::Vector{OfEntity}, case_path::String)

Extract some mesh information out of bounary file
"""
function getMeshInfo(bf_parsed::Vector{OfEntity}, case_path::String)
    bf = get_dictionary(bf_parsed, "FoamFile")

    if isnothing(bf)
        error("Cannot find file information in bounary file")
    end

    file_version = get_variable_value(bf.entries, "version")
    file_format = get_variable_value(bf.entries, "format")
    file_class = get_variable_value(bf.entries, "class")
    file_rel_location = get_variable_value(bf.entries, "location")
    file_case_path = case_path

    return OfMeshInfo(file_version, file_format, file_class,
                file_rel_location,file_case_path )

end

function parse_ascii_file(polymesh_path::String, dict_file::String)
    f = joinpath(polymesh_path,dict_file)
    # boundary file should be ascii every time ...
    f_tokens = tokenize_ascii_ofmeshfile(f)
    f_preparse = parser(f_tokens)
    f_parsed = dfs_parse(f_preparse, f_tokens)

    return f_parsed
end

function parse_binary_file()
    # TODO: Implement binary tokenizer / reader
end

function parse_mesh(meshinfo::OfMeshInfo, pasing_fct::Function)

    polymesh_path = joinpath(meshinfo.casepath, meshinfo.location)

    bf_raw = pasing_fct(polymesh_path, "boundary")
    patches = parse_boundary(bf_raw)

    f_raw = pasing_fct(polymesh_path, "faces")
    f = parse_face(f_raw)

    n_raw = pasing_fct(polymesh_path, "neighbour")
    n = parse_neighbour(n_raw)

    o_raw = pasing_fct(polymesh_path, "owner")
    o = parse_owner(o_raw)

    p_raw = pasing_fct(polymesh_path, "points")
    p = parse_point(p_raw::Vector{OfEntity})

    return OfMesh(meshinfo, p, f, patches, o, n)
end

"""
    parse_point(p_raw::Vector{OfEntity})

Parses openfoam point file (from finally parsed tokens) expects an OfArray containing OfArray.length OfLists of OfList.length = 3 with mesh Points in OfValue{Float64}.
"""
function parse_point(p_raw::Vector{OfEntity})

    if typeof(p_raw[2]) == OfArray

        n_points = p_raw[2].length
        x = zeros(Float64, n_points)
        y = zeros(Float64, n_points)
        z = zeros(Float64, n_points)

        for (i,p_list) in enumerate(p_raw[2].entries)
            x[i] = p_list.entries[1].value
            y[i] = p_list.entries[2].value
            z[i] = p_list.entries[3].value
        end
    else
        error("Reading points from points file")
    end

    return OfPoints(n_points,x,y,z)
end

function parse_face(f_raw::Vector{OfEntity})

    if typeof(f_raw[2]) == OfArray

        no_faces = f_raw[2].length
        face_list = Vector{OfFace}(undef,0)

        for (i,f_list) in enumerate(f_raw[2].entries)

            no_points = f_list.length
            point_List = zeros(Int64, no_points)

            for j = 1:length(f_list.entries)
                point_List[j] = f_list.entries[j].value
            end

            push!(face_list, OfFace(no_points,point_List))
        end

        return OfFaces(no_faces, face_list)
    else
        error("Reading faces from faces file")
    end
end

function parse_boundary(b_raw::Vector{OfEntity})

    if typeof(b_raw[2]) == OfArray

        no_patches = b_raw[2].length
        patches = Vector{OfPatch}(undef, 0)

        for boundary_dict in b_raw[2].entries

            type = get_variable_value(boundary_dict.entries, "type")
            start_face = get_variable_value(boundary_dict.entries, "startFace")
            no_faces = get_variable_value(boundary_dict.entries, "nFaces")

            if isnothing(type) || isnothing(start_face) || isnothing(no_faces)
                error("Reading boundary file from parsed tokens")
            end

            push!(patches, OfPatch(no_faces, start_face, type))
        end

        return patches
    else
        error("Expected Array in boundary file at pos 2 in parsed tokens")
    end

end

function parse_neighbour(n_raw::Vector{OfEntity})
# right now simply assuming that the boundarys begin at the higher numbers
# of the mesh

    if typeof(n_raw[2]) == OfArray

        n_internal_faces = n_raw[2].length
        face_neighbour_cells = zeros(Int64, n_internal_faces)

        for (i,n) in enumerate(n_raw[2].entries)
            face_neighbour_cells[i] = n.value
        end
    else
        error("Reading points from points file")
    end

    return face_neighbour_cells
end

function parse_owner(o_raw::Vector{OfEntity})
    if typeof(o_raw[2]) == OfArray

        n_faces = o_raw[2].length
        face_owner_cells = zeros(Int64, n_faces)

        for (i,o) in enumerate(o_raw[2].entries)
            face_owner_cells[i] = o.value
        end

    else
        error("Reading points from points file")
    end

    return face_owner_cells
end

end

max_stars_repo_name: c-schubert/linearMeshIO

max_issues_repo_name: c-schubert/linearMeshIO

max_forks_repo_name: c-schubert/linearMeshIO


text:
"""
    QuarticNonconvexPeriodicSemidiscretization(D, Di, split_form)

A semidiscretization of the quartic nonconvex conservation law
    \$\\partial_t u(t,x) + \\partial_x ( u(t,x)^4 - 10 u(t,x)^2 + 3 u(t,x) ) = 0\$
with periodic boundary conditions.

`D` is a first-derivative SBP operator, `Di` an associated dissipation operator
or `nothing`, and `split_form::Union{Val(true), Val(false)}` determines whether
the canonical split form or the conservative form is used.
"""
struct QuarticNonconvexPeriodicSemidiscretization{T,Derivative<:AbstractDerivativeOperator{T},
                                                  Dissipation,
                                                  SplitForm<:Union{Val{false}, Val{true}}} <: AbstractSemidiscretization
    derivative::Derivative
    dissipation::Dissipation
    tmp1::Vector{T}
    tmp2::Vector{T}
    split_form::SplitForm

    function QuarticNonconvexPeriodicSemidiscretization(derivative::Derivative, dissipation::Dissipation, split_form::SplitForm=Val{false}()) where {T, Derivative<:AbstractDerivativeOperator{T}, Dissipation, SplitForm<:Union{Val{false}, Val{true}}}
        if dissipation !== nothing
            @argcheck size(derivative) == size(dissipation) DimensionMismatch
            @argcheck grid(derivative) == grid(dissipation) ArgumentError
        end
        N = size(derivative, 2)
        tmp1 = Array{T}(undef, N)
        tmp2 = Array{T}(undef, N)
        new{T,Derivative,Dissipation,SplitForm}(derivative, dissipation, tmp1, tmp2, split_form)
    end
end

function Base.show(io::IO, semi::QuarticNonconvexPeriodicSemidiscretization)
    if get(io, :compact, false)
        print(io, "Semidiscretization of the quartic nonconvex conservation law (periodic)")
    else
        print(io, "Semidiscretization of the quartic nonconvex conservation law\n")
        print(io, "  \$ \\partial_t u(t,x) + \\partial_x ( u(t,x)^4 - 10 u(t,x)^2 + 3 u(t,x) ) = 0 \$ \n")
        print(io, "with periodic boundaries using")
        if semi.split_form == Val{true}()
            print(io, " a split form and: \n")
        else
            print(io, " no split form and: \n")
        end
        print(io, semi.derivative)
        print(io, semi.dissipation)
    end
end

function (disc::QuarticNonconvexPeriodicSemidiscretization)(du, u, p, t)
    @unpack tmp1, tmp2, derivative, dissipation, split_form = disc
    @boundscheck begin
        @argcheck length(u) == length(tmp1)
        @argcheck length(u) == length(tmp2)
        @argcheck length(u) == length(du)
    end

    # volume terms
    if typeof(split_form) <: Val{true}
        m_2_5 = -2*one(eltype(u)) / 5
        m_20_3 = -20*one(eltype(u)) / 3
        m_3 = -3*one(eltype(u))

        # # D u^4
        @. tmp2 = u^4
        mul!(tmp1, derivative, tmp2)
        @. du = m_2_5 * tmp1

        @. tmp2 = u^3
        mul!(tmp1, derivative, tmp2)
        @. du += m_2_5 * u * tmp1

        @. tmp2 = u^2
        mul!(tmp1, derivative, tmp2)
        @. du += m_2_5 * u^2 * tmp1

        mul!(tmp1, derivative, u)
        @. du += m_2_5 * u^3 * tmp1

        # D u^2
        @. tmp2 = u^2
        mul!(tmp1, derivative, tmp2)
        @. du -= m_20_3 * tmp1

        mul!(tmp1, derivative, u)
        @. du -= m_20_3 * u * tmp1

        # D u
        mul!(tmp1, derivative, u)
        @. du += m_3 * tmp1
    else
        @. tmp2 = u^4 - 10 * u^2 + 3 * u
        mul!(tmp1, derivative, tmp2)
        @. du = -tmp1
    end

    # dissipation
    if dissipation !== nothing
        mul!(tmp1, dissipation, u)
        @. du += tmp1
    end

    nothing
end


text:
"""
    hashfile("inputfile", keycol, valcols)

Create a dictionary from a tsv file with some
column of the file as the lookup key and values
from columns specified in an array.

"""
function hashfile(file::String, keycolumn::Int, valcols::Array; delim="\t")

    fhash = Dict{Any,Any}()
    d = Regex(delim)

    open(file) do f

        for i in eachline(f)
            line = split(i, d)
            fhash[line[keycolumn]] = line[valcols]
        end

    end

    return fhash

end

max_stars_repo_name: ScottWatkins/STRdiff

max_issues_repo_name: ScottWatkins/STRdiff

max_forks_repo_name: ScottWatkins/STRdiff